Changeset 21457

Timestamp:

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

Author:

Paul Price

Message:

Adding modification of variance map by the shutter correction.

Location:

trunk/psModules/src/detrend

Files:

• pmDetrendThreads.c (modified) (1 diff)
• pmShutterCorrection.c (modified) (5 diffs)
• pmShutterCorrection.h (modified) (3 diffs)

trunk/psModules/src/detrend/pmDetrendThreads.c

@@ -47 +47 @@
 
     {
-        psThreadTask *task = psThreadTaskAlloc("PSMODULES_DETREND_SHUTTER", 7);
+        psThreadTask *task = psThreadTaskAlloc("PSMODULES_DETREND_SHUTTER", 8);
         task->function = &pmShutterCorrectionApplyScan_Threaded;
         psThreadTaskAdd(task);

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.");

trunk/psModules/src/detrend/pmShutterCorrection.h

@@ -5 +5 @@
  * @author Paul Price, IfA
  *
- * @version $Revision: 1.22 $ $Name: not supported by cvs2svn $
- * @date $Date: 2009-01-27 06:39:38 $
+ * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $
+ * @date $Date: 2009-02-12 19:33:30 $
  * Copyright 2006 Institute for Astronomy, University of Hawaii
  */
@@ -91 +91 @@
     float offref,                       ///< Reference time offset
     int nIter,                          ///< Number of iterations
-    float rej                           ///< Rejection threshold (sigma)
+    float rej                           ///< Rejection threshold (sigma)
     );
 
@@ -131 +131 @@
 bool pmShutterCorrectionApplyScan(
     psImage *image,                     ///< Input image to correct
+    psImage *mask,                      ///< Input mask image
+    psImage *var,                       ///< Input variance image
     const psImage *shutterImage,        ///< Shutter correction image
-    psImage *mask,                      ///< Input mask image
     float exptime,                      ///< Exposure time to which to correct
     psImageMaskType blank,                   ///< Mask value to give blank pixels