diff --git a/R/RunModel_CemaNeigeGR4J.R b/R/RunModel_CemaNeigeGR4J.R
index 17cad7d180495651ea5dca3b97923b0f87b03bb4..1441261e144c1bf8edbe0933fa9c98f8c580f974 100644
--- a/R/RunModel_CemaNeigeGR4J.R
+++ b/R/RunModel_CemaNeigeGR4J.R
@@ -1,70 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the CemaNeige-GR4J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the CemaNeigeGR4J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Perrin, C., C. Michel and V. Andréassian (2003), 
-#'       Improvement of a parsimonious model for streamflow simulation, 
-#'       Journal of Hydrology, 279(1-4), 275-289, doi:10.1016/S0022-1694(03)00225-7. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": what is useful in a temperature-based snow-accounting routine? 
-#'       Part 1 - Comparison of six snow accounting routines on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.059. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": What is useful in a temperature-based snow-accounting routine? 
-#'       Part 2 - Sensitivity analysis of the Cemaneige snow accounting routine on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.058.
-#' @seealso \code{\link{RunModel_CemaNeigeGR5J}}, \code{\link{RunModel_CemaNeigeGR6J}}, \code{\link{RunModel_GR4J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_CemaNeigeGR4J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 6 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR4J X1      \tab production store capacity [mm]                                \cr
-#'                             GR4J X2      \tab intercatchment exchange coefficient [mm/d]                    \cr
-#'                             GR4J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR4J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             CemaNeige X1 \tab weighting coefficient for snow pack thermal state [-]         \cr
-#'                             CemaNeige X2 \tab degree-day melt coefficient [mm/degC/d]                       \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$CemaNeigeLayers}   \tab [list] list of CemaNeige outputs (1 list per layer)                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Pliq         }   \tab [numeric] series of liquid precip. [mm/d]                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Psol         }   \tab [numeric] series of solid precip. [mm/d]                           \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$SnowPack     }   \tab [numeric] series of snow pack [mm]                                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$ThermalState }   \tab [numeric] series of snow pack thermal state [degC]                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Gratio       }   \tab [numeric] series of Gratio [0-1]                                   \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PotMelt      }   \tab [numeric] series of potential snow melt [mm/d]                     \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Melt         }   \tab [numeric] series of actual snow melt [mm/d]                        \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PliqAndMelt  }   \tab [numeric] series of liquid precip. + actual snow melt [mm/d]       \cr
-#'          \emph{$StateEnd}                                  \tab [numeric] states at the end of the run: \cr\tab res. & HU levels [mm], CemaNeige states [mm & degC] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************
 RunModel_CemaNeigeGR4J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 6;
