Index: trunk/psphot/src/psphotPSFResiduals.c
===================================================================
--- trunk/psphot/src/psphotPSFResiduals.c	(revision 12403)
+++ trunk/psphot/src/psphotPSFResiduals.c	(revision 12403)
@@ -0,0 +1,58 @@
+# include "psphot.h"
+
+// construct an initial PSF model for each object 
+bool psphotResidPSF (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) {
+
+  psTimerStart ("psphot");
+
+  for (int i = 0; i < sources->n; i++) {
+    pmSource *source = sources->data[i];
+
+    // only use the psf stars to build the residual model
+    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;
+    } 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;
+	}
+    }
+
+    // set PSF parameters for this model (apply 2D shape model)
+    pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
+    psFree (modelEXT);
+
+    // XXX need to define the guess flux?
+    // set the fit radius based on the object flux limit and the model
+    psphotCheckRadiusPSF (readout, source, modelPSF);
+
+    // set the source PSF model
+    source->modelPSF = modelPSF;
+  }
+  psLogMsg ("psphot.models", 4, "built models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot"));
+  return true;
+}
+
+// XXX do we always know which model is supposed to be used?