From c16582003fb3a52e528f0a415e932c2e869e9952 Mon Sep 17 00:00:00 2001
From: Kazantsev <kjy41806@ws097.diamond.ac.uk>
Date: Fri, 14 Sep 2018 12:29:26 +0100
Subject: CPU modules fixed for big data

---
 Core/inpainters_CPU/Diffusion_Inpaint_core.c       | 32 ++++----
 Core/inpainters_CPU/Diffusion_Inpaint_core.h       |  9 ++-
 .../inpainters_CPU/NonlocalMarching_Inpaint_core.c |  2 +-
 Core/regularisers_CPU/Diffus4th_order_core.c       | 22 +++---
 Core/regularisers_CPU/Diffus4th_order_core.h       |  8 +-
 Core/regularisers_CPU/Diffusion_core.c             | 26 +++----
 Core/regularisers_CPU/Diffusion_core.h             |  8 +-
 Core/regularisers_CPU/FGP_TV_core.c                | 61 ++++++++--------
 Core/regularisers_CPU/FGP_TV_core.h                | 16 ++--
 Core/regularisers_CPU/FGP_dTV_core.c               | 85 +++++++++++-----------
 Core/regularisers_CPU/FGP_dTV_core.h               | 24 +++---
 Core/regularisers_CPU/LLT_ROF_core.c               | 53 +++++++-------
 Core/regularisers_CPU/LLT_ROF_core.h               | 14 ++--
 Core/regularisers_CPU/ROF_TV_core.c                | 37 +++++-----
 Core/regularisers_CPU/ROF_TV_core.h                | 11 +--
 Core/regularisers_CPU/SB_TV_core.c                 | 73 ++++++++++---------
 Core/regularisers_CPU/SB_TV_core.h                 | 18 ++---
 Core/regularisers_CPU/TGV_core.c                   | 59 +++++++--------
 Core/regularisers_CPU/TGV_core.h                   | 15 ++--
 Core/regularisers_CPU/TNV_core.c                   | 69 +++++++++---------
 Core/regularisers_CPU/TNV_core.h                   |  8 +-
 Core/regularisers_CPU/utils.c                      | 18 ++---
 Core/regularisers_CPU/utils.h                      |  2 +-
 Wrappers/Matlab/mex_compile/compileCPU_mex_Linux.m |  2 +-
 24 files changed, 339 insertions(+), 333 deletions(-)

diff --git a/Core/inpainters_CPU/Diffusion_Inpaint_core.c b/Core/inpainters_CPU/Diffusion_Inpaint_core.c
index 56e0421..08b168a 100644
--- a/Core/inpainters_CPU/Diffusion_Inpaint_core.c
+++ b/Core/inpainters_CPU/Diffusion_Inpaint_core.c
@@ -47,12 +47,12 @@ int signNDF_inc(float x) {
 
 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)
 {
-    int i,pointsone;
+    long i, pointsone;
     float sigmaPar2;
     sigmaPar2 = sigmaPar/sqrt(2.0f);
     
     /* copy into output */
-    copyIm(Input, Output, dimX, dimY, dimZ);
+    copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ));
     
     pointsone = 0;
     for (i=0; i<dimY*dimX*dimZ; i++) if (Mask[i] == 1) pointsone++;
@@ -63,17 +63,17 @@ float Diffusion_Inpaint_CPU_main(float *Input, unsigned char *Mask, float *Outpu
     if (dimZ == 1) {
     /* running 2D diffusion iterations */
     for(i=0; i < iterationsNumb; i++) {
-            if (sigmaPar == 0.0f) LinearDiff_Inp_2D(Input, Mask, Output, lambdaPar, tau, dimX, dimY); /* linear diffusion (heat equation) */
-            else NonLinearDiff_Inp_2D(Input, Mask, Output, lambdaPar, sigmaPar2, tau, penaltytype, dimX, dimY); /* nonlinear diffusion */
+            if (sigmaPar == 0.0f) LinearDiff_Inp_2D(Input, Mask, Output, lambdaPar, tau, (long)(dimX), (long)(dimY)); /* linear diffusion (heat equation) */
+            else NonLinearDiff_Inp_2D(Input, Mask, Output, lambdaPar, sigmaPar2, tau, penaltytype, (long)(dimX), (long)(dimY)); /* nonlinear diffusion */
 		}
 	}
 	else {
 	/* running 3D diffusion iterations */
     for(i=0; i < iterationsNumb; i++) {
-            if (sigmaPar == 0.0f) LinearDiff_Inp_3D(Input, Mask, Output, lambdaPar, tau, dimX, dimY, dimZ);
-            else NonLinearDiff_Inp_3D(Input, Mask, Output, lambdaPar, sigmaPar2, tau, penaltytype, dimX, dimY, dimZ);
-		}
-	}
+            if (sigmaPar == 0.0f) LinearDiff_Inp_3D(Input, Mask, Output, lambdaPar, tau, (long)(dimX), (long)(dimY), (long)(dimZ));
+            else NonLinearDiff_Inp_3D(Input, Mask, Output, lambdaPar, sigmaPar2, tau, penaltytype, (long)(dimX), (long)(dimY), (long)(dimZ));
+            }
+         }
 	}
     return *Output;
 }
@@ -81,9 +81,9 @@ float Diffusion_Inpaint_CPU_main(float *Input, unsigned char *Mask, float *Outpu
 /***************************2D Functions*****************************/
 /********************************************************************/
 /* linear diffusion (heat equation) */
-float LinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, int dimX, int dimY)
+float LinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, long dimX, long dimY)
 {
-	int i,j,i1,i2,j1,j2,index;
+	long i,j,i1,i2,j1,j2,index;
 	float e,w,n,s,e1,w1,n1,s1;
 	
 #pragma omp parallel for shared(Input,Mask) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1)
@@ -116,9 +116,9 @@ float LinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float
 }
 
 /* nonlinear diffusion */
-float NonLinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY)
+float NonLinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY)
 {
-	int i,j,i1,i2,j1,j2,index;
+	long i,j,i1,i2,j1,j2,index;
 	float e,w,n,s,e1,w1,n1,s1;
 	
 #pragma omp parallel for shared(Input,Mask) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1)
@@ -189,9 +189,9 @@ float NonLinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, flo
 /***************************3D Functions*****************************/
 /********************************************************************/
 /* linear diffusion (heat equation) */
-float LinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, int dimX, int dimY, int dimZ)
+float LinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, long dimX, long dimY, long dimZ)
 {
-	int i,j,k,i1,i2,j1,j2,k1,k2,index;
+	long i,j,k,i1,i2,j1,j2,k1,k2,index;
 	float e,w,n,s,u,d,e1,w1,n1,s1,u1,d1;
 	
 #pragma omp parallel for shared(Input,Mask) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1,k,k1,k2,u1,d1,u,d)
@@ -231,9 +231,9 @@ for(k=0; k<dimZ; k++) {
 	return *Output;
 }
 
-float NonLinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY, int dimZ)
+float NonLinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY, long dimZ)
 {
-	int i,j,k,i1,i2,j1,j2,k1,k2,index;
+	long i,j,k,i1,i2,j1,j2,k1,k2,index;
 	float e,w,n,s,u,d,e1,w1,n1,s1,u1,d1;
 	
 #pragma omp parallel for shared(Input,Mask) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1,k,k1,k2,u1,d1,u,d)
diff --git a/Core/inpainters_CPU/Diffusion_Inpaint_core.h b/Core/inpainters_CPU/Diffusion_Inpaint_core.h
index 0fc2608..a96fe79 100644
--- a/Core/inpainters_CPU/Diffusion_Inpaint_core.h
+++ b/Core/inpainters_CPU/Diffusion_Inpaint_core.h
@@ -51,10 +51,11 @@ limitations under the License.
 extern "C" {
 #endif
 CCPI_EXPORT 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);
-CCPI_EXPORT float LinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, int dimX, int dimY);
-CCPI_EXPORT float NonLinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY);
-CCPI_EXPORT float LinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float NonLinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY, int dimZ);
+
+CCPI_EXPORT float LinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, long dimX, long dimY);
+CCPI_EXPORT float NonLinearDiff_Inp_2D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY);
+CCPI_EXPORT float LinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float tau, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float NonLinearDiff_Inp_3D(float *Input, unsigned char *Mask, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/inpainters_CPU/NonlocalMarching_Inpaint_core.c b/Core/inpainters_CPU/NonlocalMarching_Inpaint_core.c
index 66be15c..b488ca4 100644
--- a/Core/inpainters_CPU/NonlocalMarching_Inpaint_core.c
+++ b/Core/inpainters_CPU/NonlocalMarching_Inpaint_core.c
@@ -35,7 +35,7 @@
  * 1. Inpainted image or a sinogram
  * 2. updated mask
  *
- * Reference: TBA
+ * Reference: D. Kazantsev (paper in preparation)
  */
 
 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)
diff --git a/Core/regularisers_CPU/Diffus4th_order_core.c b/Core/regularisers_CPU/Diffus4th_order_core.c
index 069a42e..01f4f64 100644
--- a/Core/regularisers_CPU/Diffus4th_order_core.c
+++ b/Core/regularisers_CPU/Diffus4th_order_core.c
@@ -50,7 +50,7 @@ float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sig
     W_Lapl = calloc(DimTotal, sizeof(float));
     
     /* copy into output */
-    copyIm(Input, Output, dimX, dimY, dimZ);
+    copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ));
     
     if (dimZ == 1) {
     /* running 2D diffusion iterations */
@@ -58,7 +58,7 @@ float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sig
             /* Calculating weighted Laplacian */
             Weighted_Laplc2D(W_Lapl, Output, sigmaPar2, dimX, dimY);
             /* Perform iteration step */
-            Diffusion_update_step2D(Output, Input, W_Lapl, lambdaPar, sigmaPar2, tau, dimX, dimY);
+            Diffusion_update_step2D(Output, Input, W_Lapl, lambdaPar, sigmaPar2, tau, (long)(dimX), (long)(dimY));
 		}
 	}
 	else {
@@ -67,7 +67,7 @@ float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sig
 		    /* Calculating weighted Laplacian */
             Weighted_Laplc3D(W_Lapl, Output, sigmaPar2, dimX, dimY, dimZ);
             /* Perform iteration step */
-            Diffusion_update_step3D(Output, Input, W_Lapl, lambdaPar, sigmaPar2, tau, dimX, dimY, dimZ);
+            Diffusion_update_step3D(Output, Input, W_Lapl, lambdaPar, sigmaPar2, tau, (long)(dimX), (long)(dimY), (long)(dimZ));
 		}
 	}
 	free(W_Lapl);
@@ -76,9 +76,9 @@ float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sig
 /********************************************************************/
 /***************************2D Functions*****************************/
 /********************************************************************/
-float Weighted_Laplc2D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY)
+float Weighted_Laplc2D(float *W_Lapl, float *U0, float sigma, long dimX, long dimY)
 {   
-    int i,j,i1,i2,j1,j2,index;
+    long i,j,i1,i2,j1,j2,index;
     float gradX, gradX_sq, gradY, gradY_sq, gradXX, gradYY, gradXY, xy_2, denom, V_norm, V_orth, c, c_sq;
 
         #pragma omp parallel for shared(W_Lapl) private(i,j,i1,i2,j1,j2,index,gradX, gradX_sq, gradY, gradY_sq, gradXX, gradYY, gradXY, xy_2, denom, V_norm, V_orth, c, c_sq)
@@ -125,9 +125,9 @@ float Weighted_Laplc2D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY
         return *W_Lapl;
 }
 
