Index: trunk/psModules/src/detrend/pmMaskBadPixels.c
===================================================================
--- trunk/psModules/src/detrend/pmMaskBadPixels.c	(revision 16413)
+++ trunk/psModules/src/detrend/pmMaskBadPixels.c	(revision 17228)
@@ -81,50 +81,45 @@
 
 
-psImage *pmMaskFlagSuspectPixels(psImage *out, const pmReadout *readout, float rej,
-                                 psMaskType maskVal)
+bool pmMaskFlagSuspectPixels(pmReadout *output, const pmReadout *readout, float median, float stdev,
+                             float rej, psMaskType maskVal)
 {
-    PS_ASSERT_PTR_NON_NULL(readout, NULL);
-    PS_ASSERT_FLOAT_LARGER_THAN(rej, 0.0, NULL);
-    PS_ASSERT_IMAGE_NON_NULL(readout->image, NULL);
-    PS_ASSERT_IMAGE_NON_EMPTY(readout->image, NULL);
-    PS_ASSERT_IMAGE_TYPE(readout->image, PS_TYPE_F32, NULL);
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+    PS_ASSERT_FLOAT_LARGER_THAN(rej, 0.0, false);
+    PS_ASSERT_IMAGE_NON_NULL(readout->image, false);
+    PS_ASSERT_IMAGE_NON_EMPTY(readout->image, false);
+    PS_ASSERT_IMAGE_TYPE(readout->image, PS_TYPE_F32, false);
     if (readout->mask) {
-        PS_ASSERT_IMAGE_NON_EMPTY(readout->mask, NULL);
-        PS_ASSERT_IMAGES_SIZE_EQUAL(readout->image, readout->mask, NULL);
-        PS_ASSERT_IMAGE_TYPE(readout->mask, PS_TYPE_MASK, NULL);
-    }
-    if (out) {
-        PS_ASSERT_IMAGE_NON_EMPTY(out, NULL);
-        PS_ASSERT_IMAGE_TYPE(out, PS_TYPE_S32, NULL);
-        PS_ASSERT_IMAGES_SIZE_EQUAL(readout->image, out, NULL);
-    }
-
-    bool status;
+        PS_ASSERT_IMAGE_NON_EMPTY(readout->mask, false);
+        PS_ASSERT_IMAGES_SIZE_EQUAL(readout->image, readout->mask, false);
+        PS_ASSERT_IMAGE_TYPE(readout->mask, PS_TYPE_MASK, false);
+    }
+    PS_ASSERT_PTR_NON_NULL(output, false);
+
+    bool mdok;                          // Status of MD lookup
+    psImage *suspect = psMetadataLookupPtr(&mdok, output->analysis, PM_MASK_ANALYSIS_SUSPECT); // Suspect img
+    if (suspect) {
+        PS_ASSERT_IMAGE_NON_EMPTY(suspect, false);
+        PS_ASSERT_IMAGE_TYPE(suspect, PS_TYPE_S32, false);
+        PS_ASSERT_IMAGES_SIZE_EQUAL(readout->image, suspect, false);
+        psMemIncrRefCounter(suspect);
+    } else {
+        suspect = psImageAlloc(readout->image->numCols, readout->image->numRows, PS_TYPE_S32);
+        psImageInit(suspect, 0);
+        psMetadataAddImage(output->analysis, PS_LIST_TAIL, PM_MASK_ANALYSIS_SUSPECT, PS_META_REPLACE,
+                           "Suspect pixels", suspect);
+        psMetadataAddS32(output->analysis, PS_LIST_TAIL, PM_MASK_ANALYSIS_NUM, PS_META_REPLACE,
+                         "Number of input images", 0);
+    }
+
+    if (!isfinite(median) || !isfinite(stdev)) {
+        // If we get down here and the statistics are missing, then we should go and mask the entire image
+        psWarning("Missing statistics --- flagging entire image as suspect.");
+        return (psImage*)psBinaryOp(suspect, suspect, "+", psScalarAlloc(1.0, PS_TYPE_S32));
+    }
+
     psImage *image = readout->image;    // Image of interest
     psImage *mask = readout->mask;      // Corresponding mask
 
-    if (!out) {
-        out = psImageAlloc(image->numCols, image->numRows, PS_TYPE_S32);
-        psImageInit(out, 0);
-    }
-
-    bool whole = false;                 // Mask whole image?
-    float median = psMetadataLookupF32 (&status, readout->analysis, "READOUT.MEDIAN");
-    if (!status || !isfinite(median)) {
-        whole = true;
-    }
-    float stdev  = psMetadataLookupF32 (&status, readout->analysis, "READOUT.STDEV");
-    if (!status || !isfinite(stdev)) {
-        whole = true;
-    }
-
-
     psTrace ("psModules.detrend", 3, "suspect: %f +/- %f\n", median, stdev);
