Index: /branches/eam_branches/20090715/psphot/src/psphotMakeResiduals.c
===================================================================
--- /branches/eam_branches/20090715/psphot/src/psphotMakeResiduals.c	(revision 25670)
+++ /branches/eam_branches/20090715/psphot/src/psphotMakeResiduals.c	(revision 25671)
@@ -236,4 +236,10 @@
 		}
 
+		float radius = hypot((ox - 0.5*resid->Ro->numCols), (oy - 0.5*resid->Ro->numRows));
+		if (radius > radiusMax) {
+                  resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = 1;
+		  continue;
+                }
+
                 resid->Ro->data.F32[oy][ox] = psStatsGetValue(fluxStats, statOption);
                 resid->Rx->data.F32[oy][ox] = resid->Ry->data.F32[oy][ox] = 0.0;
@@ -250,11 +256,13 @@
                   resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = 1;
                 }
+            } else {
+                assert (SPATIAL_ORDER == 1);
 
 		float radius = hypot((ox - 0.5*resid->Ro->numCols), (oy - 0.5*resid->Ro->numRows));
 		if (radius > radiusMax) {
                   resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = 1;
-                }
-            } else {
-                assert (SPATIAL_ORDER == 1);
+		  continue;
+                }
+
                 psImageInit(A, 0.0);
                 psVectorInit(B, 0.0);
@@ -279,8 +287,7 @@
 
                 if (!psMatrixGJSolve(A, B)) {
-                    psError(PSPHOT_ERR_PSF, false, "Singular matrix solving for (y,x) = (%d,%d)'s residuals",
-                            oy, ox);
-                    psFree(resid); resid = NULL;
-                    break;
+		    resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = 1;
+                    psWarning("Singular matrix solving for (y,x) = (%d,%d)'s residuals, masking", oy, ox);
+		    continue;
                 }