Changeset 21359 for trunk/psphot/src/psphotGuessModels.c

Timestamp:

Feb 5, 2009, 3:12:59 PM (17 years ago)

Author:

eugene

Message:

defining test unthreaded versions for these functions to avoid hitting the psThreadJob array; remove eventually

File:

• trunk/psphot/src/psphotGuessModels.c (modified) (4 diffs)

trunk/psphot/src/psphotGuessModels.c

@@ -22 +22 @@
 }
 
+bool psphotGuessModel_Unthreaded (pmReadout *readout, psArray *sources, pmPSF *psf, psImageMaskType maskVal, psImageMaskType markVal);
+
 // construct an initial PSF model for each object
 bool psphotGuessModels (pmConfig *config, pmReadout *readout, psArray *sources, pmPSF *psf) {
@@ -38 +40 @@
         nThreads = 0;
     }
-    // nThreads = 0; // XXX until testing is complete, do not thread this function
 
     // bit-masks to test for good/bad pixels
@@ -66 +67 @@
         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_GUESS_MODEL");
-            psArrayAdd(job->args, 1, readout);
-            psArrayAdd(job->args, 1, cells->data[j]); // sources
-            psArrayAdd(job->args, 1, psf);
-
-            // XXX change these to use abstract mask type info
-            PS_ARRAY_ADD_SCALAR(job->args, maskVal,  PS_TYPE_IMAGE_MASK);
-            PS_ARRAY_ADD_SCALAR(job->args, markVal,  PS_TYPE_IMAGE_MASK);
-
-            if (!psThreadJobAddPending(job)) {
+            if (nThreads) {
+                // allocate a job -- if threads are not defined, this just runs the job
+                psThreadJob *job = psThreadJobAlloc ("PSPHOT_GUESS_MODEL");
+                psArrayAdd(job->args, 1, readout);
+                psArrayAdd(job->args, 1, cells->data[j]); // sources
+                psArrayAdd(job->args, 1, psf);
+
+                // XXX change these to use abstract mask type info
+                PS_ARRAY_ADD_SCALAR(job->args, maskVal,  PS_TYPE_IMAGE_MASK);
+                PS_ARRAY_ADD_SCALAR(job->args, markVal,  PS_TYPE_IMAGE_MASK);
+
+                if (!psThreadJobAddPending(job)) {
+                    psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+                    psFree (job);
+                    return false;
+                }
+                psFree(job);
+            } else {
+                if (!psphotGuessModel_Unthreaded (readout, cells->data[j], psf, maskVal, markVal)) {
+                    psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+                    return false;
+                }
+            }
+        }
+
+        if (nThreads) {
+            // wait for the threads to finish and manage results
+            // wait here for the threaded jobs to finish
+            // fprintf (stderr, "wait for threads (%d, %d)\n", jx, jy);
+            if (!psThreadPoolWait (false)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                psFree (job);
                 return false;
             }
-            psFree(job);
-        }
-
-        // wait for the threads to finish and manage results
-        // wait here for the threaded jobs to finish
-        // fprintf (stderr, "wait for threads (%d, %d)\n", jx, jy);
-        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) {
-            // we have no returned data from this operation
-            if (job->args->n < 1) {
-                fprintf (stderr, "error with job\n");
-            }
-            psFree(job);
+
+            // we have only supplied one type of job, so we can assume the types here
+            psThreadJob *job = NULL;
+            while ((job = psThreadJobGetDone()) != NULL) {
+                // we have no returned data from this operation
+                if (job->args->n < 1) {
+                    fprintf (stderr, "error with job\n");
+                }
+                psFree(job);
+            }
         }
     }
@@ -211 +221 @@
 }
 
+// construct models only for sources in the specified region
+bool psphotGuessModel_Unthreaded (pmReadout *readout, psArray *sources, pmPSF *psf, psImageMaskType maskVal, psImageMaskType markVal) {
+
+    int nSrc = 0;
+
+    for (int i = 0; i < sources->n; i++) {
+        pmSource *source = sources->data[i];
+
+        // XXXX this is just for a test: use this to mark sources for which the model is measured
+        // check later that all are used.
+        source->mode |= PM_SOURCE_MODE_EXT_LIMIT;
+
+        // skip non-astronomical objects (very likely defects)
+        if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
+        if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
+        if (!source->peak) continue;
+
+        nSrc ++;
+        
+        // XXX if a source is faint, it will not have moments measured.
+        // it must be modelled as a PSF.  In this case, we need to use
+        // the peak centroid to get the coordinates and get the peak flux
+        // from the image?
+        pmModel *modelEXT;
+        if (!source->moments) {
+            modelEXT = wildGuess(source, psf);
+        } else {
+            // use the source moments, etc to guess basic model parameters
+            modelEXT = pmSourceModelGuess (source, psf->type); // ALLOC
+            if (!modelEXT) {
+                modelEXT = wildGuess(source, psf);
+            }
+            // these valuse are set in pmSourceModelGuess, should this rule be in there as well?
+            if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
+                modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->Mx;
+                modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->My;
+            } else {
+                modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
+                modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
+            }
+        }
+
+        // set PSF parameters for this model (apply 2D shape model)
+        pmModel *modelPSF = pmModelFromPSF (modelEXT, psf); // ALLOC
+        if (modelPSF == NULL) {
+            psError(PSPHOT_ERR_PSF, false,
+                    "Failed to determine PSF model at r,c = (%d,%d); trying centre of image",
+                    source->peak->y, source->peak->x);
+            //
+            // Try the centre of the image
+            //
+            modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*readout->image->numCols;
+            modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*readout->image->numRows;
+            modelPSF = pmModelFromPSF (modelEXT, psf);
+            if (modelPSF == NULL) {
+                psError(PSPHOT_ERR_PSF, false,
+                        "Failed to determine PSF model at centre of image");
+                psFree(modelEXT);
+                return false;
+            }
+
+            source->mode |= PM_SOURCE_MODE_BADPSF;
+        }
+        psFree (modelEXT);
+
+        // XXX need to define the guess flux?
+        // set the fit radius based on the object flux limit and the model
+        // this function affects the mask pixels
+        psphotCheckRadiusPSF (readout, source, modelPSF, markVal);
+
+        // set the source PSF model
+        source->modelPSF = modelPSF;
+        source->modelPSF->residuals = psf->residuals;
+
+        pmSourceCacheModel (source, maskVal);
+
+    }
+
+    return true;
+}