Changeset 21457 for trunk/psModules/src/detrend/pmShutterCorrection.c

Timestamp:

Feb 12, 2009, 9:33:30 AM (17 years ago)

Author:

Paul Price

Message:

Adding modification of variance map by the shutter correction.

File:

• trunk/psModules/src/detrend/pmShutterCorrection.c (modified) (5 diffs)

trunk/psModules/src/detrend/pmShutterCorrection.c

@@ -659 +659 @@
     PS_ASSERT_THREAD_JOB_NON_NULL(job, false);
 
-    psImage *image        = job->args->data[0];
-    const psImage *shutterImage = job->args->data[1];
-    psImage *mask         = job->args->data[2];
-
-    float exptime    = PS_SCALAR_VALUE(job->args->data[3],F32);
-    psImageMaskType blank = PS_SCALAR_VALUE(job->args->data[4],PS_TYPE_IMAGE_MASK_DATA);
-    int rowStart     = PS_SCALAR_VALUE(job->args->data[5],S32);
-    int rowStop      = PS_SCALAR_VALUE(job->args->data[6],S32);
-    return pmShutterCorrectionApplyScan (image, shutterImage, mask, exptime, blank, rowStart, rowStop);
-}
-
-bool pmShutterCorrectionApplyScan(psImage *image, const psImage *shutterImage, psImage *mask, float exptime,
+    psImage *image = job->args->data[0];
+    psImage *mask  = job->args->data[1];
+    psImage *var   = job->args->data[2];
+    const psImage *shutterImage = job->args->data[3];
+
+    float exptime    = PS_SCALAR_VALUE(job->args->data[4], F32);
+    psImageMaskType blank = PS_SCALAR_VALUE(job->args->data[5], PS_TYPE_IMAGE_MASK_DATA);
+    int rowStart     = PS_SCALAR_VALUE(job->args->data[6], S32);
+    int rowStop      = PS_SCALAR_VALUE(job->args->data[7], S32);
+    return pmShutterCorrectionApplyScan(image, mask, var, shutterImage, exptime, blank, rowStart, rowStop);
+}
+
+bool pmShutterCorrectionApplyScan(psImage *image, psImage *mask, psImage *var,
+                                  const psImage *shutterImage, float exptime,
                                   psImageMaskType blank, int rowStart, int rowStop)
 {
+    // Neglecting asserts because inputs should have been checked already
+
+    int numCols = image->numCols;       // Number of columns
+
     for (int y = rowStart; y < rowStop; y++) {
-        for (int x = 0; x < image->numCols; x++) {
+        for (int x = 0; x < numCols; x++) {
             if (mask && !isfinite(shutterImage->data.F32[y][x])) {
                 mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] |= blank;
                 image->data.F32[y][x] = NAN;
+                if (var) {
+                    var->data.F32[y][x] = NAN;
+                }
                 continue;
             }
-            image->data.F32[y][x] *= exptime / (exptime + shutterImage->data.F32[y][x]);
+            float correction = exptime / (exptime + shutterImage->data.F32[y][x]); // Correction factor
+            image->data.F32[y][x] *= correction;
+            if (var) {
+                var->data.F32[y][x] *= PS_SQR(correction);
+            }
         }
     }
@@ -688 +701 @@
 bool pmShutterCorrectionApply(pmReadout *readout, const pmReadout *shutter, psImageMaskType blank)
 {
-    PS_ASSERT_PTR_NON_NULL(readout, false);
-    PS_ASSERT_PTR_NON_NULL(shutter, false);
-    PS_ASSERT_IMAGE_NON_NULL(readout->image, false);
-    PS_ASSERT_IMAGE_NON_NULL(shutter->image, false);
-    PS_ASSERT_IMAGE_TYPE(readout->image, PS_TYPE_F32, false);
-    PS_ASSERT_IMAGE_TYPE(shutter->image, PS_TYPE_F32, false);
+    PM_ASSERT_READOUT_NON_NULL(readout, false);
+    PM_ASSERT_READOUT_NON_NULL(shutter, false);
+    PM_ASSERT_READOUT_IMAGE(readout, false);
+    PM_ASSERT_READOUT_IMAGE(shutter, false);
 
     psRegion region = psRegionSet(readout->col0, readout->col0 + readout->image->numCols,
@@ -731 +742 @@
     psImage *image = readout->image;    // Image to correct
     psImage *mask = readout->mask;      // Corresponding mask
+    psImage *var = readout->variance;   // Corresponding variance map
 
     bool threaded = true;
@@ -766 +778 @@
                 psThreadJob *job = psThreadJobAlloc("PSMODULES_DETREND_SHUTTER");
                 psArrayAdd(job->args, 1, image);
+                psArrayAdd(job->args, 1, mask);
+                psArrayAdd(job->args, 1, var);
                 psArrayAdd(job->args, 1, shutterImage);
-                psArrayAdd(job->args, 1, mask);
                 PS_ARRAY_ADD_SCALAR(job->args, exptime, PS_TYPE_F32);
                 PS_ARRAY_ADD_SCALAR(job->args, blank, PS_TYPE_IMAGE_MASK);
@@ -778 +791 @@
                 }
                 psFree(job);
-            } else if (!pmShutterCorrectionApplyScan(image, shutterImage, mask, exptime, blank,
+            } else if (!pmShutterCorrectionApplyScan(image, mask, var, shutterImage, exptime, blank,
                                                      rowStart, rowStop)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to apply shutter correction.");