Index: trunk/psModules/src/objects/pmSource.c
===================================================================
--- trunk/psModules/src/objects/pmSource.c	(revision 30621)
+++ trunk/psModules/src/objects/pmSource.c	(revision 31153)
@@ -171,9 +171,12 @@
     // peak has the same values as the original
     if (in->peak != NULL) {
-        source->peak = pmPeakAlloc (in->peak->x, in->peak->y, in->peak->value, in->peak->type);
+        source->peak = pmPeakAlloc (in->peak->x, in->peak->y, in->peak->detValue, in->peak->type);
         source->peak->xf = in->peak->xf;
         source->peak->yf = in->peak->yf;
-        source->peak->flux = in->peak->flux;
-        source->peak->SN = in->peak->SN;
+        source->peak->rawFlux         = in->peak->rawFlux;
+        source->peak->rawFluxStdev    = in->peak->rawFluxStdev;
+        source->peak->smoothFlux      = in->peak->smoothFlux;
+        source->peak->smoothFluxStdev = in->peak->smoothFluxStdev;
+        // source->peak->SN = in->peak->SN;
     }
 
@@ -328,8 +331,4 @@
 *****************************************************************************/
 
-// XXX EAM include a S/N cutoff in selecting the sources?
-// XXX this function should probably accept the values, not a recipe. wrap with a
-// psphot-specific function which applies the recipe values
-// only apply selection to sources within specified region
 pmPSFClump pmSourcePSFClump(psImage **savedImage, psRegion *region, psArray *sources, float PSF_SN_LIM, float PSF_CLUMP_GRID_SCALE, psF32 SX_MAX, psF32 SY_MAX, psF32 AR_MAX)
 {
@@ -429,4 +428,5 @@
 
 	psfClump.nSigma = stats->sampleStdev;
+	psfClump.nTotal = nValid;
 
 	if (savedImage) {
@@ -465,13 +465,13 @@
         psStats *stats  = NULL;
 
-        // select the single highest peak
-        psArraySort (peaks, pmPeaksCompareDescend);
+        // select the single highest peak (note that we only have detValue, not rawFlux, etc
+        psArraySort (peaks, pmPeaksSortByDetValueDescend);
         clump = peaks->data[0];
-        psTrace ("psModules.objects", 2, "clump is at %d, %d (%f)\n", clump->x, clump->y, clump->value);
+        psTrace ("psModules.objects", 2, "clump is at %d, %d (%f)\n", clump->x, clump->y, clump->detValue);
 
 	// XXX store the mean sigma?
 	float meanSigma = psfClump.nSigma;
-	psfClump.nStars = clump->value;
-	psfClump.nSigma = clump->value / meanSigma;
+	psfClump.nStars = clump->detValue;
+	psfClump.nSigma = clump->detValue / meanSigma;
 
         // define section window for clump
@@ -643,4 +643,22 @@
             }
 
+	    // check for insignificant sources or excessively low-surface brightness
+	    float coreSN = source->moments->KronCore / source->moments->KronCoreErr;
+	    float coreKR = source->moments->KronCore / source->moments->KronFlux;
+
+	    // XXX these values need to be in the recipe...
+	    if (false && isfinite(coreSN) && (coreSN < 5.0)) {
+                source->type = PM_SOURCE_TYPE_DEFECT;
+                source->mode |= PM_SOURCE_MODE_DEFECT;
+                Ncr ++;
+                continue;
+            }
+	    if (false && isfinite(coreKR) && (coreKR < 0.1)) {
+                source->type = PM_SOURCE_TYPE_DEFECT;
+                source->mode |= PM_SOURCE_MODE_DEFECT;
+                Ncr ++;
+                continue;
+            }
+
             // likely unsaturated extended source (too large to be stellar)
             if (sigX > clump.X + 3*clump.dX || sigY > clump.Y + 3*clump.dY) {
@@ -651,7 +669,8 @@
 
             // the rest are probable stellar objects
+	    // the vectors below are accumulated to give user feedback on the S/N ranges
             starsn_moments->data.F32[starsn_moments->n] = source->moments->SN;
             starsn_moments->n ++;
-            starsn_peaks->data.F32[starsn_peaks->n] = source->peak->SN;
+            starsn_peaks->data.F32[starsn_peaks->n] = sqrt(source->peak->detValue);
             starsn_peaks->n ++;
             Nstar ++;
@@ -1149,12 +1168,33 @@
 }
 
+// this function decides if the source position should be based on the peak or the moments.
+// this is only used if we know we should not use a model fit position (eg, no model, or no
+// model yet)
+bool pmSourcePositionUseMoments(pmSource *source) {
+
+    if (!source->moments) return false;		 // can't if there are no moments
+    if (!source->moments->nPixels) return false; // can't if the moments were not measured
+    if (source->mode && PM_SOURCE_MODE_MOMENTS_FAILURE) return false; // can't if the moments failed...
+
+    if (source->mode & PM_SOURCE_MODE_SATSTAR) return true; // moments are best for SATSTARs
+
+    float dX = source->moments->Mx - source->peak->xf;
+    float dY = source->moments->My - source->peak->yf;
+    float dR = hypot(dX, dY);
+    
+    // only use the moments position if the moment-peak offset is small or the star is saturated
+    if (dR > 1.5) return false;
+
+    return true;
+}
+
 // sort by SN (descending)
-int pmSourceSortBySN (const void **a, const void **b)
+int pmSourceSortByFlux (const void **a, const void **b)
 {
     pmSource *A = *(pmSource **)a;
     pmSource *B = *(pmSource **)b;
 
-    psF32 fA = (A->peak == NULL) ? 0 : A->peak->SN;
-    psF32 fB = (B->peak == NULL) ? 0 : B->peak->SN;
+    psF32 fA = (A->peak == NULL) ? 0 : A->peak->rawFlux;
+    psF32 fB = (B->peak == NULL) ? 0 : B->peak->rawFlux;
     if (isnan (fA)) fA = 0;
     if (isnan (fB)) fB = 0;
@@ -1166,4 +1206,23 @@
 }
 
+// sort by SN (descending)
+int pmSourceSortByParentFlux (const void **a, const void **b)
+{
+    pmSource *Ao = *(pmSource **)a;
+    pmSource *Bo = *(pmSource **)b;
+    pmSource *A  = Ao->parent;
+    pmSource *B  = Bo->parent;
+
+    psF32 fA = (A->peak == NULL) ? 0 : A->peak->rawFlux;
+    psF32 fB = (B->peak == NULL) ? 0 : B->peak->rawFlux;
+    if (isnan (fA)) fA = 0;
+    if (isnan (fB)) fB = 0;
+
+    psF32 diff = fA - fB;
+    if (diff > FLT_EPSILON) return (-1);
+    if (diff < FLT_EPSILON) return (+1);
+    return (0);
+}
+
 // sort by Y (ascending)
 int pmSourceSortByY (const void **a, const void **b)
@@ -1201,4 +1260,21 @@
     pmSource *A = *(pmSource **)a;
     pmSource *B = *(pmSource **)b;
+
+    int iA = A->seq;
+    int iB = B->seq;
+
+    int diff = iA - iB;
+    if (diff > 0) return (+1);
+    if (diff < 0) return (-1);
+    return (0);
+}
+
+// sort by Seq (ascending)
+int pmSourceSortByParentSeq (const void **a, const void **b)
+{
+    pmSource *Ao = *(pmSource **)a;
+    pmSource *Bo = *(pmSource **)b;
+    pmSource *A  = Ao->parent;
+    pmSource *B  = Bo->parent;
 
     int iA = A->seq;