-float Diffusion_update_step2D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, int dimX, int dimY)
+float Diffusion_update_step2D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, long dimX, long dimY)
 {
-	int i,j,i1,i2,j1,j2,index;
+	long i,j,i1,i2,j1,j2,index;
     float gradXXc, gradYYc;
 
             #pragma omp parallel for shared(Output, Input, W_Lapl) private(i,j,i1,i2,j1,j2,index,gradXXc,gradYYc)
@@ -152,9 +152,9 @@ float Diffusion_update_step2D(float *Output, float *Input, float *W_Lapl, float
 /********************************************************************/
 /***************************3D Functions*****************************/
 /********************************************************************/
-float Weighted_Laplc3D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY, int dimZ)
+float Weighted_Laplc3D(float *W_Lapl, float *U0, float sigma, long dimX, long dimY, long dimZ)
 {   
-    int i,j,k,i1,i2,j1,j2,k1,k2,index;
+    long i,j,k,i1,i2,j1,j2,k1,k2,index;
     float gradX, gradX_sq, gradY, gradY_sq, gradXX, gradYY, gradXY, xy_2, denom, V_norm, V_orth, c, c_sq, gradZ, gradZ_sq, gradZZ, gradXZ, gradYZ, xyz_1, xyz_2;
         
         #pragma omp parallel for shared(W_Lapl) private(i,j,k,i1,i2,j1,j2,k1,k2,index,gradX, gradX_sq, gradY, gradY_sq, gradXX, gradYY, gradXY, xy_2, denom, V_norm, V_orth, c, c_sq, gradZ, gradZ_sq, gradZZ, gradXZ, gradYZ, xyz_1, xyz_2)
@@ -216,9 +216,9 @@ float Weighted_Laplc3D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY
         return *W_Lapl;
 }
 
-float Diffusion_update_step3D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, int dimX, int dimY, int dimZ)
+float Diffusion_update_step3D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, long dimX, long dimY, long dimZ)
 {
-	int i,j,i1,i2,j1,j2,index, k,k1,k2;
+	long i,j,i1,i2,j1,j2,index,k,k1,k2;
     float gradXXc, gradYYc, gradZZc;
 
         #pragma omp parallel for shared(Output, Input, W_Lapl) private(i,j,i1,i2,j1,j2,k,k1,k2,index,gradXXc,gradYYc,gradZZc)
diff --git a/Core/regularisers_CPU/Diffus4th_order_core.h b/Core/regularisers_CPU/Diffus4th_order_core.h
index 7b8fc08..d81afcb 100644
--- a/Core/regularisers_CPU/Diffus4th_order_core.h
+++ b/Core/regularisers_CPU/Diffus4th_order_core.h
@@ -46,10 +46,10 @@ limitations under the License.
 extern "C" {
 #endif
 CCPI_EXPORT float Diffus4th_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Weighted_Laplc2D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY);
-CCPI_EXPORT float Diffusion_update_step2D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, int dimX, int dimY);
-CCPI_EXPORT float Weighted_Laplc3D(float *W_Lapl, float *U0, float sigma, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Diffusion_update_step3D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float Weighted_Laplc2D(float *W_Lapl, float *U0, float sigma, long dimX, long dimY);
+CCPI_EXPORT float Diffusion_update_step2D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, long dimX, long dimY);
+CCPI_EXPORT float Weighted_Laplc3D(float *W_Lapl, float *U0, float sigma, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Diffusion_update_step3D(float *Output, float *Input, float *W_Lapl, float lambdaPar, float sigmaPar2, float tau, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/Diffusion_core.c b/Core/regularisers_CPU/Diffusion_core.c
index 51d0a57..b765796 100644
--- a/Core/regularisers_CPU/Diffusion_core.c
+++ b/Core/regularisers_CPU/Diffusion_core.c
@@ -55,20 +55,20 @@ float Diffusion_CPU_main(float *Input, float *Output, float lambdaPar, float sig
     sigmaPar2 = sigmaPar/sqrt(2.0f);
     
     /* copy into output */
-    copyIm(Input, Output, dimX, dimY, dimZ);
+    copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ));
     
     if (dimZ == 1) {
     /* running 2D diffusion iterations */
     for(i=0; i < iterationsNumb; i++) {
-            if (sigmaPar == 0.0f) LinearDiff2D(Input, Output, lambdaPar, tau, dimX, dimY); /* linear diffusion (heat equation) */
-            else NonLinearDiff2D(Input, Output, lambdaPar, sigmaPar2, tau, penaltytype, dimX, dimY); /* nonlinear diffusion */
+            if (sigmaPar == 0.0f) LinearDiff2D(Input, Output, lambdaPar, tau, (long)(dimX), (long)(dimY)); /* linear diffusion (heat equation) */
+            else NonLinearDiff2D(Input, Output, lambdaPar, sigmaPar2, tau, penaltytype, (long)(dimX), (long)(dimY)); /* nonlinear diffusion */
 		}
 	}
 	else {
 	/* running 3D diffusion iterations */
     for(i=0; i < iterationsNumb; i++) {
-            if (sigmaPar == 0.0f) LinearDiff3D(Input, Output, lambdaPar, tau, dimX, dimY, dimZ);
-            else NonLinearDiff3D(Input, Output, lambdaPar, sigmaPar2, tau, penaltytype, dimX, dimY, dimZ);
+            if (sigmaPar == 0.0f) LinearDiff3D(Input, Output, lambdaPar, tau, (long)(dimX), (long)(dimY), (long)(dimZ));
+            else NonLinearDiff3D(Input, Output, lambdaPar, sigmaPar2, tau, penaltytype, (long)(dimX), (long)(dimY), (long)(dimZ));
 		}
 	}
     return *Output;
@@ -79,9 +79,9 @@ float Diffusion_CPU_main(float *Input, float *Output, float lambdaPar, float sig
 /***************************2D Functions*****************************/
 /********************************************************************/
 /* linear diffusion (heat equation) */
-float LinearDiff2D(float *Input, float *Output, float lambdaPar, float tau, int dimX, int dimY)
+float LinearDiff2D(float *Input, float *Output, float lambdaPar, float tau, long dimX, long dimY)
 {
-	int i,j,i1,i2,j1,j2,index;
+	long i,j,i1,i2,j1,j2,index;
 	float e,w,n,s,e1,w1,n1,s1;
 	
 #pragma omp parallel for shared(Input) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1)
@@ -111,9 +111,9 @@ float LinearDiff2D(float *Input, float *Output, float lambdaPar, float tau, int
 }
 
 /* nonlinear diffusion */
-float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY)
+float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY)
 {
-	int i,j,i1,i2,j1,j2,index;
+	long i,j,i1,i2,j1,j2,index;
 	float e,w,n,s,e1,w1,n1,s1;
 	
 #pragma omp parallel for shared(Input) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1)
@@ -181,9 +181,9 @@ float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaP
 /***************************3D Functions*****************************/
 /********************************************************************/
 /* linear diffusion (heat equation) */
-float LinearDiff3D(float *Input, float *Output, float lambdaPar, float tau, int dimX, int dimY, int dimZ)
+float LinearDiff3D(float *Input, float *Output, float lambdaPar, float tau, long dimX, long dimY, long dimZ)
 {
-	int i,j,k,i1,i2,j1,j2,k1,k2,index;
+	long i,j,k,i1,i2,j1,j2,k1,k2,index;
 	float e,w,n,s,u,d,e1,w1,n1,s1,u1,d1;
 	
 #pragma omp parallel for shared(Input) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1,k,k1,k2,u1,d1,u,d)
@@ -219,9 +219,9 @@ for(k=0; k<dimZ; k++) {
 	return *Output;
 }
 
-float NonLinearDiff3D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY, int dimZ)
+float NonLinearDiff3D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY, long dimZ)
 {
-	int i,j,k,i1,i2,j1,j2,k1,k2,index;
+	long i,j,k,i1,i2,j1,j2,k1,k2,index;
 	float e,w,n,s,u,d,e1,w1,n1,s1,u1,d1;
 	
 #pragma omp parallel for shared(Input) private(index,i,j,i1,i2,j1,j2,e,w,n,s,e1,w1,n1,s1,k,k1,k2,u1,d1,u,d)
diff --git a/Core/regularisers_CPU/Diffusion_core.h b/Core/regularisers_CPU/Diffusion_core.h
index 0b4149a..cc36dad 100644
--- a/Core/regularisers_CPU/Diffusion_core.h
+++ b/Core/regularisers_CPU/Diffusion_core.h
@@ -50,10 +50,10 @@ limitations under the License.
 extern "C" {
 #endif
 CCPI_EXPORT float Diffusion_CPU_main(float *Input, float *Output, float lambdaPar, float sigmaPar, int iterationsNumb,  float tau, int penaltytype, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float LinearDiff2D(float *Input, float *Output, float lambdaPar, float tau, int dimX, int dimY);
-CCPI_EXPORT float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY);
-CCPI_EXPORT float LinearDiff3D(float *Input, float *Output, float lambdaPar, float tau, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float NonLinearDiff3D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float LinearDiff2D(float *Input, float *Output, float lambdaPar, float tau, long dimX, long dimY);
+CCPI_EXPORT float NonLinearDiff2D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY);
+CCPI_EXPORT float LinearDiff3D(float *Input, float *Output, float lambdaPar, float tau, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float NonLinearDiff3D(float *Input, float *Output, float lambdaPar, float sigmaPar, float tau, int penaltytype, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/FGP_TV_core.c b/Core/regularisers_CPU/FGP_TV_core.c
index 15bb45c..d828d48 100644
--- a/Core/regularisers_CPU/FGP_TV_core.c
+++ b/Core/regularisers_CPU/FGP_TV_core.c
@@ -39,7 +39,8 @@ limitations under the License.
  
 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)
 {
-	int ll, j, DimTotal;
+	int ll;
+    long j, DimTotal;
 	float re, re1;
 	float tk = 1.0f;
     float tkp1=1.0f;
@@ -48,7 +49,7 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
 	if (dimZ <= 1) {
 		/*2D case */
 		float *Output_prev=NULL, *P1=NULL, *P2=NULL, *P1_prev=NULL, *P2_prev=NULL, *R1=NULL, *R2=NULL;
-		DimTotal = dimX*dimY;
+		DimTotal = (long)(dimX*dimY);
 		
         Output_prev = calloc(DimTotal, sizeof(float));
         P1 = calloc(DimTotal, sizeof(float));
@@ -62,13 +63,13 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
         for(ll=0; ll<iterationsNumb; ll++) {
             
             /* computing the gradient of the objective function */
-            Obj_func2D(Input, Output, R1, R2, lambdaPar, dimX, dimY);
+            Obj_func2D(Input, Output, R1, R2, lambdaPar, (long)(dimY), (long)(dimZ));
             
             /* apply nonnegativity */
             if (nonneg == 1) for(j=0; j<DimTotal; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;}            
             
             /*Taking a step towards minus of the gradient*/
-            Grad_func2D(P1, P2, Output, R1, R2, lambdaPar, dimX, dimY);
+            Grad_func2D(P1, P2, Output, R1, R2, lambdaPar, (long)(dimY), (long)(dimZ));
             
             /* projection step */
             Proj_func2D(P1, P2, methodTV, DimTotal);
@@ -89,9 +90,9 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
 				if (count > 4) break;
             
             /*storing old values*/
-            copyIm(Output, Output_prev, dimX, dimY, 1);
-            copyIm(P1, P1_prev, dimX, dimY, 1);
-            copyIm(P2, P2_prev, dimX, dimY, 1);
+            copyIm(Output, Output_prev, (long)(dimY), (long)(dimZ), 1l);
+            copyIm(P1, P1_prev, (long)(dimY), (long)(dimZ), 1l);
+            copyIm(P2, P2_prev, (long)(dimY), (long)(dimZ), 1l);
             tk = tkp1;
         }
         if (printM == 1) printf("FGP-TV iterations stopped at iteration %i \n", ll);   
@@ -100,7 +101,7 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
 	else {
 		/*3D case*/
 		float *Output_prev=NULL, *P1=NULL, *P2=NULL, *P3=NULL, *P1_prev=NULL, *P2_prev=NULL, *P3_prev=NULL, *R1=NULL, *R2=NULL, *R3=NULL;		
-		DimTotal = dimX*dimY*dimZ;        
+		DimTotal = (long)(dimX*dimY*dimZ);        
         
         Output_prev = calloc(DimTotal, sizeof(float));
         P1 = calloc(DimTotal, sizeof(float));
@@ -117,13 +118,13 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
         for(ll=0; ll<iterationsNumb; ll++) {
             
             /* computing the gradient of the objective function */
-            Obj_func3D(Input, Output, R1, R2, R3, lambdaPar, dimX, dimY, dimZ);
+            Obj_func3D(Input, Output, R1, R2, R3, lambdaPar, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* apply nonnegativity */
             if (nonneg == 1) for(j=0; j<DimTotal; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;}  
             
             /*Taking a step towards minus of the gradient*/
-            Grad_func3D(P1, P2, P3, Output, R1, R2, R3, lambdaPar, dimX, dimY, dimZ);
+            Grad_func3D(P1, P2, P3, Output, R1, R2, R3, lambdaPar, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* projection step */
             Proj_func3D(P1, P2, P3, methodTV, DimTotal);
@@ -145,10 +146,10 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
             if (count > 4) break;            
                         
             /*storing old values*/
-            copyIm(Output, Output_prev, dimX, dimY, dimZ);
-            copyIm(P1, P1_prev, dimX, dimY, dimZ);
-            copyIm(P2, P2_prev, dimX, dimY, dimZ);
-            copyIm(P3, P3_prev, dimX, dimY, dimZ);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P1, P1_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P2, P2_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P3, P3_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
             tk = tkp1;            
         }	
 		if (printM == 1) printf("FGP-TV iterations stopped at iteration %i \n", ll);   
@@ -157,10 +158,10 @@ float TV_FGP_CPU_main(float *Input, float *Output, float lambdaPar, int iteratio
 	return *Output;
 }
 
-float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, int dimX, int dimY)
+float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, long dimX, long dimY)
 {
     float val1, val2;
-    int i,j,index;
+    long i,j,index;
 #pragma omp parallel for shared(A,D,R1,R2) private(index,i,j,val1,val2)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -172,10 +173,10 @@ float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, int dim
         }}
     return *D;
 }
-float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float lambda, int dimX, int dimY)
+float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float lambda,  long dimX, long dimY)
 {
     float val1, val2, multip;
-    int i,j,index;
+    long i,j,index;
     multip = (1.0f/(8.0f*lambda));
 #pragma omp parallel for shared(P1,P2,D,R1,R2,multip) private(index,i,j,val1,val2)
     for(i=0; i<dimX; i++) {
@@ -189,10 +190,10 @@ float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float la
         }}
     return 1;
 }
-float Proj_func2D(float *P1, float *P2, int methTV, int DimTotal)
+float Proj_func2D(float *P1, float *P2, int methTV, long DimTotal)
 {
     float val1, val2, denom, sq_denom;
-    int i;
+    long i;
     if (methTV == 0) {
         /* isotropic TV*/
 #pragma omp parallel for shared(P1,P2) private(i,denom,sq_denom)
@@ -219,9 +220,9 @@ float Proj_func2D(float *P1, float *P2, int methTV, int DimTotal)
     }
     return 1;
 }
-float Rupd_func2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, int DimTotal)
+float Rupd_func2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, long DimTotal)
 {
-    int i;
+    long i;
     float multip;
     multip = ((tk-1.0f)/tkp1);
 #pragma omp parallel for shared(P1,P2,P1_old,P2_old,R1,R2,multip) private(i)
@@ -234,10 +235,10 @@ float Rupd_func2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1,
 
 /* 3D-case related Functions */
 /*****************************************************************/
-float Obj_func3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ)
+float Obj_func3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ)
 {
     float val1, val2, val3;
-    int i,j,k,index;
+    long i,j,k,index;
 #pragma omp parallel for shared(A,D,R1,R2,R3) private(index,i,j,k,val1,val2,val3)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -251,10 +252,10 @@ float Obj_func3D(float *A, float *D, float *R1, float *R2, float *R3, float lamb
             }}}
     return *D;
 }
