Index: /trunk/psphot/doc/notes.txt
===================================================================
--- /trunk/psphot/doc/notes.txt	(revision 14654)
+++ /trunk/psphot/doc/notes.txt	(revision 14655)
@@ -1,2 +1,49 @@
+
+2007.08.17
+
+ I am working on a number of cleanup / fixes.  I have made an overhaul
+ of the pmModel APIs, adding function pointers in the pmModel
+ structure to the class-specific utility functions (eg,
+ pmModel->modelFunc is the actual function which is evaluated). 
+
+ TO DO:
+
+ * update pmSourceFitSet to be able to include more than one model
+   type (currently it assumes the sources are all, eg, PSFs).  This is
+   now needed because the old implementation used the function lookups
+   which have been dropped.
+
+ * define a generic API set to handle 2D modelling of a scalar using
+   either polynomials (as the pmPSF code currently does) or an
+   image-based representation (as the psphot sky model currently
+   does).  
+
+   The image-based representation can automatically step down from NxM
+   super pixels to a smaller number based on the density of
+   measurements.  
+
+   Include ways to smooth and regularize the output result.
+
+   This mechanism can be applied to: psf parameters, aperture
+   residual, psf peak-to-flux variations, the psphot background
+   representation.
+
+ * generate and store the output radial profile for objects
+
+ * finish testing and incorportate the CR / EXT measurements
+
+ * adjustments to pixel center based on second derivatives : needed
+   for the sersic models.
+
+ * adjustments to pixel flux for extreme values (r ~ 0) : needed for
+   the sersic models.
+
+ * on psf stars : fall back on a Gaussian.
+
+ * OPTIMIZATIONS !!
+
+ * drop the model sky element : should only be in source->sky,dsky
+
+2006.11.16
 
 ensemble:
@@ -20,6 +67,4 @@
   * solve for source amplitudes
   * update models
-
-2006.11.16
 
   * create psSparseBorder to solve matrix equations which have a large
Index: /trunk/psphot/src/models/pmModel_STRAIL.c
===================================================================
--- /trunk/psphot/src/models/pmModel_STRAIL.c	(revision 14654)
+++ /trunk/psphot/src/models/pmModel_STRAIL.c	(revision 14655)
@@ -12,11 +12,12 @@
 *****************************************************************************/ 
  
-# define PM_MODEL_FUNC       pmModelFunc_STRAIL
-# define PM_MODEL_FLUX       pmModelFlux_STRAIL
-# define PM_MODEL_GUESS      pmModelGuess_STRAIL
-# define PM_MODEL_LIMITS     pmModelLimits_STRAIL
-# define PM_MODEL_RADIUS     pmModelRadius_STRAIL
-# define PM_MODEL_FROM_PSF   pmModelFromPSF_STRAIL
-# define PM_MODEL_FIT_STATUS pmModelFitStatus_STRAIL
+# define PM_MODEL_FUNC       	  pmModelFunc_STRAIL
+# define PM_MODEL_FLUX       	  pmModelFlux_STRAIL
+# define PM_MODEL_GUESS      	  pmModelGuess_STRAIL
+# define PM_MODEL_LIMITS     	  pmModelLimits_STRAIL
+# define PM_MODEL_RADIUS     	  pmModelRadius_STRAIL
+# define PM_MODEL_FROM_PSF   	  pmModelFromPSF_STRAIL
+# define PM_MODEL_PARAMS_FROM_PSF pmModelParamsFromPSF_STRAIL
+# define PM_MODEL_FIT_STATUS      pmModelFitStatus_STRAIL
 
 psF32 PM_MODEL_FUNC(psVector *deriv, 
@@ -474,6 +475,6 @@
  
 //fixed I think...no good way of guessing as far as I can tell 
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source) { 
- 
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+{ 
     pmMoments *Smoments = source->moments; 
     psF32     *params  = model->params->data.F32; 
@@ -516,6 +517,6 @@
  
 //fixed 
-bool PM_MODEL_FROM_PSF (pmModel *modelPSF, pmModel *modelFLT, pmPSF *psf) { 
- 
+bool PM_MODEL_FROM_PSF (pmModel *modelPSF, pmModel *modelFLT, pmPSF *psf)
+{
     psF32 *out = modelPSF->params->data.F32; 
     psF32 *in  = modelFLT->params->data.F32; 
@@ -529,5 +530,5 @@
  
     for (int i = 4; i < 7; i++) { 
-      psPolynomial2D *poly = psf->params_NEW->data[i-4]; 
+      psPolynomial2D *poly = psf->params->data[i-4]; 
 	out[i] = psPolynomial2DEval (poly, out[2], out[3]); 
     } 
@@ -535,7 +536,56 @@
 } 
  
+// construct the PSF model from the FLT model and the psf
+// XXX is this sufficiently general do be a global function, not a pmModelClass function?
+bool PM_MODEL_PARAMS_FROM_PSF (pmModel *model, pmPSF *psf, float Xo, float Yo, float Io)
+{
+    psF32 *PAR = model->params->data.F32;
+
+    // we require these two parameters to exist
+    assert (psf->params->n > PM_PAR_YPOS);
+    assert (psf->params->n > PM_PAR_XPOS);
+
+    PAR[PM_PAR_SKY]  = 0.0;
+    PAR[PM_PAR_I0]   = Io;
+    PAR[PM_PAR_XPOS] = Xo;
+    PAR[PM_PAR_YPOS] = Yo;
+    
+    // supply the model-fitted parameters, or copy from the input
+    for (int i = 0; i < psf->params->n; i++) {
+	if (i == PM_PAR_SKY) continue;
+	psPolynomial2D *poly = psf->params->data[i];
+	assert (poly);
+	PAR[i] = psPolynomial2DEval(poly, Xo, Yo);
+    }
+
+    // the 2D PSF model fits polarization terms (E0,E1,E2)
+    // convert to shape terms (SXX,SYY,SXY)
+    // XXX user-defined value for limit?
+    if (!pmPSF_FitToModel (PAR, 0.1)) {
+	psError(PM_ERR_PSF, false, "Failed to fit object at (r,c) = (%.1f,%.1f)", Xo, Yo);
+	return false;
+    }
+
+    // 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->n; i++) {
+        // apply the limits to all components or just the psf-model parameters?
+        if (psf->params->data[i] == NULL)
+            continue;
+
+	bool status = true;
+        status &= PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MIN, i, PAR, NULL);
+        status &= PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MAX, i, PAR, NULL);
+	if (!status) {
+	    psTrace ("psModules.objects", 5, "Hitting parameter limits at (r,c) = (%.1f, %.1f)", Xo, Yo);
+	    model->flags |= PM_MODEL_STATUS_LIMITS;
+	}
+    }
+    return(true);
+}
+
 //done I think 
