SBTV,LLTROF completed

3 files changed, 127 insertions, 107 deletions

diff --git a/src/Python/ccpi/filters/regularisers.py b/src/Python/ccpi/filters/regularisers.py
index 05120af..09b465a 100644
--- a/src/Python/ccpi/filters/regularisers.py
+++ b/src/Python/ccpi/filters/regularisers.py
@@ -52,21 +52,30 @@ def FGP_TV(inputData, regularisation_parameter,iterations,
         raise ValueError('Unknown device {0}. Expecting gpu or cpu'\
                          .format(device))
 def SB_TV(inputData, regularisation_parameter, iterations,
-                     tolerance_param, methodTV, printM, device='cpu'):
+                     tolerance_param, methodTV, device='cpu'):
     if device == 'cpu':
         return TV_SB_CPU(inputData,
                      regularisation_parameter,
                      iterations, 
                      tolerance_param,
-                     methodTV,
-                     printM)
+                     methodTV)
     elif device == 'gpu' and gpu_enabled:
         return TV_SB_GPU(inputData,
                      regularisation_parameter,
                      iterations, 
                      tolerance_param,
-                     methodTV,
-                     printM)
+                     methodTV)
+    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, tolerance_param, device='cpu'):
+    if device == 'cpu':
+        return LLT_ROF_CPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
+    elif device == 'gpu' and gpu_enabled:
+        return LLT_ROF_GPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
     else:
         if not gpu_enabled and device == 'gpu':
             raise ValueError ('GPU is not available')
@@ -194,17 +203,6 @@ def TGV(inputData, regularisation_parameter, alpha1, alpha0, iterations,
             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, 
diff --git a/src/Python/src/cpu_regularisers.pyx b/src/Python/src/cpu_regularisers.pyx
index aeca141..f2276bb 100644
--- a/src/Python/src/cpu_regularisers.pyx
+++ b/src/Python/src/cpu_regularisers.pyx
@@ -20,8 +20,8 @@ cimport numpy as np
 
 cdef extern float TV_ROF_CPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iterationsNumb, float tau, float epsil, int dimX, int dimY, int dimZ);
 cdef extern float TV_FGP_CPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iterationsNumb, float epsil, int methodTV, int nonneg, int dimX, int dimY, int dimZ);
-cdef extern float 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 SB_TV_CPU_main(float *Input, float *Output, float *infovector, float mu, int iter, float epsil, int methodTV, int dimX, int dimY, int dimZ);
+cdef extern float LLT_ROF_CPU_main(float *Input, float *Output, float *infovector, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, float epsil, int dimX, int dimY, int dimZ);
 cdef extern float TGV_main(float *Input, float *Output, float 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);
@@ -149,18 +149,17 @@ def TV_FGP_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
 #********************** Total-variation SB *********************#
 #***************************************************************#
 #*************** Total-variation Split Bregman (SB)*************#
-def TV_SB_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM):
+def TV_SB_CPU(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV):
     if inputData.ndim == 2:
-        return TV_SB_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM)
+        return TV_SB_2D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV)
     elif inputData.ndim == 3:
-        return TV_SB_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV, printM)
+        return TV_SB_3D(inputData, regularisation_parameter, iterationsNumb, tolerance_param, methodTV)
 
 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):
+                     int methodTV):
                          
     cdef long dims[2]
     dims[0] = inputData.shape[0]
@@ -168,23 +167,24 @@ def TV_SB_2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
     
     cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
             np.zeros([dims[0],dims[1]], dtype='float32')
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.zeros([2], dtype='float32')
                    
     #/* Run SB-TV iterations for 2D data */
-    SB_TV_CPU_main(&inputData[0,0], &outputData[0,0], regularisation_parameter, 
+    SB_TV_CPU_main(&inputData[0,0], &outputData[0,0], &infovec[0],
+                       regularisation_parameter, 
                        iterationsNumb, 
                        tolerance_param,
                        methodTV,
-                       printM,
-                       dims[1],dims[0],1)
+                       dims[1],dims[0], 1)
     
-    return outputData        
+    return (outputData,infovec)
             
 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):
+                     int methodTV):
     cdef long dims[3]
     dims[0] = inputData.shape[0]
     dims[1] = inputData.shape[1]
