Index: trunk/psModules/src/objects/pmSource.h
===================================================================
--- trunk/psModules/src/objects/pmSource.h	(revision 36710)
+++ trunk/psModules/src/objects/pmSource.h	(revision 36718)
@@ -119,5 +119,5 @@
     pmSourceExtendedPars *extpars;      ///< extended source parameters
     pmSourceDiffStats *diffStats;       ///< extra parameters for difference detections
-    pmSourceGalaxyFits *galaxyFits;     ///< fits to galaxy models (psphotFullForce only)
+    psArray *galaxyFits;                ///< fits to galaxy models (psphotFullForce only)
     pmSourceLensing *lensingOBJ;        ///< lensing moments parameters (per object)
     pmSourceLensing *lensingPSF;        ///< lensing moments parameters (psf, interpolated)
Index: trunk/psModules/src/objects/pmSourceExtendedPars.h
===================================================================
--- trunk/psModules/src/objects/pmSourceExtendedPars.h	(revision 36710)
+++ trunk/psModules/src/objects/pmSourceExtendedPars.h	(revision 36718)
@@ -83,4 +83,5 @@
 
 typedef struct {
+  int       modelType;
   psVector *Flux;
   psVector *dFlux;
Index: trunk/psModules/src/objects/pmSourceIO_CMF.c.in
===================================================================
--- trunk/psModules/src/objects/pmSourceIO_CMF.c.in	(revision 36710)
+++ trunk/psModules/src/objects/pmSourceIO_CMF.c.in	(revision 36718)
@@ -482,6 +482,8 @@
     }
 
+#ifdef SORT_EXTENSIONS_BY_FLUX
     // let's write these out in S/N order
     sources = psArraySort (sources, pmSourceSortByFlux);
+#endif
 
     table = psArrayAllocEmpty (sources->n);
@@ -821,6 +823,8 @@
     psMetadataAddStr (outhead, PS_LIST_TAIL, "EXTNAME", PS_META_REPLACE, "xsrc table extension", extname);
 
+#ifdef SORT_EXTENSIONS_BY_FLUX
     // let's write these out in S/N order
     sources = psArraySort (sources, pmSourceSortByFlux);
+#endif
 
     float magOffset; 
@@ -1225,6 +1229,8 @@
     psVector *fwhmValues = psMetadataLookupVector(&status, readout->analysis, "STACK.PSF.FWHM.VALUES");
 
+#ifdef SORT_EXTENSIONS_BY_FLUX
     // let's write these out in S/N order
     sources = psArraySort (sources, pmSourceSortByFlux);
+#endif
 
     table = psArrayAllocEmpty (sources->n);
@@ -1441,6 +1447,6 @@
     psMetadataAddStr (outhead, PS_LIST_TAIL, "HI", PS_META_REPLACE, "does this get through?", "THERE");
 
-    // let's write these out in S/N order
-    sources = psArraySort (sources, pmSourceSortByFlux);
+    // They are probably already in this order but ...
+    sources = psArraySort (sources, pmSourceSortBySeq);
 
     psArray *table = psArrayAllocEmpty (sources->n);
@@ -1459,28 +1465,36 @@
         if (source->galaxyFits == NULL) continue;
 
