Changeset 32850

Timestamp:

Dec 2, 2011, 11:42:55 AM (15 years ago)

Author:

eugene

Message:

interpolate the sky-limit radius to 0.0 flux; allow faint sources to have a somewhat larger max size; ignore DEFECT and SATURATION (garbage) in masking the neighbors for kron mags

Location:

branches/eam_branches/ipp-20111122/psphot

Files:

• . (modified) (1 prop)
• src (modified) (1 prop)
• src/psphotKronIterate.c (modified) (3 diffs)
• src/psphotRadialProfileWings.c (modified) (4 diffs)
• src/psphotStackImageLoop.c (modified) (1 prop)

branches/eam_branches/ipp-20111122/psphot/src/psphotKronIterate.c

@@ -111 +111 @@
 
         pmSource *source = sources->data[i];
-        if (!source->peak) continue; // XXX how can we have a peak-less source?
-
-        // allocate space for moments
-        if (!source->moments) continue;
 
         // set a window function for each source based on the moments
+        // (this skips really bad sources (no peak, no moments, DEFECT)
         psphotKronWindowSetSource (source, kronWindow, false, true);
     }
@@ -221 +218 @@
 
             // XXX TEST : use a window based on the radial profile numbers: max is skyRadius, min is RADIUS
-            float maxWindow = source->skyRadius;
+            // if we lack the skyRadius (eg MATCHED sources), go to the default value
+            float maxWindow = isfinite(source->skyRadius) ? source->skyRadius : RADIUS;
             float windowRadius = PS_MAX(RADIUS, maxWindow);
 
@@ -400 +398 @@
 
     if (!source) return false;
+    if (!source->peak) return false; // XXX how can we have a peak-less source?
     if (!source->moments) return false;
+    if (source->type == PM_SOURCE_TYPE_DEFECT) return false;
+    if (source->type == PM_SOURCE_TYPE_SATURATED) return false;
     psAssert(kronWindow, "need a window");
 

branches/eam_branches/ipp-20111122/psphot/src/psphotRadialProfileWings.c

@@ -80 +80 @@
     }
 
-    MIN_RADIUS = 0.25*psMetadataLookupF32 (&status, readout->analysis, "PSF_MOMENTS_RADIUS");
+    MIN_RADIUS = psMetadataLookupF32 (&status, readout->analysis, "PSF_MOMENTS_RADIUS");
     if (!status) {
-        MIN_RADIUS = 0.25*psMetadataLookupF32 (&status, recipe, "PSF_MOMENTS_RADIUS");
+        MIN_RADIUS = psMetadataLookupF32 (&status, recipe, "PSF_MOMENTS_RADIUS");
     }
 
@@ -209 +209 @@
 # define TEST_Y 3096
 
+float InterpolateValues (float X0, float Y0, float X1, float Y1, float X);
+
 // XXX use integer radius values?  the rings assume integer values, right? or do they?
 bool psphotRadialProfileWingsSource (pmSource *source, pmReadout *readout, psImageMaskType maskVal) {
@@ -223 +225 @@
     if (isfinite(source->skyRadius)) return true;
 
+    if (source->mode2 & PM_SOURCE_MODE2_MATCHED) return true; // skip matched sources (no signal)
+
     // radii will be MIN_RADIUS to MAX_RADIUS in NN log steps:
-    float NSTEP = 10.0;
+    float NSTEP = 25.0;
     float MIN_DR = 2;
     float NSIGMA = 1.0;
@@ -286 +290 @@
         // fprintf (stderr, "%f %f : %f : %f %f  :  %f\n", source->peak->xf, source->peak->yf, radius, meanFlux, meanFluxError, slope);
 
-        if (!limit) {
-            limit |= (meanFlux - NSIGMA * meanFluxError < THRESHOLD); // dropped to sky level
-            limit |= isfinite(slope) && (fabs(slope) < 3.0); // SB no longer changing.
-            limitRadius = meanRadius;
+        if (!limit && (meanFlux - NSIGMA * meanFluxError < THRESHOLD)) {
+            // dropped to sky level
+            limit = true;
+            // linearly interpolate to the radius at which we hit the sky
+            if (isfinite(lastFlux)) {
+                limitRadius = InterpolateValues(lastFlux, lastRadius, meanFlux, meanRadius, 0.0);
+            } else {
+                limitRadius = meanRadius;
+            }
+            limitFlux = meanFlux;
+            limitSlope = slope;
+        }
+        if (!limit && isfinite(slope) && (fabs(slope) < 3.0)) { 
+            // SB no longer changing.       
+            limit = true;
+            // linearly interpolate to the radius at which we hit the sky, using the last flux and the limiting slope
+            if (isfinite(lastFlux)) {
+                limitRadius = lastRadius + lastFlux / 3.0;
+            } else {
+                limitRadius = meanRadius;
+            }
             limitFlux = meanFlux;
             limitSlope = slope;