From 2de9b5dc23488ba35ad14f5b13d00db0de5e325c Mon Sep 17 00:00:00 2001
From: David Rohr <drohr@jwdt.org>
Date: Tue, 8 Jul 2014 16:44:27 +0200
Subject: [PATCH] fix sum and sumsquare calculation: for the non FFT algorithm
sum was calculated twice before (stupid mistake), for the FFT algorithm there
is actually a dependency for the sumsquare calculation on whether
NumberPixels is odd or even. Although, the effect is very small for normal
maps since high order fourier coefficients are small then. Anyhow, should be
correct now
---
bioem.cpp | 10 ++++++++--
1 file changed, 8 insertions(+), 2 deletions(-)
diff --git a/bioem.cpp b/bioem.cpp
index fe47338..347215a 100644
--- a/bioem.cpp
+++ b/bioem.cpp
@@ -761,15 +761,22 @@ int bioem::createConvolutedProjectionMap(int iMap, int iConv, mycomplex_t* lproj
sumsquareC = 0;
if (FFTAlgo)
{
+ int jloopend = param.param_device.NumberFFTPixels1D;
+ if ((param.param_device.NumberPixels & 1) == 0) jloopend--;
for(int i = 0; i < param.param_device.NumberPixels; i++)
{
- for (int j = 1;j < param.param_device.NumberFFTPixels1D;j++)
+ for (int j = 1;j < jloopend;j++)
{
int k = i * param.param_device.NumberFFTPixels1D + j;
sumsquareC += (localmultFFT[k][0] * localmultFFT[k][0] + localmultFFT[k][1] * localmultFFT[k][1]) * 2;
}
int k = i * param.param_device.NumberFFTPixels1D;
sumsquareC += localmultFFT[k][0] * localmultFFT[k][0] + localmultFFT[k][1] * localmultFFT[k][1];
+ if ((param.param_device.NumberPixels & 1) == 0)
+ {
+ k += param.param_device.NumberFFTPixels1D - 2;
+ sumsquareC += localmultFFT[k][0] * localmultFFT[k][0] + localmultFFT[k][1] * localmultFFT[k][1];
+ }
}
myfloat_t norm2 = (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels);
@@ -786,7 +793,6 @@ int bioem::createConvolutedProjectionMap(int iMap, int iConv, mycomplex_t* lproj
for(int i = 0; i < param.param_device.NumberPixels * param.param_device.NumberPixels; i++)
{
sumsquareC += Mapconv[i] * Mapconv[i];
- sumC += Mapconv[i];
}
myfloat_t norm2 = (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels);
--
GitLab