Changeset 33807 for branches/eam_branches/ipp-20120405/Ohana

Timestamp:

Apr 19, 2012, 6:22:30 AM (14 years ago)

Author:

eugene

Message:

only load images for which the full set of overlapping catalogs were loaded; generate RA range test bands to correctly include / exclude images close to the region boundary; the -reset option should not reset Mcal for any measurements from images which were not loaded (bcatalog needs to accept the current catalog number as a result); when calculating mag:ave values (setMrel), the measurements from unloaded images need to use their own Mcal values (measure.Mcal)

Location:

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

Files:

• include/relphot.h (modified) (1 diff)
• src/StarOps.c (modified) (3 diffs)
• src/bcatalog.c (modified) (2 diffs)
• src/load_catalogs.c (modified) (1 diff)
• src/select_images.c (modified) (2 diffs)
• src/setMrelFinal.c (modified) (1 diff)

branches/eam_branches/ipp-20120405/Ohana/src/relphot/include/relphot.h

@@ -198 +198 @@
 int           args                PROTO((int argc, char **argv));
 int           args_client         PROTO((int argc, char **argv));
-int           bcatalog            PROTO((Catalog *subcatalog, Catalog *catalog));
+int           bcatalog            PROTO((Catalog *subcatalog, Catalog *catalog, int Ncat));
 void          clean_images        PROTO((void));
 void          clean_measures      PROTO((Catalog *catalog, int Ncatalog, int final, FlatCorrectionTable *flatcorr));

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

