diff --git a/DESCRIPTION b/DESCRIPTION
index 454397ae87dd6a1590d78364134efc1760534780..5b1ac5990fd5b7879a4c4bc2720aab7703fe32c0 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -86,4 +86,4 @@ VignetteBuilder: knitr
 License: Apache License 2.0
 Encoding: UTF-8
 LazyData: true
-RoxygenNote: 7.0.2
+RoxygenNote: 7.2.0
diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R
index 263128aa1681a7d5f4133e1364fdf7f1d66cb562..c3097070b8a7067ac7c97951944ee98b9ad080e6 100644
--- a/R/CST_BiasCorrection.R
+++ b/R/CST_BiasCorrection.R
@@ -1,26 +1,37 @@
 #' Bias Correction based on the mean and standard deviation adjustment
 #'
 #'@author Verónica Torralba, \email{veronica.torralba@bsc.es}
-#'@description This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset.
+#'@description This function applies the simple bias adjustment technique 
+#'described in Torralba et al. (2017). The adjusted forecasts have an equivalent 
+#'standard deviation and mean to that of the reference dataset.
 #'
-#'@param exp an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonal forecast experiment data in the element named \code{$data}
-#'@param obs an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the observed data in the element named \code{$data}.
-#'@param exp_cor an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.
-#'@param na.rm a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.
+#'@param exp An object of class \code{s2dv_cube} as returned by \code{CST_Load} 
+#'  function, containing the seasonal forecast experiment data in the element 
+#'  named \code{$data}
+#'@param obs An object of class \code{s2dv_cube} as returned by \code{CST_Load} 
+#'  function, containing the observed data in the element named \code{$data}.
+#'@param exp_cor An object of class \code{s2dv_cube} as returned by 
+#'  \code{CST_Load} function, containing the seasonl forecast experiment to be 
+#'  corrected. If it is NULL, the 'exp' forecast will be corrected.
+#'@param na.rm A logical value indicating whether missing values should be 
+#'  stripped before the computation proceeds, by default it is set to FALSE.
+#'@param memb_dim A character string indicating the name of the member 
+#'  dimension. By default, it is set to 'member'.
+#'@param sdate_dim A character string indicating the name of the start date 
+#'  dimension. By default, it is set to 'sdate'.
+#'@param ncores An integer that indicates the number of cores for parallel 
+#'  computations using multiApply function. The default value is NULL.
+#'@return An object of class \code{s2dv_cube} containing the bias corrected 
+#'forecasts with the same dimensions of the experimental data.
 #'
-#'@param memb_dim is a character string indicating the name of the member dimension. By default, it is set to 'member'.
-#'@param sdate_dim is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.
-#'@param ncores is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.
-#'@return an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data.
+#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. 
+#'Davis (2017). Seasonal climate prediction: a new source of information for 
+#'the management of wind energy resources. Journal of Applied Meteorology and 
+#'Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, 
+#'EUPORIAS, NEWA, RESILIENCE, SPECS)
 #'
-#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS)
-#'
-#'@import multiApply
 #'@examples
-#'
 #'# Example
-#'# Creation of sample s2dverification objects. These are not complete
-#'# s2dverification objects though. The Load function returns complete objects.
 #'mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7)
 #'dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7)
 #'obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7)
@@ -32,118 +43,168 @@
 #'attr(exp, 'class') <- 's2dv_cube'
 #'attr(obs, 'class') <- 's2dv_cube'
 #'a <- CST_BiasCorrection(exp = exp, obs = obs)
-#'str(a)
+#'@import multiApply
 #'@export
 CST_BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, 
                                memb_dim = 'member', sdate_dim = 'sdate',
