Index: trunk/psphot/src/psphotApResid.c
===================================================================
--- trunk/psphot/src/psphotApResid.c	(revision 23978)
+++ trunk/psphot/src/psphotApResid.c	(revision 23989)
@@ -355,12 +355,15 @@
 
         if (j > 2) {
-            bool status = true;
-            status &= psVectorStats (statsS, dASubset, NULL, mkSubset, 0xff);
-            status &= psVectorStats (statsM, dMSubset, NULL, mkSubset, 0xff);
-            if (!status) { psErrorClear (); }
+            if (!psVectorStats (statsS, dASubset, NULL, mkSubset, 0xff)) {
+		psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
+		return false;
+	    }
+            if (!psVectorStats (statsM, dMSubset, NULL, mkSubset, 0xff)) {
+		psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
+		return false;
+	    }
             dSo->data.F32[i] = statsS->robustStdev;
             dMo->data.F32[i] = statsM->sampleMean;
             dRo->data.F32[i] = statsS->robustStdev / statsM->sampleMean;
-            //      fprintf (stderr, "%d (%d) : sys: %f, phot: %f, rat: %f\n", i, j, dSo->data.F32[i], dMo->data.F32[i], dRo->data.F32[i]);
         } else {
             dSo->data.F32[i] = NAN;
@@ -375,6 +378,6 @@
     psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     if (!psVectorStats (stats, dRo, NULL, NULL, 0)) {
-        // XXX better testing of raised error
-        psErrorClear();
+	psError(PS_ERR_UNKNOWN, false, "failure to measure stats");
+	return false;
     }
 
Index: trunk/psphot/src/psphotChoosePSF.c
===================================================================
--- trunk/psphot/src/psphotChoosePSF.c	(revision 23978)
+++ trunk/psphot/src/psphotChoosePSF.c	(revision 23989)
@@ -361,5 +361,9 @@
 
     psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_STDEV | PS_STAT_SAMPLE_QUARTILE);
-    psVectorStats (stats, fwhmMajor, NULL, NULL, 0);
+    if (!psVectorStats (stats, fwhmMajor, NULL, NULL, 0)) {
+        psError(PS_ERR_UNKNOWN, false, "failure to measure stats for FWHM MAJOR");
+        return false;
+    }
+
     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MAJ",   PS_META_REPLACE, "PSF FWHM Major axis (mean)", stats->sampleMean);
     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_SG",   PS_META_REPLACE, "PSF FWHM Major axis (sigma)", stats->sampleStdev);
@@ -367,5 +371,8 @@
     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_UQ",   PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", stats->sampleUQ);
 
-    psVectorStats (stats, fwhmMinor, NULL, NULL, 0);
+    if (!psVectorStats (stats, fwhmMinor, NULL, NULL, 0)) {
+        psError(PS_ERR_UNKNOWN, false, "failure to measure stats for FWHM MINOR");
+        return false;
+    }
     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MIN",   PS_META_REPLACE, "PSF FWHM Minor axis (mean)", stats->sampleMean);
     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_SG",   PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", stats->sampleStdev);
Index: trunk/psphot/src/psphotDiagnosticPlots.c
===================================================================
--- trunk/psphot/src/psphotDiagnosticPlots.c	(revision 23978)
+++ trunk/psphot/src/psphotDiagnosticPlots.c	(revision 23989)
@@ -40,5 +40,8 @@
 
     psStats *stats = psStatsAlloc (PS_STAT_MAX | PS_STAT_MIN);
-    psVectorStats (stats, values, NULL, NULL, 0);
+    if (!psVectorStats (stats, values, NULL, NULL, 0)) {
+	psError(PS_ERR_UNKNOWN, false, "failure to measure stats for histogram");
+        return false;
+    }
 
     psHistogram *histogram = psHistogramAlloc (stats->min, stats->max, 1000);
@@ -63,5 +66,8 @@
     psStatsInit (stats);
     stats->options = PS_STAT_MAX | PS_STAT_MIN;
-    psVectorStats (stats, histogram->nums, NULL, NULL, 0);
+    if (!psVectorStats (stats, histogram->nums, NULL, NULL, 0)) {
+	psError(PS_ERR_UNKNOWN, false, "failure to measure stats for histogram");
+        return false;
+    }
 
     // scale the plot to hold the histogram