@@ -178 +178 @@
 
         if (getImageEntry (m, Nc) < 0) {
-          Mcal = Mmos = Mgrid = 0;
+          // measurements without an image are either external reference photometry or
+          // data for which the associated image has not been loaded (probably because of
+          // overlaps).  Msys + measure.Mcal is our best guess of the true magnitude
+          Mmos = Mgrid = 0;
+          Mcal = catalog[Nc].measureT[m].Mcal; // check that this is zero for loaded REF value
         } else {
           Mcal  = getMcal  (m, Nc, flatcorr, catalog);
@@ -429 +433 @@
 
           if (getImageEntry (m, i) < 0) {
+          XXXX : likely problem -- compare with setMrel_catalog;
             // these detetions have no image (eg, ref values such as 2MASS)
             Mcal = Mmos = Mgrid = 0;
@@ -795 +800 @@
         if (isnan(Mgrid)) continue;
 
+        // note that measurements for which the image is not selected will not up modified
+        // (that is a good thing! -- we keep a prior calibration)
+
         // set the output calibration
         catalog[i].measure[m].Mcal = Mcal + Mmos + Mgrid;

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

@@ -3 +3 @@
 int LimitDensityCatalog_ByNmeasure (Catalog *subcatalog, Catalog *oldcatalog);
 
-int bcatalog (Catalog *subcatalog, Catalog *catalog) {
+int bcatalog (Catalog *subcatalog, Catalog *catalog, int Ncat) {
   
   off_t i, j, offset;
@@ -159 +159 @@
       subcatalog[0].measureT[Nmeasure].averef = Naverage;
       if (RESET) { 
-        if (!KEEP_UBERCAL) {
-          subcatalog[0].measureT[Nmeasure].Mcal = 0;
-          subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_PHOTOM_UBERCAL;
+        // only reset Mcal for measures with a matching image
+        // do not reset Mcal for ubercal images unless explicitly requested
+        if (subcatalog[0].measureT[Nmeasure].dbFlags & ID_MEAS_PHOTOM_UBERCAL) {
+          if (!KEEP_UBERCAL) {
+            subcatalog[0].measureT[Nmeasure].Mcal = 0.0;
+            subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_PHOTOM_UBERCAL;
+          } 
+        } else {
+          if (getImageEntry (Nmeasure, Ncat) >= 0) {
+            subcatalog[0].measureT[Nmeasure].Mcal = 0.0;
+          }
         }
         subcatalog[0].measureT[Nmeasure].dbFlags &= ~PHOTOM_FLAGS;

branches/eam_branches/ipp-20120405/Ohana/src/relphot/src/load_catalogs.c

@@ -68 +68 @@
 
     // select only the brighter stars
-    bcatalog (&catalog[i], &tcatalog);
+    bcatalog (&catalog[i], &tcatalog, i);
     dvo_catalog_unlock (&tcatalog);
     dvo_catalog_free (&tcatalog);

branches/eam_branches/ipp-20120405/Ohana/src/relphot/src/select_images.c

@@ -122 +122 @@
     RmaxSkyRegion = RmaxSkyRegionAlt;
   }
-  
   double RmidSkyRegion = 0.5*(RminSkyRegion + RmaxSkyRegion);
   MARKTIME("create sky region coords: %f sec\n", dtime);
 
-  // if we overlap 0,360, RminSkyRegion will be 180 (or so)
-
   dsortindex (RmaxSky, index, skylist[0].Nregions);
   MARKTIME("sort sky coords: %f sec\n", dtime);
+
+  // generate a set of Rmin/Rmax values using this RmidSkyRegion for DEC bands of width
+  // dDecBand (0.1 deg).  DecIndex = (Dec + 90) / dDecBand
+  double *RminBand, *RmaxBand;
+  float dDecBand = 0.1;
+  int NDecBands = 180 / dDecBand + 1;
+  ALLOCATE (RminBand, double, NDecBands);
+  ALLOCATE (RmaxBand, double, NDecBands);
+  for (i = 0; i < NDecBands; i++) {
+    RminBand[i] = +360.0;
+    RmaxBand[i] = -360.0;
+  }
+  
+   /* compare with each region file */
+  for (i = 0; i < skylist[0].Nregions; i++) { 
+    int iDecBandMin = (skylist[0].regions[i][0].Dmin + 90.0) / dDecBand;
+    int iDecBandMax = (skylist[0].regions[i][0].Dmax + 90.0) / dDecBand;
+
+    double RminAlt = ohana_normalize_angle_to_midpoint (skylist[0].regions[i][0].Rmin, RmidSkyRegion);
+    double RmaxAlt = ohana_normalize_angle_to_midpoint (skylist[0].regions[i][0].Rmax, RmidSkyRegion);
+
+    for (j = iDecBandMin; j <= iDecBandMax; j++) {
+      RminBand[j] = MIN(RminBand[j], RminAlt);
+      RmaxBand[j] = MAX(RmaxBand[j], RmaxAlt);
+    }
+  }
+
+# if (0)
+  // XXX quick test
+  FILE *fout = fopen ("dec.bands.dat", "w");
+  for (i = 0; i < NDecBands; i++) {
+    fprintf (fout, "%d %f : %f - %f\n", (int) i, dDecBand * i - 90.0, RminBand[i], RmaxBand[i]);
+  }
+  fclose (fout);
+# endif
 
   if (VERBOSE) fprintf (stderr, "finding images\n");
@@ -209 +241 @@
     if (RminImage > RmaxSkyRegion) continue;
     if (RmaxImage < RminSkyRegion) continue;
+
+    // the above checks are only valid for the outermost region.  however, at a given
+    // declination, the range of RA is limited by the actual catalog boundaries
+    // check if this image is in range for the Dmax and Dmin locations
+
+    int iDecBandMin = (DminImage + 90.0) / dDecBand;
+    int iDecBandMax = (DmaxImage + 90.0) / dDecBand;
+
+    // the sky region RA is defined to be 0 - 360.0
+    if (RminImage > RmaxBand[iDecBandMin]) { 
+      fprintf (stderr, "skip image %s (%f,%f) on boundary\n", timage[i].name, 0.5*(RminImage + RmaxImage), 0.5*(DminImage + DmaxImage));
+      continue;
+    }
+    if (RminImage > RmaxBand[iDecBandMax]) {
+      fprintf (stderr, "skip image %s (%f,%f) on boundary\n", timage[i].name, 0.5*(RminImage + RmaxImage), 0.5*(DminImage + DmaxImage));
+      continue;
+    }
+    if (RmaxImage < RminBand[iDecBandMin]) {
+      fprintf (stderr, "skip image %s (%f,%f) on boundary\n", timage[i].name, 0.5*(RminImage + RmaxImage), 0.5*(DminImage + DmaxImage));
+      continue;
+    }
+    if (RmaxImage < RminBand[iDecBandMax]) {
+      fprintf (stderr, "skip image %s (%f,%f) on boundary\n", timage[i].name, 0.5*(RminImage + RmaxImage), 0.5*(DminImage + DmaxImage));
+      continue;
+    }
 
     // image overlaps region, keep it

branches/eam_branches/ipp-20120405/Ohana/src/relphot/src/setMrelFinal.c

@@ -66 +66 @@
           }
         
-          if (!KEEP_UBERCAL) {
-            catalog[0].measure[m].Mcal = 0;
-            catalog[0].measure[m].dbFlags &= ~ID_MEAS_PHOTOM_UBERCAL;
+          // only reset Mcal for measures with a matching image
+          // do not reset Mcal for ubercal images unless explicitly requested
+
+          if (catalog[0].measure[m].dbFlags & ID_MEAS_PHOTOM_UBERCAL) {
+            if (!KEEP_UBERCAL) {
+              catalog[0].measure[m].Mcal = 0.0;
+              catalog[0].measure[m].dbFlags &= ~ID_MEAS_PHOTOM_UBERCAL;
+            } 
+          } else {
+            if (getImageEntry (m, 0) >= 0) {
+              catalog[0].measure[m].Mcal = 0.0;
+            }
           }
           catalog[0].measure[m].dbFlags &= ~PHOTOM_FLAGS;