IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 30441


Ignore:
Timestamp:
Jan 31, 2011, 2:34:15 PM (15 years ago)
Author:
eugene
Message:

update exclusions to be more flexible

Location:
branches/eam_branches/ipp-20101205/Ohana/src/photdbc
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20101205/Ohana/src/photdbc/include/photdbc.h

    r29938 r30441  
    5151double DMGAIN;
    5252double CHISQ_MAX;
     53double SIGMA_MIN_KEEP;
    5354double SIGMA_MAX;
    5455double AVE_SIGMA_LIM;
    5556int    NMEAS_MIN;
     57int    NMEAS_MIN_FILTERED;
     58int    NCODE_MIN;
    5659double ZERO_POINT;
     60
     61int ExcludeByMinSigma;
    5762
    5863int ExcludeByInstMag;
    … …  
    6570SkyRegion REGION;
    6671PhotCodeData photcodes;
     72
     73char          *PHOTCODE_DROP_LIST, *PHOTCODE_SKIP_LIST;
     74int           NphotcodesDrop,      NphotcodesSkip;
     75PhotCode     **photcodesDrop,     **photcodesSkip;
    6776
    6877# define FLAG_AREA            0X0001

  • branches/eam_branches/ipp-20101205/Ohana/src/photdbc/src/ConfigInit.c

    r28331 r30441  
    2626  if (VERBOSE) fprintf (stderr, "loaded config file: %s\n", file);
    2727
    28   WarnConfig (config, "PHOTDBC_JOIN_RADIUS",       "%lf", 0, &JOIN_RADIUS);
     28  // XXX join_stars (in photdbc.c) is currently disabled
     29  // WarnConfig (config, "PHOTDBC_JOIN_RADIUS",       "%lf", 0, &JOIN_RADIUS);
     30
     31  // XXX unique_measures (in photdbc.c) is currently disabled
    2932  // WarnConfig (config, "UNIQ_RADIUS",            "%lf", 0, &UNIQ_RADIUS);
    3033
     34  // XXX flag_measures (in photdbc.c) is currently disabled
    3135  // WarnConfig (config, "XMIN",                   "%lf", 0, &XMIN);
    3236  // WarnConfig (config, "XMAX",                   "%lf", 0, &XMAX);
    … …  
    3741  // WarnConfig (config, "MMIN",                   "%lf", 0, &tmp);  MMIN      = 1000*tmp;
    3842  // WarnConfig (config, "MMAX",                   "%lf", 0, &tmp);  MMAX      = 1000*tmp;
     43  // WarnConfig (config, "CHISQ_MAX",              "%lf", 0, &CHISQ_MAX);
     44
     45  // XXX get_mags (in photdbc.c) is currently disabled
    3946  // WarnConfig (config, "DMSYS",                  "%lf", 0, &tmp);  DMSYS     = SQ(1000*tmp);
    4047  // WarnConfig (config, "DMGAIN",                 "%lf", 0, &DMGAIN);
    41   // WarnConfig (config, "CHISQ_MAX",              "%lf", 0, &CHISQ_MAX);
    4248
    4349  ScanConfig (config, "SIGMA_MAX",              "%lf", 0, &SIGMA_MAX);
    4450  ScanConfig (config, "AVE_SIGMA_LIM",          "%lf", 0, &AVE_SIGMA_LIM);
    4551  ScanConfig (config, "NMEAS_MIN",              "%d",  0, &NMEAS_MIN);
     52  ScanConfig (config, "NMEAS_MIN_FILTERED",     "%d",  0, &NMEAS_MIN_FILTERED);
    4653
    4754  WarnConfig (config, "GSCFILE",                "%s",  0, GSCFILE);

  • branches/eam_branches/ipp-20101205/Ohana/src/photdbc/src/args.c

    r29938 r30441  
    1919  ExcludeByInstMag = FALSE;
    2020  if ((N = get_argument (argc, argv, "-instmag"))) {
    21     remove_argument (N, &argc, argv);
    2221    ExcludeByInstMag = TRUE;
    2322    remove_argument (N, &argc, argv);
    … …  
    2524    remove_argument (N, &argc, argv);
    2625    INST_MAG_MAX = atof(argv[N]);
     26    remove_argument (N, &argc, argv);
     27  }
     28
     29  ExcludeByMinSigma = FALSE;
     30  if ((N = get_argument (argc, argv, "-min-sigma"))) {
     31    ExcludeByMinSigma = TRUE;
     32    remove_argument (N, &argc, argv);
     33    SIGMA_MIN_KEEP = atof(argv[N]);
     34    remove_argument (N, &argc, argv);
    2735  }
    2836
    2937  ExcludeByMaxMinMag = FALSE;
    3038  if ((N = get_argument (argc, argv, "-maxminmag"))) {
    31     remove_argument (N, &argc, argv);
    3239    ExcludeByMaxMinMag = TRUE;
    3340    remove_argument (N, &argc, argv);
    3441    MAX_MIN_MAG = atof(argv[N]);
     42    remove_argument (N, &argc, argv);
    3543  }
    3644
    … …  
    6169  }
    6270
     71  PHOTCODE_DROP_LIST = NULL;
     72  if ((N = get_argument (argc, argv, "-photcode-drop"))) {
     73    remove_argument (N, &argc, argv);
     74    PHOTCODE_DROP_LIST = strcreate(argv[N]);
     75    remove_argument (N, &argc, argv);
     76  }
     77
     78  PHOTCODE_SKIP_LIST = NULL;
     79  if ((N = get_argument (argc, argv, "-photcode-skip"))) {
     80    remove_argument (N, &argc, argv);
     81    PHOTCODE_SKIP_LIST = strcreate(argv[N]);
     82    remove_argument (N, &argc, argv);
     83  }
     84
    6385  if (argc != 2) usage();
    6486
    … …  
    81103  fprintf (stderr, "USAGE: photdbc (output)\n\n");
    82104  fprintf (stderr, " this program takes an existing DVO database and makes a copy, applying a number of optional\n");
    83   fprintf (stderr, " filters in the process.  \n");
    84   fprintf (stderr, "\n");
    85   fprintf (stderr, " photdbc (output)\n");
    86   fprintf (stderr, "\n");
    87   fprintf (stderr, " allowed filters / restrictions include:\n");
    88   fprintf (stderr, "\n");
     105  fprintf (stderr, " filters in the process.  \n\n");
     106
     107  fprintf (stderr, " photdbc (output)\n\n");
     108
     109  fprintf (stderr, " allowed filters / restrictions include:\n\n");
     110
    89111  fprintf (stderr, " -region Rmin Rmax Dmin Dmax : limit operation to the specified region \n");
    90   fprintf (stderr, " -join                 : join measurements between stars using JOIN_RADIUS\n");
    91   fprintf (stderr, " -ccdregion X Y X Y    : only keep detections within the specified detector region\n");
    92   fprintf (stderr, "                         (can this be limited to specific photcodes? cameras? detectors?)\n");
    93   fprintf (stderr, " -photcode_limits code Mmin Mmax : allow multiples of these\n");
    94   fprintf (stderr, "\n");
    95   fprintf (stderr, " SIGMA_MAX\n");
    96   fprintf (stderr, " NMEAS_MIN\n");
    97   fprintf (stderr, " INST_MAG_MAX\n");
    98   fprintf (stderr, " INST_MAG_MIN\n");
     112 
     113  fprintf (stderr, " -photcode-drop   : remove these photcodes from the output (REF or DEP only)\n");
     114  fprintf (stderr, " -photcode-skip  : ignore these photcodes when assessing the validity (keep unless object is dropped)\n");
     115
     116  fprintf (stderr, " -instmag (min) (max) : range of valid instrumental magnitudes (or measurements are dropped)\n");
     117  fprintf (stderr, " -min-sigma (sigma)   : object must have one measurement error less than sigma or object is dropped\n");
     118
     119  fprintf (stderr, " -maxminmag (mag)   : object must have one magnitude less than this or object is dropped\n");
     120
     121  fprintf (stderr, " option options:\n");
     122  fprintf (stderr, " -v : verbose mode\n");
     123  fprintf (stderr, " -params : list the current parameters\n\n");
     124
     125  fprintf (stderr, "ptolemy.rc config values used by this program:\n");
     126  fprintf (stderr, " SIGMA_MAX : drop measurements with errors greater than this\n");
     127  fprintf (stderr, " NMEAS_MIN : drop objects with fewer measurements than this\n");
     128  fprintf (stderr, " NMEAS_MIN_FILTERED : drop objects with fewer measurements than this after filtering above\n");
     129  fprintf (stderr, " AVE_SIGMA_LIM : drop objects if all average mags are greater than this\n");
     130
    99131  exit (2);
    100132}
     133
     134// fprintf (stderr, " -join                 : join measurements between stars using JOIN_RADIUS\n");
     135// fprintf (stderr, " -ccdregion X Y X Y    : only keep detections within the specified detector region\n");
     136// fprintf (stderr, "                         (can this be limited to specific photcodes? cameras? detectors?)\n");
     137// fprintf (stderr, " -photcode_limits code Mmin Mmax : allow multiples of these\n");

  • branches/eam_branches/ipp-20101205/Ohana/src/photdbc/src/initialize.c

    r29938 r30441  
    22
    33void initialize (int argc, char **argv) {
     4
     5  int NPHOTCODES;
     6  char *codename, *ptr, *list;
    47
    58  /* are these set correctly? */
    … …  
    1114  ConfigInit (&argc, argv);
    1215  args (argc, argv);
     16
     17  NphotcodesDrop = 0;
     18  photcodesDrop = NULL;
     19  if (PHOTCODE_DROP_LIST != NULL) {
     20    NPHOTCODES = 10;
     21    ALLOCATE (photcodesDrop, PhotCode *, NPHOTCODES);
     22
     23    /* parse the comma-separated list of photcodesDrop */
     24    list = PHOTCODE_DROP_LIST;
     25    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     26      list = NULL; // pass NULL on successive strtok_r calls
     27      fprintf (stderr, "PHOTCODE_LIST: %s\n", PHOTCODE_DROP_LIST);
     28      fprintf (stderr, "codename: %s\n", codename);
     29      if ((photcodesDrop[NphotcodesDrop] = GetPhotcodebyName (codename)) == NULL) {
     30        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     31        exit (1);
     32      }
     33      NphotcodesDrop ++;
     34      CHECK_REALLOCATE (photcodesDrop, PhotCode *, NPHOTCODES, NphotcodesDrop, 10);
     35    }
     36  }
     37
     38  NphotcodesSkip = 0;
     39  photcodesSkip = NULL;
     40  if (PHOTCODE_SKIP_LIST != NULL) {
     41    NPHOTCODES = 10;
     42    ALLOCATE (photcodesSkip, PhotCode *, NPHOTCODES);
     43
     44    /* parse the comma-separated list of photcodesSkip */
     45    list = PHOTCODE_SKIP_LIST;
     46    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     47      list = NULL; // pass NULL on successive strtok_r calls
     48      fprintf (stderr, "PHOTCODE_LIST: %s\n", PHOTCODE_SKIP_LIST);
     49      fprintf (stderr, "codename: %s\n", codename);
     50      if ((photcodesSkip[NphotcodesSkip] = GetPhotcodebyName (codename)) == NULL) {
     51        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     52        exit (1);
     53      }
     54      NphotcodesSkip ++;
     55      CHECK_REALLOCATE (photcodesSkip, PhotCode *, NPHOTCODES, NphotcodesSkip, 10);
     56    }
     57  }
    1358
    1459  if (SHOW_PARAMS) {

  • branches/eam_branches/ipp-20101205/Ohana/src/photdbc/src/make_subcatalog.c

    r28331 r30441  
    55int make_subcatalog (Catalog *subcatalog, Catalog *catalog) {
    66 
    7   off_t i, j, offset;
     7  int found;
     8  off_t i, j, k, offset;
    89  off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nm, Nsecfilt;
    9   double mag, minMag;
    10   int keep;
     10  double mag, minMag, minSigma;
     11  int keep, *secSkip;
     12  PhotCode *photcode;
    1113 
    1214  Nsecfilt = GetPhotcodeNsecfilt ();
    1315  assert (catalog[0].Nsecfilt == Nsecfilt);
     16
     17  // set up a list of SEC entries to ignore when evaluating a source
     18  ALLOCATE (secSkip, int, Nsecfilt);
     19  for (i = 0; i < Nsecfilt; i++) {
     20      secSkip[i] = FALSE;
     21      photcode = GetPhotcodebyNsec(i);
     22      for (k = 0; k < NphotcodesSkip; k++) {
     23          if (photcodesSkip[k][0].code != photcode[0].code) continue;
     24          secSkip[i] = TRUE;
     25      }
     26  }
    1427
    1528  /* we are moving only the subset of measurements from catalog[0] to subcatalog[0] */
    … …  
    3851    }
    3952
     53    // exclude stars with too few measurements
     54    if (NCODE_MIN) {
     55      // drop if all of the allowed average photcodes have ncode < NCODE_MIN values
     56      keep = FALSE;
     57      for (j = 0; !keep && (j < Nsecfilt); j++) {
     58          if (secSkip[j]) continue;
     59          if (catalog[0].secfilt[Nsecfilt*i+j].Ncode >= NCODE_MIN) {
     60              keep = TRUE;
     61        }
     62      }
     63      if (!keep) continue;
     64    }
     65
    4066    /* assign average and secfilt values */
    4167    subcatalog[0].average[Naverage] = catalog[0].average[i];
    … …  
    4571    }
    4672
    47     minMag = 32;
     73    minMag   = 32;
     74    minSigma = 32;
    4875    Nm = 0;
    4976    for (j = 0; j < catalog[0].average[i].Nmeasure; j++) {
    … …  
    5279
    5380      # if 0
    54       if (DropPhotcode) {
    55         ecode = GetPhotcodeEquivCodebyCode (catalog[0].measure[offset].photcode);
    56         if (ecode == photcode[0].code) continue;
    57       }
    5881      if (DropNonStellar && catalog[0].measure[offset].dophot != 1) continue;
    5982      # endif
    6083
    61       /* exclude measurements by measurement error */
     84      // remove certain photcodes from the output measurements
     85      if (NphotcodesDrop > 0) {
     86          found = FALSE;
     87          for (k = 0; (k < NphotcodesSkip) && !found; k++) {
     88              if (photcodesSkip[k][0].code == catalog[0].measure[offset].photcode) found = TRUE;
     89              if (photcodesSkip[k][0].code == GetPhotcodeEquivCodebyCode(catalog[0].measure[offset].photcode)) found = TRUE;
     90          }
     91          if (found) continue;
     92      } 
     93 
     94      // ignore certain photcodes to assess the measurements
     95      if (NphotcodesSkip > 0) {
     96          found = FALSE;
     97          for (k = 0; (k < NphotcodesSkip) && !found; k++) {
     98              if (photcodesSkip[k][0].code == catalog[0].measure[offset].photcode) found = TRUE;
     99              if (photcodesSkip[k][0].code == GetPhotcodeEquivCodebyCode(catalog[0].measure[offset].photcode)) found = TRUE;
     100          }
     101          if (found) goto keep;
     102      } 
     103 
     104      // exclude measurements by measurement error -- drop exactly this measurement
    62105      if (SIGMA_MAX && (catalog[0].measure[offset].dM > SIGMA_MAX)) continue;
    63106
    64       /* select measurements by mag limit */
     107      // select measurements by mag limit -- drop exactly this measurement 
    65108      if (ExcludeByInstMag) {
    66109        mag = PhotInst (&catalog[0].measure[offset]);
    … …  
    69112      }
    70113
     114      // check measurements for this object -- drop object if no measurements pass
     115      if (ExcludeByMinSigma) {
     116          minSigma = MIN (minSigma, catalog[0].measure[offset].dM);
     117      }
     118
     119      // check measurements for this object -- drop object if no measurements pass
    71120      if (ExcludeByMaxMinMag) {
    72121        mag = PhotSys (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]);
    73122        minMag = MIN (minMag, mag);
    74123      }
     124
     125    keep:
    75126
    76127      subcatalog[0].measure[Nmeasure]        = catalog[0].measure[offset];
    … …  
    90141    }
    91142
    92     // after measurement exclusion, exclude stars with too few measurements
    93     if (Nm < NMEAS_MIN) {
     143    // exclude faint objects
     144    if (ExcludeByMinSigma && (minSigma > SIGMA_MIN_KEEP)) {
    94145      Nmeasure -= Nm;
    95146      continue;
    96147    }
     148
     149    // after measurement exclusion, exclude stars with too few measurements
     150    if (NMEAS_MIN_FILTERED && (Nm < NMEAS_MIN_FILTERED)) {
     151      Nmeasure -= Nm;
     152      continue;
     153    }
     154
    97155    subcatalog[0].average[Naverage].Nmissing = 0;
    98156    subcatalog[0].average[Naverage].Nmeasure = Nm;
    … …  
    122180  return (TRUE);
    123181}
     182
     183/**
     184    the purpose of this function is to create a subset database.  there are several ways this could be done:
     185
     186    * reject all measurements from all objects that fail to meet some criteria (objects may lose measurements and/or be dropped)
     187
     188    * only reject an object if all measurements fail to meet some criteria (ie, keep all measurements if any measurement passes)
     189
     190    * only apply the keep / reject criteria to certain photcodes
     191
     192
     193    ***
     194
     195    * options:
     196
     197    minimum value for nmeas : NMEAS_MIN
     198    minimum value for nmeas : NMEAS_MIN_FILTERED
     199    minimum value for ncode : NCODE_MIN (drop if all ncode < NCODE_MIN) (respect photcodes)
     200    reject faint meas       : SIGMA_MAX
     201    reject faint source     : SIGMA_MIN_KEEP (keep source source if its minimum dMag < this limit)
     202
     203    **/
Note: See TracChangeset for help on using the changeset viewer.