Index: trunk/psModules/src/objects/Makefile.am =================================================================== --- trunk/psModules/src/objects/Makefile.am (revision 20010) +++ trunk/psModules/src/objects/Makefile.am (revision 20076) @@ -44,4 +44,5 @@ pmPSFtry.c \ pmTrend2D.c \ + pmGrowthCurveGenerate.c \ pmGrowthCurve.c Index: trunk/psModules/src/objects/pmGrowthCurveGenerate.c =================================================================== --- trunk/psModules/src/objects/pmGrowthCurveGenerate.c (revision 20076) +++ trunk/psModules/src/objects/pmGrowthCurveGenerate.c (revision 20076) @@ -0,0 +1,157 @@ +/** @file pmGrowthCurveGenerate.c + * + * Generate the curve-of-growth + * + * @author EAM, IfA + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2008-10-13 01:56:42 $ + * + * Copyright 2004 Institute for Astronomy, University of Hawaii + * + */ + +#ifdef HAVE_CONFIG_H +#include +#endif + +/*****************************************************************************/ +/* INCLUDE FILES */ +/*****************************************************************************/ + +#include // for strcasecmp +#include +#include "pmHDU.h" +#include "pmFPA.h" +#include "pmPeaks.h" +#include "pmMoments.h" +#include "pmResiduals.h" +#include "pmGrowthCurve.h" +#include "pmTrend2D.h" +#include "pmPSF.h" +#include "pmModel.h" +#include "pmSource.h" +#include "pmModelClass.h" +#include "pmModelUtils.h" +#include "pmSourcePhotometry.h" +#include "pmFPAMaskWeight.h" +#include "psVectorBracket.h" +#include "pmErrorCodes.h" + +/*****************************************************************************/ +/* FUNCTION IMPLEMENTATION - PUBLIC */ +/*****************************************************************************/ + +// we generate the growth curve for the center of the image with the specified psf model +bool pmGrowthCurveGenerate (pmReadout *readout, pmPSF *psf, bool ignore, psMaskType maskVal, psMaskType markVal) +{ + PS_ASSERT_PTR_NON_NULL(readout, false); + PS_ASSERT_PTR_NON_NULL(readout->image, false); + + // bool status; + float xc, yc, dx, dy; + float fitMag, apMag; + float radius; + + // create template model + pmModel *modelRef = pmModelAlloc(psf->type); + + // maskVal is used to test for rejected pixels, and must include markVal + maskVal |= markVal; + + // XXXX A Serious hack: the psf might not be determined in the field center. + // for the moment, just try a few locations until it is defined! + + pmModel *model = NULL; + for (int ix = -1; ix <= +1; ix ++) { + for (int iy = -1; iy <= +1; iy ++) { + + // use the center of the center pixel of the image + xc = (0.5 + 0.3*ix)*readout->image->numCols + readout->image->col0 + 0.5; + yc = (0.5 + 0.3*ix)*readout->image->numRows + readout->image->row0 + 0.5; + dx = psf->growth->maxRadius + 1; + dy = psf->growth->maxRadius + 1; + + // assign the x and y coords to the image center + // create an object with center intensity of 1000 + modelRef->params->data.F32[PM_PAR_SKY] = 0; + modelRef->params->data.F32[PM_PAR_I0] = 1000; + modelRef->params->data.F32[PM_PAR_XPOS] = xc; + modelRef->params->data.F32[PM_PAR_YPOS] = yc; + + // create modelPSF from this model + model = pmModelFromPSF (modelRef, psf); + if (model) goto got_model; + } + } + psAssert (model, "cannot build growth curve (psf model is invalid everywhere)"); + +got_model: + // measure the fitMag for this model + pmSourcePhotometryModel (&fitMag, model); + psf->growth->fitMag = fitMag; + + // generate working image for this source + psRegion region = {xc - dx, xc + dx, yc - dy, yc + dy}; + + // force region to stop at dimensions of image + region = psRegionForImage (readout->image, region); + + // the view, image, and mask retain col0,row0 + psImage *view = psImageSubset (readout->image, region); + psImage *image = psImageCopy (NULL, view, PS_TYPE_F32); + psImage *mask = psImageCopy (NULL, view, PS_TYPE_U8); + + psImageInit (image, 0.0); + psImageInit (mask, 0); + + // place the reference object in the image center + // no need to mask the source here + // XXX should we measure this for the analytical model only or the full model? + pmModelAdd (image, NULL, model, PM_MODEL_OP_FULL, maskVal); + + // loop over a range of source fluxes + // no need to interpolate since we have forced the object center + // to 0.5, 0.5 above + bool completeGrowthCurve = true; // do we have a complete curve of growth? + for (int i = 0; i < psf->growth->radius->n; i++) { + + radius = psf->growth->radius->data.F32[i]; + + // mask the given aperture and measure the apMag + psImageKeepCircle (mask, xc, yc, radius, "OR", markVal); + if (!pmSourcePhotometryAper (&apMag, model, image, mask, maskVal)) { + psError(PM_ERR_PHOTOM, false, "Measuring apMag for radius == %g", radius); + psErrorStackPrint (stderr, "failure to get growth curve"); + psErrorClear(); + completeGrowthCurve = false; + break; + } + psImageKeepCircle (mask, xc, yc, radius, "AND", PS_NOT_U8(markVal)); + + // the 'ignore' mode is for testing + if (ignore) { + psf->growth->apMag->data.F32[i] = fitMag; + } else { + psf->growth->apMag->data.F32[i] = apMag; + } + } + + if (completeGrowthCurve) { + psf->growth->apRef = psVectorInterpolate (psf->growth->radius, psf->growth->apMag, psf->growth->refRadius); + psf->growth->apLoss = psf->growth->fitMag - psf->growth->apRef; + } else { + psf->growth->apRef = NAN; + psf->growth->apLoss = 0; + } + + psLogMsg ("psphot.growth", 4, "GrowthCurve : apLoss : %f\n", psf->growth->apLoss); + + psFree (view); + psFree (image); + psFree (mask); + psFree (model); + psFree (modelRef); + + return completeGrowthCurve ? true : false; +} Index: trunk/psModules/src/objects/pmPSF.c =================================================================== --- trunk/psModules/src/objects/pmPSF.c (revision 20010) +++ trunk/psModules/src/objects/pmPSF.c (revision 20076) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.37 $ $Name: not supported by cvs2svn $ - * @date $Date: 2008-10-08 21:53:29 $ + * @version $Revision: 1.38 $ $Name: not supported by cvs2svn $ + * @date $Date: 2008-10-13 01:56:42 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -401,116 +401,2 @@ return psf; } - -// we generate the growth curve for the center of the image with the specified psf model -bool pmGrowthCurveGenerate (pmReadout *readout, pmPSF *psf, bool ignore, psMaskType maskVal, psMaskType markVal) -{ - PS_ASSERT_PTR_NON_NULL(readout, false); - PS_ASSERT_PTR_NON_NULL(readout->image, false); - - // bool status; - float xc, yc, dx, dy; - float fitMag, apMag; - float radius; - - // create template model - pmModel *modelRef = pmModelAlloc(psf->type); - - // maskVal is used to test for rejected pixels, and must include markVal - maskVal |= markVal; - - // XXXX A Serious hack: the psf might not be determined in the field center. - // for the moment, just try a few locations until it is defined! - - pmModel *model = NULL; - for (int ix = -1; ix <= +1; ix ++) { - for (int iy = -1; iy <= +1; iy ++) { - - // use the center of the center pixel of the image - xc = (0.5 + 0.3*ix)*readout->image->numCols + readout->image->col0 + 0.5; - yc = (0.5 + 0.3*ix)*readout->image->numRows + readout->image->row0 + 0.5; - dx = psf->growth->maxRadius + 1; - dy = psf->growth->maxRadius + 1; - - // assign the x and y coords to the image center - // create an object with center intensity of 1000 - modelRef->params->data.F32[PM_PAR_SKY] = 0; - modelRef->params->data.F32[PM_PAR_I0] = 1000; - modelRef->params->data.F32[PM_PAR_XPOS] = xc; - modelRef->params->data.F32[PM_PAR_YPOS] = yc; - - // create modelPSF from this model - model = pmModelFromPSF (modelRef, psf); - if (model) goto got_model; - } - } - psAssert (model, "cannot build growth curve (psf model is invalid everywhere)"); - -got_model: - // measure the fitMag for this model - pmSourcePhotometryModel (&fitMag, model); - psf->growth->fitMag = fitMag; - - // generate working image for this source - psRegion region = {xc - dx, xc + dx, yc - dy, yc + dy}; - - // force region to stop at dimensions of image - region = psRegionForImage (readout->image, region); - - // the view, image, and mask retain col0,row0 - psImage *view = psImageSubset (readout->image, region); - psImage *image = psImageCopy (NULL, view, PS_TYPE_F32); - psImage *mask = psImageCopy (NULL, view, PS_TYPE_U8); - - psImageInit (image, 0.0); - psImageInit (mask, 0); - - // place the reference object in the image center - // no need to mask the source here - // XXX should we measure this for the analytical model only or the full model? - pmModelAdd (image, NULL, model, PM_MODEL_OP_FULL, maskVal); - - // loop over a range of source fluxes - // no need to interpolate since we have forced the object center - // to 0.5, 0.5 above - bool completeGrowthCurve = true; // do we have a complete curve of growth? - for (int i = 0; i < psf->growth->radius->n; i++) { - - radius = psf->growth->radius->data.F32[i]; - - // mask the given aperture and measure the apMag - psImageKeepCircle (mask, xc, yc, radius, "OR", markVal); - if (!pmSourcePhotometryAper (&apMag, model, image, mask, maskVal)) { - psError(PM_ERR_PHOTOM, false, "Measuring apMag for radius == %g", radius); - psErrorStackPrint (stderr, "failure to get growth curve"); - psErrorClear(); - completeGrowthCurve = false; - break; - } - psImageKeepCircle (mask, xc, yc, radius, "AND", PS_NOT_U8(markVal)); - - // the 'ignore' mode is for testing - if (ignore) { - psf->growth->apMag->data.F32[i] = fitMag; - } else { - psf->growth->apMag->data.F32[i] = apMag; - } - } - - if (completeGrowthCurve) { - psf->growth->apRef = psVectorInterpolate (psf->growth->radius, psf->growth->apMag, psf->growth->refRadius); - psf->growth->apLoss = psf->growth->fitMag - psf->growth->apRef; - } else { - psf->growth->apRef = NAN; - psf->growth->apLoss = 0; - } - - psLogMsg ("psphot.growth", 4, "GrowthCurve : apLoss : %f\n", psf->growth->apLoss); - - psFree (view); - psFree (image); - psFree (mask); - psFree (model); - psFree (modelRef); - - return completeGrowthCurve ? true : false; -}