Changeset 13514

Timestamp:

May 24, 2007, 2:16:03 PM (19 years ago)

Author:

magnier

Message:

adding some test output code; removed asymmetrical clipping for apresd stats

File:

• trunk/psModules/src/objects/pmPSFtry.c (modified) (6 diffs)

trunk/psModules/src/objects/pmPSFtry.c

@@ -5 +5 @@
  *  @author EAM, IfA
  *
- *  @version $Revision: 1.40 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-05-18 17:15:35 $
+ *  @version $Revision: 1.41 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-05-25 00:16:03 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -168 +168 @@
         if (!status) {
             psfTry->mask->data.U8[i] = PSFTRY_MASK_PSF_FAIL;
+            psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : failed PSF fit\n", i, source->peak->x, source->peak->y);
             continue;
         }
@@ -174 +175 @@
         if (!pmSourceMagnitudes (source, psfTry->psf, PM_SOURCE_PHOT_INTERP)) {
             psfTry->mask->data.U8[i] = PSFTRY_MASK_BAD_PHOT;
+            psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : poor photometry\n", i, source->peak->x, source->peak->y);
             continue;
         }
@@ -261 +263 @@
     }
 
+    // XXX this analysis of the apResid statistics is only approximate.  The fitted magnitudes
+    // measure at this point (in the PSF fit) use Poisson errors, and are thus biased as a
+    // function of magnitude.  We re-measure the apResid statistics later in psphot using the
+    // linear, constant-error fitting.  Do not reject outliers with excessive vigor here.
+
     // use 3hi/1lo sigma clipping on the r2rflux vs metric fit
     psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV);
-    stats->min = 1.0;
-    stats->max = 3.0;
+    stats->clipSigma = 3.0;
     stats->clipIter = 3;
+
+    // XXX we used to include an asymmetric clipping in order to toss out contaminated stars.
+    // test this on the very crowded field data.
+    // stats->min = 1.0;
+    // stats->max = 3.0;
 
     // fit ApTrend only to r2rflux, ignore x,y,flux variations for now
     // linear clipped fit of ApResid to r2rflux
     psPolynomial1D *poly = psPolynomial1DAlloc (PS_POLYNOMIAL_ORD, 1);
-
-    // XXX test to only fit a constant offset (no SKYBIAS)
-    poly->mask[1] = 1;
+    poly->mask[1] = 1; // fit only a constant offset (no SKYBIAS)
 
     bool result = psVectorClipFitPolynomial1D (poly, stats, psfTry->mask, PSFTRY_MASK_ALL, psfTry->metric, psfTry->metricErr, r2rflux);
@@ -284 +293 @@
         return false;
     }
-    psLogMsg ("pmPSFtryMetric", 4, "fit stats: %f +/- %f\n", stats->robustMedian, stats->robustStdev);
-    psLogMsg ("pmPSFtryMetric", 4, "apresid: %f +/- %f\n", poly->coeff[0], poly->coeffErr[0]);
+    psLogMsg ("pmPSFtryMetric", 4, "apresid: %f +/- %f; keeping %ld of %ld psf stars\n", poly->coeff[0], stats->robustStdev, stats->clippedNvalues, psfTry->sources->n);
 
     // XXX test dump of fitted model (dump when tracing?)
@@ -403 +411 @@
     for (int i = 0; i < stats->clipIter; i++) {
         psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E0], stats, psfTry->mask, 0xff, e0, dz, x, y);
+        psTrace ("psModules.pmPSFtry", 4, "clipped E0 : keeping %ld of %ld\n", stats->clippedNvalues, e0->n);
         psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E1], stats, psfTry->mask, 0xff, e1, dz, x, y);
+        psTrace ("psModules.pmPSFtry", 4, "clipped E1 : keeping %ld of %ld\n", stats->clippedNvalues, e1->n);
         psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E2], stats, psfTry->mask, 0xff, e2, dz, x, y);
-    }
+        psTrace ("psModules.pmPSFtry", 4, "clipped E2 : keeping %ld of %ld\n", stats->clippedNvalues, e2->n);
+    }
+
+    // XXX temporary dump of the psf parameters
+    if (psTraceGetLevel("psModules.objects") >= 4) {
+        FILE *f = fopen ("pol.dat", "w");
+        for (int i = 0; i < e0->n; i++) {
+            fprintf (f, "%f %f  :  %f %f %f  : %d\n", 
+                     x->data.F64[i], y->data.F64[i], 
+                     e0->data.F64[i], e1->data.F64[i], e2->data.F64[i], psfTry->mask->data.U8[i]);
+        }
+        fclose (f);
+    }
+
     psFree (e0);
     psFree (e1);