Changeset 23849 for trunk/psModules/src/imcombine/pmSubtractionStamps.c

Timestamp:

Apr 14, 2009, 8:53:51 AM (17 years ago)

Author:

Paul Price

Message:

Fixing how stamps are found: ensuring NANs don't get into the stamp footprint, but trying not to make it go slow when looking for the maximum pixel in a subregion.

File:

• trunk/psModules/src/imcombine/pmSubtractionStamps.c (modified) (14 diffs)

trunk/psModules/src/imcombine/pmSubtractionStamps.c

@@ -80 +80 @@
 // Is this position unmasked?
 static bool checkStampMask(int x, int y, // Coordinates of stamp
+                           int numCols, int numRows, // Size of image
                            const psImage *mask, // Mask
                            pmSubtractionMode mode, // Mode for subtraction
-                           int footprint // Footprint to check for Bad Stuff
+                           int footprint, // Size of stamp
+                           int border // Size of border
                            )
 {
+    // Check the footprint bounds
+    if (x < border || x >= numCols - border || y < border || y >= numRows - border) {
+        if (mask) {
+            goto BAD_STAMP;
+        }
+        return false;
+    }
+
     if (!mask) {
         return true;
     }
 
-    bool clean = true;                  // Is the footprint clean?
-    int numCols = mask->numCols, numRows = mask->numRows; // Size of image
-
-    int xMin = PS_MAX(x - footprint, 0), xMax = PS_MIN(x + footprint, numCols - 1); // Bounds in x
-    int yMin = PS_MAX(y - footprint, 0), yMax = PS_MIN(y + footprint, numRows - 1); // Bounds in y
-
-    // Check the footprint bounds
-    if (x < footprint || x >= numCols - footprint || y < footprint || y >= numRows - footprint) {
-        clean = false;
-        goto CHECK_STAMP_MASK_DONE;
-
-    }
-
     // Determine mask value
-    psImageMaskType maskVal = PM_SUBTRACTION_MASK_BORDER | PM_SUBTRACTION_MASK_BAD_1 | PM_SUBTRACTION_MASK_BAD_2;
+    psImageMaskType maskVal = PM_SUBTRACTION_MASK_BAD_1 | PM_SUBTRACTION_MASK_BAD_2 |
+        PM_SUBTRACTION_MASK_BORDER;     // Mask value to check for
     switch (mode) {
       case PM_SUBTRACTION_MODE_1:
-        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_1;
+        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_1 | PM_SUBTRACTION_MASK_CONVOLVE_BAD_1;
         break;
       case PM_SUBTRACTION_MODE_2:
-        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_2;
+        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_2 | PM_SUBTRACTION_MASK_CONVOLVE_BAD_2;
         break;
       case PM_SUBTRACTION_MODE_UNSURE:
       case PM_SUBTRACTION_MODE_DUAL:
-        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_1 | PM_SUBTRACTION_MASK_CONVOLVE_2;
+        maskVal |= PM_SUBTRACTION_MASK_CONVOLVE_1 | PM_SUBTRACTION_MASK_CONVOLVE_BAD_1 |
+            PM_SUBTRACTION_MASK_CONVOLVE_2 | PM_SUBTRACTION_MASK_CONVOLVE_BAD_2;
         break;
       default:
@@ -119 +118 @@
     }
 
-    // Check the immediate pixel
-    if (clean && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & (maskVal | PM_SUBTRACTION_MASK_REJ))) {
-        goto CHECK_STAMP_MASK_DONE;
-        clean = false;
+    // Check the central pixel
+    if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & (maskVal | PM_SUBTRACTION_MASK_REJ)) {
+        goto BAD_STAMP;
     }
 
     // Check the footprint
-    if (clean) {
-        for (int j = yMin; j <= yMax; j++) {
-            for (int i = xMin; i <= xMax; i++) {
-                if (mask->data.PS_TYPE_IMAGE_MASK_DATA[j][i] & maskVal) {
-                    clean = false;
-                    goto CHECK_STAMP_MASK_DONE;
-                }
+    int xMin = PS_MAX(x - footprint, 0), xMax = PS_MIN(x + footprint, numCols - 1); // Bounds in x
+    int yMin = PS_MAX(y - footprint, 0), yMax = PS_MIN(y + footprint, numRows - 1); // Bounds in y
+    for (int j = yMin; j <= yMax; j++) {
+        for (int i = xMin; i <= xMax; i++) {
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[j][i] & maskVal) {
+                goto BAD_STAMP;
             }
         }
     }
 
- CHECK_STAMP_MASK_DONE:
-    if (!clean) {
-        // Mask the footprint, so we don't select something near it again
-        for (int j = yMin; j <= yMax; j++) {
-            for (int i = xMin; i <= xMax; i++) {
-                mask->data.PS_TYPE_IMAGE_MASK_DATA[j][i] |= PM_SUBTRACTION_MASK_REJ;
-            }
-        }
-    }
-
-    return clean;
+    return true;
+
+BAD_STAMP:
+    // Stamp is bad, so mark it as such so we don't find it again
+    mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] |= PM_SUBTRACTION_MASK_REJ;
+    return false;
 }
 
