Index: /tags/ipp-20110505/psModules/src/imcombine/pmSubtractionStamps.c
===================================================================
--- /tags/ipp-20110505/psModules/src/imcombine/pmSubtractionStamps.c	(revision 31538)
+++ /tags/ipp-20110505/psModules/src/imcombine/pmSubtractionStamps.c	(revision 31539)
@@ -834,15 +834,21 @@
             }
 
-            psStatsInit (stats);
-            if (!psVectorStats (stats, flux1, NULL, NULL, 0)) {
-                psAbort ("failed to generate stats");
-            }
-            float f1 = stats->sampleMedian;
-
-            psStatsInit (stats);
-            if (!psVectorStats (stats, flux2, NULL, NULL, 0)) {
-                psAbort ("failed to generate stats");
-            }
-            float f2 = stats->sampleMedian;
+            float f1 = NAN;
+            if (flux1->n > 0) {
+                psStatsInit (stats);
+                if (!psVectorStats (stats, flux1, NULL, NULL, 0)) {
+                    psAbort ("failed to generate stats");
+                }
+                f1 = stats->sampleMedian;
+            }
+
+            float f2 = NAN;
+            if (flux2->n > 0) {
+                psStatsInit (stats);
+                if (!psVectorStats (stats, flux2, NULL, NULL, 0)) {
+                    psAbort ("failed to generate stats");
+                }
+                f2 = stats->sampleMedian;
+            }
 
             stamps->window1->kernel[y][x] = f1;
Index: /tags/ipp-20110505/psModules/src/objects/models/pmModel_PS1_V1.c
===================================================================
--- /tags/ipp-20110505/psModules/src/objects/models/pmModel_PS1_V1.c	(revision 31538)
+++ /tags/ipp-20110505/psModules/src/objects/models/pmModel_PS1_V1.c	(revision 31539)
@@ -288,5 +288,4 @@
     if (PAR[PM_PAR_I0] <= 0) return 1.0;
     if (flux >= PAR[PM_PAR_I0]) return 1.0;
-    if (PAR[PM_PAR_7] == 0.0) return powf(PAR[PM_PAR_I0] / flux - 1.0, 1.0 / ALPHA);
 
     shape.sx  = PAR[PM_PAR_SXX] / M_SQRT2;
@@ -296,4 +295,9 @@
     psEllipseAxes axes = psEllipseShapeToAxes (shape, 20.0);
     psF64 sigma = axes.major;
+
+    if (PAR[PM_PAR_7] == 0.0) { 
+        psF32 z =  powf(PAR[PM_PAR_I0] / flux - 1.0, 1.0 / ALPHA);
+        return ( sigma * sqrt (2.0 * z) );
+    }
 
     psF64 limit = flux / PAR[PM_PAR_I0];