Changeset 18752 for branches/eam_branch_20080719/pswarp/src/pswarpTransformReadout_Threaded.c

Timestamp:

Jul 27, 2008, 8:00:23 AM (18 years ago)

Author:

eugene

Message:

threaded and unthreaded equivalent versions now compile

File:

• branches/eam_branch_20080719/pswarp/src/pswarpTransformReadout_Threaded.c (modified) (6 diffs)

branches/eam_branch_20080719/pswarp/src/pswarpTransformReadout_Threaded.c

@@ -15 +15 @@
     int nGridX = psMetadataLookupS32(NULL, config->arguments, "GRID.NX");
     int nGridY = psMetadataLookupS32(NULL, config->arguments, "GRID.NY");
-    psImageInterpolateMode interpolationMode = psMetadataLookupS32(NULL, config->arguments,
-                                                                   "INTERPOLATION.MODE");
+    psImageInterpolateMode interpolationMode = psMetadataLookupS32(NULL, config->arguments, "INTERPOLATION.MODE");
+
     // load the recipe
     psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSWARP_RECIPE);
-    if (!recipe) {
-        psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSWARP_RECIPE);
-        return false;
-    }
+    psAssert (recipe, "missing recipe %s", PSWARP_RECIPE);
 
     // output mask bits
@@ -30 +27 @@
     psAssert (mdok, "MASK.INPUT was not defined");
 
-    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
-    // psImagePixelInterpolate below
-    int inCol0 = input->image->col0;
-    int inRow0 = input->image->row0;
-    int outCol0 = output->image->col0;
-    int outRow0 = output->image->row0;
-
-    // we might want to do the rectangular regions outside of the selection independently
-    // psImageInit (output->image, NAN);
+    // Flux fraction for "poor"
+    float poorFrac = psMetadataLookupF32(NULL, config->arguments, "POOR.FRAC"); 
 
     // find the output pixel range
@@ -54 +42 @@
     psLogMsg ("pswarp", 3, "maximum error using this grid sampling: %f\n", maxError);
 
-    int gridX, gridY, nextGridX, nextGridY;
+    int gridX, gridY;
     pswarpMapGridSetGrid (grid, minX, minY, &gridX, &gridY);
-    nextGridY = pswarpMapGridNextGrid_Y (grid, gridY);
-    pswarpMap *map = grid->maps[gridX][gridY];
-
+# if (0)
+    // XXX these asserts probably belong in the function (don't have outCol0,Row0 anyway)
     assert ((int)(minX - outCol0) >= 0);
     assert ((int)(maxX - outCol0) <= output->image->numCols);
     assert ((int)(minY - outRow0) >= 0);
     assert ((int)(maxY - outRow0) <= output->image->numRows);
-
-    int gridXo = gridX;
-    int nextGridXo = pswarpMapGridNextGrid_X (grid, gridX);
+# endif
 
     psImage *inImage = input->image;    // Input image
@@ -76 +61 @@
                                                                        maskBad, maskPoor, poorFrac);
 
-    psPlane *inPix = psPlaneAlloc();    // Coordinates on the input detector
-    psF32 **outImageData = output->image->data.F32; // Output image pixels
-    psF32 **outVarData = NULL;          // Output variance pixels
-    psMaskType **outMaskData = NULL;    // Output mask pixels
-    if (inVar) {
-        if (!output->weight) {
-            output->weight = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_F32);
-            psImageInit(output->weight, NAN);
-        }
-        outVarData = output->weight->data.F32;
-    }
-    if (inMask || maskPoor || maskBad) {
-        if (!output->mask) {
-            output->mask = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_MASK);
-            psImageInit(output->mask, maskBad);
-        }
-        outMaskData = output->mask->data.PS_TYPE_MASK_DATA;
-    }
-
-    // create N threads
-    for (gridY = gridYmin; gridY < gridYmax; gridY++) {
-        for (gridX = gridXmin; gridX < gridXmax; gridX++) {
-            pswarpTransformTile (gridX, gridY);
-            // send function and arguments to thread manager
-            // threadFunction (pswarpTransformTile)
-            // while (threadFull()) {
-            // usleep (50000);
-            // }
+    if (inVar && !output->weight) {
+        output->weight = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_F32);
+        psImageInit(output->weight, NAN);
+    }
+    if ((inMask || maskPoor || maskBad) && !output->mask) {
+        output->mask = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_MASK);
+        psImageInit(output->mask, maskBad);
+    }
+
+    // total number of good pixels across all tiles (summed below)
+    int goodPixels = 0;
+
+    // create jobs and supply them to the threads
+    for (gridY = 0; gridY < nGridX; gridY++) {
+        for (gridX = 0; gridX < nGridY; gridX++) {
+
+            pswarpTransformTileArgs *args = pswarpTransformTileArgsAlloc();
+
+            // these items are just views to the data; they are not freed with args
+            args->input = input;
+            args->output = output;
+            args->grid = grid;
+            args->interp = interp;
+            args->region = region;
+
+            args->gridX = gridX;
+            args->gridY = gridY;
+            args->goodPixels = 0;
+
+            // allocate a job
+            psThreadJob *job = psThreadJobAlloc ("PSWARP_TRANSFORM_TILE", 0);
+
+            // construct the arguments for this job
+            // job is pswarpTransformTile (gridX, gridY);
+            psArrayAdd (job->args, 1, args);
+            psThreadJobAddPending (job);
+            psFree (args);
         }
     }