@@ -192,16 +192,71 @@ def TV_SB_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
     
     cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
             np.zeros([dims[0], dims[1], dims[2]], dtype='float32')
-           
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.zeros([2], dtype='float32')
+            
     #/* Run SB-TV iterations for 3D data */
-    SB_TV_CPU_main(&inputData[0,0,0], &outputData[0,0,0], regularisation_parameter,
+    SB_TV_CPU_main(&inputData[0,0,0], &outputData[0,0,0], &infovec[0],
+                       regularisation_parameter,
                        iterationsNumb, 
                        tolerance_param,
                        methodTV,
-                       printM,
                        dims[2], dims[1], dims[0])
-    return outputData 
+    return (outputData,infovec) 
+#***************************************************************#
+#******************* ROF - LLT regularisation ******************#
+#***************************************************************#
+def LLT_ROF_CPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param):
+    if inputData.ndim == 2:
+        return LLT_ROF_2D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
+    elif inputData.ndim == 3:
+        return LLT_ROF_3D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
 
+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,
+                     float tolerance_param):
+                         
+    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.float32_t, ndim=1, mode="c"] infovec = \
+            np.zeros([2], dtype='float32')
+            
+    #/* Run ROF-LLT iterations for 2D data */
+    LLT_ROF_CPU_main(&inputData[0,0], &outputData[0,0], &infovec[0], regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, 
+                     tolerance_param,
+                     dims[1],dims[0],1)
+    return (outputData,infovec) 
+
+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,
+                     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')
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.zeros([2], dtype='float32')
+            
+    #/* Run ROF-LLT iterations for 3D data */
+    LLT_ROF_CPU_main(&inputData[0,0,0], &outputData[0,0,0], &infovec[0], regularisation_parameterROF, regularisation_parameterLLT, iterations, 
+                     time_marching_parameter,
+                     tolerance_param, 
+                     dims[2], dims[1], dims[0])
+    return (outputData,infovec)
 #***************************************************************#
 #***************** Total Generalised Variation *****************#
 #***************************************************************#
@@ -259,49 +314,6 @@ def TGV_3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
                        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 ******************#
diff --git a/src/Python/src/gpu_regularisers.pyx b/src/Python/src/gpu_regularisers.pyx
index db4fa67..6fc4135 100644
--- a/src/Python/src/gpu_regularisers.pyx
+++ b/src/Python/src/gpu_regularisers.pyx
@@ -22,9 +22,9 @@ CUDAErrorMessage = 'CUDA error'
 
 cdef extern int TV_ROF_GPU_main(float* Input, float* Output, float *infovector, float lambdaPar, int iter, float tau, float epsil, int N, int M, int Z);
 cdef extern int TV_FGP_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, int methodTV, int nonneg, int N, int M, int Z);
-cdef extern int TV_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 TV_SB_GPU_main(float *Input, float *Output, float *infovector, float lambdaPar, int iter, float epsil, int methodTV, int N, int M, int Z);
+cdef extern int LLT_ROF_GPU_main(float *Input, float *Output, float *infovector, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau,  float epsil, int N, int M, int Z);
 cdef extern int TGV_GPU_main(float *Input, float *Output, float 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);
@@ -75,28 +75,25 @@ def TV_SB_GPU(inputData,
                      regularisation_parameter,
                      iterations, 
                      tolerance_param,
-                     methodTV,
-                     printM):
+                     methodTV):
     if inputData.ndim == 2:
         return SBTV2D(inputData,
                      regularisation_parameter,
                      iterations, 
                      tolerance_param,
-                     methodTV,
-                     printM)
+                     methodTV)
     elif inputData.ndim == 3:
         return SBTV3D(inputData,
                      regularisation_parameter,
                      iterations, 
                      tolerance_param,
-                     methodTV,
-                     printM)
+                     methodTV)
 # LLT-ROF model
-def LLT_ROF_GPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter):
+def LLT_ROF_GPU(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param):
     if inputData.ndim == 2:
-        return LLT_ROF_GPU2D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+        return LLT_ROF_GPU2D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
     elif inputData.ndim == 3:
