Index: trunk/psphot/src/psphotMergeSources.c
===================================================================
--- trunk/psphot/src/psphotMergeSources.c	(revision 38389)
+++ trunk/psphot/src/psphotMergeSources.c	(revision 38515)
@@ -853,22 +853,37 @@
 // array containing the child sources.  XXX currently, this is only used by psphotStackReadout
 // (sources go on allSources so that psphotChoosePSF can be called repeatedly)
-psArray *psphotSourceChildrenByObject (pmConfig *config, const pmFPAview *view, const char *filerule, psArray *objectsSrc, bool sourcesSubtracted) {
+psArray *psphotSourceChildrenByObject (pmConfig *config, const pmFPAview *view, const char *fileruleOut, const char *fileruleSrc, psArray *objectsSrc, bool sourcesSubtracted) {
 
     bool status;
 
-    int nImages = psphotFileruleCount(config, filerule);
+    int nImages = psphotFileruleCount(config, fileruleOut);
 
     // generate look-up arrays for detections and readouts
     psArray *detArrays = psArrayAlloc(nImages);
     psArray *readouts = psArrayAlloc(nImages);
+    psArray *fitOptionsArray = psArrayAlloc(nImages);
+
+    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
+    assert (recipe);
+    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT");
+    assert (maskVal);
+    int psfSize  = psMetadataLookupS32 (&status, recipe, "PCM_BOX_SIZE");
+    assert (status);
 
     for (int i = 0; i < nImages; i++) {
 
 	// find the currently selected readout
-	pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest
+	pmFPAfile *file = pmFPAfileSelectSingle(config->files, fileruleOut, i); // File of interest
 	psAssert (file, "missing file?");
 
 	pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
 	psAssert (readout, "missing readout?");
+
+	pmFPAfile *fileSrc = pmFPAfileSelectSingle(config->files, fileruleSrc, i); // File of interest
+	psAssert (file, "missing file?");
+
+	pmReadout *readoutSrc = pmFPAviewThisReadout(view, fileSrc->fpa);
+	psAssert (readoutSrc, "missing readout?");
+
 
 	// create DETECTIONS containers for each image, in case one lacks it
@@ -886,8 +901,12 @@
 	    psAssert (detections, "missing detections?");
 	}
+        pmSourceFitOptions *fitOptions = psMetadataLookupPtr (&status, readoutSrc->analysis, "PCM_FIT_OPTIONS");
+        psAssert (fitOptions, "missing pcm fit options");
+        psMetadataAddPtr (readout->analysis, PS_LIST_TAIL, "PCM_FIT_OPTIONS", PS_DATA_UNKNOWN | PS_META_REPLACE, "pcm fit options", fitOptions);
 
 	// we need to save the new sources on the detection arrays of the appropriate image
 	detArrays->data[i] = psMemIncrRefCounter(detections);
 	readouts->data[i] = psMemIncrRefCounter(readout);
+        fitOptionsArray->data[i] = psMemIncrRefCounter(fitOptions);
     }
 
@@ -932,12 +951,22 @@
 	    // does this copy all model data? (NO)
 	    sourceOut->modelPSF = pmModelCopy(sourceSrc->modelPSF);
-	    sourceOut->modelEXT = pmModelCopy(sourceSrc->modelEXT);
-
+
+            bool foundModelEXT = false;
 	    if (sourceSrc->modelFits) {
 		sourceOut->modelFits = psArrayAlloc(sourceSrc->modelFits->n);
 		for (int j = 0; j < sourceSrc->modelFits->n; j++) {
-		    sourceOut->modelFits->data[j] = pmModelCopy(sourceSrc->modelFits->data[j]);
-		}
+                    pmModel *modelSrc = sourceSrc->modelFits->data[j];
+		    pmModel *modelOut = sourceOut->modelFits->data[j] = pmModelCopy(modelSrc);
+                    if (modelSrc == sourceSrc->modelEXT) {
+                        foundModelEXT = true;
+                        sourceOut->modelEXT = psMemIncrRefCounter (modelOut);
+                    }
+                    modelOut->isPCM = modelSrc->isPCM;
+                }
 	    }
+            if (!foundModelEXT && sourceSrc->modelEXT) {
+                // Will this ever happen?
+                sourceOut->modelEXT = pmModelCopy(sourceSrc->modelEXT);
+            }
 
 	    // drop the references to the original image pixels:
@@ -949,8 +978,10 @@
 	    pmReadout *readout = readouts->data[index];
 
+            pmSourceFitOptions *fitOptions = fitOptionsArray->data[index];
+
 	    // allocate image, weight, mask for the new image for each peak
 	    if (sourceOut->modelPSF) {
                 pmSourceRedefinePixels (sourceOut, readout, sourceOut->peak->x, sourceOut->peak->y, 
-                                                                        sourceOut->modelPSF->fitRadius);
+                                                                        sourceSrc->windowRadius);
 	    } else {
                 // if we have no pixels we can't use it to determine the psf so make sure this bit is off
@@ -965,5 +996,23 @@
 	    if (!sourcesSubtracted) {
 		sourceOut->tmpFlags &= ~PM_SOURCE_TMPF_SUBTRACTED;
-	    }
+	    } else {
+                if (sourceSrc->modelFlux) {
+                    bool isPSF = false;
+                    pmModel *model = pmSourceGetModel (&isPSF, sourceOut);
+                    if (model->isPCM) {
+                        pmPCMdata *pcm = pmPCMinit (sourceOut, fitOptions, model, maskVal, psfSize);
+                        if (pcm) {
+                            // pmPCMMakeModel (sourceOut, model, pcm->nsigma, maskVal, psfSize);
+                            pmPCMCacheModel (sourceOut, maskVal, psfSize, pcm->nsigma);
+                            psFree(pcm);
+                        } else {
+                            // What to do here? 
+                            psAssert (pcm, "pmPCMinit failed!");
+                        }
+                    } else {
+                        pmSourceCacheModel (sourceOut, maskVal);
+                    }
+                }
+            }
 
 	    // set the output detections:
@@ -981,4 +1030,5 @@
     psFree (detArrays);
     psFree (readouts);
+    psFree (fitOptionsArray);
 
     return objectsOut;