-                               ncores = 1) {
+                               ncores = NULL) {
   if (!inherits(exp, 's2dv_cube') || !inherits(obs, 's2dv_cube')) {
     stop("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ",
          "as output by CSTools::CST_Load.")
   }
-  if (dim(obs$data)['member'] != 1) {
-    stop("The length of the dimension 'member' in the component 'data' ",
-         "of the parameter 'obs' must be equal to 1.")
+  if (!is.null(exp_cor)) {
+    if (!inherits(exp_cor, 's2dv_cube')) {
+      stop("Parameter 'exp_cor' must be of the class 's2dv_cube', ",
+         "as output by CSTools::CST_Load.")
+    }
+    dimnames <- names(dim(exp_cor$data))
+  } else {
+    dimnames <- names(dim(exp$data))
   }
+
+  BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data, exp_cor = exp_cor$data, 
+                                  memb_dim = memb_dim, sdate_dim = sdate_dim,
+                                  na.rm = na.rm, ncores = ncores)
+  
+  pos <- match(dimnames, names(dim(BiasCorrected)))
+  BiasCorrected <- aperm(BiasCorrected, pos)
+
   if (is.null(exp_cor)) {
-    dimnames <- names(dim(exp$data))
-    BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data, na.rm = na.rm,
-                                    memb_dim = memb_dim, sdate_dim = sdate_dim,
-                                    ncores = ncores)
-    pos <- match(dimnames, names(dim(BiasCorrected)))
-    BiasCorrected <- aperm(BiasCorrected, pos)
-    names(dim(BiasCorrected)) <- dimnames
     exp$data <- BiasCorrected
     exp$Datasets <- c(exp$Datasets, obs$Datasets)
     exp$source_files <- c(exp$source_files, obs$source_files)
+    
     return(exp)
-  } else { 
-    if (!inherits(exp_cor, 's2dv_cube')) {
-      stop("Parameter 'var_cor' must be of the class 's2dv_cube'.")
-    }
-    dimnames <- names(dim(exp_cor$data))
-    BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data,
-                                    exp_cor = exp_cor$data, na.rm = na.rm)
-    pos <- match(dimnames, names(dim(BiasCorrected)))
-    BiasCorrected <- aperm(BiasCorrected, pos)
-    names(dim(BiasCorrected)) <- dimnames
+
+  } else {
     exp_cor$data <- BiasCorrected
     exp_cor$Datasets <- c(exp_cor$Datasets, exp$Datasets, obs$Datasets)
-    exp_cor$source_files <- c(exp_cor$source_files,
-                          exp$source_files, obs$source_files)
+    exp_cor$source_files <- c(exp_cor$source_files, exp$source_files, obs$source_files)
+
     return(exp_cor)
   }
 }
+
 #' Bias Correction based on the mean and standard deviation adjustment
 #'
 #'@author Verónica Torralba, \email{veronica.torralba@bsc.es}
-#'@description This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset.
+#'@description This function applies the simple bias adjustment technique 
+#'described in Torralba et al. (2017). The adjusted forecasts have an equivalent 
+#'standard deviation and mean to that of the reference dataset.
 #'
-#'@param exp a multidimensional array with named dimensions containing the seasonal forecast experiment data with at least 'member' and 'sdate' dimensions.
-#'@param obs a multidimensional array with named dimensions containing the observed data with at least 'sdate' dimension.
-#'@param exp_cor a multidimensional array with named dimensions containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.
-#'@param na.rm a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.
-#'@param memb_dim is a character string indicating the name of the member dimension. By default, it is set to 'member'.
-#'@param sdate_dim is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.
-#'@param ncores is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.
+#'@param exp A multidimensional array with named dimensions containing the 
+#'  seasonal forecast experiment data with at least 'member' and 'sdate' 
+#'  dimensions.
+#'@param obs A multidimensional array with named dimensions containing the 
+#'  observed data with at least 'sdate' dimension.
+#'@param exp_cor A multidimensional array with named dimensions containing the 
+#'  seasonl forecast experiment to be corrected. If it is NULL, the 'exp' 
+#'  forecast will be corrected.
+#'@param na.rm A logical value indicating whether missing values should be 
+#'  stripped before the computation proceeds, by default it is set to FALSE.
+#'@param memb_dim A character string indicating the name of the member 
+#'  dimension. By default, it is set to 'member'.
+#'@param sdate_dim A character string indicating the name of the start date 
+#'  dimension. By default, it is set to 'sdate'.
+#'@param ncores An integer that indicates the number of cores for parallel 
+#'  computations using multiApply function. The default value is NULL.
 #'
