Index: trunk/psLib/src/math/psStats.c
===================================================================
--- trunk/psLib/src/math/psStats.c	(revision 41494)
+++ trunk/psLib/src/math/psStats.c	(revision 41502)
@@ -324,4 +324,6 @@
     psF32 *output = outVector->data.F32; // Dereference the vector
 
+    if (maskData) psAssert (maskVector-> == inVector->n, "oops");
+
     long count = 0;                     // Number of valid entries
 
@@ -329,12 +331,13 @@
     // into the temporary vectors.
     for (long i = 0; i < inVector->n; i++) {
-        if (!isfinite(input[i]))
-            continue;
-        if (useRange && (input[i] < stats->min))
-            continue;
-        if (useRange && (input[i] > stats->max))
-            continue;
-        if (maskData && (maskData[i] & maskVal))
-            continue;
+	psAssert (count >= 0, "oops");
+	psAssert (count < outVector->n, "oops");
+	psAssert (i >= 0, "oops");
+	psAssert (i < inVector->n, "oops");
+
+        if (!isfinite(input[i])) continue;
+        if (useRange && (input[i] < stats->min)) continue;
+        if (useRange && (input[i] > stats->max)) continue;
+        if (maskData && (maskData[i] & maskVal)) continue;
 
         output[count] = input[i];
@@ -363,9 +366,21 @@
 
     // Calculate the median.  Use the average if the number of samples if even.
+    int midPt = (count/2);
+    psAssert (midPt >=           0, "oops");
+    psAssert (midPt < outVector->n, "oops");
     if (count % 2 == 0) {
-        stats->sampleMedian = 0.5 * (output[(count/2) - 1] + output[count/2]);
+	psAssert ((midPt - 1) >=           0, "oops");
+	psAssert ((midPt - 1) < outVector->n, "oops");
+        stats->sampleMedian = 0.5 * (output[midPt - 1] + output[midPt]);
     } else {
-        stats->sampleMedian = output[count/2];
-    }
+        stats->sampleMedian = output[midPt];
+    }
+
+    int Qmin = (int)0.25*count;
+    int Qmax = (int)0.75*count;
+    psAssert (Qmin >=           0, "oops");
+    psAssert (Qmin < outVector->n, "oops");
+    psAssert (Qmax >=           0, "oops");
+    psAssert (Qmax < outVector->n, "oops");
 
     // Calculate the quartile points exactly.