Changeset 1802 for trunk/psModules/src/pmFlatField.c

Timestamp:

Sep 13, 2004, 4:29:40 PM (22 years ago)

Author:

harman

Message:

Added initial module code

File:

• trunk/psModules/src/pmFlatField.c (modified) (2 diffs)

trunk/psModules/src/pmFlatField.c

@@ -1 @@
-/** @file  pmFlatField.c
+/** @file  pmFlatField.cS
  *
- *  @brief Contains the definitions for ...
+ *  @brief Given an input image and a flat-field image, pmFlatField shall divide the input image by the flat-field image.
  *
- *  @author Robert DeSonia, MHPCC
+ *  The input image, in, and the flat-field image, flat, need not be the same size, since the input image may already have
+ *  been trimmed (following overscan subtraction), but the function shall use the offsets in the image (in->x0
+ *  and in->y0) to determine the appropriate offsets to obtain the correct pixel on the flat-field. In the event that the flat
+ *  image is too small (i.e., pixels on the input image refer to pixels outside the range of the flat image), the function shall
+ *  generate an error. Pixels which are negative or zero in the flat shall be masked in the input image with the value PM_MASK_FLAT
+ *  Negative pixels in the flat may be set to zero so that they are treated identically to zeroes. Any pixels masked in the flat
+ *  shall be masked with corresponding values in the output. The function shall not normalize the flat; this responsibility is
+ *  left to the caller. This function is basically equivalent to a divide (with psImageOp), but with care for the region that is
+ *  divided, checking for negative pixels, and copying of the mask from the flat to the output.
  *
- *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-09-10 23:42:18 $
+ *  @author Ross Harman, MHPCC
+ *
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-09-14 02:29:27 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -16 +26 @@
 
 #include<stdio.h>
+#include<math.h>
 
 #include "pmFlatField.h"
 
-void pmFlatField(psVector* vector)
+
+#define CHECK_NULL(STRUCT)                                                                                   \
+if (STRUCT == NULL) {                                                                                        \
+    psAbort(__func__, " : Line %d - Null not allowed for %s", __LINE__, #STRUCT);                            \
+}                                                                                                            \
+
+psReadout *pmFlatField(psReadout *in, psReadout *inMask, const psReadout *flat, const psReadout *flatMask)
 {
-    if (vector == NULL) {
-        printf("Nothing to see here... move along.\n");
-        return;
+    int i = 0;
+    int j = 0;
+    psElemType inType;
+    psElemType flatType;
+    psImage *inImage = NULL;
+    psImage *flatImage = NULL;
+    psImage *inMaskImage = NULL;
+    psImage *flatMaskImage = NULL;
+
+    CHECK_NULL(in);
+    CHECK_NULL(flat);
+    CHECK_NULL(inMask);
+    CHECK_NULL(flatMask);
+
+    inImage = in->image;
+    flatImage = flat->image;
+    inMaskImage = inMask->image;
+    flatMaskImage = flatMask->image;
+
+    CHECK_NULL(inImage);
+    CHECK_NULL(flatImage);
+    CHECK_NULL(inMaskImage);
+    CHECK_NULL(flatMaskImage);
+
+    if(inImage->numRows > flatImage->numRows) {
+        psError(__func__, " : Line %d - Input image exceeds flat image in number of rows. (%d vs %d).",
+                __LINE__, inImage->numRows, flatImage->numRows);
+        return NULL;
     }
 
-    for (int lcv=0;lcv<vector->n;lcv++) {
-        printf("%f\n",vector->data.F32[lcv]);
+    if(inImage->numCols > flatImage->numCols) {
+        psError(__func__, " : Line %d - Input image exceeds flat image in number of columns. (%d vs %d).",
+                __LINE__, inImage->numCols, flatImage->numCols);
+        return NULL;
     }
 
-    return;
+    if(inImage->row0 > flatImage->numRows-1) {
+        psError(__func__, " : Line %d - Input image row offset exceeds number of flat image rows. (%d vs %d).",
+                __LINE__, inImage->row0, flatImage->numRows);
+        return NULL;
+    }
 
+    if(inImage->col0 > flatImage->numCols-1) {
+        psError(__func__, " : Line %d - Input image column offset exceeds number of flat image columns. (%d vs %d).",
+                __LINE__, inImage->col0, flatImage->numCols);
+        return NULL;
+    }
+
+    inType = in->image->type.type;
+    flatType = flat->image->type.type;
+    if(inType != flatType) {
+        psError(__func__, " : Line %d - Input and flat image types don't agree. (%d vs %d).",
+                __LINE__, inType, flatType);
+        return NULL;
+    }
+
+    // Macro for all PS types
+    #define PM_FLAT_DIVISION(TYPE)                                                                           \
+case PS_TYPE_##TYPE:                                                                                 \
+    for(j = inImage->row0; j < inImage->numRows; j++) {                                              \
+        for(i = inImage->col0; i < inImage->numCols; i++) {                                          \
+            if(creal(flatImage->data.TYPE[j][i]) <= 0.0) {                                           \
+                /* Negative or zero flat pixels shall be masked in input image as  PM_MASK_FLAT */   \
+                inMaskImage->data.F64[j][i] = PM_MASK_FLAT;                                          \
+            } else if(flatMaskImage->data.TYPE[j][i]){                                               \
+                /* Pixels masked in the flat shall be masked with same values in the output */       \
+                inMaskImage->data.TYPE[j][i] = flatMaskImage->data.TYPE[j][i];                       \
+            }                                                                                        \
+            if(creal(flatImage->data.TYPE[j][i]) < 0.0) {                                            \
+                /* Negative pixels in the flat may be set to zero */                                 \
+                inImage->data.TYPE[j][i] = 0;                                                        \
+            } else if(fabs(flatImage->data.TYPE[j][i]) < FLT_EPSILON) {                              \
+                psError(__func__, " : Line %d - Divide by zero. Row: %d Col: %d",__LINE__, j, i);    \
+            } else {                                                                                 \
+                /* Module shall divide the input image by the flat-fielded image */                  \
+                inImage->data.TYPE[j][i] = inImage->data.TYPE[j][i]/flatImage->data.TYPE[j][i];      \
+            }                                                                                        \
+        }                                                                                            \
+    }                                                                                                \
+    break;
+
+    switch(inType) {
+        PM_FLAT_DIVISION(U8);
+        PM_FLAT_DIVISION(U16);
+        PM_FLAT_DIVISION(U32);
+        PM_FLAT_DIVISION(U64);
+        PM_FLAT_DIVISION(S8);
+        PM_FLAT_DIVISION(S16);
+        PM_FLAT_DIVISION(S32);
+        PM_FLAT_DIVISION(S64);
+        PM_FLAT_DIVISION(F32);
+        PM_FLAT_DIVISION(F64);
+        PM_FLAT_DIVISION(C32);
+        PM_FLAT_DIVISION(C64);
+    default:
+        psError(__func__, "Unsupported image datatype (%d)", inType);
+    }
+
+    return in;
 }