Changeset 28405

Timestamp:

Jun 18, 2010, 2:25:38 PM (16 years ago)

Author:

Paul Price

Message:

Discovered a race condition in the threading process. The thread process was: psThreadJobAlloc, populate job with arguments, psThreadJobAddPending, psFree the job. This is not thread-safe, because we're decrementing the job's reference counter while the reference count is also being fiddled within the thread system (psThread.c; by pulling the job off the 'pending' queue and putting it on the 'done' queue). The reference count fiddling within the thread system was protected by a mutex, but the external free was not. Therefore, I'm pulling the free into the thread system so it can be protected (without the user having to worry about locks to do the simple case). This means that the user should no longer touch the job once he hands it over to the thread system (unless it's protected by calls to psThreadLock/psThreadUnlock, or it's on the other side of psThreadPoolWait from all the adds so that the threads aren't doing anything).

Location:

trunk

Files:

• ppStack/src/ppStackCombineFinal.c (modified) (1 diff)
• ppStack/src/ppStackCombineInitial.c (modified) (1 diff)
• ppStack/src/ppStackReject.c (modified) (1 diff)
• psLib/src/imageops/psImageConvolve.c (modified) (7 diffs)
• psLib/src/imageops/psImageCovariance.c (modified) (2 diffs)
• psLib/src/sys/psThread.c (modified) (4 diffs)
• psLib/src/sys/psThread.h (modified) (3 diffs)
• psModules/src/camera/pmReadoutFake.c (modified) (1 diff)
• psModules/src/detrend/pmBias.c (modified) (1 diff)
• psModules/src/detrend/pmDark.c (modified) (1 diff)
• psModules/src/detrend/pmFlatField.c (modified) (1 diff)
• psModules/src/detrend/pmShutterCorrection.c (modified) (3 diffs)
• psModules/src/imcombine/pmStackReject.c (modified) (1 diff)
• psModules/src/imcombine/pmSubtraction.c (modified) (1 diff)
• psModules/src/imcombine/pmSubtractionEquation.c (modified) (1 diff)
• psModules/src/imcombine/pmSubtractionMatch.c (modified) (1 diff)
• psphot/src/psphotApResid.c (modified) (14 diffs)
• psphot/src/psphotBlendFit.c (modified) (7 diffs)
• psphot/src/psphotGuessModels.c (modified) (2 diffs)
• psphot/src/psphotMagnitudes.c (modified) (5 diffs)
• psphot/src/psphotSourceStats.c (modified) (8 diffs)
• psphot/src/psphotTest.c (modified) (5 diffs)
• pswarp/src/pswarpTransformReadout.c (modified) (1 diff)

trunk/ppStack/src/ppStackCombineFinal.c

@@ -79 +79 @@
         PS_ARRAY_ADD_SCALAR(job->args, normalise, PS_TYPE_U8);
         if (!psThreadJobAddPending(job)) {
-            psFree(job);
             psFree(reject);
             return false;
         }
-        psFree(job);
     }
 

trunk/ppStack/src/ppStackCombineInitial.c

@@ -50 +50 @@
         psArrayAdd(job->args, 1, config);
         if (!psThreadJobAddPending(job)) {
-            psFree(job);
             return false;
         }
-        psFree(job);
     }
 

trunk/ppStack/src/ppStackReject.c

@@ -54 +54 @@
         PS_ARRAY_ADD_SCALAR(job->args, i, PS_TYPE_S32);
         if (!psThreadJobAddPending(job)) {
-            psFree(job);
             return false;
         }
-        psFree(job);
     }
     if (!psThreadPoolWait(true)) {

trunk/psLib/src/imageops/psImageConvolve.c

@@ -864 +864 @@
             PS_ARRAY_ADD_SCALAR(job->args, stop, PS_TYPE_S32);
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 psFree(gaussNorm);
                 psFree(out);
                 return NULL;
             }
-            psFree(job);
         }
         if (!psThreadPoolWait(true)) {
@@ -1213 +1211 @@
                   if (!psThreadJobAddPending(job)) {
                       psError(PS_ERR_UNKNOWN, false, "Unable to smooth image");
-                      psFree(job);
                       psFree(calculation);
                       psFree(calcMask);
@@ -1219 +1216 @@
                       return false;
                   }
-                  psFree(job);
               }
               // wait here for the threaded jobs to finish (NOP if threading is not active)
@@ -1261 +1257 @@
                   if (!psThreadJobAddPending(job)) {
                       psError(PS_ERR_UNKNOWN, false, "Unable to smooth image");
-                      psFree(job);
                       psFree(calculation);
                       psFree(calcMask);
@@ -1267 +1262 @@
                       return false;
                   }
-                  psFree(job);
               }
 
@@ -1580 +1574 @@
             PS_ARRAY_ADD_SCALAR(job->args, 0x00, PS_TYPE_U8); // specify rows
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 psFree(conv);
                 psFree(out);
                 return NULL;
             }
-            psFree(job);
         }
         if (!psThreadPoolWait(true)) {
@@ -1617 +1609 @@
             PS_ARRAY_ADD_SCALAR(job->args, 0xff, PS_TYPE_U8); // specify cols
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 psFree(conv);
                 psFree(out);
                 return NULL;
             }