@@ -279 +271 @@
 pmSubtractionStampList *pmSubtractionStampsFind(pmSubtractionStampList *stamps, const psImage *image,
                                                 const psImage *subMask, const psRegion *region,
-                                                float threshold, int footprint, float spacing,
+                                                float threshold, int size, int footprint, float spacing,
                                                 pmSubtractionMode mode)
 {
@@ -289 +281 @@
         PS_ASSERT_IMAGE_TYPE(subMask, PS_TYPE_IMAGE_MASK, NULL);
     }
+    PS_ASSERT_INT_POSITIVE(size, NULL);
     PS_ASSERT_INT_NONNEGATIVE(footprint, NULL);
     PS_ASSERT_FLOAT_LARGER_THAN(spacing, 0.0, NULL);
@@ -309 +302 @@
 
     int numRows = image->numRows, numCols = image->numCols; // Size of image
+    int border = size + footprint;      // Border size
 
     if (!stamps) {
@@ -367 +361 @@
                     for (int x = subRegion->x0; x <= subRegion->x1; x++) {
                         if (image->data.F32[y][x] > fluxStamp &&
-                            checkStampMask(x, y, subMask, mode, footprint)) {
+                            checkStampMask(x, y, numCols, numRows, subMask, mode, footprint, border)) {
                             fluxStamp = image->data.F32[y][x];
                             xStamp = x;
@@ -419 +413 @@
 pmSubtractionStampList *pmSubtractionStampsSet(const psVector *x, const psVector *y,
                                                const psImage *image, const psImage *subMask,
-                                               const psRegion *region, int footprint, float spacing,
-                                               pmSubtractionMode mode)
+                                               const psRegion *region, int size, int footprint,
+                                               float spacing, pmSubtractionMode mode)
 
 {
@@ -432 +426 @@
         PS_ASSERT_IMAGE_NON_NULL(subMask, NULL);
         PS_ASSERT_IMAGE_TYPE(subMask, PS_TYPE_IMAGE_MASK, NULL);
-        if (image) {
-            PS_ASSERT_IMAGE_NON_NULL(image, NULL);
-            PS_ASSERT_IMAGES_SIZE_EQUAL(image, subMask, NULL);
-        }
-    }
+        PS_ASSERT_IMAGES_SIZE_EQUAL(image, subMask, NULL);
+    }
+    PS_ASSERT_INT_POSITIVE(size, NULL);
+    PS_ASSERT_INT_POSITIVE(footprint, NULL);
     PS_ASSERT_FLOAT_LARGER_THAN(spacing, 0.0, NULL);
 
@@ -443 +436 @@
                                                                  region, footprint, spacing); // Stamp list
     int numStamps = stamps->num;        // Number of stamps
+    int numCols = image->numCols, numRows = image->numRows; // Size of image
+    int border = footprint + size;      // Border size
 
     psString ds9name = NULL;            // Filename for ds9 region file
@@ -471 +466 @@
             continue;
         }
-        if (!checkStampMask(xPix, yPix, subMask, mode, footprint)) {
+        if (!checkStampMask(xPix, yPix, numCols, numRows, subMask, mode, footprint, border)) {
             // Not a good stamp
             psTrace("psModules.imcombine", 9, "Rejecting input stamp (%d,%d) because bad mask",
@@ -613 +608 @@
 pmSubtractionStampList *pmSubtractionStampsSetFromSources(const psArray *sources, const psImage *image,
                                                           const psImage *subMask, const psRegion *region,
-                                                          int footprint, float spacing,
+                                                          int size, int footprint, float spacing,
                                                           pmSubtractionMode mode)
 {
@@ -642 +637 @@
     y->n = numOut;
 
-    pmSubtractionStampList *stamps = pmSubtractionStampsSet(x, y, image, subMask, region,
+    pmSubtractionStampList *stamps = pmSubtractionStampsSet(x, y, image, subMask, region, size,
                                                             footprint, spacing, mode); // Stamps to return
     psFree(x);
@@ -657 +652 @@
 pmSubtractionStampList *pmSubtractionStampsSetFromFile(const char *filename, const psImage *image,
                                                        const psImage *subMask, const psRegion *region,
-                                                       int footprint, float spacing, pmSubtractionMode mode)
+                                                       int size, int footprint, float spacing,
+                                                       pmSubtractionMode mode)
 {
     PS_ASSERT_STRING_NON_EMPTY(filename, NULL);
@@ -673 +669 @@
     psBinaryOp(y, y, "-", psScalarAlloc(1.0, PS_TYPE_F32));
 
-    pmSubtractionStampList *stamps = pmSubtractionStampsSet(x, y, image, subMask, region, footprint,
+    pmSubtractionStampList *stamps = pmSubtractionStampsSet(x, y, image, subMask, region, size, footprint,
                                                             spacing, mode);
     psFree(data);