-        return LLT_ROF_GPU3D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter)
+        return LLT_ROF_GPU3D(inputData, regularisation_parameterROF, regularisation_parameterLLT, iterations, time_marching_parameter, tolerance_param)
 # Total Generilised Variation (TGV)
 def TGV_GPU(inputData, regularisation_parameter, alpha1, alpha0, iterations, LipshitzConst):
     if inputData.ndim == 2:
@@ -300,8 +297,7 @@ def SBTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
                      float regularisation_parameter,
                      int iterations, 
                      float tolerance_param,
-                     int methodTV,
-                     int printM):
+                     int methodTV):
     
     cdef long dims[2]
     dims[0] = inputData.shape[0]
@@ -309,16 +305,17 @@ def SBTV2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
 
     cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
 		    np.zeros([dims[0],dims[1]], dtype='float32')
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.ones([2], dtype='float32')
           
-    # Running CUDA code here    
-    if (TV_SB_GPU_main(&inputData[0,0], &outputData[0,0],
+    # Running CUDA code here
+    if (TV_SB_GPU_main(&inputData[0,0], &outputData[0,0],&infovec[0],
                        regularisation_parameter, 
                        iterations, 
                        tolerance_param,
                        methodTV,
-                       printM,
                        dims[1], dims[0], 1)==0):
-        return outputData;
+        return (outputData,infovec)
     else:
         raise ValueError(CUDAErrorMessage);
 
@@ -327,8 +324,7 @@ def SBTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
                      float regularisation_parameter,
                      int iterations, 
                      float tolerance_param,
-                     int methodTV,
-                     int printM):
+                     int methodTV):
     
     cdef long dims[3]
     dims[0] = inputData.shape[0]
@@ -337,16 +333,17 @@ def SBTV3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
 
     cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
 		    np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
-          
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.ones([2], dtype='float32')
+            
     # Running CUDA code here    
-    if (TV_SB_GPU_main(&inputData[0,0,0], &outputData[0,0,0],
+    if (TV_SB_GPU_main(&inputData[0,0,0], &outputData[0,0,0],&infovec[0],
                        regularisation_parameter , 
                        iterations, 
                        tolerance_param,
                        methodTV,
-                       printM,
                        dims[2], dims[1], dims[0])==0):
-        return outputData;
+        return (outputData,infovec)
     else:
         raise ValueError(CUDAErrorMessage);
 
@@ -359,7 +356,8 @@ 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):
+                     float time_marching_parameter,
+                     float tolerance_param):
     
     cdef long dims[2]
     dims[0] = inputData.shape[0]
@@ -367,10 +365,15 @@ def LLT_ROF_GPU2D(np.ndarray[np.float32_t, ndim=2, mode="c"] inputData,
 
     cdef np.ndarray[np.float32_t, ndim=2, mode="c"] outputData = \
 		    np.zeros([dims[0],dims[1]], dtype='float32')
-          
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.ones([2], 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;
+    if (LLT_ROF_GPU_main(&inputData[0,0], &outputData[0,0],&infovec[0],regularisation_parameterROF, regularisation_parameterLLT, iterations, 
+                         time_marching_parameter, 
+                         tolerance_param,
+                         dims[1],dims[0],1)==0):
+        return (outputData,infovec)
     else:
         raise ValueError(CUDAErrorMessage);
 
@@ -379,7 +382,8 @@ 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):
+                     float time_marching_parameter,
+                     float tolerance_param):
     
     cdef long dims[3]
     dims[0] = inputData.shape[0]
@@ -388,10 +392,16 @@ def LLT_ROF_GPU3D(np.ndarray[np.float32_t, ndim=3, mode="c"] inputData,
 
     cdef np.ndarray[np.float32_t, ndim=3, mode="c"] outputData = \
 		    np.zeros([dims[0],dims[1],dims[2]], dtype='float32')
+    cdef np.ndarray[np.float32_t, ndim=1, mode="c"] infovec = \
+            np.ones([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;
+    # Running CUDA code here
+    if (LLT_ROF_GPU_main(&inputData[0,0,0], &outputData[0,0,0], &infovec[0], regularisation_parameterROF, regularisation_parameterLLT, 
+                         iterations, 
+                         time_marching_parameter, 
+                         tolerance_param,
+                         dims[2], dims[1], dims[0])==0):
+        return (outputData,infovec)
     else:
         raise ValueError(CUDAErrorMessage);