Index: trunk/psModules/src/objects/pmModel.c
===================================================================
--- trunk/psModules/src/objects/pmModel.c	(revision 13034)
+++ trunk/psModules/src/objects/pmModel.c	(revision 13038)
@@ -6,6 +6,6 @@
  *  @author EAM, IfA
  *
- *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-04-21 19:47:14 $
+ *  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-04-26 01:45:35 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -135,6 +135,5 @@
                    psImage *mask,
                    pmModel *model,
-                   bool center,
-                   bool sky,
+                   pmModelOpMode mode,
                    bool add
                       )
@@ -151,10 +150,29 @@
     psS32 imageCol;
     psS32 imageRow;
-    psF32 skyValue = params->data.F32[0];
     psF32 pixelValue;
     
-    float xCenter = model->params->data.F32[PM_PAR_XPOS];
-    float yCenter = model->params->data.F32[PM_PAR_YPOS];
-    float Io = model->params->data.F32[PM_PAR_I0];
+    // save original values; restore before returning
+    // use the true source position for the residual model
+    // the PSF model has presumably already been set for this coordinate
+    float xPos    = params->data.F32[PM_PAR_XPOS];
+    float yPos    = params->data.F32[PM_PAR_YPOS];
+    float IoSave  = params->data.F32[PM_PAR_I0];
+    float skySave = params->data.F32[PM_PAR_SKY];
+
+    if (mode & PM_MODEL_OP_NORM) {
+	params->data.F32[PM_PAR_I0] = 1.0;
+    }
+    if (!(mode & PM_MODEL_OP_SKY)) {
+	params->data.F32[PM_PAR_SKY] = 0.0;
+    }
+    if (mode & PM_MODEL_OP_CENTER) {
+	params->data.F32[PM_PAR_XPOS] = image->col0 + 0.5*image->numCols;
+	params->data.F32[PM_PAR_YPOS] = image->row0 + 0.5*image->numRows;
+    }
+
+    // use these values for this realization
+    float xCenter  = params->data.F32[PM_PAR_XPOS];
+    float yCenter  = params->data.F32[PM_PAR_YPOS];
+    float Io       = params->data.F32[PM_PAR_I0];
 
     int xBin = 1;
@@ -166,5 +184,5 @@
     psImageInterpolateOptions *Rx = NULL;
     psImageInterpolateOptions *Ry = NULL;
-    if (model->residuals) {
+    if (model->residuals && (mode & (PM_MODEL_OP_RES0 | PM_MODEL_OP_RES1))) {
 	Ro = psImageInterpolateOptionsAlloc(
 	    PS_INTERPOLATE_BILINEAR,
@@ -189,42 +207,41 @@
 
             // Convert i/j to image coord space:
-            // 'center' option changes meaning of i,j:
-            if (center) {
-                imageCol = ix - 0.5*image->numCols + xCenter;
-                imageRow = iy - 0.5*image->numRows + yCenter;
-            } else {
-                imageCol = ix + image->col0;
-                imageRow = iy + image->row0;
-            }
+	    // XXX should we use using 0.5 pixel offset?
+	    imageCol = ix + image->col0;
+	    imageRow = iy + image->row0;
 
             x->data.F32[0] = (float) imageCol;
             x->data.F32[1] = (float) imageRow;
 
-            // set the appropriate pixel value for this coordinate
-            if (sky) {
-                pixelValue = modelFunc (NULL, params, x);
-            } else {
-                pixelValue = modelFunc (NULL, params, x) - skyValue;
-            }
+	    pixelValue = 0.0;
+
+            // add in the desired components for this coordinate
+            if (mode & PM_MODEL_OP_FUNC) {
+                pixelValue += modelFunc (NULL, params, x);
+            } 
 
 	    // get the contribution from the residual model
-	    // XXX for a test, do this for all sources and all pixels
 	    if (Ro) {
 		// fractional image position
-		// this is wrong for the 'center' case
-		float ox = xBin*(ix + 0.5 + image->col0 - xCenter) + xResidCenter;
-		float oy = yBin*(iy + 0.5 + image->row0 - yCenter) + yResidCenter;
-
-		psU8 mflux = 0;
-		double Fo = 0.0;
-		double Fx = 0.0;
-		double Fy = 0.0;
-		psImageInterpolate (&Fo, NULL, &mflux, ox, oy, Ro);
-		psImageInterpolate (&Fx, NULL, &mflux, ox, oy, Rx);
-		psImageInterpolate (&Fy, NULL, &mflux, ox, oy, Ry);
-
-		if (!mflux && isfinite(Fo) && isfinite(Fx) && isfinite(Fy)) {
-		    double flux = Fo + xCenter*Fx + yCenter*Fy;
-		    pixelValue += Io*flux;
+		float ox = xBin*(imageCol + 0.5 - xCenter) + xResidCenter;
+		float oy = yBin*(imageRow + 0.5 - yCenter) + yResidCenter;
+
+		if (mode & PM_MODEL_OP_RES0) {
+		    psU8 mflux = 0;
+		    double Fo = 0.0;
+		    psImageInterpolate (&Fo, NULL, &mflux, ox, oy, Ro);
+		    if (!mflux && isfinite(Fo)) {
+			pixelValue += Io*Fo;
+		    }
+		}
+		if (mode & PM_MODEL_OP_RES1) {
+		    psU8 mflux = 0;
+		    double Fx = 0.0;
+		    double Fy = 0.0;
+		    psImageInterpolate (&Fx, NULL, &mflux, ox, oy, Rx);
+		    psImageInterpolate (&Fy, NULL, &mflux, ox, oy, Ry);
+		    if (!mflux && isfinite(Fx) && isfinite(Fy)) {
+			pixelValue += Io*(xPos*Fx + yPos*Fy);
+		    }
 		}
 	    }
@@ -238,4 +255,11 @@
         }
     }
+
+    // restore original values
+    params->data.F32[PM_PAR_I0]   = IoSave;
+    params->data.F32[PM_PAR_SKY]  = skySave;
+    params->data.F32[PM_PAR_XPOS] = xPos;
+    params->data.F32[PM_PAR_YPOS] = yPos;
+
     psFree(x);
     psFree(Ro);
@@ -251,9 +275,8 @@
                 psImage *mask,
                 pmModel *model,
-                bool center,
-                bool sky)
-{
-    psTrace("psModules.objects", 3, "---- %s() begin ----\n", __func__);
-    psBool rc = AddOrSubModel(image, mask, model, center, sky, true);
+                pmModelOpMode mode)
+{
+    psTrace("psModules.objects", 3, "---- %s() begin ----\n", __func__);
+    psBool rc = AddOrSubModel(image, mask, model, mode, true);
     psTrace("psModules.objects", 3, "---- %s(%d) end ----\n", __func__, rc);
     return(rc);
@@ -265,9 +288,8 @@
                 psImage *mask,
                 pmModel *model,
-                bool center,
-                bool sky)
-{
-    psTrace("psModules.objects", 3, "---- %s() begin ----\n", __func__);
-    psBool rc = AddOrSubModel(image, mask, model, center, sky, false);
+                pmModelOpMode mode)
+{
+    psTrace("psModules.objects", 3, "---- %s() begin ----\n", __func__);
+    psBool rc = AddOrSubModel(image, mask, model, mode, false);
     psTrace("psModules.objects", 3, "---- %s(%d) end ----\n", __func__, rc);
     return(rc);
