Changeset 13970 for trunk/ppImage/src/ppImageDetrendFringe.c

Timestamp:

Jun 25, 2007, 2:27:15 PM (19 years ago)

Author:

eugene

Message:

adding fringe residual measurement in addition to fringe amplitude; clean up ppImageLoop calls to stats and fringe

File:

• trunk/ppImage/src/ppImageDetrendFringe.c (modified) (13 diffs)

trunk/ppImage/src/ppImageDetrendFringe.c

@@ -10 +10 @@
 #include "ppImageDetrendFringe.h"
 
-bool ppImageDetrendFringeMeasure(pmReadout *readout, pmCell *fringe, const ppImageOptions *options)
+bool ppImageDetrendFringeMeasure(pmReadout *readout, pmCell *fringe, const bool isResidual, const ppImageOptions *options)
 {
     PS_ASSERT_PTR_NON_NULL(readout, false);
@@ -49 +49 @@
     psBinaryOp(measurements->df, measurements->df, "*", psScalarAlloc(expTime, PS_TYPE_F32));
 
+    char *scienceFringes = NULL;
+    if (isResidual) {
+        scienceFringes = psStringCopy ("FRINGE.RESIDUALS");
+    } else {
+        scienceFringes = psStringCopy ("FRINGE.MEASUREMENTS");
+    }
+
     // Science fringe measurements
-    pmFringeStats *previous = psMetadataLookupPtr(NULL, readout->parent->analysis, "FRINGE.MEASUREMENTS");
+    pmFringeStats *previous = psMetadataLookupPtr(NULL, readout->parent->analysis, scienceFringes);
     if (previous) {
         // Multiple readouts: concatenate
@@ -73 +80 @@
     }
 
-    psMetadataAdd(readout->parent->analysis, PS_LIST_TAIL, "FRINGE.MEASUREMENTS",
+    psMetadataAdd(readout->parent->analysis, PS_LIST_TAIL, scienceFringes,
                   PS_DATA_UNKNOWN | PS_META_REPLACE, "Fringe measurements", measurements);
     psFree(measurements);
-
+    psFree(scienceFringes);
     psFree(references);
 
@@ -84 +91 @@
 
 // Pull the fringes out of the cell analysis FRINGE.MEASUREMENTS for a chip
-static psArray *getFringes(const pmChip *chip // Chip of interest
-    )
+// XXX need some error checks
+static psArray *getFringes(const pmChip *chip, const char *source)
 {
     psArray *cells = chip->cells;       // Component cells
@@ -91 +98 @@
     for (int i = 0; i < cells->n; i++) {
         pmCell *cell = cells->data[i];  // Cell of interest
-        fringes->data[i] = psMemIncrRefCounter(psMetadataLookupPtr(NULL, cell->analysis,
-                                                                   "FRINGE.MEASUREMENTS"));
+        fringes->data[i] = psMemIncrRefCounter(psMetadataLookupPtr(NULL, cell->analysis, source));
     }
 
@@ -101 +107 @@
 // Solve the fringe system: we have science fringe measurements for each cell, and an array of reference
 // fringe measurements for each cell.  Need to concatenate these together first, and then solve.
-bool ppImageDetrendFringeSolve(pmChip *scienceChip, const pmChip *refChip, const ppImageOptions *options)
+bool ppImageDetrendFringeSolve(pmChip *scienceChip, const pmChip *refChip, const bool isResidual, const ppImageOptions *options)
 {
     PS_ASSERT_PTR_NON_NULL(scienceChip, NULL);
@@ -107 +113 @@
     PS_ASSERT_PTR_NON_NULL(options, NULL);
 
-    psArray *science = getFringes(scienceChip); // Fringe measurements on science chip
+    psArray *science = NULL;
+    if (isResidual) {
+        science = getFringes(scienceChip, "FRINGE.RESIDUALS"); // Fringe residuals on science chip
+    } else {
+        science = getFringes(scienceChip, "FRINGE.MEASUREMENTS"); // Fringe measurements on science chip
+    }
+
     pmFringeStats *scienceCat = pmFringeStatsConcatenate(science, NULL, NULL); // Science fringes
     psFree(science);
@@ -113 +125 @@
     // Need to transform the array of cells each with an array of fringes --> array of fringes for the chip as
     // a whole
-    psArray *references = getFringes(refChip); // Fringe measurements on reference chip
+    psArray *references = getFringes(refChip, "FRINGE.MEASUREMENTS"); // Fringe measurements on reference chip
     int numRefs = ((psArray*)references->data[0])->n; // Number of reference fringes
     psArray *referencesCat = psArrayAlloc(numRefs); // Reference fringes
@@ -132 +144 @@
                                                    options->fringeIter, options->fringeKeep);
 
-    psMetadataAdd(scienceChip->analysis, PS_LIST_TAIL, "FRINGE.SOLUTION", PS_DATA_UNKNOWN,
-                  "Fringe solution", solution);
+    if (isResidual) {
+        psMetadataAdd(scienceChip->analysis, PS_LIST_TAIL, "FRINGE.RESIDUAL.SOLUTION", PS_DATA_UNKNOWN, "Fringe solution", solution);
+    } else {
+        psMetadataAdd(scienceChip->analysis, PS_LIST_TAIL, "FRINGE.SOLUTION", PS_DATA_UNKNOWN, "Fringe solution", solution);
+    }
+
+# if (0)
+    // write the fringe amplitude or residual amplitude to the header
+    // XXX this is measured per cell, but we only have headers per chip
+    pmHDU *hdu = pmHDUFromCell(science);// HDU  of interest
+    for (int i = 0; i < solution->nFringeFrames; i++) {
+        // write metadata header value
+        psString keyword = NULL;
+        if (isResidual) {
+            psStringAppend (&keyword, "FRES_%02dV", i);
+        } else {
+            psStringAppend (&keyword, "FRNG_%02dV", i);
+        }
+        psMetadataAddF32(hdu->header, PS_LIST_TAIL, keyword, PS_META_REPLACE, "Fringe Amplitude", solution->coeff->data.F32[i + 1]);
+        psFree (keyword);
+
+        keyword = NULL;
+        if (isResidual) {
+            psStringAppend (&keyword, "FRES_%02dE", i);
+        } else {
+            psStringAppend (&keyword, "FRNG_%02dE", i);
+        }
+        psMetadataAddF32(hdu->header, PS_LIST_TAIL, keyword, PS_META_REPLACE, "Fringe Amplitude error", solution->coeffErr->data.F32[i + 1]);
+        psFree (keyword);
+    }
+# endif
+
     psFree(solution);
     psFree(scienceCat);
@@ -142 +184 @@
 
 
-bool ppImageDetrendFringeGenerate(pmCell *science, pmCell *fringes)
+bool ppImageDetrendFringeGenerate(pmCell *science, pmCell *fringes, const ppImageOptions *options)
 {
     PS_ASSERT_PTR_NON_NULL(science, false);
@@ -191 +233 @@
         psFree(md5string);
 
+# if (0)
         // write metadata header value
+        // XXX this is measured per cell, but we only have headers per chip
         psString keyword = NULL;
-        psStringAppend (&keyword, "FRNG_%02d", i);
+        psStringAppend (&keyword, "FRNG_%02dV", i);
         psMetadataAddF32(hdu->header, PS_LIST_TAIL, keyword, PS_META_REPLACE, "Fringe Amplitude", solution->coeff->data.F32[i + 1]);
         psFree (keyword);
 
         keyword = NULL;
-        psStringAppend (&keyword, "FRNG_%02dD", i);
+        psStringAppend (&keyword, "FRNG_%02dE", i);
         psMetadataAddF32(hdu->header, PS_LIST_TAIL, keyword, PS_META_REPLACE, "Fringe Amplitude error", solution->coeffErr->data.F32[i + 1]);
         psFree (keyword);
+# endif
+
     }
     if (expTime != 1.0) {
@@ -219 +265 @@
             return false;
         }
+
+        // measure residual fringe amplitude. results go to FRINGE.RESIDUALS
+        ppImageDetrendFringeMeasure (readout, fringes, true, options);
     }
     psFree(sumFringe);
@@ -233 +282 @@
 }
 
+bool ppImageDetrendFringeApply (pmConfig *config, pmChip *chip, const pmFPAview *inputView, const ppImageOptions *options) {
+
+    pmCell *cell = NULL;
+
+    assert (options->doFringe); // do not call if not needed
+    assert (inputView->chip != -1);
+    assert (inputView->cell == -1);
+    assert (inputView->readout == -1);
+
+    pmFPAview *view = pmFPAviewAlloc(0); // View for local processing
+    *view = *inputView;
+
+    // select the reference chip
+    pmChip *fringe = pmFPAfileThisChip(config->files, view, "PPIMAGE.FRINGE");
+    if (!fringe) {
+        psError(PS_ERR_UNKNOWN, false, "missing fringe reference data.\n");
+        psFree (view);
+        return false;
+    }
+
+    // Solve the fringe system
+    if (!ppImageDetrendFringeSolve(chip, fringe, false, options)) {
+        psError(PS_ERR_UNKNOWN, false, "failed to solve the fringe system.\n");
+        psFree (view);
+        return false;
+    }
+
+    // Go back over the cells to apply the fringe correction
+    view->cell = view->readout = -1;
+    while ((cell = pmFPAviewNextCell(view, chip->parent, 1)) != NULL) {
+        if (!cell->process || !cell->file_exists) {
+            continue;
+        }
+
+        // Apply the fringe correction
+        psTrace("ppImage", 3, "Applying fringe correction...\n");
+        pmCell *fringeCell = pmFPAfileThisCell(config->files, view, "PPIMAGE.FRINGE");
+        if (!fringeCell) {
+            psError(PS_ERR_UNKNOWN, false, "missing fringe reference data.\n");
+            psFree (view);
+            return false;
+        }
+
+        if (!ppImageDetrendFringeGenerate(cell, fringeCell, options)) {
+            psError(PS_ERR_UNKNOWN, false, "failed to apply fringe image.\n");
+            psFree (view);
+            return false;
+        }
+    }
+
+    // Solve the residual fringe system
+    if (!ppImageDetrendFringeSolve(chip, fringe, true, options)) {
+        psError(PS_ERR_UNKNOWN, false, "failed to solve the residual fringe system.\n");
+        psFree (view);
+        return false;
+    }
+
+    psFree (view);
+    return true;
+}