Index: trunk/psModules/src/imcombine/pmSubtractionStamps.c
===================================================================
--- trunk/psModules/src/imcombine/pmSubtractionStamps.c	(revision 14342)
+++ trunk/psModules/src/imcombine/pmSubtractionStamps.c	(revision 14455)
@@ -42,5 +42,5 @@
 
 psArray *pmSubtractionFindStamps(psArray *stamps, const psImage *image, const psImage *subMask,
-                                 float threshold, float spacing)
+                                 const psRegion *region, float threshold, float spacing)
 {
     PS_ASSERT_IMAGE_NON_NULL(image, NULL);
@@ -52,12 +52,36 @@
     }
     PS_ASSERT_FLOAT_LARGER_THAN(spacing, 0.0, NULL);
+    if (region) {
+        if (psRegionIsNaN(*region)) {
+            psString string = psRegionToString(*region);
+            psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Input region (%s) contains NAN values", string);
+            psFree(string);
+            return false;
+        }
+        if (region->x0 < 0 || region->x1 > image->numCols ||
+            region->y0 < 0 || region->y1 > image->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);
+            psFree(string);
+            return false;
+        }
+    }
 
-    // Size of image
-    int numRows = image->numRows;
-    int numCols = image->numCols;
+    int numRows = image->numRows, numCols = image->numCols; // Size of image
+
+    // Get region in which to find stamps: [xMin:xMax,yMin:yMax]
+    int xMin = 0, xMax = numCols, yMin = 0, yMax = numRows;
+    if (region) {
+        xMin = PS_MAX(region->x0, xMin);
+        xMax = PS_MIN(region->x1, xMax);
+        yMin = PS_MAX(region->y0, yMin);
+        yMax = PS_MIN(region->y1, yMax);
+    }
+    int xSize = xMax - xMin, ySize = yMax - yMin; // Size of region of interest
 
     // Number of stamps
-    int xNumStamps = numCols / spacing + 1;
-    int yNumStamps = numRows / spacing + 1;
+    int xNumStamps = xSize / spacing + 1;
+    int yNumStamps = ySize / spacing + 1;
 
     int numFound = 0;                   // Number of stamps found
@@ -92,12 +116,12 @@
 
                 // Bounds of region to search for stamp
-                int yMin = j * numRows / (yNumStamps + 1);
-                int yMax = (j + 1) * numRows / (yNumStamps + 1) - 1;
-                int xMin = i * numCols / (xNumStamps + 1);
-                int xMax = (i + 1) * numCols / (xNumStamps + 1) - 1;
-                assert(yMax < image->numRows && xMax < image->numCols && yMin >= 0 && xMin >= 0);
+                int yStart = yMin + j * ySize / (yNumStamps + 1);
+                int yStop = yMin + (j + 1) * ySize / (yNumStamps + 1) - 1;
+                int xStart = xMin + i * xSize / (xNumStamps + 1);
+                int xStop = xMin + (i + 1) * xSize / (xNumStamps + 1) - 1;
+                assert(yStop < numRows && xStop < numCols && yStart >= 0 && xStart >= 0);
 
-                for (int y = yMin; y <= yMax ; y++) {
-                    for (int x = xMin; x <= xMax ; x++) {
+                for (int y = yStart; y <= yStop ; y++) {
+                    for (int x = xStart; x <= xStop ; x++) {
                         if ((!subMask || !(subMask->data.PS_TYPE_MASK_DATA[y][x] &
                                            (PM_SUBTRACTION_MASK_BORDER | PM_SUBTRACTION_MASK_FOOTPRINT |