Index: /trunk/psModules/src/objects/pmSourcePhotometry.c
===================================================================
--- /trunk/psModules/src/objects/pmSourcePhotometry.c	(revision 9526)
+++ /trunk/psModules/src/objects/pmSourcePhotometry.c	(revision 9527)
@@ -3,6 +3,6 @@
  *  @author EAM, IfA; GLG, MHPCC
  *
- *  @version $Revision: 1.8 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2006-09-15 09:49:01 $
+ *  @version $Revision: 1.9 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-10-13 02:24:34 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -95,6 +95,6 @@
     }
 
-    if (model->dparams->data.F32[1] > 0) {
-        SN = model->params->data.F32[1] / model->dparams->data.F32[1];
+    if (model->dparams->data.F32[PM_PAR_FLUX] > 0) {
+        SN = model->params->data.F32[PM_PAR_FLUX] / model->dparams->data.F32[PM_PAR_FLUX];
         source->errMag = 1.0 / SN;
     } else {
@@ -102,6 +102,6 @@
         source->errMag = 0.0;
     }
-    x = model->params->data.F32[2];
-    y = model->params->data.F32[3];
+    x = model->params->data.F32[PM_PAR_XPOS];
+    y = model->params->data.F32[PM_PAR_YPOS];
 
     // measure object model photometry
@@ -194,5 +194,5 @@
 
     if (DO_SKY) {
-        sky = model->params->data.F32[0];
+        sky = model->params->data.F32[PM_PAR_SKY];
     } else {
         sky = 0;
@@ -233,5 +233,5 @@
     // we only care about the value of the object model, not the local sky
     if (DO_SKY) {
-        sky = model->params->data.F32[0];
+        sky = model->params->data.F32[PM_PAR_SKY];
     } else {
         sky = 0;
@@ -245,6 +245,6 @@
     psVector *params = model->params;
 
-    Xo = params->data.F32[2];
-    Yo = params->data.F32[3];
+    Xo = params->data.F32[PM_PAR_XPOS];
+    Yo = params->data.F32[PM_PAR_YPOS];
 
     dX = Xo - image->col0;
@@ -297,5 +297,7 @@
 }
 
-float pmSourceCrossProduct (pmSource *Mi, pmSource *Mj)
+double pmSourceCrossProduct (const pmSource *Mi,
+                             const pmSource *Mj,
+                             const bool unweighted_sum) // should the cross product be weighted?
 {
 
@@ -304,13 +306,13 @@
     int xIs, xJs, yIs, yJs;
     int xIe, yIe;
-    float flux, wt;
-
-    psImage *Pi = Mi->pixels;
-    psImage *Pj = Mj->pixels;
-
-    psImage *Wi = Mi->weight;
-
-    psImage *Ti = Mi->mask;
-    psImage *Tj = Mj->mask;
+    double flux, wt;
+
+    const psImage *Pi = Mi->pixels;
+    const psImage *Pj = Mj->pixels;
+
+    const psImage *Wi = Mi->weight;
+
+    const psImage *Ti = Mi->mask;
+    const psImage *Tj = Mj->mask;
 
     Xs = PS_MAX (Pi->col0, Pj->col0);
@@ -337,14 +339,21 @@
             if (Tj->data.U8[yj][xj])
                 continue;
-            wt = Wi->data.F32[yi][xi];
-            if (wt > 0) {
-                flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]) / wt;
+
+            if (unweighted_sum) {
+                flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]);
+            } else {
+                wt = Wi->data.F32[yi][xi];
+                if (wt > 0) {
+                    flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]) / wt;
+                }
             }
         }
     }
