Index: /trunk/psphot/doc/notes.txt
===================================================================
--- /trunk/psphot/doc/notes.txt	(revision 7757)
+++ /trunk/psphot/doc/notes.txt	(revision 7758)
@@ -1,2 +1,15 @@
+
+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
Index: /trunk/psphot/src/psphot.h
===================================================================
--- /trunk/psphot/src/psphot.h	(revision 7757)
+++ /trunk/psphot/src/psphot.h	(revision 7758)
@@ -21,5 +21,4 @@
 
 psArray        *psphotFakeSources ();
-int             psphotSaveImage (psMetadata *header, psImage *image, char *filename);
 
 // psphotReadout functions
@@ -31,4 +30,5 @@
 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,4 +37,5 @@
 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,4 +59,5 @@
 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
Index: /trunk/psphot/src/psphotBasicDeblend.c
===================================================================
--- /trunk/psphot/src/psphotBasicDeblend.c	(revision 7757)
+++ /trunk/psphot/src/psphotBasicDeblend.c	(revision 7758)
@@ -125,7 +125,4 @@
     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);
Index: /trunk/psphot/src/psphotChoosePSF.c
===================================================================
--- /trunk/psphot/src/psphotChoosePSF.c	(revision 7757)
+++ /trunk/psphot/src/psphotChoosePSF.c	(revision 7758)
@@ -29,5 +29,4 @@
 
     stars = psArrayAlloc (sources->n);
-    // DROP stars->n = 0;
 
     // select the candidate PSF stars (pointers to original sources)
@@ -105,7 +104,4 @@
     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,2 +147,44 @@
     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;
+}
Index: /trunk/psphot/src/psphotFindPeaks.c
===================================================================
--- /trunk/psphot/src/psphotFindPeaks.c	(revision 7757)
+++ /trunk/psphot/src/psphotFindPeaks.c	(revision 7758)
@@ -28,4 +28,5 @@
     psImage *mask = readout->mask;
 
+    // XXX make this a user-option?
     // psphotSaveImage (NULL, smooth_im, "imsmooth.fits");
     // psphotSaveImage (NULL, smooth_wt, "wtsmooth.fits");
@@ -64,7 +65,4 @@
     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);
 }
Index: /trunk/psphot/src/psphotReadout.c
===================================================================
--- /trunk/psphot/src/psphotReadout.c	(revision 7757)
+++ /trunk/psphot/src/psphotReadout.c	(revision 7758)
@@ -8,4 +8,7 @@
     // 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,9 +37,14 @@
     // 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,4 +52,8 @@
         pmFPAfileDropInternal (config->files, "PSPHOT.BACKGND");
         return false;
+    }
+    if (!strcasecmp (breakPt, "MOMENTS")) {
+	psphotReadoutCleanup (config, readout, recipe, NULL, sources);
+	return true;
     }
 
@@ -51,8 +63,15 @@
     // 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,7 +99,4 @@
     psphotMagnitudes (sources, recipe, psf);
 
-    // create an output header with stats results
-    psMetadata *header = psphotDefineHeader (recipe);
-
     // replace background in residual image
     psphotSkyReplace (config, view);
@@ -89,4 +105,21 @@
     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,9 +132,8 @@
     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;
 }
Index: /trunk/psphot/src/psphotRoughClass.c
===================================================================
--- /trunk/psphot/src/psphotRoughClass.c	(revision 7757)
+++ /trunk/psphot/src/psphotRoughClass.c	(revision 7758)
@@ -4,5 +4,4 @@
 bool psphotRoughClass (psArray *sources, psMetadata *recipe) {
 
-    bool status;
     pmPSFClump   psfClump;
 
@@ -20,7 +19,4 @@
     psphotDumpMoments (recipe, sources);
 
-    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
-    if (!strcasecmp (breakPt, "CLASS")) exit (0);
-
     return true;
 }
Index: /trunk/psphot/src/psphotSourceStats.c
===================================================================
--- /trunk/psphot/src/psphotSourceStats.c	(revision 7757)
+++ /trunk/psphot/src/psphotSourceStats.c	(revision 7758)
@@ -14,7 +14,7 @@
     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,4 +26,5 @@
 	// 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,7 +73,4 @@
     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);
 }
