Index: /branches/eam_branch_20080324/psModules/src/objects/pmSourceIO_PS1_DEV_1.c =================================================================== --- /branches/eam_branch_20080324/psModules/src/objects/pmSourceIO_PS1_DEV_1.c (revision 17308) +++ /branches/eam_branch_20080324/psModules/src/objects/pmSourceIO_PS1_DEV_1.c (revision 17309) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ - * @date $Date: 2008-03-05 01:08:08 $ + * @version $Revision: 1.9.2.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2008-04-03 06:01:38 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -256,7 +256,5 @@ psArray *table; psMetadata *row; - int i; psF32 *PAR, *dPAR; - psEllipseAxes axes; psF32 xPos, yPos; psF32 xErr, yErr; @@ -274,5 +272,5 @@ // we write out all sources, regardless of quality. the source flags tell us the state - for (i = 0; i < sources->n; i++) { + 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: @@ -283,13 +281,10 @@ pmSource *source = sources->data[i]; - // choose the convolved EXT model, if available, otherwise the simple one - pmModel *model = source->modelConv; - if (model == NULL) { - model = source->modelEXT; - } - if (model == NULL) continue; - - // XXX do we need to do this? + // skip sources without measurements if (source->extpars == NULL) continue; + + // we require a PSF model fit (ignore the real crud) + pmModel *model = source->modelPSF; + if (model == NULL) continue; // XXX I need to split the extended models from the extended aperture measurements @@ -301,8 +296,6 @@ yErr = dPAR[PM_PAR_YPOS]; - // XXX for the aperture values, I probably can / should just refer to the psf position and not write any of this - axes = pmPSF_ModelToAxes (PAR, 20.0); - row = psMetadataAlloc (); + // XXX we are not writing out the mode (flags) or the type (psf, ext, etc) psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET", PS_DATA_U32, "IPP detection identifier index", source->seq); @@ -311,61 +304,47 @@ psMetadataAdd (row, PS_LIST_TAIL, "X_EXT_SIG", PS_DATA_F32, "Sigma in EXT x coordinate", xErr); psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT_SIG", PS_DATA_F32, "Sigma in EXT y coordinate", yErr); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG", PS_DATA_F32, "EXT fit instrumental magnitude", source->extMag); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", source->errMag); - - // XXX these should be major and minor, not 'x' and 'y' - psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", PS_DATA_F32, "EXT width in x coordinate", axes.major); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", PS_DATA_F32, "EXT width in y coordinate", axes.minor); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_THETA", PS_DATA_F32, "EXT orientation angle", axes.theta); - - // XXX at the moment, this will be a fixed sized table, but perhaps have multiple output formats depending on what is measured? - - // other values that I need to report: - // R, Mag Petrosian, ..... - if (source->extpars) { - - // Petrosian measurements - pmSourcePetrosianValues *petrosian = source->extpars->petrosian; - if (petrosian) { - psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG", PS_DATA_F32, "Petrosian Magnitude", petrosian->mag); - psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG_ERR", PS_DATA_F32, "Petrosian Magnitude Error", petrosian->magErr); - psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS", PS_DATA_F32, "Petrosian Radius", petrosian->rad); - psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_ERR", PS_DATA_F32, "Petrosian Radius Error", petrosian->radErr); - } - - // Kron measurements - 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); - } - - // Isophot measurements - // XXX insert header data: isophotal level - 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); - } - - pmSourceAnnuli *annuli = source->extpars->annuli; - if (annuli) { - psVector *fluxVal = annuli->flux; - psVector *fluxErr = annuli->fluxErr; - psVector *fluxVar = annuli->fluxVar; - - for (int j = 0; j < fluxVal->n; j++) { - char name[32]; - sprintf (name, "FLUX_VAL_R_%02d", j); - psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux value in annulus", fluxVal->data.F32[j]); - sprintf (name, "FLUX_ERR_R_%02d", j); - psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux error in annulus", fluxErr->data.F32[j]); - sprintf (name, "FLUX_VAR_R_%02d", j); - psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux stdev in annulus", fluxVar->data.F32[j]); - } + + // Petrosian measurements + pmSourcePetrosianValues *petrosian = source->extpars->petrosian; + if (petrosian) { + psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG", PS_DATA_F32, "Petrosian Magnitude", petrosian->mag); + psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG_ERR", PS_DATA_F32, "Petrosian Magnitude Error", petrosian->magErr); + psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS", PS_DATA_F32, "Petrosian Radius", petrosian->rad); + psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_ERR", PS_DATA_F32, "Petrosian Radius Error", petrosian->radErr); + } + + // Kron measurements + 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); + } + + // Isophot measurements + // XXX insert header data: isophotal level + 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); + } + + pmSourceAnnuli *annuli = source->extpars->annuli; + if (annuli) { + psVector *fluxVal = annuli->flux; + psVector *fluxErr = annuli->fluxErr; + psVector *fluxVar = annuli->fluxVar; + + for (int j = 0; j < fluxVal->n; j++) { + char name[32]; + sprintf (name, "FLUX_VAL_R_%02d", j); + psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux value in annulus", fluxVal->data.F32[j]); + sprintf (name, "FLUX_ERR_R_%02d", j); + psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux error in annulus", fluxErr->data.F32[j]); + sprintf (name, "FLUX_VAR_R_%02d", j); + psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux stdev in annulus", fluxVar->data.F32[j]); } } @@ -397,9 +376,9 @@ psArray *table; psMetadata *row; - int i; psF32 *PAR, *dPAR; psEllipseAxes axes; psF32 xPos, yPos; psF32 xErr, yErr; + char name[64]; // create a header to hold the output data @@ -415,53 +394,63 @@ // we write out all sources, regardless of quality. the source flags tell us the state - for (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? + for (int i = 0; i < sources->n; i++) { pmSource *source = sources->data[i]; - // XXX need to have an array of model fits; loop over the array and write out all - // which we calculated - - // choose the convolved EXT model, if available, otherwise the simple one - // XXX should not need to choose: write both out - pmModel *model = source->modelConv; - if (model == NULL) { - model = source->modelEXT; + // XXX if no model fits are saved, write out modelEXT? + if (source->modelFits == NULL) continue; + + // We have multiple sources : need to flag the one used to subtract the light (the 'best' model) + for (int j = 0; j < source->modelFits->n; j++) { + + // choose the convolved EXT model, if available, otherwise the simple one + pmModel *model = source->modelFits->data[j]; + assert (model); + + PAR = model->params->data.F32; + dPAR = model->dparams->data.F32; + xPos = PAR[PM_PAR_XPOS]; + yPos = PAR[PM_PAR_YPOS]; + xErr = dPAR[PM_PAR_XPOS]; + yErr = dPAR[PM_PAR_YPOS]; + + axes = pmPSF_ModelToAxes (PAR, 20.0); + + row = psMetadataAlloc (); + + // XXX we are not writing out the mode (flags) or the type (psf, ext, etc) + psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET", PS_DATA_U32, "IPP detection identifier index", source->seq); + psMetadataAdd (row, PS_LIST_TAIL, "X_EXT", PS_DATA_F32, "EXT model x coordinate", xPos); + psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT", PS_DATA_F32, "EXT model y coordinate", yPos); + psMetadataAdd (row, PS_LIST_TAIL, "X_EXT_SIG", PS_DATA_F32, "Sigma in EXT x coordinate", xErr); + psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT_SIG", PS_DATA_F32, "Sigma in EXT y coordinate", yErr); + psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG", PS_DATA_F32, "EXT fit instrumental magnitude", model->mag); + psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", model->magErr); + + // XXX these should be major and minor, not 'x' and 'y' + psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", PS_DATA_F32, "EXT width in x coordinate", axes.major); + psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", PS_DATA_F32, "EXT width in y coordinate", axes.minor); + psMetadataAdd (row, PS_LIST_TAIL, "EXT_THETA", PS_DATA_F32, "EXT orientation angle", axes.theta); + + // write out the other generic parameters + for (int k = 0; k < model->params->n; k++) { + if (k == PM_PAR_I0) continue; + if (k == PM_PAR_XPOS) continue; + if (k == PM_PAR_YPOS) continue; + if (k == PM_PAR_SXX) continue; + if (k == PM_PAR_SXY) continue; + if (k == PM_PAR_SYY) continue; + + snprintf (name, 64, "EXT_PAR_%02d", k); + psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "", model->params->data.F32[k]); + } + + // XXX other parameters which may be set. + // XXX flag / value to define the model + // XXX write out the model type, fit status flags + + psArrayAdd (table, 100, row); + psFree (row); } - if (model == NULL) continue; - - PAR = model->params->data.F32; - dPAR = model->dparams->data.F32; - xPos = PAR[PM_PAR_XPOS]; - yPos = PAR[PM_PAR_YPOS]; - xErr = dPAR[PM_PAR_XPOS]; - yErr = dPAR[PM_PAR_YPOS]; - - axes = pmPSF_ModelToAxes (PAR, 20.0); - - row = psMetadataAlloc (); - // XXX we are not writing out the mode (flags) or the type (psf, ext, etc) - psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET", PS_DATA_U32, "IPP detection identifier index", source->seq); - psMetadataAdd (row, PS_LIST_TAIL, "X_EXT", PS_DATA_F32, "EXT model x coordinate", xPos); - psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT", PS_DATA_F32, "EXT model y coordinate", yPos); - psMetadataAdd (row, PS_LIST_TAIL, "X_EXT_SIG", PS_DATA_F32, "Sigma in EXT x coordinate", xErr); - psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT_SIG", PS_DATA_F32, "Sigma in EXT y coordinate", yErr); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG", PS_DATA_F32, "EXT fit instrumental magnitude", source->extMag); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", source->errMag); - - // XXX these should be major and minor, not 'x' and 'y' - psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", PS_DATA_F32, "EXT width in x coordinate", axes.major); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", PS_DATA_F32, "EXT width in y coordinate", axes.minor); - psMetadataAdd (row, PS_LIST_TAIL, "EXT_THETA", PS_DATA_F32, "EXT orientation angle", axes.theta); - - // XXX other parameters which may be set. - // XXX flag / value to define the model - - psArrayAdd (table, 100, row); - psFree (row); }