-            psFree(job);
         }
         if (!psThreadPoolWait(true)) {

trunk/psLib/src/imageops/psImageCovariance.c

@@ -169 +169 @@
                     return NULL;
                 }
-                psFree(job);
             } else {
                 out->kernel[y][x] = imageCovarianceCalculate(covar, kernel, x, y);
@@ -334 +333 @@
                     return NULL;
                 }
-                psFree(job);
             } else {
                 out->kernel[y][x] = imageCovarianceBin(covar, bin, binVal, x, y);

trunk/psLib/src/sys/psThread.c

@@ -99 +99 @@
 bool psThreadJobAddPending(psThreadJob *job)
 {
-    PS_ASSERT_THREAD_JOB_NON_NULL(job, false);
+    if (!job) {
+        // Asking for a no-op
+        return true;
+    }
 
     psThreadTask *task = psHashLookup(tasks, job->type); // Task to execute job
@@ -114 +117 @@
         }
         psListAdd(done, PS_LIST_TAIL, job);
-
+        psFree(job);
         return task->function(job);
     }
@@ -123 +126 @@
     }
     psListAdd(pending, PS_LIST_TAIL, job);
+    psFree(job);
     psThreadUnlock();
 
@@ -279 +283 @@
     for (int i = 0; i < nThreads; i++) {
         psThread *thread = pool->data[i] = psThreadAlloc(); // Thread for pool
-        if (!pthread_create(&threads[i], NULL, psThreadLauncher, thread)) {
+        if (pthread_create(&threads[i], NULL, psThreadLauncher, thread)) {
             psAbort("Unable to create thread");
         }

trunk/psLib/src/sys/psThread.h

@@ -73 +73 @@
 
 /// Add a pending job to the queue
+///
+/// This function swallows the provided job, so that the user no longer owns it.  This is because freeing the
+/// job is not thread-safe (its reference count is being changed within the threads) so we handle it ourselves
+/// and absolve the user from all responsibility.  If the user stores the job, he should only access it while
+/// threads are processing in code protected by psThreadLock/psThreadUnlock.
 bool psThreadJobAddPending(psThreadJob *job);
 
 /// Get a job off the queue of pending jobs
+///
+/// This function is not thread-safe.  Protect with psThreadLock/psThreadUnlock if threads are running.
 psThreadJob *psThreadJobGetPending(void);
 
 /// Get a job off the queue of done jobs
+///
+/// This function is not thread-safe.  Protect with psThreadLock/psThreadUnlock if threads are running.
 psThreadJob *psThreadJobGetDone(void);
 
@@ -115 +124 @@
 bool psThreadPoolFinalize(void);
 
-
+#if 0
 /// Add thread-specific data
 ///
@@ -134 +143 @@
 bool psThreadDataRemove(const char *name // Name of data
     );
-
+#endif
 
 /// @}

trunk/psModules/src/camera/pmReadoutFake.c

@@ -303 +303 @@
 
                 if (!psThreadJobAddPending(job)) {
-                    psFree(job);
                     psFree(groups);
                     return false;
                 }
-                psFree(job);
             }
             if (!psThreadPoolWait(true)) {

trunk/psModules/src/detrend/pmBias.c

@@ -144 +144 @@
 
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 return false;
             }
-            psFree(job);
         } else if (!pmBiasSubtractScan(in, sub, scale, xOffset, yOffset, rowStart, rowStop)) {
             psError(PS_ERR_UNKNOWN, false, "Unable to apply bias correction.");

trunk/psModules/src/detrend/pmDark.c

@@ -587 +587 @@
 
             if (!psThreadJobAddPending (job)) {
-                psFree(job);
                 psFree(orders);
                 psFree(values);
                 return false;
             }
-            psFree(job);
         } else if (!pmDarkApplyScan(readout, dark, orders, values, bad, doNorm, norm, rowStart, rowStop)) {
             psError(PS_ERR_UNKNOWN, false, "Unable to apply dark.");

trunk/psModules/src/detrend/pmFlatField.c

@@ -150 +150 @@
 
           if (!psThreadJobAddPending(job)) {
-              psFree(job);
               return false;
           }
-          psFree(job);
       } else if (!pmFlatFieldScan(inImage, inMask, inVar, flatImage, flatMask, badFlat,
                                   xOffset, yOffset, rowStart, rowStop)) {

trunk/psModules/src/detrend/pmShutterCorrection.c

@@ -351 +351 @@
     psStats *resStats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN | PS_STAT_SAMPLE_STDEV);
     if (!psVectorStats (rawStats, counts, NULL, NULL, 0)) {
-        psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
-        return NULL;
+        psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
+        return NULL;
     }
     if (!psVectorStats (resStats, resid, NULL, NULL, 0)) {
-        psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
-        return NULL;
+        psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
+        return NULL;
     }
 
@@ -794 +794 @@
 
                 if (!psThreadJobAddPending(job)) {
-                    psFree(job);
                     return false;
                 }
