Index: /trunk/psLib/src/imageops/psImageInterpolate.c
===================================================================
--- /trunk/psLib/src/imageops/psImageInterpolate.c	(revision 14825)
+++ /trunk/psLib/src/imageops/psImageInterpolate.c	(revision 14826)
@@ -7,6 +7,6 @@
  *  @author Paul Price, IfA
  *
- *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-06-05 22:04:28 $
+ *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-09-13 02:43:55 $
  *
  *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
@@ -312,9 +312,13 @@
         #define KERNEL_VARIANCE_CASE(TYPE) \
             case PS_TYPE_##TYPE: { \
+                double sumKernel2 = 0.0; /* Sum of kernel squares */ \
                 for (int j = 0, yPix = yCentral - (yNum - 1) / 2; j < yNum; j++, yPix++) { \
                     for (int i = 0, xPix = xCentral - (xNum - 1) / 2; i < xNum; i++, xPix++) { \
-                        *varianceValue += PS_SQR(kernel[j][i]) * variance->data.TYPE[yPix][xPix]; \
+                        double kernel2 = PS_SQR(kernel[j][i]); /* Kernel squared */ \
+                        sumKernel2 += kernel2; \
+                        *varianceValue += kernel2 * variance->data.TYPE[yPix][xPix]; \
                     } \
                 } \
+                *varianceValue /= sumKernel2; /* Normalise so that sum of kernel squares is unity */ \
                 break; \
             }
@@ -378,6 +382,6 @@
 // Interpolation engine for separable interpolation kernels (either for good reasons or for practical reasons)
 static psImageInterpolateStatus interpolateSeparate(double *imageValue, double *varianceValue,
-                                               psMaskType *maskValue, float x, float y,
-                                               const psImageInterpolateOptions *options)
+                                                    psMaskType *maskValue, float x, float y,
+                                                    const psImageInterpolateOptions *options)
 {
     // Parameters have been checked by psImageInterpolate()
@@ -542,11 +546,18 @@
         #define SEPARATE_VARIANCE_CASE(TYPE) \
           case PS_TYPE_##TYPE: { \
+            double ySumKernel2 = 0.0; /* Sum of kernel squared in y */ \
             for (int j = 0, yPix = yCentral - (yNum - 1) / 2; j < yNum; j++, yPix++) { \
                 double xInterpValue = 0.0; /* Interpolation in x */ \
+                double xSumKernel2 = 0.0; /* Sum of kernel squared in x */ \
                 for (int i = 0, xPix = xCentral - (xNum - 1) / 2; i < xNum; i++, xPix++) { \
-                    xInterpValue += PS_SQR(xKernel[i]) * variance->data.TYPE[yPix][xPix]; \
+                    double kernel2 = PS_SQR(xKernel[i]); /* Kernel squared */ \
+                    xSumKernel2 += kernel2; \
+                    xInterpValue += kernel2 * variance->data.TYPE[yPix][xPix]; \
                 } \
-                *varianceValue += xInterpValue * PS_SQR(yKernel[j]); /* Interpolating in y */ \
+                double kernel2 = PS_SQR(yKernel[j]); /* Kernel squared */ \
+                ySumKernel2 += xSumKernel2 * kernel2; \
+                *varianceValue += xInterpValue * kernel2; /* Interpolating in y */ \
             } \
+            *varianceValue /= ySumKernel2; /* Normalise so that sum of kernel squares is unity */ \
             break; \
           }