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

Timestamp:

May 11, 2007, 5:38:32 PM (19 years ago)

Author:

Paul Price

Message:

Had left out raising numbers to the power zero.

File:

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

trunk/psModules/src/imcombine/pmSubtractionKernels.c

@@ -25 +25 @@
     )
 {
+    if (exp == 0) {
+        return 1.0;
+    }
     float result = value;               // Result to return
     for (int i = 2; i < exp; i++) {
@@ -82 +85 @@
     }
 
-    kernels->subIndex = PS_SQR(size + 1);
+    kernels->subIndex = num / 2;
     assert(kernels->u->data.S32[kernels->subIndex] == 0 &&
            kernels->v->data.S32[kernels->subIndex] == 0 &&
@@ -106 +109 @@
     int num = 0;                        // Number of basis functions
     for (int i = 0; i < numGaussians; i++) {
-        num += sigmas->data.F32[i] * orders->data.F32[i];
+        int gaussOrder = orders->data.S32[i]; // Polynomial order to apply to Gaussian
+        num += (gaussOrder + 1) * (gaussOrder + 2) / 2;
     }
     num *= (spatialOrder + 1) * (spatialOrder + 2) / 2;
@@ -124 +128 @@
                 psKernel *preCalc = psKernelAlloc(-size, size, -size, size);
                 for (int v = -size; v <= size; v++) {
-                    for (int u = -size; u < size; u++) {
+                    for (int u = -size; u <= size; u++) {
                         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]));
+                            expf(-0.5 * (PS_SQR(u) + PS_SQR(v)) / PS_SQR(sigmas->data.F32[i]));
                     }
                 }
@@ -147 +151 @@
     }
 
+    kernels->subIndex = (spatialOrder + 1) * (spatialOrder + 2) / 4;
+    assert(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);
+
     return kernels;
 }