-                psFree(job);
             } else if (!pmShutterCorrectionApplyScan(image, mask, var, shutterImage, exptime, blank,
                                                      rowStart, rowStop)) {
@@ -1017 +1015 @@
 
         if (corr) {
-            psTrace("psModules.detrend", 5, "Shutter correction fit: scale: %f, offset: %f, offref: %f\n", corr->scale, corr->offset, corr->offref);
-            if (isfinite(corr->offref) && corr->valid) {
-                psTrace("psModules.detrend", 5, "Sample reference value: %f\n", corr->offref);
-                meanRef += corr->offref;
-                numGood++;
-            }
+            psTrace("psModules.detrend", 5, "Shutter correction fit: scale: %f, offset: %f, offref: %f\n", corr->scale, corr->offset, corr->offref);
+            if (isfinite(corr->offref) && corr->valid) {
+                psTrace("psModules.detrend", 5, "Sample reference value: %f\n", corr->offref);
+                meanRef += corr->offref;
+                numGood++;
+            }
         } else {
-            psTrace("psModules.detrend", 5, "failed Shutter correction fit\n");
-        }
+            psTrace("psModules.detrend", 5, "failed Shutter correction fit\n");
+        }
 
         psFree(corr);

trunk/psModules/src/imcombine/pmStackReject.c

@@ -289 +289 @@
                     PS_ARRAY_ADD_SCALAR(args, poorFrac, PS_TYPE_F32);
                     if (!psThreadJobAddPending(job)) {
-                        psFree(job);
                         psFree(source);
                         psFree(target);
                         return NULL;
                     }
-                    psFree(job);
                 } else if (!stackRejectGrow(target, source, kernels, numCols, numRows,
                                             i, xSubMax, j, ySubMax, poorFrac)) {

trunk/psModules/src/imcombine/pmSubtraction.c

@@ -1291 +1291 @@
 
                 if (!psThreadJobAddPending(job)) {
-                    psFree(job);
                     return false;
                 }
-                psFree(job);
             } else {
                 subtractionConvolvePatch(numCols, numRows, x0, y0, out1, out2, convMask, ro1, ro2,

trunk/psModules/src/imcombine/pmSubtractionEquation.c

@@ -860 +860 @@
             PS_ARRAY_ADD_SCALAR(job->args, mode, PS_TYPE_S32);
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 return false;
             }
-            psFree(job);
         } else {
             pmSubtractionCalculateEquationStamp(stamps, kernels, i, mode);

trunk/psModules/src/imcombine/pmSubtractionMatch.c

@@ -1060 +1060 @@
             PS_ARRAY_ADD_SCALAR(job->args, bg2, PS_TYPE_F32);
             if (!psThreadJobAddPending(job)) {
-                psFree(job);
                 return false;
             }
-            psFree(job);
         } else {
             if (!pmSubtractionOrderStamp(ratios, mask, stamps, models, modelSums, i, bg1, bg2)) {

trunk/psphot/src/psphotApResid.c

@@ -22 +22 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        if (i == chisqNum) continue; // skip chisq image
-        if (!psphotApResidReadout (config, view, filerule, i, recipe)) {
+        if (i == chisqNum) continue; // skip chisq image
+        if (!psphotApResidReadout (config, view, filerule, i, recipe)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed to measure aperture residual for %s entry %d", filerule, i);
-            return false;
-        }
+            return false;
+        }
     }
     return true;
@@ -56 +56 @@
 
     if (!sources->n) {
-        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping ap resid");
-        return true;
+        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping ap resid");
+        return true;
     }
 
@@ -66 +66 @@
     int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
     if (!status) {
-        nThreads = 0;
+        nThreads = 0;
     }
 
@@ -128 +128 @@
     for (int i = 0; i < cellGroups->n; i++) {
 
-        psArray *cells = cellGroups->data[i];
-
-        for (int j = 0; j < cells->n; j++) {
-
-            // allocate a job -- if threads are not defined, this just runs the job
-            psThreadJob *job = psThreadJobAlloc ("PSPHOT_APRESID_MAGS");
-
-            psArrayAdd(job->args, 1, cells->data[j]); // sources
-            psArrayAdd(job->args, 1, psf);
-            PS_ARRAY_ADD_SCALAR(job->args, photMode, PS_TYPE_S32);
-            PS_ARRAY_ADD_SCALAR(job->args, maskVal,  PS_TYPE_IMAGE_MASK);
-            PS_ARRAY_ADD_SCALAR(job->args, markVal,  PS_TYPE_IMAGE_MASK);
-
-            PS_ARRAY_ADD_SCALAR(job->args, 0,        PS_TYPE_S32); // this is used as a return value for Nskip
-            PS_ARRAY_ADD_SCALAR(job->args, 0,        PS_TYPE_S32); // this is used as a return value for Nfail
-
-            if (!psThreadJobAddPending(job)) {
-                psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
-                return false;
-            }
-            psFree(job);
-        }
-
-        // wait for the threads to finish and manage results
-        if (!psThreadPoolWait (false)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-            return false;
-        }
-
-        // we have only supplied one type of job, so we can assume the types here
-        psThreadJob *job = NULL;
-        while ((job = psThreadJobGetDone()) != NULL) {
-            if (job->args->n < 1) {
-                fprintf (stderr, "error with job\n");
-            } else {
-                psScalar *scalar = NULL;
-                scalar = job->args->data[5];
-                Nskip += scalar->data.S32;
-                scalar = job->args->data[6];
-                Nfail += scalar->data.S32;
-            }
-            psFree(job);
-        }
+        psArray *cells = cellGroups->data[i];
+
+        for (int j = 0; j < cells->n; j++) {
+
+            // allocate a job -- if threads are not defined, this just runs the job
+            psThreadJob *job = psThreadJobAlloc ("PSPHOT_APRESID_MAGS");
+
+            psArrayAdd(job->args, 1, cells->data[j]); // sources
+            psArrayAdd(job->args, 1, psf);
+            PS_ARRAY_ADD_SCALAR(job->args, photMode, PS_TYPE_S32);
+            PS_ARRAY_ADD_SCALAR(job->args, maskVal,  PS_TYPE_IMAGE_MASK);
+            PS_ARRAY_ADD_SCALAR(job->args, markVal,  PS_TYPE_IMAGE_MASK);
+
+            PS_ARRAY_ADD_SCALAR(job->args, 0,        PS_TYPE_S32); // this is used as a return value for Nskip
+            PS_ARRAY_ADD_SCALAR(job->args, 0,        PS_TYPE_S32); // this is used as a return value for Nfail
+
+            if (!psThreadJobAddPending(job)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+                return false;
+            }
+        }
+
+        // wait for the threads to finish and manage results
+        if (!psThreadPoolWait (false)) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+            return false;
+        }
+
+        // we have only supplied one type of job, so we can assume the types here
+        psThreadJob *job = NULL;
+        while ((job = psThreadJobGetDone()) != NULL) {
+            if (job->args->n < 1) {
+                fprintf (stderr, "error with job\n");
+            } else {
+                psScalar *scalar = NULL;
+                scalar = job->args->data[5];
+                Nskip += scalar->data.S32;
+                scalar = job->args->data[6];
+                Nfail += scalar->data.S32;
+            }
+            psFree(job);
+        }
     }
 