-    return (flux);
-}
-
-float pmSourceCrossWeight (pmSource *Mi, pmSource *Mj)
+    return flux;
+}
+
+double pmSourceCrossWeight(const pmSource *Mi,
+                           const pmSource *Mj,
+                           const bool unweighted_sum) // should the cross product be weighted?
 {
 
@@ -353,13 +362,13 @@
     int xIs, xJs, yIs, yJs;
     int xIe, yIe;
-    float flux, wt;
-
-    psImage *Pi = Mi->pixels;
-    psImage *Pj = Mj->pixels;
-
-    psImage *Wi = Mi->weight;
-
-    psImage *Ti = Mi->mask;
-    psImage *Tj = Mj->mask;
+    double flux, wt;
+
+    const psImage *Pi = Mi->pixels;
+    const psImage *Pj = Mj->pixels;
+
+    const psImage *Wi = Mi->weight;
+
+    const psImage *Ti = Mi->mask;
+    const psImage *Tj = Mj->mask;
 
     Xs = PS_MAX (Pi->col0, Pj->col0);
@@ -386,11 +395,16 @@
             if (Tj->data.U8[yj][xj])
                 continue;
-            wt = Wi->data.F32[yi][xi];
-            if (wt > 0) {
-                flux += 1.0 / wt;
+
+            if (unweighted_sum) {
+                flux++;
+            } else {
+                wt = Wi->data.F32[yi][xi];
+                if (wt > 0) {
+                    flux += 1.0 / wt;
+                }
             }
         }
     }
-    return (flux);
-}
-
+    return flux;
+}
+
Index: /trunk/psModules/src/objects/pmSourcePhotometry.h
===================================================================
--- /trunk/psModules/src/objects/pmSourcePhotometry.h	(revision 9526)
+++ /trunk/psModules/src/objects/pmSourcePhotometry.h	(revision 9527)
@@ -3,6 +3,6 @@
  *  @author EAM, IfA; GLG, MHPCC
  *
- *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2006-04-17 18:01:05 $
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-10-13 02:24:34 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -47,6 +47,6 @@
 bool pmSourceMagnitudesInit (psMetadata *config);
 bool pmSourceMagnitudes (pmSource *source, pmPSF *psf, pmSourcePhotometryMode mode);
-float pmSourceCrossProduct (pmSource *Mi, pmSource *Mj);
-float pmSourceCrossWeight (pmSource *Mi, pmSource *Mj);
+double pmSourceCrossProduct(const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum);
+double pmSourceCrossWeight(const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum);
 bool pmSourcePixelWeight (float *pixWeight, pmModel *model, psImage *image, psImage *mask);
 
Index: /trunk/psphot/src/psphotEnsemblePSF.c
===================================================================
--- /trunk/psphot/src/psphotEnsemblePSF.c	(revision 9526)
+++ /trunk/psphot/src/psphotEnsemblePSF.c	(revision 9527)
@@ -37,4 +37,10 @@
     // psRegion AnalysisRegion = psRegionForImage (readout->image, psRegionFromString (region));
     if (psRegionIsNaN (AnalysisRegion)) psAbort ("psphot", "analysis region mis-defined");
