Index: trunk/psModules/src/objects/pmPSFtry.c
===================================================================
--- trunk/psModules/src/objects/pmPSFtry.c	(revision 20937)
+++ trunk/psModules/src/objects/pmPSFtry.c	(revision 21173)
@@ -5,6 +5,6 @@
  *  @author EAM, IfA
  *
- *  @version $Revision: 1.67 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2008-12-08 02:51:14 $
+ *  @version $Revision: 1.68 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2009-01-27 00:01:33 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -175,9 +175,18 @@
 
         pmSource *source = psfTry->sources->data[i];
+	if (!source->moments) {
+            psfTry->mask->data.U8[i] = PSFTRY_MASK_EXT_FAIL;
+	    continue;
+	}
+	if (!source->moments->nPixels) {
+            psfTry->mask->data.U8[i] = PSFTRY_MASK_EXT_FAIL;
+	    continue;
+	}
+
         source->modelEXT = pmSourceModelGuess (source, psfTry->psf->type);
         if (source->modelEXT == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "failed to build model");
-            psFree(psfTry);
-            return NULL;
+            psfTry->mask->data.U8[i] = PSFTRY_MASK_EXT_FAIL;
+            psTrace ("psModules.objects", 4, "masking %d (%d,%d) : failed to generate model guess\n", i, source->peak->x, source->peak->y);
+            continue;
         }
 
@@ -591,9 +600,11 @@
     psVector *e2   = psVectorAlloc (sources->n, PS_TYPE_F32);
     psVector *mag  = psVectorAlloc (sources->n, PS_TYPE_F32);
+
     for (int i = 0; i < sources->n; i++) {
+	// skip any masked sources (failed to fit one of the model steps or get a magnitude)
+	if (srcMask->data.U8[i]) continue;
+	
         pmSource *source = sources->data[i];
-        if (source->modelEXT == NULL) continue;
-        // XXX I am relying on the fact that none of the masked sources
-        // have modelEXT set here.  perhaps use the value of psfTry->mask instead?
+        assert (source->modelEXT); // all unmasked sources should have modelEXT
 
         psEllipsePol pol = pmPSF_ModelToFit (source->modelEXT->params->data.F32);
@@ -619,5 +630,8 @@
         for (int i = 1; i <= PS_MAX (psf->trendNx, psf->trendNy); i++) {
 
-            psVectorInit (mask, 0);
+	    // copy srcMask to mask (we do not want the mask values set in pmPSFFitShapeParamsMap to be sticky)
+	    for (int i = 0; i < mask->n; i++) {
+		mask->data.U8[i] = srcMask->data.U8[i];
+	    }
             if (!pmPSFFitShapeParamsMap (psf, i, &scatterTotal, mask, x, y, mag, e0, e1, e2, dz)) {
                 break;
@@ -635,6 +649,8 @@
         }
 
-        // XXX supply the resulting mask values back to srcMask
-        psVectorInit (mask, 0);
+	// copy srcMask to mask (we do not want the mask values set in pmPSFFitShapeParamsMap to be sticky)
+	for (int i = 0; i < mask->n; i++) {
+	    mask->data.U8[i] = srcMask->data.U8[i];
+	}
         if (!pmPSFFitShapeParamsMap (psf, entryMin, &scatterTotal, mask, x, y, mag, e0, e1, e2, dz)) {
             psAbort ("failed pmPSFFitShapeParamsMap on second pass?");
@@ -836,5 +852,8 @@
     // the mask marks the values not used to calculate the ApTrend
     psVector *fitMask = psVectorAlloc (x_fit->n, PS_TYPE_U8);
-    psVectorInit (fitMask, 0);
+    // copy mask values to fitMask as a starting point
+    for (int i = 0; i < fitMask->n; i++) {
+	fitMask->data.U8[i] = mask->data.U8[i];
+    }
 
     // we run 'clipIter' cycles clipping in each of x and y, with only one iteration each.
@@ -884,7 +903,7 @@
 
     // construct the fitted values and the residuals
-    psVector *e0fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E0], x_tst, y_tst);
-    psVector *e1fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E1], x_tst, y_tst);
-    psVector *e2fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E2], x_tst, y_tst);
+    psVector *e0fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E0], fitMask, 0xff, x_tst, y_tst);
+    psVector *e1fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E1], fitMask, 0xff, x_tst, y_tst);
+    psVector *e2fit = pmTrend2DEvalVector (psf->params->data[PM_PAR_E2], fitMask, 0xff, x_tst, y_tst);
 
     psVector *e0res = (psVector *) psBinaryOp (NULL, (void *) e0obs_tst, "-", (void *) e0fit);
@@ -895,5 +914,5 @@
     // psTraceSetLevel ("psLib.math.vectorSampleStdev", 10);
     // psTraceSetLevel ("psLib.math.vectorClippedStats", 10);
-    pmPSFShapeParamsScatter (scatterTotal, e0res, e1res, e2res, psStatsStdevOption(psf->psfTrendStats->options));
+    pmPSFShapeParamsScatter (scatterTotal, e0res, e1res, e2res, fitMask, 0xff, psStatsStdevOption(psf->psfTrendStats->options));
     // psTraceSetLevel ("psLib.math.vectorSampleStdev", 0);
     // psTraceSetLevel ("psLib.math.vectorClippedStats", 0);
@@ -932,5 +951,5 @@
 
 // calculate the scatter of the parameters
-bool pmPSFShapeParamsScatter(float *scatterTotal, psVector *e0res, psVector *e1res, psVector *e2res, psStatsOptions stdevOpt)
+bool pmPSFShapeParamsScatter(float *scatterTotal, psVector *e0res, psVector *e1res, psVector *e2res, psVector *mask, psMaskType maskValue, psStatsOptions stdevOpt)
 {
 
@@ -941,13 +960,13 @@
     float dEsquare = 0.0;
     psStatsInit (stats);
-    psVectorStats (stats, e0res, NULL, NULL, 0xff);
+    psVectorStats (stats, e0res, NULL, mask, maskValue);
     dEsquare += PS_SQR(psStatsGetValue(stats, stdevOpt));
 
     psStatsInit (stats);
-    psVectorStats (stats, e1res, NULL, NULL, 0xff);
+    psVectorStats (stats, e1res, NULL, mask, maskValue);
     dEsquare += PS_SQR(psStatsGetValue(stats, stdevOpt));
 
     psStatsInit (stats);
-    psVectorStats (stats, e2res, NULL, NULL, 0xff);
+    psVectorStats (stats, e2res, NULL, mask, maskValue);
     dEsquare += PS_SQR(psStatsGetValue(stats, stdevOpt));