diff options
author | Daniil Kazantsev <dkazanc3@googlemail.com> | 2018-04-09 15:17:24 +0100 |
---|---|---|
committer | GitHub <noreply@github.com> | 2018-04-09 15:17:24 +0100 |
commit | 62635199f4e5a464a267ffce070ecec68bfdcfe8 (patch) | |
tree | cdc7c4469e210a52cb416b2747ca2d954da073cc /Wrappers | |
parent | a5b5872b76bf00023a7e7cee97e028003ccbc45e (diff) | |
parent | b9fafd363d1d181a4a8b42ea4038924097207913 (diff) | |
download | regularization-62635199f4e5a464a267ffce070ecec68bfdcfe8.tar.gz regularization-62635199f4e5a464a267ffce070ecec68bfdcfe8.tar.bz2 regularization-62635199f4e5a464a267ffce070ecec68bfdcfe8.tar.xz regularization-62635199f4e5a464a267ffce070ecec68bfdcfe8.zip |
Merge pull request #47 from vais-ral/add3Dtests
major renaming and new 3D demos for Matlab
Diffstat (limited to 'Wrappers')
27 files changed, 229 insertions, 829 deletions
diff --git a/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m b/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m new file mode 100644 index 0000000..f5c3ad1 --- /dev/null +++ b/Wrappers/Matlab/demos/demoMatlab_3Ddenoise.m @@ -0,0 +1,44 @@ +% Volume (3D) denoising demo using CCPi-RGL + +addpath('../mex_compile/installed'); +addpath('../../../data/'); + +N = 256; +slices = 30; +vol3D = zeros(N,N,slices, 'single'); +Im = double(imread('lena_gray_256.tif'))/255; % loading image +for i = 1:slices +vol3D(:,:,i) = Im + .05*randn(size(Im)); +end +vol3D(vol3D < 0) = 0; +figure; imshow(vol3D(:,:,15), [0 1]); title('Noisy image'); + +%% +fprintf('Denoise using ROF-TV model (CPU) \n'); +lambda_rof = 0.03; % regularisation parameter +tau_rof = 0.0025; % time-marching constant +iter_rof = 1000; % number of ROF iterations +tic; u_rof = ROF_TV(single(vol3D), lambda_rof, iter_rof, tau_rof); toc; +figure; imshow(u_rof(:,:,15), [0 1]); title('ROF-TV denoised volume (CPU)'); +%% +% fprintf('Denoise using ROF-TV model (GPU) \n'); +% lambda_rof = 0.03; % regularisation parameter +% tau_rof = 0.0025; % time-marching constant +% iter_rof = 1000; % number of ROF iterations +% tic; u_rofG = ROF_TV_GPU(single(vol3D), lambda_rof, iter_rof, tau_rof); toc; +% figure; imshow(u_rofG(:,:,15), [0 1]); title('ROF-TV denoised volume (GPU)'); +%% +fprintf('Denoise using FGP-TV model (CPU) \n'); +lambda_fgp = 0.03; % regularisation parameter +iter_fgp = 500; % number of FGP iterations +epsil_tol = 1.0e-05; % tolerance +tic; u_fgp = FGP_TV(single(vol3D), lambda_fgp, iter_fgp, epsil_tol); toc; +figure; imshow(u_fgp(:,:,15), [0 1]); title('FGP-TV denoised volume (CPU)'); +%% +% fprintf('Denoise using FGP-TV model (GPU) \n'); +% lambda_fgp = 0.03; % regularisation parameter +% iter_fgp = 500; % number of FGP iterations +% epsil_tol = 1.0e-05; % tolerance +% tic; u_fgpG = FGP_TV_GPU(single(vol3D), lambda_fgp, iter_fgp, epsil_tol); toc; +% figure; imshow(u_fgpG(:,:,15), [0 1]); title('FGP-TV denoised volume (GPU)'); +%%
\ No newline at end of file diff --git a/Wrappers/Matlab/demos/demoMatlab_denoise.m b/Wrappers/Matlab/demos/demoMatlab_denoise.m index 7258e5e..ab4e95d 100644 --- a/Wrappers/Matlab/demos/demoMatlab_denoise.m +++ b/Wrappers/Matlab/demos/demoMatlab_denoise.m @@ -9,28 +9,28 @@ figure; imshow(u0, [0 1]); title('Noisy image'); %% fprintf('Denoise using ROF-TV model (CPU) \n'); -lambda_rof = 0.03; % regularization parameter +lambda_rof = 0.03; % regularisation parameter tau_rof = 0.0025; % time-marching constant iter_rof = 2000; % number of ROF iterations tic; u_rof = ROF_TV(single(u0), lambda_rof, iter_rof, tau_rof); toc; figure; imshow(u_rof, [0 1]); title('ROF-TV denoised image (CPU)'); %% % fprintf('Denoise using ROF-TV model (GPU) \n'); -% lambda_rof = 0.03; % regularization parameter +% lambda_rof = 0.03; % regularisation parameter % tau_rof = 0.0025; % time-marching constant % iter_rof = 2000; % number of ROF iterations % tic; u_rofG = ROF_TV_GPU(single(u0), lambda_rof, iter_rof, tau_rof); toc; % figure; imshow(u_rofG, [0 1]); title('ROF-TV denoised image (GPU)'); %% fprintf('Denoise using FGP-TV model (CPU) \n'); -lambda_fgp = 0.03; % regularization parameter +lambda_fgp = 0.03; % regularisation parameter iter_fgp = 1000; % number of FGP iterations epsil_tol = 1.0e-05; % tolerance tic; u_fgp = FGP_TV(single(u0), lambda_fgp, iter_fgp, epsil_tol); toc; figure; imshow(u_fgp, [0 1]); title('FGP-TV denoised image (CPU)'); %% % fprintf('Denoise using FGP-TV model (GPU) \n'); -% lambda_fgp = 0.03; % regularization parameter +% lambda_fgp = 0.03; % regularisation parameter % iter_fgp = 1000; % number of FGP iterations % epsil_tol = 1.0e-05; % tolerance % tic; u_fgpG = FGP_TV_GPU(single(u0), lambda_fgp, iter_fgp, epsil_tol); toc; diff --git a/Wrappers/Matlab/mex_compile/compileCPU_mex.m b/Wrappers/Matlab/mex_compile/compileCPU_mex.m index fcee53a..8da81ad 100644 --- a/Wrappers/Matlab/mex_compile/compileCPU_mex.m +++ b/Wrappers/Matlab/mex_compile/compileCPU_mex.m @@ -1,10 +1,10 @@ % execute this mex file in Matlab once -copyfile ../../../Core/regularizers_CPU/ regularizers_CPU/ -copyfile ../../../Core/CCPiDefines.h regularizers_CPU/ +copyfile ../../../Core/regularisers_CPU/ regularisers_CPU/ +copyfile ../../../Core/CCPiDefines.h regularisers_CPU/ -cd regularizers_CPU/ +cd regularisers_CPU/ -fprintf('%s \n', 'Compiling CPU regularizers...'); +fprintf('%s \n', 'Compiling CPU regularisers...'); mex ROF_TV.c ROF_TV_core.c utils.c CFLAGS="\$CFLAGS -fopenmp -Wall -std=c99" LDFLAGS="\$LDFLAGS -fopenmp" movefile ROF_TV.mex* ../installed/ diff --git a/Wrappers/Matlab/mex_compile/compileGPU_mex.m b/Wrappers/Matlab/mex_compile/compileGPU_mex.m index df29a3e..45236fa 100644 --- a/Wrappers/Matlab/mex_compile/compileGPU_mex.m +++ b/Wrappers/Matlab/mex_compile/compileGPU_mex.m @@ -9,12 +9,12 @@ % tested on Ubuntu 16.04/MATLAB 2016b -copyfile ../../../Core/regularizers_GPU/ regularizers_GPU/ -copyfile ../../../Core/CCPiDefines.h regularizers_GPU/ +copyfile ../../../Core/regularisers_GPU/ regularisers_GPU/ +copyfile ../../../Core/CCPiDefines.h regularisers_GPU/ -cd regularizers_GPU/ +cd regularisers_GPU/ -fprintf('%s \n', 'Compiling GPU regularizers (CUDA)...'); +fprintf('%s \n', 'Compiling GPU regularisers (CUDA)...'); !/usr/local/cuda/bin/nvcc -O0 -c TV_ROF_GPU_core.cu -Xcompiler -fPIC -I~/SOFT/MATLAB9/extern/include/ mex -g -I/usr/local/cuda-7.5/include -L/usr/local/cuda-7.5/lib64 -lcudart -lcufft -lmwgpu ROF_TV_GPU.cpp TV_ROF_GPU_core.o movefile ROF_TV_GPU.mex* ../installed/ diff --git a/Wrappers/Matlab/mex_compile/regularizers_CPU/FGP_TV.c b/Wrappers/Matlab/mex_compile/regularisers_CPU/FGP_TV.c index ba06cc7..ba06cc7 100644 --- a/Wrappers/Matlab/mex_compile/regularizers_CPU/FGP_TV.c +++ b/Wrappers/Matlab/mex_compile/regularisers_CPU/FGP_TV.c diff --git a/Wrappers/Matlab/mex_compile/regularisers_CPU/FGP_TV.c~ b/Wrappers/Matlab/mex_compile/regularisers_CPU/FGP_TV.c~ new file mode 100644 index 0000000..30d61cd --- /dev/null +++ b/Wrappers/Matlab/mex_compile/regularisers_CPU/FGP_TV.c~ @@ -0,0 +1,91 @@ +/* + * This work is part of the Core Imaging Library developed by + * Visual Analytics and Imaging System Group of the Science Technology + * Facilities Council, STFC + * + * Copyright 2017 Daniil Kazantsev + * Copyright 2017 Srikanth Nagella, Edoardo Pasca + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * http://www.apache.org/licenses/LICENSE-2.0 + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ +#include "matrix.h" +#include "mex.h" +#include "FGP_TV_core.h" + +/* C-OMP implementation of FGP-TV [1] denoising/regularization model (2D/3D case) + * + * Input Parameters: + * 1. Noisy image/volume + * 2. lambdaPar - regularization parameter + * 3. Number of iterations + * 4. eplsilon: tolerance constant + * 5. TV-type: methodTV - 'iso' (0) or 'l1' (1) + * 6. nonneg: 'nonnegativity (0 is OFF by default) + * 7. print information: 0 (off) or 1 (on) + * + * Output: + * [1] Filtered/regularized image + * + * This function is based on the Matlab's code and paper by + * [1] Amir Beck and Marc Teboulle, "Fast Gradient-Based Algorithms for Constrained Total Variation Image Denoising and Deblurring Problems" + */ + + +void mexFunction( + int nlhs, mxArray *plhs[], + int nrhs, const mxArray *prhs[]) + +{ + int number_of_dims, iter, dimX, dimY, dimZ, methTV, printswitch; + const int *dim_array; + float *Input, *Output, lambda, epsil; + + number_of_dims = mxGetNumberOfDimensions(prhs[0]); + dim_array = mxGetDimensions(prhs[0]); + + /*Handling Matlab input data*/ + if ((nrhs < 2) || (nrhs > 6)) mexErrMsgTxt("At least 2 parameters is required: Image(2D/3D), Regularization parameter. The full list of parameters: Image(2D/3D), Regularization parameter, iterations number, tolerance, penalty type ('iso' or 'l1'), print switch"); + + Input = (float *) mxGetData(prhs[0]); /*noisy image (2D/3D) */ + lambda = (float) mxGetScalar(prhs[1]); /* regularization parameter */ + iter = 300; /* default iterations number */ + epsil = 0.0001; /* default tolerance constant */ + methTV = 0; /* default isotropic TV penalty */ + printswitch = 0; /*default print is switched off - 0 */ + + 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)) iter = (int) mxGetScalar(prhs[2]); /* iterations number */ + if ((nrhs == 4) || (nrhs == 5) || (nrhs == 6)) epsil = (float) mxGetScalar(prhs[3]); /* tolerance constant */ + if ((nrhs == 5) || (nrhs == 6)) { + 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) { + printswitch = (int) mxGetScalar(prhs[5]); + if ((printswitch != 0) || (printswitch != 1)) {mexErrMsgTxt("Print can be enabled by choosing 1 or off - 0"); } + } + + /*Handling Matlab output data*/ + 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)); + } + if (number_of_dims == 3) Output = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(3, dim_array, mxSINGLE_CLASS, mxREAL)); + + + TV_FGP_CPU_main(Input, Output, lambda, iter, epsil, methTV, nonneg, printswitch, dimX, dimY, dimZ) +} diff --git a/Wrappers/Matlab/mex_compile/regularizers_CPU/ROF_TV.c b/Wrappers/Matlab/mex_compile/regularisers_CPU/ROF_TV.c index 6b9e1ea..6b9e1ea 100644 --- a/Wrappers/Matlab/mex_compile/regularizers_CPU/ROF_TV.c +++ b/Wrappers/Matlab/mex_compile/regularisers_CPU/ROF_TV.c diff --git a/Wrappers/Matlab/mex_compile/regularizers_GPU/FGP_TV_GPU.cpp b/Wrappers/Matlab/mex_compile/regularisers_GPU/FGP_TV_GPU.cpp index 9ed9ae0..9ed9ae0 100644 --- a/Wrappers/Matlab/mex_compile/regularizers_GPU/FGP_TV_GPU.cpp +++ b/Wrappers/Matlab/mex_compile/regularisers_GPU/FGP_TV_GPU.cpp diff --git a/Wrappers/Matlab/mex_compile/regularizers_GPU/ROF_TV_GPU.cpp b/Wrappers/Matlab/mex_compile/regularisers_GPU/ROF_TV_GPU.cpp index 7bbe3af..7bbe3af 100644 --- a/Wrappers/Matlab/mex_compile/regularizers_GPU/ROF_TV_GPU.cpp +++ b/Wrappers/Matlab/mex_compile/regularisers_GPU/ROF_TV_GPU.cpp diff --git a/Wrappers/Python/CMakeLists.txt b/Wrappers/Python/CMakeLists.txt index eb3bac7..fb00706 100644 --- a/Wrappers/Python/CMakeLists.txt +++ b/Wrappers/Python/CMakeLists.txt @@ -1,7 +1,7 @@ # Copyright 2018 Edoardo Pasca cmake_minimum_required (VERSION 3.0) -project(RegularizerPython) +project(regulariserPython) #https://stackoverflow.com/questions/13298504/using-cmake-with-setup-py # The version number. @@ -25,8 +25,8 @@ if (PYTHONINTERP_FOUND) endif() -## Build the regularizers package as a library -message("Creating Regularizers as shared library") +## Build the regularisers package as a library +message("Creating Regularisers as shared library") message("CMAKE_CXX_FLAGS ${CMAKE_CXX_FLAGS}") @@ -52,7 +52,7 @@ elseif(UNIX) ) endif() -# GPU Regularizers +# GPU regularisers find_package(CUDA) if (CUDA_FOUND) @@ -60,12 +60,12 @@ if (CUDA_FOUND) set (SETUP_GPU_WRAPPERS "extra_libraries += ['cilregcuda']\n\ setup( \n\ name='ccpi', \n\ - description='CCPi Core Imaging Library - Image Regularizers GPU',\n\ + description='CCPi Core Imaging Library - Image regularisers GPU',\n\ version=cil_version,\n\ cmdclass = {'build_ext': build_ext},\n\ - ext_modules = [Extension('ccpi.filters.gpu_regularizers',\n\ + ext_modules = [Extension('ccpi.filters.gpu_regularisers',\n\ sources=[ \n\ - os.path.join('.' , 'src', 'gpu_regularizers.pyx' ),\n\ + os.path.join('.' , 'src', 'gpu_regularisers.pyx' ),\n\ ],\n\ include_dirs=extra_include_dirs, \n\ library_dirs=extra_library_dirs, \n\ @@ -80,9 +80,9 @@ else() set(SETUP_GPU_WRAPPERS "#CUDA NOT FOUND") endif() -configure_file("setup-regularizers.py.in" "setup-regularizers.py") +configure_file("setup-regularisers.py.in" "setup-regularisers.py") -#add_executable(regularizer_test ${CMAKE_CURRENT_SOURCE_DIR}/test/test_regularizer.cpp) +#add_executable(regulariser_test ${CMAKE_CURRENT_SOURCE_DIR}/test/test_regulariser.cpp) -#target_link_libraries (regularizer_test LINK_PUBLIC regularizers_lib) +#target_link_libraries (regulariser_test LINK_PUBLIC regularisers_lib) diff --git a/Wrappers/Python/ccpi/filters/regularizers.py b/Wrappers/Python/ccpi/filters/regularisers.py index d6dfa8c..039daab 100644 --- a/Wrappers/Python/ccpi/filters/regularizers.py +++ b/Wrappers/Python/ccpi/filters/regularisers.py @@ -2,30 +2,30 @@ script which assigns a proper device core function based on a flag ('cpu' or 'gpu') """ -from ccpi.filters.cpu_regularizers_cython import TV_ROF_CPU, TV_FGP_CPU -from ccpi.filters.gpu_regularizers import TV_ROF_GPU, TV_FGP_GPU +from ccpi.filters.cpu_regularisers_cython import TV_ROF_CPU, TV_FGP_CPU +from ccpi.filters.gpu_regularisers import TV_ROF_GPU, TV_FGP_GPU -def ROF_TV(inputData, regularization_parameter, iterations, +def ROF_TV(inputData, regularisation_parameter, iterations, time_marching_parameter,device='cpu'): if device == 'cpu': return TV_ROF_CPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, time_marching_parameter) elif device == 'gpu': return TV_ROF_GPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, time_marching_parameter) else: raise ValueError('Unknown device {0}. Expecting gpu or cpu'\ .format(device)) -def FGP_TV(inputData, regularization_parameter,iterations, +def FGP_TV(inputData, regularisation_parameter,iterations, tolerance_param, methodTV, nonneg, printM, device='cpu'): if device == 'cpu': return TV_FGP_CPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -33,7 +33,7 @@ def FGP_TV(inputData, regularization_parameter,iterations, printM) elif device == 'gpu': return TV_FGP_GPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -41,4 +41,4 @@ def FGP_TV(inputData, regularization_parameter,iterations, printM) else: raise ValueError('Unknown device {0}. Expecting gpu or cpu'\ - .format(device))
\ No newline at end of file + .format(device)) diff --git a/Wrappers/Python/conda-recipe/bld.bat b/Wrappers/Python/conda-recipe/bld.bat index 850905c..e47f8d9 100644 --- a/Wrappers/Python/conda-recipe/bld.bat +++ b/Wrappers/Python/conda-recipe/bld.bat @@ -11,7 +11,7 @@ cd %SRC_DIR%\ccpi\Python :: issue cmake to create setup.py cmake . -%PYTHON% setup-regularizers.py build_ext +%PYTHON% setup-regularisers.py build_ext if errorlevel 1 exit 1 -%PYTHON% setup-regularizers.py install +%PYTHON% setup-regularisers.py install if errorlevel 1 exit 1 diff --git a/Wrappers/Python/conda-recipe/build.sh b/Wrappers/Python/conda-recipe/build.sh index 9ea4161..8b05663 100644 --- a/Wrappers/Python/conda-recipe/build.sh +++ b/Wrappers/Python/conda-recipe/build.sh @@ -13,7 +13,7 @@ echo "$SRC_DIR/ccpi/Python" cmake . -$PYTHON setup-regularizers.py build_ext -$PYTHON setup-regularizers.py install +$PYTHON setup-regularisers.py build_ext +$PYTHON setup-regularisers.py install diff --git a/Wrappers/Python/conda-recipe/meta.yaml b/Wrappers/Python/conda-recipe/meta.yaml index f4cb471..5336d14 100644 --- a/Wrappers/Python/conda-recipe/meta.yaml +++ b/Wrappers/Python/conda-recipe/meta.yaml @@ -1,5 +1,5 @@ package: - name: ccpi-regularizer + name: ccpi-regulariser version: {{ environ['CIL_VERSION'] }} @@ -17,7 +17,7 @@ requirements: #- boost ==1.64.0 #- boost-cpp ==1.64.0 - cython - - cil_regularizer + - cil_regulariser - vc 14 # [win and py35] - vc 9 # [win and py27] - cmake @@ -26,7 +26,7 @@ requirements: - python - numpy x.x #- boost ==1.64 - - cil_regularizer + - cil_regulariser - vc 14 # [win and py35] - vc 9 # [win and py27] diff --git a/Wrappers/Python/demo/test_cpu_regularizers.py b/Wrappers/Python/demo/test_cpu_regularisers.py index 76b9de7..4e4a2dd 100644 --- a/Wrappers/Python/demo/test_cpu_regularizers.py +++ b/Wrappers/Python/demo/test_cpu_regularisers.py @@ -5,24 +5,13 @@ Created on Fri Aug 4 11:10:05 2017 @author: ofn77899 """ - import matplotlib.pyplot as plt import numpy as np import os from enum import Enum import timeit -from ccpi.filters.cpu_regularizers_boost import SplitBregman_TV, LLT_model, PatchBased_Regul, TGV_PD -from ccpi.filters.regularizers import ROF_TV, FGP_TV +from ccpi.filters.regularisers import ROF_TV, FGP_TV ############################################################################### -#https://stackoverflow.com/questions/13875989/comparing-image-in-url-to-image-in-filesystem-in-python/13884956#13884956 -#NRMSE a normalization of the root of the mean squared error -#NRMSE is simply 1 - [RMSE / (maxval - minval)]. Where maxval is the maximum -# intensity from the two images being compared, and respectively the same for -# minval. RMSE is given by the square root of MSE: -# sqrt[(sum(A - B) ** 2) / |A|], -# where |A| means the number of elements in A. By doing this, the maximum value -# given by RMSE is maxval. - def nrmse(im1, im2): a, b = im1.shape rmse = np.sqrt(np.sum((im2 - im1) ** 2) / float(a * b)) @@ -51,17 +40,8 @@ def printParametersToString(pars): # 2D Regularizers # ############################################################################### -#Example: -# figure; -# Im = double(imread('lena_gray_256.tif'))/255; % loading image -# u0 = Im + .05*randn(size(Im)); u0(u0 < 0) = 0; -# u = SplitBregman_TV(single(u0), 10, 30, 1e-04); - # assumes the script is launched from the test directory filename = os.path.join(".." , ".." , ".." , "data" ,"lena_gray_512.tif") -#filename = r"C:\Users\ofn77899\Documents\GitHub\CCPi-FISTA_reconstruction\data\lena_gray_512.tif" -#filename = r"/home/ofn77899/Reconstruction/CCPi-FISTA_Reconstruction/data/lena_gray_512.tif" -#filename = r'/home/algol/Documents/Python/STD_test_images/lena_gray_512.tif' Im = plt.imread(filename) Im = np.asarray(Im, dtype='float32') @@ -86,48 +66,14 @@ imgplot = plt.imshow(u0,cmap="gray" reg_output = [] ############################################################################## -# Call regularizer - -####################### SplitBregman_TV ##################################### -# u = SplitBregman_TV(single(u0), 10, 30, 1e-04); - -start_time = timeit.default_timer() -pars = {'algorithm' : SplitBregman_TV , \ - 'input' : u0, - 'regularization_parameter':15. , \ - 'number_of_iterations' :40 ,\ - 'tolerance_constant':0.0001 , \ - 'TV_penalty': 0 -} - -out = SplitBregman_TV (pars['input'], pars['regularization_parameter'], - pars['number_of_iterations'], - pars['tolerance_constant'], - pars['TV_penalty']) -splitbregman = out[0] -rms = rmse(Im, splitbregman) -pars['rmse'] = rms -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) - -a=fig.add_subplot(2,4,2) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(splitbregman,\ - cmap="gray" - ) +# Call regularisers ###################### FGP_TV ######################################### start_time = timeit.default_timer() pars = {'algorithm' : FGP_TV , \ 'input' : u0,\ - 'regularization_parameter':0.07, \ + 'regularisation_parameter':0.07, \ 'number_of_iterations' :300 ,\ 'tolerance_constant':0.00001,\ 'methodTV': 0 ,\ @@ -136,7 +82,7 @@ pars = {'algorithm' : FGP_TV , \ } fgp = FGP_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['tolerance_constant'], pars['methodTV'], @@ -163,129 +109,18 @@ imgplot = plt.imshow(fgp, \ a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, verticalalignment='top', bbox=props) - -###################### LLT_model ######################################### -start_time = timeit.default_timer() - -pars = {'algorithm': LLT_model , \ - 'input' : u0, - 'regularization_parameter': 5,\ - 'time_step':0.00035, \ - 'number_of_iterations' :350,\ - 'tolerance_constant':0.0001,\ - 'restrictive_Z_smoothing': 0 -} -out = LLT_model(pars['input'], - pars['regularization_parameter'], - pars['time_step'] , - pars['number_of_iterations'], - pars['tolerance_constant'], - pars['restrictive_Z_smoothing'] ) - -llt = out[0] -rms = rmse(Im, out[0]) -pars['rmse'] = rms - -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,4,4) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(llt,\ - cmap="gray" - ) -# ###################### PatchBased_Regul ######################################### -# # Quick 2D denoising example in Matlab: -# # Im = double(imread('lena_gray_256.tif'))/255; % loading image -# # u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -# # ImDen = PB_Regul_CPU(single(u0), 3, 1, 0.08, 0.05); - -start_time = timeit.default_timer() - -pars = {'algorithm': PatchBased_Regul , \ - 'input' : u0, - 'regularization_parameter': 0.05,\ - 'searching_window_ratio':3, \ - 'similarity_window_ratio':1,\ - 'PB_filtering_parameter': 0.06 -} -out = PatchBased_Regul(pars['input'], - pars['regularization_parameter'], - pars['searching_window_ratio'] , - pars['similarity_window_ratio'] , - pars['PB_filtering_parameter']) -pbr = out[0] -rms = rmse(Im, out[0]) -pars['rmse'] = rms - -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) - -a=fig.add_subplot(2,4,5) - - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(pbr ,cmap="gray") - - -# ###################### TGV_PD ######################################### -# # Quick 2D denoising example in Matlab: -# # Im = double(imread('lena_gray_256.tif'))/255; % loading image -# # u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -# # u = PrimalDual_TGV(single(u0), 0.02, 1.3, 1, 550); - -start_time = timeit.default_timer() - -pars = {'algorithm': TGV_PD , \ - 'input' : u0,\ - 'regularization_parameter':0.07,\ - 'first_order_term': 1.3,\ - 'second_order_term': 1, \ - 'number_of_iterations': 550 - } -out = TGV_PD(pars['input'], - pars['regularization_parameter'], - pars['first_order_term'] , - pars['second_order_term'] , - pars['number_of_iterations']) -tgv = out[0] -rms = rmse(Im, out[0]) -pars['rmse'] = rms - -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,4,6) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(tgv, cmap="gray") - # ###################### ROF_TV ######################################### start_time = timeit.default_timer() pars = {'algorithm': ROF_TV , \ 'input' : u0,\ - 'regularization_parameter':0.07,\ + 'regularisation_parameter':0.07,\ 'marching_step': 0.0025,\ 'number_of_iterations': 300 } rof = ROF_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['marching_step'], 'cpu') diff --git a/Wrappers/Python/fista-recipe/bld.bat b/Wrappers/Python/fista-recipe/bld.bat deleted file mode 100644 index 32c1bc6..0000000 --- a/Wrappers/Python/fista-recipe/bld.bat +++ /dev/null @@ -1,13 +0,0 @@ -IF NOT DEFINED CIL_VERSION ( -ECHO CIL_VERSION Not Defined. -exit 1 -) - -mkdir "%SRC_DIR%\ccpifista" -xcopy /e "%RECIPE_DIR%\.." "%SRC_DIR%\ccpifista" - -cd "%SRC_DIR%\ccpifista" -::%PYTHON% setup-fista.py -q bdist_egg -:: %PYTHON% setup.py install --single-version-externally-managed --record=record.txt -%PYTHON% setup-fista.py install -if errorlevel 1 exit 1 diff --git a/Wrappers/Python/fista-recipe/build.sh b/Wrappers/Python/fista-recipe/build.sh deleted file mode 100644 index e3f3552..0000000 --- a/Wrappers/Python/fista-recipe/build.sh +++ /dev/null @@ -1,10 +0,0 @@ -if [ -z "$CIL_VERSION" ]; then - echo "Need to set CIL_VERSION" - exit 1 -fi -mkdir "$SRC_DIR/ccpifista" -cp -r "$RECIPE_DIR/.." "$SRC_DIR/ccpifista" - -cd $SRC_DIR/ccpifista - -$PYTHON setup-fista.py install diff --git a/Wrappers/Python/fista-recipe/meta.yaml b/Wrappers/Python/fista-recipe/meta.yaml deleted file mode 100644 index 4e5cba6..0000000 --- a/Wrappers/Python/fista-recipe/meta.yaml +++ /dev/null @@ -1,29 +0,0 @@ -package: - name: ccpi-fista - version: {{ environ['CIL_VERSION'] }} - - -build: - preserve_egg_dir: False - script_env: - - CIL_VERSION -# number: 0 - -requirements: - build: - - python - - numpy - - setuptools - - run: - - python - - numpy - #- astra-toolbox - - ccpi-regularizer - - - -about: - home: http://www.ccpi.ac.uk - license: Apache v.2.0 license - summary: 'CCPi Core Imaging Library (Viewer)' diff --git a/Wrappers/Python/setup-fista.py b/Wrappers/Python/setup-fista.py deleted file mode 100644 index c5c9f4d..0000000 --- a/Wrappers/Python/setup-fista.py +++ /dev/null @@ -1,27 +0,0 @@ -from distutils.core import setup -#from setuptools import setup, find_packages -import os - -cil_version=os.environ['CIL_VERSION'] -if cil_version == '': - print("Please set the environmental variable CIL_VERSION") - sys.exit(1) - -setup( - name="ccpi-fista", - version=cil_version, - packages=['ccpi','ccpi.reconstruction'], - install_requires=['numpy'], - - zip_safe = False, - - # metadata for upload to PyPI - author="Edoardo Pasca", - author_email="edo.paskino@gmail.com", - description='CCPi Core Imaging Library - FISTA Reconstructor module', - license="Apache v2.0", - keywords="tomography interative reconstruction", - url="http://www.ccpi.ac.uk", # project home page, if any - - # could also include long_description, download_url, classifiers, etc. -) diff --git a/Wrappers/Python/setup-regularizers.py.in b/Wrappers/Python/setup-regularisers.py.in index 0811372..a1c1ab6 100644 --- a/Wrappers/Python/setup-regularizers.py.in +++ b/Wrappers/Python/setup-regularisers.py.in @@ -33,9 +33,9 @@ extra_link_args = [] extra_libraries = ['cilreg'] extra_include_dirs += [os.path.join(".." , ".." , "Core"), - os.path.join(".." , ".." , "Core", "regularizers_CPU"), - os.path.join(".." , ".." , "Core", "regularizers_GPU" , "TV_FGP" ) , - os.path.join(".." , ".." , "Core", "regularizers_GPU" , "TV_ROF" ) , + os.path.join(".." , ".." , "Core", "regularisers_CPU"), + os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_FGP" ) , + os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_ROF" ) , "."] if platform.system() == 'Windows': @@ -46,11 +46,11 @@ else: setup( name='ccpi', - description='CCPi Core Imaging Library - Image Regularizers', + description='CCPi Core Imaging Library - Image regularisers', version=cil_version, cmdclass = {'build_ext': build_ext}, - ext_modules = [Extension("ccpi.filters.cpu_regularizers_cython", - sources=[os.path.join("." , "src", "cpu_regularizers.pyx" ) ], + ext_modules = [Extension("ccpi.filters.cpu_regularisers_cython", + sources=[os.path.join("." , "src", "cpu_regularisers.pyx" ) ], include_dirs=extra_include_dirs, library_dirs=extra_library_dirs, extra_compile_args=extra_compile_args, diff --git a/Wrappers/Python/src/cpu_regularizers.pyx b/Wrappers/Python/src/cpu_regularisers.pyx index f993e54..248bad1 100644 --- a/Wrappers/Python/src/cpu_regularizers.pyx +++ b/Wrappers/Python/src/cpu_regularisers.pyx @@ -25,14 +25,14 @@ cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, #****************************************************************# #********************** Total-variation ROF *********************# #****************************************************************# -def TV_ROF_CPU(inputData, regularization_parameter, iterationsNumb, marching_step_parameter): +def TV_ROF_CPU(inputData, regularisation_parameter, iterationsNumb, marching_step_parameter): if inputData.ndim == 2: - return TV_ROF_2D(inputData, regularization_parameter, iterationsNumb, marching_step_parameter) + return TV_ROF_2D(inputData, regularisation_parameter, iterationsNumb, marching_step_parameter) elif inputData.ndim == 3: - return TV_ROF_3D(inputData, regularization_parameter, iterationsNumb, marching_step_parameter) + return TV_ROF_3D(inputData, regularisation_parameter, iterationsNumb, marching_step_parameter) def TV_ROF_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterationsNumb, float marching_step_parameter): cdef long dims[2] @@ -43,13 +43,13 @@ 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_main(&inputData[0,0], &outputData[0,0], regularization_parameter, iterationsNumb, marching_step_parameter, dims[0], dims[1], 1) + TV_ROF_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, iterationsNumb, marching_step_parameter, dims[0], dims[1], 1) return outputData def TV_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, int iterationsNumb, - float regularization_parameter, + float regularisation_parameter, float marching_step_parameter): cdef long dims[3] dims[0] = inputData.shape[0] @@ -60,7 +60,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_ROF_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularization_parameter, iterationsNumb, marching_step_parameter, dims[0], dims[1], dims[2]) + TV_ROF_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, iterationsNumb, marching_step_parameter, dims[0], dims[1], dims[2]) return outputData @@ -68,14 +68,14 @@ def TV_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, #********************** Total-variation FGP *********************# #****************************************************************# #******** Total-variation Fast-Gradient-Projection (FGP)*********# -def TV_FGP_CPU(inputData, regularization_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM): +def TV_FGP_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM): if inputData.ndim == 2: - return TV_FGP_2D(inputData, regularization_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM) + return TV_FGP_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM) elif inputData.ndim == 3: - return TV_FGP_3D(inputData, regularization_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM) + return TV_FGP_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM) def TV_FGP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterationsNumb, float tolerance_param, int methodTV, @@ -90,7 +90,7 @@ def TV_FGP_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_FGP_CPU_main(&inputData[0,0], &outputData[0,0], regularization_parameter, + TV_FGP_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, iterationsNumb, tolerance_param, methodTV, @@ -101,7 +101,7 @@ def TV_FGP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, return outputData def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterationsNumb, float tolerance_param, int methodTV, @@ -116,7 +116,7 @@ def TV_FGP_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_main(&inputData[0,0,0], &outputData[0,0,0], regularization_parameter, + TV_FGP_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, iterationsNumb, tolerance_param, methodTV, diff --git a/Wrappers/Python/src/gpu_regularizers.pyx b/Wrappers/Python/src/gpu_regularisers.pyx index a44bd1d..7ebd011 100644 --- a/Wrappers/Python/src/gpu_regularizers.pyx +++ b/Wrappers/Python/src/gpu_regularisers.pyx @@ -23,23 +23,23 @@ cdef extern void TV_FGP_GPU_main(float *Input, float *Output, float lambdaPar, i # Total-variation Rudin-Osher-Fatemi (ROF) def TV_ROF_GPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, time_marching_parameter): if inputData.ndim == 2: return ROFTV2D(inputData, - regularization_parameter, + regularisation_parameter, iterations, time_marching_parameter) elif inputData.ndim == 3: return ROFTV3D(inputData, - regularization_parameter, + regularisation_parameter, iterations, time_marching_parameter) # Total-variation Fast-Gradient-Projection (FGP) def TV_FGP_GPU(inputData, - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -47,7 +47,7 @@ def TV_FGP_GPU(inputData, printM): if inputData.ndim == 2: return FGPTV2D(inputData, - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -55,7 +55,7 @@ def TV_FGP_GPU(inputData, printM) elif inputData.ndim == 3: return FGPTV3D(inputData, - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -66,7 +66,7 @@ def TV_FGP_GPU(inputData, #********************** Total-variation ROF *********************# #****************************************************************# def ROFTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterations, float time_marching_parameter): @@ -80,7 +80,7 @@ def ROFTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, # Running CUDA code here TV_ROF_GPU_main( &inputData[0,0], &outputData[0,0], - regularization_parameter, + regularisation_parameter, iterations , time_marching_parameter, dims[0], dims[1], 1); @@ -88,7 +88,7 @@ def ROFTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, return outputData def ROFTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterations, float time_marching_parameter): @@ -103,7 +103,7 @@ def ROFTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, # Running CUDA code here TV_ROF_GPU_main( &inputData[0,0,0], &outputData[0,0,0], - regularization_parameter, + regularisation_parameter, iterations , time_marching_parameter, dims[0], dims[1], dims[2]); @@ -114,7 +114,7 @@ def ROFTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, #****************************************************************# #******** Total-variation Fast-Gradient-Projection (FGP)*********# def FGPTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterations, float tolerance_param, int methodTV, @@ -130,7 +130,7 @@ def FGPTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, # Running CUDA code here TV_FGP_GPU_main(&inputData[0,0], &outputData[0,0], - regularization_parameter, + regularisation_parameter, iterations, tolerance_param, methodTV, @@ -141,7 +141,7 @@ def FGPTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData, return outputData def FGPTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, - float regularization_parameter, + float regularisation_parameter, int iterations, float tolerance_param, int methodTV, @@ -159,7 +159,7 @@ def FGPTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData, # Running CUDA code here TV_FGP_GPU_main( &inputData[0,0,0], &outputData[0,0,0], - regularization_parameter , + regularisation_parameter , iterations, tolerance_param, methodTV, diff --git a/Wrappers/Python/test/test_cpu_regularizers.py b/Wrappers/Python/test/test_cpu_regularisers.py index 9713baa..9713baa 100644 --- a/Wrappers/Python/test/test_cpu_regularizers.py +++ b/Wrappers/Python/test/test_cpu_regularisers.py diff --git a/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py b/Wrappers/Python/test/test_cpu_vs_gpu_regularisers.py index 63be1a0..15e9042 100644 --- a/Wrappers/Python/test/test_cpu_vs_gpu_regularizers.py +++ b/Wrappers/Python/test/test_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.regularizers import ROF_TV, FGP_TV +from ccpi.filters.regularisers import ROF_TV, FGP_TV ############################################################################### def printParametersToString(pars): @@ -54,7 +54,7 @@ print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot fig = plt.figure(1) -plt.suptitle('Comparison of ROF-TV regularizer using CPU and GPU implementations') +plt.suptitle('Comparison of ROF-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") @@ -62,14 +62,14 @@ imgplot = plt.imshow(u0,cmap="gray") # set parameters pars = {'algorithm': ROF_TV, \ 'input' : u0,\ - 'regularization_parameter':0.04,\ + 'regularisation_parameter':0.04,\ 'number_of_iterations': 1200,\ 'time_marching_parameter': 0.0025 } print ("#############ROF TV CPU####################") start_time = timeit.default_timer() rof_cpu = ROF_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['time_marching_parameter'],'cpu') rms = rmse(Im, rof_cpu) @@ -92,7 +92,7 @@ plt.title('{}'.format('CPU results')) print ("##############ROF TV GPU##################") start_time = timeit.default_timer() rof_gpu = ROF_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['time_marching_parameter'],'gpu') @@ -132,7 +132,7 @@ print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%") ## plot fig = plt.figure(2) -plt.suptitle('Comparison of FGP-TV regularizer using CPU and GPU implementations') +plt.suptitle('Comparison of FGP-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") @@ -140,7 +140,7 @@ imgplot = plt.imshow(u0,cmap="gray") # set parameters pars = {'algorithm' : FGP_TV, \ 'input' : u0,\ - 'regularization_parameter':0.04, \ + 'regularisation_parameter':0.04, \ 'number_of_iterations' :1200 ,\ 'tolerance_constant':0.00001,\ 'methodTV': 0 ,\ @@ -151,7 +151,7 @@ pars = {'algorithm' : FGP_TV, \ print ("#############FGP TV CPU####################") start_time = timeit.default_timer() fgp_cpu = FGP_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['tolerance_constant'], pars['methodTV'], @@ -179,7 +179,7 @@ plt.title('{}'.format('CPU results')) print ("##############FGP TV GPU##################") start_time = timeit.default_timer() fgp_gpu = FGP_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['tolerance_constant'], pars['methodTV'], diff --git a/Wrappers/Python/test/test_gpu_regularizers.py b/Wrappers/Python/test/test_gpu_regularisers.py index 640b3f9..2103c0e 100644 --- a/Wrappers/Python/test/test_gpu_regularizers.py +++ b/Wrappers/Python/test/test_gpu_regularisers.py @@ -11,8 +11,7 @@ import numpy as np import os from enum import Enum import timeit -from ccpi.filters.gpu_regularizers import Diff4thHajiaboli, NML -from ccpi.filters.regularizers import ROF_TV, FGP_TV +from ccpi.filters.regularisers import ROF_TV, FGP_TV ############################################################################### def printParametersToString(pars): txt = r'' @@ -32,9 +31,6 @@ def rmse(im1, im2): return rmse filename = os.path.join(".." , ".." , ".." , "data" ,"lena_gray_512.tif") -#filename = r"C:\Users\ofn77899\Documents\GitHub\CCPi-FISTA_reconstruction\data\lena_gray_512.tif" -#filename = r"/home/ofn77899/Reconstruction/CCPi-FISTA_Reconstruction/data/lena_gray_512.tif" -#filename = r'/home/algol/Documents/Python/STD_test_images/lena_gray_512.tif' Im = plt.imread(filename) Im = np.asarray(Im, dtype='float32') @@ -56,112 +52,20 @@ a.set_title('noise') imgplot = plt.imshow(u0,cmap="gray") -## Diff4thHajiaboli -start_time = timeit.default_timer() -pars = { -'algorithm' : Diff4thHajiaboli , \ - 'input' : u0, - 'edge_preserv_parameter':0.1 , \ -'number_of_iterations' :250 ,\ -'time_marching_parameter':0.03 ,\ -'regularization_parameter':0.7 -} - - -d4h = Diff4thHajiaboli(pars['input'], - pars['edge_preserv_parameter'], - pars['number_of_iterations'], - pars['time_marching_parameter'], - pars['regularization_parameter']) -rms = rmse(Im, d4h) -pars['rmse'] = rms -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,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=12, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(d4h, cmap="gray") - -a=fig.add_subplot(2,4,6) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, 'd4h - u0', transform=a.transAxes, fontsize=12, - verticalalignment='top', bbox=props) -imgplot = plt.imshow((d4h - u0)**2, vmin=0, vmax=0.03, cmap="gray") -plt.colorbar(ticks=[0, 0.03], orientation='vertical') - - -## Patch Based Regul NML -start_time = timeit.default_timer() -""" -pars = {'algorithm' : NML , \ - 'input' : u0, - 'SearchW_real':3 , \ -'SimilW' :1,\ -'h':0.05 ,# -'lambda' : 0.08 -} -""" -pars = { -'algorithm' : NML , \ - 'input' : u0, - 'regularization_parameter': 0.01,\ - 'searching_window_ratio':3, \ - 'similarity_window_ratio':1,\ - 'PB_filtering_parameter': 0.2 -} - -nml = NML(pars['input'], - pars['searching_window_ratio'], - pars['similarity_window_ratio'], - pars['PB_filtering_parameter'], - pars['regularization_parameter']) -rms = rmse(Im, nml) -pars['rmse'] = rms -txtstr = printParametersToString(pars) -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,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=12, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(nml, cmap="gray") - -a=fig.add_subplot(2,4,7) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, 'nml - u0', transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow((nml - u0)**2, vmin=0, vmax=0.03, cmap="gray") -plt.colorbar(ticks=[0, 0.03], orientation='vertical') - - -## Rudin-Osher-Fatemi (ROF) TV regularization +## Rudin-Osher-Fatemi (ROF) TV regularisation start_time = timeit.default_timer() pars = { 'algorithm' : ROF_TV , \ 'input' : u0, - 'regularization_parameter': 0.04,\ + 'regularisation_parameter': 0.04,\ 'number_of_iterations':300,\ 'time_marching_parameter': 0.0025 } rof_tv = TV_ROF_GPU(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['time_marching_parameter'],'gpu') @@ -190,13 +94,13 @@ imgplot = plt.imshow((rof_tv - u0)**2, vmin=0, vmax=0.03, cmap="gray") plt.colorbar(ticks=[0, 0.03], orientation='vertical') plt.show() -## Fast-Gradient Projection TV regularization +## Fast-Gradient Projection TV regularisation """ start_time = timeit.default_timer() pars = {'algorithm' : FGP_TV, \ 'input' : u0,\ - 'regularization_parameter':0.04, \ + 'regularisation_parameter':0.04, \ 'number_of_iterations' :1200 ,\ 'tolerance_constant':0.00001,\ 'methodTV': 0 ,\ @@ -205,7 +109,7 @@ pars = {'algorithm' : FGP_TV, \ } fgp_gpu = FGP_TV(pars['input'], - pars['regularization_parameter'], + pars['regularisation_parameter'], pars['number_of_iterations'], pars['tolerance_constant'], pars['methodTV'], diff --git a/Wrappers/Python/test/test_regularizers_3d.py b/Wrappers/Python/test/test_regularisers_3d.py index 2d11a7e..2d11a7e 100644 --- a/Wrappers/Python/test/test_regularizers_3d.py +++ b/Wrappers/Python/test/test_regularisers_3d.py diff --git a/Wrappers/Python/test/test_regularizers.py b/Wrappers/Python/test/test_regularizers.py deleted file mode 100644 index cf5da2b..0000000 --- a/Wrappers/Python/test/test_regularizers.py +++ /dev/null @@ -1,395 +0,0 @@ -# -*- coding: utf-8 -*- -""" -Created on Fri Aug 4 11:10:05 2017 - -@author: ofn77899 -""" - -#from ccpi.viewer.CILViewer2D import Converter -#import vtk - -import matplotlib.pyplot as plt -import numpy as np -import os -from enum import Enum -import timeit -#from PIL import Image -#from Regularizer import Regularizer -from ccpi.filters.Regularizer import Regularizer - -############################################################################### -#https://stackoverflow.com/questions/13875989/comparing-image-in-url-to-image-in-filesystem-in-python/13884956#13884956 -#NRMSE a normalization of the root of the mean squared error -#NRMSE is simply 1 - [RMSE / (maxval - minval)]. Where maxval is the maximum -# intensity from the two images being compared, and respectively the same for -# minval. RMSE is given by the square root of MSE: -# sqrt[(sum(A - B) ** 2) / |A|], -# where |A| means the number of elements in A. By doing this, the maximum value -# given by RMSE is maxval. - -def nrmse(im1, im2): - a, b = im1.shape - rmse = np.sqrt(np.sum((im2 - im1) ** 2) / float(a * b)) - max_val = max(np.max(im1), np.max(im2)) - min_val = min(np.min(im1), np.min(im2)) - return 1 - (rmse / (max_val - min_val)) -############################################################################### - -############################################################################### -# -# 2D Regularizers -# -############################################################################### -#Example: -# figure; -# Im = double(imread('lena_gray_256.tif'))/255; % loading image -# u0 = Im + .05*randn(size(Im)); u0(u0 < 0) = 0; -# u = SplitBregman_TV(single(u0), 10, 30, 1e-04); - - -filename = r"C:\Users\ofn77899\Documents\GitHub\CCPi-FISTA_reconstruction\data\lena_gray_512.tif" -#filename = r"/home/ofn77899/Reconstruction/CCPi-FISTA_Reconstruction/data/lena_gray_512.tif" -#filename = r'/home/algol/Documents/Python/STD_test_images/lena_gray_512.tif' - -#reader = vtk.vtkTIFFReader() -#reader.SetFileName(os.path.normpath(filename)) -#reader.Update() -Im = plt.imread(filename) -#Im = Image.open('/home/algol/Documents/Python/STD_test_images/lena_gray_512.tif')/255 -#img.show() -Im = np.asarray(Im, dtype='float32') - - - - -#imgplot = plt.imshow(Im) -perc = 0.05 -u0 = Im + (perc* np.random.normal(size=np.shape(Im))) -# map the u0 u0->u0>0 -f = np.frompyfunc(lambda x: 0 if x < 0 else x, 1,1) -u0 = f(u0).astype('float32') - -## plot -fig = plt.figure() -#a=fig.add_subplot(3,3,1) -#a.set_title('Original') -#imgplot = plt.imshow(Im) - -a=fig.add_subplot(2,3,1) -a.set_title('noise') -imgplot = plt.imshow(u0,cmap="gray") - -reg_output = [] -############################################################################## -# Call regularizer - -####################### SplitBregman_TV ##################################### -# u = SplitBregman_TV(single(u0), 10, 30, 1e-04); - -start_time = timeit.default_timer() -reg = Regularizer(Regularizer.Algorithm.SplitBregman_TV) -print (reg.pars) -reg.setParameter(input=u0) -reg.setParameter(regularization_parameter=10.) -# or -# reg.setParameter(input=u0, regularization_parameter=10., #number_of_iterations=30, - #tolerance_constant=1e-4, - #TV_Penalty=Regularizer.TotalVariationPenalty.l1) -plotme = reg(output_all=True) [0] -pars = reg.pars -txtstr = reg.printParametersToString() -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) - - -a=fig.add_subplot(2,3,2) - - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(plotme,cmap="gray") - -###################### FGP_TV ######################################### -# u = FGP_TV(single(u0), 0.05, 100, 1e-04); -start_time = timeit.default_timer() -reg = Regularizer(Regularizer.Algorithm.FGP_TV) -out2 = reg(input=u0, regularization_parameter=5e-4, - number_of_iterations=10) -txtstr = reg.printParametersToString() -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) - - -a=fig.add_subplot(2,3,3) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -imgplot = plt.imshow(out2,cmap="gray") -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) - -###################### LLT_model ######################################### -# * u0 = Im + .03*randn(size(Im)); % adding noise -# [Den] = LLT_model(single(u0), 10, 0.1, 1); -#Den = LLT_model(single(u0), 25, 0.0003, 300, 0.0001, 0); -#input, regularization_parameter , time_step, number_of_iterations, -# tolerance_constant, restrictive_Z_smoothing=0 - -del out2 -start_time = timeit.default_timer() -reg = Regularizer(Regularizer.Algorithm.LLT_model) -out2 = reg(input=u0, regularization_parameter=25, - time_step=0.0003, - tolerance_constant=0.001, - number_of_iterations=300) -txtstr = reg.printParametersToString() -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,3,4) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(out2,cmap="gray") - - -# ###################### PatchBased_Regul ######################################### -# # Quick 2D denoising example in Matlab: -# # Im = double(imread('lena_gray_256.tif'))/255; % loading image -# # u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -# # ImDen = PB_Regul_CPU(single(u0), 3, 1, 0.08, 0.05); -start_time = timeit.default_timer() -reg = Regularizer(Regularizer.Algorithm.PatchBased_Regul) -out2 = reg(input=u0, regularization_parameter=0.05, - searching_window_ratio=3, - similarity_window_ratio=1, - PB_filtering_parameter=0.08) -txtstr = reg.printParametersToString() -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) - -a=fig.add_subplot(2,3,5) - - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(out2,cmap="gray") - - -# ###################### TGV_PD ######################################### -# # Quick 2D denoising example in Matlab: -# # Im = double(imread('lena_gray_256.tif'))/255; % loading image -# # u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -# # u = PrimalDual_TGV(single(u0), 0.02, 1.3, 1, 550); - -start_time = timeit.default_timer() -reg = Regularizer(Regularizer.Algorithm.TGV_PD) -out2 = reg(input=u0, regularization_parameter=0.05, - first_order_term=1.3, - second_order_term=1, - number_of_iterations=550) -txtstr = reg.printParametersToString() -txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time) -print (txtstr) -a=fig.add_subplot(2,3,6) - -# these are matplotlib.patch.Patch properties -props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -# place a text box in upper left in axes coords -a.text(0.05, 0.95, txtstr, transform=a.transAxes, fontsize=14, - verticalalignment='top', bbox=props) -imgplot = plt.imshow(out2,cmap="gray") - - -plt.show() - -################################################################################ -## -## 3D Regularizers -## -################################################################################ -##Example: -## figure; -## Im = double(imread('lena_gray_256.tif'))/255; % loading image -## u0 = Im + .05*randn(size(Im)); u0(u0 < 0) = 0; -## u = SplitBregman_TV(single(u0), 10, 30, 1e-04); -# -##filename = r"C:\Users\ofn77899\Documents\GitHub\CCPi-Reconstruction\python\test\reconstruction_example.mha" -#filename = r"C:\Users\ofn77899\Documents\GitHub\CCPi-Simpleflex\data\head.mha" -# -#reader = vtk.vtkMetaImageReader() -#reader.SetFileName(os.path.normpath(filename)) -#reader.Update() -##vtk returns 3D images, let's take just the one slice there is as 2D -#Im = Converter.vtk2numpy(reader.GetOutput()) -#Im = Im.astype('float32') -##imgplot = plt.imshow(Im) -#perc = 0.05 -#u0 = Im + (perc* np.random.normal(size=np.shape(Im))) -## map the u0 u0->u0>0 -#f = np.frompyfunc(lambda x: 0 if x < 0 else x, 1,1) -#u0 = f(u0).astype('float32') -#converter = Converter.numpy2vtkImporter(u0, reader.GetOutput().GetSpacing(), -# reader.GetOutput().GetOrigin()) -#converter.Update() -#writer = vtk.vtkMetaImageWriter() -#writer.SetInputData(converter.GetOutput()) -#writer.SetFileName(r"C:\Users\ofn77899\Documents\GitHub\CCPi-FISTA_reconstruction\data\noisy_head.mha") -##writer.Write() -# -# -### plot -#fig3D = plt.figure() -##a=fig.add_subplot(3,3,1) -##a.set_title('Original') -##imgplot = plt.imshow(Im) -#sliceNo = 32 -# -#a=fig3D.add_subplot(2,3,1) -#a.set_title('noise') -#imgplot = plt.imshow(u0.T[sliceNo]) -# -#reg_output3d = [] -# -############################################################################### -## Call regularizer -# -######################## SplitBregman_TV ##################################### -## u = SplitBregman_TV(single(u0), 10, 30, 1e-04); -# -##reg = Regularizer(Regularizer.Algorithm.SplitBregman_TV) -# -##out = reg(input=u0, regularization_parameter=10., #number_of_iterations=30, -## #tolerance_constant=1e-4, -## TV_Penalty=Regularizer.TotalVariationPenalty.l1) -# -#out2 = Regularizer.SplitBregman_TV(input=u0, regularization_parameter=10., number_of_iterations=30, -# tolerance_constant=1e-4, -# TV_Penalty=Regularizer.TotalVariationPenalty.l1) -# -# -#pars = out2[-2] -#reg_output3d.append(out2) -# -#a=fig3D.add_subplot(2,3,2) -# -# -#textstr = out2[-1] -# -# -## these are matplotlib.patch.Patch properties -#props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -## place a text box in upper left in axes coords -#a.text(0.05, 0.95, textstr, transform=a.transAxes, fontsize=14, -# verticalalignment='top', bbox=props) -#imgplot = plt.imshow(reg_output3d[-1][0].T[sliceNo]) -# -####################### FGP_TV ######################################### -## u = FGP_TV(single(u0), 0.05, 100, 1e-04); -#out2 = Regularizer.FGP_TV(input=u0, regularization_parameter=0.005, -# number_of_iterations=200) -#pars = out2[-2] -#reg_output3d.append(out2) -# -#a=fig3D.add_subplot(2,3,2) -# -# -#textstr = out2[-1] -# -# -## these are matplotlib.patch.Patch properties -#props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -## place a text box in upper left in axes coords -#a.text(0.05, 0.95, textstr, transform=a.transAxes, fontsize=14, -# verticalalignment='top', bbox=props) -#imgplot = plt.imshow(reg_output3d[-1][0].T[sliceNo]) -# -####################### LLT_model ######################################### -## * u0 = Im + .03*randn(size(Im)); % adding noise -## [Den] = LLT_model(single(u0), 10, 0.1, 1); -##Den = LLT_model(single(u0), 25, 0.0003, 300, 0.0001, 0); -##input, regularization_parameter , time_step, number_of_iterations, -## tolerance_constant, restrictive_Z_smoothing=0 -#out2 = Regularizer.LLT_model(input=u0, regularization_parameter=25, -# time_step=0.0003, -# tolerance_constant=0.0001, -# number_of_iterations=300) -#pars = out2[-2] -#reg_output3d.append(out2) -# -#a=fig3D.add_subplot(2,3,2) -# -# -#textstr = out2[-1] -# -# -## these are matplotlib.patch.Patch properties -#props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -## place a text box in upper left in axes coords -#a.text(0.05, 0.95, textstr, transform=a.transAxes, fontsize=14, -# verticalalignment='top', bbox=props) -#imgplot = plt.imshow(reg_output3d[-1][0].T[sliceNo]) -# -####################### PatchBased_Regul ######################################### -## Quick 2D denoising example in Matlab: -## Im = double(imread('lena_gray_256.tif'))/255; % loading image -## u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -## ImDen = PB_Regul_CPU(single(u0), 3, 1, 0.08, 0.05); -# -#out2 = Regularizer.PatchBased_Regul(input=u0, regularization_parameter=0.05, -# searching_window_ratio=3, -# similarity_window_ratio=1, -# PB_filtering_parameter=0.08) -#pars = out2[-2] -#reg_output3d.append(out2) -# -#a=fig3D.add_subplot(2,3,2) -# -# -#textstr = out2[-1] -# -# -## these are matplotlib.patch.Patch properties -#props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -## place a text box in upper left in axes coords -#a.text(0.05, 0.95, textstr, transform=a.transAxes, fontsize=14, -# verticalalignment='top', bbox=props) -#imgplot = plt.imshow(reg_output3d[-1][0].T[sliceNo]) -# - -###################### TGV_PD ######################################### -# Quick 2D denoising example in Matlab: -# Im = double(imread('lena_gray_256.tif'))/255; % loading image -# u0 = Im + .03*randn(size(Im)); u0(u0<0) = 0; % adding noise -# u = PrimalDual_TGV(single(u0), 0.02, 1.3, 1, 550); - - -#out2 = Regularizer.TGV_PD(input=u0, regularization_parameter=0.05, -# first_order_term=1.3, -# second_order_term=1, -# number_of_iterations=550) -#pars = out2[-2] -#reg_output3d.append(out2) -# -#a=fig3D.add_subplot(2,3,2) -# -# -#textstr = out2[-1] -# -# -## these are matplotlib.patch.Patch properties -#props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) -## place a text box in upper left in axes coords -#a.text(0.05, 0.95, textstr, transform=a.transAxes, fontsize=14, -# verticalalignment='top', bbox=props) -#imgplot = plt.imshow(reg_output3d[-1][0].T[sliceNo]) |