-#'@return an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data.
+#'@return An array containing the bias corrected forecasts with the same 
+#'dimensions of the experimental data.
 #'
-#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS)
+#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. 
+#'Davis (2017). Seasonal climate prediction: a new source of information for the 
+#'management of wind energy resources. Journal of Applied Meteorology and 
+#'Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, 
+#'EUPORIAS, NEWA, RESILIENCE, SPECS)
 #'
-#'@import multiApply
 #'@examples
-#'
 #'# Example
-#'# Creation of sample s2dverification objects. These are not complete
-#'# s2dverification objects though. The Load function returns complete objects.
 #'mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7)
 #'dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7)
 #'obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7)
 #'dim(obs1) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, lat = 6, lon = 7)
 #'a <- BiasCorrection(exp = mod1, obs = obs1)
-#'str(a)
+#'@import multiApply
 #'@export
-BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE,
-                            memb_dim = 'member', sdate_dim = 'sdate',
-                            ncores = 1) {
- 
-   if (!all(c(memb_dim, sdate_dim) %in% names(dim(exp)))) {
-    stop(paste("Parameter 'exp' must have the dimensions set up in 'memb_dim'", 
-               "and 'sdate_dim' parameters."))
+BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE,
+                           memb_dim = 'member', sdate_dim = 'sdate',
+                           ncores = NULL) {
+  # Check inputs
+  ## exp, obs
+  if (!is.array(exp) || !is.numeric(exp)) {
+    stop("Parameter 'exp' must be a numeric array.")
   }
-
-  if (!all(sdate_dim %in% names(dim(obs)))) {
-    stop("Parameter 'obs' must have the dimension specified in 'sdate_dim'.")
+  if (!is.array(obs) || !is.numeric(obs)) {
+    stop("Parameter 'obs' must be a numeric array.")
+  }
+  obsdims <- names(dim(obs))
+  expdims <- names(dim(exp))
+  if (is.null(expdims)) {
+    stop("Parameter 'exp' must have dimension names.")
+  }
+  if (is.null(obsdims)) {
+    stop("Parameter 'obs' must have dimension names.")
   }
-
   if (any(is.na(exp)))  {
     warning("Parameter 'exp' contains NA values.")
   }
-  
   if (any(is.na(obs)))  {
     warning("Parameter 'obs' contains NA values.")
   }
-
+  ## exp_cor
+  if (!is.null(exp_cor)) {
+    exp_cordims <- names(dim(exp_cor))
+    if (is.null(exp_cordims)) {
+      stop("Parameter 'exp_cor' must have dimension names.")
+    }
+  }
+  ## sdate_dim, memb_dim
+  if (!is.character(sdate_dim) || length(sdate_dim) != 1) {
+    stop("Parameter 'sdate_dim' must be a character string.")
+  }
+  if (!sdate_dim %in% expdims || !sdate_dim %in% obsdims) {
+    stop("Parameter 'sdate_dim' is not found in 'exp' or 'obs' dimension.")
+  }
+  if (dim(exp)[sdate_dim] == 1) {
+    stop("Parameter 'exp' must have dimension length of 'sdate_dim' bigger than 1.")
+  }
+  if (!all(c(memb_dim, sdate_dim) %in% expdims)) {
+    stop("Parameter 'exp' requires 'sdate_dim' and 'memb_dim' dimensions.")  
+  }
+  if (memb_dim %in% obsdims) {
+    if (dim(obs)[memb_dim] != 1) {
+      stop("If parameter 'obs' has dimension 'memb_dim' its length must be equal to 1.")
+    }
+  }
+  ## na.rm
   if (!is.logical(na.rm))  {
     na.rm <- FALSE
     warning("Paramater 'na.rm' must be a logical, it has been set to FALSE.")
   }
-
-  if (length(na.rm)>1) {
+  if (length(na.rm) > 1) {
     na.rm <- na.rm[1]
     warning("Paramter 'na.rm' has length greater than 1, and only the fist element is used.")
   }
+  ## ncores
+  if (!is.null(ncores)) {
+    if (!is.numeric(ncores) | ncores %% 1 != 0 | ncores <= 0 |
+      length(ncores) > 1) {
+      stop("Parameter 'ncores' must be either NULL or a positive integer.")
+    }
+  }
 
   target_dims_obs <- sdate_dim
   if (memb_dim %in% names(dim(obs))) {
     target_dims_obs <- c(memb_dim, target_dims_obs)
   }
+  
   if (is.null(exp_cor)) {
     BiasCorrected <- Apply(data = list(var_obs = obs, var_exp = exp),
                            target_dims = list(target_dims_obs,
                                               c(memb_dim, sdate_dim)),
-                      fun = .sbc , na.rm = na.rm, ncores = ncores)$output1
+                           fun = .sbc, 
+                           na.rm = na.rm, ncores = ncores)$output1
   } else {
     BiasCorrected <- Apply(data = list(var_obs = obs,
                                        var_exp = exp,
@@ -151,29 +212,28 @@ BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE,
                            target_dims = list(target_dims_obs,
                                               c(memb_dim, sdate_dim),
                                               c(memb_dim, sdate_dim)),
-                      fun = .sbc , output_dims = c(memb_dim,sdate_dim), 
-                      na.rm = na.rm, ncores = ncores)$output1
+                           fun = .sbc, 
+                           output_dims = c(memb_dim, sdate_dim),
+                           na.rm = na.rm, ncores = ncores)$output1
   }
   return(BiasCorrected)
 }
 
 .sbc <- function(var_obs, var_exp , var_cor = NULL, na.rm = FALSE) {
-  nmembers <- dim(var_exp)[1]
-  ntime <- dim(var_exp)[2][]
-  
+
+  ntime <- dim(var_exp)[2]
   corrected <- NA * var_exp
 
-  if (is.null(var_cor)){
-    for (t in 1 : ntime) {
+  if (is.null(var_cor)) {
+    for (t in 1:ntime) {
       # parameters
       sd_obs <- sd(var_obs[-t], na.rm = na.rm)
-      sd_exp <- sd(var_exp[ , -t], na.rm = na.rm)
-      clim_exp <- mean(var_exp[ , -t], na.rm = na.rm)
+      sd_exp <- sd(var_exp[, -t], na.rm = na.rm)
+      clim_exp <- mean(var_exp[, -t], na.rm = na.rm)
       clim_obs <- mean(var_obs[-t], na.rm = na.rm)
       
       # bias corrected forecast
-      corrected[ , t] <- ((var_exp[ , t] - clim_exp) * (sd_obs / sd_exp)) + clim_obs
-      names(dim(corrected)) <- c('member', 'sdate')
+      corrected[, t] <- ((var_exp[, t] - clim_exp) * (sd_obs / sd_exp)) + clim_obs
     }
   } else {
     # parameters
@@ -184,7 +244,7 @@ BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE,
     
     # bias corrected forecast
     corrected <- ((var_cor - clim_exp) * (sd_obs / sd_exp)) + clim_obs
-    names(dim(corrected)) <- c('member')
   }
+
   return(corrected)
 }
diff --git a/man/BiasCorrection.Rd b/man/BiasCorrection.Rd
index 3f8d72ce279759be58c6a7fc2b3fff12b615fe62..a00c3461af6aa44dff200f06070536966cd707a0 100644
--- a/man/BiasCorrection.Rd
+++ b/man/BiasCorrection.Rd
@@ -11,32 +11,43 @@ BiasCorrection(
   na.rm = FALSE,
   memb_dim = "member",
   sdate_dim = "sdate",
-  ncores = 1
+  ncores = NULL
 )
 }
 \arguments{
-\item{exp}{a multidimensional array with named dimensions containing the seasonal forecast experiment data with at least 'member' and 'sdate' dimensions.}
+\item{exp}{A multidimensional array with named dimensions containing the 
+seasonal forecast experiment data with at least 'member' and 'sdate' 
+dimensions.}
 
-\item{obs}{a multidimensional array with named dimensions containing the observed data with at least 'sdate' dimension.}
+\item{obs}{A multidimensional array with named dimensions containing the 
+observed data with at least 'sdate' dimension.}
 
-\item{exp_cor}{a multidimensional array with named dimensions containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.}
+\item{exp_cor}{A multidimensional array with named dimensions containing the 
+seasonl forecast experiment to be corrected. If it is NULL, the 'exp' 
+forecast will be corrected.}
 
-\item{na.rm}{a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.}
+\item{na.rm}{A logical value indicating whether missing values should be 
+stripped before the computation proceeds, by default it is set to FALSE.}
 
-\item{memb_dim}{is a character string indicating the name of the member dimension. By default, it is set to 'member'.}
+\item{memb_dim}{A character string indicating the name of the member 
+dimension. By default, it is set to 'member'.}
 
-\item{sdate_dim}{is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.}
+\item{sdate_dim}{A character string indicating the name of the start date 
+dimension. By default, it is set to 'sdate'.}
 
-\item{ncores}{is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.}
+\item{ncores}{An integer that indicates the number of cores for parallel 
+computations using multiApply function. The default value is NULL.}
 }
 \value{
-an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data.
+An array containing the bias corrected forecasts with the same 
+dimensions of the experimental data.
 }
 \description{
-This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset.
+This function applies the simple bias adjustment technique 
+described in Torralba et al. (2017). The adjusted forecasts have an equivalent 
+standard deviation and mean to that of the reference dataset.
 }
 \examples{
-
 # Example
 # Creation of sample s2dverification objects. These are not complete
 # s2dverification objects though. The Load function returns complete objects.
@@ -48,7 +59,11 @@ a <- BiasCorrection(exp = mod1, obs = obs1)
 str(a)
 }
 \references{
-Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS)
+Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. 
+Davis (2017). Seasonal climate prediction: a new source of information for the 
+management of wind energy resources. Journal of Applied Meteorology and 
+Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, 
+EUPORIAS, NEWA, RESILIENCE, SPECS)
 }
 \author{
 Verónica Torralba, \email{veronica.torralba@bsc.es}
diff --git a/man/CST_BiasCorrection.Rd b/man/CST_BiasCorrection.Rd
index 202e7651d661797db14514aa5e69835e4e6c39fd..45014f7d57095f53e9ae0dd5cbf6e3ea94c6b68e 100644
--- a/man/CST_BiasCorrection.Rd
+++ b/man/CST_BiasCorrection.Rd
@@ -11,32 +11,43 @@ CST_BiasCorrection(
   na.rm = FALSE,
   memb_dim = "member",
   sdate_dim = "sdate",
-  ncores = 1
+  ncores = NULL
 )
 }
 \arguments{
-\item{exp}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonal forecast experiment data in the element named \code{$data}}
+\item{exp}{An object of class \code{s2dv_cube} as returned by \code{CST_Load} 
+function, containing the seasonal forecast experiment data in the element 
+named \code{$data}}
 
-\item{obs}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the observed data in the element named \code{$data}.}
+\item{obs}{An object of class \code{s2dv_cube} as returned by \code{CST_Load} 
+function, containing the observed data in the element named \code{$data}.}
 