@@ -184 +182 @@
     Npsf = 0;
 
-# ifdef DEBUG    
+# ifdef DEBUG
     FILE *f = fopen ("apresid.dat", "w");
     psAssert (f, "failed open");
@@ -199 +197 @@
         if (source->mode &  PM_SOURCE_MODE_POOR) SKIPSTAR ("POOR STAR");
 
-        if (source->mode &  PM_SOURCE_MODE_EXT_LIMIT) SKIPSTAR ("EXTENDED");
-        if (source->mode &  PM_SOURCE_MODE_CR_LIMIT) SKIPSTAR ("COSMIC RAY");
-        if (source->mode &  PM_SOURCE_MODE_DEFECT) SKIPSTAR ("DEFECT");
-            
+        if (source->mode &  PM_SOURCE_MODE_EXT_LIMIT) SKIPSTAR ("EXTENDED");
+        if (source->mode &  PM_SOURCE_MODE_CR_LIMIT) SKIPSTAR ("COSMIC RAY");
+        if (source->mode &  PM_SOURCE_MODE_DEFECT) SKIPSTAR ("DEFECT");
+
         if (!isfinite(source->apMag) || !isfinite(source->psfMag)) {
             continue;
         }
 
-        // XXX make this user-configurable?
-        if (source->errMag > 0.01) continue;
+        // XXX make this user-configurable?
+        if (source->errMag > 0.01) continue;
 
         // aperture residual for this source
@@ -221 +219 @@
 
 # ifdef DEBUG
-        fprintf (f, "%6.1f %6.1f : %6.1f %6.1f : %8.3f %8.3f %8.3f : %f : %f %f %f : %f\n",
-                 source->peak->xf, source->peak->yf, 
-                 source->modelPSF->params->data.F32[PM_PAR_XPOS], source->modelPSF->params->data.F32[PM_PAR_YPOS], 
-                 source->psfMag, source->apMag, source->errMag,
-                 source->modelPSF->params->data.F32[PM_PAR_I0], 
-                 source->modelPSF->params->data.F32[PM_PAR_SXX], source->modelPSF->params->data.F32[PM_PAR_SXY], source->modelPSF->params->data.F32[PM_PAR_SYY], 
-                 source->modelPSF->params->data.F32[PM_PAR_7]);
+        fprintf (f, "%6.1f %6.1f : %6.1f %6.1f : %8.3f %8.3f %8.3f : %f : %f %f %f : %f\n",
+                 source->peak->xf, source->peak->yf,
+                 source->modelPSF->params->data.F32[PM_PAR_XPOS], source->modelPSF->params->data.F32[PM_PAR_YPOS],
+                 source->psfMag, source->apMag, source->errMag,
+                 source->modelPSF->params->data.F32[PM_PAR_I0],
+                 source->modelPSF->params->data.F32[PM_PAR_SXX], source->modelPSF->params->data.F32[PM_PAR_SXY], source->modelPSF->params->data.F32[PM_PAR_SYY],
+                 source->modelPSF->params->data.F32[PM_PAR_7]);
 # endif