-bool PM_MODEL_FIT_STATUS (pmModel *model) { 
- 
+bool PM_MODEL_FIT_STATUS (pmModel *model)
+{
     psF32 dP; 
     bool  status; 
@@ -566,3 +616,4 @@
 # undef PM_MODEL_RADIUS
 # undef PM_MODEL_FROM_PSF
+# undef PM_MODEL_PARAMS_FROM_PSF
 # undef PM_MODEL_FIT_STATUS
Index: /trunk/psphot/src/models/pmModel_TEST1.c
===================================================================
--- /trunk/psphot/src/models/pmModel_TEST1.c	(revision 14654)
+++ /trunk/psphot/src/models/pmModel_TEST1.c	(revision 14655)
@@ -16,19 +16,17 @@
  *****************************************************************************/
 
-# define PM_MODEL_FUNC       pmModelFunc_TEST1
-# define PM_MODEL_FLUX       pmModelFlux_TEST1
-# define PM_MODEL_GUESS      pmModelGuess_TEST1
-# define PM_MODEL_LIMITS     pmModelLimits_TEST1
-# define PM_MODEL_RADIUS     pmModelRadius_TEST1
-# define PM_MODEL_FROM_PSF   pmModelFromPSF_TEST1
-# define PM_MODEL_FIT_STATUS pmModelFitStatus_TEST1
-
-// XXX consider changing to form PAR[SXY]*(1/PAR[SXX]^2 + 1/PAR[SYY]^2)*X*Y
-// this would provide natural limits on PAR[SXY] of -0.25 : +0.25
+# define PM_MODEL_FUNC       	  pmModelFunc_TEST1
+# define PM_MODEL_FLUX       	  pmModelFlux_TEST1
+# define PM_MODEL_GUESS      	  pmModelGuess_TEST1
+# define PM_MODEL_LIMITS     	  pmModelLimits_TEST1
+# define PM_MODEL_RADIUS     	  pmModelRadius_TEST1
+# define PM_MODEL_FROM_PSF   	  pmModelFromPSF_TEST1
+# define PM_MODEL_PARAMS_FROM_PSF pmModelParamsFromPSF_TEST1
+# define PM_MODEL_FIT_STATUS      pmModelFitStatus_TEST1
 
 // the model is a function of the pixel coordinate (pixcoord[0,1] = x,y)
 psF32 PM_MODEL_FUNC(psVector *deriv,
-                         const psVector *params,
-                         const psVector *pixcoord)
+		    const psVector *params,
+		    const psVector *pixcoord)
 {
     psF32 *PAR = params->data.F32;
@@ -127,5 +125,4 @@
 bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
 {
-
     pmMoments *moments = source->moments;
     psF32     *PAR  = model->params->data.F32;
@@ -160,5 +157,5 @@
     norm = 0.0;
 
-    # define DZ 0.25
+# define DZ 0.25
 
     float f0 = 1.0;
@@ -209,17 +206,16 @@
 bool PM_MODEL_FROM_PSF (pmModel *modelPSF, pmModel *modelFLT, pmPSF *psf)
 {
-
     psF32 *out = modelPSF->params->data.F32;
     psF32 *in  = modelFLT->params->data.F32;
 
     // we require these two parameters to exist
-    assert (psf->params_NEW->n > PM_PAR_YPOS);
-    assert (psf->params_NEW->n > PM_PAR_XPOS);
-
-    for (int i = 0; i < psf->params_NEW->n; i++) {
-	if (psf->params_NEW->data[i] == NULL) {
+    assert (psf->params->n > PM_PAR_YPOS);
+    assert (psf->params->n > PM_PAR_XPOS);
+
+    for (int i = 0; i < psf->params->n; i++) {
+	if (psf->params->data[i] == NULL) {
 	    out[i] = in[i];
 	} else {	    
-	    psPolynomial2D *poly = psf->params_NEW->data[i];
+	    psPolynomial2D *poly = psf->params->data[i];
 	    out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]);
 	}
@@ -229,4 +225,53 @@
     out[PM_PAR_SXY] = pmPSF_SXYtoModel (out);
 
+    return(true);
+}
+
+// construct the PSF model from the FLT model and the psf
+// XXX is this sufficiently general do be a global function, not a pmModelClass function?
+bool PM_MODEL_PARAMS_FROM_PSF (pmModel *model, pmPSF *psf, float Xo, float Yo, float Io)
+{
+    psF32 *PAR = model->params->data.F32;
+
+    // we require these two parameters to exist
+    assert (psf->params->n > PM_PAR_YPOS);
+    assert (psf->params->n > PM_PAR_XPOS);
+
+    PAR[PM_PAR_SKY]  = 0.0;
+    PAR[PM_PAR_I0]   = Io;
+    PAR[PM_PAR_XPOS] = Xo;
+    PAR[PM_PAR_YPOS] = Yo;
+    
+    // supply the model-fitted parameters, or copy from the input
+    for (int i = 0; i < psf->params->n; i++) {
+	if (i == PM_PAR_SKY) continue;
+	psPolynomial2D *poly = psf->params->data[i];
+	assert (poly);
+	PAR[i] = psPolynomial2DEval(poly, Xo, Yo);
+    }
+
+    // the 2D PSF model fits polarization terms (E0,E1,E2)
+    // convert to shape terms (SXX,SYY,SXY)
+    // XXX user-defined value for limit?
+    if (!pmPSF_FitToModel (PAR, 0.1)) {
+	psError(PM_ERR_PSF, false, "Failed to fit object at (r,c) = (%.1f,%.1f)", Xo, Yo);
+	return false;
+    }
+
+    // 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->n; i++) {
+        // apply the limits to all components or just the psf-model parameters?
+        if (psf->params->data[i] == NULL)
+            continue;
+
+	bool status = true;
+        status &= PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MIN, i, PAR, NULL);
+        status &= PM_MODEL_LIMITS (PS_MINIMIZE_PARAM_MAX, i, PAR, NULL);
+	if (!status) {
+	    psTrace ("psModules.objects", 5, "Hitting parameter limits at (r,c) = (%.1f, %.1f)", Xo, Yo);
+	    model->flags |= PM_MODEL_STATUS_LIMITS;
+	}
+    }
     return(true);
 }
@@ -237,5 +282,4 @@
 bool PM_MODEL_FIT_STATUS (pmModel *model)
 {
-
     psF32 dP;
     bool  status;
@@ -265,3 +309,4 @@
 # undef PM_MODEL_RADIUS
 # undef PM_MODEL_FROM_PSF
+# undef PM_MODEL_PARAMS_FROM_PSF
 # undef PM_MODEL_FIT_STATUS
Index: /trunk/psphot/src/pmFootprint.c
===================================================================
--- /trunk/psphot/src/pmFootprint.c	(revision 14654)
+++ /trunk/psphot/src/pmFootprint.c	(revision 14655)
@@ -29,5 +29,5 @@
 }
 
-bool pmIsSpan(const psPtr ptr)
+bool pmSpanTest(const psPtr ptr)
 {
     return (psMemGetDeallocator(ptr) == (psFreeFunc)spanFree);
@@ -115,5 +115,5 @@
 }
 
-bool pmIsFootprint(const psPtr ptr) {
+bool pmFootprintTest(const psPtr ptr) {
     return (psMemGetDeallocator(ptr) == (psFreeFunc)footprintFree);
 }
@@ -762,5 +762,5 @@
    // Set stop bits from peaks list
    //
-   assert (peaks == NULL || peaks->n == 0 || pmIsPeak(peaks->data[0]));
+   assert (peaks == NULL || peaks->n == 0 || pmPeakTest(peaks->data[0]));
    if (peaks != NULL) {
        for (int i = 0; i < peaks->n; i++) {
@@ -859,5 +859,5 @@
    }
    const pmFootprint *fp = footprints->data[0];
-   assert(pmIsFootprint((const psPtr)fp));
+   assert(pmFootprintTest((const psPtr)fp));
    const int numCols = fp->region.x1 - fp->region.x0 + 1;
    const int numRows = fp->region.y1 - fp->region.y0 + 1;
@@ -884,5 +884,5 @@
 psImage *pmSetFootprintID(const pmFootprint *fp, // the footprint to insert
 			  const int id) {	// the desired ID
-   assert(fp != NULL && pmIsFootprint((const psPtr)fp));
+   assert(fp != NULL && pmFootprintTest((const psPtr)fp));
    const int numCols = fp->region.x1 - fp->region.x0 + 1;
    const int numRows = fp->region.y1 - fp->region.y0 + 1;
@@ -910,5 +910,5 @@
 psArray *pmGrowFootprintArray(const psArray *footprints, // footprints to grow
 			      int r) {	// how much to grow each footprint
-    assert (footprints->n == 0 || pmIsFootprint(footprints->data[0]));
+    assert (footprints->n == 0 || pmFootprintTest(footprints->data[0]));
 
     if (footprints->n == 0) {		// we don't know the size of the footprint's region
@@ -964,6 +964,6 @@
 				const psArray *footprints2, // the other set
 				const int includePeaks) { // which peaks to set? 0x1 => footprints1, 0x2 => 2
-    assert (footprints1->n == 0 || pmIsFootprint(footprints1->data[0]));
-    assert (footprints2->n == 0 || pmIsFootprint(footprints2->data[0]));
+    assert (footprints1->n == 0 || pmFootprintTest(footprints1->data[0]));
+    assert (footprints2->n == 0 || pmFootprintTest(footprints2->data[0]));
 
     if (footprints1->n == 0 || footprints2->n == 0) {		// nothing to do but put copies on merged
@@ -1032,7 +1032,7 @@
 			  const psArray *peaks) { // the pmPeaks
     assert (footprints != NULL);
-    assert (footprints->n == 0 || pmIsFootprint(footprints->data[0]));
+    assert (footprints->n == 0 || pmFootprintTest(footprints->data[0]));
     assert (peaks != NULL);
-    assert (peaks->n == 0 || pmIsPeak(peaks->data[0]));
+    assert (peaks->n == 0 || pmPeakTest(peaks->data[0]));
     
     if (footprints->n == 0) {
@@ -1222,5 +1222,5 @@
 psArray *pmFootprintArrayToPeaks(const psArray *footprints) {
    assert(footprints != NULL);
-   assert(footprints->n == 0 || pmIsFootprint(footprints->data[0]));
+   assert(footprints->n == 0 || pmFootprintTest(footprints->data[0]));
 
    int npeak = 0;
Index: /trunk/psphot/src/pmFootprint.h
===================================================================
--- /trunk/psphot/src/pmFootprint.h	(revision 14654)
+++ /trunk/psphot/src/pmFootprint.h	(revision 14655)
@@ -11,5 +11,5 @@
 
 pmSpan *pmSpanAlloc(int y, int x1, int x2);
-bool pmIsSpan(const psPtr ptr);
+bool pmSpanTest(const psPtr ptr);
 int pmSpanSortByYX (const void **a, const void **b);
 
@@ -25,5 +25,5 @@
 
 pmFootprint *pmFootprintAlloc(int nspan, const psImage *img);
-bool pmIsFootprint(const psPtr ptr);
+bool pmFootprintTest(const psPtr ptr);
 
 pmFootprint *pmFootprintNormalize(pmFootprint *fp);
Index: /trunk/psphot/src/psphot.c
===================================================================
--- /trunk/psphot/src/psphot.c	(revision 14654)
+++ /trunk/psphot/src/psphot.c	(revision 14655)
@@ -11,5 +11,5 @@
     pmErrorRegister();			// register psModule's error codes/messages
     psphotErrorRegister();              // register our error codes/messages
-    psphotModelGroupInit ();            // load implementation-specific models
+    psphotModelClassInit ();            // load implementation-specific models
 
     // load command-line arguments, options, and system config data
Index: /trunk/psphot/src/psphot.h
===================================================================
--- /trunk/psphot/src/psphot.h	(revision 14654)
+++ /trunk/psphot/src/psphot.h	(revision 14655)
@@ -51,5 +51,5 @@
 
 // basic support functions
-void            psphotModelGroupInit (void);
+void            psphotModelClassInit (void);
 int             pmPeakSortBySN (const void **a, const void **b);
 int             pmPeakSortByY (const void **a, const void **b);
Index: /trunk/psphot/src/psphotAddNoise.c
===================================================================
--- /trunk/psphot/src/psphotAddNoise.c	(revision 14654)
+++ /trunk/psphot/src/psphotAddNoise.c	(revision 14655)
@@ -68,5 +68,5 @@
 
         // XXX if we use pmSourceOp, the size (and possibly Io) will not be respected
-        pmSourceOp (source, PM_MODEL_OP_FULL | PM_MODEL_OP_NOISE, add, maskVal);
+        pmSourceOp (source, PM_MODEL_OP_FULL | PM_MODEL_OP_NOISE, add, maskVal, 0, 0);
 
         // restore original values
Index: /trunk/psphot/src/psphotChoosePSF.c
===================================================================
--- /trunk/psphot/src/psphotChoosePSF.c	(revision 14654)
+++ /trunk/psphot/src/psphotChoosePSF.c	(revision 14655)
@@ -127,6 +127,6 @@
     // print/dump psf parameters
     if (psTraceGetLevel("psphot") >= 5) {
-        for (int i = PM_PAR_SXX; i < try->psf->params_NEW->n; i++) {
-            psPolynomial2D *poly = try->psf->params_NEW->data[i];
+        for (int i = PM_PAR_SXX; i < try->psf->params->n; i++) {
+            psPolynomial2D *poly = try->psf->params->data[i];
             for (int nx = 0; nx <= poly->nX; nx++) {
                 for (int ny = 0; ny <= poly->nY; ny++) {
@@ -295,5 +295,5 @@
     }
 
-    char *modelName = pmModelGetType (psf->type);
+    char *modelName = pmModelClassGetName (psf->type);
     psLogMsg ("psphot.pspsf", PS_LOG_INFO, "select psf model: %f sec\n", psTimerMark ("psphot"));
     psLogMsg ("psphot.pspsf", PS_LOG_INFO, "psf model %s, ApResid: %f +/- %f\n", modelName, psf->ApResid, psf->dApResid);
@@ -315,26 +315,13 @@
     PS_ASSERT_PTR_NON_NULL(image, false);
 
-    pmModel *modelEXT = pmModelAlloc (psf->type);
-    PS_ASSERT_PTR_NON_NULL(modelEXT, false);
-
-    // make a model with unit central intensity at the image center
-    modelEXT->params->data.F32[PM_PAR_SKY] = 0;
-    modelEXT->params->data.F32[PM_PAR_I0] = 1;
-    modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*image->numCols;
-    modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*image->numRows;
-
-    // construct a PSF model at this coordinate
-    pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
+    // construct a normalized PSF model at this coordinate (Io = 1.0)
+    pmModel *modelPSF = pmModelFromPSFforXY (psf, 0.5*image->numCols, 0.5*image->numRows, 1.0);
     if (modelPSF == NULL) {
         psError(PSPHOT_ERR_PSF, false, "Failed to estimate PSF model at image centre");
-        psFree(modelEXT);
         return false;
     }
 
-    // get the correct model-radius function
-    pmModelRadius modelRadiusFunc = pmModelRadius_GetFunction (psf->type);
-
     // get the model full-width at half-max
-    psF64 FWHM_X = 2*modelRadiusFunc (modelPSF->params, 0.5);
+    psF64 FWHM_X = 2*modelPSF->modelRadius (modelPSF->params, 0.5);
 
     // XXX make sure this is consistent with the re-definition of PM_PAR_SXX
@@ -351,6 +338,4 @@
     psMetadataAddS32 (recipe, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", psf->nPSFstars);
     psMetadataAddBool(recipe, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "Valid PSF Model?", true);
-
-    psFree (modelEXT);
     psFree (modelPSF);
 
Index: /trunk/psphot/src/psphotCleanup.c
===================================================================
--- /trunk/psphot/src/psphotCleanup.c	(revision 14654)
+++ /trunk/psphot/src/psphotCleanup.c	(revision 14655)
@@ -7,5 +7,5 @@
     psTimerStop ();
     psMemCheckCorruption (stderr, true);
-    pmModelGroupCleanup ();
+    pmModelClassCleanup ();
     psTimeFinalize ();
     pmConceptsDone ();
Index: /trunk/psphot/src/psphotEvalFLT.c
===================================================================
--- /trunk/psphot/src/psphotEvalFLT.c	(revision 14654)
+++ /trunk/psphot/src/psphotEvalFLT.c	(revision 14655)
@@ -37,5 +37,5 @@
     } 
 
-    keep = pmModelFitStatus (model);
+    keep = model->modelFitStatus(model);
     if (keep) return true;
 
Index: /trunk/psphot/src/psphotFakeSources.c
===================================================================
--- /trunk/psphot/src/psphotFakeSources.c	(revision 14654)
+++ /trunk/psphot/src/psphotFakeSources.c	(revision 14655)
@@ -23,5 +23,5 @@
         source->type = PM_SOURCE_TYPE_STAR;
 
-        pmModelType modelType = pmModelSetType ("PS_MODEL_QGAUSS");
+        pmModelType modelType = pmModelClassGetType ("PS_MODEL_QGAUSS");
         source->modelPSF = pmSourceModelGuess (source, modelType);
         sources->data[i] = source;
Index: /trunk/psphot/src/psphotFindPeaks.c
===================================================================
--- /trunk/psphot/src/psphotFindPeaks.c	(revision 14654)
+++ /trunk/psphot/src/psphotFindPeaks.c	(revision 14655)
@@ -95,5 +95,5 @@
 
     // find the peaks in the smoothed image
-    psArray *peaks = pmFindImagePeaks (smooth_im, threshold);
+    psArray *peaks = pmPeaksInImage (smooth_im, threshold);
     if (peaks == NULL) {
         // XXX this may also be due to a programming or config error
Index: /trunk/psphot/src/psphotGrowthCurve.c
===================================================================
--- /trunk/psphot/src/psphotGrowthCurve.c	(revision 14654)
+++ /trunk/psphot/src/psphotGrowthCurve.c	(revision 14655)
@@ -14,7 +14,4 @@
     float radius;
 
-    // create template model
-    pmModel *modelRef = pmModelAlloc(psf->type);
-
     // use the center of the center pixel of the image
     xc = 0.5*readout->image->numCols + readout->image->col0 + 0.5;
@@ -23,13 +20,6 @@
     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
-    pmModel *model = pmModelFromPSF (modelRef, psf);
+    // create normalized model object at xc,yc 
+    pmModel *model = pmModelFromPSFforXY (psf, xc, yc, 1.0);
 
     // measure the fitMag for this model
@@ -76,5 +66,4 @@
     psFree (mask);
     psFree (model);
-    psFree (modelRef);
 
     return true;
Index: /trunk/psphot/src/psphotModelGroupInit.c
===================================================================
--- /trunk/psphot/src/psphotModelGroupInit.c	(revision 14654)
+++ /trunk/psphot/src/psphotModelGroupInit.c	(revision 14655)
@@ -7,18 +7,18 @@
 # include "models/pmModel_STRAIL.c"
 
-static pmModelGroup userModels[] = {
-    {"PS_MODEL_TEST1", 7, pmModelFunc_TEST1,  pmModelFlux_TEST1,  pmModelRadius_TEST1,  pmModelLimits_TEST1,  pmModelGuess_TEST1, pmModelFromPSF_TEST1, pmModelFitStatus_TEST1},
-    {"PS_MODEL_STRAIL", 9, pmModelFunc_STRAIL,  pmModelFlux_STRAIL,  pmModelRadius_STRAIL,  pmModelLimits_STRAIL,  pmModelGuess_STRAIL, pmModelFromPSF_STRAIL, pmModelFitStatus_STRAIL},
+static pmModelClass userModels[] = {
+    {"PS_MODEL_TEST1", 7, pmModelFunc_TEST1,  pmModelFlux_TEST1,  pmModelRadius_TEST1,  pmModelLimits_TEST1,  pmModelGuess_TEST1, pmModelFromPSF_TEST1, pmModelParamsFromPSF_TEST1, pmModelFitStatus_TEST1},
+    {"PS_MODEL_STRAIL", 9, pmModelFunc_STRAIL,  pmModelFlux_STRAIL,  pmModelRadius_STRAIL,  pmModelLimits_STRAIL,  pmModelGuess_STRAIL, pmModelFromPSF_STRAIL, pmModelParamsFromPSF_STRAIL, pmModelFitStatus_STRAIL},
 };
 
-void psphotModelGroupInit (void) 
+void psphotModelClassInit (void) 
 { 
 
-    // if pmModelGroupInit returns false, we have already init'ed
-    if (!pmModelGroupInit ()) return;
+    // if pmModelClassInit returns false, we have already init'ed
+    if (!pmModelClassInit ()) return;
 
-    int Nmodels = sizeof (userModels) / sizeof (pmModelGroup);
+    int Nmodels = sizeof (userModels) / sizeof (pmModelClass);
     for (int i = 0; i < Nmodels; i++) {
-	pmModelGroupAdd (&userModels[i]);
+	pmModelClassAdd (&userModels[i]);
     }
     return;
Index: /trunk/psphot/src/psphotModelTest.c
===================================================================
--- /trunk/psphot/src/psphotModelTest.c	(revision 14654)
+++ /trunk/psphot/src/psphotModelTest.c	(revision 14655)
@@ -89,5 +89,5 @@
             modelName = item->data.V;
         }
-        modelType = pmModelSetType (modelName);
+        modelType = pmModelClassGetType (modelName);
         if (modelType < 0) psAbort("unknown model %s", modelName);
 	source->type = PM_SOURCE_TYPE_EXTENDED;
@@ -118,5 +118,5 @@
             modelName = item->data.V;
         }
-        modelType = pmModelSetType (modelName);
+        modelType = pmModelClassGetType (modelName);
         if (modelType < 0) psAbort("unknown model %s", modelName);
 	source->type = PM_SOURCE_TYPE_EXTENDED;
@@ -148,5 +148,5 @@
 
     // if any parameters are defined by the user, take those values
-    int nParams = pmModelParameterCount (modelType);
+    int nParams = pmModelClassParameterCount (modelType);
     psF32 *params = model->params->data.F32;
     params[PM_PAR_XPOS] = xObj; // XXX use the user-supplied value, 
Index: /trunk/psphot/src/psphotPSFConvModel.c
===================================================================
--- /trunk/psphot/src/psphotPSFConvModel.c	(revision 14654)
+++ /trunk/psphot/src/psphotPSFConvModel.c	(revision 14655)
@@ -58,18 +58,8 @@
     psVector *dparams = modelConv->dparams;
 
-    // get the model function for this model
-    pmModelFunc modelFunc = pmModelFunc_GetFunction (modelConv->type);
-    if (!modelFunc)
-	psAbort("invalid model function");
-
-    // get the limits function for this model
-    pmModelLimits checkLimits = pmModelLimits_GetFunction (modelConv->type);
-    if (!checkLimits)
-	psAbort("invalid model limits function");
-
     // create the minimization constraints
     psMinConstraint *constraint = psMinConstraintAlloc();
     constraint->paramMask = psVectorAlloc (params->n, PS_TYPE_U8);
-    constraint->checkLimits = checkLimits;
+    constraint->checkLimits = modelConv->modelLimits;
 
     // set parameter mask based on fitting mode
@@ -81,6 +71,6 @@
     // force the floating parameters to fall within the contraint ranges
     for (int i = 0; i < params->n; i++) {
-	checkLimits (PS_MINIMIZE_PARAM_MIN, i, params->data.F32, NULL);
-	checkLimits (PS_MINIMIZE_PARAM_MAX, i, params->data.F32, NULL);
+	modelConv->modelLimits (PS_MINIMIZE_PARAM_MIN, i, params->data.F32, NULL);
+	modelConv->modelLimits (PS_MINIMIZE_PARAM_MAX, i, params->data.F32, NULL);
     }
 
@@ -90,5 +80,5 @@
     psImage *covar = psImageAlloc (params->n, params->n, PS_TYPE_F32);
 
-    bool fitStatus = psphotModelWithPSF_LMM (myMin, covar, params, constraint, source, psf, modelFunc);
+    bool fitStatus = psphotModelWithPSF_LMM (myMin, covar, params, constraint, source, psf, modelConv->modelFunc);
     for (int i = 0; i < dparams->n; i++) {
         if (psTraceGetLevel("psphot") >= 4) {
Index: /trunk/psphot/src/psphotPSFResiduals.c
===================================================================
--- /trunk/psphot/src/psphotPSFResiduals.c	(revision 14654)
+++ /trunk/psphot/src/psphotPSFResiduals.c	(revision 14655)
@@ -12,36 +12,16 @@
     if (!(source->mode & PM_SOURCE_MODE_PSFSTAR)) continue;
 
-    // XXX if a source is faint, it will not have moments measured.
-    // it must be modelled as a PSF.  In this case, we need to use 
-    // the peak centroid to get the coordinates and get the peak flux 
-    // from the image?
-
-    // use the source moments, etc to guess basic model parameters
-    pmModel *modelEXT = pmSourceModelGuess (source, psf->type);
-    
-    // XXX put this in a function of its own..
-    if (modelEXT == NULL) {
-	psErrorClear (); // XXX need to clear the error from failing the model
-	modelEXT = pmModelAlloc(psf->type);
-	psF32 *PAR = modelEXT->params->data.F32;
-	PAR[PM_PAR_SKY]  = 0;
-	// XXX get this from the image pixels
-	PAR[PM_PAR_I0]   = source->peak->flux;
-	PAR[PM_PAR_XPOS] = source->peak->xf;
-	PAR[PM_PAR_YPOS] = source->peak->yf;
+    if (source->mode & PM_SOURCE_MODE_SATSTAR) {
+	Xo = source->moments->x;
+	Yo = source->moments->y;
+	Io = source->peak->flux;
     } else {
-	// these valuse are set in pmSourceModelGuess, should this rule be in there as well?
-	if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
-	    modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->x;
-	    modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->y;
-	} else {
-	    modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
-	    modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
-	}
+	Xo = source->peak->xf;
+	Yo = source->peak->yf;
+	Io = source->peak->flux;
     }
 
     // set PSF parameters for this model (apply 2D shape model)
-    pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
-    psFree (modelEXT);
+    pmModel *modelPSF = pmModelFromPSFforXY (psf, Xo, Yo, Io);
 
     // XXX need to define the guess flux?
Index: /trunk/psphot/src/psphotRadiusChecks.c
===================================================================
--- /trunk/psphot/src/psphotRadiusChecks.c	(revision 14654)
+++ /trunk/psphot/src/psphotRadiusChecks.c	(revision 14655)
@@ -6,5 +6,4 @@
 static float PSF_FIT_RADIUS = 0;	// radius to use in fitting (ignored if <= 0,
 					// and a per-object radius is calculated)
-static pmModelRadius modelRadiusPSF;
 
 bool psphotInitRadiusPSF(const psMetadata *recipe,
@@ -17,6 +16,4 @@
     PSF_FIT_RADIUS =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
 
-    // this function specifies the radius at this the model hits the given flux
-    modelRadiusPSF       = pmModelRadius_GetFunction (type);
     return true;
 }
@@ -34,7 +31,7 @@
     if (radiusFit <= 0) {		// use fixed radius
 	if (moments == NULL) {
-	    radiusFit = modelRadiusPSF(model->params, PSF_FIT_NSIGMA*moments->dSky);
+	    radiusFit = model->modelRadius(model->params, PSF_FIT_NSIGMA*moments->dSky);
 	} else {
-	    radiusFit = modelRadiusPSF(model->params, 1.0);
+	    radiusFit = model->modelRadius(model->params, 1.0);
 	}
     }
@@ -62,5 +59,5 @@
 
     // set the fit radius based on the object flux limit and the model
-    model->radiusFit = (RADIUS_TYPE) (modelRadiusPSF (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
+    model->radiusFit = (RADIUS_TYPE) (model->modelRadius (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
     if (isnan(model->radiusFit)) psAbort("error in radius");
 	
@@ -78,5 +75,4 @@
 static float EXT_FIT_NSIGMA;
 static float EXT_FIT_PADDING;
-static pmModelRadius modelRadiusEXT;
 
 bool psphotInitRadiusEXT (psMetadata *recipe, pmModelType type) {
@@ -87,6 +83,4 @@
     EXT_FIT_PADDING  = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");
 
-    // this function specifies the radius at this the model hits the given flux
-    modelRadiusEXT       = pmModelRadius_GetFunction (type);
     return true;
 }
@@ -101,5 +95,5 @@
 
     // set the fit radius based on the object flux limit and the model
-    model->radiusFit = (RADIUS_TYPE) (modelRadiusEXT (model->params, EXT_FIT_NSIGMA*moments->dSky) + EXT_FIT_PADDING);
+    model->radiusFit = (RADIUS_TYPE) (model->modelRadius (model->params, EXT_FIT_NSIGMA*moments->dSky) + EXT_FIT_PADDING);
     if (isnan(model->radiusFit)) psAbort("error in radius");
 
Index: /trunk/psphot/src/psphotReadout.c
===================================================================
--- /trunk/psphot/src/psphotReadout.c	(revision 14654)
+++ /trunk/psphot/src/psphotReadout.c	(revision 14655)
@@ -251,5 +251,5 @@
 
     // plot positive sources
-    // psphotSourcePlots (readout, sources, recipe);
+    psphotSourcePlots (readout, sources, recipe, maskVal);
 
     // measure aperture photometry corrections
Index: /trunk/psphot/src/psphotReadoutCleanup.c
===================================================================
--- /trunk/psphot/src/psphotReadoutCleanup.c	(revision 14654)
+++ /trunk/psphot/src/psphotReadoutCleanup.c	(revision 14655)
@@ -57,6 +57,10 @@
     if (psf) {
 	// save the psf for possible output.  if there was already an entry, it was loaded from external sources
-	// the new one may have been updated or modified, so replace the existing entry
-	psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSPHOT.PSF", PS_DATA_UNKNOWN | PS_META_REPLACE,  "psphot psf", psf);
+	// the new one may have been updated or modified, so replace the existing entry.  We
+	// are required to save it on the chip, but this will cause problems if we ever want to
+	// run psphot on an unmosaiced image
+	pmCell *cell = readout->parent;
+	pmChip *chip = cell->parent;
+	psMetadataAdd (chip->analysis, PS_LIST_TAIL, "PSPHOT.PSF", PS_DATA_UNKNOWN | PS_META_REPLACE,  "psphot psf", psf);
     }
 
Index: /trunk/psphot/src/psphotSourceFits.c
===================================================================
--- /trunk/psphot/src/psphotSourceFits.c	(revision 14654)
+++ /trunk/psphot/src/psphotSourceFits.c	(revision 14655)
@@ -195,5 +195,5 @@
     // extended source model descriptions
     char *modelNameEXT = psMetadataLookupStr (&status, recipe, "EXT_MODEL");
-    modelTypeEXT = pmModelSetType (modelNameEXT);
+    modelTypeEXT = pmModelClassGetType (modelNameEXT);
     psphotInitRadiusEXT (recipe, modelTypeEXT);
 
Index: /trunk/psphot/src/psphotTestSourceOutput.c
===================================================================
--- /trunk/psphot/src/psphotTestSourceOutput.c	(revision 14654)
+++ /trunk/psphot/src/psphotTestSourceOutput.c	(revision 14655)
@@ -21,5 +21,4 @@
     psVector *x = psVectorAlloc(2, PS_TYPE_F32);
     psVector *params = model->params;
-    pmModelFunc modelFunc = pmModelFunc_GetFunction (model->type);
     psS32 imageCol;
     psS32 imageRow;
@@ -68,5 +67,5 @@
             // set the appropriate pixel value for this coordinate
 	    if (mode & PSPHOT_ADD_MODEL) {
-		pixelValue = modelFunc (NULL, params, x) - skyValue;
+		pixelValue = model->modelFunc (NULL, params, x) - skyValue;
 	    } else {
 		pixelValue = 0.0;