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

Timestamp:

Sep 24, 2007, 11:27:58 AM (19 years ago)

Author:

eugene

Message:

update from eam_branch_20070921

File:

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

trunk/psphot/src/psphotApResid.c

@@ -1 +1 @@
 # include "psphotInternal.h"
 
+# define SKIPSTAR(MSG) { psTrace ("psphot", 3, "invalid : %s", MSG); continue; }
 // measure the aperture residual statistics and 2D variations
 bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal, psMaskType mark)
@@ -63 +64 @@
         model = source->modelPSF;
 
-        if (source->type != PM_SOURCE_TYPE_STAR) continue;
-        if (source->mode &  PM_SOURCE_MODE_SATSTAR) continue;
-        if (source->mode &  PM_SOURCE_MODE_BLEND) continue;
-        if (source->mode &  PM_SOURCE_MODE_FAIL) continue;
-        if (source->mode &  PM_SOURCE_MODE_POOR) continue;
+        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");
 
         // get growth-corrected, apTrend-uncorrected magnitudes in scaled apertures
@@ -73 +74 @@
         if (!pmSourceMagnitudes (source, psf, photMode, maskVal, mark)) {
             Nskip ++;
+            psTrace ("psphot", 3, "skip : bad source mag");
             continue;
         }
 
-        if (!isfinite(source->apMag) || !isfinite(source->apMag)) {
+        if (!isfinite(source->apMag) || !isfinite(source->psfMag)) {
             Nfail ++;
+            psTrace ("psphot", 3, "fail : nan mags : %f %f", source->apMag, source->psfMag);
             continue;
         }
@@ -87 +90 @@
         if ((MAX_AP_OFFSET > 0) && (fabs(dap) > MAX_AP_OFFSET)) {
             Nfail ++;
+            psTrace ("psphot", 3, "fail : bad dap %f %f", dap, MAX_AP_OFFSET);
             continue;
         }
@@ -99 +103 @@
         dMag->data.F32[Npsf] = model->dparams->data.F32[PM_PAR_I0] / model->params->data.F32[PM_PAR_I0];
 
+        psVectorExtend (mag,     100, 1);
         psVectorExtend (mask,    100, 1);
         psVectorExtend (xPos,    100, 1);
@@ -146 +151 @@
     }
 
+    // XXX catch error condition 
     psphotApResidTrend (readout, psf, Npsf, entryMin, &errorScale, &errorFloor, mask, xPos, yPos, apResid, dMag);
 
@@ -152 +158 @@
     psVector *apResidRes = (psVector *) psBinaryOp (NULL, (void *) apResid, "-", (void *) apResidFit);
     psVector *dMagSys = (psVector *) psBinaryOp (NULL, (void *) dMag, "*", (void *) psScalarAlloc(errorScale, PS_TYPE_F32));
-
-    psphotSaveImage (NULL, psf->ApTrend->map->map, "apresid.map.fits");
 
     if (psTraceGetLevel("psphot") >= 5) {
@@ -170 +174 @@
     float xc = 0.5*readout->image->numCols + readout->image->col0 + 0.5;
     float yc = 0.5*readout->image->numRows + readout->image->row0 + 0.5;
+
     psf->ApResid  = pmTrend2DEval (psf->ApTrend, xc, yc); // ap-fit at chip center
     psf->dApResid = errorFloor;
@@ -206 +211 @@
 */
 
-bool psphotMagErrorScale (float *errorScale, float *errorFloor, psVector *dMag, psVector *dap, int nGroup) {
+bool psphotMagErrorScale (float *errorScale, float *errorFloor, psVector *dMag, psVector *dap, psVector *mask, int nGroup) {
 
     psStats *statsS = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV);
@@ -225 +230 @@
     psVector *dMSubset = psVectorAllocEmpty (nGroup, PS_TYPE_F32);
     psVector *dASubset = psVectorAllocEmpty (nGroup, PS_TYPE_F32);
+    psVector *mkSubset = psVectorAllocEmpty (nGroup, PS_TYPE_U8);
 
     int n = 0;
@@ -233 +239 @@
             dMSubset->data.F32[j] = dMag->data.F32[N];
             dASubset->data.F32[j] = dap->data.F32[N];
+            mkSubset->data.U8[j]  = mask->data.U8[N];
         }
         dMSubset->n = j;
         dASubset->n = j;
+        mkSubset->n = j;
 
         psStatsInit (statsS);
         psStatsInit (statsM);
 
-        psVectorStats (statsS, dASubset, NULL, NULL, 0xff);
-        psVectorStats (statsM, dMSubset, NULL, NULL, 0xff);
+        psVectorStats (statsS, dASubset, NULL, mkSubset, 0xff);
+        psVectorStats (statsM, dMSubset, NULL, mkSubset, 0xff);
 
         dSo->data.F32[i] = statsS->robustStdev;
@@ -250 +258 @@
     psFree (dMSubset);
     psFree (dASubset);
+    psFree (mkSubset);
     
     psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
@@ -286 +295 @@
     }
 
+    // the mask marks the values not used to calculate the ApTrend
+    psVectorInit (mask, 0);
+
     // XXX stats structure for use by ApTrend : make parameters user setable
     psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV);
@@ -305 +317 @@
     // XXX this is a bit arbitrary, but it forces ~3 stars from the bright bin per spatial bin
     int nGroup = PS_MAX (3*Nx*Ny, 10);
-    psphotMagErrorScale (errorScale, errorFloor, dMag, apResidRes, nGroup);
+    psphotMagErrorScale (errorScale, errorFloor, dMag, apResidRes, mask, nGroup);
 
     psLogMsg ("psphot.apresid", PS_LOG_INFO, "result of %d x %d grid (%d stars per bin)\n", Nx, Ny, nGroup);