IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 9633


Ignore:
Timestamp:
Oct 18, 2006, 8:54:26 AM (20 years ago)
Author:
eugene
Message:

work to allow REF measurements in the solution

Location:
trunk/Ohana/src/relphot
Files:
1 added
2 deleted
13 edited

Legend:

Unmodified
Added
Removed
  • trunk/Ohana/src/relphot/Makefile

    r8389 r9633  
    2121RELPHOT = \
    2222$(SRC)/ConfigInit.$(ARCH).o      \
    23 $(SRC)/GridOps.v2.$(ARCH).o      \
     23$(SRC)/GridOps.$(ARCH).o         \
    2424$(SRC)/ImageOps.$(ARCH).o        \
    2525$(SRC)/MosaicOps.$(ARCH).o       \
  • trunk/Ohana/src/relphot/doc/notes.txt

    r7080 r9633  
     1
     22006.10.17
     3
     4  I am working on a partial upgrade of relphot to handle more general
     5  problems than I have treated in the past.  In particular, I would
     6  like to be able to have it assign the average magnitudes, regardless
     7  of whether any or all images and/or objects are well-treated for
     8  relative photometry.  In fact, I would like it to be able to
     9  determine average magnitudes even for sources which have only
     10  externally supplied data, and thus have no matching images.  I would
     11  eventually like to have more flexibility about the filtering which
     12  is performed on the measurements to determine the average
     13  magnitudes.  Some of this work will require a bit of a more careful
     14  treatment of the process. 
     15
     16  At the moment, we can consider the relphot process to consist of the
     17  following steps:
     18
     19  1) load the complete catalog data
     20  2) select a subset of the sources:
     21     - appropriate photcode (equiv == selected photcode)
     22     - in time range, if specified
     23     - appropriate dophot type, if specified
     24     - bright enough:
     25       - mag < MAG_LIM 
     26       - dMag > SIGMA_LIM
     27       - ImagMin < iMag < ImagMax
     28     - other restrictions
     29       - within valid image region for chip
     30  3) iterate to a solution for Mcal (image offsets)
     31     - mark specific images as bad
     32     - mark specific sources as bad
     33     - mark specific measurements as bad
     34     - include desired restrictions on Mcal (eg, all chips of mosaic
     35       matched)
     36     - only keep the image offsets
     37
     38  4) if the results are to be kept, recalculate the average magnitudes
     39     using the image exclusions and Mcal values determined above
    140
    2412006.05.03
  • trunk/Ohana/src/relphot/include/relphot.h

    r7390 r9633  
    66/* # define GRID_V1 */
    77# define GRID_V2
     8# define NO_IMAGE -100
    89
    910# if (0)
     
    9192
    9293int ImagSelect, ImagMin, ImagMax;
     94
     95int DophotSelect, DophotValue;
    9396
    9497double  PlotMmin, PlotMmax, PlotdMmin, PlotdMmax;
  • trunk/Ohana/src/relphot/src/ImageOps.c

    r6683 r9633  
    142142
    143143  i = bin[cat][meas];
    144   if (i == -1) return (NO_MAG);
     144  if (i == -1) return (NO_IMAGE);
    145145
    146146  if (image[i].code & IMAGE_BAD)  return (NO_MAG); 
     
    194194     
    195195      if (catalog[c].measure[m].flags & MEAS_BAD) continue;
    196       if ((Mmos  = getMmos  (m, c)) == NO_MAG) continue;
     196      Mmos  = getMmos  (m, c);
     197      if (Mmos == NO_MAG) continue;
     198      if (Mmos == NO_IMAGE) continue;
    197199      if ((Mgrid = getMgrid (m, c)) == NO_MAG) continue;
    198200      if ((Mrel  = getMrel  (catalog, m, c)) == NO_MAG) continue;
     
    334336  int i, j, m, c, n, N;
    335337  double *list, *dlist;
     338  float Mcal, Mmos;
    336339  StatType stats;
    337340
     
    352355      c = clist[i][j];
    353356
    354       if (getMcal  (m, c) == NO_MAG) continue;
    355       if (getMmos  (m, c) == NO_MAG) continue;
     357      Mcal = getMcal  (m, c);
     358      if (Mcal == NO_MAG) continue;
     359      if (Mcal == NO_IMAGE) continue;
     360      Mmos = getMmos  (m, c);
     361      if (Mmos == NO_MAG) continue;
     362      if (Mmos == NO_IMAGE) continue;
    356363      if (getMgrid (m, c) == NO_MAG) continue;
    357364      N++;
  • trunk/Ohana/src/relphot/src/MosaicOps.c

    r6683 r9633  
    299299  if (!MOSAICNAME[0]) return (0);
    300300  i = bin[cat][meas];
    301   if (i == -1) return (NO_MAG);
     301  if (i == -1) return (NO_IMAGE);
    302302
    303303  if (mosaic[i].code & IMAGE_BAD) return (NO_MAG); 
     
    345345      c = clist[i][j];
    346346     
     347      Mcal  = getMcal  (m, c);
     348      if (Mcal == NO_MAG) continue;
     349      if (Mcal == NO_IMAGE) continue;
    347350      if (catalog[c].measure[m].flags & MEAS_BAD) continue;
    348       if ((Mcal  = getMcal  (m, c)) == NO_MAG) continue;
    349351      if ((Mgrid = getMgrid (m, c)) == NO_MAG) continue;
    350352      if ((Mrel  = getMrel  (catalog, m, c)) == NO_MAG) continue;
     
    447449  int i, j, m, c, n, N;
    448450  double *list, *dlist;
     451  float Mcal;
    449452  StatType stats;
    450453
     
    465468      c = clist[i][j];
    466469
    467       if (getMcal  (m, c) == NO_MAG) continue;
     470      Mcal = getMcal  (m, c);
     471      if (Mcal == NO_MAG) continue;
     472      if (Mcal == NO_IMAGE) continue;
    468473      if (getMgrid (m, c) == NO_MAG) continue;
    469474      if (getMrel  (catalog, m, c) == NO_MAG) continue;
  • trunk/Ohana/src/relphot/src/StarOps.c

    r6683 r9633  
    4545      m = catalog[i].average[j].offset;
    4646
     47      // XXX allow REF stars to be included in the calculation
     48      // this should be optionally set, and should allow for
     49      // REF stars to be downweighted by more than their reported
     50      // errors.  how such info is carried is unclear...
    4751      N = 0;
    4852      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    4953        if (catalog[i].measure[m].flags & MEAS_BAD) continue;
    50         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    51         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
    52         if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
     54        Mcal  = getMcal  (m, i);
     55        if (Mcal == NO_MAG) continue;
     56        if (Mcal == NO_IMAGE) {
     57          Mcal = Mmos = Mgrid = 0;
     58        } else {
     59          Mmos  = getMmos  (m, i);
     60          if (Mmos == NO_MAG) continue;
     61          if (Mmos == NO_IMAGE) continue;
     62          if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
     63        }
    5364
    5465        Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
     
    106117      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    107118        if (catalog[i].measure[m].flags & MEAS_BAD) continue;
    108         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    109         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
    110         if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
     119        Mcal  = getMcal  (m, i);
     120        if (Mcal == NO_MAG) continue;
     121        if (Mcal == NO_IMAGE) {
     122          Mcal = Mmos = Mgrid = 0;
     123        } else {
     124          Mmos  = getMmos  (m, i);
     125          if (Mmos == NO_MAG) continue;
     126          if (Mmos == NO_IMAGE) continue;
     127          if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
     128        }
    111129
    112130        Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
     
    153171      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    154172        if (catalog[i].measure[m].flags & MEAS_BAD) continue;
    155         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    156         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
     173        Mcal  = getMcal  (m, i);
     174        if (Mcal == NO_MAG) continue;
     175        if (Mcal == NO_IMAGE) continue;
     176        Mmos  = getMmos  (m, i);
     177        if (Mmos == NO_MAG) continue;
     178        if (Mmos == NO_IMAGE) continue;
    157179        if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
    158180        catalog[i].measure[m].Mcal_PS = Mcal + Mmos + Mgrid;
     
    263285      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    264286        /* if (catalog[i].measure[m].flags & MEAS_BAD) continue; */
    265         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    266         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
     287        Mcal  = getMcal  (m, i);
     288        if (Mcal == NO_MAG) continue;
     289        if (Mcal == NO_IMAGE) continue;
     290        Mmos  = getMmos  (m, i);
     291        if (Mmos == NO_MAG) continue;
     292        if (Mmos == NO_IMAGE) continue;
    267293        if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
    268294
     
    295321      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    296322        /* if (catalog[i].measure[m].flags & MEAS_BAD) continue; */
    297         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    298         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
     323        Mcal  = getMcal  (m, i);
     324        if (Mcal == NO_MAG) continue;
     325        if (Mcal == NO_IMAGE) continue;
     326        Mmos  = getMmos  (m, i);
     327        if (Mmos == NO_MAG) continue;
     328        if (Mmos == NO_IMAGE) continue;
    299329        if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
    300330
     
    328358  int i, j, k, m, n, N, Ntot;
    329359  double *list, *dlist;
     360  float Mcal, Mmos;
    330361  StatType stats;
    331362
     
    348379      N = 0;
    349380      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    350         if (getMcal  (m, i) == NO_MAG) continue;
    351         if (getMmos  (m, i) == NO_MAG) continue;
     381        Mcal = getMcal  (m, i);
     382        if (Mcal == NO_MAG) continue;
     383        if (Mcal == NO_IMAGE) continue;
     384        Mmos = getMmos  (m, i);
     385        if (Mmos == NO_MAG) continue;
     386        if (Mmos == NO_IMAGE) continue;
    352387        if (getMgrid (m, i) == NO_MAG) continue;
    353388        N++;
  • trunk/Ohana/src/relphot/src/args.c

    r7080 r9633  
    171171    remove_argument (N, &argc, argv);
    172172    ImagSelect = TRUE;
     173  }
     174
     175  DophotSelect = FALSE;
     176  if ((N = get_argument (argc, argv, "-dophot"))) {
     177    remove_argument (N, &argc, argv);
     178    DophotValue = atof (argv[N]);
     179    remove_argument (N, &argc, argv);
     180    DophotSelect = TRUE;
    173181  }
    174182
  • trunk/Ohana/src/relphot/src/bcatalog.c

    r8645 r9633  
    1818  for (i = 0; i < catalog[0].Naverage; i++) {
    1919    if (catalog[0].average[i].Nm < 2) continue;
    20 
    21     /* XXX this limitation is absurd
    22     if (catalog[0].average[i].R < fullregion[0].RA[0]) continue;
    23     if (catalog[0].average[i].R > fullregion[0].RA[1]) continue;
    24     if (catalog[0].average[i].D < fullregion[0].DEC[0]) continue;
    25     if (catalog[0].average[i].D > fullregion[0].DEC[1]) continue;
    26     */
    2720
    2821    /* start with all stars good */
     
    5952
    6053      /* select measurements by quality */
    61       if (catalog[0].measure[offset].dophot != 1) continue;
     54      // XXX ignore this criterion for REF measurements?
     55      if (DophotSelect && (catalog[0].measure[offset].dophot != DophotValue)) continue;
    6256
    6357      /* select measurements by mag limit */
     
    6660
    6761      /* select measurements by measurement error */
    68       if (catalog[0].measure[offset].dM_PS > SIGMA_LIM) continue;
     62      if ((SIGMA_LIM > 0) && (catalog[0].measure[offset].dM_PS > SIGMA_LIM)) continue;
    6963
    7064      /* select measurements by mag limit */
     
    9286      }
    9387    }
     88
     89    // XXXX test : what checks do I need to make elsewhere to avoid problems here?
     90    # if 1
    9491    if (Nm < 2) { /* enough measurements in band? */
    9592      Nmeasure -= Nm;
    9693      continue;
    9794    }
     95    # endif
    9896    subcatalog[0].average[Naverage].Nm = Nm;
    9997    Naverage ++;
  • trunk/Ohana/src/relphot/src/load_catalogs.c

    r8386 r9633  
    1717    tcatalog.catformat = dvo_catalog_catformat (CATFORMAT);    // set the default catformat from config data
    1818    tcatalog.catmode   = dvo_catalog_catmode (CATMODE);        // set the default catmode from config data
    19     tcatalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_SECF; // don't need to load all data at this point
     19    tcatalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_SECF;    // don't need to load all data at this point
     20    tcatalog.Nsecfilt  = GetPhotcodeNsecfilt ();               // set the desired number in case we need to create the catalog
    2021
    2122    if (!dvo_catalog_open (&tcatalog, skylist[0].regions[i], VERBOSE, "r")) {
     
    3738  }
    3839  if (Nstar < 2) {
    39     fprintf (stderr, "insufficient stars %d\n", Nstar);
    40     exit (0);
     40    fprintf (stderr, "warning: insufficient stars %d\n", Nstar);
    4141  }
    4242
  • trunk/Ohana/src/relphot/src/plot_scatter.c

    r6683 r9633  
    2828      for (k = 0; k < catalog[i].average[j].Nm; k++, m++) {
    2929        if (catalog[i].measure[m].flags & MEAS_BAD) continue;
    30         if ((Mcal  = getMcal  (m, i)) == NO_MAG) continue;
    31         if ((Mmos  = getMmos  (m, i)) == NO_MAG) continue;
     30        Mcal  = getMcal  (m, i);
     31        if (Mcal == NO_MAG) continue;
     32        if (Mcal == NO_IMAGE) continue;
     33        Mmos  = getMmos  (m, i);
     34        if (Mmos == NO_MAG) continue;
     35        if (Mmos == NO_IMAGE) continue;
    3236        if ((Mgrid = getMgrid (m, i)) == NO_MAG) continue;
    3337
  • trunk/Ohana/src/relphot/src/reload_catalogs.c

    r8386 r9633  
    2424    catalog.catmode   = dvo_catalog_catmode (CATMODE);        // set the default catmode from config data
    2525    catalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
     26    catalog.Nsecfilt  = GetPhotcodeNsecfilt ();               // set the desired number in case we need to create the catalog
    2627
    2728    if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) {
  • trunk/Ohana/src/relphot/src/relphot.c

    r8386 r9633  
    1414  /* register database handle with shutdown procedure */
    1515  set_db (&db);
     16  db.mode   = dvo_catalog_catmode (CATMODE);
     17  db.format = dvo_catalog_catformat (CATFORMAT);
    1618
    1719  /* lock and load the image db table */
    1820  status = dvo_image_lock (&db, ImageCat, 60.0, (UPDATE ? LCK_XCLD : LCK_SOFT));
    19   if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
    20   if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
    21   if (!dvo_image_load (&db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db.filename);
     21  if (!status && UPDATE) {
     22    fprintf (stderr, "error\n");
     23    Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
     24  }
     25  // if the file is missing, db.dbstate will have a value of either:
     26  // LCK_EMPTY (if UPDATE) or LCK_MISSING (if !UPDATE)
     27  if ((db.dbstate == LCK_EMPTY) || (db.dbstate == LCK_MISSING)) {
     28    // XXX get ZERO_POINT from config
     29    dvo_image_create (&db, 25.0);
     30    // Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
     31  } else {
     32    if (!dvo_image_load (&db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db.filename);
     33  }
    2234
    2335  /* load regions and images based on specified sky patch */
  • trunk/Ohana/src/relphot/src/setMrelFinal.c

    r6683 r9633  
    110110
    111111      /* skip measurements from BAD images and mosaics */
     112      /* do NOT skip measurements without a matching image */
    112113      if ((getMcal  (m, 0)) == NO_MAG) goto skip;
    113114      if ((getMmos  (m, 0)) == NO_MAG) goto skip;
Note: See TracChangeset for help on using the changeset viewer.