Changeset 35768 for trunk/psModules/src/objects/pmPSF.c

Timestamp:

Jul 3, 2013, 2:37:22 PM (13 years ago)

Author:

eugene

Message:

deprecate KiiOpen,KiiClose (now KapaOpen,etc); major rework of psEllipse translations : use common functions pmModelAxesToParams and pmModelParamsToAxes ; use new convergence method in pmPCM_MinimizeChisq; add convergence crerition options to psMinimization; threaded versions of pmPSFtryFitEXT and pmPSFtryFitPSF

Location:

trunk/psModules

Files:

• . (modified) (1 prop)
• src/objects/pmPSF.c (modified) (6 diffs)

trunk/psModules/src/objects/pmPSF.c

@@ -249 +249 @@
 // Mxy = SXY * (SXX^-4 + SYY^-4 - 2 SXY ^2)
 
+// XXX deprecated
 // input: model->param, output: psf->param[PM_PAR_SXY]
-double pmPSF_SXYfromModel (psF32 *modelPar)
-{
-    PS_ASSERT_PTR_NON_NULL(modelPar, NAN);
-
-    double SXX = modelPar[PM_PAR_SXX];
-    double SYY = modelPar[PM_PAR_SYY];
-    double SXY = modelPar[PM_PAR_SXY];
-
-    double par = SXY / PS_SQR(1.0 / PS_SQR(SXX) + 1.0 / PS_SQR(SYY));
-    return (par);
-}
-
+// XXX double pmPSF_SXYfromModel (psF32 *modelPar)
+// XXX {
+// XXX     PS_ASSERT_PTR_NON_NULL(modelPar, NAN);
+// XXX 
+// XXX     double SXX = modelPar[PM_PAR_SXX];
+// XXX     double SYY = modelPar[PM_PAR_SYY];
+// XXX     double SXY = modelPar[PM_PAR_SXY];
+// XXX 
+// XXX     double par = SXY / PS_SQR(1.0 / PS_SQR(SXX) + 1.0 / PS_SQR(SYY));
+// XXX     return (par);
+// XXX }
+
+// XXX deprecated
 // input: fitted psf->param, output: model->param[PM_PAR_SXY]
-double pmPSF_SXYtoModel (psF32 *fittedPar)
-{
-    PS_ASSERT_PTR_NON_NULL(fittedPar, NAN);
-
-    double SXX = fittedPar[PM_PAR_SXX];
-    double SYY = fittedPar[PM_PAR_SYY];
-    double fit = fittedPar[PM_PAR_SXY];
-
-    double SXY = fit * PS_SQR(1.0 / PS_SQR(SXX) + 1.0 / PS_SQR(SYY));
-
-    assert (!isnan(SXY));
-
-    return SXY;
-}
-
-// 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
-// XXX this function may be invalid for SERSIC, DEV, EXP models (SQRT2 not used?)
-bool pmPSF_FitToModel (psF32 *fittedPar, float minMinorAxis)
+// XXX double pmPSF_SXYtoModel (psF32 *fittedPar)
+// XXX {
+// XXX     PS_ASSERT_PTR_NON_NULL(fittedPar, NAN);
+// XXX 
+// XXX     double SXX = fittedPar[PM_PAR_SXX];
+// XXX     double SYY = fittedPar[PM_PAR_SYY];
+// XXX     double fit = fittedPar[PM_PAR_SXY];
+// XXX 
+// XXX     double SXY = fit * PS_SQR(1.0 / PS_SQR(SXX) + 1.0 / PS_SQR(SYY));
+// XXX 
+// XXX     assert (!isnan(SXY));
+// XXX 
+// XXX     return SXY;
+// XXX }
+
+// The PSF modelling function fits the polarization terms e0, e1, e2:
+
+// the FIT is the 2D representation of the shape using polarization parameters for the elliptical contour
+// the MODEL is the realized psf model for a given location
+
+// convert the parameters (in situ) used in the fitted source model to the parameters used in
+// the 2D PSF model
+bool pmPSF_FitToModel (psF32 *fittedPar, float minMinorAxis, bool useReff)
 {
     PS_ASSERT_PTR_NON_NULL(fittedPar, false);
@@ -298 +302 @@
         return false;
     }
-    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;
-
+
+    pmModelAxesToParams (&fittedPar[PM_PAR_SXX], &fittedPar[PM_PAR_SXY], &fittedPar[PM_PAR_SYY], axes, useReff);
     return true;
 }
 
