Index: trunk/ppStack/src/ppStackPrepare.c
===================================================================
--- trunk/ppStack/src/ppStackPrepare.c	(revision 30620)
+++ trunk/ppStack/src/ppStackPrepare.c	(revision 31158)
@@ -1,7 +1,9 @@
-#include "ppStack.h"
-
-#define PHOT_SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_SATSTAR | PM_SOURCE_MODE_BLEND | \
-                          PM_SOURCE_MODE_BADPSF | PM_SOURCE_MODE_DEFECT | PM_SOURCE_MODE_SATURATED | \
-                          PM_SOURCE_MODE_CR_LIMIT | PM_SOURCE_MODE_EXT_LIMIT) // Mask to apply to sources
+# include "ppStack.h"
+
+# define RE_PHOTOMETER 0
+
+# define PHOT_SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_SATSTAR | PM_SOURCE_MODE_BLEND | \
+			   PM_SOURCE_MODE_BADPSF | PM_SOURCE_MODE_DEFECT | PM_SOURCE_MODE_SATURATED | \
+			   PM_SOURCE_MODE_CR_LIMIT | PM_SOURCE_MODE_EXT_LIMIT) // Mask to apply to sources
 
 bool ppStackInputPhotometer(const pmReadout *ro, const psArray *sources, const pmConfig *config)
@@ -13,4 +15,5 @@
     psAssert(recipe, "We've thrown an error on this before.");
 
+# if (RE_PHOTOMETER)
     bool mdok;                          // Status of MD lookup
 
@@ -38,7 +41,6 @@
     psImageMaskType maskVal = pmConfigMaskGet(maskValStr, config); // Bits to mask
 
+    // Measure background
     psImage *image = ro->image, *mask = ro->mask; // Image and mask from readout
-
-    // Measure background
     psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator
     psStats *stats = psStatsAlloc(PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV); // Statistics
@@ -52,7 +54,7 @@
     psFree(stats);
     psFree(rng);
-
-    // Photometer sources
-    int numCols = image->numCols, numRows = image->numRows; // Size of image
+# endif
+
+    // Examine the sources (photometer if RE_PHOTOMETER is set)
     int numSources = sources->n; // Number of sources
     int numGood = 0;            // Number of good sources
@@ -64,5 +66,13 @@
             continue;
         }
-        float x = source->peak->xf, y = source->peak->yf; // Coordinates of source
+
+// This code re-measures the magnitude in the assumption that the warp analysis failed It would
+// be better to catch and correct such cases.  In any case, we should check on the validity of
+// this assumption (eg, compare the mean psf-ap offset) and raise an error on that?
+
+# if (RE_PHOTOMETER)
+	
+	int numCols = image->numCols, numRows = image->numRows; // Size of image
+	float x = source->peak->xf, y = source->peak->yf; // Coordinates of source
         int xCentral = x + 0.5, yCentral = y + 0.5; // Central pixel
         // Range of integration
@@ -96,4 +106,8 @@
             }
         }
+# else
+	numGood++;
+	maxMag = PS_MAX(maxMag, source->psfMag);
+# endif
     }
     psLogMsg("ppStack", PS_LOG_INFO, "Photometered %d good sources in image; max mag %f", numGood, maxMag);
@@ -217,4 +231,8 @@
         }
     }
+
+    bool mdok = false;
+    float maxFWHM = psMetadataLookupF32(&mdok, recipe, "PSF.INPUT.MAX"); // max allowed input fwhm
+    float clipFWHMnSig = psMetadataLookupF32(&mdok, recipe, "PSF.INPUT.CLIP.NSIGMA"); // sigma clipping of inputs
 
     psString log = psStringCopy("Input seeing FWHMs:\n"); // Log message
@@ -251,4 +269,10 @@
         }
 
+	// reject any input images which exceed the specified max FWHM
+	if (isfinite(maxFWHM) && (options->inputSeeing->data.F32[i] > maxFWHM)) {
+            options->inputMask->data.PS_TYPE_VECTOR_MASK_DATA[i] = PPSTACK_MASK_PSF;
+            psLogMsg("ppStack", PS_LOG_INFO, "PSF FWHM for image %d is too large (%f vs %f maxFWHM) --- rejected.", i, options->inputSeeing->data.F32[i], maxFWHM);
+	}
+
         psStringAppend(&log, "Input %d: %f\n", i, options->inputSeeing->data.F32[i]);
     }
@@ -257,4 +281,25 @@
     }
     psFree(log);
+
+    // We should have the ability to filter the input list based on the seeing:
+    // * reject above some max value and/or min value
+    // * reject images out-of-line with the rest of the images
+
+    // measure stats
+    psStats *fwhmStats = psStatsAlloc(PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
+    psVectorStats (fwhmStats, options->inputSeeing, NULL, options->inputMask, 0xff);
+    psLogMsg("ppStack", PS_LOG_INFO, "Input FWHMs : %f +/- %f", fwhmStats->clippedMean, fwhmStats->clippedStdev);
+
+    if (isfinite(clipFWHMnSig)) {
+	float fwhmLimit = fwhmStats->clippedMean + clipFWHMnSig * fwhmStats->clippedStdev;
+	fwhmLimit = isfinite(maxFWHM) ? PS_MIN (maxFWHM, fwhmLimit) : fwhmLimit;
+	for (int i = 0; i < num; i++) {
+            if (options->inputMask->data.PS_TYPE_VECTOR_MASK_DATA[i]) continue;
+            if (options->inputSeeing->data.F32[i] > fwhmLimit) {
+		options->inputMask->data.PS_TYPE_VECTOR_MASK_DATA[i] = PPSTACK_MASK_PSF;
+		psLogMsg("ppStack", PS_LOG_INFO, "PSF FWHM for image %d is too large (%f vs %f fwhmLimit) --- rejected.", i, options->inputSeeing->data.F32[i], fwhmLimit);
+	    }
+        }
+    }
 
     // Generate target PSF
