Index: /branches/eam_branches/stackphot.20100406/psphot/src/Makefile.am
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/Makefile.am	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/Makefile.am	(revision 27625)
@@ -17,13 +17,13 @@
 # Force recompilation of psphotVersion.c, since it gets the version information
 psphotVersion.c: psphotVersionDefinitions.h
-psphotVersionDefinitions.h: psphotVersionDefinitions.h.in FORCE
-	-$(RM) psphotVersionDefinitions.h
-	$(SED) -e "s|@PSPHOT_VERSION@|\"$(PSPHOT_VERSION)\"|" -e "s|@PSPHOT_BRANCH@|\"$(PSPHOT_BRANCH)\"|" -e "s|@PSPHOT_SOURCE@|\"$(PSPHOT_SOURCE)\"|" psphotVersionDefinitions.h.in > psphotVersionDefinitions.h
-FORCE: ;
+# psphotVersionDefinitions.h: psphotVersionDefinitions.h.in FORCE
+# 	-$(RM) psphotVersionDefinitions.h
+# 	$(SED) -e "s|@PSPHOT_VERSION@|\"$(PSPHOT_VERSION)\"|" -e "s|@PSPHOT_BRANCH@|\"$(PSPHOT_BRANCH)\"|" -e "s|@PSPHOT_SOURCE@|\"$(PSPHOT_SOURCE)\"|" psphotVersionDefinitions.h.in > psphotVersionDefinitions.h
+# FORCE: ;
 
 libpsphot_la_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
 libpsphot_la_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
 
-bin_PROGRAMS = psphot psphotForced psphotMakePSF
+bin_PROGRAMS = psphot psphotForced psphotMakePSF psphotStack
 # bin_PROGRAMS = psphotPetrosianStudy psphotTest psphotMomentsStudy 
 
@@ -39,4 +39,8 @@
 psphotMakePSF_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
 psphotMakePSF_LDADD = libpsphot.la
+
+psphotStack_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
+psphotStack_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
+psphotStack_LDADD = libpsphot.la
 
 # psphotMomentsStudy_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
@@ -80,4 +84,17 @@
 	psphotMosaicChip.c	   \
 	psphotCleanup.c
+
+# a psphot-variant for stack photometry
+psphotStack_SOURCES = \
+        psphotStack.c              \
+	psphotStackArguments.c     \
+	psphotStackParseCamera.c   \
+	psphotStackImageLoop.c     \
+	psphotStackReadout.c	   \
+	psphotStackChisqImage.c	   \
+	psphotCleanup.c
+
+# psphotFitSourceLinearStack.c	  
+
 
 # # psphot analysis of the detectability of specified positions
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphot.h
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphot.h	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphot.h	(revision 27625)
@@ -341,3 +341,17 @@
 bool psphotMakePSFReadout(pmConfig *config, const pmFPAview *view);
 
+/**** psphotStack prototypes ****/
+
+pmConfig *psphotStackArguments(int argc, char **argv);
+bool psphotStackParseCamera (pmConfig *config);
+bool psphotStackImageLoop (pmConfig *config);
+bool psphotStackReadout (pmConfig *config, const pmFPAview *view);
+bool psphotStackChisqImage (pmConfig *config, const pmFPAview *view);
+bool psphotStackChisqImageAddReadout(const pmConfig *config, // Configuration
+				     const pmFPAview *view,
+				     pmReadout *chiReadout,
+				     char *filename, 
+				     int index);
+
+
 #endif
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotFitSourcesLinearStack.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotFitSourcesLinearStack.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotFitSourcesLinearStack.c	(revision 27625)
@@ -1,23 +1,7 @@
 # include "psphotInternal.h"
