Changeset 12771 for trunk/pswarp/src/pswarpTransformReadout_Opt.c

Timestamp:

Apr 5, 2007, 5:38:16 PM (19 years ago)

Author:

Paul Price

Message:

Adding variance and mask output.

File:

• trunk/pswarp/src/pswarpTransformReadout_Opt.c (modified) (7 diffs)

trunk/pswarp/src/pswarpTransformReadout_Opt.c

@@ -2 +2 @@
 
 // NOTE: in this function, the coordinates are transformed from the OUTPUT to the INPUT
-bool pswarpTransformReadout_Opt (pmReadout *output, pmReadout *input, pmConfig *config) {
-
-    bool status;
-    int minX, minY, maxX, maxY;
-    int gridXo, gridX, gridY, nextGridXo, nextGridX, nextGridY;
-    pswarpMap *map = NULL;
-    psImageInterpolateMode interpolationMode = PS_INTERPOLATE_NONE;
-
+bool pswarpTransformReadout_Opt (pmReadout *output, pmReadout *input, pmConfig *config)
+{
     // XXX this implementation currently ignores the use of the region
     psImage *region = NULL;
@@ -15 +9 @@
     psTimerStart ("warp");
 
-    // select the current recipe
-    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSWARP_RECIPE);
-    if (!recipe) {
-        psError(PSWARP_ERR_CONFIG, true, "Can't find PSWARP recipe!\n");
-        return false;
-    }
-
-    int nGridX = psMetadataLookupS32 (&status, recipe, "GRID.NX");
-    if (!status) nGridX = 128;
-    int nGridY = psMetadataLookupS32 (&status, recipe, "GRID.NY");
-    if (!status) nGridY = 128;
-
-    char *name = psMetadataLookupStr (&status, recipe, "INTERPOLATION.MODE");
-    if (!name) {
-        interpolationMode = PS_INTERPOLATE_BILINEAR;
-        psLogMsg ("pswarp", 3, "defaulting to bilinear interpolation\n");
-    } else {
-        interpolationMode = psImageInterpolateModeFromString (name);
-        if (interpolationMode == PS_INTERPOLATE_NONE) {
-            interpolationMode = PS_INTERPOLATE_BILINEAR;
-            psLogMsg ("pswarp", 3, "unknown interpolation mode %s, defaulting to bilinear interpolation\n", name);
-        }
-    }
-
-    psPlane *inPix = psPlaneAlloc();    // Coordinates on the input detector
-    psImage *inImage = input->image;
-    psF32 **outData = output->image->data.F32;
+    // Get grid size
+    int nGridX = psMetadataLookupS32(NULL, config->arguments, "GRID.NX");
+    int nGridY = psMetadataLookupS32(NULL, config->arguments, "GRID.NY");
+    psImageInterpolateMode interpolationMode = psMetadataLookupS32(NULL, config->arguments,
+                                                                   "INTERPOLATION.MODE");
+    psMaskType maskIn = psMetadataLookupU8(NULL, config->arguments, "MASK.IN"); // Mask for input data
+    psMaskType maskPoor = psMetadataLookupU8(NULL, config->arguments, "MASK.POOR"); // Mask for "poor" data
+    psMaskType maskBad = psMetadataLookupU8(NULL, config->arguments, "MASK.IN"); // Mask for bad data
+    float poorFrac = psMetadataLookupF32(NULL, config->arguments, "POOR.FRAC"); // Flux fraction for "poor"
 
     // we need to apply the offset to convert parent coordinates to child coordinates for
@@ -54 +30 @@
 
     // find the output pixel range
+    int minX, minY, maxX, maxY;
     pswarpMatchRange (&minX, &minY, &maxX, &maxY, input, output);
 
@@ -64 +41 @@
     psLogMsg ("pswarp", 3, "maximum error using this grid sampling: %f\n", maxError);
 
+    int gridX, gridY, nextGridX, nextGridY;
     pswarpMapGridSetGrid (grid, minX, minY, &gridX, &gridY);
     nextGridY = pswarpMapGridNextGrid_Y (grid, gridY);
-    map = grid->maps[gridX][gridY];
+    pswarpMap *map = grid->maps[gridX][gridY];
 
     assert ((int)(minX - outCol0) >= 0);
@@ -73 +51 @@
     assert ((int)(maxY - outRow0) <= output->image->numRows);
 
-    gridXo = gridX;
-    nextGridXo = pswarpMapGridNextGrid_X (grid, gridX);
+    int gridXo = gridX;
+    int nextGridXo = pswarpMapGridNextGrid_X (grid, gridX);
 
-    psImageInterpolateOptions *interp = psImageInterpolateOptionsAlloc(PS_INTERPOLATE_BILINEAR, inImage,
-                                                                       NULL, NULL, 0, NAN, NAN, 0, 0, 0.0);
+    psImage *inImage = input->image;    // Input image
+    psImage *inVar   = input->weight;   // Input weight map
+    psImage *inMask  = input->mask;     // Input mask
+
+    // Interpolation options
+    psImageInterpolateOptions *interp = psImageInterpolateOptionsAlloc(interpolationMode, inImage,
+                                                                       inVar, inMask, maskIn, NAN, NAN,
+                                                                       maskPoor, maskBad, poorFrac);
+
+    psPlane *inPix = psPlaneAlloc();    // Coordinates on the input detector
+    psF32 **outImageData = output->image->data.F32; // Output image pixels
+    psF32 **outVarData;                 // Output variance pixels
+    psMaskType **outMaskData;           // Output mask pixels
+    if (inVar) {
+        if (!output->weight) {
+            output->weight = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_F32);
+        }
+        outVarData = output->weight->data.F32;
+    }
+    if (inMask) {
+        if (!output->mask) {
+            output->mask = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_MASK);
+        }
+        outMaskData = output->mask->data.PS_TYPE_MASK_DATA;
+    }
 
     // Iterate over the output image pixels (parent frame)
@@ -112 +113 @@
             // XXX apply scale and offset?
             // psImagePixelInterpolate determines the value at pixel coordinate (x,y) in child coordinates
-            double value;
-            if (!psImageInterpolate(&value, NULL, NULL, inPix->x - inCol0, inPix->y - inRow0, interp)) {
+            double imageValue, varValue;
+            psMaskType maskValue;
+            if (!psImageInterpolate(&imageValue, &varValue, &maskValue,
+                                    inPix->x - inCol0, inPix->y - inRow0, interp)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to interpolate image.");
                 psFree(interp);
@@ -120 +123 @@
                 return false;
             }
-            outData[y-outRow0][x-outCol0] = value;
+            outImageData[y-outRow0][x-outCol0] = imageValue;
+            if (inVar) {
+                outVarData[y-outRow0][x-outCol0] = varValue;
+            }
+            if (inMask) {
+                outMaskData[y-outRow0][x-outCol0] = maskValue;
+            }
         }
     }