-\item{exp_cor}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.}
+\item{exp_cor}{An object of class \code{s2dv_cube} as returned by 
+\code{CST_Load} function, containing the seasonl forecast experiment to be 
+corrected. If it is NULL, the 'exp' forecast will be corrected.}
 
-\item{na.rm}{a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.}
+\item{na.rm}{A logical value indicating whether missing values should be 
+stripped before the computation proceeds, by default it is set to FALSE.}
 
-\item{memb_dim}{is a character string indicating the name of the member dimension. By default, it is set to 'member'.}
+\item{memb_dim}{A character string indicating the name of the member 
+dimension. By default, it is set to 'member'.}
 
-\item{sdate_dim}{is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.}
+\item{sdate_dim}{A character string indicating the name of the start date 
+dimension. By default, it is set to 'sdate'.}
 
-\item{ncores}{is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.}
+\item{ncores}{An integer that indicates the number of cores for parallel 
+computations using multiApply function. The default value is NULL.}
 }
 \value{
-an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data.
+An object of class \code{s2dv_cube} containing the bias corrected 
+forecasts with the same dimensions of the experimental data.
 }
 \description{
-This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset.
+This function applies the simple bias adjustment technique 
+described in Torralba et al. (2017). The adjusted forecasts have an equivalent 
+standard deviation and mean to that of the reference dataset.
 }
 \examples{
-
 # Example
 # Creation of sample s2dverification objects. These are not complete
 # s2dverification objects though. The Load function returns complete objects.
@@ -54,7 +65,11 @@ a <- CST_BiasCorrection(exp = exp, obs = obs)
 str(a)
 }
 \references{
-Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS)
+Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. 
+Davis (2017). Seasonal climate prediction: a new source of information for 
+the management of wind energy resources. Journal of Applied Meteorology and 
+Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, 
+EUPORIAS, NEWA, RESILIENCE, SPECS)
 }
 \author{
 Verónica Torralba, \email{veronica.torralba@bsc.es}
diff --git a/tests/testthat/test-CST_BiasCorrection.R b/tests/testthat/test-CST_BiasCorrection.R
index 205aef445b17bcbab8457e209961e98e776a7720..826fcf127806fb77413e7d5cd7e78c49d7c09359 100644
--- a/tests/testthat/test-CST_BiasCorrection.R
+++ b/tests/testthat/test-CST_BiasCorrection.R
@@ -1,60 +1,184 @@
-context("Generic tests")
-test_that("Sanity checks", {
+context("CSTools::CST_BiasCorrection tests")
+
+##############################################
+
+# dat1
+mod <- 1 : (1 * 3 * 4 * 5 * 6 * 7)
+obs <- 1 : (1 * 1 * 4 * 5 * 6 * 7)
+dim(mod) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, 
+                lat = 6, lon = 7)
+dim(obs) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, 
+                lat = 6, lon = 7)
+lon <- seq(0, 30, 5)
+lat <- seq(0, 25, 5)
+exp <- list(data = mod, lat = lat, lon = lon)
+obs <- list(data = obs, lat = lat, lon = lon)
+attr(exp, 'class') <- 's2dv_cube'
+attr(obs, 'class') <- 's2dv_cube'
+
+exp1 <- list(data = array(1:20, dim = c(time = 20)))
+class(exp1) <- 's2dv_cube'
+
+obs1 <- list(data = array(1:20, dim = c(time = 20)))
+class(obs1) <- 's2dv_cube'
+
+exp1_2 <- list(data = array(1:20, dim = c(20)))
+class(exp1_2) <- 's2dv_cube'
+
+obs1_2 <- list(data = array(1:20, dim = c(20)))
+class(obs1_2) <- 's2dv_cube'
+
+exp_cor1 <- list(data = array(1:20, dim = c(20)))
+class(exp_cor1) <- 's2dv_cube'
+
+# dat2
+exp2 <- exp
+exp2$data[1, 2, 1, 1, 1, 1] <- NA
+obs2 <- obs
+obs2$data[1, 1, 2, 1, 1, 1] <- NA
+
+# dat3
+exp3 <- array(1:6, c(sdate = 3, member = 2))
+obs3 <- array(3:6, c(sdate = 3, member = 1))
+obs3_2 <- array(3:6, c(sdate = 3))
+obs3_3 <- array(3:6, c(sdate = 3, member = 2))
+
+# dat4
+exp4 <- array(1:100, dim = c(time = 5, members = 5, lat = 2, lon = 5))
+obs4 <- array(1:200, dim = c(time = 5, members = 1, lat = 2, lon = 5))
+obs4_1 <- obs4
+obs4_1[1,1,1,1] <- NA
+
+##############################################
+
+test_that("1. Inpput checks", {
+  # s2dv_cube
   expect_error(
     CST_BiasCorrection(exp = 1),
     paste0("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ",
-         "as output by CSTools::CST_Load."))
-  
-  mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7)
-  obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7)
-  dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, 
-                 lat = 6, lon = 7)
-  dim(obs1) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, 
-                 lat = 6, lon = 7)
-  lon <- seq(0, 30, 5)
-  lat <- seq(0, 25, 5)
-  exp <- list(data = mod1, lat = lat, lon = lon)
-  obs <- list(data = obs1, lat = lat, lon = lon)
-  attr(exp, 'class') <- 's2dv_cube'
-  attr(obs, 'class') <- 's2dv_cube'
-
-  bc <- CST_BiasCorrection(exp = exp, obs = obs)  
-  expect_equal(length(bc), 3)
-  expect_equal(dim(bc$data), 
-               c(dataset = 1, member = 3, sdate = 4, ftime = 5, 
-                 lat = 6, lon = 7))
-  expect_equal(bc$lat, lat)
-  expect_equal(bc$lon, lon)
-  
-  expect_error(CST_BiasCorrection(exp = exp, obs = exp),
-        paste0("The length of the dimension 'member' in the component 'data' ",
-         "of the parameter 'obs' must be equal to 1."))
-  
-  exp2 <- exp
-  exp2$data[1, 2, 1, 1, 1, 1] <- NA
-  expect_warning(CST_BiasCorrection(exp = exp2, obs = obs),
-    "Parameter 'exp' contains NA values.")
-  
-  obs2 <- obs
-  obs2$data[1, 1, 2, 1, 1, 1] <- NA
-  expect_warning(CST_BiasCorrection(exp = exp, obs = obs2),
-    "Parameter 'obs' contains NA values.")
-  
-  expect_warning(CST_BiasCorrection(exp = exp2, obs = obs2),
-    "Parameter 'obs' contains NA values", "Parameter 'exp' contains NA values.")
-
-  hinc <- array(1:6, c(sdate = 3, member = 2))
-  obs <-  array(3:6, c(sdate = 3, member = 1))
-  res <- round(BiasCorrection(exp = hinc, obs = obs, exp_cor = hinc), 2)
-  expect_equal(res, array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3)))
-
-  # if obs doesn't have memb_dim it works the same:
-  hinc <- array(1:6, c(sdate = 3, member = 2))
-  obs <-  array(3:6, c(sdate = 3))
-  res <- round(BiasCorrection(exp = hinc, obs = obs, exp_cor = hinc), 2)
-  expect_equal(res, array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3)))
-
-
+         "as output by CSTools::CST_Load.")
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp, obs = 1),
+    paste0("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ",
+           "as output by CSTools::CST_Load.")
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp1),
+    'argument "obs" is missing, with no default'
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp1, obs = obs1, exp_cor = 1),
+    paste0("Parameter 'exp_cor' must be of the class 's2dv_cube', as output ",
+    "by CSTools::CST_Load.")
+  )
+  # exp and obs
+  expect_error(
+    CST_BiasCorrection(exp = exp1_2, obs = obs1),
+    "Parameter 'exp' must have dimension names."
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp1, obs = obs1_2),
+    "Parameter 'obs' must have dimension names."
+  )
+  expect_warning(
+    CST_BiasCorrection(exp = exp2, obs = obs2),
+    "Parameter 'exp' contains NA values."
+  )
+  expect_warning(
+    CST_BiasCorrection(exp = exp, obs = obs2),
+    "Parameter 'obs' contains NA values."
+  )
+  expect_warning(
+    CST_BiasCorrection(exp = exp2, obs = obs2),
+    "Parameter 'obs' contains NA values",
+    "Parameter 'exp' contains NA values."
+    )
+  # exp_cor
+  expect_error(
+    CST_BiasCorrection(exp = exp1, obs = obs1, exp_cor = exp_cor1, sdate_dim = 'time'),
+    "Parameter 'exp_cor' must have dimension names."
+  )
+  # sdate_dim, memb_dim
+  expect_error(
+    CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 1),
+    paste0("Parameter 'sdate_dim' must be a character string.")
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp, obs = obs, sdate_dim = 'time'),
+    paste0("Parameter 'sdate_dim' is not found in 'exp' or 'obs' dimension.")
+  )
+  expect_error(
+    BiasCorrection(exp = array(1:20, dim = c(time = 1, member = 1)), 
+                   obs = array(1:20, dim = c(time = 2, member = 1)), sdate_dim = 'time'),
+    paste0("Parameter 'exp' must have dimension length of 'sdate_dim' bigger than 1.")
+  )
+  expect_error(
+    CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 'time'),
+    paste0("Parameter 'exp' requires 'sdate_dim' and 'memb_dim' dimensions.")
+  )
+  expect_error(
+    BiasCorrection(exp = exp3, obs = obs3_3),
+    paste0("If parameter 'obs' has dimension 'memb_dim' its length must be equal to 1.")
+  )
+  ## na.rm
+  expect_warning(
+    CST_BiasCorrection(exp = exp, obs = obs, na.rm = 1),
+    "Paramater 'na.rm' must be a logical, it has been set to FALSE."
+  )
+  expect_warning(
+    CST_BiasCorrection(exp = exp, obs = obs, na.rm = c(T,F)),
+    "Paramter 'na.rm' has length greater than 1, and only the fist element is used."
+  )
+  # ncores
+  expect_error(
+    CST_BiasCorrection(exp = exp, obs = obs, ncores = TRUE),
+    "Parameter 'ncores' must be either NULL or a positive integer."
+  )
+})
 
