Index: trunk/psModules/src/objects/models/pmModel_GAUSS.c =================================================================== --- trunk/psModules/src/objects/models/pmModel_GAUSS.c (revision 13038) +++ trunk/psModules/src/objects/models/pmModel_GAUSS.c (revision 13064) @@ -208,6 +208,6 @@ psF32 *PAR = params->data.F32; - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -237,6 +237,6 @@ return (1.0); - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -267,8 +267,14 @@ } - // the 2D model for SXY actually fits SXY / (SXX^-2 + SYY^-2); correct here - out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); + // the OLD 2D model for SXY actually fitted SXY / (SXX^-2 + SYY^-2); correct here + // out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); + + // the 2D PSF model fits polarization terms (E0,E1,E2) + // convert to shape terms (SXX,SYY,SXY) + // XXX user-defined value for limit? + pmPSF_FitToModel (out, 0.1); // apply the model limits here: this truncates excessive extrapolation + // XXX do we need to do this still? should we put in asserts to test? for (int i = 0; i < psf->params_NEW->n; i++) { // apply the limits to all components or just the psf-model parameters? Index: trunk/psModules/src/objects/models/pmModel_PGAUSS.c =================================================================== --- trunk/psModules/src/objects/models/pmModel_PGAUSS.c (revision 13038) +++ trunk/psModules/src/objects/models/pmModel_PGAUSS.c (revision 13064) @@ -209,6 +209,6 @@ psF32 *PAR = params->data.F32; - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -253,6 +253,6 @@ return (1.0); - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -288,7 +288,20 @@ } - // the 2D model for SXY actually fits SXY / (SXX^-2 + SYY^-2); correct here - out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); - + // the OLD 2D model for SXY actually fitted SXY / (SXX^-2 + SYY^-2); correct here + // out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); + + // the 2D PSF model fits polarization terms (E0,E1,E2) + // convert to shape terms (SXX,SYY,SXY) + pmPSF_FitToModel (out, 0.1); + + // apply the model limits here: this truncates excessive extrapolation + // XXX do we need to do this still? should we put in asserts to test? + for (int i = 0; i < psf->params_NEW->n; i++) { + // apply the limits to all components or just the psf-model parameters? + if (psf->params_NEW->data[i] == NULL) + continue; + PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MIN, i, out, NULL); + PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MAX, i, out, NULL); + } return(true); } Index: trunk/psModules/src/objects/models/pmModel_QGAUSS.c =================================================================== --- trunk/psModules/src/objects/models/pmModel_QGAUSS.c (revision 13038) +++ trunk/psModules/src/objects/models/pmModel_QGAUSS.c (revision 13064) @@ -238,6 +238,6 @@ psF32 *PAR = params->data.F32; - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -284,6 +284,6 @@ return (1.0); - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); + shape.sx = PAR[PM_PAR_SXX] / M_SQRT2; + shape.sy = PAR[PM_PAR_SYY] / M_SQRT2; shape.sxy = PAR[PM_PAR_SXY]; @@ -345,8 +345,13 @@ } - // the 2D model for SXY actually fits SXY / (SXX^-2 + SYY^-2); correct here - out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); + // the OLD 2D model for SXY actually fitted SXY / (SXX^-2 + SYY^-2); correct here + // out[PM_PAR_SXY] = pmPSF_SXYtoModel (out); + + // the 2D PSF model fits polarization terms (E0,E1,E2) + // convert to shape terms (SXX,SYY,SXY) + pmPSF_FitToModel (out, 0.1); // apply the model limits here: this truncates excessive extrapolation + // XXX do we need to do this still? should we put in asserts to test? for (int i = 0; i < psf->params_NEW->n; i++) { // apply the limits to all components or just the psf-model parameters? Index: trunk/psModules/src/objects/pmPSF.c =================================================================== --- trunk/psModules/src/objects/pmPSF.c (revision 13038) +++ trunk/psModules/src/objects/pmPSF.c (revision 13064) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.21 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.22 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -236,4 +236,84 @@ } +// New Concept: the PSF modelling function fits the polarization terms e0, e1, e2: + +// convert the parameters used in the fitted source model +// to the parameters used in the 2D PSF model +bool pmPSF_FitToModel (psF32 *fittedPar, float minMinorAxis) +{ + psEllipsePol pol; + + pol.e0 = fittedPar[PM_PAR_E0]; + pol.e1 = fittedPar[PM_PAR_E1]; + pol.e2 = fittedPar[PM_PAR_E2]; + + psEllipseAxes axes = psEllipsePolToAxes (pol, minMinorAxis); + psEllipseShape shape = psEllipseAxesToShape (axes); + + fittedPar[PM_PAR_SXX] = shape.sx * M_SQRT2; + fittedPar[PM_PAR_SYY] = shape.sy * M_SQRT2; + fittedPar[PM_PAR_SXY] = shape.sxy; + + return true; +} + +// convert the PSF parameters used in the 2D PSF model fit into the +// parameters used in the source model +psEllipsePol pmPSF_ModelToFit (psF32 *modelPar) +{ + psEllipseShape shape; + + shape.sx = modelPar[PM_PAR_SXX] / M_SQRT2; + shape.sy = modelPar[PM_PAR_SYY] / M_SQRT2; + shape.sxy = modelPar[PM_PAR_SXY]; + + psEllipsePol pol = psEllipseShapeToPol (shape); + + return pol; +} + +// convert the parameters used in the fitted source model to the psEllipseAxes representation +// (major,minor,theta) +psEllipseAxes pmPSF_ModelToAxes (psF32 *modelPar, double maxAR) +{ + psEllipseShape shape; + psEllipseAxes axes; + + shape.sx = modelPar[PM_PAR_SXX] / M_SQRT2; + shape.sy = modelPar[PM_PAR_SYY] / M_SQRT2; + shape.sxy = modelPar[PM_PAR_SXY]; + + if ((shape.sx == 0) || (shape.sy == 0)) { + axes.major = 0.0; + axes.minor = 0.0; + axes.theta = 0.0; + } else { + // XXX this is not really consistent with the model fit range above + axes = psEllipseShapeToAxes (shape, maxAR); + } + + return axes; +} + +// convert the psEllipseAxes representation (major,minor,theta) to the parameters used in the +// fitted source model +bool pmPSF_AxesToModel (psF32 *modelPar, psEllipseAxes axes) +{ + if ((axes.major <= 0) || (axes.minor <= 0)) { + modelPar[PM_PAR_SXX] = 0.0; + modelPar[PM_PAR_SYY] = 0.0; + modelPar[PM_PAR_SXY] = 0.0; + return true; + } + + psEllipseShape shape = psEllipseAxesToShape (axes); + + modelPar[PM_PAR_SXX] = shape.sx * M_SQRT2; + modelPar[PM_PAR_SYY] = shape.sy * M_SQRT2; + modelPar[PM_PAR_SXY] = shape.sxy; + + return true; +} + // generate a psf model of the requested type, with fixed shape pmPSF *pmPSFBuildSimple (char *typeName, float sxx, float syy, float sxy, ...) @@ -247,13 +327,21 @@ pmPSF *psf = pmPSFAlloc (type, true, psfTrend); + psVector *par = psVectorAlloc (psf->params_NEW->n, PS_TYPE_F32); + par->data.F32[PM_PAR_SXX] = sxx; + par->data.F32[PM_PAR_SYY] = syy; + par->data.F32[PM_PAR_SXY] = sxy; + + psEllipsePol pol = pmPSF_ModelToFit(par->data.F32); + // set the psf shape parameters - psPolynomial2D *poly = psf->params_NEW->data[PM_PAR_SXX]; - poly->coeff[0][0] = M_SQRT2*sxx; - - poly = psf->params_NEW->data[PM_PAR_SYY]; - poly->coeff[0][0] = M_SQRT2*syy; - - poly = psf->params_NEW->data[PM_PAR_SXY]; - poly->coeff[0][0] = sxy; + psPolynomial2D *poly = NULL; + poly = psf->params_NEW->data[PM_PAR_E0]; + poly->coeff[0][0] = pol.e0; + + poly = psf->params_NEW->data[PM_PAR_E1]; + poly->coeff[0][0] = pol.e1; + + poly = psf->params_NEW->data[PM_PAR_E2]; + poly->coeff[0][0] = pol.e2; for (int i = PM_PAR_SXY + 1; i < psf->params_NEW->n; i++) { @@ -263,4 +351,5 @@ va_end(ap); + psFree (par); psFree (psfTrend); return psf; Index: trunk/psModules/src/objects/pmPSF.h =================================================================== --- trunk/psModules/src/objects/pmPSF.h (revision 13038) +++ trunk/psModules/src/objects/pmPSF.h (revision 13064) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.12 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-21 19:47:14 $ + * @version $Revision: 1.13 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * Copyright 2004 Maui High Performance Computing Center, University of Hawaii */ @@ -64,4 +64,8 @@ pmPSF; +# define PM_PAR_E0 PM_PAR_SXX +# define PM_PAR_E1 PM_PAR_SYY +# define PM_PAR_E2 PM_PAR_SXY + /** * @@ -99,4 +103,9 @@ pmPSF *pmPSFBuildSimple (char *typeName, float sxx, float syy, float sxy, ...); +bool pmPSF_AxesToModel (psF32 *modelPar, psEllipseAxes axes); +psEllipseAxes pmPSF_ModelToAxes (psF32 *modelPar, double maxAR); +bool pmPSF_FitToModel (psF32 *fittedPar, float minMinorAxis); +psEllipsePol pmPSF_ModelToFit (psF32 *modelPar); + /// @} # endif Index: trunk/psModules/src/objects/pmPSFtry.c =================================================================== --- trunk/psModules/src/objects/pmPSFtry.c (revision 13038) +++ trunk/psModules/src/objects/pmPSFtry.c (revision 13064) @@ -5,6 +5,6 @@ * @author EAM, IfA * - * @version $Revision: 1.36 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.37 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -352,4 +352,9 @@ x->data.F64[i] = source->modelEXT->params->data.F32[PM_PAR_XPOS]; y->data.F64[i] = source->modelEXT->params->data.F32[PM_PAR_YPOS]; + + // weight by the error on the source flux + if (dz) { + dz->data.F64[i] = source->modelEXT->dparams->data.F32[PM_PAR_I0]; + } } @@ -360,14 +365,50 @@ psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV); - // skip the unfitted parameters (X, Y, Io, Sky) + // The shape parameters (SXX, SXY, SYY) are strongly coupled. We have to handle them very + // carefully. First, we convert them to the Ellipse Polarization terms (E0, E1, E2) for + // each source and fit this set of parameters. These values are less tightly coupled, but + // are still inter-related. The fitted values do a good job of constraining the major axis + // and the position angle, but the minor axis is weakly measured. When we apply the PSF + // model to construct a source model, we convert the fitted values of E0,E1,E2 to the shape + // parameters, with the constraint that the minor axis must be greater than a minimum + // threshold. + + // convert the measured source shape paramters to polarization terms + psVector *e0 = psVectorAlloc (psfTry->sources->n, PS_TYPE_F64); + psVector *e1 = psVectorAlloc (psfTry->sources->n, PS_TYPE_F64); + psVector *e2 = psVectorAlloc (psfTry->sources->n, PS_TYPE_F64); + for (int i = 0; i < psfTry->sources->n; i++) { + pmSource *source = psfTry->sources->data[i]; + if (source->modelEXT == NULL) continue; + + psEllipsePol pol = pmPSF_ModelToFit (source->modelEXT->params->data.F32); + + e0->data.F64[i] = pol.e0; + e1->data.F64[i] = pol.e1; + e2->data.F64[i] = pol.e2; + } + + // we run 'clipIter' cycles clipping in each of x and y, with only one iteration each. + // This way, the parameters masked by one of the fits will be applied to the others + for (int i = 0; i < stats->clipIter; i++) { + psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E0], stats, psfTry->mask, 0xff, e0, dz, x, y); + psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E1], stats, psfTry->mask, 0xff, e1, dz, x, y); + psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E2], stats, psfTry->mask, 0xff, e2, dz, x, y); + } + + // skip the unfitted parameters (X, Y, Io, Sky) and the shape parameters (SXX, SYY, SXY) for (int i = 0; i < psf->params_NEW->n; i++) { - if (i == PM_PAR_SKY) - continue; - if (i == PM_PAR_I0) - continue; - if (i == PM_PAR_XPOS) - continue; - if (i == PM_PAR_YPOS) - continue; + switch (i) { + case PM_PAR_SKY: + case PM_PAR_I0: + case PM_PAR_XPOS: + case PM_PAR_YPOS: + case PM_PAR_SXX: + case PM_PAR_SYY: + case PM_PAR_SXY: + continue; + default: + break; + } // select the per-object fitted data for this parameter @@ -376,19 +417,10 @@ if (source->modelEXT == NULL) continue; - z->data.F64[j] = source->modelEXT->params->data.F32[i]; - if (applyWeight) { - dz->data.F64[j] = source->modelEXT->dparams->data.F32[i]; - } - - // for SXY, we actually fit SXY / (SXX^-2 + SYY^-2) - if (i == PM_PAR_SXY) { - z->data.F64[j] = pmPSF_SXYfromModel (source->modelEXT->params->data.F32); - } } // fit the collection of measured parameters to the PSF 2D model // the mask is carried from previous steps and updated with this operation - // the weight is either the parameter error or NULL, depending on 'applyWeights' + // the weight is either the flux error or NULL, depending on 'applyWeights' if (!psVectorClipFitPolynomial2D(psf->params_NEW->data[i], stats, psfTry->mask, 0xff, z, NULL, x, y)) { psError(PS_ERR_UNKNOWN, false, "failed to build psf model for parameter %d", i); Index: trunk/psModules/src/objects/pmSourceIO_CMF.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_CMF.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_CMF.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.20 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.21 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -68,8 +68,5 @@ lsky = (source->sky < 1.0) ? 0.0 : log10(source->sky); - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - axes = psEllipseShapeToAxes (shape, 20.0); + axes = pmPSF_ModelToAxes (PAR, 20.0); row = psMetadataAlloc (); @@ -158,9 +155,5 @@ axes.theta = psMetadataLookupF32 (&status, row, "THETA"); - shape = psEllipseAxesToShape (axes); - - PAR[PM_PAR_SXX] = shape.sx; - PAR[PM_PAR_SYY] = shape.sy; - PAR[PM_PAR_SXY] = shape.sxy; + pmPSF_AxesToModel (PAR, axes); source->psfMag = psMetadataLookupF32 (&status, row, "MAG_RAW") - ZERO_POINT; Index: trunk/psModules/src/objects/pmSourceIO_CMP.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_CMP.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_CMP.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.27 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.28 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -53,5 +53,4 @@ float lsky = 0; bool status; - psEllipseShape shape; psEllipseAxes axes; @@ -119,8 +118,5 @@ lsky = (source->sky < 1.0) ? 0.0 : log10(source->sky); - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - axes = psEllipseShapeToAxes (shape, 20.0); + axes = pmPSF_ModelToAxes (PAR, 20.0); float errMag = isfinite(source->errMag) ? source->errMag : 999; @@ -155,5 +151,4 @@ int Ninstar; psF32 *PAR, *dPAR; - psEllipseShape shape; psEllipseAxes axes; @@ -274,9 +269,5 @@ goto skip_source; - shape = psEllipseAxesToShape (axes); - - PAR[PM_PAR_SXX] = shape.sx; - PAR[PM_PAR_SYY] = shape.sy; - PAR[PM_PAR_SXY] = shape.sxy; + pmPSF_AxesToModel (PAR, axes); psArrayAdd (sources, 100, source); Index: trunk/psModules/src/objects/pmSourceIO_OBJ.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_OBJ.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_OBJ.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.11 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.12 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -45,5 +45,4 @@ psF32 *PAR, *dPAR; float dmag, apResid; - psEllipseShape shape; psEllipseAxes axes; @@ -78,8 +77,5 @@ } - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - axes = psEllipseShapeToAxes (shape, 20.0); + axes = pmPSF_ModelToAxes (PAR, 20.0); psLineInit (line); Index: trunk/psModules/src/objects/pmSourceIO_PS1_DEV_0.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_PS1_DEV_0.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_PS1_DEV_0.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.8 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -53,5 +53,4 @@ int i; psF32 *PAR, *dPAR; - psEllipseShape shape; psEllipseAxes axes; psF32 xPos, yPos; @@ -75,14 +74,5 @@ yErr = dPAR[PM_PAR_YPOS]; - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - if ((shape.sx == 0) || (shape.sy == 0)) { - axes.major = 0.0; - axes.minor = 0.0; - axes.theta = 0.0; - } else { - axes = psEllipseShapeToAxes (shape, 20.0); - } + axes = pmPSF_ModelToAxes (PAR, 20.0); } else { xPos = source->peak->xf; @@ -159,5 +149,4 @@ bool status; psF32 *PAR, *dPAR; - psEllipseShape shape; psEllipseAxes axes; @@ -203,17 +192,5 @@ source->errMag = psMetadataLookupF32 (&status, row, "PSF_INST_MAG_ERR"); - bool skip = false; - skip |= (axes.major <= 0); - skip |= (axes.minor <= 0); - if (!skip) { - shape = psEllipseAxesToShape (axes); - PAR[PM_PAR_SXX] = shape.sx; - PAR[PM_PAR_SYY] = shape.sy; - PAR[PM_PAR_SXY] = shape.sxy; - } else { - PAR[PM_PAR_SXX] = 0.0; - PAR[PM_PAR_SYY] = 0.0; - PAR[PM_PAR_SXY] = 0.0; - } + pmPSF_AxesToModel (PAR, axes); float lflux = psMetadataLookupF32 (&status, row, "PEAK_FLUX_AS_MAG"); Index: trunk/psModules/src/objects/pmSourceIO_SMPDATA.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_SMPDATA.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_SMPDATA.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.7 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.8 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -49,5 +49,4 @@ psF32 *PAR, *dPAR; bool status; - psEllipseShape shape; psEllipseAxes axes; psF32 xPos, yPos; @@ -74,14 +73,6 @@ lsky = (source->sky < 1.0) ? 0.0 : log10(source->sky); - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - if ((shape.sx == 0) || (shape.sy == 0)) { - axes.major = 0.0; - axes.minor = 0.0; - axes.theta = 0.0; - } else { - axes = psEllipseShapeToAxes (shape, 20.0); - } + axes = pmPSF_ModelToAxes (PAR, 20.0); + } else { xPos = source->peak->xf; @@ -134,5 +125,4 @@ bool status; psF32 *PAR, *dPAR; - psEllipseShape shape; psEllipseAxes axes; float lsky; @@ -177,17 +167,5 @@ axes.theta = psMetadataLookupF32 (&status, row, "THETA"); - bool skip = false; - skip |= (axes.major <= 0); - skip |= (axes.minor <= 0); - if (!skip) { - shape = psEllipseAxesToShape (axes); - PAR[PM_PAR_SXX] = shape.sx; - PAR[PM_PAR_SYY] = shape.sy; - PAR[PM_PAR_SXY] = shape.sxy; - } else { - PAR[PM_PAR_SXX] = 0.0; - PAR[PM_PAR_SYY] = 0.0; - PAR[PM_PAR_SXY] = 0.0; - } + pmPSF_AxesToModel (PAR, axes); source->psfMag = psMetadataLookupF32 (&status, row, "MAG_RAW") - ZERO_POINT; Index: trunk/psModules/src/objects/pmSourceIO_SX.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_SX.c (revision 13038) +++ trunk/psModules/src/objects/pmSourceIO_SX.c (revision 13064) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.10 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-26 01:20:29 $ + * @version $Revision: 1.11 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -43,5 +43,4 @@ psF32 *PAR, *dPAR; - psEllipseShape shape; psEllipseAxes axes; @@ -68,8 +67,5 @@ // pmSourceSextractType (source, &type, &flags); - shape.sx = PAR[PM_PAR_SXX]; - shape.sy = PAR[PM_PAR_SYY]; - shape.sxy = PAR[PM_PAR_SXY]; - axes = psEllipseShapeToAxes (shape, 20.0); + axes = pmPSF_ModelToAxes (PAR, 20.0); psLineInit (line); Index: trunk/psModules/src/objects/pmSourcePlotPSFModel.c =================================================================== --- trunk/psModules/src/objects/pmSourcePlotPSFModel.c (revision 13038) +++ trunk/psModules/src/objects/pmSourcePlotPSFModel.c (revision 13064) @@ -5,6 +5,6 @@ * @author EAM, IfA * - * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-21 19:47:14 $ + * @version $Revision: 1.7 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-04-27 22:14:08 $ * * Copyright 2006 IfA, University of Hawaii @@ -120,12 +120,7 @@ moments.xy = source->moments->Sxy; - psEllipseShape shape; - shape.sx = PAR[PM_PAR_SXX] / sqrt(2.0); - shape.sy = PAR[PM_PAR_SYY] / sqrt(2.0); - shape.sxy = PAR[PM_PAR_SXY]; - // force the axis ratio to be < 20.0 psEllipseAxes axes_mnt = psEllipseMomentsToAxes (moments, 20.0); - psEllipseAxes axes_psf = psEllipseShapeToAxes (shape, 20.0); + psEllipseAxes axes_psf = pmPSF_ModelToAxes (PAR, 20.0); // moments major axis