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

Timestamp:

Apr 20, 2009, 2:48:42 PM (17 years ago)

Author:

Paul Price

Message:

Implement stamp thresholds on both images. This is to avoid transients being selected as stamps.

File:

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

trunk/psModules/src/imcombine/pmSubtractionStamps.c

@@ -232 +232 @@
 
 
-pmSubtractionStampList *pmSubtractionStampsFind(pmSubtractionStampList *stamps, const psImage *image,
-                                                const psImage *subMask, const psRegion *region,
-                                                float threshold, int size, int footprint, float spacing,
+pmSubtractionStampList *pmSubtractionStampsFind(pmSubtractionStampList *stamps, const psImage *image1,
+                                                const psImage *image2, const psImage *subMask,
+                                                const psRegion *region, float thresh1, float thresh2,
+                                                int size, int footprint, float spacing,
                                                 pmSubtractionMode mode)
 {
-    PS_ASSERT_IMAGE_NON_NULL(image, NULL);
-    PS_ASSERT_IMAGE_TYPE(image, PS_TYPE_F32, NULL);
+    if (!image1 && !image2) {
+        psError(PS_ERR_UNEXPECTED_NULL, true, "Must specify an image");
+        return NULL;
+    }
+    int numCols = 0, numRows = 0;       // Size of images
+    if (image1) {
+        PS_ASSERT_IMAGE_NON_NULL(image1, NULL);
+        PS_ASSERT_IMAGE_TYPE(image1, PS_TYPE_F32, NULL);
+        if (subMask) {
+            PS_ASSERT_IMAGES_SIZE_EQUAL(image1, subMask, NULL);
+        }
+        numCols = image1->numCols;
+        numRows = image1->numRows;
+    }
+    if (image2) {
+        PS_ASSERT_IMAGE_NON_NULL(image2, NULL);
+        PS_ASSERT_IMAGE_TYPE(image2, PS_TYPE_F32, NULL);
+        if (subMask) {
+            PS_ASSERT_IMAGES_SIZE_EQUAL(image2, subMask, NULL);
+        }
+        numCols = image2->numCols;
+        numRows = image2->numRows;
+    }
+    if (image1 && image2) {
+        PS_ASSERT_IMAGES_SIZE_EQUAL(image1, image2, NULL);
+    }
     if (subMask) {
         PS_ASSERT_IMAGE_NON_NULL(subMask, NULL);
-        PS_ASSERT_IMAGES_SIZE_EQUAL(image, subMask, NULL);
         PS_ASSERT_IMAGE_TYPE(subMask, PS_TYPE_IMAGE_MASK, NULL);
+        PS_ASSERT_IMAGE_SIZE(subMask, numCols, numRows, NULL);
     }
     PS_ASSERT_INT_POSITIVE(size, NULL);
@@ -254 +279 @@
             return false;
         }
-        if (region->x0 < 0 || region->x1 > image->numCols ||
-            region->y0 < 0 || region->y1 > image->numRows) {
+        if (region->x0 < 0 || region->x1 > numCols ||
+            region->y0 < 0 || region->y1 > numRows) {
             psString string = psRegionToString(*region);
             psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Input region (%s) does not fit in image (%dx%d)",
-                    string, image->numCols, image->numRows);
+                    string, numCols, numRows);
             psFree(string);
             return false;
@@ -264 +289 @@
     }
 
-    int numRows = image->numRows, numCols = image->numCols; // Size of image
     int border = size + footprint;      // Border size
 
@@ -312 +336 @@
                     fluxList->n = index;
 
-                    goodStamp = (fluxStamp > threshold) ? true : false;
+                    goodStamp = true;
+                    int x = xStamp + 0.5, y = yStamp + 0.5; // Pixel location
+                    if (image1 && image1->data.F32[y][x] < thresh1) {
+                        goodStamp = false;
+                    } else if (image2 && image2->data.F32[y][x] < thresh2) {
+                        goodStamp = false;
+                    }
                 } else {
                     psTrace("psModules.imcombine", 9, "No sources in subregion %d", i);
@@ -319 +349 @@
                 // Use a simple method of automatically finding stamps --- take the highest pixel in the
                 // subregion
-                fluxStamp = threshold;
+                fluxStamp = 0.0;
                 psRegion *subRegion = stamps->regions->data[i]; // Sub-region of interest
+
+                const psImage *image = NULL; // Image for flux of stamp
+                switch (mode) {
+                  case PM_SUBTRACTION_MODE_1:
+                    image = image2 ? image2 : image1;
+                    break;
+                  case PM_SUBTRACTION_MODE_2:
+                    image = image1 ? image1 : image2;
+                    break;
+                  case PM_SUBTRACTION_MODE_UNSURE:
+                  case PM_SUBTRACTION_MODE_DUAL:
+                    // If both are available, we'll take the MIN value below in the hope of rejecting
+                    // transients in both
+                    image = (image1 && image2) ? NULL : (image1 ? image1 : image2);
+                  default:
+                    psAbort("Unrecognised mode: %x", mode);
+                }
+
                 for (int y = subRegion->y0; y <= subRegion->y1; y++) {
                     for (int x = subRegion->x0; x <= subRegion->x1; x++) {
-                        if (image->data.F32[y][x] > fluxStamp &&
+                        if ((image1 && image1->data.F32[y][x] < thresh1) ||
+                            (image2 && image2->data.F32[y][x] < thresh2)) {
+                            continue;
+                        }
+                        float value = image ? image->data.F32[y][x] :
+                            PS_MIN(image1->data.F32[y][x], image2->data.F32[y][x]); // Value of image
+                        if (value > fluxStamp &&
                             checkStampMask(x, y, numCols, numRows, subMask, mode, footprint, border)) {
                             fluxStamp = image->data.F32[y][x];