Changeset 14420 for trunk/psModules/src/imcombine/pmSubtractionKernels.c

Timestamp:

Aug 7, 2007, 9:02:55 AM (19 years ago)

Author:

Paul Price

Message:

Need to normalise sum of kernel components to unity in order to preserve photometric scaling.

File:

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

trunk/psModules/src/imcombine/pmSubtractionKernels.c

@@ -147 +147 @@
                 // Set the pre-calculated kernel
                 psKernel *preCalc = psKernelAlloc(-size, size, -size, size);
+                double sum = 0.0;       // Normalisation
                 for (int v = -size; v <= size; v++) {
                     for (int u = -size; u <= size; u++) {
-                        preCalc->kernel[v][u] = power(u, uOrder) * power(v, vOrder) *
+                        sum += preCalc->kernel[v][u] = power(u, uOrder) * power(v, vOrder) *
                             expf(-0.5 * (PS_SQR(u) + PS_SQR(v)) / PS_SQR(sigmas->data.F32[i]));
                     }
                 }
+                // Normalise sum of kernel component to unity
+                psBinaryOp(preCalc->image, preCalc->image, "*", psScalarAlloc(1.0/sum, PS_TYPE_F32));
 
                 // Iterate over spatial order.  This loop creates the terms for
@@ -450 +453 @@
                 // Set the pre-calculated kernel
                 psKernel *preCalc = psKernelAlloc(-size, size, -size, size);
+                double sum = 0.0;       // Normalisation
                 for (int v = -size; v <= size; v++) {
                     for (int u = -size; u <= size; u++) {
-                        preCalc->kernel[v][u] = power(u, uOrder) * power(v, vOrder) *
+                        sum += preCalc->kernel[v][u] = power(u, uOrder) * power(v, vOrder) *
                             expf(-0.5 * (PS_SQR(u) + PS_SQR(v)) / PS_SQR(sigmas->data.F32[i]));
                     }
                 }
+                // Normalise sum of kernel component to unity
+                psBinaryOp(preCalc->image, preCalc->image, "*", psScalarAlloc(1.0/sum, PS_TYPE_F32));
 
                 // Iterate over spatial order.  This loop creates the terms for