Index: /trunk/psModules/src/imcombine/pmStack.c
===================================================================
--- /trunk/psModules/src/imcombine/pmStack.c	(revision 19481)
+++ /trunk/psModules/src/imcombine/pmStack.c	(revision 19482)
@@ -8,6 +8,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.38 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2008-09-03 19:42:54 $
+ *  @version $Revision: 1.39 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2008-09-11 03:37:58 $
  *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
  *
@@ -309,8 +309,9 @@
               }
           }
-          if (useVariance && safe && numIter > 0) {
+          if (useVariance && numIter > 0) {
               // Use variance to check that the two are consistent
               float diff = pixelData->data.F32[0] - pixelData->data.F32[1];
-              float sigma2 = pixelVariances->data.F32[0] + pixelVariances->data.F32[1];
+              float sigma2 = pixelVariances->data.F32[0] * varFactors->data.F32[pixelSources->data.U16[0]] +
+                  pixelVariances->data.F32[1] * varFactors->data.F32[pixelSources->data.U16[1]];
               if (PS_SQR(diff) > PS_SQR(rej) * sigma2) {
                   // Not consistent: mark both for inspection
@@ -341,5 +342,5 @@
                   float rej2 = PS_SQR(rej); // Rejection level squared
                   for (int i = 0; i < num; i++) {
-                      pixelVariances->data.F32[i] *= rej2 * varFactors->data.F32[i];
+                      pixelVariances->data.F32[i] *= rej2 * varFactors->data.F32[pixelSources->data.U16[i]];
                   }
               }
@@ -348,19 +349,17 @@
           // The clipping that follows is solely to identify suspect pixels.
           // These suspect pixels will be inspected in more detail by other functions.
-          int numClipped = INT_MAX;       // Number of pixels clipped per iteration
-          int totalClipped = 0;           // Total number of pixels clipped
+          int numClipped = INT_MAX;     // Number of pixels clipped per iteration
+          int totalClipped = 0;         // Total number of pixels clipped
           for (int i = 0; i < numIter && numClipped > 0 && num - totalClipped > 2; i++) {
               numClipped = 0;
               float median, stdev;    // Median and stdev of the combination, for rejection
 
-              if (!combinationMedianStdev(&median, useVariance ? NULL : &stdev, pixelData, pixelMasks, sort)) {
+              if (!combinationMedianStdev(&median, useVariance ? NULL : &stdev,
+                                          pixelData, pixelMasks, sort)) {
                   psWarning("Bad median/stdev at %d,%d", x, y);
                   break;
               }
 
-              float limit = NAN;        // Rejection limit, when rejecting based on data properties
-              if (!useVariance) {
-                  limit = rej * stdev;
-              }
+              float limit = rej * stdev; // Rejection limit, when rejecting based on data properties
 
 // Mask a pixel for inspection
@@ -378,4 +377,5 @@
                   if (useVariance) {
                       // Comparing squares --- cheaper than lots of sqrts
+                      // pixelVariances includes the variance factor and the rejection limit, from above
                       if (PS_SQR(diff) > pixelVariances->data.F32[j]) {
                           MASK_PIXEL_FOR_INSPECTION();
Index: /trunk/psModules/src/imcombine/pmStackReject.c
===================================================================
--- /trunk/psModules/src/imcombine/pmStackReject.c	(revision 19481)
+++ /trunk/psModules/src/imcombine/pmStackReject.c	(revision 19482)
@@ -295,7 +295,14 @@
     }
 
-    if (threaded) {
+    // Harvest the jobs
+    if (threaded) {
+        psThreadJob *job;                   // Job to destroy
+        while ((job = psThreadJobGetDone())) {
+            psFree(job);
+        }
+
         psMutexDestroy(source);
     }
+
     psFree(source);
     bad = psPixelsFromMask(bad, target, PM_STACK_MASK_ALL);
Index: /trunk/psModules/src/imcombine/pmSubtraction.c
===================================================================
--- /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 19481)
+++ /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 19482)
@@ -1141,5 +1141,5 @@
                     psMutexLock(subMask);
                 }
-                psArrayAdd(args, 1, (psImage*)subMask); // Casting away const
+                psArrayAdd(args, 1, subMask);
                 if (subMask) {
                     psMutexUnlock(subMask);
@@ -1166,11 +1166,22 @@
     }
 
-    if (!psThreadPoolWait(true)) {
+    if (!psThreadPoolWait(false)) {
         psError(PS_ERR_UNKNOWN, false, "Error waiting for threads.");
         return false;
     }
 
-    if (subMask && threaded) {
-        psMutexDestroy(subMask);
+    // We don't rely on psThreadPoolWait to harvest the jobs because the job contains a reference to the
+    // subMask, which is being changed on a thread, and psThreadPoolWait doesn't know that it needs to be
+    // locked before freeing.  After psThreadPoolWait, however, the jobs are completed, the threads are idle,
+    // and so there's no need to lock the subMask when we're blowing away the jobs.
+    if (threaded) {
+        psThreadJob *job;               // Completed job
+        while ((job = psThreadJobGetDone())) {
+            psFree(job);
+        }
+
+        if (subMask) {
+            psMutexDestroy(subMask);
+        }
     }
 
Index: /trunk/psModules/src/imcombine/pmSubtractionMatch.c
===================================================================
--- /trunk/psModules/src/imcombine/pmSubtractionMatch.c	(revision 19481)
+++ /trunk/psModules/src/imcombine/pmSubtractionMatch.c	(revision 19482)
@@ -98,5 +98,5 @@
 {
     if (region) {
-        varFactor *= (region->x1 - region->x1 + 1) * (region->y1 * region->y1 + 1) /
+        varFactor *= (region->x1 - region->x0 + 1) * (region->y1 - region->y0 + 1) /
             (ro->image->numCols * ro->image->numRows);
     }