Index: /trunk/psModules/src/imcombine/pmSubtraction.c
===================================================================
--- /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 13397)
+++ /trunk/psModules/src/imcombine/pmSubtraction.c	(revision 13398)
@@ -4,6 +4,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.9 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-05-16 03:53:36 $
+ *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-05-16 21:13:51 $
  *
  *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
@@ -54,9 +54,9 @@
 
 // Generate an image of the solved kernel
-static psKernel *solvedKernel(psKernel *kernel, // Kernel, to return
-                              const psVector *solution, // Solution to the least-squares problem
-                              const pmSubtractionKernels *kernels, // Kernel basis functions
-                              const psImage *polyValues, // Spatial polynomial values
-                              double (*weightFunc)(double value) // Function for weighting
+static inline psKernel *solvedKernel(psKernel *kernel, // Kernel, to return
+                                     const psVector *solution, // Solution to the least-squares problem
+                                     const pmSubtractionKernels *kernels, // Kernel basis functions
+                                     const psImage *polyValues, // Spatial polynomial values
+                                     double (*weightFunc)(double value) // Function for weighting
     )
 {
@@ -82,13 +82,15 @@
         int xOrder = kernels->xOrder->data.S32[i]; // Polynomial order in x
         int yOrder = kernels->yOrder->data.S32[i]; // Polynomial order in y
-        double polyValue = polyValues->data.F64[yOrder][xOrder]; // Value of spatial polynomial
+        float value = weightFunc(polyValues->data.F64[yOrder][xOrder] *
+                                 solution->data.F64[i]); // Value to sum
+        float subValue = weightFunc(-solution->data.F64[i]); // Value to subtract (actually add, because of -)
         switch (kernels->type) {
           case PM_SUBTRACTION_KERNEL_POIS: {
               int u = kernels->u->data.S32[i]; // Offset in x
               int v = kernels->v->data.S32[i]; // Offset in y
-              kernel->kernel[v][u] += weightFunc(polyValue * solution->data.F64[i]);
+              kernel->kernel[v][u] += value;
               if (kernels->spatialOrder > 0 && 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]);
+                  kernel->kernel[0][0] += subValue;
               }
               break;
@@ -99,10 +101,8 @@
               for (int v = -size; v <= size; v++) {
                   for (int u = -size; u <= size; u++) {
-                      kernel->kernel[v][u] += weightFunc(solution->data.F64[i] * polyValue *
-                                                         preCalc->kernel[v][u]);
+                      kernel->kernel[v][u] += value * weightFunc(preCalc->kernel[v][u]);
                       // The (0,0) kernel is subtracted from other kernels to preserve photometric scaling
                       if (kernels->spatialOrder > 0 && i != kernels->subIndex) {
-                          kernel->kernel[v][u] += weightFunc(-solution->data.F64[i] *
-                                                             subKernel->kernel[v][u]);
+                          kernel->kernel[v][u] += subValue * weightFunc(subKernel->kernel[v][u]);
                       }
                   }