Changeset 13873 for trunk/ppMerge/src/ppMergeScaleZero.c

Timestamp:

Jun 18, 2007, 5:59:32 PM (19 years ago)

Author:

Paul Price

Message:

Changes to allow shutter correction to be processed piece by piece.

File:

• trunk/ppMerge/src/ppMergeScaleZero.c (modified) (14 diffs)

trunk/ppMerge/src/ppMergeScaleZero.c

@@ -15 +15 @@
 bool ppMergeScaleZero(psImage **scales, // The scales for each integration/cell
                       psImage **zeros, // The zeroes for each integration/cell
+                      psArray **shutters, // The shutter correction data for each cell
                       ppMergeData *data,// The data
                       const ppMergeOptions *options, // The options
@@ -24 +25 @@
     assert(config);
 
-    if (!options->scale && !options->zero && !options->darktime) {
+    if (!options->scale && !options->zero && !options->darktime && !options->shutter) {
         return true;                    // We did everything we were asked for
     }
@@ -43 +44 @@
                                  (*zeros)->numCols == data->numCells &&
                                  (*zeros)->numRows == filenames->n));
+    assert(!options->shutter || shutters);
+    assert(!shutters || !*shutters || (*shutters)->n == data->numCells);
 
     // Allocate the outputs
@@ -57 +60 @@
         *zeros = psImageAlloc(data->numCells, filenames->n, PS_TYPE_F32);
     }
-
-    bool fromConcepts = false;          // Do we get the scale and zero points from the concepts
+    psRandom *rng = NULL;               // Random number generator
+    if (options->shutter) {
+        if (*shutters) {
+            psFree(*shutters);
+        }
+        *shutters = psArrayAlloc(data->numCells);
+        rng = psRandomAlloc(PS_RANDOM_TAUS, 0);
+    }
+
+    bool fromConcepts = false;          // Do we get the scale and zero points from the concepts?
     bool first = true;                  // Are we on the first cell (that sets the standard for the rest)?
     bool done = false;                  // Are we done going through the list?
@@ -86 +97 @@
                     if (mdok && !isnan(scale)) {
                         if (!first && !fromConcepts) {
-                            psLogMsg(__func__, PS_LOG_WARN, "PPMERGE.SCALE and PPMERGE.ZERO have been set "
+                            psWarning("PPMERGE.SCALE and PPMERGE.ZERO have been set "
                                      "for some, but not all cells --- we will re-measure it for all cells.");
                             done = true;
@@ -93 +104 @@
                         fromConcepts = true;
                         (*scales)->data.F32[i][cellNum] = scale;
-                        psTrace("ppMerge", 9, "Scale for input %ld, chip %ld, cell %ld: %f\n", i, j, k, scale);
+                        psTrace("ppMerge", 9, "Scale for input %ld, chip %ld, cell %ld: %f\n",
+                                i, j, k, scale);
                     } else if (!first && fromConcepts) {
-                        psLogMsg(__func__, PS_LOG_WARN, "PPMERGE.SCALE and PPMERGE.ZERO have been set "
+                        psWarning("PPMERGE.SCALE and PPMERGE.ZERO have been set "
                                  "for some, but not all cells --- we will re-measure it for all cells.");
                         fromConcepts = false;
@@ -107 +119 @@
                     if (mdok && !isnan(zero)) {
                         if (!first && !fromConcepts) {
-                            psLogMsg(__func__, PS_LOG_WARN, "PPMERGE.SCALE and PPMERGE.ZERO have been set "
+                            psWarning("PPMERGE.SCALE and PPMERGE.ZERO have been set "
                                      "for some, but not all cells --- we will re-measure it for all cells.");
                             done = true;
@@ -116 +128 @@
                         psTrace("ppMerge", 9, "Zero for input %ld, chip %ld, cell %ld: %f\n", i, j, k, zero);
                     } else if (!first && fromConcepts) {
-                        psLogMsg(__func__, PS_LOG_WARN, "PPMERGE.SCALE and PPMERGE.ZERO have been set "
+                        psWarning("PPMERGE.SCALE and PPMERGE.ZERO have been set "
                                  "for some, but not all cells --- we will re-measure it for all cells.");
                         fromConcepts = false;
@@ -131 +143 @@
     if (fromConcepts) {
         // We've already done everything we need to
+        psFree(rng);
         return true;
     }
@@ -205 +218 @@
 
                     if (cell->readouts->n > 1) {
-                        psLogMsg(__func__, PS_LOG_WARN, "File %s chip %d cell %d contains more than one "
+                        psWarning("File %s chip %d cell %d contains more than one "
                                  "readout --- ignoring all but the first.\n", name, j, k);
                         status = false;
@@ -222 +235 @@
                         gains->data.F32[cellNum] = psMetadataLookupF32(&mdok, cell->concepts, "CELL.GAIN");
                         if (!mdok || isnan(gains->data.F32[cellNum])) {
-                            psLogMsg(__func__, PS_LOG_WARN, "CELL.GAIN for file %s chip %d cell %d is not "
+                            psWarning("CELL.GAIN for file %s chip %d cell %d is not "
                                      "set.\n", name, j, k);
                             gains->data.F32[cellNum] = NAN;
@@ -259 +272 @@
                 }
 
+                // Shutter correction
+                if (options->shutter) {
+                    if (!pmCellRead(cell, inFile, config->database)) {
+                        // Nothing here
+                        pmCellFreeData(cell);
+                        continue;
+                    }
+
+                    if (cell->readouts->n > 1) {
+                        psWarning("File %s chip %d cell %d contains more than one "
+                                  "readout --- ignoring all but the first.\n", name, j, k);
+                        status = false;
+                    }
+                    pmReadout *readout = cell->readouts->data[0]; // The readout of interest
+
+                    pmShutterCorrectionData *shutter = (*shutters)->data[cellNum]; // Shutter correction data
+                    if (!shutter) {
+                        shutter = pmShutterCorrectionDataAlloc(readout->image->numCols,
+                                                               readout->image->numRows,
+                                                               options->shutterSize);
+                        (*shutters)->data[cellNum] = shutter;
+                    }
+                    if (!pmShutterCorrectionAddReadout(shutter, readout, options->mean, options->stdev,
+                                                       options->combine->maskVal, rng)) {
+                        psWarning("Can't add file %s chip %d cell %d to shutter correction --- ignored.",
+                                  name, j, k);
+                        status = false;
+                    }
+                }
+
+
                 pmCellFreeData(cell);
             }
@@ -265 +309 @@
         pmFPAFreeData(fpa);
     }
+    psFree(rng);
     psFree(bgStats);
     psFree(view);
@@ -277 +322 @@
         psVector *fluxes = NULL;        // Solution to fluxes
         if (!pmFlatNormalize(&fluxes, &gains, background)) {
-            psLogMsg(__func__, PS_LOG_WARN, "Normalisation failed to converge --- continuing anyway.\n");
+            psWarning("Normalisation failed to converge --- continuing anyway.\n");
             status = false;
         }