Index: /trunk/psModules/src/imcombine/pmSubtraction.c
===================================================================
--- /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 13388)
+++ /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 13389)
@@ -4,6 +4,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-05-15 04:33:58 $
+ *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-05-16 00:01:03 $
  *
  *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
@@ -66,4 +66,10 @@
     assert(kernels->xOrder->n == numKernels);
     assert(kernels->yOrder->n == numKernels);
+
+    // Ensure the subIndex for POIS kernels is what is expected
+    assert(kernels->type != PM_SUBTRACTION_KERNEL_POIS ||
+           (kernels->u->data.S32[kernels->subIndex] == 0 && kernels->v->data.S32[kernels->subIndex] == 0 &&
+            kernels->xOrder->data.S32[kernels->subIndex] == 0 &&
+            kernels->yOrder->data.S32[kernels->subIndex] == 0));
 
     int size = kernels->size;           // Kernel half-size
@@ -82,4 +88,8 @@
               int v = kernels->v->data.S32[i]; // Offset in y
               kernel->kernel[v][u] += weightFunc(polyValue * solution->data.F64[i]);
+              if (i != kernels->subIndex) {
+                  // The (0,0) element is subtracted from most kernels to preserve photometric scaling
+                  kernel->kernel[0][0] -= weightFunc(solution->data.F64[i]);
+              }
               break;
           }
@@ -119,5 +129,10 @@
           int u = kernels->u->data.S32[index]; // Offset in x
           int v = kernels->v->data.S32[index]; // Offset in y
-          return weightFunc(polyValue) * image->data.F32[y + v][x + u];
+          double value = weightFunc(polyValue) * image->data.F32[y + v][x + u]; // Value of convolution
+          if (index != kernels->subIndex) {
+              // The (0,0) element is subtracted from most kernels to preserve photometric scaling
+              value -= image->data.F32[y][x];
+          }
+          return value;
       }
       case PM_SUBTRACTION_KERNEL_ISIS: {
@@ -243,15 +258,4 @@
                         convolutions->data.F64[i] = convolvePixel(kernels, i, x, y, reference, polyValues,
                                                                   imageWeighting);
-                    }
-
-                    if (spatialOrder > 0) {
-                        // Subtract a nominated convolution, to force flux conservation
-                        int subIndex = kernels->subIndex; // Index of convolution to subtract
-                        for (int i = 0; i < subIndex; i++) {
-                            convolutions->data.F64[i] -= convolutions->data.F64[subIndex];
-                        }
-                        for (int i = subIndex + 1; i < numKernels; i++) {
-                            convolutions->data.F64[i] -= convolutions->data.F64[subIndex];
-                        }
                     }
 
Index: /trunk/psModules/src/imcombine/pmSubtractionKernels.c
===================================================================
--- /trunk/psModules/src/imcombine/pmSubtractionKernels.c	(revision 13388)
+++ /trunk/psModules/src/imcombine/pmSubtractionKernels.c	(revision 13389)
@@ -168,4 +168,5 @@
     }
 
+    // Subtract a particular kernel in order to preserve photometric calibration across image
     kernels->subIndex = 0;
     assert(kernels->u->data.S32[kernels->subIndex] == 0 &&
@@ -173,4 +174,14 @@
            kernels->xOrder->data.S32[kernels->subIndex] == 0 &&
            kernels->yOrder->data.S32[kernels->subIndex] == 0);
+
+    psKernel *subKernel = kernels->preCalc->data[kernels->subIndex]; // Kernel to subtract
+    for (int i = 0; i < kernels->subIndex; i++) {
+        psKernel *kernel = kernels->preCalc->data[i]; // Kernel of interest
+        (void)psBinaryOp(kernel->image, kernel->image, "-", subKernel->image);
+    }
+    for (int i = kernels->subIndex + 1; i < num; i++) {
+        psKernel *kernel = kernels->preCalc->data[i]; // Kernel of interest
+        (void)psBinaryOp(kernel->image, kernel->image, "-", subKernel->image);
+    }