diff --git a/DESCRIPTION b/DESCRIPTION
index ca33631094054ddd06b5d8a9d081cdbad81c0c9d..c8e6049b6981951a41ab7738b9db7bd6dafd2e7c 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.2.14.3
+Version: 1.2.14.4
Date: 2019-04-16
Authors@R: c(
person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")),
diff --git a/NEWS.rmd b/NEWS.rmd
index 059a294da78e9872c88658e53452b12d8c95b90d..d55dcb006d647be36c676493bd68e39c369a0dc4 100644
--- a/NEWS.rmd
+++ b/NEWS.rmd
@@ -14,7 +14,7 @@ output:
-### 1.2.14.3 Release Notes (2019-04-16)
+### 1.2.14.4 Release Notes (2019-04-16)
#### New features
diff --git a/R/ErrorCrit_NSE.R b/R/ErrorCrit_NSE.R
index 599c6f138a2aa312fafd55fa83ed3d63e550bd27..72c2205168f5341a720bf2648df4dcd70c5e1187 100644
--- a/R/ErrorCrit_NSE.R
+++ b/R/ErrorCrit_NSE.R
@@ -1,145 +1,40 @@
ErrorCrit_NSE <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
-
- ##Arguments_check________________________________
- if (!inherits(InputsCrit, "InputsCrit")) {
- stop("'InputsCrit' must be of class 'InputsCrit'")
- }
- if (inherits(InputsCrit, "Multi") | inherits(InputsCrit, "Compo")) {
- stop("'InputsCrit' must be of class 'Single'. Use the 'ErrorCrit' function on objects of class 'Multi' or 'Compo' with NSE")
- }
+ ## Arguments check
if (!inherits(OutputsModel, "OutputsModel")) {
stop("'OutputsModel' must be of class 'OutputsModel'")
}
+ OutputsCritCheck <- .ErrorCrit(InputsCrit = InputsCrit, crit = "NSE")
- ##Initialisation_________________________________
- CritName <- NA
- CritVar <- InputsCrit$VarObs
- if (InputsCrit$transfo == "") {
- CritName <- "NSE[CritVar]"
- }
- if (InputsCrit$transfo == "sqrt") {
- CritName <- "NSE[sqrt(CritVar)]"
- }
- if (InputsCrit$transfo == "log") {
- CritName <- "NSE[log(CritVar)]"
- }
- if (InputsCrit$transfo == "inv") {
- CritName <- "NSE[1/CritVar]"
- }
- if (InputsCrit$transfo == "sort") {
- CritName <- "NSE[sort(CritVar)]"
- }
- CritName <- gsub(pattern = "CritVar", replacement = CritVar, x = CritName)
- CritValue <- NA
- CritBestValue <- +1
- Multiplier <- -1
- ### must be equal to -1 or +1 only
-
-
- ##Data_preparation_______________________________
- VarObs <- InputsCrit$Obs
- VarObs[!InputsCrit$BoolCrit] <- NA
- if (InputsCrit$VarObs == "Q") {
- VarSim <- OutputsModel$Qsim
- }
- if (InputsCrit$VarObs == "SCA") {
- VarSim <- rowMeans(sapply(OutputsModel$CemaNeigeLayers[InputsCrit$idLayer], FUN = "[[", "Gratio"))
- }
- if (InputsCrit$VarObs == "SWE") {
- VarSim <- rowMeans(sapply(OutputsModel$CemaNeigeLayers[InputsCrit$idLayer], FUN = "[[", "SnowPack"))
- }
- VarSim[!InputsCrit$BoolCrit] <- NA
-
- ##Data_transformation
- if (InputsCrit$transfo %in% c("log", "inv") & is.null(InputsCrit$epsilon) & warnings) {
- if (any(VarObs %in% 0)) {
- warning("zeroes detected in 'Qobs': the corresponding time-steps will be excluded from the criteria computation if the epsilon argument of 'CreateInputsCrit' = NULL")
+ if (!OutputsCritCheck$CritCompute) {
+ CritValue <- NA
+ } else {
+ ## Other variables preparation
+ meanVarObs <- mean(VarObs[!TS_ignore])
+ meanVarSim <- mean(VarSim[!TS_ignore])
+
+ ## ErrorCrit
+ Emod <- sum((VarSim[!TS_ignore] - VarObs[!TS_ignore])^2)
+ Eref <- sum((VarObs[!TS_ignore] - mean(VarObs[!TS_ignore]))^2)
+
+ if (Emod == 0 & Eref == 0) {
+ Crit <- 0
+ } else {
+ Crit <- (1 - Emod / Eref)
+ }
+ if (is.numeric(Crit) & is.finite(Crit)) {
+ CritValue <- Crit
}
- if (any(VarSim %in% 0)) {
- warning("zeroes detected in 'Qsim': the corresponding time-steps will be excluded from the criteria computation if the epsilon argument of 'CreateInputsCrit' = NULL")
- }
- }
- if ("epsilon" %in% names(InputsCrit) & !is.null(InputsCrit$epsilon)) {
- VarObs <- VarObs + InputsCrit$epsilon
- VarSim <- VarSim + InputsCrit$epsilon
- }
- if (InputsCrit$transfo == "sqrt") {
- VarObs <- sqrt(VarObs)
- VarSim <- sqrt(VarSim)
- }
- if (InputsCrit$transfo == "log") {
- VarObs <- log(VarObs)
- VarSim <- log(VarSim)
- VarSim[VarSim < -1e100] <- NA
- }
- if (InputsCrit$transfo == "inv") {
- VarObs <- 1 / VarObs
- VarSim <- 1 / VarSim
- VarSim[abs(VarSim) > 1e+100] <- NA
- }
- if (InputsCrit$transfo == "sort") {
- VarSim[is.na(VarObs)] <- NA
- VarSim <- sort(VarSim, na.last = TRUE)
- VarObs <- sort(VarObs, na.last = TRUE)
- InputsCrit$BoolCrit <- sort(InputsCrit$BoolCrit, decreasing = TRUE)
- }
-
- ##TS_ignore
- TS_ignore <- !is.finite(VarObs) | !is.finite(VarSim) | !InputsCrit$BoolCrit
- Ind_TS_ignore <- which(TS_ignore)
- if (length(Ind_TS_ignore) == 0) {
- Ind_TS_ignore <- NULL
- }
- if (sum(!TS_ignore) == 0) {
- OutputsCrit <- list(NA)
- names(OutputsCrit) <- c("CritValue")
- return(OutputsCrit)
- }
- if (inherits(OutputsModel, "hourly")) {
- WarningTS <- 365
- }
- if (inherits(OutputsModel, "daily")) {
- WarningTS <- 365
- }
- if (inherits(OutputsModel, "monthly")) {
- WarningTS <- 12
- }
- if (inherits(OutputsModel, "yearly")) {
- WarningTS <- 3
- }
- if (sum(!TS_ignore) < WarningTS & warnings) {
- warning("\t criterion computed on less than ", WarningTS, " time-steps")
- }
-
- ##Other_variables_preparation
- meanVarObs <- mean(VarObs[!TS_ignore])
- meanVarSim <- mean(VarSim[!TS_ignore])
-
-
- ##ErrorCrit______________________________________
- Emod <- sum((VarSim[!TS_ignore] - VarObs[!TS_ignore])^2)
- Eref <- sum((VarObs[!TS_ignore] - mean(VarObs[!TS_ignore]))^2)
-
- if (Emod == 0 & Eref == 0) {
- Crit <- 0
- } else {
- Crit <- (1 - Emod / Eref)
- }
- if (is.numeric(Crit) & is.finite(Crit)) {
- CritValue <- Crit
- }
-
-
- ##Verbose______________________________________
- if (verbose) {
- message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue), "\n")
+ ## Verbose
+ if (verbose) {
+ message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue), "\n")
+ }
}
- ##Output_________________________________________
+ ## Output
OutputsCrit <- list(CritValue = CritValue,
CritName = CritName,
CritBestValue = CritBestValue,
@@ -148,4 +43,5 @@ ErrorCrit_NSE <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = T
class(OutputsCrit) <- c("NSE", "ErrorCrit")
return(OutputsCrit)
+
}