From 80c5a5e5de2aca8d5c7b96f0adc91b5738cc9025 Mon Sep 17 00:00:00 2001 From: Daniil Kazantsev Date: Mon, 16 Apr 2018 13:38:40 +0100 Subject: SB TV method CPU/GPU added --- Wrappers/Python/ccpi/filters/regularisers.py | 19 ++++ Wrappers/Python/demos/demo_cpu_regularisers.py | 103 ++++++++++++++++++++- .../Python/demos/demo_cpu_vs_gpu_regularisers.py | 90 +++++++++++++++++- Wrappers/Python/demos/demo_gpu_regularisers.py | 102 +++++++++++++++++++- Wrappers/Python/setup-regularisers.py.in | 1 + Wrappers/Python/src/cpu_regularisers.pyx | 58 ++++++++++++ Wrappers/Python/src/gpu_regularisers.pyx | 75 +++++++++++++++ 7 files changed, 436 insertions(+), 12 deletions(-) (limited to 'Wrappers/Python') diff --git a/Wrappers/Python/ccpi/filters/regularisers.py b/Wrappers/Python/ccpi/filters/regularisers.py index 376cc9c..53623c0 100644 --- a/Wrappers/Python/ccpi/filters/regularisers.py +++ b/Wrappers/Python/ccpi/filters/regularisers.py @@ -42,6 +42,25 @@ def FGP_TV(inputData, regularisation_parameter,iterations, else: raise ValueError('Unknown device {0}. Expecting gpu or cpu'\ .format(device)) +def SB_TV(inputData, regularisation_parameter, iterations, + tolerance_param, methodTV, printM, device='cpu'): + if device == 'cpu': + return TV_SB_CPU(inputData, + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM) + elif device == 'gpu': + return TV_SB_GPU(inputData, + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM) + else: + raise ValueError('Unknown device {0}. Expecting gpu or cpu'\ + .format(device)) def FGP_dTV(inputData, refdata, regularisation_parameter, iterations, tolerance_param, eta_const, methodTV, nonneg, printM, device='cpu'): if device == 'cpu': diff --git a/Wrappers/Python/demos/demo_cpu_regularisers.py b/Wrappers/Python/demos/demo_cpu_regularisers.py index 00beb0b..0e4355b 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, FGP_dTV +from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, FGP_dTV from qualitymetrics import rmse ############################################################################### def printParametersToString(pars): @@ -141,13 +141,60 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14, imgplot = plt.imshow(fgp_cpu, cmap="gray") plt.title('{}'.format('CPU results')) +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") +print ("_______________SB-TV (2D)__________________") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") + +## plot +fig = plt.figure(3) +plt.suptitle('Performance of SB-TV 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' : SB_TV, \ + 'input' : u0,\ + 'regularisation_parameter':0.04, \ + 'number_of_iterations' :150 ,\ + 'tolerance_constant':1e-06,\ + 'methodTV': 0 ,\ + 'printingOut': 0 + } + +print ("#############SB TV CPU####################") +start_time = timeit.default_timer() +sb_cpu = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'cpu') + + +rms = rmse(Im, sb_cpu) +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(sb_cpu, cmap="gray") +plt.title('{}'.format('CPU results')) + print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") print ("_____________FGP-dTV (2D)__________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(3) +fig = plt.figure(4) plt.suptitle('Performance of FGP-dTV regulariser using the CPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') @@ -223,7 +270,7 @@ print ("_______________ROF-TV (3D)_________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(4) +fig = plt.figure(5) plt.suptitle('Performance of ROF-TV regulariser using the CPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy 15th slice of a volume') @@ -263,7 +310,7 @@ print ("_______________FGP-TV (3D)__________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(5) +fig = plt.figure(6) plt.suptitle('Performance of FGP-TV regulariser using the CPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') @@ -307,13 +354,59 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14, imgplot = plt.imshow(fgp_cpu3D[10,:,:], cmap="gray") plt.title('{}'.format('Recovered volume on the CPU using FGP-TV')) +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") +print ("_______________SB-TV (3D)_________________") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") + +## plot +fig = plt.figure(7) +plt.suptitle('Performance of SB-TV 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' : SB_TV, \ + 'input' : noisyVol,\ + 'regularisation_parameter':0.04, \ + 'number_of_iterations' :150 ,\ + 'tolerance_constant':0.00001,\ + 'methodTV': 0 ,\ + 'printingOut': 0 + } + +print ("#############SB TV CPU####################") +start_time = timeit.default_timer() +sb_cpu3D = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'cpu') + +rms = rmse(idealVol, sb_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(sb_cpu3D[10,:,:], cmap="gray") +plt.title('{}'.format('Recovered volume on the CPU using SB-TV')) + print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") print ("_______________FGP-dTV (3D)__________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(6) +fig = plt.figure(8) plt.suptitle('Performance of FGP-dTV regulariser using the CPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') diff --git a/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py b/Wrappers/Python/demos/demo_cpu_vs_gpu_regularisers.py index 310cf75..d8e2da7 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, FGP_dTV +from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, FGP_dTV from qualitymetrics import rmse ############################################################################### def printParametersToString(pars): @@ -218,13 +218,99 @@ if (diff_im.sum() > 1): else: print ("Arrays match") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") +print ("____________SB-TV bench___________________") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") + +## plot +fig = plt.figure(3) +plt.suptitle('Comparison of SB-TV regulariser using CPU and GPU implementations') +a=fig.add_subplot(1,4,1) +a.set_title('Noisy Image') +imgplot = plt.imshow(u0,cmap="gray") + +# set parameters +pars = {'algorithm' : SB_TV, \ + 'input' : u0,\ + 'regularisation_parameter':0.04, \ + 'number_of_iterations' :150 ,\ + 'tolerance_constant':1e-05,\ + 'methodTV': 0 ,\ + 'printingOut': 0 + } + +print ("#############SB-TV CPU####################") +start_time = timeit.default_timer() +sb_cpu = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'cpu') + + +rms = rmse(Im, sb_cpu) +pars['rmse'] = rms + +txtstr = printParametersToString(pars) +txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) +print (txtstr) +a=fig.add_subplot(1,4,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(sb_cpu, cmap="gray") +plt.title('{}'.format('CPU results')) + + +print ("##############SB TV GPU##################") +start_time = timeit.default_timer() +sb_gpu = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'gpu') + +rms = rmse(Im, sb_gpu) +pars['rmse'] = rms +pars['algorithm'] = SB_TV +txtstr = printParametersToString(pars) +txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) +print (txtstr) +a=fig.add_subplot(1,4,3) + +# 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(sb_gpu, cmap="gray") +plt.title('{}'.format('GPU results')) + +print ("--------Compare the results--------") +tolerance = 1e-05 +diff_im = np.zeros(np.shape(rof_cpu)) +diff_im = abs(sb_cpu - sb_gpu) +diff_im[diff_im > tolerance] = 1 +a=fig.add_subplot(1,4,4) +imgplot = plt.imshow(diff_im, vmin=0, vmax=1, cmap="gray") +plt.title('{}'.format('Pixels larger threshold difference')) +if (diff_im.sum() > 1): + print ("Arrays do not match!") +else: + print ("Arrays match") + print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") print ("____________FGP-dTV bench___________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(3) +fig = plt.figure(4) plt.suptitle('Comparison of FGP-dTV regulariser using CPU and GPU implementations') a=fig.add_subplot(1,4,1) a.set_title('Noisy Image') diff --git a/Wrappers/Python/demos/demo_gpu_regularisers.py b/Wrappers/Python/demos/demo_gpu_regularisers.py index 24a3c88..25d8d85 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, FGP_dTV +from ccpi.filters.regularisers import ROF_TV, FGP_TV, SB_TV, FGP_dTV from qualitymetrics import rmse ############################################################################### def printParametersToString(pars): @@ -139,12 +139,59 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14, imgplot = plt.imshow(fgp_gpu, cmap="gray") plt.title('{}'.format('GPU results')) + print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") -print ("____________FGP-dTV bench___________________") +print ("____________SB-TV bench___________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot fig = plt.figure(3) +plt.suptitle('Performance of the SB-TV 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' : SB_TV, \ + 'input' : u0,\ + 'regularisation_parameter':0.04, \ + 'number_of_iterations' :150 ,\ + 'tolerance_constant':1e-06,\ + 'methodTV': 0 ,\ + 'printingOut': 0 + } + +print ("##############SB TV GPU##################") +start_time = timeit.default_timer() +sb_gpu = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'gpu') + +rms = rmse(Im, sb_gpu) +pars['rmse'] = rms +pars['algorithm'] = SB_TV +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(sb_gpu, cmap="gray") +plt.title('{}'.format('GPU results')) + +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") +print ("____________FGP-dTV bench___________________") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") + +## plot +fig = plt.figure(4) plt.suptitle('Performance of the FGP-dTV regulariser using the GPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') @@ -219,7 +266,7 @@ print ("_______________ROF-TV (3D)_________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(4) +fig = plt.figure(5) plt.suptitle('Performance of ROF-TV regulariser using the GPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy 15th slice of a volume') @@ -259,7 +306,7 @@ print ("_______________FGP-TV (3D)__________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(5) +fig = plt.figure(6) plt.suptitle('Performance of FGP-TV regulariser using the GPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') @@ -302,13 +349,58 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14, imgplot = plt.imshow(fgp_gpu3D[10,:,:], cmap="gray") plt.title('{}'.format('Recovered volume on the GPU using FGP-TV')) +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") +print ("_______________SB-TV (3D)__________________") +print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") + +## plot +fig = plt.figure(7) +plt.suptitle('Performance of SB-TV 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' : SB_TV, \ + 'input' : noisyVol,\ + 'regularisation_parameter':0.04, \ + 'number_of_iterations' :100 ,\ + 'tolerance_constant':1e-05,\ + 'methodTV': 0 ,\ + 'printingOut': 0 + } + +print ("#############SB TV GPU####################") +start_time = timeit.default_timer() +sb_gpu3D = SB_TV(pars['input'], + pars['regularisation_parameter'], + pars['number_of_iterations'], + pars['tolerance_constant'], + pars['methodTV'], + pars['printingOut'],'gpu') + +rms = rmse(idealVol, sb_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(sb_gpu3D[10,:,:], cmap="gray") +plt.title('{}'.format('Recovered volume on the GPU using SB-TV')) print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") print ("_______________FGP-dTV (3D)________________") print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot -fig = plt.figure(6) +fig = plt.figure(8) plt.suptitle('Performance of FGP-dTV regulariser using the GPU') a=fig.add_subplot(1,2,1) a.set_title('Noisy Image') diff --git a/Wrappers/Python/setup-regularisers.py.in b/Wrappers/Python/setup-regularisers.py.in index c7ebb5c..0681cc4 100644 --- a/Wrappers/Python/setup-regularisers.py.in +++ b/Wrappers/Python/setup-regularisers.py.in @@ -36,6 +36,7 @@ extra_include_dirs += [os.path.join(".." , ".." , "Core"), os.path.join(".." , ".." , "Core", "regularisers_CPU"), os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_FGP" ) , os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_ROF" ) , + os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_SB" ) , os.path.join(".." , ".." , "Core", "regularisers_GPU" , "dTV_FGP" ) , "."] diff --git a/Wrappers/Python/src/cpu_regularisers.pyx b/Wrappers/Python/src/cpu_regularisers.pyx index 1661375..b8d2523 100644 --- a/Wrappers/Python/src/cpu_regularisers.pyx +++ b/Wrappers/Python/src/cpu_regularisers.pyx @@ -20,6 +20,7 @@ cimport numpy as np cdef extern float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ); 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 TV_SB_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 dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambdaPar, int iterationsNumb, float epsil, float eta, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); @@ -125,6 +126,63 @@ def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, printM, dims[2], dims[1], dims[0]) return outputData + +#***************************************************************# +#********************** Total-variation SB *********************# +#***************************************************************# +#*************** Total-variation Split Bregman (SB)*************# +def TV_SB_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM): + if inputData.ndim == 2: + return TV_SB_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM) + elif inputData.ndim == 3: + return TV_SB_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM) + +def TV_SB_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, + float regularisation_parameter, + int iterationsNumb, + float tolerance_param, + int methodTV, + int printM): + + cdef long dims[2] + dims[0] = inputData.shape[0] + dims[1] = inputData.shape[1] + + cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \ + np.zeros([dims[0],dims[1]], dtype='float32') + + #/* Run SB-TV iterations for 2D data */ + TV_SB_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, + iterationsNumb, + tolerance_param, + methodTV, + printM, + dims[0], dims[1], 1) + + return outputData + +def TV_SB_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, + float regularisation_parameter, + int iterationsNumb, + float tolerance_param, + int methodTV, + int printM): + 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 SB-TV iterations for 3D data */ + TV_SB_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, + iterationsNumb, + tolerance_param, + methodTV, + printM, + dims[2], dims[1], dims[0]) + return outputData #****************************************************************# #**************Directional Total-variation FGP ******************# #****************************************************************# diff --git a/Wrappers/Python/src/gpu_regularisers.pyx b/Wrappers/Python/src/gpu_regularisers.pyx index 18efdcd..36eec95 100644 --- a/Wrappers/Python/src/gpu_regularisers.pyx +++ b/Wrappers/Python/src/gpu_regularisers.pyx @@ -20,6 +20,7 @@ cimport numpy as np cdef extern void TV_ROF_GPU_main(float* Input, float* Output, float lambdaPar, int iter, float tau, int N, int M, int Z); cdef extern void 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 void 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 void 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); # Total-variation Rudin-Osher-Fatemi (ROF) @@ -62,6 +63,27 @@ def TV_FGP_GPU(inputData, methodTV, nonneg, printM) +# Total-variation Split Bregman (SB) +def TV_SB_GPU(inputData, + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM): + if inputData.ndim == 2: + return SBTV2D(inputData, + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM) + elif inputData.ndim == 3: + return SBTV3D(inputData, + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM) # Directional Total-variation Fast-Gradient-Projection (FGP) def dTV_FGP_GPU(inputData, refdata, @@ -197,7 +219,60 @@ def FGPTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, dims[2], dims[1], dims[0]); return outputData +#***************************************************************# +#********************** Total-variation SB *********************# +#***************************************************************# +#*************** Total-variation Split Bregman (SB)*************# +def SBTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, + float regularisation_parameter, + int iterations, + float tolerance_param, + int methodTV, + int printM): + + cdef long dims[2] + dims[0] = inputData.shape[0] + dims[1] = inputData.shape[1] + + cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \ + np.zeros([dims[0],dims[1]], dtype='float32') + + # Running CUDA code here + TV_SB_GPU_main(&inputData[0,0], &outputData[0,0], + regularisation_parameter, + iterations, + tolerance_param, + methodTV, + printM, + dims[0], dims[1], 1); + + return outputData +def SBTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, + float regularisation_parameter, + int iterations, + float tolerance_param, + int methodTV, + int printM): + + 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 + TV_SB_GPU_main(&inputData[0,0,0], &outputData[0,0,0], + regularisation_parameter , + iterations, + tolerance_param, + methodTV, + printM, + dims[2], dims[1], dims[0]); + + return outputData #****************************************************************# #**************Directional Total-variation FGP ******************# #****************************************************************# -- cgit v1.2.3 From 57727584760e6b1a980071587e1f1e8910c7d6a3 Mon Sep 17 00:00:00 2001 From: algol Date: Mon, 16 Apr 2018 15:16:32 +0100 Subject: SB TV demos all work --- Wrappers/Python/ccpi/filters/regularisers.py | 4 ++-- Wrappers/Python/src/cpu_regularisers.pyx | 6 +++--- 2 files changed, 5 insertions(+), 5 deletions(-) (limited to 'Wrappers/Python') diff --git a/Wrappers/Python/ccpi/filters/regularisers.py b/Wrappers/Python/ccpi/filters/regularisers.py index 53623c0..50c4374 100644 --- a/Wrappers/Python/ccpi/filters/regularisers.py +++ b/Wrappers/Python/ccpi/filters/regularisers.py @@ -2,8 +2,8 @@ script which assigns a proper device core function based on a flag ('cpu' or 'gpu') """ -from ccpi.filters.cpu_regularisers_cython import TV_ROF_CPU, TV_FGP_CPU, dTV_FGP_CPU -from ccpi.filters.gpu_regularisers import TV_ROF_GPU, TV_FGP_GPU, dTV_FGP_GPU +from ccpi.filters.cpu_regularisers_cython import TV_ROF_CPU, TV_FGP_CPU, TV_SB_CPU, dTV_FGP_CPU +from ccpi.filters.gpu_regularisers import TV_ROF_GPU, TV_FGP_GPU, TV_SB_GPU, dTV_FGP_GPU def ROF_TV(inputData, regularisation_parameter, iterations, time_marching_parameter,device='cpu'): diff --git a/Wrappers/Python/src/cpu_regularisers.pyx b/Wrappers/Python/src/cpu_regularisers.pyx index b8d2523..417670d 100644 --- a/Wrappers/Python/src/cpu_regularisers.pyx +++ b/Wrappers/Python/src/cpu_regularisers.pyx @@ -20,7 +20,7 @@ cimport numpy as np cdef extern float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ); 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 TV_SB_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 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 dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambdaPar, int iterationsNumb, float epsil, float eta, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); @@ -152,7 +152,7 @@ def TV_SB_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, np.zeros([dims[0],dims[1]], dtype='float32') #/* Run SB-TV iterations for 2D data */ - TV_SB_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, + SB_TV_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, iterationsNumb, tolerance_param, methodTV, @@ -176,7 +176,7 @@ def TV_SB_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, np.zeros([dims[0], dims[1], dims[2]], dtype='float32') #/* Run SB-TV iterations for 3D data */ - TV_SB_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, + SB_TV_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, iterationsNumb, tolerance_param, methodTV, -- cgit v1.2.3