-float Grad_func3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ)
+float Grad_func3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ)
 {
     float val1, val2, val3, multip;
-    int i,j,k, index;
+    long i,j,k, index;
     multip = (1.0f/(26.0f*lambda));
 #pragma omp parallel for shared(P1,P2,P3,D,R1,R2,R3,multip) private(index,i,j,k,val1,val2,val3)
     for(i=0; i<dimX; i++) {
@@ -271,10 +272,10 @@ float Grad_func3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R
             }}}
     return 1;
 }
-float Proj_func3D(float *P1, float *P2, float *P3, int methTV, int DimTotal)
+float Proj_func3D(float *P1, float *P2, float *P3, int methTV, long DimTotal)
 {		
     float val1, val2, val3, denom, sq_denom;
-    int i;
+    long i;
     if (methTV == 0) {
 	/* isotropic TV*/
 	#pragma omp parallel for shared(P1,P2,P3) private(i,val1,val2,val3,sq_denom)
@@ -305,9 +306,9 @@ float Proj_func3D(float *P1, float *P2, float *P3, int methTV, int DimTotal)
 		}
     return 1;
 }
-float Rupd_func3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, int DimTotal)
+float Rupd_func3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, long DimTotal)
 {
-    int i;
+    long i;
     float multip;
     multip = ((tk-1.0f)/tkp1);
 #pragma omp parallel for shared(P1,P2,P3,P1_old,P2_old,P3_old,R1,R2,R3,multip) private(i)
diff --git a/Core/regularisers_CPU/FGP_TV_core.h b/Core/regularisers_CPU/FGP_TV_core.h
index 37e32f7..3418604 100644
--- a/Core/regularisers_CPU/FGP_TV_core.h
+++ b/Core/regularisers_CPU/FGP_TV_core.h
@@ -49,15 +49,15 @@ extern "C" {
 #endif
 CCPI_EXPORT 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);
 
-CCPI_EXPORT float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, int dimX, int dimY);
-CCPI_EXPORT float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float lambda, int dimX, int dimY);
-CCPI_EXPORT float Proj_func2D(float *P1, float *P2, int methTV, int DimTotal);
-CCPI_EXPORT float Rupd_func2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, int DimTotal);
+CCPI_EXPORT float Obj_func2D(float *A, float *D, float *R1, float *R2, float lambda, long dimX, long dimY);
+CCPI_EXPORT float Grad_func2D(float *P1, float *P2, float *D, float *R1, float *R2, float lambda, long dimX, long dimY);
+CCPI_EXPORT float Proj_func2D(float *P1, float *P2, int methTV, long DimTotal);
+CCPI_EXPORT float Rupd_func2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, long DimTotal);
 
-CCPI_EXPORT float Obj_func3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Grad_func3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Proj_func3D(float *P1, float *P2, float *P3, int methTV, int DimTotal);
-CCPI_EXPORT float Rupd_func3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, int DimTotal);
+CCPI_EXPORT float Obj_func3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Grad_func3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Proj_func3D(float *P1, float *P2, float *P3, int methTV, long DimTotal);
+CCPI_EXPORT float Rupd_func3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, long DimTotal);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/FGP_dTV_core.c b/Core/regularisers_CPU/FGP_dTV_core.c
index f6b4f79..17b75ff 100644
--- a/Core/regularisers_CPU/FGP_dTV_core.c
+++ b/Core/regularisers_CPU/FGP_dTV_core.c
@@ -44,7 +44,8 @@ limitations under the License.
  
 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)
 {
-	int ll, j, DimTotal;
+	int ll;
+    long j, DimTotal;
 	float re, re1;
 	float tk = 1.0f;
     float tkp1=1.0f;
@@ -53,7 +54,7 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
 	if (dimZ <= 1) {
 		/*2D case */
 		float *Output_prev=NULL, *P1=NULL, *P2=NULL, *P1_prev=NULL, *P2_prev=NULL, *R1=NULL, *R2=NULL, *InputRef_x=NULL, *InputRef_y=NULL;
-		DimTotal = dimX*dimY;
+		DimTotal = (long)(dimX*dimY);
 		
         Output_prev = calloc(DimTotal, sizeof(float));
         P1 = calloc(DimTotal, sizeof(float));
@@ -66,22 +67,22 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
         InputRef_y = calloc(DimTotal, sizeof(float)); 
 
 		/* calculate gradient field (smoothed) for the reference image */
-		GradNorm_func2D(InputRef, InputRef_x, InputRef_y, eta, dimX, dimY);
+		GradNorm_func2D(InputRef, InputRef_x, InputRef_y, eta, (long)(dimX), (long)(dimY));
 		
 		/* begin iterations */
         for(ll=0; ll<iterationsNumb; ll++) {
             
             /*projects a 2D vector field R-1,2 onto the orthogonal complement of another 2D vector field InputRef_xy*/                    
-            ProjectVect_func2D(R1, R2, InputRef_x, InputRef_y, dimX, dimY);
+            ProjectVect_func2D(R1, R2, InputRef_x, InputRef_y, (long)(dimX), (long)(dimY));
             
             /* computing the gradient of the objective function */
-            Obj_dfunc2D(Input, Output, R1, R2, lambdaPar, dimX, dimY);
+            Obj_dfunc2D(Input, Output, R1, R2, lambdaPar, (long)(dimX), (long)(dimY));
             
             /* apply nonnegativity */
             if (nonneg == 1) for(j=0; j<DimTotal; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;}
             
             /*Taking a step towards minus of the gradient*/
-            Grad_dfunc2D(P1, P2, Output, R1, R2, InputRef_x, InputRef_y, lambdaPar, dimX, dimY);
+            Grad_dfunc2D(P1, P2, Output, R1, R2, InputRef_x, InputRef_y, lambdaPar, (long)(dimX), (long)(dimY));
             
             /* projection step */
             Proj_dfunc2D(P1, P2, methodTV, DimTotal);
@@ -102,9 +103,9 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
 				if (count > 4) break;
             
             /*storing old values*/
-            copyIm(Output, Output_prev, dimX, dimY, 1);
-            copyIm(P1, P1_prev, dimX, dimY, 1);
-            copyIm(P2, P2_prev, dimX, dimY, 1);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), 1l);
+            copyIm(P1, P1_prev, (long)(dimX), (long)(dimY), 1l);
+            copyIm(P2, P2_prev, (long)(dimX), (long)(dimY), 1l);
             tk = tkp1;
         }
         if (printM == 1) printf("FGP-dTV iterations stopped at iteration %i \n", ll);   
