diff --git a/DESCRIPTION b/DESCRIPTION
index 91026ab85edaef44f465a2e8b1ca5c0b911d0f14..bae6fed2f86d812bcff869984b75b4eb1dc62164 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: airGR
 Type: Package
 Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling
-Version: 1.6.3.64
+Version: 1.6.3.65
 Date: 2020-11-17
 Authors@R: c(
   person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")),
diff --git a/NEWS.md b/NEWS.md
index 37d17e1a1ad3ac1729e9558858fc4bf5235f879c..29cfe222b7368dc08d361050d11a9742a795c113 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -4,7 +4,7 @@
 
 
 
-### 1.6.3.64 Release Notes (2020-11-17)
+### 1.6.3.65 Release Notes (2020-11-17)
 
 #### New features
 
diff --git a/R/RunModel_CemaNeige.R b/R/RunModel_CemaNeige.R
index e1dbbf659969732b8fbeca38246e7db2d8aa8bf2..a9f35e0e189465dd24e23bb0b4904216e98b6766 100644
--- a/R/RunModel_CemaNeige.R
+++ b/R/RunModel_CemaNeige.R
@@ -8,7 +8,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
   FortranOutputsCemaNeige <- .FortranOutputs(GR = NULL, isCN = TRUE)$CN
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -39,7 +39,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
     time_mult <- 24
   }
 
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, as.integer(0))) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -52,7 +52,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   ParamCemaNeige <- Param
   NLayers        <- length(InputsModel$LayerPrecip)
   
@@ -71,7 +71,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
   NameCemaNeigeLayers <- "CemaNeigeLayers"
   
   
-  ## Call_DLL_CemaNeige_________________________
+  ## Call CemaNeige Fortran_________________________
   for (iLayer in 1:NLayers) {
     
     if (!IsHyst) {
@@ -105,7 +105,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
     
     
     
-    ## Data_storage
+    ## Data storage
     CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(CemaNeigeLayers[[iLayer]]) <- FortranOutputsCemaNeige[IndOutputsCemaNeige]
     if (ExportStateEnd) {
@@ -113,7 +113,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
     }
     rm(RESULTS)
     
-  } ### ENDFOR_iLayer
+  } ### ENDFOR iLayer
   
   names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
   
@@ -129,7 +129,7 @@ RunModel_CemaNeige <- function(InputsModel, RunOptions, Param) {
                                          verbose = FALSE)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- list(CemaNeigeLayers)
     names(OutputsModel) <- NameCemaNeigeLayers
