Index: trunk/psModules/src/detrend/pmShutterCorrection.c
===================================================================
--- trunk/psModules/src/detrend/pmShutterCorrection.c	(revision 18859)
+++ trunk/psModules/src/detrend/pmShutterCorrection.c	(revision 18860)
@@ -957,4 +957,62 @@
 }
 
+bool pmShutterCorrectionGeneratePrepare(pmReadout *shutter, pmReadout *pattern, const psArray *inputs, psMaskType maskVal)
+{
+    PS_ASSERT_PTR_NON_NULL(shutter, false);
+    PS_ASSERT_PTR_NON_NULL(pattern, false);
+    PS_ASSERT_ARRAY_NON_NULL(inputs, false);
+
+    // determine the output image size based on the input images
+    int row0, col0, numCols, numRows;
+    if (!pmReadoutStackSetOutputSize(&col0, &row0, &numCols, &numRows, inputs)) {
+        psError(PS_ERR_UNKNOWN, false, "problem setting output readout size.");
+        return false;
+    }
+
+    // generate the required output image based on the specified sizes
+    pmReadoutStackDefineOutput(shutter, col0, row0, numCols, numRows, false, false, maskVal);
+    if (pattern) {
+	pmReadoutStackDefineOutput(pattern, col0, row0, numCols, numRows, false, false, maskVal);
+    }
+
+    psImage *nums = pmReadoutSetAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_COUNT, numCols, numRows, PS_TYPE_U16, 0); // Image with number fitted per pixel
+    if (!nums) {
+        return false;
+    }
+    psImage *sigma = pmReadoutSetAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_SIGMA, numCols, numRows, PS_TYPE_F32, NAN); // Image with stdev per pixel
+    if (!sigma) {
+        return false;
+    }
+
+    // Update the "concepts"
+    psList *inputCells = psListAlloc(NULL); // List of cells
+    for (long i = 0; i < inputs->n; i++) {
+        pmReadout *readout = inputs->data[i]; // Readout of interest
+        psListAdd(inputCells, PS_LIST_TAIL, readout->parent);
+    }
+    bool success = pmConceptsAverageCells(shutter->parent, inputCells, NULL, NULL, true);
+    psFree(inputCells);
+
+    // Correct the exposure times --- they don't make sense any more.
+    psMetadataItem *item = psMetadataLookup(shutter->parent->concepts, "CELL.EXPOSURE");
+    item->data.F32 = NAN;
+    item = psMetadataLookup(shutter->parent->concepts, "CELL.DARKTIME");
+    item->data.F32 = NAN;
+
+    shutter->data_exists = true;
+    shutter->parent->data_exists = true;
+    shutter->parent->parent->data_exists = true;
+
+    pattern->data_exists = true;
+    if (pattern->parent) { 
+	pattern->parent->data_exists = true;
+	if (pattern->parent->parent) {
+	    pattern->parent->parent->data_exists = true;
+	}
+    }
+
+    return success;
+}
+
 bool pmShutterCorrectionGenerate(pmReadout *shutter, pmReadout *pattern, const psArray *inputs,
                                  float reference, const pmShutterCorrectionData *data,
@@ -962,4 +1020,5 @@
 {
     PS_ASSERT_PTR_NON_NULL(shutter, false);
+    PS_ASSERT_PTR_NON_NULL(pattern, false);
     PS_ASSERT_ARRAY_NON_NULL(inputs, false);
     PS_ASSERT_INT_EQUAL(data->num, inputs->n, false);
@@ -975,28 +1034,15 @@
     }
 
-    pmReadoutUpdateSize(shutter, minInputCols, minInputRows, xSize, ySize, false, false, maskVal);
-    if (pattern) {
-        pmReadoutUpdateSize(pattern, minInputCols, minInputRows, xSize, ySize, false, false, maskVal);
-    }
-
-    psImage *nums = pmReadoutAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_COUNT, xSize, ySize,
-                                           PS_TYPE_U16, 0); // Image with number fitted per pixel
+    psImage *nums = pmReadoutGetAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_COUNT);
     if (!nums) {
         return false;
     }
-    psImage *sigma = pmReadoutAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_SIGMA, xSize, ySize,
-                                            PS_TYPE_F32, NAN); // Image with stdev per pixel
+    psImage *sigma = pmReadoutGetAnalysisImage(shutter, PM_READOUT_STACK_ANALYSIS_SIGMA);
     if (!sigma) {
-        psFree(nums);
         return false;
     }
 
     psImage *shutterImage = shutter->image; // Shutter correction image
-    psImage *patternImage; // Illumination pattern
-    if (pattern) {
-        patternImage = psMemIncrRefCounter(pattern->image);
-    } else {
-        patternImage = psImageAlloc(xSize, ySize, PS_TYPE_F32);
-    }
+    psImage *patternImage = pattern->image; // Illumination pattern
 
     int num = data->num;                // Number of images
@@ -1024,5 +1070,6 @@
                     errors->data.F32[r] = sqrtf(readout->weight->data.F32[yIn][xIn]) * ref;
                 } else {
-                    errors->data.F32[r] = sqrtf(image->data.F32[yIn][xIn]) * ref;
+		    // XXX guess that the input data is Poisson distributed; if we go negative, force high
+                    errors->data.F32[r] = sqrtf(fabs(image->data.F32[yIn][xIn])) * ref;
                 }
             }
@@ -1049,32 +1096,5 @@
     psFree(errors);
     psFree(counts);
-    psFree(nums);
-    psFree(sigma);
-
-    // Update the "concepts"
-    psList *inputCells = psListAlloc(NULL); // List of cells
-    for (long i = 0; i < inputs->n; i++) {
-        pmReadout *readout = inputs->data[i]; // Readout of interest
-        psListAdd(inputCells, PS_LIST_TAIL, readout->parent);
-    }
-    bool success = pmConceptsAverageCells(shutter->parent, inputCells, NULL, NULL, true);
-    psFree(inputCells);
-
-    // Correct the exposure times --- they don't make sense any more.
-    psMetadataItem *item = psMetadataLookup(shutter->parent->concepts, "CELL.EXPOSURE");
-    item->data.F32 = NAN;
-    item = psMetadataLookup(shutter->parent->concepts, "CELL.DARKTIME");
-    item->data.F32 = NAN;
-
-    shutter->data_exists = true;
-    shutter->parent->data_exists = true;
-    shutter->parent->parent->data_exists = true;
-
-    if (pattern) {
-        pattern->data_exists = true;
-        pattern->parent->data_exists = true;
-        pattern->parent->parent->data_exists = true;
-    }
-
-    return success;
-}
+
+    return true;
+}
Index: trunk/psModules/src/detrend/pmShutterCorrection.h
===================================================================
--- trunk/psModules/src/detrend/pmShutterCorrection.h	(revision 18859)
+++ trunk/psModules/src/detrend/pmShutterCorrection.h	(revision 18860)
@@ -5,6 +5,6 @@
  * @author Paul Price, IfA
  *
- * @version $Revision: 1.16 $ $Name: not supported by cvs2svn $
- * @date $Date: 2008-06-09 00:39:37 $
+ * @version $Revision: 1.17 $ $Name: not supported by cvs2svn $
+ * @date $Date: 2008-08-01 23:58:37 $
  * Copyright 2006 Institute for Astronomy, University of Hawaii
  */
@@ -189,5 +189,6 @@
     );
 
-
+// prepare outputs for shutter correction
+bool pmShutterCorrectionGeneratePrepare(pmReadout *shutter, pmReadout *pattern, const psArray *inputs, psMaskType maskVal);
 
 /// @}
