diff options
-rw-r--r-- | Core/CMakeLists.txt | 3 | ||||
-rw-r--r-- | Core/regularizers_CPU/FGP_TV_core.c | 12 | ||||
-rw-r--r-- | Core/regularizers_CPU/FGP_TV_core.h | 2 | ||||
-rw-r--r-- | Core/regularizers_CPU/ROF_TV_core.c | 4 | ||||
-rw-r--r-- | Core/regularizers_CPU/ROF_TV_core.h | 2 | ||||
-rwxr-xr-x | Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu | 12 | ||||
-rwxr-xr-x | Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.h | 2 | ||||
-rw-r--r-- | Wrappers/Python/demo/test_cpu_regularizers.py | 25 | ||||
-rw-r--r-- | Wrappers/Python/src/cpu_regularizers.cpp | 1 | ||||
-rw-r--r-- | Wrappers/Python/src/cpu_regularizers.pyx | 25 | ||||
-rw-r--r-- | Wrappers/Python/src/gpu_regularizers.pyx | 2 | ||||
-rw-r--r-- | Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py | 9 |
12 files changed, 51 insertions, 48 deletions
diff --git a/Core/CMakeLists.txt b/Core/CMakeLists.txt index 3e3f89e..68b7111 100644 --- a/Core/CMakeLists.txt +++ b/Core/CMakeLists.txt @@ -129,7 +129,8 @@ if (CUDA_FOUND) CUDA_ADD_LIBRARY(cilregcuda SHARED ${CMAKE_CURRENT_SOURCE_DIR}/regularizers_GPU/Diffus_HO/Diff4th_GPU_kernel.cu ${CMAKE_CURRENT_SOURCE_DIR}/regularizers_GPU/NL_Regul/NLM_GPU_kernel.cu - ${CMAKE_CURRENT_SOURCE_DIR}/regularizers_GPU/TV_ROF/TV_ROF_GPU.cu + ${CMAKE_CURRENT_SOURCE_DIR}/regularizers_GPU/TV_ROF/TV_ROF_GPU_core.cu + ${CMAKE_CURRENT_SOURCE_DIR}/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu ) if (UNIX) message ("I'd install into ${CMAKE_INSTALL_PREFIX}/lib") diff --git a/Core/regularizers_CPU/FGP_TV_core.c b/Core/regularizers_CPU/FGP_TV_core.c index 2f1439d..5365631 100644 --- a/Core/regularizers_CPU/FGP_TV_core.c +++ b/Core/regularizers_CPU/FGP_TV_core.c @@ -23,7 +23,7 @@ limitations under the License. * * Input Parameters: * 1. Noisy image/volume - * 2. lambda - regularization parameter + * 2. lambdaPar - regularization parameter * 3. Number of iterations * 4. eplsilon: tolerance constant * 5. TV-type: methodTV - 'iso' (0) or 'l1' (1) @@ -37,7 +37,7 @@ limitations under the License. * [1] Amir Beck and Marc Teboulle, "Fast Gradient-Based Algorithms for Constrained Total Variation Image Denoising and Deblurring Problems" */ -float TV_FGP_CPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ) +float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ) { int ll, j, DimTotal; float re, re1; @@ -62,13 +62,13 @@ float TV_FGP_CPU(float *Input, float *Output, float lambda, int iter, float epsi for(ll=0; ll<iter; ll++) { /* computing the gradient of the objective function */ - Obj_func2D(Input, Output, R1, R2, lambda, dimX, dimY); + Obj_func2D(Input, Output, R1, R2, lambdaPar, dimX, dimY); /* apply nonnegativity */ if (nonneg == 1) for(j=0; j<dimX*dimY; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;} /*Taking a step towards minus of the gradient*/ - Grad_func2D(P1, P2, Output, R1, R2, lambda, dimX, dimY); + Grad_func2D(P1, P2, Output, R1, R2, lambdaPar, dimX, dimY); /* projection step */ Proj_func2D(P1, P2, methodTV, dimX, dimY); @@ -117,13 +117,13 @@ float TV_FGP_CPU(float *Input, float *Output, float lambda, int iter, float epsi for(ll=0; ll<iter; ll++) { /* computing the gradient of the objective function */ - Obj_func3D(Input, Output, R1, R2, R3, lambda, dimX, dimY, dimZ); + Obj_func3D(Input, Output, R1, R2, R3, lambdaPar, dimX, dimY, dimZ); /* apply nonnegativity */ if (nonneg == 1) for(j=0; j<dimX*dimY*dimZ; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;} /*Taking a step towards minus of the gradient*/ - Grad_func3D(P1, P2, P3, Output, R1, R2, R3, lambda, dimX, dimY, dimZ); + Grad_func3D(P1, P2, P3, Output, R1, R2, R3, lambdaPar, dimX, dimY, dimZ); /* projection step */ Proj_func3D(P1, P2, P3, methodTV, dimX, dimY, dimZ); diff --git a/Core/regularizers_CPU/FGP_TV_core.h b/Core/regularizers_CPU/FGP_TV_core.h index 98ceaec..360251d 100644 --- a/Core/regularizers_CPU/FGP_TV_core.h +++ b/Core/regularizers_CPU/FGP_TV_core.h @@ -47,7 +47,7 @@ limitations under the License. #ifdef __cplusplus extern "C" { #endif -CCPI_EXPORT float TV_FGP_CPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); +CCPI_EXPORT float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); CCPI_EXPORT float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, int dimX, int dimY); CCPI_EXPORT float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float lambda, int dimX, int dimY); diff --git a/Core/regularizers_CPU/ROF_TV_core.c b/Core/regularizers_CPU/ROF_TV_core.c index b2c6f00..a058f54 100644 --- a/Core/regularizers_CPU/ROF_TV_core.c +++ b/Core/regularizers_CPU/ROF_TV_core.c @@ -19,7 +19,7 @@ #include "ROF_TV_core.h" -#define EPS 1.0e-4 +#define EPS 1.0e-5 #define MAX(x, y) (((x) > (y)) ? (x) : (y)) #define MIN(x, y) (((x) < (y)) ? (x) : (y)) @@ -46,7 +46,7 @@ int sign(float x) { */ /* Running iterations of TV-ROF function */ -float TV_ROF_CPU(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float lambda) +float TV_ROF_CPU_main(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float lambda) { float *D1, *D2, *D3; int i, DimTotal; diff --git a/Core/regularizers_CPU/ROF_TV_core.h b/Core/regularizers_CPU/ROF_TV_core.h index b32d0d5..63c0419 100644 --- a/Core/regularizers_CPU/ROF_TV_core.h +++ b/Core/regularizers_CPU/ROF_TV_core.h @@ -47,7 +47,7 @@ limitations under the License. extern "C" { #endif CCPI_EXPORT float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambda, float tau, int dimY, int dimX, int dimZ); -CCPI_EXPORT float TV_ROF_CPU(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float lambda); +CCPI_EXPORT float TV_ROF_CPU_main(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float lambda); CCPI_EXPORT float D1_func(float *A, float *D1, int dimY, int dimX, int dimZ); CCPI_EXPORT float D2_func(float *A, float *D2, int dimY, int dimX, int dimZ); CCPI_EXPORT float D3_func(float *A, float *D3, int dimY, int dimX, int dimZ); diff --git a/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu b/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu index 0533a85..cbe9b5e 100755 --- a/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu +++ b/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.cu @@ -25,7 +25,7 @@ limitations under the License. * * Input Parameters: * 1. Noisy image/volume - * 2. lambda - regularization parameter + * 2. lambdaPar - regularization parameter * 3. Number of iterations * 4. eplsilon: tolerance constant * 5. TV-type: methodTV - 'iso' (0) or 'l1' (1) @@ -341,7 +341,7 @@ __global__ void copy_kernel3D(float *Input, float* Output, int N, int M, int Z, /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/ ////////////MAIN HOST FUNCTION /////////////// -extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ) +extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ) { int deviceCount = -1; // number of devices cudaGetDeviceCount(&deviceCount); @@ -383,13 +383,13 @@ extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambda, int iter, cudaMemset(R2, 0, ImSize*sizeof(float)); /********************** Run CUDA 2D kernel here ********************/ - multip = (1.0f/(8.0f*lambda)); + multip = (1.0f/(8.0f*lambdaPar)); /* The main kernel */ for (i = 0; i < iter; i++) { /* computing the gradient of the objective function */ - Obj_func2D_kernel<<<dimGrid,dimBlock>>>(d_input, d_update, R1, R2, dimX, dimY, ImSize, lambda); + Obj_func2D_kernel<<<dimGrid,dimBlock>>>(d_input, d_update, R1, R2, dimX, dimY, ImSize, lambdaPar); checkCudaErrors( cudaDeviceSynchronize() ); checkCudaErrors(cudaPeekAtLastError() ); @@ -493,13 +493,13 @@ extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambda, int iter, cudaMemset(R2, 0, ImSize*sizeof(float)); cudaMemset(R3, 0, ImSize*sizeof(float)); /********************** Run CUDA 3D kernel here ********************/ - multip = (1.0f/(8.0f*lambda)); + multip = (1.0f/(8.0f*lambdaPar)); /* The main kernel */ for (i = 0; i < iter; i++) { /* computing the gradient of the objective function */ - Obj_func3D_kernel<<<dimGrid,dimBlock>>>(d_input, d_update, R1, R2, R3, dimX, dimY, dimZ, ImSize, lambda); + Obj_func3D_kernel<<<dimGrid,dimBlock>>>(d_input, d_update, R1, R2, R3, dimX, dimY, dimZ, ImSize, lambdaPar); checkCudaErrors( cudaDeviceSynchronize() ); checkCudaErrors(cudaPeekAtLastError() ); diff --git a/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.h b/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.h index 15c7120..2b8fdcf 100755 --- a/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.h +++ b/Core/regularizers_GPU/TV_FGP/TV_FGP_GPU_core.h @@ -5,6 +5,6 @@ #ifndef _TV_FGP_GPU_ #define _TV_FGP_GPU_ -extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); +extern "C" void TV_FGP_GPU(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); #endif diff --git a/Wrappers/Python/demo/test_cpu_regularizers.py b/Wrappers/Python/demo/test_cpu_regularizers.py index b08c418..53b8538 100644 --- a/Wrappers/Python/demo/test_cpu_regularizers.py +++ b/Wrappers/Python/demo/test_cpu_regularizers.py @@ -11,10 +11,9 @@ import numpy as np import os from enum import Enum import timeit -from ccpi.filters.cpu_regularizers_boost import SplitBregman_TV , FGP_TV ,\ - LLT_model, PatchBased_Regul ,\ +from ccpi.filters.cpu_regularizers_boost import SplitBregman_TV, LLT_model, PatchBased_Regul ,\ TGV_PD -from ccpi.filters.cpu_regularizers_cython import ROF_TV +from ccpi.filters.cpu_regularizers_cython import TV_ROF_CPU, TV_FGP_CPU ############################################################################### #https://stackoverflow.com/questions/13875989/comparing-image-in-url-to-image-in-filesystem-in-python/13884956#13884956 @@ -128,21 +127,25 @@ imgplot = plt.imshow(splitbregman,\ ) ###################### FGP_TV ######################################### -# u = FGP_TV(single(u0), 0.05, 100, 1e-04); + start_time = timeit.default_timer() -pars = {'algorithm' : FGP_TV , \ +pars = {'algorithm' : TV_FGP_CPU , \ 'input' : u0, 'regularization_parameter':0.05, \ 'number_of_iterations' :200 ,\ - 'tolerance_constant':1e-4,\ - 'TV_penalty': 0 + 'tolerance_constant':1e-5,\ + 'methodTV': 0 ,\ + 'nonneg': 0 ,\ + 'printingOut': 0 } -out = FGP_TV (pars['input'], +out = TV_FGP_CPU (pars['input'], pars['regularization_parameter'], pars['number_of_iterations'], pars['tolerance_constant'], - pars['TV_penalty']) + pars['methodTV'], + pars['nonneg'], + pars['printingOut']) fgp = out[0] rms = rmse(Im, fgp) @@ -282,13 +285,13 @@ imgplot = plt.imshow(tgv, cmap="gray") start_time = timeit.default_timer() -pars = {'algorithm': ROF_TV , \ +pars = {'algorithm': TV_ROF_CPU , \ 'input' : u0,\ 'regularization_parameter':0.04,\ 'marching_step': 0.0025,\ 'number_of_iterations': 300 } -rof = ROF_TV(pars['input'], +rof = TV_ROF_CPU(pars['input'], pars['number_of_iterations'], pars['regularization_parameter'], pars['marching_step'] diff --git a/Wrappers/Python/src/cpu_regularizers.cpp b/Wrappers/Python/src/cpu_regularizers.cpp index 43d5d11..b8156d5 100644 --- a/Wrappers/Python/src/cpu_regularizers.cpp +++ b/Wrappers/Python/src/cpu_regularizers.cpp @@ -1040,7 +1040,6 @@ BOOST_PYTHON_MODULE(cpu_regularizers_boost) np::dtype dt2 = np::dtype::get_builtin<uint16_t>(); def("SplitBregman_TV", SplitBregman_TV); - def("FGP_TV", FGP_TV); def("LLT_model", LLT_model); def("PatchBased_Regul", PatchBased_Regul); def("TGV_PD", TGV_PD); diff --git a/Wrappers/Python/src/cpu_regularizers.pyx b/Wrappers/Python/src/cpu_regularizers.pyx index b8089a8..448da31 100644 --- a/Wrappers/Python/src/cpu_regularizers.pyx +++ b/Wrappers/Python/src/cpu_regularizers.pyx @@ -18,8 +18,8 @@ import cython import numpy as np cimport numpy as np -cdef extern float TV_ROF_CPU(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float flambda); -cdef extern float TV_FGP_CPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); +cdef extern float TV_ROF_CPU_main(float *Input, float *Output, int dimX, int dimY, int dimZ, int iterationsNumb, float tau, float flambda); +cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ); # Can we use the same name here in "def" as the C function? def TV_ROF_CPU(inputData, iterations, regularization_parameter, @@ -45,11 +45,10 @@ def TV_ROF_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, np.zeros([dims[0],dims[1]], dtype='float32') #/* Run ROF iterations for 2D data */ - TV_ROF_CPU(&inputData[0,0], &B[0,0], dims[0], dims[1], 1, iterations, + TV_ROF_CPU_main(&inputData[0,0], &B[0,0], dims[0], dims[1], 1, iterations, marching_step_parameter, regularization_parameter) - return B - + return B def TV_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, int iterations, @@ -65,7 +64,7 @@ def TV_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, np.zeros([dims[0],dims[1],dims[2]], dtype='float32') #/* Run ROF iterations for 3D data */ - TV_FGP_CPU(&inputData[0,0,0], &B[0,0,0], dims[0], dims[1], dims[2], iterations, + TV_ROF_CPU_main(&inputData[0,0,0], &B[0,0,0], dims[0], dims[1], dims[2], iterations, marching_step_parameter, regularization_parameter) return B @@ -88,11 +87,11 @@ def TV_FGP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, dims[0] = inputData.shape[0] dims[1] = inputData.shape[1] - cdef np.ndarray[np.float32_t, ndim=2, mode="c"] B = \ + cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \ np.zeros([dims[0],dims[1]], dtype='float32') #/* Run ROF iterations for 2D data */ - TV_FGP_CPU(&inputData[0,0], &B[0,0], regularization_parameter, + TV_FGP_CPU_main(&inputData[0,0], &outputData[0,0], regularization_parameter, iterations, tolerance_param, methodTV, @@ -100,8 +99,7 @@ def TV_FGP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, printM, dims[0], dims[1], 1) - return B - + return outputData def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, float regularization_parameter, @@ -115,15 +113,16 @@ def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, dims[1] = inputData.shape[1] dims[2] = inputData.shape[2] - cdef np.ndarray[np.float32_t, ndim=3, mode="c"] B = \ + cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \ np.zeros([dims[0],dims[1],dims[2]], dtype='float32') #/* Run ROF iterations for 3D data */ - TV_FGP_CPU(&inputData[0,0, 0], &B[0,0, 0], iterations, + TV_FGP_CPU_main(&inputData[0,0,0], &outputData[0,0, 0], regularization_parameter, + iterations, tolerance_param, methodTV, nonneg, printM, dims[0], dims[1], dims[2]) - return B + return outputData diff --git a/Wrappers/Python/src/gpu_regularizers.pyx b/Wrappers/Python/src/gpu_regularizers.pyx index a14a20d..e99bfa7 100644 --- a/Wrappers/Python/src/gpu_regularizers.pyx +++ b/Wrappers/Python/src/gpu_regularizers.pyx @@ -26,7 +26,7 @@ cdef extern void NLM_GPU_kernel(float *A, float* B, float *Eucl_Vec, int SearchW, int SimilW, int SearchW_real, float denh2, float lambdaf); cdef extern void TV_ROF_GPU(float* Input, float* Output, int N, int M, int Z, int iter, float tau, float lambdaf); -cdef extern void TV_FGP_GPU(float *Input, float *Output, float lambda, int iter, float epsil, int methodTV, int nonneg, int printM, int N, int M, int Z); +cdef extern void TV_FGP_GPU(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int N, int M, int Z); cdef extern float pad_crop(float *A, float *Ap, int OldSizeX, int OldSizeY, int OldSizeZ, diff --git a/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py b/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py index 6344021..e162afa 100644 --- a/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py +++ b/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py @@ -12,7 +12,7 @@ import matplotlib.pyplot as plt import numpy as np import os import timeit -from ccpi.filters.gpu_regularizers import Diff4thHajiaboli, NML, TV_ROF_GPU +from ccpi.filters.gpu_regularizers import Diff4thHajiaboli, NML, GPU_ROF_TV from ccpi.filters.cpu_regularizers_cython import TV_ROF_CPU ############################################################################### def printParametersToString(pars): @@ -56,11 +56,11 @@ imgplot = plt.imshow(u0,cmap="gray") # set parameters -pars = {'algorithm': ROF_TV , \ +pars = {'algorithm': TV_ROF_CPU , \ 'input' : u0,\ 'regularization_parameter':0.04,\ 'time_marching_parameter': 0.0025,\ - 'number_of_iterations': 600 + 'number_of_iterations': 1200 } print ("#################ROF TV CPU#####################") start_time = timeit.default_timer() @@ -89,13 +89,14 @@ plt.title('{}'.format('CPU results')) print ("#################ROF TV GPU#####################") start_time = timeit.default_timer() -rof_gpu = TV_ROF_GPU(pars['input'], +rof_gpu = GPU_ROF_TV(pars['input'], pars['number_of_iterations'], pars['time_marching_parameter'], pars['regularization_parameter']) rms = rmse(Im, rof_gpu) pars['rmse'] = rms +pars['algorithm'] = GPU_ROF_TV txtstr = printParametersToString(pars) txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) print (txtstr) |