summaryrefslogtreecommitdiffstats
path: root/src/Python
diff options
context:
space:
mode:
Diffstat (limited to 'src/Python')
-rw-r--r--src/Python/CMakeLists.txt141
-rw-r--r--src/Python/ccpi/__init__.py0
-rw-r--r--src/Python/ccpi/filters/__init__.py0
-rw-r--r--src/Python/ccpi/filters/regularisers.py214
-rw-r--r--src/Python/setup-regularisers.py.in75
-rw-r--r--src/Python/src/cpu_regularisers.pyx685
-rw-r--r--src/Python/src/gpu_regularisers.pyx640
7 files changed, 1755 insertions, 0 deletions
diff --git a/src/Python/CMakeLists.txt b/src/Python/CMakeLists.txt
new file mode 100644
index 0000000..c2ef855
--- /dev/null
+++ b/src/Python/CMakeLists.txt
@@ -0,0 +1,141 @@
+# Copyright 2018 Edoardo Pasca
+cmake_minimum_required (VERSION 3.0)
+
+project(regulariserPython)
+#https://stackoverflow.com/questions/13298504/using-cmake-with-setup-py
+
+# The version number.
+
+#set (CIL_VERSION $ENV{CIL_VERSION} CACHE INTERNAL "Core Imaging Library version" FORCE)
+
+# conda orchestrated build
+message("CIL_VERSION: ${CIL_VERSION}")
+#include (GenerateExportHeader)
+
+find_package(PythonInterp REQUIRED)
+if (PYTHONINTERP_FOUND)
+ message ("Current Python " ${PYTHON_VERSION_STRING} " found " ${PYTHON_EXECUTABLE})
+endif()
+
+
+## Build the regularisers package as a library
+message("Creating Regularisers as shared library")
+
+message("CMAKE_CXX_FLAGS ${CMAKE_CXX_FLAGS}")
+
+set(CMAKE_BUILD_TYPE "Release")
+
+if(WIN32)
+ set (FLAGS "/DWIN32 /EHsc /openmp /DCCPiCore_EXPORTS")
+ set (CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} /NODEFAULTLIB:MSVCRT.lib")
+
+ set (EXTRA_LIBRARIES)
+
+ message("library lib: ${LIBRARY_LIB}")
+
+elseif(UNIX)
+ set (FLAGS "-fopenmp -O2 -funsigned-char -Wall -Wl,--no-undefined -DCCPiReconstructionIterative_EXPORTS -std=c++0x")
+ set (EXTRA_LIBRARIES
+ "gomp"
+ )
+endif()
+
+# GPU regularisers
+if (BUILD_CUDA)
+ find_package(CUDA)
+ if (CUDA_FOUND)
+ message("CUDA FOUND")
+ set (SETUP_GPU_WRAPPERS "extra_libraries += ['cilregcuda']\n\
+setup( \n\
+ name='ccpi', \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_regularisers',\n\
+ sources=[ \n\
+ os.path.join('.' , 'src', 'gpu_regularisers.pyx' ),\n\
+ ],\n\
+ include_dirs=extra_include_dirs, \n\
+ library_dirs=extra_library_dirs, \n\
+ extra_compile_args=extra_compile_args, \n\
+ libraries=extra_libraries ), \n\
+ ],\n\
+ zip_safe = False, \n\
+ packages = {'ccpi','ccpi.filters'},\n\
+ )")
+ else()
+ message("CUDA NOT FOUND")
+ set(SETUP_GPU_WRAPPERS "#CUDA NOT FOUND")
+ endif()
+endif()
+configure_file("${CMAKE_CURRENT_SOURCE_DIR}/setup-regularisers.py.in" "${CMAKE_CURRENT_BINARY_DIR}/setup-regularisers.py")
+
+
+find_package(PythonInterp)
+find_package(PythonLibs)
+if (PYTHONINTERP_FOUND)
+ message(STATUS "Found PYTHON_EXECUTABLE=${PYTHON_EXECUTABLE}")
+ message(STATUS "Python version ${PYTHON_VERSION_STRING}")
+endif()
+if (PYTHONLIBS_FOUND)
+ message(STATUS "Found PYTHON_INCLUDE_DIRS=${PYTHON_INCLUDE_DIRS}")
+ message(STATUS "Found PYTHON_LIBRARIES=${PYTHON_LIBRARIES}")
+endif()
+
+if (PYTHONINTERP_FOUND)
+ message("Python found " ${PYTHON_EXECUTABLE})
+ set(SETUP_PY_IN "${CMAKE_CURRENT_SOURCE_DIR}/setup-regularisers.py.in")
+ set(SETUP_PY "${CMAKE_CURRENT_BINARY_DIR}/setup-regularisers.py")
+ #set(DEPS "${CMAKE_CURRENT_SOURCE_DIR}/module/__init__.py")
+ set (DEPS "${CMAKE_BINARY_DIR}/Core/")
+ set(OUTPUT "${CMAKE_CURRENT_BINARY_DIR}/build/timestamp")
+
+ configure_file(${SETUP_PY_IN} ${SETUP_PY})
+
+ message("Core binary dir " ${CMAKE_BINARY_DIR}/Core/${CMAKE_BUILD_TYPE})
+
+ if (CONDA_BUILD)
+ add_custom_command(OUTPUT ${OUTPUT}
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/src ${CMAKE_CURRENT_BINARY_DIR}/src
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/ccpi ${CMAKE_CURRENT_BINARY_DIR}/ccpi
+ COMMAND ${CMAKE_COMMAND} -E env CIL_VERSION=${CIL_VERSION}
+ PREFIX=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_INC=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_LIB=${CMAKE_BINARY_DIR}/Core
+ ${PYTHON_EXECUTABLE} ${SETUP_PY} install
+ COMMAND ${CMAKE_COMMAND} -E touch ${OUTPUT}
+ DEPENDS cilreg)
+
+ else()
+ if (WIN32)
+ add_custom_command(OUTPUT ${OUTPUT}
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/src ${CMAKE_CURRENT_BINARY_DIR}/src
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/ccpi ${CMAKE_CURRENT_BINARY_DIR}/ccpi
+ COMMAND ${CMAKE_COMMAND} -E env CIL_VERSION=${CIL_VERSION}
+ PREFIX=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_INC=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_LIB=${CMAKE_BINARY_DIR}/Core/${CMAKE_BUILD_TYPE}
+ ${PYTHON_EXECUTABLE} ${SETUP_PY} build_ext --inplace
+ COMMAND ${CMAKE_COMMAND} -E touch ${OUTPUT}
+ DEPENDS cilreg)
+ else()
+ add_custom_command(OUTPUT ${OUTPUT}
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/src ${CMAKE_CURRENT_BINARY_DIR}/src
+ COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/ccpi ${CMAKE_CURRENT_BINARY_DIR}/ccpi
+ COMMAND ${CMAKE_COMMAND} -E env CIL_VERSION=${CIL_VERSION}
+ PREFIX=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_INC=${CMAKE_SOURCE_DIR}/Core
+ LIBRARY_LIB=${CMAKE_BINARY_DIR}/Core
+ ${PYTHON_EXECUTABLE} ${SETUP_PY} build_ext --inplace
+ COMMAND ${CMAKE_COMMAND} -E touch ${OUTPUT}
+ DEPENDS cilreg)
+ endif()
+ install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/ccpi
+ DESTINATION ${PYTHON_DEST})
+ endif()
+
+
+ add_custom_target(PythonWrapper ALL DEPENDS ${OUTPUT})
+
+ #install(CODE "execute_process(COMMAND ${PYTHON} ${SETUP_PY} install)")
+endif()
diff --git a/src/Python/ccpi/__init__.py b/src/Python/ccpi/__init__.py
new file mode 100644
index 0000000..e69de29
--- /dev/null
+++ b/src/Python/ccpi/__init__.py
diff --git a/src/Python/ccpi/filters/__init__.py b/src/Python/ccpi/filters/__init__.py
new file mode 100644
index 0000000..e69de29
--- /dev/null
+++ b/src/Python/ccpi/filters/__init__.py
diff --git a/src/Python/ccpi/filters/regularisers.py b/src/Python/ccpi/filters/regularisers.py
new file mode 100644
index 0000000..588ea32
--- /dev/null
+++ b/src/Python/ccpi/filters/regularisers.py
@@ -0,0 +1,214 @@
+"""
+script which assigns a proper device core function based on a flag ('cpu' or 'gpu')
+"""
+
+from ccpi.filters.cpu_regularisers import TV_ROF_CPU, TV_FGP_CPU, TV_SB_CPU, dTV_FGP_CPU, TNV_CPU, NDF_CPU, Diff4th_CPU, TGV_CPU, LLT_ROF_CPU, PATCHSEL_CPU, NLTV_CPU
+try:
+ from ccpi.filters.gpu_regularisers import TV_ROF_GPU, TV_FGP_GPU, TV_SB_GPU, dTV_FGP_GPU, NDF_GPU, Diff4th_GPU, TGV_GPU, LLT_ROF_GPU, PATCHSEL_GPU
+ gpu_enabled = True
+except ImportError:
+ gpu_enabled = False
+from ccpi.filters.cpu_regularisers import NDF_INPAINT_CPU, NVM_INPAINT_CPU
+
+def ROF_TV(inputData, regularisation_parameter, iterations,
+ time_marching_parameter,device='cpu'):
+ if device == 'cpu':
+ return TV_ROF_CPU(inputData,
+ regularisation_parameter,
+ iterations,
+ time_marching_parameter)
+ elif device == 'gpu' and gpu_enabled:
+ return TV_ROF_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ time_marching_parameter)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+
+def FGP_TV(inputData, regularisation_parameter,iterations,
+ tolerance_param, methodTV, nonneg, printM, device='cpu'):
+ if device == 'cpu':
+ return TV_FGP_CPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM)
+ elif device == 'gpu' and gpu_enabled:
+ return TV_FGP_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def SB_TV(inputData, regularisation_parameter, iterations,
+ tolerance_param, methodTV, printM, device='cpu'):
+ if device == 'cpu':
+ return TV_SB_CPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM)
+ elif device == 'gpu' and gpu_enabled:
+ return TV_SB_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def FGP_dTV(inputData, refdata, regularisation_parameter, iterations,
+ tolerance_param, eta_const, methodTV, nonneg, printM, device='cpu'):
+ if device == 'cpu':
+ return dTV_FGP_CPU(inputData,
+ refdata,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM)
+ elif device == 'gpu' and gpu_enabled:
+ return dTV_FGP_GPU(inputData,
+ refdata,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def TNV(inputData, regularisation_parameter, iterations, tolerance_param):
+ return TNV_CPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param)
+def NDF(inputData, regularisation_parameter, edge_parameter, iterations,
+ time_marching_parameter, penalty_type, device='cpu'):
+ if device == 'cpu':
+ return NDF_CPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter,
+ penalty_type)
+ elif device == 'gpu' and gpu_enabled:
+ return NDF_GPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter,
+ penalty_type)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def Diff4th(inputData, regularisation_parameter, edge_parameter, iterations,
+ time_marching_parameter, device='cpu'):
+ if device == 'cpu':
+ return Diff4th_CPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter)
+ elif device == 'gpu' and gpu_enabled:
+ return Diff4th_GPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+
+def PatchSelect(inputData, searchwindow, patchwindow, neighbours, edge_parameter, device='cpu'):
+ if device == 'cpu':
+ return PATCHSEL_CPU(inputData,
+ searchwindow,
+ patchwindow,
+ neighbours,
+ edge_parameter)
+ elif device == 'gpu' and gpu_enabled:
+ return PATCHSEL_GPU(inputData,
+ searchwindow,
+ patchwindow,
+ neighbours,
+ edge_parameter)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+
+def NLTV(inputData, H_i, H_j, H_k, Weights, regularisation_parameter, iterations):
+ return NLTV_CPU(inputData,
+ H_i,
+ H_j,
+ H_k,
+ Weights,
+ regularisation_parameter,
+ iterations)
+
+def TGV(inputData, regularisation_parameter, alpha1, alpha0, iterations,
+ LipshitzConst, device='cpu'):
+ if device == 'cpu':
+ return TGV_CPU(inputData,
+ regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterations,
+ LipshitzConst)
+ elif device == 'gpu' and gpu_enabled:
+ return TGV_GPU(inputData,
+ regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterations,
+ LipshitzConst)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def LLT_ROF(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations,
+ time_marching_parameter, device='cpu'):
+ if device == 'cpu':
+ return LLT_ROF_CPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+ elif device == 'gpu' and gpu_enabled:
+ return LLT_ROF_GPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+ else:
+ if not gpu_enabled and device == 'gpu':
+ raise ValueError ('GPU is not available')
+ raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
+ .format(device))
+def NDF_INP(inputData, maskData, regularisation_parameter, edge_parameter, iterations,
+ time_marching_parameter, penalty_type):
+ return NDF_INPAINT_CPU(inputData, maskData, regularisation_parameter,
+ edge_parameter, iterations, time_marching_parameter, penalty_type)
+
+def NVM_INP(inputData, maskData, SW_increment, iterations):
+ return NVM_INPAINT_CPU(inputData, maskData, SW_increment, iterations)
diff --git a/src/Python/setup-regularisers.py.in b/src/Python/setup-regularisers.py.in
new file mode 100644
index 0000000..462edda
--- /dev/null
+++ b/src/Python/setup-regularisers.py.in
@@ -0,0 +1,75 @@
+#!/usr/bin/env python
+
+import setuptools
+from distutils.core import setup
+from distutils.extension import Extension
+from Cython.Distutils import build_ext
+
+import os
+import sys
+import numpy
+import platform
+
+cil_version=os.environ['CIL_VERSION']
+if cil_version == '':
+ print("Please set the environmental variable CIL_VERSION")
+ sys.exit(1)
+
+library_include_path = ""
+library_lib_path = ""
+try:
+ library_include_path = os.environ['LIBRARY_INC']
+ library_lib_path = os.environ['LIBRARY_LIB']
+except:
+ library_include_path = os.environ['PREFIX']+'/include'
+ pass
+
+extra_include_dirs = [numpy.get_include(), library_include_path]
+#extra_library_dirs = [os.path.join(library_include_path, "..", "lib")]
+extra_compile_args = []
+extra_library_dirs = [library_lib_path]
+extra_compile_args = []
+extra_link_args = []
+extra_libraries = ['cilreg']
+
+print ("extra_library_dirs " , extra_library_dirs)
+
+extra_include_dirs += [os.path.join(".." , ".." , "Core"),
+ os.path.join(".." , ".." , "Core", "regularisers_CPU"),
+ os.path.join(".." , ".." , "Core", "inpainters_CPU"),
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_FGP" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_ROF" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TV_SB" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "TGV" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "LLTROF" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "NDF" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "dTV_FGP" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "DIFF4th" ) ,
+ os.path.join(".." , ".." , "Core", "regularisers_GPU" , "PatchSelect" ) ,
+ "."]
+
+if platform.system() == 'Windows':
+ extra_compile_args[0:] = ['/DWIN32','/EHsc','/DBOOST_ALL_NO_LIB' , '/openmp' ]
+else:
+ extra_compile_args = ['-fopenmp','-O2', '-funsigned-char', '-Wall', '-std=c++0x']
+ extra_libraries += [@EXTRA_OMP_LIB@]
+
+setup(
+ name='ccpi',
+ description='CCPi Core Imaging Library - Image regularisers',
+ version=cil_version,
+ cmdclass = {'build_ext': build_ext},
+ ext_modules = [Extension("ccpi.filters.cpu_regularisers",
+ sources=[os.path.join("." , "src", "cpu_regularisers.pyx" ) ],
+ include_dirs=extra_include_dirs,
+ library_dirs=extra_library_dirs,
+ extra_compile_args=extra_compile_args,
+ libraries=extra_libraries ),
+
+ ],
+ zip_safe = False,
+ packages = {'ccpi','ccpi.filters'},
+)
+
+
+@SETUP_GPU_WRAPPERS@
diff --git a/src/Python/src/cpu_regularisers.pyx b/src/Python/src/cpu_regularisers.pyx
new file mode 100644
index 0000000..11a0617
--- /dev/null
+++ b/src/Python/src/cpu_regularisers.pyx
@@ -0,0 +1,685 @@
+# distutils: language=c++
+"""
+Copyright 2018 CCPi
+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.
+
+Author: Edoardo Pasca, Daniil Kazantsev
+"""
+
+import cython
+import numpy as np
+cimport numpy as np
+
+cdef extern float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
+cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ);
+cdef extern float SB_TV_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int printM, int dimX, int dimY, int dimZ);
+cdef extern float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
+cdef extern float TGV_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY, int dimZ);
+cdef extern float Diffusion_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int penaltytype, int dimX, int dimY, int dimZ);
+cdef extern float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
+cdef extern float TNV_CPU_main(float *Input, float *u, float lambdaPar, int maxIter, float tol, int dimX, int dimY, int dimZ);
+cdef extern float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambdaPar, int iterationsNumb, float epsil, float eta, int methodTV, int nonneg, int printM, int dimX, int dimY, int dimZ);
+cdef extern float PatchSelect_CPU_main(float *Input, unsigned short *H_i, unsigned short *H_j, unsigned short *H_k, float *Weights, int dimX, int dimY, int dimZ, int SearchWindow, int SimilarWin, int NumNeighb, float h, int switchM);
+cdef extern float Nonlocal_TV_CPU_main(float *A_orig, float *Output, unsigned short *H_i, unsigned short *H_j, unsigned short *H_k, float *Weights, int dimX, int dimY, int dimZ, int NumNeighb, float lambdaReg, int IterNumb);
+
+cdef extern float Diffusion_Inpaint_CPU_main(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int penaltytype, int dimX, int dimY, int dimZ);
+cdef extern float NonlocalMarching_Inpaint_main(float *Input, unsigned char *M, float *Output, unsigned char *M_upd, int SW_increment, int iterationsNumb, int trigger, int dimX, int dimY, int dimZ);
+cdef extern float TV_energy2D(float *U, float *U0, float *E_val, float lambdaPar, int type, int dimX, int dimY);
+cdef extern float TV_energy3D(float *U, float *U0, float *E_val, float lambdaPar, int type, int dimX, int dimY, int dimZ);
+#****************************************************************#
+#********************** Total-variation ROF *********************#
+#****************************************************************#
+def TV_ROF_CPU(inputData, regularisation_parameter, iterationsNumb, marching_step_parameter):
+ if inputData.ndim == 2:
+ return TV_ROF_2D(inputData, regularisation_parameter, iterationsNumb, marching_step_parameter)
+ elif inputData.ndim == 3:
+ 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 regularisation_parameter,
+ int iterationsNumb,
+ float marching_step_parameter):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run ROF iterations for 2D data
+ TV_ROF_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, iterationsNumb, marching_step_parameter, dims[1], dims[0], 1)
+
+ return outputData
+
+def TV_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float marching_step_parameter):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run ROF iterations for 3D data
+ TV_ROF_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, iterationsNumb, marching_step_parameter, dims[2], dims[1], dims[0])
+
+ return outputData
+
+#****************************************************************#
+#********************** Total-variation FGP *********************#
+#****************************************************************#
+#******** Total-variation Fast-Gradient-Projection (FGP)*********#
+def TV_FGP_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM):
+ if inputData.ndim == 2:
+ return TV_FGP_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, nonneg, printM)
+ elif inputData.ndim == 3:
+ 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 regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run FGP-TV iterations for 2D data */
+ TV_FGP_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM,
+ dims[1],dims[0],1)
+
+ return outputData
+
+def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ int methodTV,
+ int nonneg,
+ int printM):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run FGP-TV iterations for 3D data */
+ TV_FGP_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM,
+ dims[2], dims[1], dims[0])
+ return outputData
+
+#***************************************************************#
+#********************** Total-variation SB *********************#
+#***************************************************************#
+#*************** Total-variation Split Bregman (SB)*************#
+def TV_SB_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM):
+ if inputData.ndim == 2:
+ return TV_SB_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM)
+ elif inputData.ndim == 3:
+ return TV_SB_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM)
+
+def TV_SB_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ int methodTV,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run SB-TV iterations for 2D data */
+ SB_TV_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ methodTV,
+ printM,
+ dims[1],dims[0],1)
+
+ return outputData
+
+def TV_SB_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ int methodTV,
+ int printM):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run SB-TV iterations for 3D data */
+ SB_TV_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ methodTV,
+ printM,
+ dims[2], dims[1], dims[0])
+ return outputData
+
+#***************************************************************#
+#***************** Total Generalised Variation *****************#
+#***************************************************************#
+def TGV_CPU(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst):
+ if inputData.ndim == 2:
+ return TGV_2D(inputData, regularisation_parameter, alpha1, alpha0,
+ iterations, LipshitzConst)
+ elif inputData.ndim == 3:
+ return TGV_3D(inputData, regularisation_parameter, alpha1, alpha0,
+ iterations, LipshitzConst)
+
+def TGV_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run TGV iterations for 2D data */
+ TGV_main(&inputData[0,0], &outputData[0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[1],dims[0],1)
+ return outputData
+def TGV_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run TGV iterations for 3D data */
+ TGV_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[2], dims[1], dims[0])
+ return outputData
+
+#***************************************************************#
+#******************* ROF - LLT regularisation ******************#
+#***************************************************************#
+def LLT_ROF_CPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter):
+ if inputData.ndim == 2:
+ return LLT_ROF_2D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+ elif inputData.ndim == 3:
+ return LLT_ROF_3D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+
+def LLT_ROF_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameterROF,
+ float regularisation_parameterLLT,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run ROF-LLT iterations for 2D data */
+ LLT_ROF_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, dims[1],dims[0],1)
+ return outputData
+
+def LLT_ROF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameterROF,
+ float regularisation_parameterLLT,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run ROF-LLT iterations for 3D data */
+ LLT_ROF_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, dims[2], dims[1], dims[0])
+ return outputData
+
+#****************************************************************#
+#**************Directional Total-variation FGP ******************#
+#****************************************************************#
+#******** Directional TV Fast-Gradient-Projection (FGP)*********#
+def dTV_FGP_CPU(inputData, refdata, regularisation_parameter, iterationsNumb, tolerance_param, eta_const, methodTV, nonneg, printM):
+ if inputData.ndim == 2:
+ return dTV_FGP_2D(inputData, refdata, regularisation_parameter, iterationsNumb, tolerance_param, eta_const, methodTV, nonneg, printM)
+ elif inputData.ndim == 3:
+ return dTV_FGP_3D(inputData, refdata, regularisation_parameter, iterationsNumb, tolerance_param, eta_const, methodTV, nonneg, printM)
+
+def dTV_FGP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=2, mode="c"] refdata,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ float eta_const,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run FGP-dTV iterations for 2D data */
+ dTV_FGP_CPU_main(&inputData[0,0], &refdata[0,0], &outputData[0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM,
+ dims[1], dims[0], 1)
+
+ return outputData
+
+def dTV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=3, mode="c"] refdata,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param,
+ float eta_const,
+ int methodTV,
+ int nonneg,
+ int printM):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
+
+ #/* Run FGP-dTV iterations for 3D data */
+ dTV_FGP_CPU_main(&inputData[0,0,0], &refdata[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ iterationsNumb,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM,
+ dims[2], dims[1], dims[0])
+ return outputData
+
+#****************************************************************#
+#*********************Total Nuclear Variation********************#
+#****************************************************************#
+def TNV_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param):
+ if inputData.ndim == 2:
+ return
+ elif inputData.ndim == 3:
+ return TNV_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param)
+
+def TNV_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterationsNumb,
+ float tolerance_param):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run TNV iterations for 3D (X,Y,Channels) data
+ TNV_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, iterationsNumb, tolerance_param, dims[2], dims[1], dims[0])
+ return outputData
+#****************************************************************#
+#***************Nonlinear (Isotropic) Diffusion******************#
+#****************************************************************#
+def NDF_CPU(inputData, regularisation_parameter, edge_parameter, iterationsNumb,time_marching_parameter, penalty_type):
+ if inputData.ndim == 2:
+ return NDF_2D(inputData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type)
+ elif inputData.ndim == 3:
+ return NDF_3D(inputData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type)
+
+def NDF_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run Nonlinear Diffusion iterations for 2D data
+ Diffusion_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[1], dims[0], 1)
+ return outputData
+
+def NDF_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run Nonlinear Diffusion iterations for 3D data
+ Diffusion_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[2], dims[1], dims[0])
+
+ return outputData
+
+#****************************************************************#
+#*************Anisotropic Fourth-Order diffusion*****************#
+#****************************************************************#
+def Diff4th_CPU(inputData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter):
+ if inputData.ndim == 2:
+ return Diff4th_2D(inputData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter)
+ elif inputData.ndim == 3:
+ return Diff4th_3D(inputData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter)
+
+def Diff4th_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run Anisotropic Fourth-Order diffusion for 2D data
+ Diffus4th_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, dims[1], dims[0], 1)
+ return outputData
+
+def Diff4th_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run Anisotropic Fourth-Order diffusion for 3D data
+ Diffus4th_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, dims[2], dims[1], dims[0])
+
+ return outputData
+
+#****************************************************************#
+#***************Patch-based weights calculation******************#
+#****************************************************************#
+def PATCHSEL_CPU(inputData, searchwindow, patchwindow, neighbours, edge_parameter):
+ if inputData.ndim == 2:
+ return PatchSel_2D(inputData, searchwindow, patchwindow, neighbours, edge_parameter)
+ elif inputData.ndim == 3:
+ return 1
+def PatchSel_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ int searchwindow,
+ int patchwindow,
+ int neighbours,
+ float edge_parameter):
+ cdef long dims[3]
+ dims[0] = neighbours
+ dims[1] = inputData.shape[0]
+ dims[2] = inputData.shape[1]
+
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] Weights = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='float32')
+
+ cdef np.ndarray[np.uint16_t, ndim=3, mode="c"] H_i = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='uint16')
+
+ cdef np.ndarray[np.uint16_t, ndim=3, mode="c"] H_j = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='uint16')
+
+ # Run patch-based weight selection function
+ PatchSelect_CPU_main(&inputData[0,0], &H_j[0,0,0], &H_i[0,0,0], &H_i[0,0,0], &Weights[0,0,0], dims[2], dims[1], 0, searchwindow, patchwindow, neighbours, edge_parameter, 1)
+ return H_i, H_j, Weights
+"""
+def PatchSel_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ int searchwindow,
+ int patchwindow,
+ int neighbours,
+ float edge_parameter):
+ cdef long dims[4]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+ dims[3] = neighbours
+
+ cdef np.ndarray[np.float32_t, ndim=4, mode="c"] Weights = \
+ np.zeros([dims[3],dims[0],dims[1],dims[2]], dtype='float32')
+
+ cdef np.ndarray[np.uint16_t, ndim=4, mode="c"] H_i = \
+ np.zeros([dims[3],dims[0],dims[1],dims[2]], dtype='uint16')
+
+ cdef np.ndarray[np.uint16_t, ndim=4, mode="c"] H_j = \
+ np.zeros([dims[3],dims[0],dims[1],dims[2]], dtype='uint16')
+
+ cdef np.ndarray[np.uint16_t, ndim=4, mode="c"] H_k = \
+ np.zeros([dims[3],dims[0],dims[1],dims[2]], dtype='uint16')
+
+ # Run patch-based weight selection function
+ PatchSelect_CPU_main(&inputData[0,0,0], &H_i[0,0,0,0], &H_j[0,0,0,0], &H_k[0,0,0,0], &Weights[0,0,0,0], dims[2], dims[1], dims[0], searchwindow, patchwindow, neighbours, edge_parameter, 1)
+ return H_i, H_j, H_k, Weights
+"""
+
+#****************************************************************#
+#***************Non-local Total Variation******************#
+#****************************************************************#
+def NLTV_CPU(inputData, H_i, H_j, H_k, Weights, regularisation_parameter, iterations):
+ if inputData.ndim == 2:
+ return NLTV_2D(inputData, H_i, H_j, Weights, regularisation_parameter, iterations)
+ elif inputData.ndim == 3:
+ return 1
+def NLTV_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.uint16_t, ndim=3, mode="c"] H_i,
+ np.ndarray[np.uint16_t, ndim=3, mode="c"] H_j,
+ np.ndarray[np.float32_t, ndim=3, mode="c"] Weights,
+ float regularisation_parameter,
+ int iterations):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ neighbours = H_i.shape[0]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run nonlocal TV regularisation
+ Nonlocal_TV_CPU_main(&inputData[0,0], &outputData[0,0], &H_i[0,0,0], &H_j[0,0,0], &H_i[0,0,0], &Weights[0,0,0], dims[1], dims[0], 0, neighbours, regularisation_parameter, iterations)
+ return outputData
+
+#*********************Inpainting WITH****************************#
+#***************Nonlinear (Isotropic) Diffusion******************#
+#****************************************************************#
+def NDF_INPAINT_CPU(inputData, maskData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type):
+ if inputData.ndim == 2:
+ return NDF_INP_2D(inputData, maskData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type)
+ elif inputData.ndim == 3:
+ return NDF_INP_3D(inputData, maskData, regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type)
+
+def NDF_INP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.uint8_t, ndim=2, mode="c"] maskData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run Inpaiting by Diffusion iterations for 2D data
+ Diffusion_Inpaint_CPU_main(&inputData[0,0], &maskData[0,0], &outputData[0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[1], dims[0], 1)
+ return outputData
+
+def NDF_INP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ np.ndarray[np.uint8_t, ndim=3, mode="c"] maskData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run Inpaiting by Diffusion iterations for 3D data
+ Diffusion_Inpaint_CPU_main(&inputData[0,0,0], &maskData[0,0,0], &outputData[0,0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[2], dims[1], dims[0])
+
+ return outputData
+#*********************Inpainting WITH****************************#
+#***************Nonlocal Vertical Marching method****************#
+#****************************************************************#
+def NVM_INPAINT_CPU(inputData, maskData, SW_increment, iterationsNumb):
+ if inputData.ndim == 2:
+ return NVM_INP_2D(inputData, maskData, SW_increment, iterationsNumb)
+ elif inputData.ndim == 3:
+ return
+
+def NVM_INP_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.uint8_t, ndim=2, mode="c"] maskData,
+ int SW_increment,
+ int iterationsNumb):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ cdef np.ndarray[np.uint8_t, ndim=2, mode="c"] maskData_upd = \
+ np.zeros([dims[0],dims[1]], dtype='uint8')
+
+ # Run Inpaiting by Nonlocal vertical marching method for 2D data
+ NonlocalMarching_Inpaint_main(&inputData[0,0], &maskData[0,0], &outputData[0,0],
+ &maskData_upd[0,0],
+ SW_increment, iterationsNumb, 1, dims[1], dims[0], 1)
+
+ return (outputData, maskData_upd)
+
+
+#****************************************************************#
+#***************Calculation of TV-energy functional**************#
+#****************************************************************#
+def TV_ENERGY(inputData, inputData0, regularisation_parameter, typeFunctional):
+ if inputData.ndim == 2:
+ return TV_ENERGY_2D(inputData, inputData0, regularisation_parameter, typeFunctional)
+ elif inputData.ndim == 3:
+ return TV_ENERGY_3D(inputData, inputData0, regularisation_parameter, typeFunctional)
+
+def TV_ENERGY_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=2, mode="c"] inputData0,
+ float regularisation_parameter,
+ int typeFunctional):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=1, mode="c"] outputData = \
+ np.zeros([1], dtype='float32')
+
+ # run function
+ TV_energy2D(&inputData[0,0], &inputData0[0,0], &outputData[0], regularisation_parameter, typeFunctional, dims[1], dims[0])
+
+ return outputData
+
+def TV_ENERGY_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=3, mode="c"] inputData0,
+ float regularisation_parameter,
+ int typeFunctional):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=1, mode="c"] outputData = \
+ np.zeros([1], dtype='float32')
+
+ # Run function
+ TV_energy3D(&inputData[0,0,0], &inputData0[0,0,0], &outputData[0], regularisation_parameter, typeFunctional, dims[2], dims[1], dims[0])
+
+ return outputData
diff --git a/src/Python/src/gpu_regularisers.pyx b/src/Python/src/gpu_regularisers.pyx
new file mode 100644
index 0000000..b52f669
--- /dev/null
+++ b/src/Python/src/gpu_regularisers.pyx
@@ -0,0 +1,640 @@
+# distutils: language=c++
+"""
+Copyright 2018 CCPi
+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.
+
+Author: Edoardo Pasca, Daniil Kazantsev
+"""
+
+import cython
+import numpy as np
+cimport numpy as np
+
+CUDAErrorMessage = 'CUDA error'
+
+cdef extern int TV_ROF_GPU_main(float* Input, float* Output, float lambdaPar, int iter, float tau, int N, int M, int Z);
+cdef extern int TV_FGP_GPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int printM, int N, int M, int Z);
+cdef extern int TV_SB_GPU_main(float *Input, float *Output, float lambdaPar, int iter, float epsil, int methodTV, int printM, int N, int M, int Z);
+cdef extern int TGV_GPU_main(float *Input, float *Output, float lambdaPar, float alpha1, float alpha0, int iterationsNumb, float L2, int dimX, int dimY, int dimZ);
+cdef extern int LLT_ROF_GPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int N, int M, int Z);
+cdef extern int NonlDiff_GPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int penaltytype, int N, int M, int Z);
+cdef extern int dTV_FGP_GPU_main(float *Input, float *InputRef, float *Output, float lambdaPar, int iterationsNumb, float epsil, float eta, int methodTV, int nonneg, int printM, int N, int M, int Z);
+cdef extern int Diffus4th_GPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int N, int M, int Z);
+cdef extern int PatchSelect_GPU_main(float *Input, unsigned short *H_i, unsigned short *H_j, float *Weights, int N, int M, int SearchWindow, int SimilarWin, int NumNeighb, float h);
+
+# Total-variation Rudin-Osher-Fatemi (ROF)
+def TV_ROF_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ time_marching_parameter):
+ if inputData.ndim == 2:
+ return ROFTV2D(inputData,
+ regularisation_parameter,
+ iterations,
+ time_marching_parameter)
+ elif inputData.ndim == 3:
+ return ROFTV3D(inputData,
+ regularisation_parameter,
+ iterations,
+ time_marching_parameter)
+
+# Total-variation Fast-Gradient-Projection (FGP)
+def TV_FGP_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM):
+ if inputData.ndim == 2:
+ return FGPTV2D(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM)
+ elif inputData.ndim == 3:
+ return FGPTV3D(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM)
+# Total-variation Split Bregman (SB)
+def TV_SB_GPU(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM):
+ if inputData.ndim == 2:
+ return SBTV2D(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM)
+ elif inputData.ndim == 3:
+ return SBTV3D(inputData,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM)
+# LLT-ROF model
+def LLT_ROF_GPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter):
+ if inputData.ndim == 2:
+ return LLT_ROF_GPU2D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+ elif inputData.ndim == 3:
+ return LLT_ROF_GPU3D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+# Total Generilised Variation (TGV)
+def TGV_GPU(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst):
+ if inputData.ndim == 2:
+ return TGV2D(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst)
+ elif inputData.ndim == 3:
+ return TGV3D(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst)
+# Directional Total-variation Fast-Gradient-Projection (FGP)
+def dTV_FGP_GPU(inputData,
+ refdata,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM):
+ if inputData.ndim == 2:
+ return FGPdTV2D(inputData,
+ refdata,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM)
+ elif inputData.ndim == 3:
+ return FGPdTV3D(inputData,
+ refdata,
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM)
+# Nonlocal Isotropic Diffusion (NDF)
+def NDF_GPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter,
+ penalty_type):
+ if inputData.ndim == 2:
+ return NDF_GPU_2D(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter,
+ penalty_type)
+ elif inputData.ndim == 3:
+ return NDF_GPU_3D(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter,
+ penalty_type)
+# Anisotropic Fourth-Order diffusion
+def Diff4th_GPU(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter):
+ if inputData.ndim == 2:
+ return Diff4th_2D(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter)
+ elif inputData.ndim == 3:
+ return Diff4th_3D(inputData,
+ regularisation_parameter,
+ edge_parameter,
+ iterations,
+ time_marching_parameter)
+
+#****************************************************************#
+#********************** Total-variation ROF *********************#
+#****************************************************************#
+def ROFTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_ROF_GPU_main(
+ &inputData[0,0], &outputData[0,0],
+ regularisation_parameter,
+ iterations ,
+ time_marching_parameter,
+ dims[1], dims[0], 1)==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+def ROFTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_ROF_GPU_main(
+ &inputData[0,0,0], &outputData[0,0,0],
+ regularisation_parameter,
+ iterations ,
+ time_marching_parameter,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+#****************************************************************#
+#********************** Total-variation FGP *********************#
+#****************************************************************#
+#******** Total-variation Fast-Gradient-Projection (FGP)*********#
+def FGPTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_FGP_GPU_main(&inputData[0,0], &outputData[0,0],
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM,
+ dims[1], dims[0], 1)==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def FGPTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_FGP_GPU_main(&inputData[0,0,0], &outputData[0,0,0],
+ regularisation_parameter ,
+ iterations,
+ tolerance_param,
+ methodTV,
+ nonneg,
+ printM,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+#***************************************************************#
+#********************** Total-variation SB *********************#
+#***************************************************************#
+#*************** Total-variation Split Bregman (SB)*************#
+def SBTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ int methodTV,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_SB_GPU_main(&inputData[0,0], &outputData[0,0],
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM,
+ dims[1], dims[0], 1)==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def SBTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ int methodTV,
+ int printM):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (TV_SB_GPU_main(&inputData[0,0,0], &outputData[0,0,0],
+ regularisation_parameter ,
+ iterations,
+ tolerance_param,
+ methodTV,
+ printM,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+#***************************************************************#
+#************************ LLT-ROF model ************************#
+#***************************************************************#
+#************Joint LLT-ROF model for higher order **************#
+def LLT_ROF_GPU2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameterROF,
+ float regularisation_parameterLLT,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Running CUDA code here
+ if (LLT_ROF_GPU_main(&inputData[0,0], &outputData[0,0],regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, dims[1],dims[0],1)==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def LLT_ROF_GPU3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameterROF,
+ float regularisation_parameterLLT,
+ int iterations,
+ float time_marching_parameter):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (LLT_ROF_GPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+#***************************************************************#
+#***************** Total Generalised Variation *****************#
+#***************************************************************#
+def TGV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #/* Run TGV iterations for 2D data */
+ if (TGV_GPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[1],dims[0], 1)==0):
+ return outputData
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+def TGV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float alpha1,
+ float alpha0,
+ int iterationsNumb,
+ float LipshitzConst):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (TGV_GPU_main(
+ &inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+ alpha1,
+ alpha0,
+ iterationsNumb,
+ LipshitzConst,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+#****************************************************************#
+#**************Directional Total-variation FGP ******************#
+#****************************************************************#
+#******** Directional TV Fast-Gradient-Projection (FGP)*********#
+def FGPdTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=2, mode="c"] refdata,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ float eta_const,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Running CUDA code here
+ if (dTV_FGP_GPU_main(&inputData[0,0], &refdata[0,0], &outputData[0,0],
+ regularisation_parameter,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM,
+ dims[1], dims[0], 1)==0):
+ return outputData
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def FGPdTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ np.ndarray[np.float32_t, ndim=3, mode="c"] refdata,
+ float regularisation_parameter,
+ int iterations,
+ float tolerance_param,
+ float eta_const,
+ int methodTV,
+ int nonneg,
+ int printM):
+
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Running CUDA code here
+ if (dTV_FGP_GPU_main(&inputData[0,0,0], &refdata[0,0,0], &outputData[0,0,0],
+ regularisation_parameter ,
+ iterations,
+ tolerance_param,
+ eta_const,
+ methodTV,
+ nonneg,
+ printM,
+ dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+#****************************************************************#
+#***************Nonlinear (Isotropic) Diffusion******************#
+#****************************************************************#
+def NDF_GPU_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ #rangecheck = penalty_type < 1 and penalty_type > 3
+ #if not rangecheck:
+# raise ValueError('Choose penalty type as 1 for Huber, 2 - Perona-Malik, 3 - Tukey Biweight')
+
+ # Run Nonlinear Diffusion iterations for 2D data
+ # Running CUDA code here
+ if (NonlDiff_GPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[1], dims[0], 1)==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def NDF_GPU_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter,
+ int penalty_type):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run Nonlinear Diffusion iterations for 3D data
+ # Running CUDA code here
+ if (NonlDiff_GPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, penalty_type, dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+#****************************************************************#
+#************Anisotropic Fourth-Order diffusion******************#
+#****************************************************************#
+def Diff4th_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter):
+ cdef long dims[2]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1]], dtype='float32')
+
+ # Run Anisotropic Fourth-Order diffusion for 2D data
+ # Running CUDA code here
+ if (Diffus4th_GPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, dims[1], dims[0], 1)==0):
+ return outputData
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+
+def Diff4th_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
+ float regularisation_parameter,
+ float edge_parameter,
+ int iterationsNumb,
+ float time_marching_parameter):
+ cdef long dims[3]
+ dims[0] = inputData.shape[0]
+ dims[1] = inputData.shape[1]
+ dims[2] = inputData.shape[2]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
+ np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+
+ # Run Anisotropic Fourth-Order diffusion for 3D data
+ # Running CUDA code here
+ if (Diffus4th_GPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter, edge_parameter, iterationsNumb, time_marching_parameter, dims[2], dims[1], dims[0])==0):
+ return outputData;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
+#****************************************************************#
+#************Patch-based weights pre-selection******************#
+#****************************************************************#
+def PATCHSEL_GPU(inputData, searchwindow, patchwindow, neighbours, edge_parameter):
+ if inputData.ndim == 2:
+ return PatchSel_2D(inputData, searchwindow, patchwindow, neighbours, edge_parameter)
+ elif inputData.ndim == 3:
+ return 1
+def PatchSel_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
+ int searchwindow,
+ int patchwindow,
+ int neighbours,
+ float edge_parameter):
+ cdef long dims[3]
+ dims[0] = neighbours
+ dims[1] = inputData.shape[0]
+ dims[2] = inputData.shape[1]
+
+ cdef np.ndarray[np.float32_t, ndim=3, mode="c"] Weights = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='float32')
+
+ cdef np.ndarray[np.uint16_t, ndim=3, mode="c"] H_i = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='uint16')
+
+ cdef np.ndarray[np.uint16_t, ndim=3, mode="c"] H_j = \
+ np.zeros([dims[0], dims[1],dims[2]], dtype='uint16')
+
+ # Run patch-based weight selection function
+ if (PatchSelect_GPU_main(&inputData[0,0], &H_j[0,0,0], &H_i[0,0,0], &Weights[0,0,0], dims[2], dims[1], searchwindow, patchwindow, neighbours, edge_parameter)==0):
+ return H_i, H_j, Weights;
+ else:
+ raise ValueError(CUDAErrorMessage);
+
generated by cgit v1.2.3 (git 2.25.1) at 2025-11-04 11:28:57 +0000