Changeset 21166 for trunk/psphot/src/psphotApResid.c

Timestamp:

Jan 24, 2009, 10:54:29 AM (17 years ago)

Author:

eugene

Message:

clean up threading model for psphotGuessModel (defined standard way to split sources into regions); extending the threading to psphotMagnitudes, psphotApReset, psphotBlendFit, psphotSourceStats

File:

• trunk/psphot/src/psphotApResid.c (modified) (7 diffs)

trunk/psphot/src/psphotApResid.c

@@ -4 +4 @@
 // measure the aperture residual statistics and 2D variations
 
-bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf)
+bool psphotApResid (pmConfig *config, pmReadout *readout, psArray *sources, pmPSF *psf)
 {
     int Nfail = 0;
@@ -13 +13 @@
     pmSource *source;
 
-    PS_ASSERT_PTR_NON_NULL(psf, false);
+    PS_ASSERT_PTR_NON_NULL(config, false);
     PS_ASSERT_PTR_NON_NULL(readout, false);
     PS_ASSERT_PTR_NON_NULL(sources, false);
-    PS_ASSERT_PTR_NON_NULL(recipe, false);
+    PS_ASSERT_PTR_NON_NULL(psf, false);
 
     psTimerStart ("psphot.apresid");
+
+    // select the appropriate recipe information
+    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
+    assert (recipe);
+
+    // determine the number of allowed threads
+    int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
+    if (!status) {
+        nThreads = 0;
+    }
+    // nThreads = 0; // XXX until testing is complete, do not thread this function
 
     bool measureAptrend = psMetadataLookupBool (&status, recipe, "MEASURE.APTREND");
@@ -70 +81 @@
     pmSourceMagnitudesInit (recipe);
 
+    // threaded measurement of the source magnitudes
+    // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)
+    int Cx = 1, Cy = 1;
+    psphotChooseCellSizes (&Cx, &Cy, readout, nThreads);
+
+    psArray *cellGroups = psphotAssignSources (Cx, Cy, sources);
+
+    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_U8); // XXX change this to use abstract mask type info
+
+            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[4];
+                Nskip += scalar->data.S32;
+                scalar = job->args->data[5];
+                Nfail += scalar->data.S32;
+            }
+            psFree(job);
+        }
+    }
+
+    psFree (cellGroups);
+
+    // gather the stats to assess the aperture residuals
     psVector *mask    = psVectorAllocEmpty (300, PS_TYPE_U8);
     psVector *mag     = psVectorAllocEmpty (300, PS_TYPE_F32);
@@ -78 +147 @@
     Npsf = 0;
 
-    // select all good PM_SOURCE_TYPE_STAR entries
     for (int i = 0; i < sources->n; i++) {
         source = sources->data[i];
@@ -89 +157 @@
         if (source->mode &  PM_SOURCE_MODE_POOR) SKIPSTAR ("POOR STAR");
 
-        // get growth-corrected, apTrend-uncorrected magnitudes in scaled apertures
-        // will fail if below S/N threshold or model is missing
-        if (!pmSourceMagnitudes (source, psf, photMode, maskVal)) {
-            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;
         }
@@ -206 +264 @@
 
     psLogMsg ("psphot.apresid", PS_LOG_DETAIL, "aperture residual: %f +/- %f\n", psf->ApResid, psf->dApResid);
-    psLogMsg ("psphot.apresid", PS_LOG_INFO, "measure full-frame aperture residuals for: %f sec\n", psTimerMark ("psphot.apresid"));
+    psLogMsg ("psphot.apresid", PS_LOG_INFO, "measure full-frame aperture residuals for %d sources: %f sec\n", Npsf, psTimerMark ("psphot.apresid"));
 
     psFree (mag);
@@ -392 +450 @@
 }
 
+bool psphotApResidMags_Threaded (psThreadJob *job) {
+
+    int Nskip = 0;
+    int Nfail = 0;
+
+    psScalar *scalar = NULL;
+
+    psArray *sources                = job->args->data[0];
+    pmPSF *psf                      = job->args->data[1];
+    pmSourcePhotometryMode photMode = PS_SCALAR_VALUE(job->args->data[2],S32);
+    psMaskType maskVal              = PS_SCALAR_VALUE(job->args->data[3],U8);
+
+    for (int i = 0; i < sources->n; i++) {
+        pmSource *source = (pmSource *) sources->data[i];
+
+        if (source->type != PM_SOURCE_TYPE_STAR) SKIPSTAR ("NOT STAR");
+        if (source->mode &  PM_SOURCE_MODE_SATSTAR) SKIPSTAR ("SATSTAR");
+        if (source->mode &  PM_SOURCE_MODE_BLEND) SKIPSTAR ("BLEND");
+        if (source->mode &  PM_SOURCE_MODE_FAIL) SKIPSTAR ("FAIL STAR");
+        if (source->mode &  PM_SOURCE_MODE_POOR) SKIPSTAR ("POOR STAR");
+
+        if (!pmSourceMagnitudes (source, psf, photMode, maskVal)) {
+            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;
+        }
+    }
+
+    // change the value of a scalar on the array (wrap this and put it in psArray.h)
+    scalar = job->args->data[4];
+    scalar->data.S32 = Nskip;
+
+    scalar = job->args->data[5];
+    scalar->data.S32 = Nfail;
+
+    return true;
+}