Changeset 23218

Timestamp:

Mar 6, 2009, 3:11:24 PM (17 years ago)

Author:

Paul Price

Message:

Add additional diagnostic: number of sources rejected as too faint.

Location:

trunk/psphot/src

Files:

• psphotSetThreads.c (modified) (1 diff)
• psphotSourceStats.c (modified) (14 diffs)

trunk/psphot/src/psphotSetThreads.c

@@ -25 +25 @@
     psFree(task);
 
-    task = psThreadTaskAlloc("PSPHOT_SOURCE_STATS", 4);
+    task = psThreadTaskAlloc("PSPHOT_SOURCE_STATS", 5);
     task->function = &psphotSourceStats_Threaded;
     psThreadTaskAdd(task);

trunk/psphot/src/psphotSourceStats.c

@@ -15 +15 @@
     int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
     if (!status) {
-        nThreads = 0;
+        nThreads = 0;
     }
 
@@ -41 +41 @@
         pmSource *source = pmSourceAlloc();
 
-        // add the peak
+        // add the peak
         source->peak = psMemIncrRefCounter(peak);
 
-        // allocate space for moments
+        // allocate space for moments
         source->moments = pmMomentsAlloc();
 
@@ -50 +50 @@
         pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
 
-        peak->assigned = true;
-        psArrayAdd (sources, 100, source);
-        psFree (source);
+        peak->assigned = true;
+        psArrayAdd (sources, 100, source);
+        psFree (source);
     }
 
     if (!strcasecmp (breakPt, "PEAKS")) {
-        psLogMsg ("psphot", PS_LOG_INFO, "%ld sources : %f sec\n", sources->n, psTimerMark ("psphot.stats"));
-        psLogMsg ("psphot", PS_LOG_INFO, "break point PEAKS, skipping MOMENTS\n");
-        psphotVisualShowMoments (sources);
-        return sources;
+        psLogMsg ("psphot", PS_LOG_INFO, "%ld sources : %f sec\n", sources->n, psTimerMark ("psphot.stats"));
+        psLogMsg ("psphot", PS_LOG_INFO, "break point PEAKS, skipping MOMENTS\n");
+        psphotVisualShowMoments (sources);
+        return sources;
     }
 
@@ -65 +65 @@
     int Nfail = 0;
     int Nmoments = 0;
+    int Nfaint = 0;
 
     // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)
@@ -74 +75 @@
     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_SOURCE_STATS");
-
-            psArrayAdd(job->args, 1, cells->data[j]); // sources
-            psArrayAdd(job->args, 1, recipe);
-            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments
-            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);
+        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_SOURCE_STATS");
+
+            psArrayAdd(job->args, 1, cells->data[j]); // sources
+            psArrayAdd(job->args, 1, recipe);
+            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments
+            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);
 
 # if (0)
-                int nfail = 0;
-                int nmoments = 0;
-                if (!psphotSourceStats_Unthreaded (&nfail, &nmoments, cells->data[j], recipe)) {
-                    psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
-                    return NULL;
-                }
-                Nfail += nfail;
-                Nmoments += nmoments;
+                int nfail = 0;
+                int nmoments = 0;
+                if (!psphotSourceStats_Unthreaded (&nfail, &nmoments, cells->data[j], recipe)) {
+                    psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
+                    return NULL;
+                }
+                Nfail += nfail;
+                Nmoments += nmoments;
 # endif
-        }
-
-        // 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[2];
-                Nmoments += scalar->data.S32;
-                scalar = job->args->data[3];
-                Nfail += scalar->data.S32;
-            }
-            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[2];
+                Nmoments += scalar->data.S32;
+                scalar = job->args->data[3];
+                Nfail += scalar->data.S32;
+                scalar = job->args->data[4];
+                Nfaint += scalar->data.S32;
+            }
+            psFree(job);
+        }
     }
 
     psFree (cellGroups);
 
-    psLogMsg ("psphot", PS_LOG_INFO, "%ld sources, %d moments, %d failed: %f sec\n", sources->n, Nmoments, Nfail, psTimerMark ("psphot.stats"));
+    psLogMsg ("psphot", PS_LOG_INFO, "%ld sources, %d moments, %d faint, %d failed: %f sec\n", sources->n, Nmoments, Nfaint, Nfail, psTimerMark ("psphot.stats"));
 
     psphotVisualShowMoments (sources);
@@ -166 +169 @@
     int Nfail = 0;
     int Nmoments = 0;
+    int Nfaint = 0;
     for (int i = 0; i < sources->n; i++) {
         pmSource *source = sources->data[i];
@@ -171 +175 @@
         // skip faint sources for moments measurement
         if (source->peak->SN < MIN_SN) {
-            source->mode |= PM_SOURCE_MODE_BELOW_MOMENTS_SN;
+            source->mode |= PM_SOURCE_MODE_BELOW_MOMENTS_SN;
+            Nfaint++;
             continue;
         }
@@ -179 +184 @@
         status = pmSourceLocalSky (source, PS_STAT_SAMPLE_MEDIAN, INNER, maskVal, markVal);
         if (!status) {
-            source->mode |= PM_SOURCE_MODE_SKY_FAILURE;
-            psErrorClear(); // XXX re-consider the errors raised here
-            Nfail ++;
-            continue;
+            source->mode |= PM_SOURCE_MODE_SKY_FAILURE;
+            psErrorClear(); // XXX re-consider the errors raised here
+            Nfail ++;
+            continue;
         }
 
@@ -189 +194 @@
         status = pmSourceLocalSkyVariance (source, PS_STAT_SAMPLE_MEDIAN, INNER, maskVal, markVal);
         if (!status) {
-            source->mode |= PM_SOURCE_MODE_SKYVAR_FAILURE;
-            Nfail ++;
-            psErrorClear();
-            continue;
+            source->mode |= PM_SOURCE_MODE_SKYVAR_FAILURE;
+            Nfail ++;
+            psErrorClear();
+            continue;
         }
 
@@ -208 +213 @@
         status = pmSourceMoments (source, BIG_RADIUS);
         if (status) {
-            source->mode |= PM_SOURCE_MODE_BIG_RADIUS;
+            source->mode |= PM_SOURCE_MODE_BIG_RADIUS;
             Nmoments ++;
             continue;
         }
 
-        source->mode |= PM_SOURCE_MODE_MOMENTS_FAILURE;
+        source->mode |= PM_SOURCE_MODE_MOMENTS_FAILURE;
         Nfail ++;
         psErrorClear();
@@ -225 +230 @@
     scalar = job->args->data[3];
     scalar->data.S32 = Nfail;
-    
+
+    scalar = job->args->data[4];
+    scalar->data.S32 = Nfaint;
+
     return true;
 }
@@ -268 +276 @@
         status = pmSourceLocalSky (source, PS_STAT_SAMPLE_MEDIAN, INNER, maskVal, markVal);
         if (!status) {
-            psErrorClear(); // XXX re-consider the errors raised here
-            Nfail ++;
-            continue;
+            psErrorClear(); // XXX re-consider the errors raised here
+            Nfail ++;
+            continue;
         }
 
@@ -277 +285 @@
         status = pmSourceLocalSkyVariance (source, PS_STAT_SAMPLE_MEDIAN, INNER, maskVal, markVal);
         if (!status) {
-            Nfail ++;
-            psErrorClear();
-            continue;
+            Nfail ++;
+            psErrorClear();
+            continue;
         }
 
@@ -307 +315 @@
     *nmoments = Nmoments;
     *nfail = Nfail;
-    
+
     return true;
 }
-# endif 
+# endif