Index: /branches/eam_rel9_b0/psModules/src/objects/pmObjects.c
===================================================================
--- /branches/eam_rel9_b0/psModules/src/objects/pmObjects.c	(revision 5726)
+++ /branches/eam_rel9_b0/psModules/src/objects/pmObjects.c	(revision 5727)
@@ -6,6 +6,6 @@
  *  @author EAM, IfA: significant modifications.
  *
- *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-10-20 23:06:24 $
+ *  @version $Revision: 1.3.4.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-07 20:01:05 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -81,4 +81,7 @@
     tmp->chisq = 0.0;
     tmp->nIter = 0;
+    tmp->radius = 0;
+    tmp->status = PM_MODEL_UNTRIED;
+
     psS32 Nparams = pmModelParameterCount(type);
     if (Nparams == 0) {
@@ -769,6 +772,7 @@
     for (psS32 row = 0; row < source->pixels->numRows ; row++) {
         for (psS32 col = 0; col < source->pixels->numCols ; col++) {
-            if ((source->mask != NULL) && (source->mask->data.U8[row][col]))
+            if ((source->mask != NULL) && (source->mask->data.U8[row][col])) {
                 continue;
+            }
 
             psF32 xDiff = col + source->pixels->col0 - xPeak;
@@ -777,6 +781,7 @@
             // XXX EAM : calculate xDiff, yDiff up front;
             //           radius is just a function of (xDiff, yDiff)
-            if (!VALID_RADIUS(xDiff, yDiff, R2))
+            if (!VALID_RADIUS(xDiff, yDiff, R2)) {
                 continue;
+            }
 
             psF32 pDiff = source->pixels->data.F32[row][col] - sky;
@@ -784,6 +789,7 @@
             // XXX EAM : check for valid S/N in pixel
             // XXX EAM : should this limit be user-defined?
-            if (pDiff / sqrt(source->weight->data.F32[row][col]) < 1)
+            if (pDiff / sqrt(source->weight->data.F32[row][col]) < 1) {
                 continue;
+            }
 
             Sum += pDiff;
@@ -913,8 +919,16 @@
         psImage *splane = NULL;
         int binX, binY;
+        bool status;
+
+        psF32 SX_MAX = psMetadataLookupF32 (&status, metadata, "MOMENTS_SX_MAX");
+        if (!status)
+            SX_MAX = 10.0;
+        psF32 SY_MAX = psMetadataLookupF32 (&status, metadata, "MOMENTS_SY_MAX");
+        if (!status)
+            SY_MAX = 10.0;
 
         // construct a sigma-plane image
         // psImageAlloc does zero the data
-        splane = psImageAlloc (NPIX/SCALE, NPIX/SCALE, PS_TYPE_F32);
+        splane = psImageAlloc (SX_MAX/SCALE, SY_MAX/SCALE, PS_TYPE_F32);
         for (int i = 0; i < splane->numRows; i++)
         {
@@ -1065,4 +1079,5 @@
     int Ncr      = 0;
     int Nsatstar = 0;
+    psRegion allArray = psRegionSet (0, 0, 0, 0);
 
     psVector *starsn = psVectorAlloc (sources->n, PS_TYPE_F32);
@@ -1096,10 +1111,5 @@
 
         // XXX EAM : can we use the value of SATURATE if mask is NULL?
-        //
-        // XXX: Must verify this region (the region argument was added to psImageCountPixelMask()
-        // after EAM wrote this code.
-        //
-        psRegion tmpRegion = psRegionSet(0, tmpSrc->mask->numCols, 0, tmpSrc->mask->numRows);
-        int Nsatpix = psImageCountPixelMask(tmpSrc->mask, tmpRegion, PSPHOT_MASK_SATURATED);
+        int Nsatpix = psImageCountPixelMask (tmpSrc->mask, allArray, PSPHOT_MASK_SATURATED);
 
         // saturated star (size consistent with PSF or larger)
@@ -1107,5 +1117,5 @@
         bool big = (sigX > (clump.X - clump.dX)) && (sigY > (clump.Y - clump.dY));
         if ((Nsatpix > 1) && big) {
-            tmpSrc->type |= PM_SOURCE_SATSTAR;
+            tmpSrc->type = PM_SOURCE_SATSTAR;
             Nsatstar ++;
             continue;
@@ -1114,5 +1124,5 @@
         // saturated object (not a star, eg bleed trails, hot pixels)
         if (Nsatpix > 1) {
-            tmpSrc->type |= PM_SOURCE_SATURATED;
+            tmpSrc->type = PM_SOURCE_SATURATED;
             Nsat ++;
             continue;
@@ -1123,5 +1133,5 @@
         // only set candidate defects if
         if ((sigX < 0.05) || (sigY < 0.05)) {
-            tmpSrc->type |= PM_SOURCE_DEFECT;
+            tmpSrc->type = PM_SOURCE_DEFECT;
             Ncr ++;
             continue;
@@ -1130,5 +1140,5 @@
         // likely unsaturated galaxy (too large to be stellar)
         if ((sigX > (clump.X + 3*clump.dX)) || (sigY > (clump.Y + 3*clump.dY))) {
-            tmpSrc->type |= PM_SOURCE_GALAXY;
+            tmpSrc->type = PM_SOURCE_GALAXY;
             Ngal ++;
             continue;
@@ -1143,5 +1153,5 @@
         psF32 radius = hypot ((sigX-clump.X)/clump.dX, (sigY-clump.Y)/clump.dY);
         if ((SN > PSF_SN_LIM) && (radius < 1.5)) {
-            tmpSrc->type |= PM_SOURCE_PSFSTAR;
+            tmpSrc->type = PM_SOURCE_PSFSTAR;
             Npsf ++;
             continue;
@@ -1149,5 +1159,5 @@
 
         // random type of star
-        tmpSrc->type |= PM_SOURCE_OTHER;
+        tmpSrc->type = PM_SOURCE_OTHER;
     }
 
@@ -1603,4 +1613,5 @@
         psTrace (".pmObjects.pmSourceFitModel", 4, "insufficient valid pixels\n");
         psTrace(__func__, 3, "---- %s(false) end ----\n", __func__);
+        model->status = PM_MODEL_BADARGS;
         return(false);
     }
@@ -1654,19 +1665,10 @@
     }
 
-    // XXX EAM: we need to do something (give an error?) if rc is false
-    // XXX EAM: psMinimizeLMChi2 does not check convergence
-
-    // XXX models can go insane: reject these
-    onPic &= (params->data.F32[2] >= source->pixels->col0);
-    onPic &= (params->data.F32[2] <  source->pixels->col0 + source->pixels->numCols);
-    onPic &= (params->data.F32[3] >= source->pixels->row0);
-    onPic &= (params->data.F32[3] <  source->pixels->row0 + source->pixels->numRows);
-
-    // XXX EAM: save the resulting chisq, nDOF, nIter
+    // save the resulting chisq, nDOF, nIter
     model->chisq = myMin->value;
     model->nIter = myMin->iter;
     model->nDOF  = y->n - nParams;
 
-    // XXX EAM get the Gauss-Newton distance for fixed model parameters
+    // get the Gauss-Newton distance for fixed model parameters
     if (paramMask != NULL) {
         psVector *delta = psVectorAlloc (params->n, PS_TYPE_F64);
@@ -1678,4 +1680,20 @@
         }
         psFree (delta);
+    }
+
+    // set the model success or failure status
+    if (!fitStatus) {
+        model->status = PM_MODEL_NONCONVERGE;
+    } else {
+        model->status = PM_MODEL_SUCCESS;
+    }
+
+    // models can go insane: reject these
+    onPic &= (params->data.F32[2] >= source->pixels->col0);
+    onPic &= (params->data.F32[2] <  source->pixels->col0 + source->pixels->numCols);
+    onPic &= (params->data.F32[3] >= source->pixels->row0);
+    onPic &= (params->data.F32[3] <  source->pixels->row0 + source->pixels->numRows);
+    if (!onPic) {
+        model->status = PM_MODEL_OFFIMAGE;
     }
 
@@ -1737,4 +1755,5 @@
         psTrace (".pmObjects.pmSourceFitModel", 4, "insufficient valid pixels\n");
         psTrace(__func__, 3, "---- %s(false) end ----\n", __func__);
+        model->status = PM_MODEL_BADARGS;
         return(false);
     }
@@ -1815,10 +1834,4 @@
     // XXX EAM: psMinimizeLMChi2 does not check convergence
 
-    // XXX models can go insane: reject these
-    onPic &= (params->data.F32[2] >= source->pixels->col0);
-    onPic &= (params->data.F32[2] <  source->pixels->col0 + source->pixels->numCols);
-    onPic &= (params->data.F32[3] >= source->pixels->row0);
-    onPic &= (params->data.F32[3] <  source->pixels->row0 + source->pixels->numRows);
-
     // XXX EAM: save the resulting chisq, nDOF, nIter
     model->chisq = myMin->value;
@@ -1837,8 +1850,26 @@
     }
 
-    psFree(paramMask);
+    // set the model success or failure status
+    if (!fitStatus) {
+        model->status = PM_MODEL_NONCONVERGE;
+    } else {
+        model->status = PM_MODEL_SUCCESS;
+    }
+
+    // XXX models can go insane: reject these
+    onPic &= (params->data.F32[2] >= source->pixels->col0);
+    onPic &= (params->data.F32[2] <  source->pixels->col0 + source->pixels->numCols);
+    onPic &= (params->data.F32[3] >= source->pixels->row0);
+    onPic &= (params->data.F32[3] <  source->pixels->row0 + source->pixels->numRows);
+    if (!onPic) {
+        model->status = PM_MODEL_OFFIMAGE;
+    }
+
     psFree(x);
     psFree(y);
+    psFree(yErr);
     psFree(myMin);
+    psFree(covar);
+    psFree(paramMask);
 
     rc = (onPic && fitStatus);
@@ -1851,5 +1882,7 @@
                              pmModel *model,
                              bool center,
-                             psS32 flag)
+                             bool sky,
+                             bool add
+                                )
 {
     psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
@@ -1864,4 +1897,6 @@
     psS32 imageCol;
     psS32 imageRow;
+    psF32 skyValue = params->data.F32[0];
+    psF32 pixelValue;
 
     for (psS32 i = 0; i < image->numRows; i++) {
@@ -1869,5 +1904,5 @@
             if ((mask != NULL) && mask->data.U8[i][j])
                 continue;
-            psF32 pixelValue;
+
             // XXX: Should you be adding the pixels for the entire subImage,
             // or a radius of pixels around it?
@@ -1887,15 +1922,21 @@
             x->data.F32[0] = (float) imageCol;
             x->data.F32[1] = (float) imageRow;
-            pixelValue = modelFunc (NULL, params, x);
-            // fprintf (stderr, "%f %f  %d %d  %f\n", x->data.F32[0], x->data.F32[1], i, j, pixelValue);
-
-            if (flag == 1) {
-                pixelValue = -pixelValue;
-            }
-
-            // XXX: Must figure out how to calculate the image coordinates and
-            // how to use the boolean "center" flag.
-
-            image->data.F32[i][j]+= pixelValue;
+
+            // set the appropriate pixel value for this coordinate
+            if (sky) {
+                pixelValue = modelFunc (NULL, params, x);
+            } else {
+                pixelValue = modelFunc (NULL, params, x) - skyValue;
+            }
+
+
+            // add or subtract the value
+            if (add
+               ) {
+                image->data.F32[i][j] += pixelValue;
+            }
+            else {
+                image->data.F32[i][j] -= pixelValue;
+            }
         }
     }
@@ -1912,8 +1953,9 @@
                       psImage *mask,
                       pmModel *model,
-                      bool center)
-{
-    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
-    psBool rc = p_pmSourceAddOrSubModel(image, mask, model, center, 0);
+                      bool center,
+                      bool sky)
+{
+    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
+    psBool rc = p_pmSourceAddOrSubModel(image, mask, model, center, sky, true);
     psTrace(__func__, 3, "---- %s(%d) end ----\n", __func__, rc);
     return(rc);
@@ -1925,8 +1967,9 @@
                       psImage *mask,
                       pmModel *model,
-                      bool center)
-{
-    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
-    psBool rc = p_pmSourceAddOrSubModel(image, mask, model, center, 1);
+                      bool center,
+                      bool sky)
+{
+    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
+    psBool rc = p_pmSourceAddOrSubModel(image, mask, model, center, sky, false);
     psTrace(__func__, 3, "---- %s(%d) end ----\n", __func__, rc);
     return(rc);
