diff --git a/RCaN/R/RcppExports.R b/RCaN/R/RcppExports.R
index 2491311fe19e8d4fbd39b8f43bb12a83b867fb60..30b0a263a858ab33eca6ac486baa89365deb2ac4 100644
--- a/RCaN/R/RcppExports.R
+++ b/RCaN/R/RcppExports.R
@@ -5,7 +5,7 @@
 #' @useDynLib RCaN
 NULL
 
-fitCaN <- function(N, A, b, C, v, L, x0, thin) {
-    .Call('_RCaN_fitCaN', PACKAGE = 'RCaN', N, A, b, C, v, L, x0, thin)
+fitCaN <- function(N, A, b, C, v, L, x0, thin, test) {
+    .Call('_RCaN_fitCaN', PACKAGE = 'RCaN', N, A, b, C, v, L, x0, thin, test)
 }
 
diff --git a/RCaN/R/fitmyCaNmod.R b/RCaN/R/fitmyCaNmod.R
index b4d19f87cbc3d15319c02ababd039f6889fcaba1..af6725675e9aa5f9066754632242571eb388b56b 100644
--- a/RCaN/R/fitmyCaNmod.R
+++ b/RCaN/R/fitmyCaNmod.R
@@ -7,6 +7,7 @@
 #' @param nchain the number of mcmc chains
 #' @param ncore number of cores to use
 #' @param thin thinning interval
+#' @param test if TRUE try a new method aiming at reducing autocorrelation
 #' @return a \code{\link[coda]{mcmc.list}}
 #' @export
 #'
@@ -41,7 +42,8 @@ fitmyCaNmod <- function(myCaNmod,
                         N,
                         nchain = 1,
                         ncore = 1,
-                        thin = 1) {
+                        thin = 1,
+                        test = FALSE) {
   ncore <- min(min(detectCores() - 1, ncore), nchain)
   `%myinfix%` <- `%do%`
 
@@ -91,7 +93,8 @@ fitmyCaNmod <- function(myCaNmod,
         as.matrix(myCaNmod$C),
         myCaNmod$v,
         as.matrix(myCaNmod$L),
-        x0
+        x0,
+        test
       )
     names(res) <- c("F", "B")
     res$F <- res$F[, -seq_len(length(myCaNmod$species))]
diff --git a/RCaN/man/fitmyCaNmod.Rd b/RCaN/man/fitmyCaNmod.Rd
index 4549a5a874eb4b49f2937e915647c6ed389b0d8c..c0f5a822610eda63e7789c5f5c8e8999e19049cf 100644
--- a/RCaN/man/fitmyCaNmod.Rd
+++ b/RCaN/man/fitmyCaNmod.Rd
@@ -4,7 +4,8 @@
 \alias{fitmyCaNmod}
 \title{fitmyCaNmod}
 \usage{
-fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1)
+fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1,
+  test = FALSE)
 }
 \arguments{
 \item{myCaNmod}{a CaNmod object with following elements}
@@ -16,6 +17,8 @@ fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1)
 \item{ncore}{number of cores to use}
 
 \item{thin}{thinning interval}
+
+\item{test}{if TRUE try a new method aiming at reducing autocorrelation}
 }
 \value{
 a \code{\link[coda]{mcmc.list}}
diff --git a/RCaN/src/RcppExports.cpp b/RCaN/src/RcppExports.cpp
index 256bb3d7c89d5398d76ba92e29c6393513ff7b8c..af2aa6d5343c035bf7a5b3325ab8e3ff288e8af3 100644
--- a/RCaN/src/RcppExports.cpp
+++ b/RCaN/src/RcppExports.cpp
@@ -7,8 +7,8 @@
 using namespace Rcpp;
 
 // fitCaN
-List fitCaN(const int N, const Eigen::MatrixXd& A, const Eigen::VectorXd& b, const Eigen::MatrixXd& C, const Eigen::VectorXd& v, const Eigen::MatrixXd& L, const Eigen::VectorXd& x0, const int thin);
-RcppExport SEXP _RCaN_fitCaN(SEXP NSEXP, SEXP ASEXP, SEXP bSEXP, SEXP CSEXP, SEXP vSEXP, SEXP LSEXP, SEXP x0SEXP, SEXP thinSEXP) {
+List fitCaN(const int N, const Eigen::MatrixXd& A, const Eigen::VectorXd& b, const Eigen::MatrixXd& C, const Eigen::VectorXd& v, const Eigen::MatrixXd& L, const Eigen::VectorXd& x0, const int thin, const bool test);
+RcppExport SEXP _RCaN_fitCaN(SEXP NSEXP, SEXP ASEXP, SEXP bSEXP, SEXP CSEXP, SEXP vSEXP, SEXP LSEXP, SEXP x0SEXP, SEXP thinSEXP, SEXP testSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
@@ -20,13 +20,14 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< const Eigen::MatrixXd& >::type L(LSEXP);
     Rcpp::traits::input_parameter< const Eigen::VectorXd& >::type x0(x0SEXP);
     Rcpp::traits::input_parameter< const int >::type thin(thinSEXP);
-    rcpp_result_gen = Rcpp::wrap(fitCaN(N, A, b, C, v, L, x0, thin));
+    Rcpp::traits::input_parameter< const bool >::type test(testSEXP);
+    rcpp_result_gen = Rcpp::wrap(fitCaN(N, A, b, C, v, L, x0, thin, test));
     return rcpp_result_gen;
 END_RCPP
 }
 
 static const R_CallMethodDef CallEntries[] = {
-    {"_RCaN_fitCaN", (DL_FUNC) &_RCaN_fitCaN, 8},
+    {"_RCaN_fitCaN", (DL_FUNC) &_RCaN_fitCaN, 9},
     {NULL, NULL, 0}
 };
 
diff --git a/RCaN/src/fitCan.cpp b/RCaN/src/fitCan.cpp
index 74319975ea0f99a974f15b9e0d783f28521f77f0..bc28f408d3fed3503f679091e7ebc4b64a079597 100644
--- a/RCaN/src/fitCan.cpp
+++ b/RCaN/src/fitCan.cpp
@@ -18,15 +18,16 @@ using Eigen::VectorXd;
 List fitCaN(const int N, const Eigen::MatrixXd &A ,const Eigen::VectorXd &b,
             const Eigen::MatrixXd &C ,const Eigen::VectorXd &v,
             const Eigen::MatrixXd &L,
-            const Eigen::VectorXd &x0, const int thin) {
+            const Eigen::VectorXd &x0, const int thin,
+            const bool test) {
   int p=A.cols();
   int m2=C.rows();
   MatrixXd F(N, p);
   MatrixXd B(N, L.rows());
   if(m2>0){ //there are equality constraints
-    F=cpgsR::cpgsEquality(N, A, b, C, v, x0, thin);
+    F=cpgsR::cpgsEquality(N, A, b, C, v, x0, thin, test);
   } else{
-    F=cpgsR::cpgs(N, A, b, x0, thin);
+    F=cpgsR::cpgs(N, A, b, x0, thin, test);
   }
   for (int i=0;i<N;++i){
     B.row(i)=L*F.row(i).transpose();