Changeset 25860 for branches/pap/psModules/src/imcombine/pmSubtractionKernels.c

Timestamp:

Oct 15, 2009, 4:58:42 PM (17 years ago)

Author:

Paul Price

Message:

Fixed implementation of chi^{2 penalties. Using penalties going as r}4 (following Yuan & Akerlof 2008).

File:

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

branches/pap/psModules/src/imcombine/pmSubtractionKernels.c

@@ -94 +94 @@
             kernels->v->data.S32[index] = v;
             kernels->preCalc->data[index] = NULL;
-            kernels->penalties->data.F32[index] = kernels->penalty * (PS_SQR(u) + PS_SQR(v));
+            kernels->penalties->data.F32[index] = kernels->penalty * PS_SQR(PS_SQR(u) + PS_SQR(v));
 
             psTrace("psModules.imcombine", 7, "Kernel %d: %d %d\n", index, u, v);
@@ -151 +151 @@
                         double value = xKernel->data.F32[x] * yKernel->data.F32[y]; // Value of kernel
                         kernel->kernel[v][u] = value;
-                        moment += value * (PS_SQR(u) + PS_SQR(v));
+                        moment += value * PS_SQR((PS_SQR(u) + PS_SQR(v)));
                     }
                 }
@@ -190 +190 @@
                 kernels->preCalc->data[index] = preCalc;
                 kernels->penalties->data.F32[index] = kernels->penalty * fabsf(moment);
+                fprintf(stderr, "Penalty %d: %f\n", index, kernels->penalty * fabsf(moment));
 
                 psTrace("psModules.imcombine", 7, "Kernel %d: %f %d %d %f\n", index,
@@ -592 +593 @@
                                     poly->data.F32[j] = polyVal;
                                     norm += polyVal;
-                                    moment += polyVal * (PS_SQR(u) + PS_SQR(v));
+                                    moment += polyVal * PS_SQR(PS_SQR(u) + PS_SQR(v));
 
                                     psVectorExtend(uCoords, RINGS_BUFFER, 1);
@@ -618 +619 @@
                         psBinaryOp(poly, poly, "*", psScalarAlloc(1.0 / norm, PS_TYPE_F32));
                     }
-//                    moment /= norm;
+                    moment /= norm;
                 }