Changeset 7758

Timestamp:

Jun 29, 2006, 10:50:09 AM (20 years ago)

Author:

eugene

Message:

added breakpoints, better exit handling on psphotReadout

Location:

trunk/psphot

Files:

• doc/notes.txt (modified) (1 diff)
• src/psphot.h (modified) (4 diffs)
• src/psphotBasicDeblend.c (modified) (1 diff)
• src/psphotChoosePSF.c (modified) (3 diffs)
• src/psphotFindPeaks.c (modified) (2 diffs)
• src/psphotReadout.c (modified) (7 diffs)
• src/psphotRoughClass.c (modified) (2 diffs)
• src/psphotSourceStats.c (modified) (3 diffs)

trunk/psphot/doc/notes.txt

@@ +1 @@
+
+2006.06.28
+
+  I am adding a few minor features.  First up is the ability to break
+  the processing at interesting stages.  These will be:
+
+  - PEAKS: after sources are generated, but before moments are
+    measured
+  - MOMENTS: after moments are measured, before the PSF is measured
+  - PSFMODEL: after the psf model is measured, before sources are
+    fitted in bulk
+  - ENSEMBLE: after the linear ensemble fitting
+  - DEFAULT: complete processing
 
 2006.04.22

trunk/psphot/src/psphot.h

@@ -21 +21 @@
 
 psArray        *psphotFakeSources ();
-int             psphotSaveImage (psMetadata *header, psImage *image, char *filename);
 
 // psphotReadout functions
@@ -31 +30 @@
 pmPSF          *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe);
 bool            psphotPSFstats (pmReadout *readout, psMetadata *recipe, pmPSF *psf);
+bool            psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources);
 bool            psphotEnsemblePSF (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final);
 bool            psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf);
@@ -37 +37 @@
 bool            psphotMagnitudes (psArray *sources, psMetadata *recipe, pmPSF *psf);
 bool            psphotSkyReplace (pmConfig *config, pmFPAview *view);
+bool            psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmPSF *psf, psArray *sources);
 
 // basic support functions
@@ -58 +59 @@
 psMetadata     *psphotDefineHeader (psMetadata *recipe);
 bool            psphotWeightBias (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf);
+int             psphotSaveImage (psMetadata *header, psImage *image, char *filename);
 
 // PSF / DBL / EXT evaluation functions

trunk/psphot/src/psphotBasicDeblend.c

@@ -125 +125 @@
     psLogMsg ("psphot.deblend", 3, "identified %d blended objects (%f sec)\n", Nblend, psTimerMark ("psphot"));
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "DEBLEND")) exit (0);
-
     psFree (SN);
     psFree (index);

trunk/psphot/src/psphotChoosePSF.c

@@ -29 +29 @@
 
     stars = psArrayAlloc (sources->n);
-    // DROP stars->n = 0;
 
     // select the candidate PSF stars (pointers to original sources)
@@ -105 +104 @@
     psLogMsg ("psphot.pspsf", 3, "selected psf model %s, ApResid: %f +/- %f\n", modelName, psf->ApResid, psf->dApResid);
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "PSFFIT")) exit (0);
-
     return (psf);
 }
@@ -151 +147 @@
     return true;
 }
+
+bool psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources) {
+
+    // without the PSF model, we can only rely on the source->moments 
+    // as a measure of the FWHM values
+
+    double FWHM_X, FWHM_Y, FWHM_T;
+    int FWHM_N;
+
+    psEllipseMoments moments;
+    psEllipseAxes axes;
+
+    FWHM_N = 0;
+    FWHM_X = FWHM_Y = 0.0;
+
+    for (int i = 0; i < sources->n; i++) {
+        pmSource *source = sources->data[i];
+        if (!(source->mode & PM_SOURCE_MODE_PSFSTAR)) continue;
+        
+        // moments->Sx,Sy are roughly sigma_x,y
+        moments.x2 = PS_SQR(source->moments->Sx);
+        moments.y2 = PS_SQR(source->moments->Sy);
+        moments.xy = source->moments->Sxy;
+
+        axes = psEllipseMomentsToAxes (moments);
+
+        FWHM_X += axes.major * 2.35;
+        FWHM_Y += axes.minor * 2.35;
+        FWHM_T += axes.theta;
+        FWHM_N ++;
+    }
+
+    FWHM_X /= FWHM_N;
+    FWHM_Y /= FWHM_N;
+    FWHM_T /= FWHM_N;
+
+    psMetadataAdd (recipe, PS_LIST_TAIL, "FWHM_X", PS_DATA_F32 | PS_META_REPLACE, "PSF FWHM Major axis", FWHM_X);
+    psMetadataAdd (recipe, PS_LIST_TAIL, "FWHM_Y", PS_DATA_F32 | PS_META_REPLACE, "PSF FWHM Minor axis", FWHM_Y);
+    psMetadataAdd (recipe, PS_LIST_TAIL, "ANGLE",  PS_DATA_F32 | PS_META_REPLACE, "PSF angle",           FWHM_T);
+    
+    return true;
+}