@@ -113,7 +114,7 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
 	else {
 		/*3D case*/
 		float *Output_prev=NULL, *P1=NULL, *P2=NULL, *P3=NULL, *P1_prev=NULL, *P2_prev=NULL, *P3_prev=NULL, *R1=NULL, *R2=NULL, *R3=NULL, *InputRef_x=NULL, *InputRef_y=NULL, *InputRef_z=NULL; 
-		DimTotal = dimX*dimY*dimZ;
+		DimTotal = (long)(dimX*dimY*dimZ);
         
         Output_prev = calloc(DimTotal, sizeof(float));
         P1 = calloc(DimTotal, sizeof(float));
@@ -130,22 +131,22 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
         InputRef_z = calloc(DimTotal, sizeof(float)); 
 
 		/* calculate gradient field (smoothed) for the reference volume */
-		GradNorm_func3D(InputRef, InputRef_x, InputRef_y, InputRef_z, eta, dimX, dimY, dimZ);
+		GradNorm_func3D(InputRef, InputRef_x, InputRef_y, InputRef_z, eta, (long)(dimX), (long)(dimY), (long)(dimZ));
 		
 		/* begin iterations */
         for(ll=0; ll<iterationsNumb; ll++) {
 
 			 /*projects a 3D vector field R-1,2,3 onto the orthogonal complement of another 3D vector field InputRef_xyz*/
-            ProjectVect_func3D(R1, R2, R3, InputRef_x, InputRef_y, InputRef_z, dimX, dimY, dimZ);
+            ProjectVect_func3D(R1, R2, R3, InputRef_x, InputRef_y, InputRef_z, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* computing the gradient of the objective function */
-            Obj_dfunc3D(Input, Output, R1, R2, R3, lambdaPar, dimX, dimY, dimZ);
+            Obj_dfunc3D(Input, Output, R1, R2, R3, lambdaPar, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* apply nonnegativity */
             if (nonneg == 1) for(j=0; j<DimTotal; j++) {if (Output[j] < 0.0f) Output[j] = 0.0f;}  
             
             /*Taking a step towards minus of the gradient*/
-            Grad_dfunc3D(P1, P2, P3, Output, R1, R2, R3, InputRef_x, InputRef_y, InputRef_z, lambdaPar, dimX, dimY, dimZ);
+            Grad_dfunc3D(P1, P2, P3, Output, R1, R2, R3, InputRef_x, InputRef_y, InputRef_z, lambdaPar, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* projection step */
             Proj_dfunc3D(P1, P2, P3, methodTV, DimTotal);
@@ -167,10 +168,10 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
             if (count > 4) break;            
                         
             /*storing old values*/
-            copyIm(Output, Output_prev, dimX, dimY, dimZ);
-            copyIm(P1, P1_prev, dimX, dimY, dimZ);
-            copyIm(P2, P2_prev, dimX, dimY, dimZ);
-            copyIm(P3, P3_prev, dimX, dimY, dimZ);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P1, P1_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P2, P2_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
+            copyIm(P3, P3_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
             tk = tkp1;            
         }	
 		if (printM == 1) printf("FGP-dTV iterations stopped at iteration %i \n", ll);   
@@ -184,9 +185,9 @@ float dTV_FGP_CPU_main(float *Input, float *InputRef, float *Output, float lambd
 /***************************2D Functions*****************************/
 /********************************************************************/
 
-float GradNorm_func2D(float *B, float *B_x, float *B_y, float eta, int dimX, int dimY)
+float GradNorm_func2D(float *B, float *B_x, float *B_y, float eta, long dimX, long dimY)
 {
-    int i,j,index;
+    long i,j,index;
     float val1, val2, gradX, gradY, magn;
 #pragma omp parallel for shared(B, B_x, B_y) private(i,j,index,val1,val2,gradX,gradY,magn)
     for(i=0; i<dimX; i++) {
@@ -205,9 +206,9 @@ float GradNorm_func2D(float *B, float *B_x, float *B_y, float eta, int dimX, int
     return 1;
 }
 
-float ProjectVect_func2D(float *R1, float *R2, float *B_x, float *B_y, int dimX, int dimY)
+float ProjectVect_func2D(float *R1, float *R2, float *B_x, float *B_y, long dimX, long dimY)
 {
-    int i,j,index;
+    long i,j,index;
     float in_prod;
 #pragma omp parallel for shared(R1, R2, B_x, B_y) private(index,i,j,in_prod)
     for(i=0; i<dimX; i++) {
@@ -220,10 +221,10 @@ float ProjectVect_func2D(float *R1, float *R2, float *B_x, float *B_y, int dimX,
     return 1;
 }
 
-float Obj_dfunc2D(float *A, float *D, float *R1, float *R2, float lambda, int dimX, int dimY)
+float Obj_dfunc2D(float *A, float *D, float *R1, float *R2, float lambda, long dimX, long dimY)
 {
     float val1, val2;
-    int i,j,index;
+    long i,j,index;
 #pragma omp parallel for shared(A,D,R1,R2) private(index,i,j,val1,val2)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -235,10 +236,10 @@ float Obj_dfunc2D(float *A, float *D, float *R1, float *R2, float lambda, int di
         }}
     return *D;
 }
-float Grad_dfunc2D(float *P1, float *P2, float *D, float *R1, float *R2, float *B_x, float *B_y, float lambda, int dimX, int dimY)
+float Grad_dfunc2D(float *P1, float *P2, float *D, float *R1, float *R2, float *B_x, float *B_y, float lambda, long dimX, long dimY)
 {
     float val1, val2, multip, in_prod;
-    int i,j,index;
+    long i,j,index;
     multip = (1.0f/(8.0f*lambda));
 #pragma omp parallel for shared(P1,P2,D,R1,R2,B_x,B_y,multip) private(i,j,index,val1,val2,in_prod)
     for(i=0; i<dimX; i++) {
@@ -258,10 +259,10 @@ float Grad_dfunc2D(float *P1, float *P2, float *D, float *R1, float *R2, float *
         }}
     return 1;
 }
-float Proj_dfunc2D(float *P1, float *P2, int methTV, int DimTotal)
+float Proj_dfunc2D(float *P1, float *P2, int methTV, long DimTotal)
 {
     float val1, val2, denom, sq_denom;
-    int i;
+    long i;
     if (methTV == 0) {
         /* isotropic TV*/
 #pragma omp parallel for shared(P1,P2) private(i,denom,sq_denom)
@@ -288,9 +289,9 @@ float Proj_dfunc2D(float *P1, float *P2, int methTV, int DimTotal)
     }
     return 1;
 }
-float Rupd_dfunc2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, int DimTotal)
+float Rupd_dfunc2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, long DimTotal)
 {
-    int i;
+    long i;
     float multip;
     multip = ((tk-1.0f)/tkp1);
 #pragma omp parallel for shared(P1,P2,P1_old,P2_old,R1,R2,multip) private(i)
@@ -304,9 +305,9 @@ float Rupd_dfunc2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1
 /********************************************************************/
 /***************************3D Functions*****************************/
 /********************************************************************/
-float GradNorm_func3D(float *B, float *B_x, float *B_y, float *B_z, float eta, int dimX, int dimY, int dimZ)
+float GradNorm_func3D(float *B, float *B_x, float *B_y, float *B_z, float eta, long dimX, long dimY, long dimZ)
 {
-    int i, j, k, index;
+    long i, j, k, index;
     float val1, val2, val3, gradX, gradY, gradZ, magn;
 #pragma omp parallel for shared(B, B_x, B_y, B_z) private(i,j,k,index,val1,val2,val3,gradX,gradY,gradZ,magn)
     for(i=0; i<dimX; i++) {
@@ -331,9 +332,9 @@ float GradNorm_func3D(float *B, float *B_x, float *B_y, float *B_z, float eta, i
     return 1;
 }
 
-float ProjectVect_func3D(float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, int dimX, int dimY, int dimZ)
+float ProjectVect_func3D(float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, long dimX, long dimY, long dimZ)
 {
-    int i,j,k,index;
+    long i,j,k,index;
     float in_prod;
 #pragma omp parallel for shared(R1, R2, R3, B_x, B_y, B_z) private(index,i,j,k,in_prod)
     for(i=0; i<dimX; i++) {
@@ -348,10 +349,10 @@ float ProjectVect_func3D(float *R1, float *R2, float *R3, float *B_x, float *B_y
     return 1;
 }
 
-float Obj_dfunc3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ)
+float Obj_dfunc3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ)
 {
     float val1, val2, val3;
-    int i,j,k,index;
+    long i,j,k,index;
 #pragma omp parallel for shared(A,D,R1,R2,R3) private(index,i,j,k,val1,val2,val3)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -365,10 +366,10 @@ float Obj_dfunc3D(float *A, float *D, float *R1, float *R2, float *R3, float lam
             }}}
     return *D;
 }
-float Grad_dfunc3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, float lambda, int dimX, int dimY, int dimZ)
+float Grad_dfunc3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, float lambda, long dimX, long dimY, long dimZ)
 {
     float val1, val2, val3, multip, in_prod;
-    int i,j,k, index;
+    long i,j,k, index;
     multip = (1.0f/(26.0f*lambda));
 #pragma omp parallel for shared(P1,P2,P3,D,R1,R2,R3,multip) private(index,i,j,k,val1,val2,val3,in_prod)
     for(i=0; i<dimX; i++) {
@@ -391,10 +392,10 @@ float Grad_dfunc3D(float *P1, float *P2, float *P3, float *D, float *R1, float *
             }}}
     return 1;
 }
-float Proj_dfunc3D(float *P1, float *P2, float *P3, int methTV, int DimTotal)
+float Proj_dfunc3D(float *P1, float *P2, float *P3, int methTV, long DimTotal)
 {		
     float val1, val2, val3, denom, sq_denom;
-    int i;
+    long i;
     if (methTV == 0) {
 	/* isotropic TV*/
 	#pragma omp parallel for shared(P1,P2,P3) private(i,val1,val2,val3,sq_denom)
@@ -425,9 +426,9 @@ float Proj_dfunc3D(float *P1, float *P2, float *P3, int methTV, int DimTotal)
 		}
     return 1;
 }
-float Rupd_dfunc3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, int DimTotal)
+float Rupd_dfunc3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, long DimTotal)
 {
-    int i;
+    long i;
     float multip;
     multip = ((tk-1.0f)/tkp1);
 #pragma omp parallel for shared(P1,P2,P3,P1_old,P2_old,P3_old,R1,R2,R3,multip) private(i)
diff --git a/Core/regularisers_CPU/FGP_dTV_core.h b/Core/regularisers_CPU/FGP_dTV_core.h
index 8d721d5..442dd30 100644
--- a/Core/regularisers_CPU/FGP_dTV_core.h
+++ b/Core/regularisers_CPU/FGP_dTV_core.h
@@ -54,19 +54,19 @@ extern "C" {
 #endif
 CCPI_EXPORT 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);
 
-CCPI_EXPORT float GradNorm_func2D(float *B, float *B_x, float *B_y, float eta, int dimX, int dimY);
-CCPI_EXPORT float ProjectVect_func2D(float *R1, float *R2, float *B_x, float *B_y, int dimX, int dimY);
-CCPI_EXPORT float Obj_dfunc2D(float *A, float *D, float *R1, float *R2, float lambda, int dimX, int dimY);
-CCPI_EXPORT float Grad_dfunc2D(float *P1, float *P2, float *D, float *R1, float *R2, float *B_x, float *B_y, float lambda, int dimX, int dimY);
-CCPI_EXPORT float Proj_dfunc2D(float *P1, float *P2, int methTV, int DimTotal);
-CCPI_EXPORT float Rupd_dfunc2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, int DimTotal);
+CCPI_EXPORT float GradNorm_func2D(float *B, float *B_x, float *B_y, float eta, long dimX, long dimY);
+CCPI_EXPORT float ProjectVect_func2D(float *R1, float *R2, float *B_x, float *B_y, long dimX, long dimY);
+CCPI_EXPORT float Obj_dfunc2D(float *A, float *D, float *R1, float *R2, float lambda, long dimX, long dimY);
+CCPI_EXPORT float Grad_dfunc2D(float *P1, float *P2, float *D, float *R1, float *R2, float *B_x, float *B_y, float lambda, long dimX, long dimY);
+CCPI_EXPORT float Proj_dfunc2D(float *P1, float *P2, int methTV, long DimTotal);
+CCPI_EXPORT float Rupd_dfunc2D(float *P1, float *P1_old, float *P2, float *P2_old, float *R1, float *R2, float tkp1, float tk, long DimTotal);
 
-CCPI_EXPORT float GradNorm_func3D(float *B, float *B_x, float *B_y, float *B_z, float eta, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float ProjectVect_func3D(float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Obj_dfunc3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Grad_dfunc3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, float lambda, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Proj_dfunc3D(float *P1, float *P2, float *P3, int methTV, int DimTotal);
-CCPI_EXPORT float Rupd_dfunc3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, int DimTotal);
+CCPI_EXPORT float GradNorm_func3D(float *B, float *B_x, float *B_y, float *B_z, float eta, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float ProjectVect_func3D(float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Obj_dfunc3D(float *A, float *D, float *R1, float *R2, float *R3, float lambda, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Grad_dfunc3D(float *P1, float *P2, float *P3, float *D, float *R1, float *R2, float *R3, float *B_x, float *B_y, float *B_z, float lambda, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Proj_dfunc3D(float *P1, float *P2, float *P3, int methTV, long DimTotal);
+CCPI_EXPORT float Rupd_dfunc3D(float *P1, float *P1_old, float *P2, float *P2_old, float *P3, float *P3_old, float *R1, float *R2, float *R3, float tkp1, float tk, long DimTotal);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/LLT_ROF_core.c b/Core/regularisers_CPU/LLT_ROF_core.c
index 1584a29..8416a14 100644
--- a/Core/regularisers_CPU/LLT_ROF_core.c
+++ b/Core/regularisers_CPU/LLT_ROF_core.c
@@ -51,10 +51,11 @@ int signLLT(float x) {
 
 float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int dimX, int dimY, int dimZ)
 {
-		int DimTotal, ll;
+		long DimTotal;
+        int ll;
 		float *D1_LLT=NULL, *D2_LLT=NULL, *D3_LLT=NULL, *D1_ROF=NULL, *D2_ROF=NULL, *D3_ROF=NULL;
 		
-		DimTotal = dimX*dimY*dimZ;
+		DimTotal = (long)(dimX*dimY*dimZ);
         
         D1_ROF = calloc(DimTotal, sizeof(float));
         D2_ROF = calloc(DimTotal, sizeof(float));
@@ -64,32 +65,32 @@ float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambd
         D2_LLT = calloc(DimTotal, sizeof(float));
         D3_LLT = calloc(DimTotal, sizeof(float));
         
-        copyIm(Input, Output, dimX, dimY, dimZ); /* initialize  */
+        copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ)); /* initialize  */
        
 		for(ll = 0; ll < iterationsNumb; ll++) {            
             if (dimZ == 1) {
 			/* 2D case */
 			/****************ROF******************/
 			 /* calculate first-order differences */
-            D1_func_ROF(Output, D1_ROF, dimX, dimY, dimZ);
-            D2_func_ROF(Output, D2_ROF, dimX, dimY, dimZ);
+            D1_func_ROF(Output, D1_ROF, (long)(dimX), (long)(dimY), 1l);
+            D2_func_ROF(Output, D2_ROF, (long)(dimX), (long)(dimY), 1l);
             /****************LLT******************/
             /* estimate second-order derrivatives */
-            der2D_LLT(Output, D1_LLT, D2_LLT, dimX, dimY, dimZ);
+            der2D_LLT(Output, D1_LLT, D2_LLT, (long)(dimX), (long)(dimY), 1l);
             /* Joint update for ROF and LLT models */
-            Update2D_LLT_ROF(Input, Output, D1_LLT, D2_LLT, D1_ROF, D2_ROF, lambdaROF, lambdaLLT, tau, dimX, dimY, dimZ);
+            Update2D_LLT_ROF(Input, Output, D1_LLT, D2_LLT, D1_ROF, D2_ROF, lambdaROF, lambdaLLT, tau, (long)(dimX), (long)(dimY), 1l);
             }
             else {
 			/* 3D case */
 			/* calculate first-order differences */
-            D1_func_ROF(Output, D1_ROF, dimX, dimY, dimZ);
-            D2_func_ROF(Output, D2_ROF, dimX, dimY, dimZ);
-            D3_func_ROF(Output, D3_ROF, dimX, dimY, dimZ); 
+            D1_func_ROF(Output, D1_ROF, (long)(dimX), (long)(dimY), (long)(dimZ));
+            D2_func_ROF(Output, D2_ROF, (long)(dimX), (long)(dimY), (long)(dimZ));
+            D3_func_ROF(Output, D3_ROF, (long)(dimX), (long)(dimY), (long)(dimZ)); 
             /****************LLT******************/
             /* estimate second-order derrivatives */
-            der3D_LLT(Output, D1_LLT, D2_LLT, D3_LLT, dimX, dimY, dimZ);
+            der3D_LLT(Output, D1_LLT, D2_LLT, D3_LLT,(long)(dimX), (long)(dimY), (long)(dimZ));
             /* Joint update for ROF and LLT models */
-            Update3D_LLT_ROF(Input, Output, D1_LLT, D2_LLT, D3_LLT, D1_ROF, D2_ROF, D3_ROF, lambdaROF, lambdaLLT, tau, dimX, dimY, dimZ);
+            Update3D_LLT_ROF(Input, Output, D1_LLT, D2_LLT, D3_LLT, D1_ROF, D2_ROF, D3_ROF, lambdaROF, lambdaLLT, tau, (long)(dimX), (long)(dimY), (long)(dimZ));
 			}
         } /*end of iterations*/
     free(D1_LLT);free(D2_LLT);free(D3_LLT);
@@ -100,9 +101,9 @@ float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambd
 /*************************************************************************/
 /**********************LLT-related functions *****************************/
 /*************************************************************************/
-float der2D_LLT(float *U, float *D1, float *D2, int dimX, int dimY, int dimZ)
+float der2D_LLT(float *U, float *D1, float *D2, long dimX, long dimY, long dimZ)
 {
-	int i, j, index, i_p, i_m, j_m, j_p;
+	long i, j, index, i_p, i_m, j_m, j_p;
 	float dxx, dyy, denom_xx, denom_yy;
 #pragma omp parallel for shared(U,D1,D2) private(i, j, index, i_p, i_m, j_m, j_p, denom_xx, denom_yy, dxx, dyy)
 	for (i = 0; i<dimX; i++) {
@@ -127,9 +128,9 @@ float der2D_LLT(float *U, float *D1, float *D2, int dimX, int dimY, int dimZ)
 	return 1;
 }
 
-float der3D_LLT(float *U, float *D1, float *D2, float *D3, int dimX, int dimY, int dimZ)
+float der3D_LLT(float *U, float *D1, float *D2, float *D3, long dimX, long dimY, long dimZ)
  {
- 	int i, j, k, i_p, i_m, j_m, j_p, k_p, k_m, index;
+ 	long i, j, k, i_p, i_m, j_m, j_p, k_p, k_m, index;
  	float dxx, dyy, dzz, denom_xx, denom_yy, denom_zz;
  #pragma omp parallel for shared(U,D1,D2,D3) private(i, j, index, k, i_p, i_m, j_m, j_p, k_p, k_m, denom_xx, denom_yy, denom_zz, dxx, dyy, dzz)
  	for (i = 0; i<dimX; i++) {
@@ -167,10 +168,10 @@ float der3D_LLT(float *U, float *D1, float *D2, float *D3, int dimX, int dimY, i
 /*************************************************************************/
 
 /* calculate differences 1 */
-float D1_func_ROF(float *A, float *D1, int dimX, int dimY, int dimZ)
+float D1_func_ROF(float *A, float *D1, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMy_0, NOMz_1, NOMz_0, denom1, denom2,denom3, T1;
-    int i,j,k,i1,i2,k1,j1,j2,k2,index;
+    long i,j,k,i1,i2,k1,j1,j2,k2,index;
     
     if (dimZ > 1) {
 #pragma omp parallel for shared (A, D1, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2, NOMx_1,NOMy_1,NOMy_0,NOMz_1,NOMz_0,denom1,denom2,denom3,T1)
@@ -232,10 +233,10 @@ float D1_func_ROF(float *A, float *D1, int dimX, int dimY, int dimZ)
     return *D1;
 }
 /* calculate differences 2 */
-float D2_func_ROF(float *A, float *D2, int dimX, int dimY, int dimZ)
+float D2_func_ROF(float *A, float *D2, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMx_0, NOMz_1, NOMz_0, denom1, denom2, denom3, T2;
-    int i,j,k,i1,i2,k1,j1,j2,k2,index;
+    long i,j,k,i1,i2,k1,j1,j2,k2,index;
     
     if (dimZ > 1) {
 #pragma omp parallel for shared (A, D2, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2,  NOMx_1, NOMy_1, NOMx_0, NOMz_1, NOMz_0, denom1, denom2, denom3, T2)
@@ -297,10 +298,10 @@ float D2_func_ROF(float *A, float *D2, int dimX, int dimY, int dimZ)
 }
 
 /* calculate differences 3 */
-float D3_func_ROF(float *A, float *D3, int dimY, int dimX, int dimZ)
+float D3_func_ROF(float *A, float *D3, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMx_0, NOMy_0, NOMz_1, denom1, denom2, denom3, T3;
-    int index,i,j,k,i1,i2,k1,j1,j2,k2;
+    long index,i,j,k,i1,i2,k1,j1,j2,k2;
     
 #pragma omp parallel for shared (A, D3, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2,  NOMx_1, NOMy_1, NOMy_0, NOMx_0, NOMz_1, denom1, denom2, denom3, T3)
     for(j=0; j<dimY; j++) {
@@ -338,9 +339,9 @@ float D3_func_ROF(float *A, float *D3, int dimY, int dimX, int dimZ)
 /**********************ROF-LLT-related functions *************************/
 /*************************************************************************/
 
-float Update2D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D1_ROF, float *D2_ROF, float lambdaROF, float lambdaLLT, float tau, int dimX, int dimY, int dimZ)
+float Update2D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D1_ROF, float *D2_ROF, float lambdaROF, float lambdaLLT, float tau, long dimX, long dimY, long dimZ)
 {
-	int i, j, index, i_p, i_m, j_m, j_p;
+	long i, j, index, i_p, i_m, j_m, j_p;
 	float div, laplc, dxx, dyy, dv1, dv2;
 #pragma omp parallel for shared(U,U0) private(i, j, index, i_p, i_m, j_m, j_p, laplc, div, dxx, dyy, dv1, dv2)
 	for (i = 0; i<dimX; i++) {
@@ -369,9 +370,9 @@ float Update2D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float
 	return *U;
 }
 
-float Update3D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D3_LLT, float *D1_ROF, float *D2_ROF, float *D3_ROF, float lambdaROF, float lambdaLLT, float tau, int dimX, int dimY, int dimZ)
+float Update3D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D3_LLT, float *D1_ROF, float *D2_ROF, float *D3_ROF, float lambdaROF, float lambdaLLT, float tau, long dimX, long dimY, long dimZ)
 {
-	int i, j, k, i_p, i_m, j_m, j_p, k_p, k_m, index;
+	long i, j, k, i_p, i_m, j_m, j_p, k_p, k_m, index;
 	float div, laplc, dxx, dyy, dzz, dv1, dv2, dv3;
 #pragma omp parallel for shared(U,U0) private(i, j, k, index, i_p, i_m, j_m, j_p, k_p, k_m, laplc, div, dxx, dyy, dzz, dv1, dv2, dv3)
  	for (i = 0; i<dimX; i++) {
diff --git a/Core/regularisers_CPU/LLT_ROF_core.h b/Core/regularisers_CPU/LLT_ROF_core.h
index 65d25cd..8e6591e 100644
--- a/Core/regularisers_CPU/LLT_ROF_core.h
+++ b/Core/regularisers_CPU/LLT_ROF_core.h
@@ -51,15 +51,15 @@ extern "C" {
 #endif
 CCPI_EXPORT float LLT_ROF_CPU_main(float *Input, float *Output, float lambdaROF, float lambdaLLT, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
 
-CCPI_EXPORT float der2D_LLT(float *U, float *D1, float *D2, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float der3D_LLT(float *U, float *D1, float *D2, float *D3, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float der2D_LLT(float *U, float *D1, float *D2, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float der3D_LLT(float *U, float *D1, float *D2, float *D3, long dimX, long dimY, long dimZ);
 
-CCPI_EXPORT float D1_func_ROF(float *A, float *D1, int dimY, int dimX, int dimZ);
-CCPI_EXPORT float D2_func_ROF(float *A, float *D2, int dimY, int dimX, int dimZ);
-CCPI_EXPORT float D3_func_ROF(float *A, float *D3, int dimY, int dimX, int dimZ);
+CCPI_EXPORT float D1_func_ROF(float *A, float *D1, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float D2_func_ROF(float *A, float *D2, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float D3_func_ROF(float *A, float *D3, long dimX, long dimY, long dimZ);
 
-CCPI_EXPORT float Update2D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D1_ROF, float *D2_ROF, float lambdaROF, float lambdaLLT, float tau, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float Update3D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D3_LLT, float *D1_ROF, float *D2_ROF, float *D3_ROF, float lambdaROF, float lambdaLLT, float tau, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float Update2D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D1_ROF, float *D2_ROF, float lambdaROF, float lambdaLLT, float tau, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float Update3D_LLT_ROF(float *U0, float *U, float *D1_LLT, float *D2_LLT, float *D3_LLT, float *D1_ROF, float *D2_ROF, float *D3_ROF, float lambdaROF, float lambdaLLT, float tau, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/ROF_TV_core.c b/Core/regularisers_CPU/ROF_TV_core.c
index fb3bc7c..1858442 100644
--- a/Core/regularisers_CPU/ROF_TV_core.c
+++ b/Core/regularisers_CPU/ROF_TV_core.c
@@ -49,34 +49,34 @@ int sign(float x) {
 float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ)
 {
     float *D1, *D2, *D3;
-    int i, DimTotal;
-    DimTotal = dimX*dimY*dimZ;    
+    int i; 
+    long DimTotal;
+    DimTotal = (long)(dimX*dimY*dimZ);    
     
     D1 = calloc(DimTotal, sizeof(float));
     D2 = calloc(DimTotal, sizeof(float));
     D3 = calloc(DimTotal, sizeof(float));
 	   
     /* copy into output */
-    copyIm(Input, Output, dimX, dimY, dimZ);
+    copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ));
         
     /* start TV iterations */
-    for(i=0; i < iterationsNumb; i++) {
-            
+    for(i=0; i < iterationsNumb; i++) {            
             /* calculate differences */
-            D1_func(Output, D1, dimX, dimY, dimZ);
-            D2_func(Output, D2, dimX, dimY, dimZ);
-            if (dimZ > 1) D3_func(Output, D3, dimX, dimY, dimZ); 
-            TV_kernel(D1, D2, D3, Output, Input, lambdaPar, tau, dimX, dimY, dimZ);
+            D1_func(Output, D1, (long)(dimX), (long)(dimY), (long)(dimZ));
+            D2_func(Output, D2, (long)(dimX), (long)(dimY), (long)(dimZ));
+            if (dimZ > 1) D3_func(Output, D3, (long)(dimX), (long)(dimY), (long)(dimZ)); 
+            TV_kernel(D1, D2, D3, Output, Input, lambdaPar, tau, (long)(dimX), (long)(dimY), (long)(dimZ));
 		}           
     free(D1);free(D2); free(D3);
     return *Output;
 }
 
 /* calculate differences 1 */
-float D1_func(float *A, float *D1, int dimX, int dimY, int dimZ)
+float D1_func(float *A, float *D1, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMy_0, NOMz_1, NOMz_0, denom1, denom2,denom3, T1;
-    int i,j,k,i1,i2,k1,j1,j2,k2,index;
+    long i,j,k,i1,i2,k1,j1,j2,k2,index;
     
     if (dimZ > 1) {
 #pragma omp parallel for shared (A, D1, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2, NOMx_1,NOMy_1,NOMy_0,NOMz_1,NOMz_0,denom1,denom2,denom3,T1)
@@ -138,10 +138,10 @@ float D1_func(float *A, float *D1, int dimX, int dimY, int dimZ)
     return *D1;
 }
 /* calculate differences 2 */
-float D2_func(float *A, float *D2, int dimX, int dimY, int dimZ)
+float D2_func(float *A, float *D2, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMx_0, NOMz_1, NOMz_0, denom1, denom2, denom3, T2;
-    int i,j,k,i1,i2,k1,j1,j2,k2,index;
+    long i,j,k,i1,i2,k1,j1,j2,k2,index;
     
     if (dimZ > 1) {
 #pragma omp parallel for shared (A, D2, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2,  NOMx_1, NOMy_1, NOMx_0, NOMz_1, NOMz_0, denom1, denom2, denom3, T2)
@@ -155,8 +155,7 @@ float D2_func(float *A, float *D2, int dimX, int dimY, int dimZ)
                     j1 = j + 1; if (j1 >= dimY) j1 = j-1;
                     j2 = j - 1; if (j2 < 0) j2 = j+1;
                     k1 = k + 1; if (k1 >= dimZ) k1 = k-1;
-                    k2 = k - 1; if (k2 < 0) k2 = k+1;
-                    
+                    k2 = k - 1; if (k2 < 0) k2 = k+1;                    
                     
                     /* Forward-backward differences */
                     NOMx_1 = A[(dimX*dimY)*k + (j1)*dimX + i] - A[index]; /* x+ */
@@ -203,10 +202,10 @@ float D2_func(float *A, float *D2, int dimX, int dimY, int dimZ)
 }
 
 /* calculate differences 3 */
-float D3_func(float *A, float *D3, int dimY, int dimX, int dimZ)
+float D3_func(float *A, float *D3, long dimX, long dimY, long dimZ)
 {
     float NOMx_1, NOMy_1, NOMx_0, NOMy_0, NOMz_1, denom1, denom2, denom3, T3;
-    int index,i,j,k,i1,i2,k1,j1,j2,k2;
+    long index,i,j,k,i1,i2,k1,j1,j2,k2;
     
 #pragma omp parallel for shared (A, D3, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2,  NOMx_1, NOMy_1, NOMy_0, NOMx_0, NOMz_1, denom1, denom2, denom3, T3)
     for(j=0; j<dimY; j++) {
@@ -241,10 +240,10 @@ float D3_func(float *A, float *D3, int dimY, int dimX, int dimZ)
 }
 
 /* calculate divergence */
-float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambda, float tau, int dimX, int dimY, int dimZ)
+float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambda, float tau, long dimX, long dimY, long dimZ)
 {
     float dv1, dv2, dv3;
-    int index,i,j,k,i1,i2,k1,j1,j2,k2;
+    long index,i,j,k,i1,i2,k1,j1,j2,k2;
     
     if (dimZ > 1) {
 #pragma omp parallel for shared (D1, D2, D3, B, dimX, dimY, dimZ) private(index, i, j, k, i1, j1, k1, i2, j2, k2, dv1,dv2,dv3)
diff --git a/Core/regularisers_CPU/ROF_TV_core.h b/Core/regularisers_CPU/ROF_TV_core.h
index 14daf58..4e320e9 100644
--- a/Core/regularisers_CPU/ROF_TV_core.h
+++ b/Core/regularisers_CPU/ROF_TV_core.h
@@ -46,11 +46,12 @@ limitations under the License.
 #ifdef __cplusplus
 extern "C" {
 #endif
-CCPI_EXPORT float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambda, float tau, int dimY, int dimX, int dimZ);
 CCPI_EXPORT float TV_ROF_CPU_main(float *Input, float *Output, float lambdaPar, int iterationsNumb, float tau, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float D1_func(float *A, float *D1, int dimY, int dimX, int dimZ);
-CCPI_EXPORT float D2_func(float *A, float *D2, int dimY, int dimX, int dimZ);
-CCPI_EXPORT float D3_func(float *A, float *D3, int dimY, int dimX, int dimZ);
+
+CCPI_EXPORT float TV_kernel(float *D1, float *D2, float *D3, float *B, float *A, float lambda, float tau, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float D1_func(float *A, float *D1, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float D2_func(float *A, float *D2, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float D3_func(float *A, float *D3, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
-#endif
+#endif
\ No newline at end of file
diff --git a/Core/regularisers_CPU/SB_TV_core.c b/Core/regularisers_CPU/SB_TV_core.c
index 06d3de3..769ea67 100755
--- a/Core/regularisers_CPU/SB_TV_core.c
+++ b/Core/regularisers_CPU/SB_TV_core.c
@@ -37,16 +37,17 @@ limitations under the License.
  
 float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsil, int methodTV, int printM, int dimX, int dimY, int dimZ)
 {
-	int ll, j, DimTotal;
+	int ll;
+    long j, DimTotal;    
 	float re, re1, lambda;
     int count = 0;
     mu = 1.0f/mu;
-	lambda = 2.0f*mu;
+    lambda = 2.0f*mu;
 
 	if (dimZ <= 1) {
 		/* 2D case */
 		float *Output_prev=NULL, *Dx=NULL, *Dy=NULL, *Bx=NULL, *By=NULL;
-		DimTotal = dimX*dimY;
+		DimTotal = (long)(dimX*dimY);
 		
 		Output_prev = calloc(DimTotal, sizeof(float));
 		Dx = calloc(DimTotal, sizeof(float));
@@ -54,26 +55,26 @@ float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsi
 		Bx = calloc(DimTotal, sizeof(float));
 		By = calloc(DimTotal, sizeof(float));
         
-        copyIm(Input, Output, dimX, dimY, 1); /*initialize */
+        copyIm(Input, Output, (long)(dimX), (long)(dimY), 1l); /*initialize */
         
         /* begin outer SB iterations */
         for(ll=0; ll<iter; ll++) {
             
             /* storing old estimate */
-            copyIm(Output, Output_prev, dimX, dimY, 1);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), 1l);
             
             /* perform two GS iterations (normally 2 is enough for the convergence) */
-            gauss_seidel2D(Output, Input, Output_prev, Dx, Dy, Bx, By, dimX, dimY, lambda, mu);
-            copyIm(Output, Output_prev, dimX, dimY, 1);
+            gauss_seidel2D(Output, Input, Output_prev, Dx, Dy, Bx, By, (long)(dimX), (long)(dimY), lambda, mu);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), 1l);
             /*GS iteration */
-            gauss_seidel2D(Output, Input, Output_prev, Dx, Dy, Bx, By, dimX, dimY, lambda, mu);
+            gauss_seidel2D(Output, Input, Output_prev, Dx, Dy, Bx, By, (long)(dimX), (long)(dimY), lambda, mu);
             
             /* TV-related step */
-            if (methodTV == 1)  updDxDy_shrinkAniso2D(Output, Dx, Dy, Bx, By, dimX, dimY, lambda);
-            else updDxDy_shrinkIso2D(Output, Dx, Dy, Bx, By, dimX, dimY, lambda);
+            if (methodTV == 1)  updDxDy_shrinkAniso2D(Output, Dx, Dy, Bx, By, (long)(dimX), (long)(dimY), lambda);
+            else updDxDy_shrinkIso2D(Output, Dx, Dy, Bx, By, (long)(dimX), (long)(dimY), lambda);
             
             /* update for Bregman variables */
-            updBxBy2D(Output, Dx, Dy, Bx, By, dimX, dimY);
+            updBxBy2D(Output, Dx, Dy, Bx, By, (long)(dimX), (long)(dimY));
             
             /* check early stopping criteria if epsilon not equal zero */
             if (epsil != 0) {
@@ -94,7 +95,7 @@ float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsi
 	else {
 		/* 3D case */
 		float *Output_prev=NULL, *Dx=NULL, *Dy=NULL, *Dz=NULL, *Bx=NULL, *By=NULL, *Bz=NULL;
-		DimTotal = dimX*dimY*dimZ;
+		DimTotal = (long)(dimX*dimY*dimZ);
 		
 		Output_prev = calloc(DimTotal, sizeof(float));
 		Dx = calloc(DimTotal, sizeof(float));
@@ -104,26 +105,26 @@ float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsi
 		By = calloc(DimTotal, sizeof(float));
 		Bz = calloc(DimTotal, sizeof(float));
         
-        copyIm(Input, Output, dimX, dimY, dimZ); /*initialize */
+        copyIm(Input, Output, (long)(dimX), (long)(dimY), (long)(dimZ)); /*initialize */
         
         /* begin outer SB iterations */
         for(ll=0; ll<iter; ll++) {
             
             /* storing old estimate */
-            copyIm(Output, Output_prev, dimX, dimY, dimZ);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
             
              /* perform two GS iterations (normally 2 is enough for the convergence) */
-            gauss_seidel3D(Output, Input, Output_prev, Dx, Dy, Dz, Bx, By, Bz, dimX, dimY, dimZ, lambda, mu);
-            copyIm(Output, Output_prev, dimX, dimY, dimZ);
+            gauss_seidel3D(Output, Input, Output_prev, Dx, Dy, Dz, Bx, By, Bz, (long)(dimX), (long)(dimY), (long)(dimZ), lambda, mu);
+            copyIm(Output, Output_prev, (long)(dimX), (long)(dimY), (long)(dimZ));
             /*GS iteration */
-            gauss_seidel3D(Output, Input, Output_prev, Dx, Dy, Dz, Bx, By, Bz, dimX, dimY, dimZ, lambda, mu);
+            gauss_seidel3D(Output, Input, Output_prev, Dx, Dy, Dz, Bx, By, Bz, (long)(dimX), (long)(dimY), (long)(dimZ), lambda, mu);
             
             /* TV-related step */
-            if (methodTV == 1)  updDxDyDz_shrinkAniso3D(Output, Dx, Dy, Dz, Bx, By, Bz, dimX, dimY, dimZ, lambda);
-            else updDxDyDz_shrinkIso3D(Output, Dx, Dy, Dz, Bx, By, Bz, dimX, dimY, dimZ, lambda);
+            if (methodTV == 1)  updDxDyDz_shrinkAniso3D(Output, Dx, Dy, Dz, Bx, By, Bz, (long)(dimX), (long)(dimY), (long)(dimZ), lambda);
+            else updDxDyDz_shrinkIso3D(Output, Dx, Dy, Dz, Bx, By, Bz, (long)(dimX), (long)(dimY), (long)(dimZ), lambda);
             
             /* update for Bregman variables */
-            updBxByBz3D(Output, Dx, Dy, Dz, Bx, By, Bz, dimX, dimY, dimZ);
+            updBxByBz3D(Output, Dx, Dy, Dz, Bx, By, Bz, (long)(dimX), (long)(dimY), (long)(dimZ));
             
             /* check early stopping criteria if epsilon not equal zero */
             if (epsil != 0) {
@@ -147,10 +148,10 @@ float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsi
 /********************************************************************/
 /***************************2D Functions*****************************/
 /********************************************************************/
-float gauss_seidel2D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda, float mu)
+float gauss_seidel2D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda, float mu)
 {
     float sum, normConst;
-    int i,j,i1,i2,j1,j2,index;
+    long i,j,i1,i2,j1,j2,index;
     normConst = 1.0f/(mu + 4.0f*lambda);
     
 #pragma omp parallel for shared(U) private(index,i,j,i1,i2,j1,j2,sum)
@@ -173,9 +174,9 @@ float gauss_seidel2D(float *U, float *A, float *U_prev, float *Dx, float *Dy, fl
     return *U;
 }
 
-float updDxDy_shrinkAniso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda)
+float updDxDy_shrinkAniso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda)
 {
-    int i,j,i1,j1,index;
+    long i,j,i1,j1,index;
     float val1, val11, val2, val22, denom_lam;
     denom_lam = 1.0f/lambda;
 #pragma omp parallel for shared(U,denom_lam) private(index,i,j,i1,j1,val1,val11,val2,val22)
@@ -198,9 +199,9 @@ float updDxDy_shrinkAniso2D(float *U, float *Dx, float *Dy, float *Bx, float *By
         }}
     return 1;
 }
-float updDxDy_shrinkIso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda)
+float updDxDy_shrinkIso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda)
 {
-    int i,j,i1,j1,index;
+    long i,j,i1,j1,index;
     float val1, val11, val2, denom, denom_lam;
     denom_lam = 1.0f/lambda;
     
@@ -230,9 +231,9 @@ float updDxDy_shrinkIso2D(float *U, float *Dx, float *Dy, float *Bx, float *By,
         }}
     return 1;
 }
-float updBxBy2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY)
+float updBxBy2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY)
 {
-    int i,j,i1,j1,index;
+    long i,j,i1,j1,index;
 #pragma omp parallel for shared(U) private(index,i,j,i1,j1)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -251,10 +252,10 @@ float updBxBy2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX,
 /***************************3D Functions*****************************/
 /********************************************************************/
 /*****************************************************************/
-float gauss_seidel3D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda, float mu)
+float gauss_seidel3D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda, float mu)
 {
     float normConst, d_val, b_val, sum;
-    int i,j,i1,i2,j1,j2,k,k1,k2,index;
+    long i,j,i1,i2,j1,j2,k,k1,k2,index;
     normConst = 1.0f/(mu + 6.0f*lambda);
 #pragma omp parallel for shared(U) private(index,i,j,i1,i2,j1,j2,k,k1,k2,d_val,b_val,sum)
     for(i=0; i<dimX; i++) {
@@ -280,9 +281,9 @@ float gauss_seidel3D(float *U, float *A, float *U_prev, float *Dx, float *Dy, fl
     return *U;
 }
 
-float updDxDyDz_shrinkAniso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda)
+float updDxDyDz_shrinkAniso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda)
 {
-    int i,j,i1,j1,k,k1,index;
+    long i,j,i1,j1,k,k1,index;
     float val1, val11, val2, val22, val3, val33, denom_lam;
     denom_lam = 1.0f/lambda;
 #pragma omp parallel for shared(U,denom_lam) private(index,i,j,i1,j1,k,k1,val1,val11,val2,val22,val3,val33)
@@ -310,9 +311,9 @@ float updDxDyDz_shrinkAniso3D(float *U, float *Dx, float *Dy, float *Dz, float *
             }}}
     return 1;
 }
-float updDxDyDz_shrinkIso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda)
+float updDxDyDz_shrinkIso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda)
 {
-    int i,j,i1,j1,k,k1,index;
+    long i,j,i1,j1,k,k1,index;
     float val1, val11, val2, val3, denom, denom_lam;
     denom_lam = 1.0f/lambda;
 #pragma omp parallel for shared(U,denom_lam) private(index,denom,i,j,i1,j1,k,k1,val1,val11,val2,val3)
@@ -346,9 +347,9 @@ float updDxDyDz_shrinkIso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx
             }}}
     return 1;
 }
-float updBxByBz3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ)
+float updBxByBz3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ)
 {
-    int i,j,k,i1,j1,k1,index;
+    long i,j,k,i1,j1,k1,index;
 #pragma omp parallel for shared(U) private(index,i,j,k,i1,j1,k1)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
diff --git a/Core/regularisers_CPU/SB_TV_core.h b/Core/regularisers_CPU/SB_TV_core.h
index ad4f529..7485e3b 100644
--- a/Core/regularisers_CPU/SB_TV_core.h
+++ b/Core/regularisers_CPU/SB_TV_core.h
@@ -47,15 +47,15 @@ extern "C" {
 #endif
 CCPI_EXPORT float SB_TV_CPU_main(float *Input, float *Output, float mu, int iter, float epsil, int methodTV, int printM, int dimX, int dimY, int dimZ);
 
-CCPI_EXPORT float gauss_seidel2D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda, float mu);
-CCPI_EXPORT float updDxDy_shrinkAniso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda);
-CCPI_EXPORT float updDxDy_shrinkIso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY, float lambda);
-CCPI_EXPORT float updBxBy2D(float *U, float *Dx, float *Dy, float *Bx, float *By, int dimX, int dimY);
-
-CCPI_EXPORT float gauss_seidel3D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda, float mu);
-CCPI_EXPORT float updDxDyDz_shrinkAniso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda);
-CCPI_EXPORT float updDxDyDz_shrinkIso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ, float lambda);
-CCPI_EXPORT float updBxByBz3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float gauss_seidel2D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda, float mu);
+CCPI_EXPORT float updDxDy_shrinkAniso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda);
+CCPI_EXPORT float updDxDy_shrinkIso2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY, float lambda);
+CCPI_EXPORT float updBxBy2D(float *U, float *Dx, float *Dy, float *Bx, float *By, long dimX, long dimY);
+
+CCPI_EXPORT float gauss_seidel3D(float *U, float *A, float *U_prev, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda, float mu);
+CCPI_EXPORT float updDxDyDz_shrinkAniso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda);
+CCPI_EXPORT float updDxDyDz_shrinkIso3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ, float lambda);
+CCPI_EXPORT float updBxByBz3D(float *U, float *Dx, float *Dy, float *Dz, float *Bx, float *By, float *Bz, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/TGV_core.c b/Core/regularisers_CPU/TGV_core.c
index db88614..edd77b2 100644
--- a/Core/regularisers_CPU/TGV_core.c
+++ b/Core/regularisers_CPU/TGV_core.c
@@ -39,10 +39,11 @@ limitations under the License.
  
 float TGV_main(float *U0, float *U, float lambda, float alpha1, float alpha0, int iter, float L2, int dimX, int dimY)
 {
-	int DimTotal, ll;
+	long DimTotal;
+    int ll;
 	float *U_old, *P1, *P2, *Q1, *Q2, *Q3, *V1, *V1_old, *V2, *V2_old, tau, sigma;
 		
-		DimTotal = dimX*dimY;
+		DimTotal = (long)(dimX*dimY);
         
         /* dual variables */
         P1 = calloc(DimTotal, sizeof(float));
@@ -59,7 +60,7 @@ float TGV_main(float *U0, float *U, float lambda, float alpha1, float alpha0, in
         V2 = calloc(DimTotal, sizeof(float));
         V2_old = calloc(DimTotal, sizeof(float));
         
-        copyIm(U0, U, dimX, dimY, 1); /* initialize  */
+        copyIm(U0, U, (long)(dimX), (long)(dimY), 1l); /* initialize  */
        
         tau = pow(L2,-0.5);
         sigma = pow(L2,-0.5);
@@ -68,36 +69,36 @@ float TGV_main(float *U0, float *U, float lambda, float alpha1, float alpha0, in
         for(ll = 0; ll < iter; ll++) {
             
             /* Calculate Dual Variable P */
-            DualP_2D(U, V1, V2, P1, P2, dimX, dimY, sigma);
+            DualP_2D(U, V1, V2, P1, P2, (long)(dimX), (long)(dimY), sigma);
             
             /*Projection onto convex set for P*/
-            ProjP_2D(P1, P2, dimX, dimY, alpha1);
+            ProjP_2D(P1, P2, (long)(dimX), (long)(dimY), alpha1);
             
             /* Calculate Dual Variable Q */
-            DualQ_2D(V1, V2, Q1, Q2, Q3, dimX, dimY, sigma);
+            DualQ_2D(V1, V2, Q1, Q2, Q3, (long)(dimX), (long)(dimY), sigma);
             
             /*Projection onto convex set for Q*/
-            ProjQ_2D(Q1, Q2, Q3, dimX, dimY, alpha0);
+            ProjQ_2D(Q1, Q2, Q3, (long)(dimX), (long)(dimY), alpha0);
             
             /*saving U into U_old*/
-            copyIm(U, U_old, dimX, dimY, 1);
+            copyIm(U, U_old, (long)(dimX), (long)(dimY), 1l);
             
             /*adjoint operation  -> divergence and projection of P*/
-            DivProjP_2D(U, U0, P1, P2, dimX, dimY, lambda, tau);
+            DivProjP_2D(U, U0, P1, P2, (long)(dimX), (long)(dimY), lambda, tau);
             
             /*get updated solution U*/
-            newU(U, U_old, dimX, dimY);
+            newU(U, U_old, (long)(dimX), (long)(dimY));
             
             /*saving V into V_old*/
-            copyIm(V1, V1_old, dimX, dimY, 1);
-            copyIm(V2, V2_old, dimX, dimY, 1);
+            copyIm(V1, V1_old, (long)(dimX), (long)(dimY), 1l);
+            copyIm(V2, V2_old, (long)(dimX), (long)(dimY), 1l);
             
             /* upd V*/
-            UpdV_2D(V1, V2, P1, P2, Q1, Q2, Q3, dimX, dimY, tau);
+            UpdV_2D(V1, V2, P1, P2, Q1, Q2, Q3, (long)(dimX), (long)(dimY), tau);
             
             /*get new V*/
-            newU(V1, V1_old, dimX, dimY);
-            newU(V2, V2_old, dimX, dimY);
+            newU(V1, V1_old, (long)(dimX), (long)(dimY));
+            newU(V2, V2_old, (long)(dimX), (long)(dimY));
         } /*end of iterations*/
 	return *U;
 }
@@ -107,9 +108,9 @@ float TGV_main(float *U0, float *U, float lambda, float alpha1, float alpha0, in
 /********************************************************************/
 
 /*Calculating dual variable P (using forward differences)*/
-float DualP_2D(float *U, float *V1, float *V2, float *P1, float *P2, int dimX, int dimY, float sigma)
+float DualP_2D(float *U, float *V1, float *V2, float *P1, float *P2, long dimX, long dimY, float sigma)
 {
-    int i,j, index;
+    long i,j, index;
 #pragma omp parallel for shared(U,V1,V2,P1,P2) private(i,j,index)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -123,10 +124,10 @@ float DualP_2D(float *U, float *V1, float *V2, float *P1, float *P2, int dimX, i
     return 1;
 }
 /*Projection onto convex set for P*/
-float ProjP_2D(float *P1, float *P2, int dimX, int dimY, float alpha1)
+float ProjP_2D(float *P1, float *P2, long dimX, long dimY, float alpha1)
 {
     float grad_magn;
-    int i,j,index;
+    long i,j,index;
 #pragma omp parallel for shared(P1,P2) private(i,j,index,grad_magn)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -141,9 +142,9 @@ float ProjP_2D(float *P1, float *P2, int dimX, int dimY, float alpha1)
     return 1;
 }
 /*Calculating dual variable Q (using forward differences)*/
-float DualQ_2D(float *V1, float *V2, float *Q1, float *Q2, float *Q3, int dimX, int dimY, float sigma)
+float DualQ_2D(float *V1, float *V2, float *Q1, float *Q2, float *Q3, long dimX, long dimY, float sigma)
 {
-    int i,j,index;
+    long i,j,index;
     float q1, q2, q11, q22;
 #pragma omp parallel for shared(Q1,Q2,Q3,V1,V2) private(i,j,index,q1,q2,q11,q22)
     for(i=0; i<dimX; i++) {
@@ -172,10 +173,10 @@ float DualQ_2D(float *V1, float *V2, float *Q1, float *Q2, float *Q3, int dimX,
         }}
     return 1;
 }
-float ProjQ_2D(float *Q1, float *Q2, float *Q3, int dimX, int dimY, float alpha0)
+float ProjQ_2D(float *Q1, float *Q2, float *Q3, long dimX, long dimY, float alpha0)
 {
     float grad_magn;
-    int i,j,index;
+    long i,j,index;
 #pragma omp parallel for shared(Q1,Q2,Q3) private(i,j,index,grad_magn)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -191,9 +192,9 @@ float ProjQ_2D(float *Q1, float *Q2, float *Q3, int dimX, int dimY, float alpha0
     return 1;
 }
 /* Divergence and projection for P*/
-float DivProjP_2D(float *U, float *U0, float *P1, float *P2, int dimX, int dimY, float lambda, float tau)
+float DivProjP_2D(float *U, float *U0, float *P1, float *P2, long dimX, long dimY, float lambda, float tau)
 {
-    int i,j,index;
+    long i,j,index;
     float P_v1, P_v2, div;
 #pragma omp parallel for shared(U,U0,P1,P2) private(i,j,index,P_v1,P_v2,div)
     for(i=0; i<dimX; i++) {
@@ -209,17 +210,17 @@ float DivProjP_2D(float *U, float *U0, float *P1, float *P2, int dimX, int dimY,
     return *U;
 }
 /*get updated solution U*/
-float newU(float *U, float *U_old, int dimX, int dimY)
+float newU(float *U, float *U_old, long dimX, long dimY)
 {
-    int i;
+    long i;
 #pragma omp parallel for shared(U,U_old) private(i)
     for(i=0; i<dimX*dimY; i++) U[i] = 2*U[i] - U_old[i];
     return *U;
 }
 /*get update for V*/
-float UpdV_2D(float *V1, float *V2, float *P1, float *P2, float *Q1, float *Q2, float *Q3, int dimX, int dimY, float tau)
+float UpdV_2D(float *V1, float *V2, float *P1, float *P2, float *Q1, float *Q2, float *Q3, long dimX, long dimY, float tau)
 {
-    int i, j, index;
+    long i, j, index;
     float q1, q11, q2, q22, div1, div2;
 #pragma omp parallel for shared(V1,V2,P1,P2,Q1,Q2,Q3) private(i, j, index, q1, q11, q2, q22, div1, div2)
     for(i=0; i<dimX; i++) {
diff --git a/Core/regularisers_CPU/TGV_core.h b/Core/regularisers_CPU/TGV_core.h
index d5632b9..29482ad 100644
--- a/Core/regularisers_CPU/TGV_core.h
+++ b/Core/regularisers_CPU/TGV_core.h
@@ -49,13 +49,14 @@ extern "C" {
 #endif
 /* 2D functions */
 CCPI_EXPORT float TGV_main(float *U0, float *U, float lambda, float alpha1, float alpha0, int iter, float L2, int dimX, int dimY);
-CCPI_EXPORT float DualP_2D(float *U, float *V1, float *V2, float *P1, float *P2, int dimX, int dimY, float sigma);
-CCPI_EXPORT float ProjP_2D(float *P1, float *P2, int dimX, int dimY, float alpha1);
-CCPI_EXPORT float DualQ_2D(float *V1, float *V2, float *Q1, float *Q2, float *Q3, int dimX, int dimY, float sigma);
-CCPI_EXPORT float ProjQ_2D(float *Q1, float *Q2, float *Q3, int dimX, int dimY, float alpha0);
-CCPI_EXPORT float DivProjP_2D(float *U, float *U0, float *P1, float *P2, int dimX, int dimY, float lambda, float tau);
-CCPI_EXPORT float UpdV_2D(float *V1, float *V2, float *P1, float *P2, float *Q1, float *Q2, float *Q3, int dimX, int dimY, float tau);
-CCPI_EXPORT float newU(float *U, float *U_old, int dimX, int dimY);
+
+CCPI_EXPORT float DualP_2D(float *U, float *V1, float *V2, float *P1, float *P2, long dimX, long dimY, float sigma);
+CCPI_EXPORT float ProjP_2D(float *P1, float *P2, long dimX, long dimY, float alpha1);
+CCPI_EXPORT float DualQ_2D(float *V1, float *V2, float *Q1, float *Q2, float *Q3, long dimX, long dimY, float sigma);
+CCPI_EXPORT float ProjQ_2D(float *Q1, float *Q2, float *Q3, long dimX, long dimY, float alpha0);
+CCPI_EXPORT float DivProjP_2D(float *U, float *U0, float *P1, float *P2, long dimX, long dimY, float lambda, float tau);
+CCPI_EXPORT float UpdV_2D(float *V1, float *V2, float *P1, float *P2, float *Q1, float *Q2, float *Q3, long dimX, long dimY, float tau);
+CCPI_EXPORT float newU(float *U, float *U_old, long dimX, long dimY);
 #ifdef __cplusplus
 }
 #endif
diff --git a/Core/regularisers_CPU/TNV_core.c b/Core/regularisers_CPU/TNV_core.c
index c51d6cd..753cc5f 100755
--- a/Core/regularisers_CPU/TNV_core.c
+++ b/Core/regularisers_CPU/TNV_core.c
@@ -39,16 +39,16 @@
 
 float TNV_CPU_main(float *Input, float *u, float lambda, int maxIter, float tol, int dimX, int dimY, int dimZ)
 {
-    int k, p, q, r, DimTotal;
+    long k, p, q, r, DimTotal;
     float taulambda;
     float *u_upd, *gx, *gy, *gx_upd, *gy_upd, *qx, *qy, *qx_upd, *qy_upd, *v, *vx, *vy, *gradx, *grady, *gradx_upd, *grady_upd, *gradx_ubar, *grady_ubar, *div, *div_upd;
     
-    p = 1;
-    q = 1;
-    r = 0;
+    p = 1l;
+    q = 1l;
+    r = 0l;
     
     lambda = 1.0f/(2.0f*lambda);
-    DimTotal = dimX*dimY*dimZ;
+    DimTotal = (long)(dimX*dimY*dimZ);
     /* PDHG algorithm parameters*/
     float tau = 0.5f;
     float sigma = 0.5f;
@@ -106,10 +106,10 @@ float TNV_CPU_main(float *Input, float *u, float lambda, int maxIter, float tol,
         for(k=0; k<dimX*dimY*dimZ; k++)  {v[k] = u[k] + tau*div[k];}
 
 // Proximal solution of fidelity term
-proxG(u_upd, v, Input, taulambda, dimX, dimY, dimZ);
+proxG(u_upd, v, Input, taulambda, (long)(dimX), (long)(dimY), (long)(dimZ));
 
 // Gradient of updated primal variable
-gradient(u_upd, gradx_upd, grady_upd, dimX, dimY, dimZ);
+gradient(u_upd, gradx_upd, grady_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
 
 // Argument of proximal mapping of regularization term
 #pragma omp parallel for shared(gradx_upd, grady_upd, gradx, grady) private(k, ubarx, ubary)
@@ -122,7 +122,7 @@ for(k=0; k<dimX*dimY*dimZ; k++) {
     grady_ubar[k] = ubary;
 }
 
-proxF(gx_upd, gy_upd, vx, vy, sigma, p, q, r, dimX, dimY, dimZ);
+proxF(gx_upd, gy_upd, vx, vy, sigma, p, q, r, (long)(dimX), (long)(dimY), (long)(dimZ));
 
 // Update dual variable
 #pragma omp parallel for shared(qx_upd, qy_upd) private(k)
@@ -174,14 +174,14 @@ if(b > 1)
     sigma = (beta * sigma) / b;
     alpha = alpha0;
     
-    copyIm(u, u_upd, dimX, dimY, dimZ);
-    copyIm(gx, gx_upd, dimX, dimY, dimZ);
-    copyIm(gy, gy_upd, dimX, dimY, dimZ);
-    copyIm(qx, qx_upd, dimX, dimY, dimZ);
-    copyIm(qy, qy_upd, dimX, dimY, dimZ);
-    copyIm(gradx, gradx_upd, dimX, dimY, dimZ);
-    copyIm(grady, grady_upd, dimX, dimY, dimZ);
-    copyIm(div, div_upd, dimX, dimY, dimZ);
+    copyIm(u, u_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(gx, gx_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(gy, gy_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(qx, qx_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(qy, qy_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(gradx, gradx_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(grady, grady_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
+    copyIm(div, div_upd, (long)(dimX), (long)(dimY), (long)(dimZ));
     
 } else if(resprimal > dual_dot_delta)
 {
@@ -205,14 +205,14 @@ taulambda = tau * lambda;
 divsigma = 1.0f / sigma;
 divtau = 1.0f / tau;
 
-copyIm(u_upd, u, dimX, dimY, dimZ);
-copyIm(gx_upd, gx, dimX, dimY, dimZ);
-copyIm(gy_upd, gy, dimX, dimY, dimZ);
-copyIm(qx_upd, qx, dimX, dimY, dimZ);
-copyIm(qy_upd, qy, dimX, dimY, dimZ);
-copyIm(gradx_upd, gradx, dimX, dimY, dimZ);
-copyIm(grady_upd, grady, dimX, dimY, dimZ);
-copyIm(div_upd, div, dimX, dimY, dimZ);
+copyIm(u_upd, u, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(gx_upd, gx, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(gy_upd, gy, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(qx_upd, qx, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(qy_upd, qy, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(gradx_upd, gradx, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(grady_upd, grady, (long)(dimX), (long)(dimY), (long)(dimZ));
+copyIm(div_upd, div, (long)(dimX), (long)(dimY), (long)(dimZ));
 
 // Compute residual at current iteration
 residual = (resprimal + resdual) / ((float) (dimX*dimY*dimZ));
@@ -227,15 +227,14 @@ if (residual < tol) {
     free (u_upd); free(gx); free(gy); free(gx_upd); free(gy_upd);
     free(qx); free(qy); free(qx_upd); free(qy_upd); free(v); free(vx); free(vy);
     free(gradx); free(grady); free(gradx_upd); free(grady_upd); free(gradx_ubar);
-    free(grady_ubar); free(div); free(div_upd);
-    
+    free(grady_ubar); free(div); free(div_upd);    
     return *u;
 }
 
-float proxG(float *u_upd, float *v, float *f, float taulambda, int dimX, int dimY, int dimZ)
+float proxG(float *u_upd, float *v, float *f, float taulambda, long dimX, long dimY, long dimZ)
 {
     float constant;
-    int k;
+    long k;
     constant = 1.0f + taulambda;
 #pragma omp parallel for shared(v, f, u_upd) private(k)
     for(k=0; k<dimZ*dimX*dimY; k++) {
@@ -245,9 +244,9 @@ float proxG(float *u_upd, float *v, float *f, float taulambda, int dimX, int dim
     return *u_upd;
 }
 
-float gradient(float *u_upd, float *gradx_upd, float *grady_upd, int dimX, int dimY, int dimZ)
+float gradient(float *u_upd, float *gradx_upd, float *grady_upd, long dimX, long dimY, long dimZ)
 {
-    int i, j, k, l;
+    long i, j, k, l;
     // Compute discrete gradient using forward differences
 #pragma omp parallel for shared(gradx_upd,grady_upd,u_upd) private(i, j, k, l)
     for(k = 0; k < dimZ; k++)   {
@@ -270,7 +269,7 @@ float gradient(float *u_upd, float *gradx_upd, float *grady_upd, int dimX, int d
     return 1;
 }
 
-float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int q, int r, int dimX, int dimY, int dimZ)
+float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int q, int r, long dimX, long dimY, long dimZ)
 {
     // (S^p, \ell^1) norm decouples at each pixel
 //   Spl1(gx, gy, vx, vy, sigma, p, num_channels, dim);
@@ -289,7 +288,7 @@ float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int
     // Auxiliar vector
     float *proj, sum, shrinkfactor ;
     float M1,M2,M3,valuex,valuey,T,D,det,eig1,eig2,sig1,sig2,V1, V2, V3, V4, v0,v1,v2, mu1,mu2,sig1_upd,sig2_upd,t1,t2,t3;
-    int i,j,k, ii, num;
+    long i,j,k, ii, num;
 #pragma omp parallel for shared (gx,gy,vx,vy,p) private(i,ii,j,k,proj,num, sum, shrinkfactor, M1,M2,M3,valuex,valuey,T,D,det,eig1,eig2,sig1,sig2,V1, V2, V3, V4,v0,v1,v2,mu1,mu2,sig1_upd,sig2_upd,t1,t2,t3)
     for(i=0; i<dimX; i++) {
         for(j=0; j<dimY; j++) {
@@ -372,7 +371,7 @@ float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int
                 
                 /*l1 projection part */
                 sum = fLarge;
-                num = 0;
+                num = 0l;
                 shrinkfactor = 0.0f;
                 while(sum > 1.0f)
                 {
@@ -424,9 +423,9 @@ float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int
     return 1;
 }
 
-float divergence(float *qx_upd, float *qy_upd, float *div_upd, int dimX, int dimY, int dimZ)
+float divergence(float *qx_upd, float *qy_upd, float *div_upd, long dimX, long dimY, long dimZ)
 {
-    int i, j, k, l;
+    long i, j, k, l;
 #pragma omp parallel for shared(qx_upd,qy_upd,div_upd) private(i, j, k, l)
     for(k = 0; k < dimZ; k++)   {
         for(j = 0; j < dimY; j++)   {
diff --git a/Core/regularisers_CPU/TNV_core.h b/Core/regularisers_CPU/TNV_core.h
index c082694..aa050a4 100644
--- a/Core/regularisers_CPU/TNV_core.h
+++ b/Core/regularisers_CPU/TNV_core.h
@@ -38,10 +38,10 @@ extern "C" {
 CCPI_EXPORT float TNV_CPU_main(float *Input, float *u, float lambda, int maxIter, float tol, int dimX, int dimY, int dimZ);
 
 /*float PDHG(float *A, float *B, float tau, float sigma, float theta, float lambda, int p, int q, int r, float tol, int maxIter, int d_c, int d_w, int d_h);*/
-CCPI_EXPORT float proxG(float *u_upd, float *v, float *f, float taulambda, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float gradient(float *u_upd, float *gradx_upd, float *grady_upd, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int q, int r, int dimX, int dimY, int dimZ);
-CCPI_EXPORT float divergence(float *qx_upd, float *qy_upd, float *div_upd, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float proxG(float *u_upd, float *v, float *f, float taulambda, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float gradient(float *u_upd, float *gradx_upd, float *grady_upd, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float proxF(float *gx, float *gy, float *vx, float *vy, float sigma, int p, int q, int r, long dimX, long dimY, long dimZ);
+CCPI_EXPORT float divergence(float *qx_upd, float *qy_upd, float *div_upd, long dimX, long dimY, long dimZ);
 #ifdef __cplusplus
 }
 #endif
\ No newline at end of file
diff --git a/Core/regularisers_CPU/utils.c b/Core/regularisers_CPU/utils.c
index 3eebdaa..7a4e80b 100644
--- a/Core/regularisers_CPU/utils.c
+++ b/Core/regularisers_CPU/utils.c
@@ -21,9 +21,9 @@ limitations under the License.
 #include <math.h>
 
 /* Copy Image (float) */
-float copyIm(float *A, float *U, int dimX, int dimY, int dimZ)
+float copyIm(float *A, float *U, long dimX, long dimY, long dimZ)
 {
-	int j;
+	long j;
 #pragma omp parallel for shared(A, U) private(j)
 	for (j = 0; j<dimX*dimY*dimZ; j++)  U[j] = A[j];
 	return *U;
@@ -88,18 +88,18 @@ float TV_energy2D(float *U, float *U0, float *E_val, float lambda, int type, int
 
 float TV_energy3D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY, int dimZ)
 {
-	int i, j, k, i1, j1, k1, index;
+	long i, j, k, i1, j1, k1, index;
 	float NOMx_2, NOMy_2, NOMz_2, E_Grad=0.0f, E_Data=0.0f;
 	
 	/* first calculate \grad U_xy*/	
-    for(j=0; j<dimY; j++) {
-        for(i=0; i<dimX; i++) {
-            for(k=0; k<dimZ; k++) {
+    for(j=0; j<(long)(dimY); j++) {
+        for(i=0; i<(long)(dimX); i++) {
+            for(k=0; k<(long)(dimZ); k++) {
 				index = (dimX*dimY)*k + j*dimX+i;
                 /* boundary conditions */
-                i1 = i + 1; if (i == dimX-1) i1 = i;
-                j1 = j + 1; if (j == dimY-1) j1 = j;
-                k1 = k + 1; if (k == dimZ-1) k1 = k;
+                i1 = i + 1; if (i == (long)(dimX-1)) i1 = i;
+                j1 = j + 1; if (j == (long)(dimY-1)) j1 = j;
+                k1 = k + 1; if (k == (long)(dimZ-1)) k1 = k;
                 
                 /* Forward differences */                
                 NOMx_2 = powf((float)(U[(dimX*dimY)*k + j1*dimX+i] - U[index]),2); /* x+ */
diff --git a/Core/regularisers_CPU/utils.h b/Core/regularisers_CPU/utils.h
index 7ad83ab..cfaf6d7 100644
--- a/Core/regularisers_CPU/utils.h
+++ b/Core/regularisers_CPU/utils.h
@@ -24,7 +24,7 @@ limitations under the License.
 #ifdef __cplusplus
 extern "C" {
 #endif
-CCPI_EXPORT float copyIm(float *A, float *U, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float copyIm(float *A, float *U, long dimX, long dimY, long dimZ);
 CCPI_EXPORT unsigned char copyIm_unchar(unsigned char *A, unsigned char *U, int dimX, int dimY, int dimZ);
 CCPI_EXPORT float copyIm_roll(float *A, float *U, int dimX, int dimY, int roll_value, int switcher);
 CCPI_EXPORT float TV_energy2D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY);
diff --git a/Wrappers/Matlab/mex_compile/compileCPU_mex_Linux.m b/Wrappers/Matlab/mex_compile/compileCPU_mex_Linux.m
index 064b416..767d811 100644
--- a/Wrappers/Matlab/mex_compile/compileCPU_mex_Linux.m
+++ b/Wrappers/Matlab/mex_compile/compileCPU_mex_Linux.m
@@ -61,7 +61,7 @@ fprintf('%s \n', 'Compiling Nonlinear/Linear diffusion inpainting...');
 mex NonlDiff_Inp.c Diffusion_Inpaint_core.c utils.c CFLAGS="\$CFLAGS -fopenmp -Wall -std=c99" LDFLAGS="\$LDFLAGS -fopenmp"
 movefile('NonlDiff_Inp.mex*',Pathmove);
 
-fprintf('%s \n', 'Compiling Nonlocal marching method for inpaiting...');
+fprintf('%s \n', 'Compiling Nonlocal marching method for inpainting...');
 mex NonlocalMarching_Inpaint.c NonlocalMarching_Inpaint_core.c utils.c CFLAGS="\$CFLAGS -fopenmp -Wall -std=c99" LDFLAGS="\$LDFLAGS -fopenmp"
 movefile('NonlocalMarching_Inpaint.mex*',Pathmove);
 
-- 
cgit v1.2.3