-        if (!isfinite(source->psfMag)) psAbort ("nan in psfMag");
-        if (!isfinite(source->errMag)) psAbort ("nan in errMag");
-        if (!isfinite(source->apMag)) psAbort ("nan in apMag");
-        if (!isfinite(model->params->data.F32[PM_PAR_XPOS])) psAbort ("nan in xPos");
-        if (!isfinite(model->params->data.F32[PM_PAR_YPOS])) psAbort ("nan in yPos");
+        if (!isfinite(source->psfMag)) psAbort ("nan in psfMag");
+        if (!isfinite(source->errMag)) psAbort ("nan in errMag");
+        if (!isfinite(source->apMag)) psAbort ("nan in apMag");
+        if (!isfinite(model->params->data.F32[PM_PAR_XPOS])) psAbort ("nan in xPos");
+        if (!isfinite(model->params->data.F32[PM_PAR_YPOS])) psAbort ("nan in yPos");
 
         psVectorAppend (mag, source->psfMag);
@@ -253 +251 @@
     if (Npsf < APTREND_NSTAR_MIN) {
         psWarning("Only %d valid aperture residual sources (need %d), giving up", Npsf, APTREND_NSTAR_MIN);
-        goto escape;
-    }
-
-    // this is a bit tricky, because we have two cases (MAP vs POLY), and they have a different 
-    // definition for 'order' (order_MAP = order_POLY + 1).  in addition, we have a 
+        goto escape;
+    }
+
+    // this is a bit tricky, because we have two cases (MAP vs POLY), and they have a different
+    // definition for 'order' (order_MAP = order_POLY + 1).  in addition, we have a
     // user-specified MAX order, which we should respect, regardless of the mode
 
@@ -270 +268 @@
     pmTrend2DMode mode = PM_TREND_MAP;
     if (mode == PM_TREND_MAP) {
-        MaxOrderForStars ++;
-    } 
+        MaxOrderForStars ++;
+    }
     APTREND_ORDER_MAX = PS_MIN (APTREND_ORDER_MAX, MaxOrderForStars);
 
@@ -283 +281 @@
     int NY = readout->image->numRows;
     for (int i = 1; i <= APTREND_ORDER_MAX; i++) {
-        
-        int Nx, Ny;
-        if (NX > NY) {
-            Nx = i;
-            Ny = PS_MAX (1, (int)(i * (NY / (float)(NX)) + 0.5));
-        } else {
-            Ny = i;
-            Nx = PS_MAX (1, (int)(i * (NX / (float)(NY)) + 0.5));
-        }
-
-        float errorFloor;
+
+        int Nx, Ny;
+        if (NX > NY) {
+            Nx = i;
+            Ny = PS_MAX (1, (int)(i * (NY / (float)(NX)) + 0.5));
+        } else {
+            Ny = i;
+            Nx = PS_MAX (1, (int)(i * (NX / (float)(NY)) + 0.5));
+        }
+
+        float errorFloor;
         pmTrend2D *apTrend = psphotApResidTrend (&errorFloor, readout, Nx, Ny, xPos, yPos, apResid, dMag);
-        if (!apTrend) {
-            continue;
-        }
-
-        // apply ApTrend results
-        // float xc = 0.5*readout->image->numCols + readout->image->col0 + 0.5;
-        // float yc = 0.5*readout->image->numRows + readout->image->row0 + 0.5;
-        // float ApResid = pmTrend2DEval (psf->ApTrend, xc, yc); // ap-fit at chip center
-        // if (!isfinite(ApResid)) psAbort("nan apresid @ center");
+        if (!apTrend) {
+            continue;
+        }
+
+        // apply ApTrend results
+        // float xc = 0.5*readout->image->numCols + readout->image->col0 + 0.5;
+        // float yc = 0.5*readout->image->numRows + readout->image->row0 + 0.5;
+        // float ApResid = pmTrend2DEval (psf->ApTrend, xc, yc); // ap-fit at chip center
+        // if (!isfinite(ApResid)) psAbort("nan apresid @ center");
 
         // store the minimum errorFloor and best ApTrend to keep
         if (errorFloor < errorFloorMin) {
             errorFloorMin = errorFloor;
-            psFree (psf->ApTrend);
-            psf->ApTrend = psMemIncrRefCounter(apTrend);
-        }
-        psFree (apTrend);
+            psFree (psf->ApTrend);
+            psf->ApTrend = psMemIncrRefCounter(apTrend);
+        }
+        psFree (apTrend);
     }
     if (psf->ApTrend == NULL) {
         psWarning("Failed to find a valid aperture residual value");
-        goto escape;
+        goto escape;
     }
 
@@ -383 +381 @@
     if (!pmTrend2DFit (apTrend, mask, 0xff, xPos, yPos, apResid, dMagSoft)) {
         psWarning("Failed to fit trend for %d x %d map", Nx, Ny);
-        psFree (apTrend);
-        return NULL;
+        psFree (apTrend);
+        return NULL;
     }
     if (apTrend->mode == PM_TREND_MAP) {
-        // p_psImagePrint (2, apTrend->map->map, "ApTrend Before"); // XXX TEST:
-        psImageMapRepair (apTrend->map->map);
-        // p_psImagePrint (2, apTrend->map->map, "ApTrend After"); // XXX TEST:
+        // p_psImagePrint (2, apTrend->map->map, "ApTrend Before"); // XXX TEST:
+        psImageMapRepair (apTrend->map->map);
+        // p_psImagePrint (2, apTrend->map->map, "ApTrend After"); // XXX TEST:
     }
 
