Index: trunk/psModules/src/imcombine/pmSubtractionMatch.c
===================================================================
--- trunk/psModules/src/imcombine/pmSubtractionMatch.c	(revision 23296)
+++ trunk/psModules/src/imcombine/pmSubtractionMatch.c	(revision 23307)
@@ -90,4 +90,105 @@
 }
 
+// Check input arguments
+static bool subtractionMatchCheck(pmReadout *conv1, pmReadout *conv2, // Convolved images
+                                  const pmReadout *ro1, const pmReadout *ro2, // Input images
+                                  int stride, // Size for convolution patches
+                                  float sysError, // Relative systematic error
+                                  psImageMaskType maskVal, // Value to mask for input
+                                  psImageMaskType maskBad, // Mask for output bad pixels
+                                  psImageMaskType maskPoor, // Mask for output poor pixels
+                                  float poorFrac, // Fraction for "poor"
+                                  float badFrac,   // Maximum fraction of bad input pixels to accept
+                                  pmSubtractionMode subMode // Mode of subtraction
+                                  )
+{
+    if (subMode != PM_SUBTRACTION_MODE_2) {
+        PM_ASSERT_READOUT_NON_NULL(conv1, false);
+        if (conv1->image) {
+            psFree(conv1->image);
+            conv1->image = NULL;
+        }
+        if (conv1->mask) {
+            psFree(conv1->mask);
+            conv1->mask = NULL;
+        }
+        if (conv1->variance) {
+            psFree(conv1->variance);
+            conv1->variance = NULL;
+        }
+    }
+    if (subMode != PM_SUBTRACTION_MODE_1) {
+        PM_ASSERT_READOUT_NON_NULL(conv2, false);
+        if (conv2->image) {
+            psFree(conv2->image);
+            conv2->image = NULL;
+        }
+        if (conv2->mask) {
+            psFree(conv2->mask);
+            conv2->mask = NULL;
+        }
+        if (conv2->variance) {
+            psFree(conv2->variance);
+            conv2->variance = NULL;
+        }
+    }
+
+    PM_ASSERT_READOUT_NON_NULL(ro1, false);
+    PM_ASSERT_READOUT_NON_NULL(ro2, false);
+    PM_ASSERT_READOUT_IMAGE(ro1, false);
+    PM_ASSERT_READOUT_IMAGE(ro2, false);
+    PS_ASSERT_IMAGES_SIZE_EQUAL(ro1->image, ro2->image, false);
+    PS_ASSERT_INT_NONNEGATIVE(stride, false);
+    if (isfinite(sysError)) {
+        PS_ASSERT_FLOAT_LARGER_THAN_OR_EQUAL(sysError, 0.0, false);
+        PS_ASSERT_FLOAT_LESS_THAN(sysError, 1.0, false);
+    }
+    // Don't care about maskVal
+    // Don't care about maskBad
+    // Don't care about maskPoor
+    PS_ASSERT_FLOAT_LARGER_THAN(poorFrac, 0.0, false);
+    PS_ASSERT_FLOAT_LESS_THAN_OR_EQUAL(poorFrac, 1.0, false);
+    if (isfinite(badFrac)) {
+        PS_ASSERT_FLOAT_LARGER_THAN(badFrac, 0.0, false);
+        PS_ASSERT_FLOAT_LESS_THAN_OR_EQUAL(badFrac, 1.0, false);
+    }
+
+    return true;
+}
+
+
+
+bool pmSubtractionMatchPrecalc(pmReadout *conv1, pmReadout *conv2, const pmReadout *ro1, const pmReadout *ro2,
+                               pmSubtractionKernels *kernel, psRegion *region, int stride, float sysError,
+                               psImageMaskType maskVal, psImageMaskType maskBad, psImageMaskType maskPoor,
+                               float poorFrac, float badFrac)
+{
+    PM_ASSERT_SUBTRACTION_KERNELS_NON_NULL(kernel, false);
+    pmSubtractionMode subMode = kernel->mode; // Subtraction mode: which input to convolve
+    if (!subtractionMatchCheck(conv1, conv2, ro1, ro2, stride, sysError, maskVal, maskBad, maskPoor,
+                               poorFrac, badFrac, subMode)) {
+        return false;
+    }
+
+    psImage *subMask = pmSubtractionMask(ro1->mask, ro2 ? ro2->mask : NULL, maskVal, kernel->size, 0,
+                                         badFrac, useFFT); // Subtraction mask
+    if (!subMask) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to generate subtraction mask.");
+        return false;
+    }
+
+    psTrace("psModules.imcombine", 2, "Convolving...\n");
+    if (!pmSubtractionConvolve(conv1, conv2, ro1, ro2, subMask, stride, maskBad, maskPoor, poorFrac,
+                               sysError, region, kernel, true, useFFT)) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to convolve image.");
+        psFree(subMask);
+        return false;
+    }
+
+    psFree(subMask);
+
+    return true;
+}
+
 
 bool pmSubtractionMatch(pmReadout *conv1, pmReadout *conv2, const pmReadout *ro1, const pmReadout *ro2,
@@ -101,43 +202,10 @@
                         psImageMaskType maskPoor, float poorFrac, float badFrac, pmSubtractionMode subMode)
 {
-    if (subMode != PM_SUBTRACTION_MODE_2) {
-        PM_ASSERT_READOUT_NON_NULL(conv1, false);
-        if (conv1->image) {
-            psFree(conv1->image);
-            conv1->image = NULL;
-        }
-        if (conv1->mask) {
-            psFree(conv1->mask);
-            conv1->mask = NULL;
-        }
-        if (conv1->variance) {
-            psFree(conv1->variance);
-            conv1->variance = NULL;
-        }
-    }
-    if (subMode != PM_SUBTRACTION_MODE_1) {
-        PM_ASSERT_READOUT_NON_NULL(conv2, false);
-        if (conv2->image) {
-            psFree(conv2->image);
-            conv2->image = NULL;
-        }
-        if (conv2->mask) {
-            psFree(conv2->mask);
-            conv2->mask = NULL;
-        }
-        if (conv2->variance) {
-            psFree(conv2->variance);
-            conv2->variance = NULL;
-        }
-    }
-
-    PM_ASSERT_READOUT_NON_NULL(ro1, false);
-    PM_ASSERT_READOUT_NON_NULL(ro2, false);
-    PM_ASSERT_READOUT_IMAGE(ro1, false);
-    PM_ASSERT_READOUT_IMAGE(ro2, false);
-    PS_ASSERT_IMAGES_SIZE_EQUAL(ro1->image, ro2->image, false);
+    if (!subtractionMatchCheck(conv1, conv2, ro1, ro2, stride, sysError, maskVal, maskBad, maskPoor,
+                               poorFrac, badFrac, subMode)) {
+        return false;
+    }
 
     PS_ASSERT_INT_NONNEGATIVE(footprint, false);
-    PS_ASSERT_INT_NONNEGATIVE(stride, false);
     // regionSize can be just about anything (except maybe negative, but it can be NAN)
     PS_ASSERT_FLOAT_LARGER_THAN(stampSpacing, 0.0, false);
@@ -168,17 +236,4 @@
     PS_ASSERT_INT_POSITIVE(iter, false);
     PS_ASSERT_FLOAT_LARGER_THAN(rej, 0.0, false);
-    if (isfinite(sysError)) {
-        PS_ASSERT_FLOAT_LARGER_THAN_OR_EQUAL(sysError, 0.0, NULL);
-        PS_ASSERT_FLOAT_LESS_THAN(sysError, 1.0, NULL);
-    }
-    // Don't care about maskVal
-    // Don't care about maskBad
-    // Don't care about maskPoor
-    PS_ASSERT_FLOAT_LARGER_THAN(poorFrac, 0.0, NULL);
-    PS_ASSERT_FLOAT_LESS_THAN_OR_EQUAL(poorFrac, 1.0, NULL);
-    if (isfinite(badFrac)) {
-        PS_ASSERT_FLOAT_LARGER_THAN(badFrac, 0.0, NULL);
-        PS_ASSERT_FLOAT_LESS_THAN_OR_EQUAL(badFrac, 1.0, NULL);
-    }
 
     // If the stamp footprint is smaller than the kernel size, then we won't get much signal in the outer
Index: trunk/psModules/src/imcombine/pmSubtractionMatch.h
===================================================================
--- trunk/psModules/src/imcombine/pmSubtractionMatch.h	(revision 23296)
+++ trunk/psModules/src/imcombine/pmSubtractionMatch.h	(revision 23307)
@@ -48,4 +48,20 @@
     );
 
+/// Match two images using precalculated kernel
+bool pmSubtractionMatchPrecalc(pmReadout *conv1, ///< Output convolved data for image 1
+                               pmReadout *conv2, ///< Output convolved data for image 2
+                               const pmReadout *ro1, ///< Image 1
+                               const pmReadout *ro2, ///< Image 2
+                               pmSubtractionKernels *kernels, ///< Pre-calculated matching kernel
+                               psRegion *region, ///< Region of applicability for kernel
+                               int stride, ///< Size for convolution patches
+                               float sysError, ///< Relative systematic error
+                               psImageMaskType maskVal, ///< Value to mask for input
+                               psImageMaskType maskBad, ///< Mask for output bad pixels
+                               psImageMaskType maskPoor, ///< Mask for output poor pixels
+                               float poorFrac, ///< Fraction for "poor"
+                               float badFrac ///< Maximum fraction of bad input pixels to accept
+    );
+
 /// Execute a thread job to measure the PSF width ratios
 bool pmSubtractionOrderThread(psThreadJob *job ///< Job to execute