@@ -100,5 +106,8 @@
     psStatsInit (stats);
     stats->options = PS_STAT_MAX | PS_STAT_MIN;
-    psVectorStats (stats, histogram->nums, NULL, NULL, 0);
+    if (!psVectorStats (stats, histogram->nums, NULL, NULL, 0)) {
+	psError(PS_ERR_UNKNOWN, false, "failure to measure stats for histogram");
+        return false;
+    }
 
     // scale the plot to hold the histogram
Index: trunk/psphot/src/psphotImageQuality.c
===================================================================
--- trunk/psphot/src/psphotImageQuality.c	(revision 23978)
+++ trunk/psphot/src/psphotImageQuality.c	(revision 23989)
@@ -88,4 +88,5 @@
                      "Number of stars used for IQ measurements", M2->n);
 
+// XXX make this a recipe option
 #if (USE_SAMPLE)
     psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_STDEV | PS_STAT_SAMPLE_QUARTILE);
Index: trunk/psphot/src/psphotMakeResiduals.c
===================================================================
--- trunk/psphot/src/psphotMakeResiduals.c	(revision 23978)
+++ trunk/psphot/src/psphotMakeResiduals.c	(revision 23989)
@@ -200,6 +200,15 @@
             // measure the robust median to determine a baseline reference value
             *fluxClip = *fluxClipDef;
-            psVectorStats (fluxClip, fluxes, NULL, fmasks, fmaskVal);
-            psErrorClear();             // clear (ignore) any outstanding errors
+            if (!psVectorStats (fluxClip, fluxes, NULL, fmasks, fmaskVal)) {
+		psError(PSPHOT_ERR_CONFIG, false, "Error calculating residual stats");
+		return false;
+	    }
+	    if (isnan(fluxClip->robustMedian)) {
+                resid->Ro->data.F32[oy][ox] = 0.0;
+                resid->Rx->data.F32[oy][ox] = 0.0;
+                resid->Ry->data.F32[oy][ox] = 0.0;
+                resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = badMask;
+                continue;
+	    }
 
             // mark input pixels which are more than N sigma from the median
@@ -220,10 +229,19 @@
                 // measure the desired statistic on the unclipped pixels
                 *fluxStats = *fluxStatsDef;
-                psVectorStats (fluxStats, fluxes, NULL, fmasks, fmaskVal);
-                psErrorClear();         // clear (ignore) any outstanding errors
+                if (!psVectorStats (fluxStats, fluxes, NULL, fmasks, fmaskVal)) {
+		    psError(PSPHOT_ERR_CONFIG, false, "Error calculating residual stats");
+		    return false;
+		}
 
                 resid->Ro->data.F32[oy][ox] = psStatsGetValue(fluxStats, statOption);
                 resid->Rx->data.F32[oy][ox] = resid->Ry->data.F32[oy][ox] = 0.0;
-                //resid->variance->data.F32[oy][ox] = fluxStats->sampleStdev;
+
+		if (isnan(resid->Ro->data.F32[oy][ox])) {
+		    resid->Ro->data.F32[oy][ox] = 0.0;
+		    resid->Rx->data.F32[oy][ox] = 0.0;
+		    resid->Ry->data.F32[oy][ox] = 0.0;
+		    resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = badMask;
+		    continue;
+		}
 
                 if (fabs(resid->Ro->data.F32[oy][ox]) < pixelSN*fluxStats->sampleStdev/sqrt(nKeep)) {
Index: trunk/psphot/src/psphotRoughClass.c
===================================================================
--- trunk/psphot/src/psphotRoughClass.c	(revision 23978)
+++ trunk/psphot/src/psphotRoughClass.c	(revision 23989)
@@ -91,5 +91,5 @@
 	return false;
     }
-    if (!psfClump.X || !psfClump.Y) {
+    if (!psfClump.X || !psfClump.Y || isnan(psfClump.X) || isnan(psfClump.Y)) {
 	psLogMsg ("psphot", 4, "Failed to find a valid PSF clump for region %f,%f - %f,%f\n", region->x0, region->y0, region->x1, region->y1);
 	return false;