@@ -400 +398 @@
     if (!isfinite(*apResidSysErr)) {
         psWarning("Failed to find systematic error for %d x %d map", Nx, Ny);
-        psFree (apTrend);
-        return NULL;
+        psFree (apTrend);
+        return NULL;
     }
 
@@ -408 +406 @@
 
     if (psTraceGetLevel("psphot") >= 4) {
-        char filename[64];
-        snprintf (filename, 64, "apresid.%dx%d.dat", Nx, Ny);
+        char filename[64];
+        snprintf (filename, 64, "apresid.%dx%d.dat", Nx, Ny);
         FILE *dumpFile = fopen (filename, "w");
         for (int i = 0; i < xPos->n; i++) {
@@ -457 +455 @@
         }
 
-        // clear the mask bit and set the circular mask pixels
-        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
-        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, source->apRadius, "OR", markVal);
-
-        bool status = pmSourceMagnitudes (source, psf, photMode, maskVal);
-
-        // clear the mask bit 
-        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
+        // clear the mask bit and set the circular mask pixels
+        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
+        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, source->apRadius, "OR", markVal);
+
+        bool status = pmSourceMagnitudes (source, psf, photMode, maskVal);
+
+        // clear the mask bit
+        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
 
         // re-subtract the object, leave local sky
         pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
 
-        if (!status) {
-            Nskip ++;
-            psTrace ("psphot", 3, "skip : bad source mag");
-            continue;
-        }
-    
-        if (!isfinite(source->apMag) || !isfinite(source->psfMag)) {
-            Nfail ++;
-            psTrace ("psphot", 3, "fail : nan mags : %f %f", source->apMag, source->psfMag);
-            continue;
-        }
-        source->mode |= PM_SOURCE_MODE_AP_MAGS;
+        if (!status) {
+            Nskip ++;
+            psTrace ("psphot", 3, "skip : bad source mag");
+            continue;
+        }
+
+        if (!isfinite(source->apMag) || !isfinite(source->psfMag)) {
+            Nfail ++;
+            psTrace ("psphot", 3, "fail : nan mags : %f %f", source->apMag, source->psfMag);
+            continue;
+        }
+        source->mode |= PM_SOURCE_MODE_AP_MAGS;
     }
 

trunk/psphot/src/psphotBlendFit.c

@@ -15 +15 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        if (!psphotBlendFitReadout (config, view, filerule, i, recipe)) {
+        if (!psphotBlendFitReadout (config, view, filerule, i, recipe)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (non-linear) for %s entry %d", filerule, i);
-            return false;
-        }
+            return false;
+        }
     }
     return true;
@@ -48 +48 @@
 
     if (!sources->n) {
-        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping blend fit");
-        return true;
+        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping blend fit");
+        return true;
     }
 
@@ -87 +87 @@
     sources = psArraySort (sources, pmSourceSortBySN);
     if (!sources->n) {
-        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping blend");
-        return true;
+        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping blend");
+        return true;
     }
 
@@ -103 +103 @@
         for (int j = 0; j < cells->n; j++) {
 
-            // allocate a job -- if threads are not defined, this just runs the job
-            psThreadJob *job = psThreadJobAlloc ("PSPHOT_BLEND_FIT");
-            psArray *newSources = psArrayAllocEmpty(16);
-
-            psArrayAdd(job->args, 1, readout);
-            psArrayAdd(job->args, 1, recipe);
-            psArrayAdd(job->args, 1, cells->data[j]); // sources
-            psArrayAdd(job->args, 1, psf);
-            psArrayAdd(job->args, 1, newSources); // return for new sources
-            psFree (newSources);
-
-            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfit
-            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Npsf
-            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Next
-            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail
-
-            if (!psThreadJobAddPending(job)) {
-                psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
-                return NULL;
-            }
-            psFree(job);
+            // allocate a job -- if threads are not defined, this just runs the job
+            psThreadJob *job = psThreadJobAlloc ("PSPHOT_BLEND_FIT");
+            psArray *newSources = psArrayAllocEmpty(16);
+
+            psArrayAdd(job->args, 1, readout);
+            psArrayAdd(job->args, 1, recipe);
+            psArrayAdd(job->args, 1, cells->data[j]); // sources
+            psArrayAdd(job->args, 1, psf);
+            psArrayAdd(job->args, 1, newSources); // return for new sources
+            psFree (newSources);
+
+            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfit
+            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Npsf
+            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Next
+            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail
+
+            if (!psThreadJobAddPending(job)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+                return NULL;
+            }
 
 # if (0)
@@ -152 +150 @@
         }
 