-// convert the PSF parameters used in the 2D PSF model fit into the
-// parameters used in the source model
-// XXX this function may be invalid for SERSIC, DEV, EXP models (SQRT2 not used?)
-psEllipsePol pmPSF_ModelToFit (psF32 *modelPar)
+// convert the parameters (in situ) used in the 2D PSF model fit into the parameters used in
+// the source model
+psEllipsePol pmPSF_ModelToFit (psF32 *modelPar, bool useReff)
 {
     // must assert non-NULL input parameter
@@ -319 +318 @@
     PS_ASSERT_PTR_NON_NULL(modelPar, pol);
 
-    psEllipseShape shape;
-
-    shape.sx  = modelPar[PM_PAR_SXX] / M_SQRT2;
-    shape.sy  = modelPar[PM_PAR_SYY] / M_SQRT2;
-    shape.sxy = modelPar[PM_PAR_SXY];
-
-    pol = psEllipseShapeToPol (shape);
+    psEllipseAxes axes;
+    pmModelParamsToAxes (&axes, modelPar[PM_PAR_SXX], modelPar[PM_PAR_SXY], modelPar[PM_PAR_SYY], useReff);
+
+    pol = psEllipseAxesToPol (axes);
 
     return pol;
@@ -332 +328 @@
 // convert the parameters used in the fitted source model to the psEllipseAxes representation
 // (major,minor,theta)
-psEllipseAxes pmPSF_ModelToAxes (psF32 *modelPar, double maxAR, pmModelType type)
-{
-    psEllipseShape shape;
+psEllipseAxes pmPSF_ModelToAxes (psF32 *modelPar, pmModelType type)
+{
     psEllipseAxes axes;
     axes.major = NAN;
     axes.minor = NAN;
     axes.theta = NAN;
-    //   XXX: must assert non-NULL input parameter
+
     PS_ASSERT_PTR_NON_NULL(modelPar, axes);
 
-    bool useReff = false;
-    useReff |= (type == pmModelClassGetType ("PS_MODEL_SERSIC"));
-    useReff |= (type == pmModelClassGetType ("PS_MODEL_DEV"));
-    useReff |= (type == pmModelClassGetType ("PS_MODEL_EXP"));
-
-    if (useReff) {
-        shape.sx  = modelPar[PM_PAR_SXX];
-        shape.sy  = modelPar[PM_PAR_SYY];
-        shape.sxy = modelPar[PM_PAR_SXY] / 2.0;
-        // XXX I *think* dividing by 2.0 is the right direction, but this 
-        // needs to be checked with a real test
-    } else {
-        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);
-    }
-
+    bool useReff = pmModelUseReff (type);
+    pmModelParamsToAxes (&axes, modelPar[PM_PAR_SXX], modelPar[PM_PAR_SXY], modelPar[PM_PAR_SYY], useReff);
     return axes;
 }
@@ -377 +348 @@
     PS_ASSERT_PTR_NON_NULL(modelPar, false);
 
+    modelPar[PM_PAR_SXX] = 0.0;
+    modelPar[PM_PAR_SYY] = 0.0;
+    modelPar[PM_PAR_SXY] = 0.0;
+    
     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);
-
-    bool useReff = false;
-    useReff |= ( type == pmModelClassGetType ("PS_MODEL_SERSIC"));
-    useReff |= ( type == pmModelClassGetType ("PS_MODEL_DEV"));
-    useReff |= ( type == pmModelClassGetType ("PS_MODEL_EXP"));
-
-    if (useReff) {
-        modelPar[PM_PAR_SXX] = shape.sx;
-        modelPar[PM_PAR_SYY] = shape.sy;
-        modelPar[PM_PAR_SXY] = shape.sxy * 2.0; // XXX NEED factor of 2 here for correct angle conversion
-    } else {
-        modelPar[PM_PAR_SXX] = shape.sx * M_SQRT2;
-        modelPar[PM_PAR_SYY] = shape.sy * M_SQRT2;
-        modelPar[PM_PAR_SXY] = shape.sxy;
-    }
+    
+    bool useReff = pmModelUseReff (type);
+    pmModelAxesToParams (&modelPar[PM_PAR_SXX], &modelPar[PM_PAR_SXY], &modelPar[PM_PAR_SYY], axes, useReff);
     return true;
 }
@@ -423 +381 @@
     par->data.F32[PM_PAR_SXY] = sxy;
 
-    psEllipsePol pol = pmPSF_ModelToFit(par->data.F32);
+    bool useReff = pmModelUseReff (options->type);
+    psEllipsePol pol = pmPSF_ModelToFit(par->data.F32, useReff);
 
     pmTrend2D *trend = NULL;