trunk/psphot/src/psphotFindPeaks.c

@@ -28 +28 @@
     psImage *mask = readout->mask;
 
+    // XXX make this a user-option?
     // psphotSaveImage (NULL, smooth_im, "imsmooth.fits");
     // psphotSaveImage (NULL, smooth_wt, "wtsmooth.fits");
@@ -64 +65 @@
     psLogMsg ("psphot", 3, "%d peaks: %f sec\n", peaks->n, psTimerMark ("psphot"));
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "PEAKS")) exit (0);
-
     return (peaks);
 }

trunk/psphot/src/psphotReadout.c

@@ -8 +8 @@
     // find the currently selected readout
     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
+
+    // optional break-point for processing
+    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
 
     // XXX does this need to invoke I/O?
@@ -34 +37 @@
     // limit moments analysis by S/N?
     psArray *sources = psphotSourceStats (readout, recipe, peaks);
+    psFree (peaks);
+
+    if (!strcasecmp (breakPt, "PEAKS")) {
+        psphotReadoutCleanup (config, readout, recipe, NULL, sources);
+        return true;
+    }
 
     // classify sources based on moments, brightness
     // faint sources not classified?
     if (!psphotRoughClass (sources, recipe)) {
-        psFree (peaks);
         psFree (sources);
         psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
@@ -44 +52 @@
         pmFPAfileDropInternal (config->files, "PSPHOT.BACKGND");
         return false;
+    }
+    if (!strcasecmp (breakPt, "MOMENTS")) {
+        psphotReadoutCleanup (config, readout, recipe, NULL, sources);
+        return true;
     }
 
@@ -51 +63 @@
     // use bright stellar objects to measure PSF
     pmPSF *psf = psphotChoosePSF (readout, sources, recipe);
-    psphotPSFstats (readout, recipe, psf);
+    if (!strcasecmp (breakPt, "PSFMODEL")) {
+        psphotReadoutCleanup (config, readout, recipe, psf, sources);
+        return true;
+    }
 
     // linear PSF fit to peaks
     psphotEnsemblePSF (readout, sources, recipe, psf, FALSE);
+    if (!strcasecmp (breakPt, "ENSEMBLE")) {
+        psphotReadoutCleanup (config, readout, recipe, psf, sources);
+        return true;
+    }
 
     // non-linear PSF and EXT fit to brighter sources
@@ -80 +99 @@
     psphotMagnitudes (sources, recipe, psf);
 
-    // create an output header with stats results
-    psMetadata *header = psphotDefineHeader (recipe);
-
     // replace background in residual image
     psphotSkyReplace (config, view);
@@ -89 +105 @@
     psphotSourceFreePixels (sources);
     // psphotSaveImage (NULL, readout->image, "resid.fits");
+
+    // create the exported-metadata and free local data
+    psphotReadoutCleanup (config, readout, recipe, psf, sources);
+    return true;
+}
+
+bool psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmPSF *psf, psArray *sources) {
+
+    // use the psf-model to measure FWHM stats
+    if (psf) {
+        psphotPSFstats (readout, recipe, psf);
+    } else {
+        psphotMomentsStats (readout, recipe, sources);
+    }
+
+    // create an output header with stats results
+    psMetadata *header = psphotDefineHeader (recipe);
 
     // save the results of the analysis
@@ -99 +132 @@
     pmFPAfileDropInternal (config->files, "PSPHOT.BACKGND");
 
-    // free up the local copies of the data
-    psFree (peaks);
-    psFree (sources);
     psFree (psf);
     psFree (header);
+    psFree (sources);
+
     return true;
 }

trunk/psphot/src/psphotRoughClass.c

@@ -4 +4 @@
 bool psphotRoughClass (psArray *sources, psMetadata *recipe) {
 
-    bool status;
     pmPSFClump   psfClump;
 
@@ -20 +19 @@
     psphotDumpMoments (recipe, sources);
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "CLASS")) exit (0);
-
     return true;
 }

trunk/psphot/src/psphotSourceStats.c

@@ -14 +14 @@
     float OUTER    = psMetadataLookupF32 (&status, recipe, "SKY_OUTER_RADIUS");
     float RADIUS   = psMetadataLookupF32 (&status, recipe, "PSF_MOMENTS_RADIUS");
+    char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
 
     sources = psArrayAlloc (peaks->n);
-    // DROP sources->n = 0;
 
     for (int i = 0; i < peaks->n; i++) {
@@ -26 +26 @@
         // allocate image, weight, mask arrays for each peak (square of radius OUTER)
         pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
+        if (!strcasecmp (breakPt, "PEAKS")) continue;
 
         // XXX skip faint sources?
@@ -72 +73 @@
     psLogMsg ("psphot", 3, "%d moments: %f sec\n", sources->n, psTimerMark ("psphot"));
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "STATS")) exit (0);
-
     return (sources);
 }