diff --git a/R/RunModel_CemaNeigeGR4H.R b/R/RunModel_CemaNeigeGR4H.R
index dced731bdf80843a13bfdf3598976a6d6f7f06dd..fe29e9096e863e14283b54635006c341d13ab9dd 100644
--- a/R/RunModel_CemaNeigeGR4H.R
+++ b/R/RunModel_CemaNeigeGR4H.R
@@ -9,7 +9,7 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR4H", isCN = TRUE)
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -55,7 +55,7 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }      
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -71,7 +71,7 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   if (inherits(RunOptions, "CemaNeige")) {
     if ("all" %in% RunOptions$Outputs_Sim) {
       IndOutputsCemaNeige <- as.integer(1:length(FortranOutputs$CN))
@@ -83,7 +83,7 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
     NameCemaNeigeLayers <- "CemaNeigeLayers"
     
     
-    ## Call_DLL_CemaNeige_________________________
+    ## Call CemaNeige Fortran_________________________
     for (iLayer in 1:NLayers) {
       
       if (!IsHyst) {
@@ -112,7 +112,7 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
       RESULTS$Outputs[ round(RESULTS$Outputs , 3) == -999.999] <- NA
       RESULTS$StateEnd[round(RESULTS$StateEnd, 3) == -999.999] <- NA
       
-      ## Data_storage
+      ## Data storage
       CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
       names(CemaNeigeLayers[[iLayer]]) <- FortranOutputs$CN[IndOutputsCemaNeige]
       IndPliqAndMelt <- which(names(CemaNeigeLayers[[iLayer]]) == "PliqAndMelt")
@@ -126,9 +126,9 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
         CemaNeigeStateEnd <- c(CemaNeigeStateEnd, RESULTS$StateEnd)
       }
       rm(RESULTS)
-    } ### ENDFOR_iLayer
+    } ### ENDFOR iLayer
     names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
-  } ### ENDIF_RunSnowModule
+  } ### ENDIF RunSnowModule
   if (!inherits(RunOptions, "CemaNeige")) {
     CemaNeigeLayers <- list()
     CemaNeigeStateEnd <- NULL
@@ -138,20 +138,20 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## MODEL
+  ## GR model
   if ("all" %in% RunOptions$Outputs_Sim) {
     IndOutputsMod <- as.integer(1:length(FortranOutputs$GR))
   } else {
     IndOutputsMod <- which(FortranOutputs$GR %in% RunOptions$Outputs_Sim)
   }
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr4h", PACKAGE = "airGR",
                       ## inputs
                       LInputs = LInputSeries,                          ### length of input and output series
@@ -186,28 +186,28 @@ RunModel_CemaNeigeGR4H <- function(InputsModel, RunOptions, Param) {
     RESULTS$Outputs[, which(FortranOutputs$GR[IndOutputsMod] ==  "Precip")] <- InputsModel$Precip[IndPeriod1]
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if ( ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c("DatesR", FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers),
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers, "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ( ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
diff --git a/R/RunModel_CemaNeigeGR4J.R b/R/RunModel_CemaNeigeGR4J.R
index d525bdbbcdd0235ca7dd517b9783db04c7049088..9156165856ae34a26358c231d9466c8174823bee 100644
--- a/R/RunModel_CemaNeigeGR4J.R
+++ b/R/RunModel_CemaNeigeGR4J.R
@@ -9,7 +9,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR4J", isCN = TRUE)
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -54,7 +54,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }      
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -70,7 +70,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   if (inherits(RunOptions, "CemaNeige")) {
     if ("all" %in% RunOptions$Outputs_Sim) {
       IndOutputsCemaNeige <- as.integer(1:length(FortranOutputs$CN)) 
@@ -82,7 +82,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
     NameCemaNeigeLayers <- "CemaNeigeLayers"
     
     
-    ## Call_DLL_CemaNeige_________________________
+    ## Call CemaNeige Fortran_________________________
     for(iLayer in 1:NLayers) {
       if (!IsHyst) {
         StateStartCemaNeige <- RunOptions$IniStates[(7 + 20 + 40) + c(iLayer, iLayer+NLayers)]
@@ -110,7 +110,7 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
       RESULTS$Outputs[ round(RESULTS$Outputs , 3) == -999.999] <- NA
       RESULTS$StateEnd[round(RESULTS$StateEnd, 3) == -999.999] <- NA
       
-      ## Data_storage
+      ## Data storage
       CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
       names(CemaNeigeLayers[[iLayer]]) <- FortranOutputs$CN[IndOutputsCemaNeige]
       IndPliqAndMelt <- which(names(CemaNeigeLayers[[iLayer]]) == "PliqAndMelt")
@@ -124,9 +124,9 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
         CemaNeigeStateEnd <- c(CemaNeigeStateEnd, RESULTS$StateEnd)
       }
       rm(RESULTS)
-    } ### ENDFOR_iLayer
+    } ### ENDFOR iLayer
     names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
-  } ### ENDIF_RunSnowModule
+  } ### ENDIF RunSnowModule
   if (!inherits(RunOptions, "CemaNeige")) {
     CemaNeigeLayers <- list()
     CemaNeigeStateEnd <- NULL
@@ -136,20 +136,20 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## MODEL______________________________________________________________________________________
+  ## GR model______________________________________________________________________________________
   if ("all" %in% RunOptions$Outputs_Sim) {
     IndOutputsMod <- as.integer(1:length(FortranOutputs$GR))
   } else {
     IndOutputsMod <- which(FortranOutputs$GR %in% RunOptions$Outputs_Sim)
   }
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr4j", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                          ### length of input and output series
@@ -184,28 +184,28 @@ RunModel_CemaNeigeGR4J <- function(InputsModel, RunOptions, Param) {
     RESULTS$Outputs[, which(FortranOutputs$GR[IndOutputsMod] == "Precip")] <- InputsModel$Precip[IndPeriod1]
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if ( ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c("DatesR", FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(CemaNeigeLayers), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers, "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_Sate
+  ## DatesR and OutputsModel and Sate
   if ( ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_CemaNeigeGR5H.R b/R/RunModel_CemaNeigeGR5H.R
index 5bb10771fe004613755dcfa8e349b6b19d55c2d6..6cfcc4180c93eb3fa361830d55b0f382df949b57 100644
--- a/R/RunModel_CemaNeigeGR5H.R
+++ b/R/RunModel_CemaNeigeGR5H.R
@@ -15,7 +15,7 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
   }
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -61,7 +61,7 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }     
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -81,12 +81,12 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
   NStatesMod     <- as.integer(length(RunOptions$IniStates) - NStates * NLayers)
   
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   if (inherits(RunOptions, "CemaNeige")) {
     if ("all" %in% RunOptions$Outputs_Sim) {
       IndOutputsCemaNeige <- as.integer(1:length(FortranOutputs$CN))
@@ -98,7 +98,7 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
     NameCemaNeigeLayers <- "CemaNeigeLayers"
     
     
-    ## Call_DLL_CemaNeige_________________________
+    ## Call CemaNeige Fortran_________________________
     for(iLayer in 1:NLayers) {
       if (!IsHyst) {
         StateStartCemaNeige <- RunOptions$IniStates[(7 + 20*24 + 40*24) + c(iLayer, iLayer+NLayers)]
@@ -126,7 +126,7 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
       RESULTS$Outputs[ round(RESULTS$Outputs , 3) == -999.999] <- NA
       RESULTS$StateEnd[round(RESULTS$StateEnd, 3) == -999.999] <- NA
       
-      ## Data_storage
+      ## Data storage
       CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
       names(CemaNeigeLayers[[iLayer]]) <- FortranOutputs$CN[IndOutputsCemaNeige]
       IndPliqAndMelt <- which(names(CemaNeigeLayers[[iLayer]]) == "PliqAndMelt")
@@ -140,9 +140,9 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
         CemaNeigeStateEnd <- c(CemaNeigeStateEnd, RESULTS$StateEnd)
       }
       rm(RESULTS)
-    } ### ENDFOR_iLayer
+    } ### ENDFOR iLayer
     names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
-  } ### ENDIF_RunSnowModule
+  } ### ENDIF RunSnowModule
   if (!inherits(RunOptions, "CemaNeige")) {
     CemaNeigeLayers <- list()
     CemaNeigeStateEnd <- NULL
@@ -152,14 +152,14 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## MODEL______________________________________________________________________________________ 
+  ## GR model______________________________________________________________________________________ 
   if ("all" %in% RunOptions$Outputs_Sim) {
     IndOutputsMod <- as.integer(1:length(FortranOutputs$GR))
   } else {
     IndOutputsMod <- which(FortranOutputs$GR %in% RunOptions$Outputs_Sim)
   }
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     RunOptions$IniStates[1] <- RunOptions$IniResLevels[1] * Param[1] ### production store level (mm)
     RunOptions$IniStates[2] <- RunOptions$IniResLevels[2] * Param[3] ### routing store level (mm)
@@ -168,7 +168,7 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
     }
   }
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr5h",PACKAGE="airGR",
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -205,29 +205,29 @@ RunModel_CemaNeigeGR5H <- function(InputsModel, RunOptions, Param) {
     RESULTS$Outputs[,which(FortranOutputs$GR[IndOutputsMod]=="Precip")] <- InputsModel$Precip[IndPeriod1]
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if ( ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c("DatesR", FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers),
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers, "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ( ExportDatesR & ExportStateEnd) {
     OutputsModel <- c( list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                        lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
diff --git a/R/RunModel_CemaNeigeGR5J.R b/R/RunModel_CemaNeigeGR5J.R
index 7482fe517e2d52eb4f2e24a3db203425c7f9ed18..3e3344d5e3b589bb310d15e075a528e86a2b91b2 100644
--- a/R/RunModel_CemaNeigeGR5J.R
+++ b/R/RunModel_CemaNeigeGR5J.R
@@ -9,7 +9,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR5J", isCN = TRUE)
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -55,7 +55,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }      
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -71,7 +71,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   if (inherits(RunOptions, "CemaNeige")) {
     if ("all" %in% RunOptions$Outputs_Sim) {
       IndOutputsCemaNeige <- as.integer(1:length(FortranOutputs$CN))
@@ -83,7 +83,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
     NameCemaNeigeLayers <- "CemaNeigeLayers"
     
     
-    ## Call_DLL_CemaNeige_________________________
+    ## Call CemaNeige Fortran_________________________
     for(iLayer in 1:NLayers) {
       if (!IsHyst) {
         StateStartCemaNeige <- RunOptions$IniStates[(7 + 20 + 40) + c(iLayer, iLayer+NLayers)]
@@ -111,7 +111,7 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
       RESULTS$Outputs[ round(RESULTS$Outputs , 3) == -999.999] <- NA
       RESULTS$StateEnd[round(RESULTS$StateEnd, 3) == -999.999] <- NA
       
-      ## Data_storage
+      ## Data storage
       CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
       names(CemaNeigeLayers[[iLayer]]) <- FortranOutputs$CN[IndOutputsCemaNeige]
       IndPliqAndMelt <- which(names(CemaNeigeLayers[[iLayer]]) == "PliqAndMelt")
@@ -125,9 +125,9 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
         CemaNeigeStateEnd <- c(CemaNeigeStateEnd, RESULTS$StateEnd)
       }
       rm(RESULTS)
-    } ### ENDFOR_iLayer
+    } ### ENDFOR iLayer
     names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
-  } ### ENDIF_RunSnowModule
+  } ### ENDIF RunSnowModule
   if (!inherits(RunOptions, "CemaNeige")) {
     CemaNeigeLayers <- list()
     CemaNeigeStateEnd <- NULL
@@ -137,20 +137,20 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## MODEL______________________________________________________________________________________
+  ## GR model______________________________________________________________________________________
   if ("all" %in% RunOptions$Outputs_Sim) {
     IndOutputsMod <- as.integer(1:length(FortranOutputs$GR))
   } else {
     IndOutputsMod <- which(FortranOutputs$GR %in% RunOptions$Outputs_Sim)
   }
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr5j", PACKAGE = "airGR",
                       ## inputs
                       LInputs = LInputSeries,                          ### length of input and output series
@@ -185,28 +185,28 @@ RunModel_CemaNeigeGR5J <- function(InputsModel, RunOptions, Param) {
     RESULTS$Outputs[,which(FortranOutputs$GR[IndOutputsMod] == "Precip")] <- InputsModel$Precip[IndPeriod1]
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c("DatesR", FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers),
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers, "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if (ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
diff --git a/R/RunModel_CemaNeigeGR6J.R b/R/RunModel_CemaNeigeGR6J.R
index 62e1786be80c93aab0c69c7ebd1b0a61270222b3..d82a0241ec42abb8a607d291d60f588885f8defb 100644
--- a/R/RunModel_CemaNeigeGR6J.R
+++ b/R/RunModel_CemaNeigeGR6J.R
@@ -9,7 +9,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR6J", isCN = TRUE)
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -59,7 +59,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
     Param[6L] <- Param_X1X3X6_threshold
   }         
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -75,7 +75,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
   
-  ## SNOW_MODULE________________________________________________________________________________
+  ## CemaNeige________________________________________________________________________________
   if (inherits(RunOptions, "CemaNeige")) {
     if ("all" %in% RunOptions$Outputs_Sim) {
       IndOutputsCemaNeige <- as.integer(1:length(FortranOutputs$CN))
@@ -87,7 +87,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
     NameCemaNeigeLayers <- "CemaNeigeLayers"
     
     
-    ## Call_DLL_CemaNeige_________________________
+    ## Call CemaNeige Fortran_________________________
     for(iLayer in 1:NLayers) {
       if (!IsHyst) {
         StateStartCemaNeige <- RunOptions$IniStates[(7 + 20 + 40) + c(iLayer, iLayer+NLayers)]
@@ -115,7 +115,7 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
       RESULTS$Outputs[ round(RESULTS$Outputs , 3) == -999.999] <- NA
       RESULTS$StateEnd[round(RESULTS$StateEnd, 3) == -999.999] <- NA
       
-      ## Data_storage
+      ## Data storage
       CemaNeigeLayers[[iLayer]] <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
       names(CemaNeigeLayers[[iLayer]]) <- FortranOutputs$CN[IndOutputsCemaNeige]
       IndPliqAndMelt <- which(names(CemaNeigeLayers[[iLayer]]) == "PliqAndMelt")
@@ -129,9 +129,9 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
         CemaNeigeStateEnd <- c(CemaNeigeStateEnd,RESULTS$StateEnd)
       }
       rm(RESULTS)
-    } ### ENDFOR_iLayer
+    } ### ENDFOR iLayer
     names(CemaNeigeLayers) <- sprintf("Layer%02i", seq_len(NLayers))
-  } ### ENDIF_RunSnowModule
+  } ### ENDIF RunSnowModule
   if (!inherits(RunOptions, "CemaNeige")) {
     CemaNeigeLayers <- list()
     CemaNeigeStateEnd <- NULL
@@ -141,21 +141,21 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
   
   
   
-  ## MODEL______________________________________________________________________________________
+  ## GR model______________________________________________________________________________________
   if ("all" %in% RunOptions$Outputs_Sim) {
     IndOutputsMod <- as.integer(1:length(FortranOutputs$GR))
   } else {
     IndOutputsMod <- which(FortranOutputs$GR %in% RunOptions$Outputs_Sim)
   }
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     RunOptions$IniStates[1] <- RunOptions$IniResLevels[1] * ParamMod[1] ### production store level (mm)
     RunOptions$IniStates[2] <- RunOptions$IniResLevels[2] * ParamMod[3] ### routing store level (mm)
     RunOptions$IniStates[3] <- RunOptions$IniResLevels[3]               ### exponential store level (mm)
   }
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr6j", PACKAGE = "airGR",
                       ## inputs
                       LInputs = LInputSeries,                          ### length of input and output series
@@ -190,28 +190,28 @@ RunModel_CemaNeigeGR6J <- function(InputsModel, RunOptions, Param) {
     RESULTS$Outputs[,which(FortranOutputs$GR[IndOutputsMod] == "Precip")] <- InputsModel$Precip[IndPeriod1]
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR== FALSE & !ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers))
     names(OutputsModel) <- c("DatesR", FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers)
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
                       list(CemaNeigeLayers),
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs$GR[IndOutputsMod], NameCemaNeigeLayers, "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if (ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]),
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),
diff --git a/R/RunModel_GR1A.R b/R/RunModel_GR1A.R
index 2470a7fd9b96a1f1e8a462364e087264eda5f568..58bca283713067eb4ebc0547671b0fc61303b790 100644
--- a/R/RunModel_GR1A.R
+++ b/R/RunModel_GR1A.R
@@ -6,7 +6,7 @@ RunModel_GR1A <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR1A")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }
@@ -50,7 +50,7 @@ RunModel_GR1A <- function(InputsModel, RunOptions, Param) {
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr1a", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
diff --git a/R/RunModel_GR2M.R b/R/RunModel_GR2M.R
index 0cf0fcec3b5765752be2c11555343caf7f2243f5..385c5f1a798de5ea4ac460fd3fb0d516c7c6430b 100644
--- a/R/RunModel_GR2M.R
+++ b/R/RunModel_GR2M.R
@@ -6,7 +6,7 @@ RunModel_GR2M <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR2M")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -40,7 +40,7 @@ RunModel_GR2M <- function(InputsModel, RunOptions, Param) {
     Param[2L] <- Param_X1X2_threshold
   }
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -52,18 +52,18 @@ RunModel_GR2M <- function(InputsModel, RunOptions, Param) {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1] ### production store level (mm)
     RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[2] ### routing store level (mm)
   }
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr2M", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -90,25 +90,25 @@ RunModel_GR2M <- function(InputsModel, RunOptions, Param) {
   }
   
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                        lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_GR4H.R b/R/RunModel_GR4H.R
index 7a53b891b94ab473c265505194ed83e8c2a7aac1..47c26ff968cce41eb982f897604c7d4b843d0162 100644
--- a/R/RunModel_GR4H.R
+++ b/R/RunModel_GR4H.R
@@ -6,7 +6,7 @@ RunModel_GR4H <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR4H")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -45,7 +45,7 @@ RunModel_GR4H <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }     
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -57,18 +57,18 @@ RunModel_GR4H <- function(InputsModel, RunOptions, Param) {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr4h", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -95,25 +95,25 @@ RunModel_GR4H <- function(InputsModel, RunOptions, Param) {
                                         verbose = FALSE)
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_GR4J.R b/R/RunModel_GR4J.R
index 92fd979df7618b92dab0c5e28c271ac54f3a90bd..fc647035051dc0f05fb741f0f43ae1f1945ded08 100644
--- a/R/RunModel_GR4J.R
+++ b/R/RunModel_GR4J.R
@@ -6,7 +6,7 @@ RunModel_GR4J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR4J")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -45,7 +45,7 @@ RunModel_GR4J <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -56,18 +56,18 @@ RunModel_GR4J <- function(InputsModel, RunOptions, Param) {
   } else {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
-  ## Input_data_preparation      
+  ## Input data preparation      
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr4j", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -94,25 +94,25 @@ RunModel_GR4J <- function(InputsModel, RunOptions, Param) {
                                         verbose = FALSE)
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_StateEnd_only
+  ## OutputsModel and StateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_StateEnd
+  ## DatesR and OutputsModel and StateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_GR5H.R b/R/RunModel_GR5H.R
index ca920b91744b2adeb54c829811f6cf6ba8483ce1..590a103291d3ea8e99f26681cd3eb1774ef05edb 100644
--- a/R/RunModel_GR5H.R
+++ b/R/RunModel_GR5H.R
@@ -12,7 +12,7 @@ RunModel_GR5H <- function(InputsModel, RunOptions, Param) {
   }
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -51,7 +51,7 @@ RunModel_GR5H <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }     
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -62,12 +62,12 @@ RunModel_GR5H <- function(InputsModel, RunOptions, Param) {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1] ### production store level (mm)
     RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3] ### routing store level (mm)
