Index: /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_SV1.c =================================================================== --- /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_SV1.c (revision 28969) +++ /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_SV1.c (revision 28970) @@ -71,4 +71,5 @@ psF32 errMag, chisq, apRadius; psS32 nPix, nDOF; + char keyword1[80], keyword2[80]; pmChip *chip = readout->parent->parent; @@ -104,7 +105,18 @@ // we use this just to define the output vectors (which must be present for all objects) bool status = false; + psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER"); psVector *radMax = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.UPPER"); psAssert (radMax, "this must have been defined and tested earlier!"); psAssert (radMax->n, "this must have been defined and tested earlier!"); + psAssert (radMin->n == radMax->n, "inconsistent annular bins"); + + // write the radial profile apertures to header + for (int i = 0; i < radMax->n; i++) { + sprintf (keyword1, "RMIN_%02d", i); + sprintf (keyword2, "RMAX_%02d", i); + psMetadataAddF32 (imageHeader, PS_LIST_TAIL, keyword1, PS_META_REPLACE, "min radius for SB profile", radMin->data.F32[i]); + psMetadataAddF32 (imageHeader, PS_LIST_TAIL, keyword2, PS_META_REPLACE, "min radius for SB profile", radMax->data.F32[i]); + } + // we write out PSF-fits for all sources, regardless of quality. the source flags tell us the state @@ -202,17 +214,46 @@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_MINOR", PS_DATA_F32, "PSF width (minor axis)", axes.minor); psMetadataAdd (row, PS_LIST_TAIL, "PSF_THETA", PS_DATA_F32, "PSF orientation angle", axes.theta); - psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF", PS_DATA_F32, "PSF coverage/quality factor", source->pixWeightNotBad); + psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF", PS_DATA_F32, "PSF coverage/quality factor (bad)", source->pixWeightNotBad); + psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF_PERFECT", PS_DATA_F32, "PSF coverage/quality factor (poor)", source->pixWeightNotPoor); psMetadataAdd (row, PS_LIST_TAIL, "PSF_NDOF", PS_DATA_S32, "degrees of freedom", nDOF); psMetadataAdd (row, PS_LIST_TAIL, "PSF_NPIX", PS_DATA_S32, "number of pixels in fit", nPix); // distinguish moments measure from window vs S/N > XX ?? - float mxx = source->moments ? source->moments->Mxx : NAN; - float mxy = source->moments ? source->moments->Mxy : NAN; - float myy = source->moments ? source->moments->Myy : NAN; - psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XX", PS_DATA_F32, "second moments (X^2)", mxx); - psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XY", PS_DATA_F32, "second moments (X*Y)", mxy); - psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_YY", PS_DATA_F32, "second moments (Y*Y)", myy); - - psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); + float Mxx = source->moments ? source->moments->Mxx : NAN; + float Mxy = source->moments ? source->moments->Mxy : NAN; + float Myy = source->moments ? source->moments->Myy : NAN; + + float Mrf = source->moments ? source->moments->Mrf : NAN; + float Mrh = source->moments ? source->moments->Mrh : NAN; + float Krf = source->moments ? source->moments->KronFlux : NAN; + float dKrf = source->moments ? source->moments->KronFluxErr : NAN; + + float Kinner = source->moments ? source->moments->KronFinner : NAN; + float Kouter = source->moments ? source->moments->KronFouter : NAN; + + float M_c3 = source->moments ? 1.0*source->moments->Mxxx - 3.0*source->moments->Mxyy : NAN; + float M_s3 = source->moments ? 3.0*source->moments->Mxxy - 1.0*source->moments->Myyy : NAN; + float M_c4 = source->moments ? 1.0*source->moments->Mxxxx - 6.0*source->moments->Mxxyy + 1.0*source->moments->Myyyy : NAN; + float M_s4 = source->moments ? 4.0*source->moments->Mxxxy - 4.0*source->moments->Mxyyy : NAN; + + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XX", PS_DATA_F32, "second moments (X^2)", Mxx); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XY", PS_DATA_F32, "second moments (X*Y)", Mxy); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_YY", PS_DATA_F32, "second moments (Y*Y)", Myy); + + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3C", PS_DATA_F32, "third momemt cos theta", M_c3); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3S", PS_DATA_F32, "third momemt sin theta", M_s3); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4C", PS_DATA_F32, "fourth momemt cos theta", M_c4); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4S", PS_DATA_F32, "fourth momemt sin theta", M_s4); + + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_R1", PS_DATA_F32, "first radial moment", Mrf); + psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_RH", PS_DATA_F32, "half radial moment", Mrh); + psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX", PS_DATA_F32, "Kron Flux (in 2.5 R1)", Krf); + psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_ERR", PS_DATA_F32, "Kron Flux Error", dKrf); + + psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_INNER", PS_DATA_F32, "Kron Flux (in 2.5 R1)", Kinner); + psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER", PS_DATA_F32, "Kron Flux (in 2.5 R1)", Kouter); + + psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); + psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2", PS_DATA_U32, "psphot analysis flags", source->mode2); psVector *radFlux = psVectorAlloc(radMax->n, PS_TYPE_F32); @@ -372,4 +413,5 @@ source->pixWeightNotBad = psMetadataLookupF32 (&status, row, "PSF_QF"); + source->pixWeightNotPoor = psMetadataLookupF32 (&status, row, "PSF_QF_PERFECT"); source->crNsigma = psMetadataLookupF32 (&status, row, "CR_NSIGMA"); source->extNsigma = psMetadataLookupF32 (&status, row, "EXT_NSIGMA"); @@ -386,5 +428,27 @@ source->moments->Myy = psMetadataLookupF32 (&status, row, "MOMENTS_YY"); + source->moments->Mrf = psMetadataLookupF32 (&status, row, "MOMENTS_R1"); + source->moments->Mrh = psMetadataLookupF32 (&status, row, "MOMENTS_RH"); + source->moments->KronFlux = psMetadataLookupF32 (&status, row, "KRON_FLUX"); + source->moments->KronFluxErr = psMetadataLookupF32 (&status, row, "KRON_FLUX_ERR"); + + source->moments->KronFinner = psMetadataLookupF32 (&status, row, "KRON_FLUX_INNER"); + source->moments->KronFouter = psMetadataLookupF32 (&status, row, "KRON_FLUX_OUTER"); + + // XXX we do not save all of the 3rd and 4th moment parameters. when we load in data, + // we are storing enough information so the output will be consistent with the input + source->moments->Mxxx = +1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3C"); + source->moments->Mxxy = 0.0; + source->moments->Mxyy = 0.0; + source->moments->Myyy = -1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3S"); + + source->moments->Mxxxx = +1.00 * psMetadataLookupF32 (&status, row, "MOMENTS_M4C"); + source->moments->Mxxxy = 0.0; + source->moments->Mxxyy = 0.0; + source->moments->Mxyyy = -0.25 * psMetadataLookupF32 (&status, row, "MOMENTS_M4S"); + source->moments->Myyyy = 0.0; + source->mode = psMetadataLookupU32 (&status, row, "FLAGS"); + source->mode2 = psMetadataLookupU32 (&status, row, "FLAGS2"); assert (status); @@ -406,4 +470,5 @@ psF32 xErr, yErr; int nRow = -1; + char keyword1[80], keyword2[80]; // create a header to hold the output data @@ -437,6 +502,4 @@ bool doAnnuli = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI"); bool doPetrosian = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN"); - // bool doIsophotal = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ISOPHOTAL"); - // bool doKron = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_KRON"); psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER"); @@ -444,24 +507,14 @@ psAssert (radMin->n == radMax->n, "inconsistent annular bins"); - // int nRadialBins = 0; - // if (doAnnuli) { - // // get the max count of radial bins - // for (int i = 0; i < sources->n; i++) { - // pmSource *source = sources->data[i]; - // if (!source->extpars) continue; - // if (!source->extpars->radProfile ) continue; - // if (!source->extpars->radProfile->binSB) continue; - // nRadialBins = PS_MAX(nRadialBins, source->extpars->radProfile->binSB->n); - // } - // } + // write the radial profile apertures to header + for (int i = 0; i < radMax->n; i++) { + sprintf (keyword1, "RMIN_%02d", i); + sprintf (keyword2, "RMAX_%02d", i); + psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword1, PS_META_REPLACE, "min radius for SB profile", radMin->data.F32[i]); + psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword2, PS_META_REPLACE, "min radius for SB profile", radMax->data.F32[i]); + } // we write out all sources, regardless of quality. the source flags tell us the state for (int i = 0; i < sources->n; i++) { - // skip source if it is not a ext sourc - // XXX we have two places that extended source parameters are measured: - // psphotExtendedSources, which measures the aperture-like parameters and (potentially) the psf-convolved extended source models, - // psphotFitEXT, which does the simple extended source model fit (not psf-convolved) - // should we require both? - pmSource *source = sources->data[i]; @@ -529,39 +582,4 @@ } -# if (0) - // Kron measurements - if (doKron) { - pmSourceKronValues *kron = source->extpars->kron; - if (kron) { - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG", PS_DATA_F32, "Kron Magnitude", kron->mag); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG_ERR", PS_DATA_F32, "Kron Magnitude Error", kron->magErr); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS", PS_DATA_F32, "Kron Radius", kron->rad); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS_ERR", PS_DATA_F32, "Kron Radius Error", kron->radErr); - } else { - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG", PS_DATA_F32, "Kron Magnitude", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG_ERR", PS_DATA_F32, "Kron Magnitude Error", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS", PS_DATA_F32, "Kron Radius", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS_ERR", PS_DATA_F32, "Kron Radius Error", NAN); - } - } - - // Isophot measurements - // XXX insert header data: isophotal level - if (doIsophotal) { - pmSourceIsophotalValues *isophot = source->extpars->isophot; - if (isophot) { - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG", PS_DATA_F32, "Isophot Magnitude", isophot->mag); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG_ERR", PS_DATA_F32, "Isophot Magnitude Error", isophot->magErr); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS", PS_DATA_F32, "Isophot Radius", isophot->rad); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS_ERR", PS_DATA_F32, "Isophot Radius Error", isophot->radErr); - } else { - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG", PS_DATA_F32, "Isophot Magnitude", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG_ERR", PS_DATA_F32, "Isophot Magnitude Error", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS", PS_DATA_F32, "Isophot Radius", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS_ERR", PS_DATA_F32, "Isophot Radius Error", NAN); - } - } -# endif - // Flux Annuli (if we have extended source measurements, we have these. only optionally save them) if (doAnnuli) { Index: /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_V3.c =================================================================== --- /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_V3.c (revision 28969) +++ /branches/eam_branches/ipp-20100621/psModules/src/objects/pmSourceIO_CMF_PS1_V3.c (revision 28970) @@ -174,6 +174,6 @@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG", PS_DATA_F32, "PSF fit instrumental magnitude", source->psfMag); psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", errMag); - psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX", PS_DATA_F32, "PSF fit instrumental magnitude", source->psfFlux); - psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental magnitude", source->psfFluxErr); + psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX", PS_DATA_F32, "PSF fit instrumental flux (counts)", source->psfFlux); + psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental flux", source->psfFluxErr); psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG", PS_DATA_F32, "magnitude in standard aperture", source->apMag); psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RADIUS", PS_DATA_F32, "radius used for aperture mags", apRadius); @@ -233,5 +233,6 @@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER", PS_DATA_F32, "Kron Flux (in 2.5 R1)", Kouter); - psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); + psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); + psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2", PS_DATA_U32, "psphot analysis flags", source->mode2); // XXX not sure how to get this : need to load Nimages with weight? @@ -353,8 +354,8 @@ source->peak = pmPeakAlloc(PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], peakFlux, PM_PEAK_LONE); source->peak->flux = peakFlux; + source->peak->xf = PAR[PM_PAR_XPOS]; // more accurate position + source->peak->yf = PAR[PM_PAR_YPOS]; // more accurate position source->peak->dx = dPAR[PM_PAR_XPOS]; source->peak->dy = dPAR[PM_PAR_YPOS]; - source->peak->xf = PAR[PM_PAR_XPOS]; // more accurate position - source->peak->yf = PAR[PM_PAR_YPOS]; // more accurate position if (isfinite (source->errMag) && (source->errMag > 0.0)) { source->peak->SN = 1.0 / source->errMag; @@ -401,4 +402,5 @@ source->mode = psMetadataLookupU32 (&status, row, "FLAGS"); + source->mode2 = psMetadataLookupU32 (&status, row, "FLAGS2"); assert (status); @@ -412,5 +414,4 @@ bool pmSourcesWrite_CMF_PS1_V3_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe) { - bool status; psArray *table; @@ -420,4 +421,5 @@ psF32 xErr, yErr; int nRow = -1; + char keyword1[80], keyword2[80]; // create a header to hold the output data @@ -451,6 +453,4 @@ bool doAnnuli = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI"); bool doPetrosian = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN"); - // bool doIsophotal = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ISOPHOTAL"); - // bool doKron = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_KRON"); psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER"); @@ -458,24 +458,14 @@ psAssert (radMin->n == radMax->n, "inconsistent annular bins"); - // int nRadialBins = 0; - // if (doAnnuli) { - // // get the max count of radial bins - // for (int i = 0; i < sources->n; i++) { - // pmSource *source = sources->data[i]; - // if (!source->extpars) continue; - // if (!source->extpars->radProfile ) continue; - // if (!source->extpars->radProfile->binSB) continue; - // nRadialBins = PS_MAX(nRadialBins, source->extpars->radProfile->binSB->n); - // } - // } + // write the radial profile apertures to header + for (int i = 0; i < radMax->n; i++) { + sprintf (keyword1, "RMIN_%02d", i); + sprintf (keyword2, "RMAX_%02d", i); + psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword1, PS_META_REPLACE, "min radius for SB profile", radMin->data.F32[i]); + psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword2, PS_META_REPLACE, "min radius for SB profile", radMax->data.F32[i]); + } // we write out all sources, regardless of quality. the source flags tell us the state for (int i = 0; i < sources->n; i++) { - // skip source if it is not a ext sourc - // XXX we have two places that extended source parameters are measured: - // psphotExtendedSources, which measures the aperture-like parameters and (potentially) the psf-convolved extended source models, - // psphotFitEXT, which does the simple extended source model fit (not psf-convolved) - // should we require both? - pmSource *source = sources->data[i]; @@ -543,39 +533,4 @@ } -# if (0) - // Kron measurements - if (doKron) { - pmSourceKronValues *kron = source->extpars->kron; - if (kron) { - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG", PS_DATA_F32, "Kron Magnitude", kron->mag); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG_ERR", PS_DATA_F32, "Kron Magnitude Error", kron->magErr); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS", PS_DATA_F32, "Kron Radius", kron->rad); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS_ERR", PS_DATA_F32, "Kron Radius Error", kron->radErr); - } else { - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG", PS_DATA_F32, "Kron Magnitude", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_MAG_ERR", PS_DATA_F32, "Kron Magnitude Error", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS", PS_DATA_F32, "Kron Radius", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "KRON_RADIUS_ERR", PS_DATA_F32, "Kron Radius Error", NAN); - } - } - - // Isophot measurements - // XXX insert header data: isophotal level - if (doIsophotal) { - pmSourceIsophotalValues *isophot = source->extpars->isophot; - if (isophot) { - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG", PS_DATA_F32, "Isophot Magnitude", isophot->mag); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG_ERR", PS_DATA_F32, "Isophot Magnitude Error", isophot->magErr); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS", PS_DATA_F32, "Isophot Radius", isophot->rad); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS_ERR", PS_DATA_F32, "Isophot Radius Error", isophot->radErr); - } else { - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG", PS_DATA_F32, "Isophot Magnitude", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_MAG_ERR", PS_DATA_F32, "Isophot Magnitude Error", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS", PS_DATA_F32, "Isophot Radius", NAN); - psMetadataAdd (row, PS_LIST_TAIL, "ISOPHOT_RADIUS_ERR", PS_DATA_F32, "Isophot Radius Error", NAN); - } - } -# endif - // Flux Annuli (if we have extended source measurements, we have these. only optionally save them) if (doAnnuli) {