-
-// for now, let's store the detections on the readout->analysis for each readout
-bool psphotFitSourcesLinearStack (pmConfig *config, const pmFPAview *view, bool final)
-{
-    bool status = true;
-
-    int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM");
-    psAssert (status, "programming error: must define PSPHOT.INPUT.NUM");
-
-    // loop over the available readouts
-    for (int i = 0; i < num; i++) {
-	if (!psphotFitSourcesLinearReadoutStack (config, view, "PSPHOT.INPUT", i, final)) {
-            psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (linear) for PSPHOT.INPUT entry %d", i);
-	    return false;
-	}
-    }
-    return true;
-}
-
-bool psphotFitSourcesLinearReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, bool final) {
+// XXX need array of covar factors for each image
+// XXX define the 'good' / 'bad' flags?
+
+bool psphotFitSourcesLinearStack (pmConfig *config, psArray *objects, bool final) {
 
     bool status;
@@ -25,5 +9,4 @@
     float y;
     float f;
-    // float r;
 
     // select the appropriate recipe information
@@ -31,25 +14,4 @@
     assert (recipe);
 
-    // find the currently selected readout
-    pmFPAfile *file = pmFPAfileSelectSingle(config->files, filename, index); // 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?");
-
-    if (!sources->n) {
-	psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping linear fit");
-	return true;
-    }
-
-    pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
-    psAssert (sources, "missing psf?");
-
     psTimerStart ("psphot.linear");
 
@@ -65,36 +27,17 @@
     maskVal |= markVal;
 
-    // source analysis is done in spatial order
-    sources = psArraySort (sources, pmSourceSortByY);
+    // analysis is done in spatial order (to speed up overlap search)
+    // sort by first element in each source list
+    objects = psArraySort (objects, pmPhotObjSortByY);
 
     // storage array for fitSources
-    psArray *fitSources = psArrayAllocEmpty (sources->n);
-
-    // option to limit analysis to a specific region
-    char *region = psMetadataLookupStr (&status, recipe, "ANALYSIS_REGION");
-    psRegion AnalysisRegion = psRegionFromString (region);
-    AnalysisRegion = psRegionForImage (readout->image, AnalysisRegion);
-    if (psRegionIsNaN (AnalysisRegion)) psAbort("analysis region mis-defined");
-
-    bool CONSTANT_PHOTOMETRIC_WEIGHTS =
-        psMetadataLookupBool(&status, recipe, "CONSTANT_PHOTOMETRIC_WEIGHTS");
-    if (!status) {
-        psAbort("You must provide a value for the BOOL recipe CONSTANT_PHOTOMETRIC_WEIGHTS");
-    }
-    int SKY_FIT_ORDER = psMetadataLookupS32(&status, recipe, "SKY_FIT_ORDER");
-    if (!status) {
-        SKY_FIT_ORDER = 0;
-    }
-    bool SKY_FIT_LINEAR = psMetadataLookupBool(&status, recipe, "SKY_FIT_LINEAR");
-    if (!status) {
-        SKY_FIT_LINEAR = false;
-    }
-
-    // XXX test: choose a larger-than expected radius:
-    float covarFactor = psImageCovarianceFactorForAperture(readout->covariance, 10.0); // Covariance matrix
-    psLogMsg ("psphot.ensemble", PS_LOG_MINUTIA, "covariance factor: %f\n", covarFactor);
-
-    // XXX do not apply covarFactor for the moment...
-    // covarFactor = 1.0;
+    psArray *fitSources = psArrayAllocEmpty (objects->n);
+
+    bool CONSTANT_PHOTOMETRIC_WEIGHTS = psMetadataLookupBool(&status, recipe, "CONSTANT_PHOTOMETRIC_WEIGHTS");
+    psAssert (status, "You must provide a value for the BOOL recipe CONSTANT_PHOTOMETRIC_WEIGHTS");
+
+    // XXX store a local static array of covar factors for each of the images (by image ID)
+    // float covarFactor = psImageCovarianceFactorForAperture(readout->covariance, 10.0); // Covariance matrix
+    // psLogMsg ("psphot.ensemble", PS_LOG_MINUTIA, "covariance factor: %f\n", covarFactor);
 
     // select the sources which will be used for the fitting analysis
@@ -104,6 +47,6 @@
 	if (!object->sources) continue;
 
-	// check an element of the group to see if we should use it
-	if (!object->flags & PM_PHOT_OBJ_BAD) continue;
+	// XXX check an element of the group to see if we should use it
+	// if (!object->flags & PM_PHOT_OBJ_BAD) continue;
 
 	for (int j = 0; j < object->sources->n; j++) {
@@ -162,6 +105,6 @@
             pmSource *SRCj = fitSources->data[j];
 
-	    // XXX I need to know if this source is on the same image as SRCi --
-	    if (!sameImge) { continue; }
+	    // we only need to generate dot terms for source on the same image
+	    if (SRCj->imageID != SRCi->imageID) { continue; }
 
             // skip over disjoint source images, break after last possible overlap
@@ -190,10 +133,4 @@
     psLogMsg ("psphot.ensemble", PS_LOG_MINUTIA, "solve matrix: %f sec (%d elements)\n", psTimerMark ("psphot.linear"), sparse->Nelem);
 
-    // XXXX **** philosophical question:
-    // we measure bright objects in three passes: 1) linear fit; 2) non-linear fit; 3) linear fit:
-    // should retain the chisq and errors from the intermediate non-linear fit?
-    // the non-linear fit provides better values for the position errors, and for
-    // extended sources, the shape errors
-
     // adjust I0 for fitSources and subtract
     for (int i = 0; i < fitSources->n; i++) {
@@ -208,5 +145,4 @@
         model->params->data.F32[PM_PAR_I0] = norm->data.F32[i];
         model->dparams->data.F32[PM_PAR_I0] = errors->data.F32[i];
-        // XXX is the value of 'errors' modified by the sky fit?
 
         // subtract object
@@ -229,100 +165,36 @@
     psFree (fitSources);
     psFree (norm);
-    psFree (skyfit);
     psFree (errors);
-    psFree (border);
 
     psLogMsg ("psphot.ensemble", PS_LOG_INFO, "measure ensemble of PSFs: %f sec\n", psTimerMark ("psphot.linear"));
-
-    psphotVisualShowResidualImage (readout);
-    psphotVisualPlotChisq (sources);
-    // psphotVisualShowFlags (sources);
-
-    // We have to place this visualization here because the models are not realized until
-    // psphotGuessModels or fitted until psphotFitSourcesLinear.
-    psphotVisualShowPSFStars (recipe, psf, sources);
 
     return true;
 }
 
-// XXX do we need this?
-// XXX disallow the simultaneous sky fit and remove this code...
-
-// Calculate the weight terms for the sky fit component of the matrix.  This function operates
-// on the pixels which correspond to all of the sources of interest.  These elements fill in
-// the border matrix components in the sparse matrix equation.
-static bool SetBorderMatrixElements (psSparseBorder *border, pmReadout *readout, psArray *sources, bool constant_weights, int SKY_FIT_ORDER, psImageMaskType markVal) {
-
-    // generate the image-wide weight terms
-    // turn on MARK for all image pixels
-    psRegion fullArray = psRegionSet (0, 0, 0, 0);
-    fullArray = psRegionForImage (readout->mask, fullArray);
-    psImageMaskRegion (readout->mask, fullArray, "OR", markVal);
-
-    // turn off MARK for all object pixels
-    for (int i = 0; i < sources->n; i++) {
-        pmSource *source = sources->data[i];
-        pmModel *model = pmSourceGetModel (NULL, source);
-        if (model == NULL) continue;
-        float x = model->params->data.F32[PM_PAR_XPOS];
-        float y = model->params->data.F32[PM_PAR_YPOS];
-        psImageMaskCircle (source->maskView, x, y, model->fitRadius, "AND", PS_NOT_IMAGE_MASK(markVal));
-    }
-
-    // accumulate the image statistics from the masked regions
-    psF32 **image  = readout->image->data.F32;
-    psF32 **variance = readout->variance->data.F32;
-    psImageMaskType  **mask   = readout->mask->data.PS_TYPE_IMAGE_MASK_DATA;
-
-    double w, x, y, x2, xy, y2, xc, yc, wt, f, fo, fx, fy;
-    w = x = y = x2 = xy = y2 = fo = fx = fy = 0;
-
-    int col0 = readout->image->col0;
-    int row0 = readout->image->row0;
-
-    for (int j = 0; j < readout->image->numRows; j++) {
-        for (int i = 0; i < readout->image->numCols; i++) {
-            if (mask[j][i]) continue;
-            if (constant_weights) {
-                wt = 1.0;
-            } else {
-                wt = variance[j][i];
-            }
-            f = image[j][i];
-            w   += 1/wt;
-            fo  += f/wt;
-            if (SKY_FIT_ORDER == 0) continue;
-
-            xc  = i + col0;
-            yc  = j + row0;
-            x  +=    xc/wt;
-            y  +=    yc/wt;
-            x2 += xc*xc/wt;
-            xy += xc*yc/wt;
-            y2 += yc*yc/wt;
-            fx +=  f*xc/wt;
-            fy +=  f*yc/wt;
-        }
-    }
-
-    // turn off MARK for all image pixels
-    psImageMaskRegion (readout->mask, fullArray, "AND", PS_NOT_IMAGE_MASK(markVal));
-
-    // set the Border T elements
-    psSparseBorderElementG (border, 0, fo);
-    psSparseBorderElementT (border, 0, 0, w);
-    if (SKY_FIT_ORDER > 0) {
-        psSparseBorderElementG (border, 0, fx);
-        psSparseBorderElementG (border, 0, fy);
-        psSparseBorderElementT (border, 1, 0, x);
-        psSparseBorderElementT (border, 2, 0, y);
-        psSparseBorderElementT (border, 0, 1, x);
-        psSparseBorderElementT (border, 1, 1, x2);
-        psSparseBorderElementT (border, 2, 1, xy);
-        psSparseBorderElementT (border, 0, 2, y);
-        psSparseBorderElementT (border, 1, 2, xy);
-        psSparseBorderElementT (border, 2, 2, y2);
-    }
-
-    return true;
+// sort by Y (ascending)
+int pmPhotObjSortBySN (const void **a, const void **b)
+{
+    pmPhotObj *objA = *(pmPhotObj **)a;
+    pmPhotObj *objB = *(pmPhotObj **)b;
+
+    psAssert (objA->sources, "missing sources?");
+    psAssert (objB->sources, "missing sources?");
+
+    psAssert (objA->sources->n, "missing sources");
+    psAssert (objB->sources->n, "missing sources");
+
+    psAssert (objA->sources->data[0], "missing sources");
+    psAssert (objB->sources->data[0], "missing sources");
+
+    pmSource *A = objA->sources->data[0];
+    pmSource *B = objB->sources->data[0];
+
+    psF32 fA = (A->peak == NULL) ? 0 : A->peak->y;
+    psF32 fB = (B->peak == NULL) ? 0 : B->peak->y;
+
+    psF32 diff = fA - fB;
+    if (diff > FLT_EPSILON) return (+1);
+    if (diff < FLT_EPSILON) return (-1);
+    return (0);
 }
+
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotReadout.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotReadout.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotReadout.c	(revision 27625)
@@ -34,5 +34,7 @@
 
     // Generate the mask and weight images, including the user-defined analysis region of interest
-    psphotSetMaskAndVariance (config, view);
+    if (!psphotSetMaskAndVariance (config, view)) {
+        return psphotReadoutCleanup(config, view);
+    }
     if (!strcasecmp (breakPt, "NOTHING")) {
         return psphotReadoutCleanup(config, view);
Index: anches/eam_branches/stackphot.20100406/psphot/src/psphotSourceMerge.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotSourceMerge.c	(revision 27624)
+++ 	(revision )
@@ -1,49 +1,0 @@
-
-// we need a couple of functions to distinguish coincident sources:
-// XXX these are similar (identical?) to the goals of pmSourceMatch.c
-
-# define NEXT1 { i++; continue; }
-# define NEXT2 { j++; continue; }
-
-psphotSourceMerge () {
-
-    float dx, dy;
-
-    // sort the source list by X
-    pmSourceSortByX (sources1);
-    pmSourceSortByX (sources2);
-
-    int i, j;
-    for (i = j = 0; (i < sources1->n) && (j < sources2->n); ) {
-
-	pmSource *src1 = sources1->data[i];
-	pmSource *src2 = sources2->data[j];
-
-	if (!src1) NEXT1;
-	if (!src1->peak) NEXT1;
-	if (!finite(src1->peak->xf)) NEXT1;
-	if (!finite(src1->peak->yf)) NEXT1;
-
-	if (!src2) NEXT2;
-	if (!src2->peak) NEXT2;
-	if (!finite(src2->peak->xf)) NEXT2;
-	if (!finite(src2->peak->yf)) NEXT2;
-
-	dx = src1->peak->xf - src2->peak->xf;
-	if (dx < -1.02*RADIUS) NEXT1;
-	if (dx > +1.02*RADIUS) NEXT2;
-
-	// we are within match range, look for matches:
-	for (int J = j; (dx > -1.02*radius) && (J < sources2->n); J++) {
-
-	    dx = src1->peak->xf - src2->peak->xf;
-	    dy = src1->peak->yf - src2->peak->yf;
-
-	    dr = dx*dx + dy*dy;
-	    if (dr > RADIUS2) continue;
-
-	    // add to group?
-	}
-	i++;
-    }
-}
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackArguments.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackArguments.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackArguments.c	(revision 27625)
@@ -2,14 +2,13 @@
 
 static void usage(pmConfig *config, int exitCode);
-static void writeHelpInfo(pmConfig* config, FILE* ofile);
+static void writeHelpInfo(FILE* ofile);
 
 pmConfig *psphotStackArguments(int argc, char **argv) {
 
     int N;
-    bool status;
 
     // print help info
-    if (psArgumentGet(argc, argv, "-help")) writeHelpInfo(argv[0], config, stdout);
-    if (psArgumentGet(argc, argv, "-h")) writeHelpInfo(argv[0], config, stdout);
+    if (psArgumentGet(argc, argv, "-help")) writeHelpInfo(stdout);
+    if (psArgumentGet(argc, argv, "-h")) writeHelpInfo(stdout);
 
     // load config data from default locations
@@ -29,15 +28,4 @@
     // Number of threads is handled
     PSARGUMENTS_INSTANTIATE_THREADSARG( psphot, config, argc, argv )
-
-    // photcode : used in output to supplement header data (argument or recipe?)
-    if ((N = psArgumentGet (argc, argv, "-photcode"))) {
-        if (argc <= N+1) {
-	  psErrorStackPrint(stderr, "Expected to see 1 more argument; saw %d", argc - 1);
-	  exit(PS_EXIT_CONFIG_ERROR);
-	}
-        psArgumentRemove (N, &argc, argv);
-        psMetadataAddStr (options, PS_LIST_TAIL, "PHOTCODE", PS_META_REPLACE, "", argv[N]);
-        psArgumentRemove (N, &argc, argv);
-    }
 
     // visual : interactive display mode
@@ -90,14 +78,14 @@
 }
 
-static void writeHelpInfo(pmConfig* config, FILE* ofile)
+static void writeHelpInfo(FILE* ofile)
 {
   fprintf(ofile,
-	  "Usage: psphotStack -input (input.mdc) outroot\n"
+	  "Usage: psphotStack -input (INPUTS.mdc) (OUTROOT)\n"
 	  "\n"
-	  "where:\n"
-	  "  FileNameList is a text file containing filenames, one per line\n"
-	  "  MaskFileNameList is a text file of mask filenames, one per line\n"
-	  "  VarFileNameList is a text file of variance filenames, one per line\n"
-	  "  OutFileBaseName is the 'root name' for output files\n"
+	  "where INPUTS.mdc contains various METADATAs, each with:\n"
+	  "\tIMAGE(STR):     Image filename\n"
+	  "\tMASK(STR):      Mask filename\n"
+	  "\tVARIANCE(STR):  Variance map filename\n"
+	  "OUTROOT is the 'root name' for output files\n"
 	  "\n"
 	  "additional options:\n"
@@ -136,6 +124,4 @@
 	  "  -trace-levels\n"
 	  "     print current trace levels\n");
-    psFree(config);
-    pmConfigDone();
     psLibFinalize();
     exit(PS_EXIT_SUCCESS);
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackChisqImage.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackChisqImage.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackChisqImage.c	(revision 27625)
@@ -12,6 +12,16 @@
     psTimerStart ("psphot.chisq.image");
 
+    pmFPAfile *chisqFile = pmFPAfileSelectSingle(config->files, "PSPHOT.CHISQ.OUTPUT", 0);
+    psAssert (chisqFile, "missing chisq image FPA?");
+
+    pmCell *chisqCell = pmFPAviewThisCell(view, chisqFile->fpa);
+
     // create a holder for the image
-    pmReadout *chiImage = pmReadoutAlloc();
+    pmReadout *chiReadout = pmFPAviewThisReadout(view, chisqFile->fpa);
+    if (!chiReadout) {
+      chiReadout = pmReadoutAlloc(chisqCell);
+    } else {
+      psMemIncrRefCounter(chiReadout);
+    }
 
     int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM");
@@ -20,5 +30,5 @@
     // loop over the available readouts
     for (int i = 0; i < num; i++) {
-        if (!psphotStackChisqImageAddReadout(config, view, chiImage, "PSPHOT.INPUT", i)) {
+        if (!psphotStackChisqImageAddReadout(config, view, chiReadout, "PSPHOT.INPUT", i)) {
             psError (PSPHOT_ERR_CONFIG, false, "failed to model background for PSPHOT.INPUT entry %d", i);
             return false;
@@ -26,13 +36,25 @@
     }
 
+    num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM");
+    psAssert (status, "programming error: must define PSPHOT.INPUT.NUM");
+
+    psMetadataAddS32(config->arguments, PS_LIST_TAIL, "PSPHOT.CHISQ.NUM", PS_META_REPLACE, "", num);
+    num++;
+    psMetadataAddS32(config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "", num);
+
     // need to save the resulting image somewhere (PSPHOT.INPUT?)
+    if (!psMetadataAddPtr(config->files, PS_LIST_TAIL, "PSPHOT.INPUT", PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "", chisqFile)) {
+        psError(PM_ERR_CONFIG, false, "could not add chisqFPA to config files");
+        return false;
+    }
 
     psLogMsg ("psphot", PS_LOG_INFO, "built chisq image: %f sec\n", psTimerMark ("psphot.chisq.image"));
 
+    psFree (chiReadout);
     return true;
 }
 
 bool psphotStackChisqImageAddReadout(const pmConfig *config, // Configuration
-				     pmFPAview *view,
+				     const pmFPAview *view,
 				     pmReadout *chiReadout,
 				     char *filename, 
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackImageLoop.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackImageLoop.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackImageLoop.c	(revision 27625)
@@ -7,10 +7,5 @@
 }
 
-bool psphotImageLoop (pmConfig *config) {
-
-    bool status;
-    pmChip *chip;
-    pmCell *cell;
-    pmReadout *readout;
+bool psphotStackImageLoop (pmConfig *config) {
 
     pmFPAview *view = pmFPAviewAlloc (0);
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackParseCamera.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackParseCamera.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackParseCamera.c	(revision 27625)
@@ -15,7 +15,7 @@
     }
 
-    psMetadata *item == NULL;
+    psMetadataItem *item = NULL;
     int nInputs = 0;
-    while ((item = psMetadataGetIndex(inputs, nInputs)) != NULL) {
+    while ((item = psMetadataGet(inputs, nInputs)) != NULL) {
 	if (item->type != PS_DATA_METADATA) {
 	    psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Component %s of the input metadata is not of type METADATA", item->name);
@@ -34,5 +34,5 @@
 	pmFPAfile *imageFile = defineFile(config, NULL, "PSPHOT.INPUT", image, PM_FPA_FILE_IMAGE); // File for image
 	if (!imageFile) {
-	    psError(PS_ERR_UNKNOWN, false, "Unable to define file from image %d (%s)", i, image);
+	    psError(PS_ERR_UNKNOWN, false, "Unable to define file from image %d (%s)", nInputs, image);
 	    return false;
 	}
@@ -41,13 +41,13 @@
 	if (mask && strlen(mask) > 0) {
 	    if (!defineFile(config, imageFile, "PSPHOT.INPUT.MASK", mask, PM_FPA_FILE_MASK)) {
-		psError(PS_ERR_UNKNOWN, false, "Unable to define file from mask %d (%s)", i, mask);
+		psError(PS_ERR_UNKNOWN, false, "Unable to define file from mask %d (%s)", nInputs, mask);
 		return false;
 	    }
 	}
 
-	psString variance = psMetadataLookupStr(&mdok, input, "VARIANCE"); // Name of variance map
+	psString variance = psMetadataLookupStr(&status, input, "VARIANCE"); // Name of variance map
 	if (variance && strlen(variance) > 0) {
 	    if (!defineFile(config, imageFile, "PSPHOT.INPUT.VARIANCE", variance, PM_FPA_FILE_VARIANCE)) {
-		psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", i, variance);
+		psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", nInputs, variance);
 		return false;
 	    }
@@ -58,4 +58,30 @@
     psMetadataAddS32 (config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "number of inputs", nInputs);
 
+    // generate an pmFPAimage for the chisqImage
+    pmFPAfile *chisqImage = pmFPAfileDefineOutput(config, NULL, "PSPHOT.CHISQ.IMAGE");
+    if (!chisqImage) {
+        psError(PSPHOT_ERR_CONFIG, false, "Trouble defining PSPHOT.CHISQ.IMAGE");
+        return false;
+    }
+    pmFPAfile *chisqMask = pmFPAfileDefineOutput(config, chisqImage->fpa, "PPIMAGE.CHISQ.MASK");
+    if (!chisqMask) {
+        psError(PS_ERR_IO, false, _("Unable to generate output file from PPIMAGE.CHISQ.MASK"));
+        return NULL;
+    }
+    if (chisqMask->type != PM_FPA_FILE_MASK) {
+        psError(PS_ERR_IO, true, "PPIMAGE.CHISQ.MASK is not of type MASK");
+        return NULL;
+    }
+    pmFPAfile *chisqVariance = pmFPAfileDefineOutput(config, chisqImage->fpa, "PPIMAGE.CHISQ.VARIANCE");
+    if (!chisqVariance) {
+        psError(PS_ERR_IO, false, _("Unable to generate output file from PPIMAGE.CHISQ.VARIANCE"));
+        return NULL;
+    }
+    if (chisqVariance->type != PM_FPA_FILE_VARIANCE) {
+        psError(PS_ERR_IO, true, "PPIMAGE.CHISQ.VARIANCE is not of type VARIANCE");
+        return NULL;
+    }
+
+# if (0)    
     // define the additional input/output files associated with psphot
     // XXX figure out which files are needed by psphotStack
@@ -64,4 +90,5 @@
         return false;
     }
+# endif
 
     psTrace("psphot", 1, "Done with psphotStackParseCamera...\n");
Index: /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackReadout.c
===================================================================
--- /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackReadout.c	(revision 27624)
+++ /branches/eam_branches/stackphot.20100406/psphot/src/psphotStackReadout.c	(revision 27625)
@@ -13,6 +13,9 @@
         return false;
     }
+    // optional break-point for processing
+    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
+    PS_ASSERT_PTR_NON_NULL (breakPt, false);
 
-    // set the photcode for this image (XXX currently goes into recipe, should it go into analysis?)
+    // set the photcode for each image
     if (!psphotAddPhotcode (config, view)) {
         psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode");
@@ -20,173 +23,121 @@
     }
 
-    // generate a background model (median, smoothed image)
-    if (!psphotSetMaskAndVariance (config, view, recipe)) {
-        return psphotReadoutCleanup (config, NULL, recipe, NULL, NULL, NULL);
+    // Generate the mask and weight images
+    if (!psphotSetMaskAndVariance (config, view)) {
+	return psphotReadoutCleanup (config, view);
     }
     if (!strcasecmp (breakPt, "NOTHING")) {
-	return psphotReadoutCleanup(config, NULL, recipe, NULL, NULL, NULL);
+	return psphotReadoutCleanup (config, view);
     }
 
-    // optional break-point for processing
-    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
-    PS_ASSERT_PTR_NON_NULL (breakPt, false);
-
     // generate a background model (median, smoothed image)
+    // XXX I think this is not defined correctly for an array of images.
     if (!psphotModelBackground (config, view)) {
-        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
+	return psphotReadoutCleanup (config, view);
     }
     if (!psphotSubtractBackground (config, view)) {
-        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
+	return psphotReadoutCleanup (config, view);
     }
     if (!strcasecmp (breakPt, "BACKMDL")) {
-        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
+	return psphotReadoutCleanup (config, view);
+    }
+
+    // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are determined and saved on
+    // readout->analysis XXX this function currently only works with a single PSPHOT.INPUT
+    if (!psphotLoadPSF (config, view)) {
+        psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model");
+        return psphotReadoutCleanup (config, view);
+    }
+
+    if (!psphotStackChisqImage(config, view)) {
+        psError (PSPHOT_ERR_UNKNOWN, false, "failure to generate chisq image");
+        return psphotReadoutCleanup (config, view);
+    }
+    if (!strcasecmp (breakPt, "CHISQ")) {
+	return psphotReadoutCleanup (config, view);
     }
 
     // find the detections (by peak and/or footprint) in the image.
-    pmDetections *detections = psphotFindDetections (NULL, readout, recipe);
-    if (!detections) {
-        psLogMsg ("psphot", 3, "unable to find detections in this image");
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, NULL);
+    if (!psphotFindDetections (config, view, true)) { // pass 1
+        // this only happens if we had an error in psphotFindDetections
+        psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
+        return psphotReadoutCleanup (config, view);
     }
 
-    // construct sources and measure basic stats
-    psArray *sources = psphotSourceStats (config, readout, detections, true);
-    if (!sources) return false;
+    // construct sources and measure basic stats (saved on detections->newSources)
+    if (!psphotSourceStats (config, view, true)) { // pass 1
+        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
+        return psphotReadoutCleanup (config, view);
+    }
     if (!strcasecmp (breakPt, "PEAKS")) {
-        return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
+        return psphotReadoutCleanup(config, view);
     }
 
-    // find blended neighbors of very saturated stars
-    // XXX merge this with Basic Deblend?
-    psphotDeblendSatstars (readout, sources, recipe);
+    // find blended neighbors of very saturated stars (detections->newSources)
+    if (!psphotDeblendSatstars (config, view)) {
+        psError (PSPHOT_ERR_UNKNOWN, false, "failed on satstar deblend analysis");
+        return psphotReadoutCleanup (config, view);
+    }
 
-    // mark blended peaks PS_SOURCE_BLEND
-    if (!psphotBasicDeblend (sources, recipe)) {
-        psLogMsg ("psphot", 3, "failed on deblend analysis");
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
+    // mark blended peaks PS_SOURCE_BLEND (detections->newSources)
+    if (!psphotBasicDeblend (config, view)) {
+        psError (PSPHOT_ERR_UNKNOWN, false, "failed on deblend analysis");
+        return psphotReadoutCleanup (config, view);
     }
 
     // classify sources based on moments, brightness
-    if (!psphotRoughClass (readout, sources, recipe, havePSF)) {
-        psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
+    if (!psphotRoughClass (config, view)) {
+        psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough classifications");
+	return psphotReadoutCleanup (config, view);
+    }
+    // if we were not supplied a PSF model, determine the IQ stats here (detections->newSources)
+    if (!psphotImageQuality (config, view)) { // pass 1
+        psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality");
+        return psphotReadoutCleanup(config, view);
     }
     if (!strcasecmp (breakPt, "MOMENTS")) {
-        return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
+	return psphotReadoutCleanup (config, view);
     }
 
-    if (!havePSF && !psphotImageQuality (recipe, sources)) {
-        psLogMsg("psphot", 3, "failed to measure image quality");
-        return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
+    // use bright stellar objects to measure PSF if we were supplied a PSF for any input file,
+    // this step is skipped
+    if (!psphotChoosePSF (config, view)) { // pass 1
+        psLogMsg ("psphot", 3, "failure to construct a psf model");
+        return psphotReadoutCleanup (config, view);
+    }
+    if (!strcasecmp (breakPt, "PSFMODEL")) {
+        return psphotReadoutCleanup (config, view);
     }
 
-    // if we were not supplied a PSF, choose one here
-    if (psf == NULL) {
-        // use bright stellar objects to measure PSF
-        // XXX if we do not have enough stars to generate the PSF, build one
-        // from the SEEING guess and model class
-        psf = psphotChoosePSF (readout, sources, recipe);
-        if (psf == NULL) {
-            psLogMsg ("psphot", 3, "failure to construct a psf model");
-            return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
-        }
-        havePSF = true;
-    }
-    if (!strcasecmp (breakPt, "PSFMODEL")) {
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
-    }
-    psphotVisualShowPSFModel (readout, psf);
-
     // include externally-supplied sources
-    psphotLoadExtSources (config, view, sources);
+    // XXX fix this in the new multi-input context
+    // psphotLoadExtSources (config, view); // pass 1
 
     // construct an initial model for each object, set the radius to fitRadius, set circular fit mask
-    psphotGuessModels (config, readout, sources, psf);
+    psphotGuessModels (config, view);
 
+    // merge the newly selected sources into the existing list
+    // NOTE: merge OLD and NEW
+    psphotMergeSources (config, view);
+
+    // *** generate the objects (which unify the sources from the different images)
+    // pmArray *objects = psphotMatchSources (config, view);
+    
     // linear PSF fit to source peaks, subtract the models from the image (in PSF mask)
-    psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE);
-
-    // We have to place these visualizations here because the models are not realized until
-    // psphotGuessModels or fitted until psphotFitSourcesLinear.
-    psphotVisualShowPSFStars (recipe, psf, sources);
+    // psphotFitSourcesLinearStack (config, objects, FALSE);
 
     // identify CRs and extended sources
-    psphotSourceSize (config, readout, sources, recipe, psf, 0);
-    if (!strcasecmp (breakPt, "ENSEMBLE")) {
-        goto finish;
-    }
-    psphotVisualShowSatStars (recipe, psf, sources);
-
-    // non-linear PSF and EXT fit to brighter sources
-    // replace model flux, adjust mask as needed, fit, subtract the models (full stamp)
-    psphotBlendFit (config, readout, sources, psf);
-
-    // replace all sources
-    psphotReplaceAllSources (sources, recipe);
-
-    // linear fit to include all sources (subtract again)
-    psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE);
-
-    // if we only do one pass, skip to extended source analysis
-    if (!strcasecmp (breakPt, "PASS1")) {
-        goto pass1finish;
-    }
-    // NOTE: possibly re-measure background model here with objects subtracted
-
-    // add noise for subtracted objects
-    psphotAddNoise (readout, sources, recipe);
-
-    // find fainter sources (pass 2)
-    detections = psphotFindDetections (detections, readout, recipe);
-
-    // remove noise for subtracted objects (ie, return to normal noise level)
-    psphotSubNoise (readout, sources, recipe);
-
-    // define new sources based on only the new peaks
-    psArray *newSources = psphotSourceStats (config, readout, detections, false);
-
-    // set source type
-    if (!psphotRoughClass (readout, newSources, recipe, havePSF)) {
-        psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
-    }
-
-    // create full input models, set the radius to fitRadius, set circular fit mask
-    psphotGuessModels (config, readout, newSources, psf);
-
-    // replace all sources so fit below applies to all at once
-    psphotReplaceAllSources (sources, recipe);
-
-    // merge the newly selected sources into the existing list
-    psphotMergeSources (sources, newSources);
-    psFree (newSources);
-
-    // linear fit to all sources
-    psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE);
-
-pass1finish:
-
-    // measure source size for the remaining sources
-    psphotSourceSize (config, readout, sources, recipe, psf, 0);
-
-    psphotExtendedSourceAnalysis (readout, sources, recipe);
-
-    psphotExtendedSourceFits (readout, sources, recipe);
-
-finish:
-
-    // plot positive sources
-    // psphotSourcePlots (readout, sources, recipe);
+    psphotSourceSize (config, view, TRUE);
 
     // measure aperture photometry corrections
-    if (!psphotApResid (config, readout, sources, psf)) {
+    if (!psphotApResid (config, view)) {
         psLogMsg ("psphot", 3, "failed on psphotApResid");
-        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
+	return psphotReadoutCleanup (config, view);
     }
 
     // calculate source magnitudes
-    psphotMagnitudes(config, readout, view, sources, psf);
+    psphotMagnitudes(config, view);
 
-    if (!psphotEfficiency(config, readout, view, psf, recipe, sources)) {
+    if (!psphotEfficiency(config, view)) {
         psErrorStackPrint(stderr, "Unable to determine detection efficiencies from fake sources");
         psErrorClear();
@@ -200,8 +151,8 @@
 
     // drop the references to the image pixels held by each source
-    psphotSourceFreePixels (sources);
+    psphotSourceFreePixels (config, view);
 
     // create the exported-metadata and free local data
-    return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
+    return psphotReadoutCleanup (config, view);
 }