-
-
-    // Iterate over the output image pixels (parent frame)
-    long goodPixels = 0;                // Number of input pixels landing on the output image
-    for (int y = minY; y < maxY; y++) {
-        if (y >= nextGridY) {
-            gridY ++;
-            nextGridY = pswarpMapGridNextGrid_Y (grid, gridY);
-            map = grid->maps[gridX][gridY];
-        }
-
-        gridX = gridXo;
-        nextGridX = nextGridXo;
-        map = grid->maps[gridX][gridY];
-        int yOut = y - outRow0;         // Position on image
-        for (int x = minX; x < maxX; x++) {
-            if (x >= nextGridX) {
-                gridX ++;
-                nextGridX = pswarpMapGridNextGrid_X (grid, gridX);
-                map = grid->maps[gridX][gridY];
-            }
-
-            // Only transform those pixels requested
-            if (region && region->data.U8[y][x]) continue;
-
-            // pswarpMapApply converts the output coordinate (x,y) to the input coordinate.
-            // both are in the parent frames of the input and output images.
-            pswarpMapApply (&inPix->x, &inPix->y, map, x + 0.5, y + 0.5);
-
-            if (inPix->x - inCol0 < 0) continue;
-            if (inPix->x - inCol0 >= inImage->numCols) continue;
-            if (inPix->y - inRow0 < 0) continue;
-            if (inPix->y - inRow0 >= inImage->numRows) continue;
-
-            goodPixels++;
-
-            // XXX include mask
-            // XXX apply scale and offset?
-            // psImagePixelInterpolate determines the value at pixel coordinate (x,y) in child coordinates
-            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);
-                psFree(inPix);
-                psFree(grid);
-                return false;
-            }
-            int xOut = x - outCol0;     // Position on image
-            outImageData[yOut][xOut] = imageValue;
-            if (inVar) {
-                outVarData[yOut][xOut] = varValue;
-            }
-            if (outMaskData) {
-                outMaskData[yOut][xOut] = maskValue;
-            }
-        }
+    // wait here for the threaded jobs to finish
+    if (!psThreadPoolWait ()) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to interpolate image.");
+        return false;
+    }
+    fprintf (stderr, "success for threaded jobs\n");
+
+    // each job records its own goodPixel values; sum them here
+    // we have only supplied one type of job, so we can assume the types here
+    psThreadJob *job = NULL;
+    while ((job = psThreadJobGetDone()) != NULL) {
+        pswarpTransformTileArgs *args = job->args->data[0];
+        goodPixels += args->goodPixels;
+        psFree (job);
     }
     psFree(interp);
-    psFree(inPix);
     psFree(grid);
 
     // Store the variance factor and number of good pixels
     if (goodPixels > 0) {
-        float varFactor = psImageInterpolateVarianceFactor(inImage->numCols / 2.0 + inCol0,
-                                                           inImage->numRows / 2.0 + inRow0,
+        float varFactor = psImageInterpolateVarianceFactor(inImage->numCols / 2.0 + inImage->col0,
+                                                           inImage->numRows / 2.0 + inImage->row0,
                                                            interp); // Variance factor: large --> small scale
         psMetadataItem *vfItem = psMetadataLookup(output->analysis, PSWARP_ANALYSIS_VARFACTOR);
@@ -208 +157 @@
             pmModel *model = source->modelPSF; // Model for this source
             float xIn, yIn;             // Coordinates of source
-            xIn = model->params->data.F32[PM_PAR_XPOS] - inCol0;
-            yIn = model->params->data.F32[PM_PAR_YPOS] - inRow0;
+            xIn = model->params->data.F32[PM_PAR_XPOS] - inImage->col0;
+            yIn = model->params->data.F32[PM_PAR_YPOS] - inImage->row0;
             int xGrid, yGrid;           // Grid coordinates for local map
             if (!pswarpMapGridSetGrid(sourceGrid, xIn + 0.5, yIn + 0.5, &xGrid, &yGrid)) {
@@ -232 +181 @@
                 return false;
             }
-            xOut += outCol0 - 0.5;
-            yOut += outRow0 - 0.5;
+            xOut += output->image->col0 - 0.5;
+            yOut += output->image->row0 - 0.5;
             if (xOut < minX || xOut > maxX || yOut < minY || yOut > maxY) {
                 // It's not in the output image