diff options
| -rw-r--r-- | Core/regularisers_CPU/ROF_TV_core.c | 12 | ||||
| -rwxr-xr-x | Core/regularisers_GPU/TV_ROF_GPU_core.cu | 18 |
2 files changed, 15 insertions, 15 deletions
diff --git a/Core/regularisers_CPU/ROF_TV_core.c b/Core/regularisers_CPU/ROF_TV_core.c index e89774f..100bf66 100644 --- a/Core/regularisers_CPU/ROF_TV_core.c +++ b/Core/regularisers_CPU/ROF_TV_core.c @@ -148,7 +148,7 @@ float D2_func(float *A, float *D2, int dimX, int dimY, int dimZ) for(j=0; j<dimY; j++) { for(i=0; i<dimX; i++) { for(k=0; k<dimZ; k++) { - index = (dimX*dimY)*k + j*dimX+i; + index = (dimX*dimY)*k + j*dimX+i; /* symmetric boundary conditions (Neuman) */ i1 = i + 1; if (i1 >= dimX) i1 = i-1; i2 = i - 1; if (i2 < 0) i2 = i+1; @@ -179,7 +179,7 @@ float D2_func(float *A, float *D2, int dimX, int dimY, int dimZ) #pragma omp parallel for shared (A, D2, dimX, dimY) private(i, j, i1, j1, i2, j2, NOMx_1,NOMy_1,NOMx_0,denom1,denom2,T2,index) for(j=0; j<dimY; j++) { for(i=0; i<dimX; i++) { - index = j*dimX+i; + index = j*dimX+i; /* symmetric boundary conditions (Neuman) */ i1 = i + 1; if (i1 >= dimX) i1 = i-1; i2 = i - 1; if (i2 < 0) i2 = i+1; @@ -265,8 +265,8 @@ float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambd dv2 = D2[index] - D2[(dimX*dimY)*k + j*dimX+i2]; dv3 = D3[index] - D3[(dimX*dimY)*k2 + j*dimX+i]; - B[index] = B[index] + tau*(lambda*(dv1 + dv2 + dv3) - (B[index] - A[index])); - }}} + B[index] += tau*(2.0f*lambda*(dv1 + dv2 + dv3) - (B[index] - A[index])); + }}} } else { #pragma omp parallel for shared (D1, D2, B, dimX, dimY) private(index, i, j, i1, j1, i2, j2,dv1,dv2) @@ -281,9 +281,9 @@ float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambd /* divergence components */ dv1 = D1[index] - D1[j2*dimX + i]; - dv2 = D2[index] - D2[j*dimX + i2]; + dv2 = D2[index] - D2[j*dimX + i2]; - B[index] = B[index] + tau*(lambda*(dv1 + dv2) - (B[index] - A[index])); + B[index] += tau*(2.0f*lambda*(dv1 + dv2) - (B[index] - A[index])); }} } return *B; diff --git a/Core/regularisers_GPU/TV_ROF_GPU_core.cu b/Core/regularisers_GPU/TV_ROF_GPU_core.cu index 67cdd5c..4610beb 100755 --- a/Core/regularisers_GPU/TV_ROF_GPU_core.cu +++ b/Core/regularisers_GPU/TV_ROF_GPU_core.cu @@ -94,7 +94,7 @@ __host__ __device__ int sign (float x) denom2 = 0.5f*(sign((float)NOMy_1) + sign((float)NOMy_0))*(MIN(abs((float)NOMy_1),abs((float)NOMy_0))); denom2 = denom2*denom2; T1 = sqrt(denom1 + denom2 + EPS); - D1[index] = NOMx_1/T1; + D1[index] = NOMx_1/T1; } } @@ -124,7 +124,7 @@ __host__ __device__ int sign (float x) denom2 = 0.5f*(sign((float)NOMx_1) + sign((float)NOMx_0))*(MIN(abs((float)NOMx_1),abs((float)NOMx_0))); denom2 = denom2*denom2; T2 = sqrt(denom1 + denom2 + EPS); - D2[index] = NOMy_1/T2; + D2[index] = NOMy_1/T2; } } @@ -139,15 +139,15 @@ __host__ __device__ int sign (float x) if ((i >= 0) && (i < (N)) && (j >= 0) && (j < (M))) { - /* boundary conditions (Neumann reflections) */ - i2 = i - 1; if (i2 < 0) i2 = i+1; + /* boundary conditions (Neumann reflections) */ + i2 = i - 1; if (i2 < 0) i2 = i+1; j2 = j - 1; if (j2 < 0) j2 = j+1; /* divergence components */ dv1 = D1[index] - D1[j2*N + i]; - dv2 = D2[index] - D2[j*N + i2]; + dv2 = D2[index] - D2[j*N + i2]; - Update[index] = Update[index] + tau*(lambda*(dv1 + dv2) - (Update[index] - Input[index])); + Update[index] += tau*(2.0f*lambda*(dv1 + dv2) - (Update[index] - Input[index])); } } @@ -268,7 +268,7 @@ __host__ __device__ int sign (float x) denom3 = 0.5*(sign(NOMy_1) + sign(NOMy_0))*(MIN(abs(NOMy_1),abs(NOMy_0))); denom3 = denom3*denom3; T3 = sqrt(denom1 + denom2 + denom3 + EPS); - D3[index] = NOMz_1/T3; + D3[index] = NOMz_1/T3; } } @@ -297,7 +297,7 @@ __host__ __device__ int sign (float x) dv2 = D2[index] - D2[(dimX*dimY)*k + j*dimX+i2]; dv3 = D3[index] - D3[(dimX*dimY)*k2 + j*dimX+i]; - Update[index] = Update[index] + tau*(lambda*(dv1 + dv2 + dv3) - (Update[index] - Input[index])); + Update[index] += tau*(2.0f*lambda*(dv1 + dv2 + dv3) - (Update[index] - Input[index])); } } @@ -341,7 +341,7 @@ extern "C" void TV_ROF_GPU_main(float* Input, float* Output, float lambdaPar, in CHECK(cudaDeviceSynchronize()); } - CHECK(cudaFree(d_D3)); + CHECK(cudaFree(d_D3)); } else { // TV - 2D case |