@@ -97,7 +30,6 @@ RunModel_CemaNeigeGR4J <- function(InputsModel,RunOptions,Param){
       ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim;
 
 
-
     ##SNOW_MODULE________________________________________________________________________________##
     if(RunOptions$RunSnowModule==TRUE){
       if("all" %in% RunOptions$Outputs_Sim){ IndOutputsCemaNeige <- as.integer(1:length(FortranOutputsCemaNeige)); 
@@ -107,7 +39,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel,RunOptions,Param){
     ##Call_DLL_CemaNeige_________________________
       for(iLayer in 1:NLayers){
         StateStartCemaNeige <- RunOptions$IniStates[ (NStatesMod+2*(iLayer-1)+1):(NStatesMod+2*(iLayer-1)+2) ];
-        RESULTS <- .Fortran("frun_cemaneige",PACKAGE="airGR",
+        RESULTS <- .Fortran("frun_CemaNeige",PACKAGE="airGR",
                         ##inputs
                             LInputs=LInputSeries,                                                          ### length of input and output series
                             InputsPrecip=InputsModel$LayerPrecip[[iLayer]][IndPeriod1],                    ### input series of total precipitation [mm/d]
@@ -150,12 +82,12 @@ RunModel_CemaNeigeGR4J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% RunOptions){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr4j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR4J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                         ### length of input and output series
                      InputsPrecip=CatchMeltAndPliq,                ### input series of total precipitation [mm/d]
@@ -168,7 +100,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel,RunOptions,Param){
                      IndOutputs=IndOutputsMod,                     ### indices of output series
                  ##outputs                                        
                      Outputs=matrix(as.double(-999.999),nrow=LInputSeries,ncol=length(IndOutputsMod)),  ### output series [mm]
-                     StateEnd=rep(as.double(-999.999),length(RunOptions$IniStates))                     ### state variables at the end of the model run
+                     StateEnd=rep(as.double(-999.999),NStatesMod)                     ### state variables at the end of the model run
                      )
       RESULTS$Outputs[ round(RESULTS$Outputs ,3)==(-999.999)] <- NA;
       RESULTS$StateEnd[round(RESULTS$StateEnd,3)==(-999.999)] <- NA;
diff --git a/R/RunModel_CemaNeigeGR5J.R b/R/RunModel_CemaNeigeGR5J.R
index d752403881a2c3997cbef2da564ac3a43628cd26..a65c664cead6913789076d83f294a1ea8d832b82 100644
--- a/R/RunModel_CemaNeigeGR5J.R
+++ b/R/RunModel_CemaNeigeGR5J.R
@@ -1,73 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the CemaNeige-GR5J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the CemaNeigeGR5J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Le Moine, N. (2008), Le bassin versant de surface vu par le souterrain : une voie d'amélioration des performances 
-#'       et du réalisme des modèles pluie-débit ?, PhD thesis (french), UPMC, Paris, France. \cr
-#'   Pushpalatha, R., C. Perrin, N. Le Moine, T. Mathevet and V. Andréassian (2011), 
-#'       A downward structural sensitivity analysis of hydrological models to improve low-flow simulation, 
-#'       Journal of Hydrology, 411(1-2), 66-76, doi:10.1016/j.jhydrol.2011.09.034. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": what is useful in a temperature-based snow-accounting routine? 
-#'       Part 1 - Comparison of six snow accounting routines on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.059. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": What is useful in a temperature-based snow-accounting routine? 
-#'       Part 2 - Sensitivity analysis of the Cemaneige snow accounting routine on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.058.
-#' @seealso \code{\link{RunModel_CemaNeigeGR4J}}, \code{\link{RunModel_CemaNeigeGR6J}}, \code{\link{RunModel_GR5J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_CemaNeigeGR5J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 7 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR5J X1      \tab production store capacity [mm]                                \cr
-#'                             GR5J X2      \tab intercatchment exchange coefficient 1 [mm/d]                  \cr
-#'                             GR5J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR5J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             GR5J X5      \tab intercatchment exchange coefficient 2 [-]                     \cr
-#'                             CemaNeige X1 \tab weighting coefficient for snow pack thermal state [-]         \cr
-#'                             CemaNeige X2 \tab degree-day melt coefficient [mm/degC/d]                       \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$CemaNeigeLayers}   \tab [list] list of CemaNeige outputs (1 list per layer)                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Pliq         }   \tab [numeric] series of liquid precip. [mm/d]                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Psol         }   \tab [numeric] series of solid precip. [mm/d]                           \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$SnowPack     }   \tab [numeric] series of snow pack [mm]                                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$ThermalState }   \tab [numeric] series of snow pack thermal state [degC]                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Gratio       }   \tab [numeric] series of Gratio [0-1]                                   \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PotMelt      }   \tab [numeric] series of potential snow melt [mm/d]                     \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Melt         }   \tab [numeric] series of actual snow melt [mm/d]                        \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PliqAndMelt  }   \tab [numeric] series of liquid precip. + actual snow melt [mm/d]       \cr
-#'          \emph{$StateEnd}                                  \tab [numeric] states at the end of the run: \cr\tab res. & HU levels [mm], CemaNeige states [mm & degC] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************
 RunModel_CemaNeigeGR5J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 7;
@@ -109,7 +39,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel,RunOptions,Param){
     ##Call_DLL_CemaNeige_________________________
       for(iLayer in 1:NLayers){
         StateStartCemaNeige <- RunOptions$IniStates[ (NStatesMod+2*(iLayer-1)+1):(NStatesMod+2*(iLayer-1)+2) ];
-        RESULTS <- .Fortran("frun_cemaneige",PACKAGE="airGR",
+        RESULTS <- .Fortran("frun_CemaNeige",PACKAGE="airGR",
                         ##inputs
                             LInputs=LInputSeries,                                                          ### length of input and output series
                             InputsPrecip=InputsModel$LayerPrecip[[iLayer]][IndPeriod1],                    ### input series of total precipitation [mm/d]
@@ -152,12 +82,12 @@ RunModel_CemaNeigeGR5J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% RunOptions){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr5j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR5J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                         ### length of input and output series
                      InputsPrecip=CatchMeltAndPliq,                ### input series of total precipitation [mm/d]
@@ -170,7 +100,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel,RunOptions,Param){
                      IndOutputs=IndOutputsMod,                     ### indices of output series
                  ##outputs                                        
                      Outputs=matrix(as.double(-999.999),nrow=LInputSeries,ncol=length(IndOutputsMod)),  ### output series [mm]
-                     StateEnd=rep(as.double(-999.999),length(RunOptions$IniStates))                     ### state variables at the end of the model run
+                     StateEnd=rep(as.double(-999.999),NStatesMod)                     ### state variables at the end of the model run
                      )
       RESULTS$Outputs[ round(RESULTS$Outputs ,3)==(-999.999)] <- NA;
       RESULTS$StateEnd[round(RESULTS$StateEnd,3)==(-999.999)] <- NA;
diff --git a/R/RunModel_CemaNeigeGR6J.R b/R/RunModel_CemaNeigeGR6J.R
index f32b0072ceeabddb00b867ba94f290f341324482..73d9e75706216b18a4602eead1f039e7cd6ff098 100644
--- a/R/RunModel_CemaNeigeGR6J.R
+++ b/R/RunModel_CemaNeigeGR6J.R
@@ -1,74 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the CemaNeige-GR6J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the CemaNeigeGR6J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Pushpalatha, R., C. Perrin, N. Le Moine, T. Mathevet and V. Andréassian (2011), 
-#'       A downward structural sensitivity analysis of hydrological models to improve low-flow simulation, 
-#'       Journal of Hydrology, 411(1-2), 66-76, doi:10.1016/j.jhydrol.2011.09.034. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": what is useful in a temperature-based snow-accounting routine? 
-#'       Part 1 - Comparison of six snow accounting routines on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.059. \cr
-#'   Valéry, A., V. Andréassian and C. Perrin (2014), 
-#'       "As simple as possible but not simpler": What is useful in a temperature-based snow-accounting routine? 
-#'       Part 2 - Sensitivity analysis of the Cemaneige snow accounting routine on 380 catchments, Journal of Hydrology, doi:10.1016/j.jhydrol.2014.04.058.
-#' @seealso \code{\link{RunModel_CemaNeigeGR4J}}, \code{\link{RunModel_CemaNeigeGR5J}}, \code{\link{RunModel_GR6J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_CemaNeigeGR6J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 8 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR6J X1      \tab production store capacity [mm]                                \cr
-#'                             GR6J X2      \tab intercatchment exchange coefficient 1 [mm/d]                  \cr
-#'                             GR6J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR6J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             GR6J X5      \tab intercatchment exchange coefficient 2 [-]                     \cr
-#'                             GR6J X6      \tab coefficient for emptying exponential store [-]                \cr
-#'                             CemaNeige X1 \tab weighting coefficient for snow pack thermal state [-]         \cr
-#'                             CemaNeige X2 \tab degree-day melt coefficient [mm/degC/d]                       \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QR1     }          \tab [numeric] series of exponential store outflow (QR1) [mm/d]                    \cr
-#'          \emph{$Exp     }          \tab [numeric] series of exponential store level (X(6)) (negative) [mm]            \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$CemaNeigeLayers}   \tab [list] list of CemaNeige outputs (1 list per layer)                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Pliq         }   \tab [numeric] series of liquid precip. [mm/d]                          \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Psol         }   \tab [numeric] series of solid precip. [mm/d]                           \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$SnowPack     }   \tab [numeric] series of snow pack [mm]                                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$ThermalState }   \tab [numeric] series of snow pack thermal state [degC]                 \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Gratio       }   \tab [numeric] series of Gratio [0-1]                                   \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PotMelt      }   \tab [numeric] series of potential snow melt [mm/d]                     \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$Melt         }   \tab [numeric] series of actual snow melt [mm/d]                        \cr
-#'          \emph{$CemaNeigeLayers[[iLayer]]$PliqAndMelt  }   \tab [numeric] series of liquid precip. + actual snow melt [mm/d]       \cr
-#'          \emph{$StateEnd}                                  \tab [numeric] states at the end of the run: \cr\tab res. & HU levels [mm], CemaNeige states [mm & degC] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************
 RunModel_CemaNeigeGR6J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 8;
@@ -111,7 +40,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel,RunOptions,Param){
     ##Call_DLL_CemaNeige_________________________
       for(iLayer in 1:NLayers){
         StateStartCemaNeige <- RunOptions$IniStates[ (NStatesMod+2*(iLayer-1)+1):(NStatesMod+2*(iLayer-1)+2) ];
-        RESULTS <- .Fortran("frun_cemaneige",PACKAGE="airGR",
+        RESULTS <- .Fortran("frun_CemaNeige",PACKAGE="airGR",
                         ##inputs
                             LInputs=LInputSeries,                                                          ### length of input and output series
                             InputsPrecip=InputsModel$LayerPrecip[[iLayer]][IndPeriod1],                    ### input series of total precipitation [mm/d]
@@ -154,12 +83,12 @@ RunModel_CemaNeigeGR6J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% RunOptions){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*ParamMod[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*ParamMod[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr6j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR6J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                         ### length of input and output series
                      InputsPrecip=CatchMeltAndPliq,                ### input series of total precipitation [mm/d]
@@ -172,7 +101,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel,RunOptions,Param){
                      IndOutputs=IndOutputsMod,                     ### indices of output series
                  ##outputs                                        
                      Outputs=matrix(as.double(-999.999),nrow=LInputSeries,ncol=length(IndOutputsMod)),  ### output series [mm]
-                     StateEnd=rep(as.double(-999.999),length(RunOptions$IniStates))                     ### state variables at the end of the model run
+                     StateEnd=rep(as.double(-999.999),NStatesMod)                     ### state variables at the end of the model run
                      )
       RESULTS$Outputs[ round(RESULTS$Outputs ,3)==(-999.999)] <- NA;
       RESULTS$StateEnd[round(RESULTS$StateEnd,3)==(-999.999)] <- NA;
diff --git a/R/RunModel_GR2M.R b/R/RunModel_GR2M.R
index 005a1b00739699ad4e87e05a9d61d404245e4b92..c0415938cbdf6b619cde71d86c0c467f39547c13 100644
--- a/R/RunModel_GR2M.R
+++ b/R/RunModel_GR2M.R
@@ -1,50 +1,11 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the GR2M monthly lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the GR2M hydrological model
-#' @author Laurent Coron (March 2015)
-#' @example tests/example_RunModel_GR2M.R
-#' @references
-#'   Mouelhi S. (2003), 
-#'       Vers une chaîne cohérente de modèles pluie-débit conceptuels globaux aux pas de temps pluriannuel, annuel, mensuel et journalier,
-#'       PhD thesis (in French), ENGREF, Cemagref Antony, France. \cr
-#'   Mouelhi, S., C. Michel, C. Perrin and V. Andréassian (2006),
-#'       Stepwise development of a two-parameter monthly water balance model,
-#'       Journal of Hydrology, 318(1-4), 200-214, doi:10.1016/j.jhydrol.2005.06.014.
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 2 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR2M X1      \tab production store capacity [mm]                 \cr
-#'                             GR2M X2      \tab groundwater exchange coefficient [-]    \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                    \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/month]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/month]                          \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/month]                                               \cr
-#'          \emph{$StateEnd}          \tab [numeric] states at the end of the run (res. levels, HU1 levels, HU2 levels) [mm] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************
 RunModel_GR2M <- function(InputsModel,RunOptions,Param){
 
     NParam <- 2;
-    FortranOutputs <- c("PotEvap","Precip","Prod","Rout","Qsim");
-    ### FortranOutputs <- c("PotEvap","Precip","Prod","AE","Perc","PR","Q9","Q1","Rout","Exch","AExch","QR","QD","Qsim");
+    FortranOutputs <- c("PotEvap","Precip","AE","Perc","P3","Exch","Prod","Rout","Qsim");
 
     ##Arguments_check
       if(inherits(InputsModel,"InputsModel")==FALSE){ stop("InputsModel must be of class 'InputsModel' \n"); return(NULL); }  
-      if(inherits(InputsModel,"monthly"    )==FALSE){ stop("InputsModel must be of class 'monthly'     \n"); return(NULL); }  
+      if(inherits(InputsModel,"monthly"    )==FALSE){ stop("InputsModel must be of class 'monthly'      \n"); return(NULL); }  
       if(inherits(InputsModel,"GR"         )==FALSE){ stop("InputsModel must be of class 'GR'          \n"); return(NULL); }  
       if(inherits(RunOptions,"RunOptions"  )==FALSE){ stop("RunOptions must be of class 'RunOptions'   \n"); return(NULL); }  
       if(inherits(RunOptions,"GR"          )==FALSE){ stop("RunOptions must be of class 'GR'           \n"); return(NULL); }  
@@ -65,7 +26,7 @@ RunModel_GR2M <- function(InputsModel,RunOptions,Param){
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr2m",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR2M",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                             ### length of input and output series
                      InputsPrecip=InputsModel$Precip[IndPeriod1],      ### input series of total precipitation [mm/month]
@@ -80,7 +41,7 @@ RunModel_GR2M <- function(InputsModel,RunOptions,Param){
                      Outputs=matrix(as.double(-999.999),nrow=LInputSeries,ncol=length(IndOutputs)),  ### output series [mm]
                      StateEnd=rep(as.double(-999.999),length(RunOptions$IniStates))                  ### state variables at the end of the model run
                      )
-      RESULTS$Outputs[ round(RESULTS$Outputs ,3)==(-999.999)] <- NA;
+      RESULTS$Outputs [round(RESULTS$Outputs ,3)==(-999.999)] <- NA;
       RESULTS$StateEnd[round(RESULTS$StateEnd,3)==(-999.999)] <- NA;
 
     ##Output_data_preparation
diff --git a/R/RunModel_GR4H.R b/R/RunModel_GR4H.R
index 1aa9ba31055165ef9a03757474a8cdf2a5104596..792642f6f833888b066c680a1000bbeddffa4c6c 100644
--- a/R/RunModel_GR4H.R
+++ b/R/RunModel_GR4H.R
@@ -1,53 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the GR4H hourly lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the GR4H hydrological model
-#' @author Laurent Coron (July 2014)
-#' @seealso \code{\link{RunModel_GR4J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_GR4H.R
-#' @references
-#'   Mathevet, T. (2005), 
-#'       Quels modèles pluie-débit globaux pour le pas de temps horaire ? Développement empirique et comparaison de modèles sur un large échantillon de bassins versants, 
-#'       PhD thesis (in French), ENGREF - Cemagref (Antony), Paris, France.
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 4 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR4H X1      \tab production store capacity [mm]                                \cr
-#'                             GR4H X2      \tab groundwater exchange coefficient [mm/h]                       \cr
-#'                             GR4H X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR4H X4      \tab unit hydrograph time constant [h]                             \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                    \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/h]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/h]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/h]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/h]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/h]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/h]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/h]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/h]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/h]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/h]                         \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/h]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/h]                                               \cr
-#'          \emph{$StateEnd}          \tab [numeric] states at the end of the run (res. levels, HU1 levels, HU2 levels) [mm] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************'
 RunModel_GR4H <- function(InputsModel,RunOptions,Param){
 
     NParam <- 4;
@@ -72,12 +22,12 @@ RunModel_GR4H <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% names(RunOptions)){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr4h",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR4H",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                             ### length of input and output series
                      InputsPrecip=InputsModel$Precip[IndPeriod1],      ### input series of total precipitation [mm/h]
diff --git a/R/RunModel_GR4J.R b/R/RunModel_GR4J.R
index 8ac0918b0ec99960256509b804c8fb4a93449e58..193dc91fa812ed61cd042884ae9c63d7d9b0c269 100644
--- a/R/RunModel_GR4J.R
+++ b/R/RunModel_GR4J.R
@@ -1,53 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the GR4J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the GR4J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Perrin, C., C. Michel and V. Andréassian (2003), 
-#'       Improvement of a parsimonious model for streamflow simulation, 
-#'       Journal of Hydrology, 279(1-4), 275-289, doi:10.1016/S0022-1694(03)00225-7.
-#' @seealso \code{\link{RunModel_GR5J}}, \code{\link{RunModel_GR6J}}, \code{\link{RunModel_CemaNeigeGR4J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_GR4J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 4 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR4J X1      \tab production store capacity [mm]                                \cr
-#'                             GR4J X2      \tab intercatchment exchange coefficient [mm/d]                    \cr
-#'                             GR4J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR4J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$StateEnd}          \tab [numeric] states at the end of the run (res. levels, HU1 levels, HU2 levels) [mm] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************'
 RunModel_GR4J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 4;
@@ -72,12 +22,12 @@ RunModel_GR4J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% names(RunOptions)){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr4j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR4J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                             ### length of input and output series
                      InputsPrecip=InputsModel$Precip[IndPeriod1],      ### input series of total precipitation [mm/d]
diff --git a/R/RunModel_GR5J.R b/R/RunModel_GR5J.R
index be26d20bdc6491f93e7d651c87e7cfc9ef74c684..63d0a9e58e89dfd9bdab6a03a06d5f22a18ebc5b 100644
--- a/R/RunModel_GR5J.R
+++ b/R/RunModel_GR5J.R
@@ -1,56 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the GR5J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the GR5J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Le Moine, N. (2008), Le bassin versant de surface vu par le souterrain : une voie d'amélioration des performances
-#'       et du réalisme des modèles pluie-débit ?, PhD thesis (french), UPMC, Paris, France. \cr
-#'   Pushpalatha, R., C. Perrin, N. Le Moine, T. Mathevet, and V. Andréassian (2011), 
-#'       A downward structural sensitivity analysis of hydrological models to improve low-flow simulation, 
-#'       Journal of Hydrology, 411(1-2), 66-76, doi:10.1016/j.jhydrol.2011.09.034. \cr
-#' @seealso \code{\link{RunModel_GR4J}}, \code{\link{RunModel_GR6J}}, \code{\link{RunModel_CemaNeigeGR5J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_GR5J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 5 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR5J X1      \tab production store capacity [mm]                                \cr
-#'                             GR5J X2      \tab intercatchment exchange coefficient 1 [mm/d]                  \cr
-#'                             GR5J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR5J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             GR5J X5      \tab intercatchment exchange coefficient 2 [-]                     \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$StateEnd}          \tab [numeric] states at the end of the run (res. levels, HU1 levels, HU2 levels) [mm] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************'
 RunModel_GR5J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 5;
@@ -75,12 +22,12 @@ RunModel_GR5J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% names(RunOptions)){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr5j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR5J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                             ### length of input and output series
                      InputsPrecip=InputsModel$Precip[IndPeriod1],      ### input series of total precipitation [mm/d]
diff --git a/R/RunModel_GR6J.R b/R/RunModel_GR6J.R
index a3cd68a217ccfaa8b1a99edba65f2072b1fb40e4..a3227689eb37b6e134aa1d524a627def4402804f 100644
--- a/R/RunModel_GR6J.R
+++ b/R/RunModel_GR6J.R
@@ -1,57 +1,3 @@
-#*****************************************************************************************************************
-#' Function which performs a single run for the GR6J daily lumped model.
-#'
-#' For further details on the model, see the references section.
-#' For further details on the argument structures and initialisation options, see \code{\link{CreateRunOptions}}.
-#*****************************************************************************************************************
-#' @title Run with the GR6J hydrological model
-#' @author Laurent Coron (December 2013)
-#' @references
-#'   Pushpalatha, R., C. Perrin, N. Le Moine, T. Mathevet and V. Andréassian (2011), 
-#'       A downward structural sensitivity analysis of hydrological models to improve low-flow simulation, 
-#'       Journal of Hydrology, 411(1-2), 66-76, doi:10.1016/j.jhydrol.2011.09.034. \cr
-#' @seealso \code{\link{RunModel_GR4J}}, \code{\link{RunModel_GR5J}}, \code{\link{RunModel_CemaNeigeGR6J}},
-#'          \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}.
-#' @example tests/example_RunModel_GR6J.R
-#' @useDynLib airGR
-#' @encoding UTF-8
-#' @export
-#_FunctionInputs__________________________________________________________________________________________________
-#' @param  InputsModel         [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details
-#' @param  RunOptions          [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details
-#' @param  Param               [numeric] vector of 6 parameters                                                             
-#'                             \tabular{ll}{                                                                      
-#'                             GR6J X1      \tab production store capacity [mm]                                \cr
-#'                             GR6J X2      \tab intercatchment exchange coefficient 1 [mm/d]                  \cr
-#'                             GR6J X3      \tab routing store capacity [mm]                                   \cr
-#'                             GR6J X4      \tab unit hydrograph time constant [d]                             \cr
-#'                             GR6J X5      \tab intercatchment exchange coefficient 2 [-]                     \cr
-#'                             GR6J X6      \tab coefficient for emptying exponential store [-]                \cr
-#'                             }                                                                                  
-#_FunctionOutputs_________________________________________________________________________________________________
-#' @return  [list] list containing the function outputs organised as follows:                                         
-#'          \tabular{ll}{                                                                                         
-#'          \emph{$DatesR  }          \tab [POSIXlt] series of dates                                                     \cr
-#'          \emph{$PotEvap }          \tab [numeric] series of input potential evapotranspiration [mm/d]                 \cr
-#'          \emph{$Precip  }          \tab [numeric] series of input total precipitation [mm/d]                          \cr
-#'          \emph{$Prod    }          \tab [numeric] series of production store level (X(2)) [mm]                        \cr
-#'          \emph{$AE      }          \tab [numeric] series of actual evapotranspiration [mm/d]                          \cr
-#'          \emph{$Perc    }          \tab [numeric] series of percolation (PERC) [mm/d]                                 \cr
-#'          \emph{$PR      }          \tab [numeric] series of PR=PN-PS+PERC [mm/d]                                      \cr
-#'          \emph{$Q9      }          \tab [numeric] series of HU1 outflow (Q9) [mm/d]                                   \cr
-#'          \emph{$Q1      }          \tab [numeric] series of HU2 outflow (Q1) [mm/d]                                   \cr
-#'          \emph{$Rout    }          \tab [numeric] series of routing store level (X(1)) [mm]                           \cr
-#'          \emph{$Exch    }          \tab [numeric] series of potential semi-exchange between catchments [mm/d]         \cr
-#'          \emph{$AExch   }          \tab [numeric] series of actual exchange between catchments (1+2) [mm/d]           \cr
-#'          \emph{$QR      }          \tab [numeric] series of routing store outflow (QR) [mm/d]                         \cr
-#'          \emph{$QR1     }          \tab [numeric] series of exponential store outflow (QR1) [mm/d]                    \cr
-#'          \emph{$Exp     }          \tab [numeric] series of exponential store level (X(6)) (negative) [mm]            \cr
-#'          \emph{$QD      }          \tab [numeric] series of direct flow from HU2 after exchange (QD) [mm/d]           \cr
-#'          \emph{$Qsim    }          \tab [numeric] series of Qsim [mm/d]                                               \cr
-#'          \emph{$StateEnd}          \tab [numeric] states at the end of the run (res. levels, HU1 levels, HU2 levels) [mm] \cr
-#'          }                                                                                                     
-#'          (refer to the provided references or to the package source code for further details on these model outputs)
-#*****************************************************************************************************************'
 RunModel_GR6J <- function(InputsModel,RunOptions,Param){
 
     NParam <- 6;
@@ -76,12 +22,12 @@ RunModel_GR6J <- function(InputsModel,RunOptions,Param){
 
     ##Use_of_IniResLevels
       if("IniResLevels" %in% names(RunOptions)){
-        RunOptions$IniStates[1] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
-        RunOptions$IniStates[2] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1];  ### production store level (mm)
+        RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3];  ### routing store level (mm)
       }
 
     ##Call_fortan
-      RESULTS <- .Fortran("frun_gr6j",PACKAGE="airGR",
+      RESULTS <- .Fortran("frun_GR6J",PACKAGE="airGR",
                  ##inputs
                      LInputs=LInputSeries,                             ### length of input and output series
                      InputsPrecip=InputsModel$Precip[IndPeriod1],      ### input series of total precipitation [mm/d]