-        // wait for the threads to finish and manage results
-        if (!psThreadPoolWait (false)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-            return NULL;
-        }
-
-        // we have only supplied one type of job, so we can assume the types here
-        psThreadJob *job = NULL;
-        while ((job = psThreadJobGetDone()) != NULL) {
-            if (job->args->n < 1) {
-                fprintf (stderr, "error with job\n");
-            } else {
-                psScalar *scalar = NULL;
-                scalar = job->args->data[5];
-                Nfit += scalar->data.S32;
-                scalar = job->args->data[6];
-                Npsf += scalar->data.S32;
-                scalar = job->args->data[7];
-                Next += scalar->data.S32;
-                scalar = job->args->data[8];
-                Nfail += scalar->data.S32;
-
-                // add these back onto sources
-                psArray *newSources = job->args->data[4];
-                for (int j = 0; j < newSources->n; j++) {
-                    psArrayAdd (sources, 16, newSources->data[j]);
-                }
-            }
-            psFree(job);
+        // wait for the threads to finish and manage results
+        if (!psThreadPoolWait (false)) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+            return NULL;
+        }
+
+        // we have only supplied one type of job, so we can assume the types here
+        psThreadJob *job = NULL;
+        while ((job = psThreadJobGetDone()) != NULL) {
+            if (job->args->n < 1) {
+                fprintf (stderr, "error with job\n");
+            } else {
+                psScalar *scalar = NULL;
+                scalar = job->args->data[5];
+                Nfit += scalar->data.S32;
+                scalar = job->args->data[6];
+                Npsf += scalar->data.S32;
+                scalar = job->args->data[7];
+                Next += scalar->data.S32;
+                scalar = job->args->data[8];
+                Nfail += scalar->data.S32;
+
+                // add these back onto sources
+                psArray *newSources = job->args->data[4];
+                for (int j = 0; j < newSources->n; j++) {
+                    psArrayAdd (sources, 16, newSources->data[j]);
+                }
+            }
+            psFree(job);
             }
     }
@@ -278 +276 @@
                 psTrace ("psphot", 5, "source at %7.1f, %7.1f is ext", source->peak->xf, source->peak->yf);
                 Next ++;
-                source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT;
+                source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT;
                 continue;
             }
@@ -286 +284 @@
                 psTrace ("psphot", 5, "source at %7.1f, %7.1f is psf", source->peak->xf, source->peak->yf);
                 Npsf ++;
-                source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT;
+                source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT;
                 continue;
             }

trunk/psphot/src/psphotGuessModels.c

@@ -21 +21 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        if (i == chisqNum) continue; // skip chisq image
+        if (i == chisqNum) continue; // skip chisq image
         if (!psphotGuessModelsReadout (config, view, filerule, i)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed on to guess models for %s entry %d", filerule, i);
@@ -106 +106 @@
             if (!psThreadJobAddPending(job)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
                 return false;
             }
-            psFree(job);
         }
 

trunk/psphot/src/psphotMagnitudes.c

@@ -18 +18 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        if (i == chisqNum) continue; // skip chisq image
-
-        // find the currently selected readout
-        pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest
-        psAssert (file, "missing file?");
-
-        pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
-        psAssert (readout, "missing readout?");
-
-        pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
-        psAssert (detections, "missing detections?");
-
-        psArray *sources = detections->allSources;
-        psAssert (sources, "missing sources?");
-
-        pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
-        psAssert (psf, "missing psf?");
-
-        if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf)) {
+        if (i == chisqNum) continue; // skip chisq image
+
+        // find the currently selected readout
+        pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest
+        psAssert (file, "missing file?");
+
+        pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
+        psAssert (readout, "missing readout?");
+
+        pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
+        psAssert (detections, "missing detections?");
+
+        psArray *sources = detections->allSources;
+        psAssert (sources, "missing sources?");
+
+        pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
+        psAssert (psf, "missing psf?");
+
+        if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for %s entry %d", filerule, i);
-            return false;
-        }
+            return false;
+        }
     }
     return true;
@@ -50 +50 @@
 
     if (!sources->n) {
-        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping source magnitudes");
-        return true;
-    }
-        
+        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping source magnitudes");
+        return true;
+    }
+
     psTimerStart ("psphot.mags");
 
@@ -122 +122 @@
             if (!psThreadJobAddPending(job)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
                 return false;
             }
-            psFree(job);
 
 # if (0)
@@ -184 +182 @@
         }
 
-        // clear the mask bit and set the circular mask pixels
-        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
-        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, source->apRadius, "OR", markVal);
+        // clear the mask bit and set the circular mask pixels
+        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
+        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, source->apRadius, "OR", markVal);
 
         status = pmSourceMagnitudes (source, psf, photMode, maskVal); // maskVal includes markVal
         if (status && isfinite(source->apMag)) Nap ++;
 
-        // clear the mask bit 
-        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
+        // clear the mask bit
+        psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
 
         // re-subtract the object, leave local sky
@@ -273 +271 @@
             if (!psThreadJobAddPending(job)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
                 return false;
             }
-            psFree(job);
-
         }
 

trunk/psphot/src/psphotSourceStats.c

@@ -3 +3 @@
 // convert detections to sources and measure their basic properties (moments, local sky, sky
 // variance) Note: this function only generates sources for the new peaks (peak->assigned).
-// The new sources are added to any existing sources on detections->newSources.  The sources 
+// The new sources are added to any existing sources on detections->newSources.  The sources
 // on detections->allSources are ignored.
 bool psphotSourceStats (pmConfig *config, const pmFPAview *view, const char *filerule, bool setWindow)
