Index: trunk/psModules/src/detrend/pmFlatNormalize.c
===================================================================
--- trunk/psModules/src/detrend/pmFlatNormalize.c	(revision 6999)
+++ trunk/psModules/src/detrend/pmFlatNormalize.c	(revision 7060)
@@ -7,25 +7,23 @@
 
 
-// Estimate the flat-field normalisation
-bool pmFlatNormalize(psVector *sourceFlux, // The source flux in each image; modified for return
-                     psVector *chipGains, // Initial guess of the chip gains; modified for return
-                     psImage *fluxLevels, // Fluxes for each integration (row) and chip (col); modified
-                     unsigned int maxIter, // Maximum number of iterations
-                     double tolerance   // Tolerance level before dying
-                    )
+// Estimate the flat-field normalisation; return the source flux in each integration
+psVector *pmFlatNormalize(bool *converge, // Did we converge?
+                          psVector *chipGains, // Initial guess of the chip gains; modified for return
+                          const psImage *fluxLevels, // Fluxes for each integration (row) and chip (col)
+                          unsigned int maxIter, // Maximum number of iterations
+                          double tolerance   // Tolerance level before dying
+                         )
 {
-    int numSources = sourceFlux->n;     // Number of measurements made
-    int numChips = chipGains->n;        // Number of chips with which each measurement is made
+    int numSources = fluxLevels->numRows; // Number of integrations
+    int numChips = fluxLevels->numCols; // Number of chips with which each integration is made
 
     // Sanity checks
-    assert(fluxLevels->numCols == numSources);
-    assert(fluxLevels->numRows == numChips);
-    assert(sourceFlux->type.type == PS_TYPE_F64);
+    assert(chipGains->n == numChips);
     assert(chipGains->type.type == PS_TYPE_F64);
-    assert(fluxLevels->type.type == PS_TYPE_F64);
     assert(maxIter >= 1);
     assert(tolerance > 0);
 
     // Take the logarithms
+    psImage *flux = psImageCopy(NULL, fluxLevels, PS_TYPE_F64); // Copy of the input flux levels matrix
     psImage *fluxMask = psImageAlloc(numSources, numChips, PS_TYPE_U8); // Mask for bad measurements
     psImageInit(fluxMask, 0);
@@ -44,17 +42,18 @@
 
         for (int j = 0; j < numSources; j++) {
-            if (isfinite(fluxLevels->data.F64[j][i]) && fluxLevels->data.F64[j][i] > 0) {
-                fluxLevels->data.F64[j][i] = log(fluxLevels->data.F64[j][i]);
+            if (isfinite(flux->data.F64[j][i]) && flux->data.F64[j][i] > 0) {
+                flux->data.F64[j][i] = log(flux->data.F64[j][i]);
             } else {
                 // Blank out this measurement
                 fluxMask->data.U8[j][i] = 1;
-                fluxLevels->data.F64[j][i] = NAN;
+                flux->data.F64[j][i] = NAN;
             }
         }
     }
-    // Don't need to initialise sourceFlux, since that is changed immediately
 
 
     double diff = INFINITY;             // Difference from previous iteration
+    psVector *sourceFlux = psVectorAlloc(numSources, PS_TYPE_F64); // The flux in each integration
+    // Don't need to initialise sourceFlux, since that is changed immediately
     for (int iter = 0; iter < maxIter && diff > tolerance; iter++) {
         diff = 0.0;
@@ -70,5 +69,5 @@
             for (int j = 0; j < numChips; j++) {
                 if (!gainMask->data.U8[j] && !fluxMask->data.U8[i][j]) {
-                    sum += fluxLevels->data.F64[i][j] - chipGains->data.F64[j];
+                    sum += flux->data.F64[i][j] - chipGains->data.F64[j];
                     number++;
                 }
@@ -93,5 +92,5 @@
             for (int j = 0; j < numSources; j++) {
                 if (!fluxMask->data.U8[j][i]) {
-                    sum += fluxLevels->data.F64[j][i] - sourceFlux->data.F64[j];
+                    sum += flux->data.F64[j][i] - sourceFlux->data.F64[j];
                     number++;
                 }
@@ -106,7 +105,8 @@
         }
     }
+    psFree(flux);
     psFree(fluxMask);
 
-    // Un-log everything
+    // Un-log the vectors
     for (int i = 0; i < numChips; i++) {
         if (!gainMask->data.U8[i]) {
@@ -122,4 +122,8 @@
     psFree(sourceMask);
 
-    return (diff <= tolerance);         // Did we converge?
+    if (converge) {
+        *converge = (diff < tolerance); // Did we converge?
+    }
+
+    return sourceFlux;
 }