summaryrefslogtreecommitdiffstats
path: root/Wrappers
diff options
context:
space:
mode:
Diffstat (limited to 'Wrappers')
-rw-r--r--Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m15
-rw-r--r--Wrappers/Matlab/demos/demoMatlab_denoise.m16
-rw-r--r--Wrappers/Matlab/mex_compile/regularisers_CPU/TGV.c29
-rw-r--r--Wrappers/Python/conda-recipe/conda_build_config.yaml2
-rwxr-xr-xWrappers/Python/conda-recipe/run_test.py12
-rw-r--r--Wrappers/Python/demos/demo_cpu_regularisers.py14
-rw-r--r--Wrappers/Python/demos/demo_cpu_regularisers3D.py60
-rw-r--r--Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py14
-rw-r--r--Wrappers/Python/demos/demo_gpu_regularisers.py12
-rw-r--r--Wrappers/Python/demos/demo_gpu_regularisers3D.py55
-rw-r--r--Wrappers/Python/src/cpu_regularisers.pyx35
-rw-r--r--Wrappers/Python/src/gpu_regularisers.pyx39
12 files changed, 227 insertions, 76 deletions
diff --git a/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m b/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m
index 5cc47b3..0c331a4 100644
--- a/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m
+++ b/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m
@@ -2,11 +2,13 @@
clear; close all
Path1 = sprintf(['..' filesep 'mex_compile' filesep 'installed'], 1i);
Path2 = sprintf(['..' filesep '..' filesep '..' filesep 'data' filesep], 1i);
+Path3 = sprintf(['..' filesep 'supp'], 1i);
addpath(Path1);
addpath(Path2);
+addpath(Path3);
N = 512;
-slices = 15;
+slices = 7;
vol3D = zeros(N,N,slices, 'single');
Ideal3D = zeros(N,N,slices, 'single');
Im = double(imread('lena_gray_512.tif'))/255; % loading image
@@ -131,7 +133,16 @@ figure; imshow(u_diff4(:,:,7), [0 1]); title('Diffusion 4thO denoised volume (CP
% fprintf('%s %f \n', 'RMSE error for Anis.Diff of 4th order is:', rmse_diff4);
% figure; imshow(u_diff4_g(:,:,7), [0 1]); title('Diffusion 4thO denoised volume (GPU)');
%%
-
+fprintf('Denoise using the TGV model (CPU) \n');
+lambda_TGV = 0.03; % regularisation parameter
+alpha1 = 1.0; % parameter to control the first-order term
+alpha0 = 2.0; % parameter to control the second-order term
+iter_TGV = 500; % number of Primal-Dual iterations for TGV
+tic; u_tgv = TGV(single(vol3D), lambda_TGV, alpha1, alpha0, iter_TGV); toc;
+rmseTGV = RMSE(Ideal3D(:),u_tgv(:));
+fprintf('%s %f \n', 'RMSE error for TGV is:', rmseTGV);
+figure; imshow(u_tgv(:,:,3), [0 1]); title('TGV denoised volume (CPU)');
+%%
%>>>>>>>>>>>>>> MULTI-CHANNEL priors <<<<<<<<<<<<<<< %
fprintf('Denoise a volume using the FGP-dTV model (CPU) \n');
diff --git a/Wrappers/Matlab/demos/demoMatlab_denoise.m b/Wrappers/Matlab/demos/demoMatlab_denoise.m
index 3506cca..14d3096 100644
--- a/Wrappers/Matlab/demos/demoMatlab_denoise.m
+++ b/Wrappers/Matlab/demos/demoMatlab_denoise.m
@@ -60,20 +60,20 @@ figure; imshow(u_sb, [0 1]); title('SB-TV denoised image (CPU)');
% figure; imshow(u_sbG, [0 1]); title('SB-TV denoised image (GPU)');
%%
fprintf('Denoise using the TGV model (CPU) \n');
-lambda_TGV = 0.04; % regularisation parameter
-alpha1 = 1; % parameter to control the first-order term
-alpha0 = 0.7; % parameter to control the second-order term
-iter_TGV = 500; % number of Primal-Dual iterations for TGV
+lambda_TGV = 0.045; % regularisation parameter
+alpha1 = 1.0; % parameter to control the first-order term
+alpha0 = 2.0; % parameter to control the second-order term
+iter_TGV = 2000; % number of Primal-Dual iterations for TGV
tic; u_tgv = TGV(single(u0), lambda_TGV, alpha1, alpha0, iter_TGV); toc;
rmseTGV = (RMSE(u_tgv(:),Im(:)));
fprintf('%s %f \n', 'RMSE error for TGV is:', rmseTGV);
figure; imshow(u_tgv, [0 1]); title('TGV denoised image (CPU)');
%%
% fprintf('Denoise using the TGV model (GPU) \n');
-% lambda_TGV = 0.04; % regularisation parameter
-% alpha1 = 1; % parameter to control the first-order term
-% alpha0 = 0.7; % parameter to control the second-order term
-% iter_TGV = 500; % number of Primal-Dual iterations for TGV
+% lambda_TGV = 0.045; % regularisation parameter
+% alpha1 = 1.0; % parameter to control the first-order term
+% alpha0 = 2.0; % parameter to control the second-order term
+% iter_TGV = 2000; % number of Primal-Dual iterations for TGV
% tic; u_tgv_gpu = TGV_GPU(single(u0), lambda_TGV, alpha1, alpha0, iter_TGV); toc;
% rmseTGV_gpu = (RMSE(u_tgv_gpu(:),Im(:)));
% fprintf('%s %f \n', 'RMSE error for TGV is:', rmseTGV_gpu);
diff --git a/Wrappers/Matlab/mex_compile/regularisers_CPU/TGV.c b/Wrappers/Matlab/mex_compile/regularisers_CPU/TGV.c
index 5459bf5..aa4eed4 100644
--- a/Wrappers/Matlab/mex_compile/regularisers_CPU/TGV.c
+++ b/Wrappers/Matlab/mex_compile/regularisers_CPU/TGV.c
@@ -21,14 +21,14 @@ limitations under the License.
#include "TGV_core.h"
/* C-OMP implementation of Primal-Dual denoising method for
- * Total Generilized Variation (TGV)-L2 model [1] (2D case only)
+ * Total Generilized Variation (TGV)-L2 model [1] (2D/3D)
*
* Input Parameters:
- * 1. Noisy image (2D) (required)
- * 2. lambda - regularisation parameter (required)
- * 3. parameter to control the first-order term (alpha1) (default - 1)
- * 4. parameter to control the second-order term (alpha0) (default - 0.5)
- * 5. Number of Chambolle-Pock (Primal-Dual) iterations (default is 300)
+ * 1. Noisy image/volume (2D/3D)
+ * 2. lambda - regularisation parameter
+ * 3. parameter to control the first-order term (alpha1)
+ * 4. parameter to control the second-order term (alpha0)
+ * 5. Number of Chambolle-Pock (Primal-Dual) iterations
* 6. Lipshitz constant (default is 12)
*
* Output:
@@ -44,7 +44,7 @@ void mexFunction(
{
int number_of_dims, iter;
- mwSize dimX, dimY;
+ mwSize dimX, dimY, dimZ;
const mwSize *dim_array;
float *Input, *Output=NULL, lambda, alpha0, alpha1, L2;
@@ -55,7 +55,7 @@ void mexFunction(
/*Handling Matlab input data*/
if ((nrhs < 2) || (nrhs > 6)) mexErrMsgTxt("At least 2 parameters is required, all parameters are: Image(2D), Regularisation parameter, alpha0, alpha1, iterations number, Lipshitz Constant");
- Input = (float *) mxGetData(prhs[0]); /*noisy image (2D) */
+ Input = (float *) mxGetData(prhs[0]); /*noisy image/volume */
lambda = (float) mxGetScalar(prhs[1]); /* regularisation parameter */
alpha1 = 1.0f; /* parameter to control the first-order term */
alpha0 = 0.5f; /* parameter to control the second-order term */
@@ -69,12 +69,15 @@ void mexFunction(
if (nrhs == 6) L2 = (float) mxGetScalar(prhs[5]); /* Lipshitz constant */
/*Handling Matlab output data*/
- dimX = dim_array[0]; dimY = dim_array[1];
+ dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2];
if (number_of_dims == 2) {
- Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL));
- /* running the function */
- TGV_main(Input, Output, lambda, alpha1, alpha0, iter, L2, dimX, dimY);
+ dimZ = 1; /*2D case*/
+ Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL));
}
- if (number_of_dims == 3) {mexErrMsgTxt("Only 2D images accepted");}
+ if (number_of_dims == 3) {
+ Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(3, dim_array, mxSINGLE_CLASS, mxREAL));
+ }
+ /* running the function */
+ TGV_main(Input, Output, lambda, alpha1, alpha0, iter, L2, dimX, dimY, dimZ);
}
diff --git a/Wrappers/Python/conda-recipe/conda_build_config.yaml b/Wrappers/Python/conda-recipe/conda_build_config.yaml
index b7977f3..d6a3915 100644
--- a/Wrappers/Python/conda-recipe/conda_build_config.yaml
+++ b/Wrappers/Python/conda-recipe/conda_build_config.yaml
@@ -2,6 +2,8 @@ python:
- 2.7 # [not win]
- 3.5
- 3.6
+ - 3.7
numpy:
- 1.12
+ - 1.14
- 1.15
diff --git a/Wrappers/Python/conda-recipe/run_test.py b/Wrappers/Python/conda-recipe/run_test.py
index cfb3f53..21f3216 100755
--- a/Wrappers/Python/conda-recipe/run_test.py
+++ b/Wrappers/Python/conda-recipe/run_test.py
@@ -2,7 +2,7 @@ import unittest
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, Diff4th
from PIL import Image
class TiffReader(object):
@@ -303,7 +303,7 @@ class TestRegularisers(unittest.TestCase):
'input' : u0,\
'regularisation_parameter':0.04, \
'alpha1':1.0,\
- 'alpha0':0.7,\
+ 'alpha0':2.0,\
'number_of_iterations' :250 ,\
'LipshitzConstant' :12 ,\
}
@@ -530,7 +530,7 @@ class TestRegularisers(unittest.TestCase):
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
# set parameters
- pars = {'algorithm' : DIFF4th, \
+ pars = {'algorithm' : Diff4th, \
'input' : u0,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -540,7 +540,7 @@ class TestRegularisers(unittest.TestCase):
print ("#############Diff4th CPU####################")
start_time = timeit.default_timer()
- diff4th_cpu = DIFF4th(pars['input'],
+ diff4th_cpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
@@ -555,7 +555,7 @@ class TestRegularisers(unittest.TestCase):
print ("##############Diff4th GPU##################")
start_time = timeit.default_timer()
try:
- diff4th_gpu = DIFF4th(pars['input'],
+ diff4th_gpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
@@ -565,7 +565,7 @@ class TestRegularisers(unittest.TestCase):
self.skipTest("Results not comparable. GPU computing error.")
rms = rmse(Im, diff4th_gpu)
pars['rmse'] = rms
- pars['algorithm'] = DIFF4th
+ pars['algorithm'] = Diff4th
txtstr = printParametersToString(pars)
txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
print (txtstr)
diff --git a/Wrappers/Python/demos/demo_cpu_regularisers.py b/Wrappers/Python/demos/demo_cpu_regularisers.py
index 78e9aff..e6befa9 100644
--- a/Wrappers/Python/demos/demo_cpu_regularisers.py
+++ b/Wrappers/Python/demos/demo_cpu_regularisers.py
@@ -12,7 +12,7 @@ import matplotlib.pyplot as plt
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, TNV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, TNV, NDF, Diff4th
from ccpi.filters.regularisers import PatchSelect, NLTV
from qualitymetrics import rmse
###############################################################################
@@ -225,8 +225,8 @@ pars = {'algorithm' : TGV, \
'input' : u0,\
'regularisation_parameter':0.04, \
'alpha1':1.0,\
- 'alpha0':0.7,\
- 'number_of_iterations' :250 ,\
+ 'alpha0':2.0,\
+ 'number_of_iterations' :1350 ,\
'LipshitzConstant' :12 ,\
}
@@ -358,13 +358,13 @@ print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
## plot
fig = plt.figure()
-plt.suptitle('Performance of DIFF4th regulariser using the CPU')
+plt.suptitle('Performance of Diff4th regulariser using the CPU')
a=fig.add_subplot(1,2,1)
a.set_title('Noisy Image')
imgplot = plt.imshow(u0,cmap="gray")
# set parameters
-pars = {'algorithm' : DIFF4th, \
+pars = {'algorithm' : Diff4th, \
'input' : u0,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -372,9 +372,9 @@ pars = {'algorithm' : DIFF4th, \
'time_marching_parameter':0.0015
}
-print ("#############DIFF4th CPU################")
+print ("#############Diff4th CPU################")
start_time = timeit.default_timer()
-diff4_cpu = DIFF4th(pars['input'],
+diff4_cpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
diff --git a/Wrappers/Python/demos/demo_cpu_regularisers3D.py b/Wrappers/Python/demos/demo_cpu_regularisers3D.py
index 9c28de1..2d2fc22 100644
--- a/Wrappers/Python/demos/demo_cpu_regularisers3D.py
+++ b/Wrappers/Python/demos/demo_cpu_regularisers3D.py
@@ -12,7 +12,7 @@ import matplotlib.pyplot as plt
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, LLT_ROF, FGP_dTV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, Diff4th
from qualitymetrics import rmse
###############################################################################
def printParametersToString(pars):
@@ -68,7 +68,7 @@ Im2[:,0:M] = Im[:,0:M]
Im = Im2
del Im2
"""
-slices = 20
+slices = 15
noisyVol = np.zeros((slices,N,M),dtype='float32')
noisyRef = np.zeros((slices,N,M),dtype='float32')
@@ -96,7 +96,7 @@ pars = {'algorithm': ROF_TV, \
'input' : noisyVol,\
'regularisation_parameter':0.04,\
'number_of_iterations': 500,\
- 'time_marching_parameter': 0.0025
+ 'time_marching_parameter': 0.0025
}
print ("#############ROF TV CPU####################")
start_time = timeit.default_timer()
@@ -264,6 +264,54 @@ plt.title('{}'.format('Recovered volume on the CPU using LLT-ROF'))
#%%
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+print ("_______________TGV (3D)_________________")
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+
+## plot
+fig = plt.figure()
+plt.suptitle('Performance of TGV regulariser using the CPU')
+a=fig.add_subplot(1,2,1)
+a.set_title('Noisy Image')
+imgplot = plt.imshow(noisyVol[10,:,:],cmap="gray")
+
+# set parameters
+pars = {'algorithm' : TGV, \
+ 'input' : noisyVol,\
+ 'regularisation_parameter':0.04, \
+ 'alpha1':1.0,\
+ 'alpha0':2.0,\
+ 'number_of_iterations' :250 ,\
+ 'LipshitzConstant' :12 ,\
+ }
+
+print ("#############TGV CPU####################")
+start_time = timeit.default_timer()
+tgv_cpu3D = TGV(pars['input'],
+ pars['regularisation_parameter'],
+ pars['alpha1'],
+ pars['alpha0'],
+ pars['number_of_iterations'],
+ pars['LipshitzConstant'],'cpu')
+
+
+rms = rmse(idealVol, tgv_cpu3D)
+pars['rmse'] = rms
+
+txtstr = printParametersToString(pars)
+txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
+print (txtstr)
+a=fig.add_subplot(1,2,2)
+
+# these are matplotlib.patch.Patch properties
+props = dict(boxstyle='round', facecolor='wheat', alpha=0.75)
+# place a text box in upper left in axes coords
+a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
+ verticalalignment='top', bbox=props)
+imgplot = plt.imshow(tgv_cpu3D[10,:,:], cmap="gray")
+plt.title('{}'.format('Recovered volume on the CPU using TGV'))
+
+#%%
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("________________NDF (3D)___________________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
@@ -322,7 +370,7 @@ a.set_title('Noisy volume')
imgplot = plt.imshow(noisyVol[10,:,:],cmap="gray")
# set parameters
-pars = {'algorithm' : DIFF4th, \
+pars = {'algorithm' : Diff4th, \
'input' : noisyVol,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -330,9 +378,9 @@ pars = {'algorithm' : DIFF4th, \
'time_marching_parameter':0.0015
}
-print ("#############DIFF4th CPU################")
+print ("#############Diff4th CPU################")
start_time = timeit.default_timer()
-diff4th_cpu3D = DIFF4th(pars['input'],
+diff4th_cpu3D = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
diff --git a/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py b/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py
index 6529b5c..230a761 100644
--- a/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py
+++ b/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py
@@ -12,7 +12,7 @@ import matplotlib.pyplot as plt
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, Diff4th
from ccpi.filters.regularisers import PatchSelect
from qualitymetrics import rmse
###############################################################################
@@ -323,8 +323,8 @@ pars = {'algorithm' : TGV, \
'input' : u0,\
'regularisation_parameter':0.04, \
'alpha1':1.0,\
- 'alpha0':0.7,\
- 'number_of_iterations' :250 ,\
+ 'alpha0':2.0,\
+ 'number_of_iterations' :400 ,\
'LipshitzConstant' :12 ,\
}
@@ -570,7 +570,7 @@ a.set_title('Noisy Image')
imgplot = plt.imshow(u0,cmap="gray")
# set parameters
-pars = {'algorithm' : DIFF4th, \
+pars = {'algorithm' : Diff4th, \
'input' : u0,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -580,7 +580,7 @@ pars = {'algorithm' : DIFF4th, \
print ("#############Diff4th CPU####################")
start_time = timeit.default_timer()
-diff4th_cpu = DIFF4th(pars['input'],
+diff4th_cpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
@@ -604,7 +604,7 @@ plt.title('{}'.format('CPU results'))
print ("##############Diff4th GPU##################")
start_time = timeit.default_timer()
-diff4th_gpu = DIFF4th(pars['input'],
+diff4th_gpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
@@ -612,7 +612,7 @@ diff4th_gpu = DIFF4th(pars['input'],
rms = rmse(Im, diff4th_gpu)
pars['rmse'] = rms
-pars['algorithm'] = DIFF4th
+pars['algorithm'] = Diff4th
txtstr = printParametersToString(pars)
txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
print (txtstr)
diff --git a/Wrappers/Python/demos/demo_gpu_regularisers.py b/Wrappers/Python/demos/demo_gpu_regularisers.py
index 2ada559..e1c6575 100644
--- a/Wrappers/Python/demos/demo_gpu_regularisers.py
+++ b/Wrappers/Python/demos/demo_gpu_regularisers.py
@@ -12,7 +12,7 @@ import matplotlib.pyplot as plt
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, Diff4th
from ccpi.filters.regularisers import PatchSelect, NLTV
from qualitymetrics import rmse
###############################################################################
@@ -223,8 +223,8 @@ pars = {'algorithm' : TGV, \
'input' : u0,\
'regularisation_parameter':0.04, \
'alpha1':1.0,\
- 'alpha0':0.7,\
- 'number_of_iterations' :250 ,\
+ 'alpha0':2.0,\
+ 'number_of_iterations' :1250 ,\
'LipshitzConstant' :12 ,\
}
@@ -355,13 +355,13 @@ print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
## plot
fig = plt.figure()
-plt.suptitle('Performance of DIFF4th regulariser using the GPU')
+plt.suptitle('Performance of Diff4th regulariser using the GPU')
a=fig.add_subplot(1,2,1)
a.set_title('Noisy Image')
imgplot = plt.imshow(u0,cmap="gray")
# set parameters
-pars = {'algorithm' : DIFF4th, \
+pars = {'algorithm' : Diff4th, \
'input' : u0,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -371,7 +371,7 @@ pars = {'algorithm' : DIFF4th, \
print ("#############DIFF4th CPU################")
start_time = timeit.default_timer()
-diff4_gpu = DIFF4th(pars['input'],
+diff4_gpu = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
diff --git a/Wrappers/Python/demos/demo_gpu_regularisers3D.py b/Wrappers/Python/demos/demo_gpu_regularisers3D.py
index d5f9a39..b6058d2 100644
--- a/Wrappers/Python/demos/demo_gpu_regularisers3D.py
+++ b/Wrappers/Python/demos/demo_gpu_regularisers3D.py
@@ -12,7 +12,7 @@ import matplotlib.pyplot as plt
import numpy as np
import os
import timeit
-from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, LLT_ROF, FGP_dTV, NDF, DIFF4th
+from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, TGV, LLT_ROF, FGP_dTV, NDF, Diff4th
from qualitymetrics import rmse
###############################################################################
def printParametersToString(pars):
@@ -67,7 +67,7 @@ Im = Im2
del Im2
"""
-#%%
+
slices = 20
filename = os.path.join(".." , ".." , ".." , "data" ,"lena_gray_512.tif")
@@ -268,6 +268,53 @@ plt.title('{}'.format('Recovered volume on the GPU using LLT-ROF'))
#%%
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+print ("_______________TGV (3D)_________________")
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+
+## plot
+fig = plt.figure()
+plt.suptitle('Performance of TGV regulariser using the GPU')
+a=fig.add_subplot(1,2,1)
+a.set_title('Noisy Image')
+imgplot = plt.imshow(noisyVol[10,:,:],cmap="gray")
+
+# set parameters
+pars = {'algorithm' : TGV, \
+ 'input' : noisyVol,\
+ 'regularisation_parameter':0.04, \
+ 'alpha1':1.0,\
+ 'alpha0':2.0,\
+ 'number_of_iterations' :600 ,\
+ 'LipshitzConstant' :12 ,\
+ }
+
+print ("#############TGV GPU####################")
+start_time = timeit.default_timer()
+tgv_gpu3D = TGV(pars['input'],
+ pars['regularisation_parameter'],
+ pars['alpha1'],
+ pars['alpha0'],
+ pars['number_of_iterations'],
+ pars['LipshitzConstant'],'gpu')
+
+
+rms = rmse(idealVol, tgv_gpu3D)
+pars['rmse'] = rms
+
+txtstr = printParametersToString(pars)
+txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
+print (txtstr)
+a=fig.add_subplot(1,2,2)
+
+# these are matplotlib.patch.Patch properties
+props = dict(boxstyle='round', facecolor='wheat', alpha=0.75)
+# place a text box in upper left in axes coords
+a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
+ verticalalignment='top', bbox=props)
+imgplot = plt.imshow(tgv_gpu3D[10,:,:], cmap="gray")
+plt.title('{}'.format('Recovered volume on the GPU using TGV'))
+#%%
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("_______________NDF-TV (3D)_________________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
@@ -326,7 +373,7 @@ a.set_title('Noisy Image')
imgplot = plt.imshow(noisyVol[10,:,:],cmap="gray")
# set parameters
-pars = {'algorithm' : DIFF4th, \
+pars = {'algorithm' : Diff4th, \
'input' : noisyVol,\
'regularisation_parameter':3.5, \
'edge_parameter':0.02,\
@@ -336,7 +383,7 @@ pars = {'algorithm' : DIFF4th, \
print ("#############DIFF4th CPU################")
start_time = timeit.default_timer()
-diff4_gpu3D = DIFF4th(pars['input'],
+diff4_gpu3D = Diff4th(pars['input'],
pars['regularisation_parameter'],
pars['edge_parameter'],
pars['number_of_iterations'],
diff --git a/Wrappers/Python/src/cpu_regularisers.pyx b/Wrappers/Python/src/cpu_regularisers.pyx
index 7d57ed1..11a0617 100644
--- a/Wrappers/Python/src/cpu_regularisers.pyx
+++ b/Wrappers/Python/src/cpu_regularisers.pyx
@@ -22,7 +22,7 @@ cdef extern float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar,
cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ);
cdef extern float SB_TV_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int printM, int dimX, int dimY, int dimZ);
cdef extern float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
-cdef extern float TGV_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY);
+cdef extern float TGV_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY, int dimZ);
cdef extern float Diffusion_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int penaltytype, int dimX, int dimY, int dimZ);
cdef extern float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
cdef extern float TNV_CPU_main(float *Input, float *u, float lambdaPar, int maxIter, float tol, int dimX, int dimY, int dimZ);
@@ -202,12 +202,8 @@ def TGV_CPU(inputData, regularisation_parameter, alpha1, alpha0, iterations, Lip
return TGV_2D(inputData, regularisation_parameter, alpha1, alpha0,
iterations, LipshitzConst)
elif inputData.ndim == 3:
- shape = inputData.shape
- out = inputData.copy()
- for i in range(shape[0]):
- out[i,:,:] = TGV_2D(inputData[i,:,:], regularisation_parameter,
- alpha1, alpha0, iterations, LipshitzConst)
- return out
+ return TGV_3D(inputData, regularisation_parameter, alpha1, alpha0,
+ iterations, LipshitzConst)
def TGV_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
float regularisation_parameter,
@@ -229,7 +225,30 @@ def TGV_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
alpha0,
iterationsNumb,
LipshitzConst,
- dims[1],dims[0])
+ dims[1],dims[0],1)
+ return outputData
+def TGV_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run TGV iterations for 3D data */
+ TGV_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[2], dims[1], dims[0])
return outputData
#***************************************************************#
diff --git a/Wrappers/Python/src/gpu_regularisers.pyx b/Wrappers/Python/src/gpu_regularisers.pyx
index 47a6149..b52f669 100644
--- a/Wrappers/Python/src/gpu_regularisers.pyx
+++ b/Wrappers/Python/src/gpu_regularisers.pyx
@@ -23,7 +23,7 @@ CUDAErrorMessage = 'CUDA error'
cdef extern int TV_ROF_GPU_main(float* Input, float* Output, float lambdaPar, int iter, float tau, int N, int M, int Z);
cdef extern int TV_FGP_GPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int N, int M, int Z);
cdef extern int TV_SB_GPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int printM, int N, int M, int Z);
-cdef extern int TGV_GPU_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY);
+cdef extern int TGV_GPU_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY, int dimZ);
cdef extern int LLT_ROF_GPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int N, int M, int Z);
cdef extern int NonlDiff_GPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int penaltytype, int N, int M, int Z);
cdef extern int dTV_FGP_GPU_main(float *Input, float *InputRef, float *Output, float lambdaPar, int iterationsNumb, float epsil, float eta, int methodTV, int nonneg, int printM, int N, int M, int Z);
@@ -102,12 +102,7 @@ def TGV_GPU(inputData, regularisation_parameter, alpha1, alpha0, iterations, Lip
if inputData.ndim == 2:
return TGV2D(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst)
elif inputData.ndim == 3:
- shape = inputData.shape
- out = inputData.copy()
- for i in range(shape[0]):
- out[i,:,:] = TGV2D(inputData[i,:,:], regularisation_parameter,
- alpha1, alpha0, iterations, LipshitzConst)
- return out
+ return TGV3D(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst)
# Directional Total-variation Fast-Gradient-Projection (FGP)
def dTV_FGP_GPU(inputData,
refdata,
@@ -393,7 +388,6 @@ def LLT_ROF_GPU3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
raise ValueError(CUDAErrorMessage);
-
#***************************************************************#
#***************** Total Generalised Variation *****************#
#***************************************************************#
@@ -417,11 +411,38 @@ def TGV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
alpha0,
iterationsNumb,
LipshitzConst,
- dims[1],dims[0])==0):
+ dims[1],dims[0], 1)==0):
return outputData
else:
raise ValueError(CUDAErrorMessage);
+def TGV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (TGV_GPU_main(
+ &inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
#****************************************************************#
#**************Directional Total-variation FGP ******************#