Index: branches/tap_branches/psphot/src/psphotMagnitudes.c
===================================================================
--- branches/tap_branches/psphot/src/psphotMagnitudes.c	(revision 25900)
+++ branches/tap_branches/psphot/src/psphotMagnitudes.c	(revision 27838)
@@ -1,14 +1,58 @@
 # include "psphotInternal.h"
 
-bool psphotMagnitudes(pmConfig *config, pmReadout *readout, const pmFPAview *view, psArray *sources, const pmPSF *psf) {
+bool psphotMagnitudes (pmConfig *config, const pmFPAview *view)
+{
+    bool status = true;
+
+    // select the appropriate recipe information
+    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
+    psAssert (recipe, "missing recipe?");
+
+    int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM");
+    psAssert (status, "programming error: must define PSPHOT.INPUT.NUM");
+
+    // skip the chisq image (optionally?)
+    int chisqNum = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.CHISQ.NUM");
+    if (!status) chisqNum = -1;
+
+    // loop over the available readouts
+    for (int i = 0; i < num; i++) {
+	if (i == chisqNum) continue; // skip chisq image
+
+	// find the currently selected readout
+	pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", i); // File of interest
+	psAssert (file, "missing file?");
+
+	pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
+	psAssert (readout, "missing readout?");
+
+	pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
+	psAssert (detections, "missing detections?");
+
+	psArray *sources = detections->allSources;
+	psAssert (sources, "missing sources?");
+
+	pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
+	psAssert (psf, "missing psf?");
+
+	if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf)) {
+            psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for PSPHOT.INPUT entry %d", i);
+	    return false;
+	}
+    }
+    return true;
+}
+
+bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf) {
 
     bool status = false;
     int Nap = 0;
 
+    if (!sources->n) {
+	psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping source magnitudes");
+	return true;
+    }
+	
     psTimerStart ("psphot.mags");
-
-    // select the appropriate recipe information
-    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
-    assert (recipe);
 
     // determine the number of allowed threads
@@ -41,10 +85,12 @@
     }
 
-    bool IGNORE_GROWTH = psMetadataLookupBool (&status, recipe, "IGNORE_GROWTH");
+    bool IGNORE_GROWTH  = psMetadataLookupBool (&status, recipe, "IGNORE_GROWTH");
     bool INTERPOLATE_AP = psMetadataLookupBool (&status, recipe, "INTERPOLATE_AP");
+    bool DIFF_STATS     = psMetadataLookupBool (&status, recipe, "INTERPOLATE_AP");
 
     pmSourcePhotometryMode photMode = PM_SOURCE_PHOT_APCORR | PM_SOURCE_PHOT_WEIGHT;
     if (!IGNORE_GROWTH) photMode |= PM_SOURCE_PHOT_GROWTH;
     if (INTERPOLATE_AP) photMode |= PM_SOURCE_PHOT_INTERP;
+    if (DIFF_STATS)     photMode |= PM_SOURCE_PHOT_DIFFSTATS;
 
     // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)
@@ -64,5 +110,5 @@
 
             psArrayAdd(job->args, 1, cells->data[j]); // sources
-            psArrayAdd(job->args, 1, (pmPSF*)psf);    // Casting away const
+            psArrayAdd(job->args, 1, psf);
             psArrayAdd(job->args, 1, binning);
             psArrayAdd(job->args, 1, backModel);
@@ -179,43 +225,4 @@
 }
 
-# if (0)
-bool psphotMagnitudes_Unthreaded (int *nap, psArray *sources, pmPSF *psf, psImageBinning *binning, pmReadout *backModel, pmReadout *backStdev, pmSourcePhotometryMode photMode, psImageMaskType maskVal) {
-
-    bool status;
-    int Nap = 0;
-
-    for (int i = 0; i < sources->n; i++) {
-        pmSource *source = (pmSource *) sources->data[i];
-        status = pmSourceMagnitudes (source, psf, photMode, maskVal);
-        if (status && isfinite(source->apMag)) Nap ++;
-
-        if (backModel) {
-            psAssert (binning, "if backModel is defined, so should binning be");
-            source->sky = psImageUnbinPixel(source->peak->x, source->peak->y, backModel->image, binning);
-            if (isnan(source->sky) && false) {
-                psLogMsg ("psphot.magnitudes", PS_LOG_DETAIL, "error setting pmSource.sky");
-            }
-        } else {
-            source->sky = NAN;
-        }
-
-        if (backStdev) {
-            psAssert (binning, "if backStdev is defined, so should binning be");
-            source->skyErr = psImageUnbinPixel(source->peak->x, source->peak->y, backStdev->image, binning);
-            if (isnan(source->skyErr) && false) {
-                psLogMsg ("psphot.magnitudes", PS_LOG_DETAIL, "error setting pmSource.skyErr");
-            }
-        } else {
-            source->skyErr = NAN;
-        }
-    }
-
-    // change the value of a scalar on the array (wrap this and put it in psArray.h)
-    *nap = Nap;
-
-    return true;
-}
-# endif
-
 bool psphotPSFWeights(pmConfig *config, pmReadout *readout, const pmFPAview *view, psArray *sources) {