+
+    const bool CONSTANT_PHOTOMETRIC_WEIGHTS =
+	psMetadataLookupBool(&status, recipe, "CONSTANT_PHOTOMETRIC_WEIGHTS");
+    if (!status) {
+	psAbort(PS_FILE_LINE, "You must provide a value for the BOOL recipe CONSTANT_PHOTOMETRIC_WEIGHTS");
+    }
 
     for (int i = 0; i < sources->n; i++) {
@@ -68,6 +74,6 @@
 	pmModel *modelEXT = pmSourceModelGuess (inSource, psf->type); 
 	if (inSource->mode &  PM_SOURCE_MODE_SATSTAR) {
-	    modelEXT->params->data.F32[2] = inSource->moments->x;
-	    modelEXT->params->data.F32[3] = inSource->moments->y;
+	    modelEXT->params->data.F32[PM_PAR_XPOS] = inSource->moments->x;
+	    modelEXT->params->data.F32[PM_PAR_YPOS] = inSource->moments->y;
 	} else {
 	    // peak-up on peak (for non-sat objects)
@@ -81,10 +87,10 @@
 
 	    psTrace ("psphotEnsemblePSF", 5, "peak coord: %f %f -> %f %f\n", 
-		     modelEXT->params->data.F32[2], modelEXT->params->data.F32[3], min.x + ix, min.y + iy);
+		     modelEXT->params->data.F32[PM_PAR_XPOS], modelEXT->params->data.F32[PM_PAR_YPOS], min.x + ix, min.y + iy);
 	    
 	    // if min point is too deviant, keep the old value
 	    if ((fabs(min.x) < 1.5) && (fabs(min.y) < 1.5)) {
-		modelEXT->params->data.F32[2] = min.x + ix;
-		modelEXT->params->data.F32[3] = min.y + iy;
+		modelEXT->params->data.F32[PM_PAR_XPOS] = min.x + ix;
+		modelEXT->params->data.F32[PM_PAR_YPOS] = min.y + iy;
 	    }
 	    psFree (bicube);
@@ -96,6 +102,6 @@
 
 	// save the original coords
-	x = model->params->data.F32[2];
-	y = model->params->data.F32[3];
+	x = model->params->data.F32[PM_PAR_XPOS];
+	y = model->params->data.F32[PM_PAR_YPOS];
 
 	// set the fit radius based on the object flux limit and the model
@@ -112,6 +118,6 @@
 
 	// set model to unit peak, zero sky (we assume sky is subtracted)
-	model->params->data.F32[0] = 0.0;
-	model->params->data.F32[1] = 1.0;
+	model->params->data.F32[PM_PAR_SKY] = 0.0;
+	model->params->data.F32[PM_PAR_FLUX] = 1.0;
 
 	// fill in the model pixel values
@@ -143,5 +149,5 @@
 
 	// scale by diagonal element (auto-cross-product)
-	r = pmSourceCrossProduct (Mi, Mi);
+	r = pmSourceCrossProduct (Mi, Mi, CONSTANT_PHOTOMETRIC_WEIGHTS);
 	weight->data.F32[i] = r;
 
@@ -149,5 +155,5 @@
 
 	// find the image x model value
-	f = pmSourceCrossProduct (Fi, Mi);
+	f = pmSourceCrossProduct (Fi, Mi, CONSTANT_PHOTOMETRIC_WEIGHTS);
 	psSparseVectorElement (sparse, i, f / r);
 
@@ -163,5 +169,5 @@
 	    
 	    // got an overlap; calculate cross-product and add to output array
-	    f = pmSourceCrossProduct (Mi, Mj);
+	    f = pmSourceCrossProduct (Mi, Mj, CONSTANT_PHOTOMETRIC_WEIGHTS);
 	    psSparseMatrixElement (sparse, j, i, f / r); 
 	}
@@ -194,6 +200,6 @@
 	    psAbort ("psphot", "ensemble source is nan");
 	}
-	Fi->modelPSF->params->data.F32[1] = norm->data.F32[i];
-	Fi->modelPSF->dparams->data.F32[1] = sqrt(sqrt(2) * norm->data.F32[i] / (sparse->Bfj->data.F32[i] * weight->data.F32[i]));
+	Fi->modelPSF->params->data.F32[PM_PAR_FLUX] = norm->data.F32[i];
+	Fi->modelPSF->dparams->data.F32[PM_PAR_FLUX] = sqrt(sqrt(2) * norm->data.F32[i] / (sparse->Bfj->data.F32[i] * weight->data.F32[i]));
 	// XXX EAM : this factor of sqrt(2) makes the errors consistent, but I don't understand it
 
@@ -210,6 +216,6 @@
 	pmModel *model = Fi->modelPSF;
 
-	x = model->params->data.F32[2];
-	y = model->params->data.F32[3];
+	x = model->params->data.F32[PM_PAR_XPOS];
+	y = model->params->data.F32[PM_PAR_YPOS];
 
 	psImageKeepCircle (Fi->mask, x, y, model->radiusTMP, "OR", PM_MASK_MARK);
