Changeset 33372

Timestamp:

Feb 25, 2012, 8:50:53 PM (14 years ago)

Author:

eugene

Message:

fix weighting scheme for different photcodes; fix reset & ubercal collisions

Location:

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src

Files:

• StarOps.c (modified) (7 diffs)
• bcatalog.c (modified) (3 diffs)

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/StarOps.c

@@ -76 +76 @@
         N = 0;
         for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
+
           // skip measurements that do not match the current photcode
-          int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
-          if (ecode != thisCode) { continue; }
+          PhotCode *code = GetPhotcodebyCode (catalog[i].measureT[m].photcode);
+          if (!code) continue;
+          if (code->equiv != thisCode) { continue; }
+
+          // int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
+          // if (ecode != thisCode) { continue; }
 
           if (catalog[i].measureT[m].dbFlags & MEAS_BAD) {
@@ -114 +119 @@
           }
           list[N] = Msys - Mcal - Mmos - Mgrid;
-          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
+
+          // dlist gives the error, which is used as the weight in WT_MEAN.
+          // we can modify the resulting weight in a few ways:
+          // 1) MIN_ERROR guarantees a floor
+          // 2) photomErrSys is added in quadrature as a sytematic error, set per photcode
+          // 3) UBERCAL measurements can have their weight increased by a big factor to help tie down the averages
+          // 4) some reference photcode of some kind can be specified as fixed and have a high weight
+          dlist[N] = MAX (hypot(catalog[i].measureT[m].dM, code->photomErrSys), MIN_ERROR);
 
           // up-weight the ubercal values (or convergence can take a long time...)
@@ -122 +134 @@
 
           // tie down reference photometry if the -refcode (code) option is selected
-          // XXX re-think this...
+          // eg, -refcode g_SDSS
+          // this probably makes no sense in the context of multifilter analysis
           if (refPhotcode) {
-            if (GetPhotcodeEquivCodebyCode(catalog[i].measureT[m].photcode) == refPhotcode[0].equiv) {
-              // increase the weight by a factor of 100:
-              dlist[N] = 0.01*catalog[i].measureT[m].dM;
+            if (code->code == refPhotcode->code) {
+              // tiny error -> large weight
+              // dlist[N] = MAX (0.01*catalog[i].measureT[m].dM, MIN_ERROR);
+              dlist[N] = 0.0001;
             }
           }
@@ -157 +171 @@
 // * setMrelOutput skips stars for which there are too few good measurements 
 // * setMrelOutput is meant to be called repeatedly, relaxing the criteria for 'good' on each pass
+
+// setMave is also similar to the above.  but:
+// * setMave is called by relphot -update-objects
+// * setMave includes PSF_QF tests 
+// * setMave updates 2MASS average flags
+// * setMave updates average EXT flags (PS1 and 2MASS)
 int setMrelOutput (Catalog *catalog, int Ncatalog, int mark, FlatCorrectionTable *flatcorr) {
 
@@ -164 +184 @@
   double *list, *dlist;
   StatType stats;
-  int Nsec, Nsecfilt, ecode;
-
-  Nsecfilt = GetPhotcodeNsecfilt ();
+  int Nsec, ecode;
+
+  int Nsecfilt = GetPhotcodeNsecfilt ();
 
   /* Nmeasure is now different, need to reallocate */
@@ -193 +213 @@
         m = catalog[i].averageT[j].measureOffset;
         for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
+
           // skip measurements that do not match the current photcode
-          ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
-          if (ecode != thisCode) { continue; }
+          PhotCode *code = GetPhotcodebyCode (catalog[i].measureT[m].photcode);
+          if (!code) continue;
+          if (code->equiv != thisCode) { continue; }
+
+          // ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
+          // if (ecode != thisCode) { continue; }
 
           if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
@@ -215 +240 @@
           Msys = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]);
           list[N] = Msys - Mcal - Mmos - Mgrid;
-          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
+
+          // dlist gives the error, which is used as the weight in WT_MEAN.
+          // we can modify the resulting weight in a few ways:
+          // 1) MIN_ERROR guarantees a floor
+          // 2) photomErrSys is added in quadrature as a sytematic error, set per photcode
+          // 3) UBERCAL measurements can have their weight increased by a big factor to help tie down the averages
+          // 4) some reference photcode of some kind can be specified as fixed and have a high weight
+          dlist[N] = MAX (hypot(catalog[i].measureT[m].dM, code->photomErrSys), MIN_ERROR);
+
+          // up-weight the ubercal values (or convergence can take a long time...)
+          if (catalog[i].measureT[m].dbFlags & ID_MEAS_PHOTOM_UBERCAL) {
+            dlist[N] = MAX (0.1*catalog[i].measureT[m].dM, MIN_ERROR);
+          }
+
+          // tie down reference photometry if the -refcode (code) option is selected
+          // eg, -refcode g_SDSS
+          // this probably makes no sense in the context of multifilter analysis
+          if (refPhotcode) {
+            if (code->code == refPhotcode->code) {
+              // tiny error -> large weight
+              // dlist[N] = MAX (0.01*catalog[i].measureT[m].dM, MIN_ERROR);
+              dlist[N] = 0.0001;
+            }
+          }
           N++;
         }

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/bcatalog.c

@@ -62 +62 @@
   subcatalog[0].filename = catalog[0].filename;
 
+  // flags used by the photometry analysis (excluding UBERCAL)
+  unsigned int PHOTOM_FLAGS = 
+    ID_MEAS_NOCAL |       // detection ignored for this analysis (photcode, time range)
+    ID_MEAS_POOR_PHOTOM | // detection is photometry outlier 
+    ID_MEAS_SKIP_PHOTOM | // detection was ignored for photometry measurement 
+    ID_MEAS_AREA;         // detetion was outside acceptable area of device
+
   /* exclude stars not in range or with too few measurements */
   for (i = 0; i < catalog[0].Naverage; i++) {
@@ -74 +81 @@
     }
 
+    // reset the calculated average magnitudes (does not affect ubercal-tied measurements or images)
     if (RESET) {
       int Ns;
@@ -147 +155 @@
       subcatalog[0].measureT[Nmeasure].averef = Naverage;
       if (RESET) { 
-        subcatalog[0].measureT[Nmeasure].Mcal = 0;
-        subcatalog[0].measureT[Nmeasure].dbFlags &= 0xff00;
-        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_POOR_PHOTOM;
-        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
-        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_NOCAL;
+        if (!KEEP_UBERCAL) {
+          subcatalog[0].measureT[Nmeasure].Mcal = 0;
+          subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_PHOTOM_UBERCAL;
+        }
+        subcatalog[0].measureT[Nmeasure].dbFlags &= ~PHOTOM_FLAGS;
       }
       Nmeasure ++;