Changeset 6862 for trunk/psphot/src/psphotChoosePSF.c

Timestamp:

Apr 14, 2006, 11:58:44 AM (20 years ago)

Author:

eugene

Message:

renamed config to recipe in lower-level functions, added pixel weights to pmSource, other cleanups

File:

• trunk/psphot/src/psphotChoosePSF.c (modified) (6 diffs)

trunk/psphot/src/psphotChoosePSF.c

@@ -1 +1 @@
 # include "psphot.h"
 
-// 2006.02.07 : no leaks!
 // try PSF models and select best option
-pmPSF *psphotChoosePSF (psArray *sources, psMetadata *config) { 
+pmPSF *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe) { 
 
     bool            status;
     char           *modelName;
-    pmPSFtry       *try   = NULL;
+    pmPSF          *psf = NULL;
+    pmPSFtry       *try = NULL;
     psArray        *stars = NULL;
 
@@ -13 +13 @@
 
     // check if a PSF model is supplied by the user
-    char *psfFile = psMetadataLookupStr (&status, config, "PSF_INPUT_FILE");
-    if (status) {
-        int Nfail;
-        psMetadata *psfData = psMetadataConfigParse (NULL, &Nfail, psfFile, FALSE);
-        pmPSF *psf = pmPSFfromMD (psfData);
-        return psf;
-    }
+    psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
+    if (psf != NULL) return psf;
 
     // examine PSF sources in S/N order (brightest first)
@@ -25 +20 @@
 
     // array to store candidate PSF stars
-    int NSTARS = psMetadataLookupS32 (&status, config, "PSF_MAX_NSTARS");
+    int NSTARS = psMetadataLookupS32 (&status, recipe, "PSF_MAX_NSTARS");
     if (!status) NSTARS = PS_MIN (sources->n, 200);
 
     // use poissonian errors or local-sky errors
-    bool POISSON_ERRORS = psMetadataLookupBool (&status, config, "POISSON_ERRORS");
+    bool POISSON_ERRORS = psMetadataLookupBool (&status, recipe, "POISSON_ERRORS");
     if (!status) POISSON_ERRORS = true;
     pmSourceFitModelInit (15, 0.1, POISSON_ERRORS);
@@ -45 +40 @@
     // get the fixed PSF fit radius
     // XXX EAM : check that PSF_FIT_RADIUS < SKY_OUTER_RADIUS
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, recipe, "PSF_FIT_RADIUS");
 
     // get the list pointers for the PSF_MODEL entries
     psList *list = NULL;
-    psMetadataItem *mdi = psMetadataLookup (config, "PSF_MODEL");
+    psMetadataItem *mdi = psMetadataLookup (recipe, "PSF_MODEL");
     if (mdi == NULL) psAbort ("psphotChoosePSF", "missing PSF_MODEL selection");
 
@@ -100 +95 @@
 
     // save only the best model;
-    pmPSF *psf = psMemCopy(try->psf);
+    psf = psMemCopy(try->psf);
     psFree (models);
 
-    psMetadataAdd (config, PS_LIST_TAIL, "NPSFSTAR", PS_DATA_S32 | PS_META_REPLACE, "number of stars used to make PSF", psf->nPSFstars);
+    psMetadataAdd (recipe, PS_LIST_TAIL, "NPSFSTAR", PS_DATA_S32 | PS_META_REPLACE, "number of stars used to make PSF", psf->nPSFstars);
 
     modelName = pmModelGetType (psf->type);
@@ -109 +104 @@
     psLogMsg ("psphot.pspsf", 3, "selected psf model %s, ApResid: %f +/- %f\n", modelName, psf->ApResid, psf->dApResid);
 
-    char *breakPt = psMetadataLookupStr (&status, config, "BREAK_POINT");
+    char *breakPt = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
     if (!strcasecmp (breakPt, "PSFFIT")) exit (0);
 
     return (psf);
 }
+
+bool psphotPSFstats (pmReadout *readout, psMetadata *recipe, pmPSF *psf) {
+
+    psEllipseShape shape;
+    psEllipseAxes axes;
+
+    psImage *image = readout->image;
+
+    pmModel *modelEXT = pmModelAlloc (psf->type);
+
+    // make a model with unit central intensity at the image center
+    modelEXT->params->data.F32[0] = 0;
+    modelEXT->params->data.F32[1] = 1;
+    modelEXT->params->data.F32[2] = 0.5*image->numCols;
+    modelEXT->params->data.F32[3] = 0.5*image->numRows;
+
+    // construct a PSF model at this coordinate
+    pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
+
+    // get the correct model-radius function
+    pmModelRadius modelRadiusFunc = pmModelRadius_GetFunction (psf->type);
+
+    // get the model full-width at half-max
+    psF64 FWHM_X = 2*modelRadiusFunc (modelPSF->params, 0.5);
+    
+    shape.sx  = modelPSF->params->data.F32[4];
+    shape.sy  = modelPSF->params->data.F32[5];
+    shape.sxy = modelPSF->params->data.F32[6];
+    axes = psEllipseShapeToAxes (shape);
+
+    psF64 FWHM_Y = FWHM_X * (axes.minor / axes.major);
+
+    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",           axes.theta);
+    
+    return true;
+}