summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
Diffstat
-rwxr-xr-xbuild/run.sh4
-rw-r--r--src/Core/regularisers_CPU/Diffusion_core.c70
-rw-r--r--src/Core/regularisers_GPU/NonlDiff_GPU_core.cu57
3 files changed, 117 insertions, 14 deletions
diff --git a/build/run.sh b/build/run.sh
index f0b3d3e..91c5f05 100755
--- a/build/run.sh
+++ b/build/run.sh
@@ -6,11 +6,11 @@ rm -r ../build_proj
mkdir ../build_proj
cd ../build_proj/
#make clean
-export CIL_VERSION=19.10
+export CIL_VERSION=19.10.1
# install Python modules without CUDA
#cmake ../ -DBUILD_PYTHON_WRAPPER=ON -DBUILD_MATLAB_WRAPPER=OFF -DBUILD_CUDA=OFF -DCMAKE_BUILD_TYPE=Release -DCMAKE_INSTALL_PREFIX=./install
# install Python modules with CUDA
-cmake ../ -DBUILD_PYTHON_WRAPPER=ON -DBUILD_MATLAB_WRAPPER=OFF -DBUILD_CUDA=ON -DCMAKE_BUILD_TYPE=Release -DCMAKE_INSTALL_PREFIX=./install
+ cmake ../ -DBUILD_PYTHON_WRAPPER=ON -DBUILD_MATLAB_WRAPPER=OFF -DBUILD_CUDA=ON -DCMAKE_BUILD_TYPE=Release -DCMAKE_INSTALL_PREFIX=./install
# install Matlab modules without CUDA
# cmake ../ -DBUILD_PYTHON_WRAPPER=OFF -DMatlab_ROOT_DIR=/dls_sw/apps/matlab/r2014a/ -DBUILD_MATLAB_WRAPPER=ON -DBUILD_CUDA=OFF -DCMAKE_BUILD_TYPE=Release -DCMAKE_INSTALL_PREFIX=./install
# install Matlab modules with CUDA
diff --git a/src/Core/regularisers_CPU/Diffusion_core.c b/src/Core/regularisers_CPU/Diffusion_core.c
index 1a05c39..fdd60de 100644
--- a/src/Core/regularisers_CPU/Diffusion_core.c
+++ b/src/Core/regularisers_CPU/Diffusion_core.c
@@ -179,10 +179,10 @@ float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaP
}
else if (penaltytype == 2) {
/* Perona-Malik */
- e1 = (e1)/(1.0f + powf((e1/sigmaPar),2));
- w1 = (w1)/(1.0f + powf((w1/sigmaPar),2));
- n1 = (n1)/(1.0f + powf((n1/sigmaPar),2));
- s1 = (s1)/(1.0f + powf((s1/sigmaPar),2));
+ e1 /= (1.0f + powf((e1/sigmaPar),2));
+ w1 /= (1.0f + powf((w1/sigmaPar),2));
+ n1 /= (1.0f + powf((n1/sigmaPar),2));
+ s1 /= (1.0f + powf((s1/sigmaPar),2));
}
else if (penaltytype == 3) {
/* Tukey Biweight */
@@ -205,8 +205,32 @@ float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaP
if (fabs(n1) > sigmaPar) n1 = 0.0f;
if (fabs(s1) > sigmaPar) s1 = 0.0f;
}
+ else if (penaltytype == 5) {
+ /*
+ Threshold constrained Huber diffusion
+ */
+ if (fabs(e1) <= 2.0f*sigmaPar) {
+ if (fabs(e1) > sigmaPar) e1 = signNDFc(e1);
+ else e1 = e1/sigmaPar; }
+ else e1 = 0.0f;
+
+ if (fabs(w1) <= 2.0f*sigmaPar) {
+ if (fabs(w1) > sigmaPar) w1 = signNDFc(w1);
+ else w1 = w1/sigmaPar; }
+ else w1 = 0.0f;
+
+ if (fabs(n1) <= 2.0f*sigmaPar) {
+ if (fabs(n1) > sigmaPar) n1 = signNDFc(n1);
+ else n1 = n1/sigmaPar; }
+ else n1 = 0.0f;
+
+ if (fabs(s1) <= 2.0f*sigmaPar) {
+ if (fabs(s1) > sigmaPar) s1 = signNDFc(s1);
+ else s1 = s1/sigmaPar; }
+ else s1 = 0.0f;
+ }
else {
- printf("%s \n", "No penalty function selected! Use 1,2,3 or 4.");
+ printf("%s \n", "No penalty function selected! Use 1,2,3,4 or 5.");
break;
}
Output[index] += tau*(lambdaPar*(e1 + w1 + n1 + s1) - (Output[index] - Input[index]));
@@ -344,8 +368,42 @@ float NonLinearDiff3D(float *Input, float *Output, float lambdaPar, float sigmaP
if (fabs(u1) > sigmaPar) u1 = 0.0f;
if (fabs(d1) > sigmaPar) d1 = 0.0f;
}
+ else if (penaltytype == 5) {
+ /*
+ Threshold constrained Huber diffusion
+ */
+ if (fabs(e1) <= 2.0f*sigmaPar) {
+ if (fabs(e1) > sigmaPar) e1 = signNDFc(e1);
+ else e1 = e1/sigmaPar; }
+ else e1 = 0.0f;
+
+ if (fabs(w1) <= 2.0f*sigmaPar) {
+ if (fabs(w1) > sigmaPar) w1 = signNDFc(w1);
+ else w1 = w1/sigmaPar; }
+ else w1 = 0.0f;
+
+ if (fabs(n1) <= 2.0f*sigmaPar) {
+ if (fabs(n1) > sigmaPar) n1 = signNDFc(n1);
+ else n1 = n1/sigmaPar; }
+ else n1 = 0.0f;
+
+ if (fabs(s1) <= 2.0f*sigmaPar) {
+ if (fabs(s1) > sigmaPar) s1 = signNDFc(s1);
+ else s1 = s1/sigmaPar; }
+ else s1 = 0.0f;
+
+ if (fabs(u1) <= 2.0f*sigmaPar) {
+ if (fabs(u1) > sigmaPar) u1 = signNDFc(u1);
+ else u1 = u1/sigmaPar; }
+ else u1 = 0.0f;
+
+ if (fabs(d1) <= 2.0f*sigmaPar) {
+ if (fabs(d1) > sigmaPar) d1 = signNDFc(d1);
+ else d1 = d1/sigmaPar; }
+ else d1 = 0.0f;
+ }
else {
- printf("%s \n", "No penalty function selected! Use 1,2,3 or 4.");
+ printf("%s \n", "No penalty function selected! Use 1,2,3,4 or 5.");
break;
}
diff --git a/src/Core/regularisers_GPU/NonlDiff_GPU_core.cu b/src/Core/regularisers_GPU/NonlDiff_GPU_core.cu
index a8876cd..0d34cf8 100644
--- a/src/Core/regularisers_GPU/NonlDiff_GPU_core.cu
+++ b/src/Core/regularisers_GPU/NonlDiff_GPU_core.cu
@@ -160,27 +160,37 @@ __global__ void LinearDiff2D_kernel(float *Input, float *Output, float lambdaPar
This means that the linear diffusion will be performed on pixels with
absolute difference less than the threshold.
*/
- if (abs(e1) <= 1.5f*sigmaPar) {
+ if (abs(e1) > sigmaPar) e1 = 0.0f;
+ if (abs(w1) > sigmaPar) w1 = 0.0f;
+ if (abs(n1) > sigmaPar) n1 = 0.0f;
+ if (abs(s1) > sigmaPar) s1 = 0.0f;
+ }
+ else if (penaltytype == 5) {
+ /* Threshold-constrained Huber nonlinear diffusion
+ This means that the linear diffusion will be performed on pixels with
+ absolute difference less than the threshold.
+ */
+ if (abs(e1) <= 2.0f*sigmaPar) {
if (abs(e1) > sigmaPar) e1 = signNDF(e1);
else e1 = e1/sigmaPar;}
else e1 = 0.0f;
- if (abs(w1) <= 1.5f*sigmaPar) {
+ if (abs(w1) <= 2.0f*sigmaPar) {
if (abs(w1) > sigmaPar) w1 = signNDF(w1);
else w1 = w1/sigmaPar;}
else w1 = 0.0f;
- if (abs(n1) <= 1.5f*sigmaPar) {
+ if (abs(n1) <= 2.0f*sigmaPar) {
if (abs(n1) > sigmaPar) n1 = signNDF(n1);
else n1 = n1/sigmaPar; }
else n1 = 0.0f;
- if (abs(s1) <= 1.5f*sigmaPar) {
+ if (abs(s1) <= 2.0f*sigmaPar) {
if (abs(s1) > sigmaPar) s1 = signNDF(s1);
else s1 = s1/sigmaPar; }
else s1 = 0.0f;
}
- else printf("%s \n", "No penalty function selected! Use 1,2,3 or 4.");
+ else printf("%s \n", "No penalty function selected! Use 1,2,3, 4 or 5.");
Output[index] += tau*(lambdaPar*(e1 + w1 + n1 + s1) - (Output[index] - Input[index]));
}
@@ -318,7 +328,42 @@ __global__ void NonLinearDiff3D_kernel(float *Input, float *Output, float lambda
if (abs(u1) > sigmaPar) u1 = 0.0f;
if (abs(d1) > sigmaPar) d1 = 0.0f;
}
- else printf("%s \n", "No penalty function selected! Use 1,2,3 or 4.");
+ else if (penaltytype == 5) {
+ /* Threshold-constrained Huber nonlinear diffusion
+ This means that the linear diffusion will be performed on pixels with
+ absolute difference less than the threshold.
+ */
+ if (abs(e1) <= 2.0f*sigmaPar) {
+ if (abs(e1) > sigmaPar) e1 = signNDF(e1);
+ else e1 = e1/sigmaPar;}
+ else e1 = 0.0f;
+
+ if (abs(w1) <= 2.0f*sigmaPar) {
+ if (abs(w1) > sigmaPar) w1 = signNDF(w1);
+ else w1 = w1/sigmaPar;}
+ else w1 = 0.0f;
+
+ if (abs(n1) <= 2.0f*sigmaPar) {
+ if (abs(n1) > sigmaPar) n1 = signNDF(n1);
+ else n1 = n1/sigmaPar; }
+ else n1 = 0.0f;
+
+ if (abs(s1) <= 2.0f*sigmaPar) {
+ if (abs(s1) > sigmaPar) s1 = signNDF(s1);
+ else s1 = s1/sigmaPar; }
+ else s1 = 0.0f;
+
+ if (abs(u1) <= 2.0f*sigmaPar) {
+ if (abs(u1) > sigmaPar) u1 = signNDF(u1);
+ else u1 = u1/sigmaPar; }
+ else u1 = 0.0f;
+
+ if (abs(d1) <= 2.0f*sigmaPar) {
+ if (abs(d1) > sigmaPar) d1 = signNDF(d1);
+ else d1 = d1/sigmaPar; }
+ else d1 = 0.0f;
+ }
+ else printf("%s \n", "No penalty function selected! Use 1,2,3,4, or 5.");
Output[index] += tau*(lambdaPar*(e1 + w1 + n1 + s1 + u1 + d1) - (Output[index] - Input[index]));
}
}
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-08 16:40:53 +0000