-	pmModel *model = source->modelFits->data[0];
-	if (!model) return false;
-
-	// X,Y coordinates are stored with the model parameters
- 	psF32 *PAR = model->params->data.F32;
-
-	psMetadata *row = psMetadataAlloc ();
-
-	// we write out the x,y positions so people can link to the psf either way (position or ID)
-	psMetadataAddU32 (row, PS_LIST_TAIL, "IPP_IDET",         0, "IPP detection identifier index", source->seq);
-	psMetadataAddF32 (row, PS_LIST_TAIL, "X_FIT",            0, "model x coordinate",             PAR[PM_PAR_XPOS]);
-	psMetadataAddF32 (row, PS_LIST_TAIL, "Y_FIT",            0, "model y coordinate",             PAR[PM_PAR_YPOS]);
-	psMetadataAddF32 (row, PS_LIST_TAIL, "NPIX",             0, "number of pixels for fits",      source->galaxyFits->nPix);
-
-	psVector *Flux = source->galaxyFits->Flux;
-	psVector *dFlux = source->galaxyFits->dFlux;
-	psVector *chisq = source->galaxyFits->chisq;
-
-	psMetadataAddVector (row, PS_LIST_TAIL, "GAL_FLUX",     PS_META_REPLACE, "normalization for galaxy flux", Flux);
-	psMetadataAddVector (row, PS_LIST_TAIL, "GAL_FLUX_ERR", PS_META_REPLACE, "error on normalization", dFlux);
-	psMetadataAddVector (row, PS_LIST_TAIL, "GAL_CHISQ",    PS_META_REPLACE, "galaxy fit chisq", chisq);
-
-	psArrayAdd (table, 100, row);
-	psFree (row);
+        for (int iModel = 0; iModel < source->modelFits->n; iModel++) {
+            pmModel *model = source->modelFits->data[iModel];
+            pmSourceGalaxyFits *galaxyFits = source->galaxyFits->data[iModel];
+
+            if (!model || !galaxyFits) return false;
+
+            // XXX: TODO: we should search for the model that has the same type as the fits
+            psAssert(model->type == galaxyFits->modelType, "galaxyFits model type does not match model type");
+
+            // X,Y coordinates are stored with the model parameters
+            psF32 *PAR = model->params->data.F32;
+
+            psMetadata *row = psMetadataAlloc ();
+
+            // we write out the x,y positions so people can link to the psf either way (position or ID)
+            psMetadataAddU32 (row, PS_LIST_TAIL, "IPP_IDET",         0, "IPP detection identifier index", source->seq);
+            psMetadataAddS32 (row, PS_LIST_TAIL, "MODEL_TYPE",       0, "model type",                     galaxyFits->modelType);
+            psMetadataAddF32 (row, PS_LIST_TAIL, "X_FIT",            0, "model x coordinate",             PAR[PM_PAR_XPOS]);
+            psMetadataAddF32 (row, PS_LIST_TAIL, "Y_FIT",            0, "model y coordinate",             PAR[PM_PAR_YPOS]);
+            psMetadataAddF32 (row, PS_LIST_TAIL, "NPIX",             0, "number of pixels for fits",      galaxyFits->nPix);
+
+            psVector *Flux = galaxyFits->Flux;
+            psVector *dFlux = galaxyFits->dFlux;
+            psVector *chisq = galaxyFits->chisq;
+
+            psMetadataAddVector (row, PS_LIST_TAIL, "GAL_FLUX",     PS_META_REPLACE, "normalization for galaxy flux", Flux);
+            psMetadataAddVector (row, PS_LIST_TAIL, "GAL_FLUX_ERR", PS_META_REPLACE, "error on normalization", dFlux);
+            psMetadataAddVector (row, PS_LIST_TAIL, "GAL_CHISQ",    PS_META_REPLACE, "galaxy fit chisq", chisq);
+
+            psArrayAdd (table, 100, row);
+            psFree (row);
+        }
     }
 
@@ -1541,4 +1555,5 @@
         }
 
+        int modelType = psMetadataLookupS32 (&status,    row, "MODEL_TYPE");
         psVector *Flux  = psMetadataLookupVector(&status, row, "GAL_FLUX");
         psVector *dFlux = psMetadataLookupVector(&status, row, "GAL_FLUX_ERR");
@@ -1546,14 +1561,22 @@
 
         if (Flux && Flux->n > 0) {
-            psFree(source->galaxyFits);
-            source->galaxyFits = pmSourceGalaxyFitsAlloc();
-            source->galaxyFits->nPix = psMetadataLookupF32(&status, row, "NPIX");
-
-            psFree(source->galaxyFits->Flux);
-            source->galaxyFits->Flux  = psMemIncrRefCounter(Flux);
-            psFree(source->galaxyFits->dFlux);
-            source->galaxyFits->dFlux = psMemIncrRefCounter(dFlux);
-            psFree(source->galaxyFits->chisq);
-            source->galaxyFits->chisq = psMemIncrRefCounter(chisq);
+            if (!source->galaxyFits) {
+                source->galaxyFits = psArrayAllocEmpty(1);
+            }
+
+            pmSourceGalaxyFits *galaxyFits = pmSourceGalaxyFitsAlloc();
+
+            psArrayAdd(source->galaxyFits, 1, galaxyFits);
+
+            psFree(galaxyFits);
+            galaxyFits->modelType = modelType;
+            galaxyFits->nPix = psMetadataLookupF32(&status, row, "NPIX");
+
+            psFree(galaxyFits->Flux);
+            galaxyFits->Flux  = psMemIncrRefCounter(Flux);
+            psFree(galaxyFits->dFlux);
+            galaxyFits->dFlux = psMemIncrRefCounter(dFlux);
+            psFree(galaxyFits->chisq);
+            galaxyFits->chisq = psMemIncrRefCounter(chisq);
         }