-
-    if (whole) {
-        // If we get down here and the statistics are missing, then we should go and mask the entire image
-        psWarning("Missing statistics --- flagging entire image as suspect.");
-        return (psImage*)psBinaryOp(out, out, "+", psScalarAlloc(1.0, PS_TYPE_S32));
-    }
 
     for (int y = 0; y < image->numRows; y++) {
@@ -132,20 +127,38 @@
             if (fabs((image->data.F32[y][x] - median) / stdev) >= rej &&
                     (!mask || !(mask->data.PS_TYPE_MASK_DATA[y][x] & maskVal))) {
-                out->data.S32[y][x]++;
-            }
-        }
-    }
-
-    return out;
-}
-
-psImage *pmMaskIdentifyBadPixels(const psImage *suspects, psMaskType maskVal, int nTotal, float thresh, pmMaskIdentifyMode mode)
+                suspect->data.S32[y][x]++;
+            }
+        }
+    }
+    psFree(suspect);                    // Drop reference
+
+    psMetadataItem *numItem = psMetadataLookup(output->analysis, PM_MASK_ANALYSIS_NUM); // Item with number
+    assert(numItem);
+    numItem->data.S32++;
+
+    return true;
+}
+
+bool pmMaskIdentifyBadPixels(pmReadout *output, psMaskType maskVal, float thresh, pmMaskIdentifyMode mode)
 {
-    PS_ASSERT_IMAGE_NON_NULL(suspects, NULL);
-    PS_ASSERT_IMAGE_NON_EMPTY(suspects, NULL);
-    PS_ASSERT_IMAGE_TYPE(suspects, PS_TYPE_S32, NULL);
+    PS_ASSERT_PTR_NON_NULL(output, false);
+    psImage *suspects = psMetadataLookupPtr(NULL, output->analysis, PM_MASK_ANALYSIS_SUSPECT); // Suspect img
+    if (!suspects) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to find image with suspected bad pixels.");
+        return false;
+    }
+    PS_ASSERT_IMAGE_NON_EMPTY(suspects, false);
+    PS_ASSERT_IMAGE_TYPE(suspects, PS_TYPE_S32, false);
+    if (output->mask) {
+        PS_ASSERT_IMAGE_NON_EMPTY(output->mask, false);
+        PS_ASSERT_IMAGES_SIZE_EQUAL(output->mask, suspects, false);
+        PS_ASSERT_IMAGE_TYPE(output->mask, PS_TYPE_MASK, false);
+    } else {
+        output->mask = psImageAlloc(suspects->numCols, suspects->numRows, PS_TYPE_MASK);
+    }
+    int num = psMetadataLookupS32(NULL, output->analysis, PM_MASK_ANALYSIS_NUM); // Number of inputs
+    PS_ASSERT_INT_POSITIVE(num, false);
 
     float limit = NAN;                  // Limit for masking
-
     switch (mode) {
       case PM_MASK_ID_VALUE:
@@ -154,5 +167,5 @@
 
       case PM_MASK_ID_FRACTION:
-        limit = thresh * nTotal;
+        limit = thresh * num;
         break;
 
@@ -200,5 +213,6 @@
         limit = max + 1.0 - thresh * sqrtf((float)max + 1.0);
 
-        psTrace ("psModules.detrend", 3, "bad: mode: %d, stdev: %f, limit: %f\n", max, sqrtf((float)max + 1.0), limit);
+        psTrace ("psModules.detrend", 3, "bad: mode: %d, stdev: %f, limit: %f\n",
+                 max, sqrtf((float)max + 1.0), limit);
         break;
       }
@@ -207,7 +221,4 @@
         return NULL;
     }
-
-    psImage *badpix = psImageAlloc(suspects->numCols, suspects->numRows, PS_TYPE_MASK); // Bad pixel mask
-    psImageInit(badpix, 0);
 
     if (psTraceGetLevel("psModules.detrend") > 9) {
@@ -227,4 +238,7 @@
     psTrace ("psModules.detrend", 3, "bad pixel threshold: %f", limit);
 
+    psImage *badpix = output->mask;     // Bad pixel mask
+    psImageInit(badpix, 0);
+
     for (int y = 0; y < suspects->numRows; y++) {
         for (int x = 0; x < suspects->numCols; x++) {
@@ -235,19 +249,19 @@
     }
 
-    return badpix;
+    return true;
 }
 
 pmMaskIdentifyMode pmMaskIdentifyModeFromString (const char *string) {
 
-    if (!strcasecmp (string, "VALUE")) {
+    if (!strcasecmp(string, "VALUE")) {
       return PM_MASK_ID_VALUE;
     }
-    if (!strcasecmp (string, "FRACTION")) {
+    if (!strcasecmp(string, "FRACTION")) {
       return PM_MASK_ID_FRACTION;
     }
-    if (!strcasecmp (string, "SIGMA")) {
+    if (!strcasecmp(string, "SIGMA")) {
       return PM_MASK_ID_SIGMA;
     }
-    if (!strcasecmp (string, "POISSON")) {
+    if (!strcasecmp(string, "POISSON")) {
       return PM_MASK_ID_POISSON;
     }