Changeset 30554

Timestamp:

Feb 10, 2011, 10:39:48 AM (15 years ago)

Author:

eugene

Message:

speed up image selection; add dMsys, nFitPhotom info; plug a leak

Location:

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

Files:

• include/relphot.h (modified) (3 diffs)
• src/ConfigInit.c (modified) (1 diff)
• src/ImageOps.c (modified) (5 diffs)
• src/MosaicOps.c (modified) (3 diffs)
• src/args.c (modified) (1 diff)
• src/select_images.c (modified) (6 diffs)

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

@@ -13 +13 @@
   float Mcal;
   float dMcal;
+  float dMsys;
+  short nFitPhotom;
   short Xm;
   float secz;
@@ -50 +52 @@
 double STAR_CHISQ;
 double MIN_ERROR;
+double IMFIT_SYS_SIGMA_LIM;
 
 int    VERBOSE;
@@ -103 +106 @@
 SkyRegion UserPatch;
 int UserPatchSelect;
+
+int USE_BASIC_CHECK;
 
 # ifdef GRID_V1

branches/eam_branches/ipp-20101205/Ohana/src/relphot/src/ConfigInit.c

@@ -20 +20 @@
   GetConfig (config, "MAG_LIM",                "%lf", 0, &MAG_LIM);
   GetConfig (config, "SIGMA_LIM",              "%lf", 0, &SIGMA_LIM);
+
+  if (!ScanConfig (config, "RELPHOT_IMFIT_SYS_SIGMA_LIM", "%lf", 0, &IMFIT_SYS_SIGMA_LIM)) IMFIT_SYS_SIGMA_LIM = 0.01;
+
   GetConfig (config, "STAR_SCATTER",           "%lf", 0, &STAR_SCATTER);
   GetConfig (config, "IMAGE_SCATTER",          "%lf", 0, &IMAGE_SCATTER);

branches/eam_branches/ipp-20101205/Ohana/src/relphot/src/ImageOps.c

@@ -124 +124 @@
   free (clist);
   free (mlist);
+  free (Nlist);
+  free (NLIST);
 }
 
@@ -299 +301 @@
 
   off_t i, j, m, c, n, N, Nmax, mark, bad;
-  int Nfew, Nbad, Nmos, Ngrid, Nrel, Nsys;
-  float Msys, Mrel, Mmos, Mgrid;
+  int Nfew, Nbad, Nmos, Ngrid, Nrel, Nsys, Nbright;
+  float Msys, Mrel, Mmos, Mgrid, McalBright, McalBright2;
   double *list, *dlist;
   StatType stats;
@@ -318 +320 @@
 
   Nfew = Nbad = Nmos = Ngrid = Nrel = Nsys = 0;
+
+  // counters to measure the bright-end scatter
+  McalBright = McalBright2 = 0.0;
+  Nbright = 0;
 
   for (i = 0; i < Nimage; i++) {
@@ -361 +367 @@
       list[N] = Msys - Mrel - Mmos - Mgrid;
       dlist[N] = MAX (catalog[c].measure[m].dM, MIN_ERROR);
+      if (catalog[c].measure[m].dM < IMFIT_SYS_SIGMA_LIM) {
+        McalBright += list[N];
+        McalBright2 += SQ(list[N]);
+        Nbright ++;
+      }
       N++;
     }
@@ -378 +389 @@
     liststats (list, dlist, N, &stats);
     image[i].Mcal  = stats.mean;
-    image[i].dMcal = stats.sigma;
+    image[i].dMcal = stats.error;
+    image[i].dMagSys = stats.sigma;
+    image[i].nFitPhotom = N;
     image[i].Xm    = 100.0*log10(stats.chisq);
   }

branches/eam_branches/ipp-20101205/Ohana/src/relphot/src/MosaicOps.c

@@ -73 +73 @@
     mosaic[Nmosaic].Mcal  = 0.0;
     mosaic[Nmosaic].dMcal = 0.0;
+    mosaic[Nmosaic].dMsys = 0.0;
     mosaic[Nmosaic].Xm    = 0.0;
     mosaic[Nmosaic].flags  = image[i].flags;
@@ -189 +190 @@
       image[im].dMcal = mosaic[i].dMcal;
       image[im].Xm = mosaic[i].Xm;
+      image[im].dMagSys = mosaic[i].dMsys;
+      image[im].nFitPhotom = mosaic[i].nFitPhotom;
       image[im].flags |= (mosaic[i].flags & ID_IMAGE_PHOTOM_FEW);
       image[im].flags |= (mosaic[i].flags & ID_IMAGE_PHOTOM_POOR);
@@ -410 +413 @@
     if (PoorImages) fprintf (stderr, "Mmos: %f %f %d "OFF_T_FMT"\n", stats.mean, stats.sigma, stats.Nmeas,  N);
     mosaic[i].Mcal  = stats.mean;
