Changeset 16352 for trunk/psModules/src/imcombine/pmSubtraction.c

Timestamp:

Feb 6, 2008, 6:03:23 PM (18 years ago)

Author:

Paul Price

Message:

I no longer believe that I should be normalising the variance kernel.
This isn't done by Andy Becker in hotpants, and I can't see it in the
math:

If I(x,y) is the image, V(x,y) is the variance map and k(u,v) is the
kernel, then the convolved image is: sum_u,v I(x-u,y-v).k(u,v). Then
the variance of the convolved image is: sum_u,v V(x-u,y-v).k(u,v)^{2.
If you want to normalise by the sum of the kernel, then the convolved
image is: sum_u,v I(x-u,y-v).k(u,v)/sum_u,v k(u,v) and the variance
is: sum_u,v V(x-u,y-v).k(u,v)}2/(sum_u,v k(u,v))^{2. But, of course in
practise you absorb the normalisation into the kernel itself, so you
get the first form. This has no division by the sum of the square.}

See also psLib/src/imageops/psImageInterpolate.c

File:

• trunk/psModules/src/imcombine/pmSubtraction.c (modified) (2 diffs)

trunk/psModules/src/imcombine/pmSubtraction.c

@@ -4 +4 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.76 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2008-01-02 20:35:32 $
+ *  @version $Revision: 1.77 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2008-02-07 04:03:23 $
  *
  *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
@@ -50 +50 @@
 
     // Take the square of the normal kernel
-    double sumNormal = 0.0, sumVariance = 0.0; // Sum of the normal and variance kernels
     for (int v = yMin; v <= yMax; v++) {
         for (int u = xMin; u <= xMax; u++) {
-            float value = normalKernel->kernel[v][u]; // Value of interest
-            float value2 = PS_SQR(value); // Value squared
-            sumNormal += value;
-            sumVariance += value2;
-            out->kernel[v][u] = value2;
-        }
-    }
-
-    // Normalise so that the sum of the variance kernel is the square of the sum of the normal kernel
-    // This is required to keep the relative scaling between the image and the weight map
-    psBinaryOp(out->image, out->image, "*", psScalarAlloc(PS_SQR(sumNormal) / sumVariance, PS_TYPE_F32));
+            out->kernel[v][u] = PS_SQR(normalKernel->kernel[v][u]);
+        }
+    }
 
     return out;