Index: trunk/psphot/src/psphotSourceStats.c
===================================================================
--- trunk/psphot/src/psphotSourceStats.c	(revision 41358)
+++ trunk/psphot/src/psphotSourceStats.c	(revision 41359)
@@ -560,5 +560,5 @@
         }
 
-#if (0)
+#if (1)
 	if (psTraceGetLevel("psphot.moments.save")) {
 	  char name[64];
@@ -604,21 +604,48 @@
 #endif
 
-        Sout[i] = sqrt(0.5*(psfClump.X + psfClump.Y)) / sigma[i];
+	// a clump with no stars is not a valid clump
 	Nout[i] = psfClump.nStars;
+        Sout[i] = (Nout[i] == 0) ? NAN : sqrt(0.5*(psfClump.X + psfClump.Y)) / sigma[i];
     }
 
     // we are looking for sigma for which Sout = 0.65 (or some other value)
+    // EAM 2020.05.19 : this analysis fails if any of the measurements above have no valid stars
+    // A result with Nout[] = 0 should be skipped in the analysis
+    
     int Nstars = 0;
     float minKronRadius = NAN;
     float Sigma = NAN;
-    float minS = Sout[0];
-    float maxS = Sout[0];
+
+    float minS = NAN;
+    float maxS = NAN;
+    int maxN = -1;
+    int minN = -1;
     for (int i = 0; i < nsigma; i++) {
-        minS = PS_MIN(Sout[i], minS);
-        maxS = PS_MAX(Sout[i], maxS);
-    }
-    if (minS > 0.65) Sigma = sigma[nsigma-1];
-    if (maxS < 0.65) Sigma = sigma[0];
-
+      if (!isfinite(Sout[i])) continue;
+      if (!isfinite(minS)) {
+	minS = Sout[i];
+	maxS = Sout[i];
+	minN = i;
+	maxN = i;
+      }
+      if (Sout[i] < minS) {
+	  minS = Sout[i];
+	  minN = i;
+      }
+      if (Sout[i] > maxS) {
+	  maxS = Sout[i];
+	  maxN = i;
+      }
+    }
+    if (!isfinite(minS)) {
+      psLogMsg ("psphot", 3, "failed to find a valid window function\n");
+      return false;
+    }
+
+    // if all Sout values are above (or below) the threshold, saturate on the minimum (or maximum) valid sigma value
+    if (minS > 0.65) { Sigma = sigma[minN]; Nstars = Nout[minN]; minKronRadius = Rmin[minN]; }
+    if (maxS < 0.65) { Sigma = sigma[maxN]; Nstars = Nout[maxN]; minKronRadius = Rmin[maxN]; }
+
+    // this loop will be skipped if Sigma is set above
     for (int i = 0; i < nsigma - 1 && isnan(Sigma); i++) {
         if (!isfinite(Sout[i]) || !isfinite(Sout[i+1])) continue;
@@ -629,5 +656,9 @@
         minKronRadius = Rmin[i] + (0.65 - Sout[i])*(Rmin[i+1] - Rmin[i])/(Sout[i+1] - Sout[i]);
     }
-    psAssert (isfinite(Sigma), "did we miss a case?");
+    // EAM: this test will fail if only a single entry above is valid
+    if (!isfinite(Sigma)) {
+      psLogMsg ("psphot", 3, "failed to find a valid window function (case 2)\n");
+      return false;
+    }
     if (!isfinite(minKronRadius)) {
         minKronRadius = Sigma;