-    mosaic[i].dMcal = stats.sigma;
+    mosaic[i].dMcal = stats.error;
+    mosaic[i].dMsys = stats.sigma;
+    mosaic[i].nFitPhotom = N;
     mosaic[i].Xm    = 100.0*log10(stats.chisq);
   }

branches/eam_branches/ipp-20101205/Ohana/src/relphot/src/args.c

@@ -52 +52 @@
   }
 
+  USE_BASIC_CHECK = FALSE;
+  if ((N = get_argument (argc, argv, "-basic-image-search"))) {
+    remove_argument (N, &argc, argv);
+    USE_BASIC_CHECK = TRUE;
+  }
+
   VERBOSE = VERBOSE2 = FALSE;
   if ((N = get_argument (argc, argv, "-v"))) {

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

@@ -33 +33 @@
   struct timeval start, stop;
   
-  double RmaxSkyRegion, RminSkyRegion, DminSkyRegion, DmaxSkyRegion;
+  double RmaxSkyRegion, RminSkyRegion, DminSkyRegion, DmaxSkyRegion, RmidSkyRegion;
 
   double *RmaxSky;
@@ -97 +97 @@
     DmaxSkyRegion = MAX(DmaxSkyRegion, skylist[0].regions[i][0].Dmax);
   }
+  RmidSkyRegion = 0.5*(RminSkyRegion + RmaxSkyRegion);
   MARKTIME("create sky region coords: %f sec\n", dtime);
 
@@ -145 +146 @@
     for (j = 0; j < 5; j++) {
       XY_to_RD (&Ri[j], &Di[j], Xi[j], Yi[j], &timage[i].coords);
+      Ri[j] = ohana_normalize_angle_to_midpoint (Ri[j], RmidSkyRegion);
+      
       RminImage = MIN(RminImage, Ri[j]);
       RmaxImage = MAX(RmaxImage, Ri[j]);
@@ -164 +167 @@
     if (RmaxImage < RminSkyRegion) continue;
 
+    // image overlaps region, keep it
+    if (USE_BASIC_CHECK) goto found_it;
+
     // RA(nStart) is guaranteed to be < RminImage: -- costs 0.5sec for 3M images
     nStart = getRegionStartByRA (RminImage, RmaxSky, skylist[0].Nregions);
@@ -186 +192 @@
       for (j = 0; (j < 4) && !found; j++) {
         found = corner_check (&Xi[j], &Yi[j], &skycoords[m].Xc[0], &skycoords[m].Yc[0]);
+        if (found) goto found_it;
       }
       /* check if catalog corner inside image */
       for (j = 0; (j < 4) && !found; j++) {
         found = corner_check (&skycoords[m].Xc[j], &skycoords[m].Yc[j], &Xi[0], &Yi[0]);
+        if (found) goto found_it;
       }
       /* check if edges cross */
@@ -195 +203 @@
         for (k = 0; (k < 4) && !found; k++) {
           found = edge_check (&Xi[j], &Yi[j], &skycoords[m].Xc[k], &skycoords[m].Yc[k]);
+          if (found) goto found_it;
         }
       }
-      if (!found) continue;
-
-      image[nimage] = timage[i]; 
-      /* always allow 'few' images to succeed, if possible */
-      if (image[nimage].flags & ID_IMAGE_PHOTOM_FEW) { 
-        image[nimage].flags &= ~(ID_IMAGE_PHOTOM_FEW | ID_IMAGE_PHOTOM_POOR);
-      }
-      if (RESET) {
-        // XXX assignMcal (&image[nimage], (double *) NULL, -1);
-        // XXX this needs to be thought through a bit more
-        image[nimage].Mcal = 0.0;
-        image[nimage].dMcal = NAN;
-        image[nimage].flags &= ~ID_IMAGE_PHOTOM_POOR;
-      }
-      line_number[nimage] = i;
-      nimage ++;
-      if (nimage == NIMAGE) {
-        NIMAGE += 100;
-        REALLOCATE (image, Image, NIMAGE);
-        REALLOCATE (line_number, off_t, NIMAGE);
-      }
+    }
+    if (!found) continue;
+
+  found_it:
+    image[nimage] = timage[i]; 
+    /* always allow 'few' images to succeed, if possible */
+    if (image[nimage].flags & ID_IMAGE_PHOTOM_FEW) { 
+      image[nimage].flags &= ~(ID_IMAGE_PHOTOM_FEW | ID_IMAGE_PHOTOM_POOR);
+    }
+    if (RESET) {
+      // XXX assignMcal (&image[nimage], (double *) NULL, -1);
+      // XXX this needs to be thought through a bit more
+      image[nimage].Mcal = 0.0;
+      image[nimage].dMcal = NAN;
+      image[nimage].flags &= ~ID_IMAGE_PHOTOM_POOR;
+    }
+    line_number[nimage] = i;
+    nimage ++;
+    if (nimage == NIMAGE) {
+      NIMAGE += 100;
+      REALLOCATE (image, Image, NIMAGE);
+      REALLOCATE (line_number, off_t, NIMAGE);
     }
   }