From 970e82c7ee5fedab771480f42293963fdc32d17b Mon Sep 17 00:00:00 2001 From: dkazanc Date: Mon, 9 Dec 2019 14:08:21 +0000 Subject: fixes gpu issues with pdtv --- demos/demo_cpu_regularisers.py | 8 ++-- demos/demo_cpu_regularisers3D.py | 6 +-- demos/demo_cpu_vs_gpu_regularisers.py | 12 ++---- demos/demo_gpu_regularisers.py | 6 +-- demos/demo_gpu_regularisers3D.py | 6 +-- src/Core/regularisers_CPU/PD_TV_core.c | 6 +-- src/Core/regularisers_CPU/PD_TV_core.h | 2 +- src/Core/regularisers_GPU/TV_PD_GPU_core.cu | 13 +++--- src/Core/regularisers_GPU/TV_PD_GPU_core.h | 2 +- src/Matlab/mex_compile/regularisers_CPU/PD_TV.c | 50 +++++++++++----------- .../mex_compile/regularisers_GPU/PD_TV_GPU.cpp | 49 ++++++++++----------- src/Python/ccpi/filters/regularisers.py | 8 ++-- src/Python/src/cpu_regularisers.pyx | 18 +++----- src/Python/src/gpu_regularisers.pyx | 16 +++---- test/test_CPU_regularisers.py | 2 +- test/test_run_test.py | 9 ++-- 16 files changed, 90 insertions(+), 123 deletions(-) diff --git a/demos/demo_cpu_regularisers.py b/demos/demo_cpu_regularisers.py index 7d66b7f..50a5065 100644 --- a/demos/demo_cpu_regularisers.py +++ b/demos/demo_cpu_regularisers.py @@ -176,11 +176,10 @@ pars = {'algorithm' : PD_TV, \ 'input' : u0,\ 'regularisation_parameter':0.02, \ 'number_of_iterations' :1500 ,\ - 'tolerance_constant':1e-08,\ + 'tolerance_constant':1e-06,\ 'methodTV': 0 ,\ 'nonneg': 1, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV CPU####################") start_time = timeit.default_timer() @@ -190,8 +189,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'cpu') + pars['lipschitz_const'],'cpu') Qtools = QualityTools(Im, pd_cpu) pars['rmse'] = Qtools.rmse() diff --git a/demos/demo_cpu_regularisers3D.py b/demos/demo_cpu_regularisers3D.py index cfdd2d4..0e7e9be 100644 --- a/demos/demo_cpu_regularisers3D.py +++ b/demos/demo_cpu_regularisers3D.py @@ -188,8 +188,7 @@ pars = {'algorithm' : PD_TV, \ 'tolerance_constant':1e-06,\ 'methodTV': 0 ,\ 'nonneg': 0, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############FGP TV GPU####################") start_time = timeit.default_timer() @@ -199,8 +198,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'cpu') + pars['lipschitz_const'],'cpu') Qtools = QualityTools(idealVol, pd_cpu3D) pars['rmse'] = Qtools.rmse() diff --git a/demos/demo_cpu_vs_gpu_regularisers.py b/demos/demo_cpu_vs_gpu_regularisers.py index 015dfc6..a34bc19 100644 --- a/demos/demo_cpu_vs_gpu_regularisers.py +++ b/demos/demo_cpu_vs_gpu_regularisers.py @@ -241,8 +241,7 @@ pars = {'algorithm' : PD_TV, \ 'tolerance_constant':0.0,\ 'methodTV': 0 ,\ 'nonneg': 0, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV CPU####################") start_time = timeit.default_timer() @@ -252,8 +251,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'cpu') + pars['lipschitz_const'],'cpu') Qtools = QualityTools(Im, pd_cpu) pars['rmse'] = Qtools.rmse() @@ -279,8 +277,7 @@ pars = {'algorithm' : PD_TV, \ 'tolerance_constant':0.0,\ 'methodTV': 0 ,\ 'nonneg': 0, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV CPU####################") start_time = timeit.default_timer() @@ -290,8 +287,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'gpu') + pars['lipschitz_const'],'gpu') Qtools = QualityTools(Im, pd_gpu) pars['rmse'] = Qtools.rmse() diff --git a/demos/demo_gpu_regularisers.py b/demos/demo_gpu_regularisers.py index 5131847..c6114db 100644 --- a/demos/demo_gpu_regularisers.py +++ b/demos/demo_gpu_regularisers.py @@ -176,8 +176,7 @@ pars = {'algorithm' : PD_TV, \ 'tolerance_constant':1e-06,\ 'methodTV': 0 ,\ 'nonneg': 1, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV CPU####################") start_time = timeit.default_timer() @@ -187,8 +186,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'gpu') + pars['lipschitz_const'],'gpu') Qtools = QualityTools(Im, pd_gpu) pars['rmse'] = Qtools.rmse() diff --git a/demos/demo_gpu_regularisers3D.py b/demos/demo_gpu_regularisers3D.py index 2c25d01..18d97e5 100644 --- a/demos/demo_gpu_regularisers3D.py +++ b/demos/demo_gpu_regularisers3D.py @@ -192,8 +192,7 @@ pars = {'algorithm' : PD_TV, \ 'tolerance_constant':1e-06,\ 'methodTV': 0 ,\ 'nonneg': 0, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV GPU####################") start_time = timeit.default_timer() @@ -203,8 +202,7 @@ start_time = timeit.default_timer() pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'gpu') + pars['lipschitz_const'],'gpu') Qtools = QualityTools(idealVol, pd_gpu3D) pars['rmse'] = Qtools.rmse() diff --git a/src/Core/regularisers_CPU/PD_TV_core.c b/src/Core/regularisers_CPU/PD_TV_core.c index 534091b..c1b21e7 100644 --- a/src/Core/regularisers_CPU/PD_TV_core.c +++ b/src/Core/regularisers_CPU/PD_TV_core.c @@ -29,7 +29,6 @@ * 5. lipschitz_const: convergence related parameter * 6. TV-type: methodTV - 'iso' (0) or 'l1' (1) * 7. nonneg: 'nonnegativity (0 is OFF by default, 1 is ON) - * 8. tau: time marching parameter * Output: * [1] TV - Filtered/regularized image/volume @@ -38,16 +37,17 @@ * [1] Antonin Chambolle, Thomas Pock. "A First-Order Primal-Dual Algorithm for Convex Problems with Applications to Imaging", 2010 */ -float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ) +float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ) { int ll; long j, DimTotal; - float re, re1, sigma, theta, lt; + float re, re1, sigma, theta, lt, tau; re = 0.0f; re1 = 0.0f; int count = 0; //tau = 1.0/powf(lipschitz_const,0.5); //sigma = 1.0/powf(lipschitz_const,0.5); + tau = lambdaPar*0.1f; sigma = 1.0/(lipschitz_const*tau); theta = 1.0f; lt = tau/lambdaPar; diff --git a/src/Core/regularisers_CPU/PD_TV_core.h b/src/Core/regularisers_CPU/PD_TV_core.h index 97edc05..294e75c 100644 --- a/src/Core/regularisers_CPU/PD_TV_core.h +++ b/src/Core/regularisers_CPU/PD_TV_core.h @@ -47,7 +47,7 @@ limitations under the License. #ifdef __cplusplus extern "C" { #endif -float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ); +float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ); CCPI_EXPORT float DualP2D(float *U, float *P1, float *P2, long dimX, long dimY, float sigma); CCPI_EXPORT float DivProj2D(float *U, float *Input, float *P1, float *P2, long dimX, long dimY, float lt, float tau); diff --git a/src/Core/regularisers_GPU/TV_PD_GPU_core.cu b/src/Core/regularisers_GPU/TV_PD_GPU_core.cu index e57020a..01e0ab5 100644 --- a/src/Core/regularisers_GPU/TV_PD_GPU_core.cu +++ b/src/Core/regularisers_GPU/TV_PD_GPU_core.cu @@ -33,7 +33,6 @@ limitations under the License. * 5. lipschitz_const: convergence related parameter * 6. TV-type: methodTV - 'iso' (0) or 'l1' (1) * 7. nonneg: 'nonnegativity (0 is OFF by default, 1 is ON) - * 8. tau: time marching parameter * Output: * [1] TV - Filtered/regularized image/volume @@ -42,8 +41,8 @@ limitations under the License. * [1] Antonin Chambolle, Thomas Pock. "A First-Order Primal-Dual Algorithm for Convex Problems with Applications to Imaging", 2010 */ -#define BLKXSIZE2D 8 -#define BLKYSIZE2D 8 +#define BLKXSIZE2D 16 +#define BLKYSIZE2D 16 #define BLKXSIZE 8 #define BLKYSIZE 8 @@ -322,12 +321,9 @@ __global__ void PDResidCalc3D_kernel(float *Input1, float *Input2, float* Output Output[index] = Input1[index] - Input2[index]; } } - - /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/ - ////////////MAIN HOST FUNCTION /////////////// -extern "C" int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ) +extern "C" int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ) { int deviceCount = -1; // number of devices cudaGetDeviceCount(&deviceCount); @@ -336,9 +332,10 @@ extern "C" int TV_PD_GPU_main(float *Input, float *Output, float *infovector, fl return -1; } int count = 0, i; - float re, sigma, theta, lt; + float re, sigma, theta, lt, tau; re = 0.0f; + tau = lambdaPar*0.1f; sigma = 1.0/(lipschitz_const*tau); theta = 1.0f; lt = tau/lambdaPar; diff --git a/src/Core/regularisers_GPU/TV_PD_GPU_core.h b/src/Core/regularisers_GPU/TV_PD_GPU_core.h index 2b123d9..48e353e 100644 --- a/src/Core/regularisers_GPU/TV_PD_GPU_core.h +++ b/src/Core/regularisers_GPU/TV_PD_GPU_core.h @@ -4,6 +4,6 @@ #include "CCPiDefines.h" #include -extern "C" CCPI_EXPORT int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ); +extern "C" CCPI_EXPORT int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ); #endif diff --git a/src/Matlab/mex_compile/regularisers_CPU/PD_TV.c b/src/Matlab/mex_compile/regularisers_CPU/PD_TV.c index e5ab1e4..f8f5272 100644 --- a/src/Matlab/mex_compile/regularisers_CPU/PD_TV.c +++ b/src/Matlab/mex_compile/regularisers_CPU/PD_TV.c @@ -30,8 +30,7 @@ * 5. TV-type: methodTV - 'iso' (0) or 'l1' (1) * 6. nonneg: 'nonnegativity (0 is OFF by default, 1 is ON) * 7. lipschitz_const: convergence related parameter - * 8. tau: convergence related parameter - + * Output: * [1] TV - Filtered/regularized image/volume * [2] Information vector which contains [iteration no., reached tolerance] @@ -41,19 +40,19 @@ void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) - + { int number_of_dims, iter, methTV, nonneg; mwSize dimX, dimY, dimZ; const mwSize *dim_array; - float *Input, *infovec=NULL, *Output=NULL, lambda, epsil, lipschitz_const, tau; - + float *Input, *infovec=NULL, *Output=NULL, lambda, epsil, lipschitz_const; + number_of_dims = mxGetNumberOfDimensions(prhs[0]); dim_array = mxGetDimensions(prhs[0]); - + /*Handling Matlab input data*/ - if ((nrhs < 2) || (nrhs > 7)) mexErrMsgTxt("At least 2 parameters is required, all parameters are: Image(2D/3D), Regularization parameter, iterations number, tolerance, penalty type ('iso' or 'l1'), nonnegativity switch, lipschitz_const"); - + if ((nrhs < 2) || (nrhs > 6)) mexErrMsgTxt("At least 2 parameters is required, all parameters are: Image(2D/3D), Regularization parameter, iterations number, tolerance, penalty type ('iso' or 'l1'), nonnegativity switch, lipschitz_const"); + Input = (float *) mxGetData(prhs[0]); /*noisy image (2D/3D) */ lambda = (float) mxGetScalar(prhs[1]); /* regularization parameter */ iter = 500; /* default iterations number */ @@ -61,40 +60,39 @@ void mexFunction( methTV = 0; /* default isotropic TV penalty */ nonneg = 0; /* default nonnegativity switch, off - 0 */ lipschitz_const = 8.0; /* lipschitz_const */ - tau = 0.0025; /* tau convergence const */ - + + if (mxGetClassID(prhs[0]) != mxSINGLE_CLASS) {mexErrMsgTxt("The input image must be in a single precision"); } - - if ((nrhs == 3) || (nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) iter = (int) mxGetScalar(prhs[2]); /* iterations number */ - if ((nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) epsil = (float) mxGetScalar(prhs[3]); /* tolerance constant */ - if ((nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) { + + if ((nrhs == 3) || (nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7)) iter = (int) mxGetScalar(prhs[2]); /* iterations number */ + if ((nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7)) epsil = (float) mxGetScalar(prhs[3]); /* tolerance constant */ + if ((nrhs == 5) || (nrhs == 6) || (nrhs == 7)) { char *penalty_type; penalty_type = mxArrayToString(prhs[4]); /* choosing TV penalty: 'iso' or 'l1', 'iso' is the default */ if ((strcmp(penalty_type, "l1") != 0) && (strcmp(penalty_type, "iso") != 0)) mexErrMsgTxt("Choose TV type: 'iso' or 'l1',"); if (strcmp(penalty_type, "l1") == 0) methTV = 1; /* enable 'l1' penalty */ mxFree(penalty_type); } - if ((nrhs == 6) || (nrhs == 7) || (nrhs == 8)) { + if ((nrhs == 6) || (nrhs == 7)) { nonneg = (int) mxGetScalar(prhs[5]); if ((nonneg != 0) && (nonneg != 1)) mexErrMsgTxt("Nonnegativity constraint can be enabled by choosing 1 or off - 0"); } - if ((nrhs == 7) || (nrhs == 8)) lipschitz_const = (float) mxGetScalar(prhs[6]); - if (nrhs == 8) tau = (float) mxGetScalar(prhs[7]); - + if (nrhs == 7) lipschitz_const = (float) mxGetScalar(prhs[6]); + /*Handling Matlab output data*/ - dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2]; - + dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2]; + if (number_of_dims == 2) { dimZ = 1; /*2D case*/ - Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL)); + Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL)); } - if (number_of_dims == 3) { + if (number_of_dims == 3) { Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(3, dim_array, mxSINGLE_CLASS, mxREAL)); - } + } mwSize vecdim[1]; vecdim[0] = 2; infovec = (float*)mxGetPr(plhs[1] = mxCreateNumericArray(1, vecdim, mxSINGLE_CLASS, mxREAL)); - - /* running the function */ - PDTV_CPU_main(Input, Output, infovec, lambda, iter, epsil, lipschitz_const, methTV, nonneg, tau, dimX, dimY, dimZ); + + /* running the function */ + PDTV_CPU_main(Input, Output, infovec, lambda, iter, epsil, lipschitz_const, methTV, nonneg, dimX, dimY, dimZ); } diff --git a/src/Matlab/mex_compile/regularisers_GPU/PD_TV_GPU.cpp b/src/Matlab/mex_compile/regularisers_GPU/PD_TV_GPU.cpp index e853dd3..2c037a5 100644 --- a/src/Matlab/mex_compile/regularisers_GPU/PD_TV_GPU.cpp +++ b/src/Matlab/mex_compile/regularisers_GPU/PD_TV_GPU.cpp @@ -30,8 +30,7 @@ * 5. TV-type: methodTV - 'iso' (0) or 'l1' (1) * 6. nonneg: 'nonnegativity (0 is OFF by default, 1 is ON) * 7. lipschitz_const: convergence related parameter - * 8. tau: convergence related parameter - + * Output: * [1] TV - Filtered/regularized image/volume * [2] Information vector which contains [iteration no., reached tolerance] @@ -42,19 +41,19 @@ void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) - + { int number_of_dims, iter, methTV, nonneg; mwSize dimX, dimY, dimZ; const mwSize *dim_array; - float *Input, *infovec=NULL, *Output=NULL, lambda, epsil, lipschitz_const, tau; - + float *Input, *infovec=NULL, *Output=NULL, lambda, epsil, lipschitz_const; + number_of_dims = mxGetNumberOfDimensions(prhs[0]); dim_array = mxGetDimensions(prhs[0]); - + /*Handling Matlab input data*/ - if ((nrhs < 2) || (nrhs > 7)) mexErrMsgTxt("At least 2 parameters is required, all parameters are: Image(2D/3D), Regularization parameter, iterations number, tolerance, penalty type ('iso' or 'l1'), nonnegativity switch, lipschitz_const"); - + if ((nrhs < 2) || (nrhs > 6)) mexErrMsgTxt("At least 2 parameters is required, all parameters are: Image(2D/3D), Regularization parameter, iterations number, tolerance, penalty type ('iso' or 'l1'), nonnegativity switch, lipschitz_const"); + Input = (float *) mxGetData(prhs[0]); /*noisy image (2D/3D) */ lambda = (float) mxGetScalar(prhs[1]); /* regularization parameter */ iter = 500; /* default iterations number */ @@ -62,40 +61,38 @@ void mexFunction( methTV = 0; /* default isotropic TV penalty */ nonneg = 0; /* default nonnegativity switch, off - 0 */ lipschitz_const = 8.0; /* lipschitz_const */ - tau = 0.0025; /* tau convergence const */ - + if (mxGetClassID(prhs[0]) != mxSINGLE_CLASS) {mexErrMsgTxt("The input image must be in a single precision"); } - - if ((nrhs == 3) || (nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) iter = (int) mxGetScalar(prhs[2]); /* iterations number */ - if ((nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) epsil = (float) mxGetScalar(prhs[3]); /* tolerance constant */ - if ((nrhs == 5) || (nrhs == 6) || (nrhs == 7) || (nrhs == 8)) { + + if ((nrhs == 3) || (nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7)) iter = (int) mxGetScalar(prhs[2]); /* iterations number */ + if ((nrhs == 4) || (nrhs == 5) || (nrhs == 6) || (nrhs == 7)) epsil = (float) mxGetScalar(prhs[3]); /* tolerance constant */ + if ((nrhs == 5) || (nrhs == 6) || (nrhs == 7)) { char *penalty_type; penalty_type = mxArrayToString(prhs[4]); /* choosing TV penalty: 'iso' or 'l1', 'iso' is the default */ if ((strcmp(penalty_type, "l1") != 0) && (strcmp(penalty_type, "iso") != 0)) mexErrMsgTxt("Choose TV type: 'iso' or 'l1',"); if (strcmp(penalty_type, "l1") == 0) methTV = 1; /* enable 'l1' penalty */ mxFree(penalty_type); } - if ((nrhs == 6) || (nrhs == 7) || (nrhs == 8)) { + if ((nrhs == 6) || (nrhs == 7)) { nonneg = (int) mxGetScalar(prhs[5]); if ((nonneg != 0) && (nonneg != 1)) mexErrMsgTxt("Nonnegativity constraint can be enabled by choosing 1 or off - 0"); } - if ((nrhs == 7) || (nrhs == 8)) lipschitz_const = (float) mxGetScalar(prhs[6]); - if (nrhs == 8) tau = (float) mxGetScalar(prhs[7]); - + if (nrhs == 7) lipschitz_const = (float) mxGetScalar(prhs[6]); + /*Handling Matlab output data*/ - dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2]; - + dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2]; + if (number_of_dims == 2) { dimZ = 1; /*2D case*/ - Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL)); + Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array, mxSINGLE_CLASS, mxREAL)); } - if (number_of_dims == 3) { + if (number_of_dims == 3) { Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(3, dim_array, mxSINGLE_CLASS, mxREAL)); - } + } mwSize vecdim[1]; vecdim[0] = 2; infovec = (float*)mxGetPr(plhs[1] = mxCreateNumericArray(1, vecdim, mxSINGLE_CLASS, mxREAL)); - - /* running the function */ - TV_PD_GPU_main(Input, Output, infovec, lambda, iter, epsil, lipschitz_const, methTV, nonneg, tau, dimX, dimY, dimZ); + + /* running the function */ + TV_PD_GPU_main(Input, Output, infovec, lambda, iter, epsil, lipschitz_const, methTV, nonneg, dimX, dimY, dimZ); } diff --git a/src/Python/ccpi/filters/regularisers.py b/src/Python/ccpi/filters/regularisers.py index 5f4001a..e3a984e 100644 --- a/src/Python/ccpi/filters/regularisers.py +++ b/src/Python/ccpi/filters/regularisers.py @@ -53,7 +53,7 @@ def FGP_TV(inputData, regularisation_parameter,iterations, .format(device)) def PD_TV(inputData, regularisation_parameter, iterations, - tolerance_param, methodTV, nonneg, lipschitz_const, tau, device='cpu'): + tolerance_param, methodTV, nonneg, lipschitz_const, device='cpu'): if device == 'cpu': return TV_PD_CPU(inputData, regularisation_parameter, @@ -61,8 +61,7 @@ def PD_TV(inputData, regularisation_parameter, iterations, tolerance_param, methodTV, nonneg, - lipschitz_const, - tau) + lipschitz_const) elif device == 'gpu' and gpu_enabled: return TV_PD_GPU(inputData, regularisation_parameter, @@ -70,8 +69,7 @@ def PD_TV(inputData, regularisation_parameter, iterations, tolerance_param, methodTV, nonneg, - lipschitz_const, - tau) + lipschitz_const) else: if not gpu_enabled and device == 'gpu': raise ValueError ('GPU is not available') diff --git a/src/Python/src/cpu_regularisers.pyx b/src/Python/src/cpu_regularisers.pyx index 8de6aea..08e247c 100644 --- a/src/Python/src/cpu_regularisers.pyx +++ b/src/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 *infovector, float *lambdaPar, int lambda_is_arr, int iterationsNumb, float tau, float epsil, int dimX, int dimY, int dimZ); cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int nonneg, int dimX, int dimY, int dimZ); -cdef extern float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ); +cdef extern float PDTV_CPU_main(float *Input, float *U, float *infovector, float lambdaPar, int iterationsNumb, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ); cdef extern float SB_TV_CPU_main(float *Input, float *Output, float *infovector, float mu, int iter, float epsil, int methodTV, int dimX, int dimY, int dimZ); cdef extern float LLT_ROF_CPU_main(float *Input, float *Output, float *infovector, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, float epsil, int dimX, int dimY, int dimZ); cdef extern float TGV_main(float *Input, float *Output, float *infovector, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, float epsil, int dimX, int dimY, int dimZ); @@ -159,11 +159,11 @@ def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, #****************************************************************# #****************** Total-variation Primal-dual *****************# #****************************************************************# -def TV_PD_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau): +def TV_PD_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const): if inputData.ndim == 2: - return TV_PD_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau) + return TV_PD_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const) elif inputData.ndim == 3: - return TV_PD_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau) + return TV_PD_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const) def TV_PD_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, float regularisation_parameter, @@ -171,8 +171,7 @@ def TV_PD_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, float tolerance_param, int methodTV, int nonneg, - float lipschitz_const, - float tau): + float lipschitz_const): cdef long dims[2] dims[0] = inputData.shape[0] @@ -191,7 +190,6 @@ def TV_PD_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, lipschitz_const, methodTV, nonneg, - tau, dims[1],dims[0], 1) return (outputData,infovec) @@ -200,9 +198,8 @@ def TV_PD_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, int iterationsNumb, float tolerance_param, int methodTV, - int nonneg, - float lipschitz_const, - float tau): + int nonneg, + float lipschitz_const): cdef long dims[3] dims[0] = inputData.shape[0] @@ -221,7 +218,6 @@ def TV_PD_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, lipschitz_const, methodTV, nonneg, - tau, dims[2], dims[1], dims[0]) return (outputData,infovec) diff --git a/src/Python/src/gpu_regularisers.pyx b/src/Python/src/gpu_regularisers.pyx index b22d15e..8a4568e 100644 --- a/src/Python/src/gpu_regularisers.pyx +++ b/src/Python/src/gpu_regularisers.pyx @@ -22,7 +22,7 @@ CUDAErrorMessage = 'CUDA error' cdef extern int TV_ROF_GPU_main(float* Input, float* Output, float *infovector, float *lambdaPar, int lambda_is_arr, int iter, float tau, float epsil, int N, int M, int Z); cdef extern int TV_FGP_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int N, int M, int Z); -cdef extern int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, float tau, int dimX, int dimY, int dimZ); +cdef extern int TV_PD_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, float lipschitz_const, int methodTV, int nonneg, int dimX, int dimY, int dimZ); cdef extern int TV_SB_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, int methodTV, int N, int M, int Z); cdef extern int LLT_ROF_GPU_main(float *Input, float *Output, float *infovector, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, float epsil, int N, int M, int Z); cdef extern int TGV_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, float epsil, int dimX, int dimY, int dimZ); @@ -72,11 +72,11 @@ def TV_FGP_GPU(inputData, methodTV, nonneg) # Total-variation Primal-Dual (PD) -def TV_PD_GPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau): +def TV_PD_GPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const): if inputData.ndim == 2: - return TVPD2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau) + return TVPD2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const) elif inputData.ndim == 3: - return TVPD3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const, tau) + return TVPD3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, lipschitz_const) def TVPD2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, float regularisation_parameter, @@ -84,8 +84,7 @@ def TVPD2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, float tolerance_param, int methodTV, int nonneg, - float lipschitz_const, - float tau): + float lipschitz_const): cdef long dims[2] dims[0] = inputData.shape[0] @@ -103,7 +102,6 @@ def TVPD2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, lipschitz_const, methodTV, nonneg, - tau, dims[1],dims[0], 1) ==0): return (outputData,infovec) else: @@ -115,8 +113,7 @@ def TVPD3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, float tolerance_param, int methodTV, int nonneg, - float lipschitz_const, - float tau): + float lipschitz_const): cdef long dims[3] dims[0] = inputData.shape[0] @@ -134,7 +131,6 @@ def TVPD3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, lipschitz_const, methodTV, nonneg, - tau, dims[2], dims[1], dims[0]) ==0): return (outputData,infovec) else: diff --git a/test/test_CPU_regularisers.py b/test/test_CPU_regularisers.py index 266ca8a..1eba479 100644 --- a/test/test_CPU_regularisers.py +++ b/test/test_CPU_regularisers.py @@ -42,7 +42,7 @@ class TestRegularisers(unittest.TestCase): def test_PD_TV_CPU(self): Im,input,ref = self.getPars() - pd_cpu,info = PD_TV(input, 0.02, 300, 0.0, 0, 0, 8, 0.0025, 'cpu'); + pd_cpu,info = PD_TV(input, 0.02, 300, 0.0, 0, 0, 8, 'cpu'); rms = rmse(Im, pd_cpu) diff --git a/test/test_run_test.py b/test/test_run_test.py index 1707aec..f562593 100755 --- a/test/test_run_test.py +++ b/test/test_run_test.py @@ -200,8 +200,7 @@ class TestRegularisers(unittest.TestCase): 'tolerance_constant':0.0,\ 'methodTV': 0 ,\ 'nonneg': 0, - 'lipschitz_const' : 8, - 'tau' : 0.0025} + 'lipschitz_const' : 8} print ("#############PD TV CPU####################") start_time = timeit.default_timer() @@ -211,8 +210,7 @@ class TestRegularisers(unittest.TestCase): pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'cpu') + pars['lipschitz_const'],'cpu') rms = rmse(Im, pd_cpu) pars['rmse'] = rms @@ -230,8 +228,7 @@ class TestRegularisers(unittest.TestCase): pars['tolerance_constant'], pars['methodTV'], pars['nonneg'], - pars['lipschitz_const'], - pars['tau'],'gpu') + pars['lipschitz_const'],'gpu') except ValueError as ve: self.skipTest("Results not comparable. GPU computing error.") -- cgit v1.2.3