@@ -76,7 +76,7 @@ RunModel_GR5H <- function(InputsModel, RunOptions, Param) {
     }
   }
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr5h", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -105,25 +105,25 @@ RunModel_GR5H <- function(InputsModel, RunOptions, Param) {
                                         verbose = FALSE)
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_StateEnd_only
+  ## OutputsModel and StateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_StateEnd
+  ## DatesR and OutputsModel and StateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_GR5J.R b/R/RunModel_GR5J.R
index 862bb1dbe760829d08dd135d51d47cead92f4581..33210729deffa802b0e51732e0659e27ee7802be 100644
--- a/R/RunModel_GR5J.R
+++ b/R/RunModel_GR5J.R
@@ -6,7 +6,7 @@ RunModel_GR5J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR5J")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -45,7 +45,7 @@ RunModel_GR5J <- function(InputsModel, RunOptions, Param) {
     Param[4L] <- Param_X4_threshold
   }      
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -57,18 +57,18 @@ RunModel_GR5J <- function(InputsModel, RunOptions, Param) {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     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
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr5j", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -95,25 +95,25 @@ RunModel_GR5J <- function(InputsModel, RunOptions, Param) {
                                         verbose = FALSE)
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
diff --git a/R/RunModel_GR6J.R b/R/RunModel_GR6J.R
index 4f4e317bced78fc08555d8a1b67b4631a39636aa..c971056fbd21f4fd423a0b98abe0065d62c5090b 100644
--- a/R/RunModel_GR6J.R
+++ b/R/RunModel_GR6J.R
@@ -6,7 +6,7 @@ RunModel_GR6J <- function(InputsModel, RunOptions, Param) {
   FortranOutputs <- .FortranOutputs(GR = "GR6J")$GR
   
   
-  ## Arguments_check
+  ## Arguments check
   if (!inherits(InputsModel, "InputsModel")) {
     stop("'InputsModel' must be of class 'InputsModel'")
   }  
@@ -49,7 +49,7 @@ RunModel_GR6J <- function(InputsModel, RunOptions, Param) {
     Param[6L] <- Param_X1X3X6_threshold
   }      
   
-  ## Input_data_preparation
+  ## Input data preparation
   if (identical(RunOptions$IndPeriod_WarmUp, 0L)) {
     RunOptions$IndPeriod_WarmUp <- NULL
   }
@@ -61,19 +61,19 @@ RunModel_GR6J <- function(InputsModel, RunOptions, Param) {
     IndOutputs <- which(FortranOutputs %in% RunOptions$Outputs_Sim)
   }
   
-  ## Output_data_preparation
+  ## Output data preparation
   IndPeriod2     <- (length(RunOptions$IndPeriod_WarmUp)+1):LInputSeries
   ExportDatesR   <- "DatesR"   %in% RunOptions$Outputs_Sim
   ExportStateEnd <- "StateEnd" %in% RunOptions$Outputs_Sim
   
-  ## Use_of_IniResLevels
+  ## Use of IniResLevels
   if (!is.null(RunOptions$IniResLevels)) {
     RunOptions$IniStates[1] <- RunOptions$IniResLevels[1]*Param[1] ### production store level (mm)
     RunOptions$IniStates[2] <- RunOptions$IniResLevels[2]*Param[3] ### routing store level (mm)
     RunOptions$IniStates[3] <- RunOptions$IniResLevels[3]          ### exponential store level (mm)
   }
   
-  ## Call_fortan
+  ## Call GR model Fortan
   RESULTS <- .Fortran("frun_gr6j", PACKAGE = "airGR", 
                       ## inputs
                       LInputs = LInputSeries,                             ### length of input and output series
@@ -100,25 +100,25 @@ RunModel_GR6J <- function(InputsModel, RunOptions, Param) {
                                         verbose = FALSE)
   }
   
-  ## Output_data_preparation
-  ## OutputsModel_only
+  ## Output data preparation
+  ## OutputsModel only
   if (!ExportDatesR & !ExportStateEnd) {
     OutputsModel <- lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i])
     names(OutputsModel) <- FortranOutputs[IndOutputs]
   }
-  ## DatesR_and_OutputsModel_only
+  ## DatesR and OutputsModel only
   if (ExportDatesR & !ExportStateEnd) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]))
     names(OutputsModel) <- c("DatesR", FortranOutputs[IndOutputs])
   }
-  ## OutputsModel_and_SateEnd_only
+  ## OutputsModel and SateEnd only
   if (!ExportDatesR & ExportStateEnd) {
     OutputsModel <- c(lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]), 
                       list(RESULTS$StateEnd))
     names(OutputsModel) <- c(FortranOutputs[IndOutputs], "StateEnd")
   }
-  ## DatesR_and_OutputsModel_and_SateEnd
+  ## DatesR and OutputsModel and SateEnd
   if ((ExportDatesR & ExportStateEnd) | "all" %in% RunOptions$Outputs_Sim) {
     OutputsModel <- c(list(InputsModel$DatesR[RunOptions$IndPeriod_Run]), 
                       lapply(seq_len(RESULTS$NOutputs), function(i) RESULTS$Outputs[IndPeriod2, i]),