@@ -22 +22 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        // if (i == chisqNum) continue; // skip chisq image
+        // if (i == chisqNum) continue; // skip chisq image
         if (!psphotSourceStatsReadout (config, view, filerule, i, recipe, setWindow)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for %s entry %d", filerule, i);
@@ -91 +91 @@
     // generate the array of sources, define the associated pixel
     if (!detections->newSources) {
-        detections->newSources = psArrayAllocEmpty (peaks->n);
+        detections->newSources = psArrayAllocEmpty (peaks->n);
     }
     sources = detections->newSources;
@@ -107 +107 @@
         // create a new source
         pmSource *source = pmSourceAlloc();
-        source->imageID = index;
+        source->imageID = index;
 
         // add the peak
@@ -184 +184 @@
             if (!psThreadJobAddPending(job)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to launch thread job PSPHOT_SOURCE_STATS");
-                psFree (job);
                 psFree(detections->newSources);
                 return false;
             }
-            psFree(job);
         }
 
@@ -194 +192 @@
         if (!psThreadPoolWait (false)) {
             psError(PS_ERR_UNKNOWN, false, "Failure in thread job PSPHOT_SOURCE_STATS");
-            psFree(detections->newSources);
+            psFree(detections->newSources);
             return false;
         }
@@ -315 +313 @@
             if (!psThreadJobAddPending(job)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to launch thread job PSPHOT_SOURCE_STATS");
-                psFree (job);
                 return NULL;
             }
-            psFree(job);
         }
 
@@ -380 +376 @@
         pmSource *source = sources->data[i];
 
-        if (source->tmpFlags & PM_SOURCE_TMPF_MOMENTS_MEASURED) continue;
-        source->tmpFlags |= PM_SOURCE_TMPF_MOMENTS_MEASURED;
+        if (source->tmpFlags & PM_SOURCE_TMPF_MOMENTS_MEASURED) continue;
+        source->tmpFlags |= PM_SOURCE_TMPF_MOMENTS_MEASURED;
 
         // skip faint sources for moments measurement

trunk/psphot/src/psphotTest.c

@@ -25 +25 @@
 
 bool FillImage_Threaded (psThreadJob *job) {
-    
+
     psImage *image = job->args->data[0];
     int xs = PS_SCALAR_VALUE(job->args->data[1],S32);
@@ -36 +36 @@
     psRegion region = psRegionSet (xs, xs + dx, ys, ys + dy);
     for (int i = 0; i < 100; i++) {
-        psImage *subset = psImageSubset (image, region);
-        psImageInit (subset, value + i);
-        psFree (subset);
+        psImage *subset = psImageSubset (image, region);
+        psImageInit (subset, value + i);
+        psFree (subset);
     }
     return true;
@@ -46 +46 @@
 
     if (argc != 3) {
-        fprintf (stderr, "USAGE: psphotTest (output.fits) (nThreads)\n");
-        exit (2);
+        fprintf (stderr, "USAGE: psphotTest (output.fits) (nThreads)\n");
+        exit (2);
     }
 
@@ -62 +62 @@
 
     for (int ix = 0; ix < 1000; ix += 100) {
-        for (int iy = 0; iy < 1000; iy += 100) {
+        for (int iy = 0; iy < 1000; iy += 100) {
 
-            // allocate a job -- if threads are not defined, this just runs the job
-            psThreadJob *job = psThreadJobAlloc ("FILL_IMAGE");
+            // allocate a job -- if threads are not defined, this just runs the job
+            psThreadJob *job = psThreadJobAlloc ("FILL_IMAGE");
 
-            psArrayAdd(job->args, 1, image);
-            PS_ARRAY_ADD_SCALAR(job->args, ix, PS_TYPE_S32);
-            PS_ARRAY_ADD_SCALAR(job->args, iy, PS_TYPE_S32);
-            PS_ARRAY_ADD_SCALAR(job->args, 100, PS_TYPE_S32);
-            PS_ARRAY_ADD_SCALAR(job->args, 100, PS_TYPE_S32);
-            PS_ARRAY_ADD_SCALAR(job->args, ix + iy, PS_TYPE_S32);
+            psArrayAdd(job->args, 1, image);
+            PS_ARRAY_ADD_SCALAR(job->args, ix, PS_TYPE_S32);
+            PS_ARRAY_ADD_SCALAR(job->args, iy, PS_TYPE_S32);
+            PS_ARRAY_ADD_SCALAR(job->args, 100, PS_TYPE_S32);
+            PS_ARRAY_ADD_SCALAR(job->args, 100, PS_TYPE_S32);
+            PS_ARRAY_ADD_SCALAR(job->args, ix + iy, PS_TYPE_S32);
 
-            // FillImage (image, ix, iy, 100, 100, ix + iy);
+            // FillImage (image, ix, iy, 100, 100, ix + iy);
 
-            if (!psThreadJobAddPending(job)) {
-                fprintf (stderr, "failure to run FillImage(1)");
-                psFree (job);
-                exit (1);
-            }
-            psFree(job);
-        }
+            if (!psThreadJobAddPending(job)) {
+                fprintf (stderr, "failure to run FillImage(1)");
+                exit (1);
+            }
+        }
     }
 
@@ -88 +86 @@
     // wait for the threads to finish and manage results
     if (!psThreadPoolWait (true)) {
-        fprintf (stderr, "failure to run FillImage (2)");
-        exit (1);
+        fprintf (stderr, "failure to run FillImage (2)");
+        exit (1);
     }
 

trunk/pswarp/src/pswarpTransformReadout.c

@@ -132 +132 @@
                 return false;
             }
-            psFree(job);
             psFree(args);
         }