+##############################################
 
+test_that("2. Output checks: dat1", {
+  bc <- CST_BiasCorrection(exp = exp, obs = obs)
+  expect_equal(
+    length(bc),
+    3
+  )
+  expect_equal(
+    dim(bc$data),
+    c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7)
+  )
+  expect_equal(
+    bc$lat,
+    lat
+  )
+  expect_equal(
+    bc$lon,
+    lon
+  )
+  expect_equal(
+    round(BiasCorrection(exp = exp3, obs = obs3, exp_cor = exp3), 2), 
+    array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3))
+  )
+  expect_equal(
+    round(BiasCorrection(exp = exp3, obs = obs3_2, exp_cor = exp3), 2),
+    array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3))
+  )
+  expect_equal(
+    dim(BiasCorrection(exp = exp4, obs = obs4, sdate_dim = 'time', memb_dim = 'members')),
+    c(members = 5, time = 5, lat = 2, lon = 5)
+  )
+  suppressWarnings(
+    expect_equal(
+      sum(is.na(BiasCorrection(exp = exp4, obs = obs4_1, sdate_dim = 'time', memb_dim = 'members', na.rm = TRUE))),
+      0
+    )
+  )
+  suppressWarnings(
+    expect_equal(
+      sum(is.na(BiasCorrection(exp = exp4, obs = obs4_1, sdate_dim = 'time', memb_dim = 'members', na.rm = FALSE))),
+      20
+    )
+  )
 })