Changeset 30772 for branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c

Timestamp:

Mar 1, 2011, 6:21:28 PM (15 years ago)

Author:

eugene

Message:

fix threading for psphotStack; use standard PS_ARGUMENTS_ macros for ppSub, ppStack, pswarp, ppImage, psphotStack; fix FPA.EXPOSURE in ppImageAddNoise; fix accounting for psf bad and poor weights; added flags for substantial poor masks due to each poor mask bit type; fix the measurement of the fwhm in pmPSFtoFWHM

File:

• branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c (modified) (9 diffs)

branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c

@@ -23 +23 @@
 #include "pmFPAMaskWeight.h"
 
+#include "pmConfigMask.h"
 #include "pmTrend2D.h"
 #include "pmResiduals.h"
@@ -49 +50 @@
 static float AP_MIN_SN = 0.0;
 
-bool pmSourceMagnitudesInit (psMetadata *config)
-{
-    PS_ASSERT_PTR_NON_NULL(config, false);
+// make this a bit more clever and dynamic
+static psImageMaskType maskSuspect  = 0;
+static psImageMaskType maskSpike    = 0;
+static psImageMaskType maskStarCore = 0;
+static psImageMaskType maskBurntool = 0;
+static psImageMaskType maskConvPoor = 0;
+
+bool pmSourceMagnitudesInit (pmConfig *config, psMetadata *recipe)
+{
+    PS_ASSERT_PTR_NON_NULL(recipe, false);
     bool status;
 
-    float limit = psMetadataLookupF32 (&status, config, "AP_MIN_SN");
+    // we are going to test specially against these poor values
+    if (config) {
+        maskSpike    = pmConfigMaskGet("SPIKE", config);
+        maskStarCore = pmConfigMaskGet("STARCORE", config);
+        maskBurntool = pmConfigMaskGet("BURNTOOL", config);
+        maskConvPoor = pmConfigMaskGet("CONV.POOR", config);
+        maskSuspect  = maskSpike | maskStarCore | maskBurntool | maskConvPoor;
+    }
+
+    float limit = psMetadataLookupF32 (&status, recipe, "AP_MIN_SN");
     if (status) {
         AP_MIN_SN = limit;
@@ -77 +94 @@
 // XXX masked region should be (optionally) elliptical
 // if mode is PM_SOURCE_PHOT_PSFONLY, we skip all other magnitudes
-bool pmSourceMagnitudes (pmSource *source, pmPSF *psf, pmSourcePhotometryMode mode, psImageMaskType maskVal, psImageMaskType markVal)
+bool pmSourceMagnitudes (pmSource *source, pmPSF *psf, pmSourcePhotometryMode mode, psImageMaskType maskVal, psImageMaskType markVal, float radius)
 {
     PS_ASSERT_PTR_NON_NULL(source, false);
@@ -166 +183 @@
     // measure the contribution of included pixels
     if (mode & PM_SOURCE_PHOT_WEIGHT) {
-        pmSourcePixelWeight (&source->pixWeightNotBad, &source->pixWeightNotPoor, model, source->maskObj, maskVal, markVal);
+        pmSourcePixelWeight (source, model, source->maskObj, maskVal, radius);
     }
 
@@ -342 +359 @@
 
 // return source aperture magnitude
-bool pmSourcePixelWeight (float *pixWeightNotBad, float *pixWeightNotPoor, pmModel *model, psImage *mask, psImageMaskType maskVal, psImageMaskType markVal)
-{
-    PS_ASSERT_PTR_NON_NULL(pixWeightNotBad, false);
-    PS_ASSERT_PTR_NON_NULL(pixWeightNotPoor, false);
+bool pmSourcePixelWeight (pmSource *source, pmModel *model, psImage *mask, psImageMaskType maskVal, float radius)
+{
+    PS_ASSERT_PTR_NON_NULL(source, false);
+    source->pixWeightNotBad = NAN;
+    source->pixWeightNotPoor = NAN;
+
     PS_ASSERT_PTR_NON_NULL(mask, false);
     PS_ASSERT_PTR_NON_NULL(model, false);
@@ -355 +374 @@
     float value;
 
+    float spikeSum = 0;
+    float starcoreSum = 0;
+    float burntoolSum = 0;
+    float convpoorSum = 0;
+
     int Xo, Yo, dP;
     int dX, DX, NX;
     int dY, DY, NY;
 
-    *pixWeightNotBad = 0.0;
-    *pixWeightNotPoor = 0.0;
+    float radius2 = PS_SQR(radius);
 
     // we only care about the value of the object model, not the local sky
@@ -387 +410 @@
     NY = mask->numRows;
 
+    psImageMaskType maskBad = maskVal;
+    maskBad &= ~maskSuspect;
+
     // measure modelSum and validSum.  this function is applied to a sources' subimage.  the
     // value of DX is chosen (see above) to cover the full possible size of the subimage if it
     // were not by an edge; ie, if the source is cut in half by an image edge, we correctly
     // count the virtual pixels off the edge in normalizing the value of the pixWeight
+
+    // we skip any pixels [real or virtual] outside of the specified radius (nominally the aperture radius)
     for (int ix = -DX; ix < DX + 1; ix++) {
+        if (ix > radius) continue;
         int mx = ix + dX;
         for (int iy = -DY; iy < DY + 1; iy++) {
+            if (iy > radius) continue;
+            if (ix*ix + iy*iy > radius2) continue;
             int my = iy + dY;
 
@@ -409 +440 @@
             if (my >= NY) continue;
 
-            if (!(mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskVal)) {
+            // count pixels which are masked only with bad pixels
+            if (!(mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskBad)) {
                 notBadSum += value;
             }
-            if (!(mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & ~markVal)) {
+
+            // count pixels which are masked with an mask bit (bad or poor)
+            if (!(mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskVal)) {
                 notPoorSum += value;
             }
+
+            // count pixels which are masked with an mask bit (bad or poor)
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskSpike) {
+                spikeSum += value;
+            }
+            // count pixels which are masked with an mask bit (bad or poor)
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskStarCore) {
+                starcoreSum += value;
+            }
+            // count pixels which are masked with an mask bit (bad or poor)
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskBurntool) {
+                burntoolSum += value;
+            }
+            // count pixels which are masked with an mask bit (bad or poor)
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[my][mx] & maskConvPoor) {
+                convpoorSum += value;
+            }
         }
     }
     psFree (coord);
 
-    *pixWeightNotBad  = notBadSum  / modelSum;
-    *pixWeightNotPoor = notPoorSum / modelSum;
-
-    if (false && isfinite(*pixWeightNotBad) && isfinite(*pixWeightNotPoor)) {
-        psAssert (*pixWeightNotBad <= *pixWeightNotPoor, "error: all bad pixels should also be poor");
+    source->pixWeightNotBad  = notBadSum  / modelSum;
+    source->pixWeightNotPoor = notPoorSum / modelSum;
+
+    if ((spikeSum/modelSum) > 0.25) {
+        source->mode2 |= PM_SOURCE_MODE2_ON_SPIKE;
+    }
+    if ((starcoreSum/modelSum) > 0.25) {
+        source->mode2 |= PM_SOURCE_MODE2_ON_STARCORE;
+    }
+    if ((burntoolSum/modelSum) > 0.25) {
+        source->mode2 |= PM_SOURCE_MODE2_ON_BURNTOOL;
+    }
+    if ((convpoorSum/modelSum) > 0.25) {
+        source->mode2 |= PM_SOURCE_MODE2_ON_CONVPOOR;
+    }
+
+    if (isfinite(source->pixWeightNotBad) && isfinite(source->pixWeightNotPoor)) {
+        psAssert (source->pixWeightNotBad <= source->pixWeightNotPoor, "error: all bad pixels should also be poor");
     }
 
@@ -431 +495 @@
 # define FLUX_LIMIT 3.0
 
-// measure stats that may be using in difference images for distinguishing real sources from bad residuals
+// measure stats that may be used in difference images for distinguishing real sources from bad residuals
 bool pmSourceMeasureDiffStats (pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {