Source

Target

Showing with 1987 additions and 1822 deletions
+1987 -1822
This diff is collapsed.
ErrorCrit <- function(InputsCrit,OutputsModel,FUN_CRIT, warnings = TRUE, verbose = TRUE){
return( FUN_CRIT(InputsCrit,OutputsModel, warnings = warnings, verbose = verbose) )
ErrorCrit <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
## ---------- Arguments check
if (!inherits(InputsCrit, "InputsCrit")) {
stop("InputsCrit must be of class 'InputsCrit'")
}
if (!inherits(OutputsModel, "OutputsModel")) {
stop("OutputsModel must be of class 'OutputsModel'")
}
## ---------- Criterion computation
## ----- Single criterion
if (inherits(InputsCrit, "Single")) {
FUN_CRIT <- match.fun(InputsCrit$FUN_CRIT)
OutputsCrit <- FUN_CRIT(InputsCrit = InputsCrit,
OutputsModel = OutputsModel,
warnings = warnings,
verbose = verbose)
}
## ----- Multiple criteria or Composite criterion
if (inherits(InputsCrit, "Multi") | inherits(InputsCrit, "Compo")) {
listOutputsCrit <- lapply(InputsCrit, FUN = function(iInputsCrit) {
FUN_CRIT <- match.fun(iInputsCrit$FUN_CRIT)
FUN_CRIT(InputsCrit = iInputsCrit,
OutputsModel = OutputsModel,
warnings = warnings,
verbose = verbose)
})
listValCrit <- sapply(listOutputsCrit, function(x) x[["CritValue"]])
listNameCrit <- sapply(listOutputsCrit, function(x) x[["CritName"]])
listweights <- unlist(lapply(InputsCrit, function(x) x[["Weights"]]))
listweights <- listweights / sum(listweights)
if (inherits(InputsCrit, "Compo")) {
CritValue <- sum(listValCrit * listweights)
OutputsCritCompo <- list(MultiCritValues = listValCrit,
MultiCritNames = listNameCrit,
MultiCritWeights = listweights)
OutputsCrit <- list(CritValue = CritValue,
CritName = "Composite",
CritBestValue = +1,
Multiplier = -1,
Ind_notcomputed = NULL,
CritCompo = OutputsCritCompo,
MultiCrit = listOutputsCrit)
class(OutputsCrit) <- c("Compo", "ErrorCrit")
if (verbose) {
message("------------------------------------\n")
message("Crit. Composite = ", sprintf("%.4f", CritValue))
msgForm <- paste(sprintf("%.2f", listweights), listNameCrit, sep = " * ", collapse = ", ")
msgForm <- unlist(strsplit(msgForm, split = ","))
msgFormSep <- rep(c(",", ",", ",\n\t\t "), times = ceiling(length(msgForm)/3))[1:length(msgForm)]
msgForm <- paste(msgForm, msgFormSep, sep = "", collapse = "")
msgForm <- gsub("\\,\\\n\\\t\\\t $|\\,$", "", msgForm)
message("\tFormula: sum(", msgForm, ")\n")
}
} else {
OutputsCrit <- listOutputsCrit
class(OutputsCrit) <- c("Multi", "ErrorCrit")
}
}
return(OutputsCrit)
}
ErrorCrit_KGE <- function(InputsCrit,OutputsModel, warnings = TRUE, verbose = TRUE){
##Arguments_check________________________________
if(inherits(InputsCrit,"InputsCrit")==FALSE){ stop("InputsCrit must be of class 'InputsCrit' \n"); return(NULL); }
if(inherits(OutputsModel,"OutputsModel")==FALSE){ stop("OutputsModel must be of class 'OutputsModel' \n"); return(NULL); }
##Initialisation_________________________________
CritName <- NA;
if(InputsCrit$transfo=="" ){ CritName <- "KGE[Q]" ; }
if(InputsCrit$transfo=="sqrt"){ CritName <- "KGE[sqrt(Q)]"; }
if(InputsCrit$transfo=="log" ){ CritName <- "KGE[log(Q)]" ; }
if(InputsCrit$transfo=="inv" ){ CritName <- "KGE[1/Q]" ; }
if(InputsCrit$transfo=="sort"){ CritName <- "KGE[sort(Q)]"; }
CritValue <- NA;
CritBestValue <- +1;
Multiplier <- -1; ### must be equal to -1 or +1 only
##Data_preparation_______________________________
VarObs <- InputsCrit$Qobs ; VarObs[!InputsCrit$BoolCrit] <- NA;
VarSim <- OutputsModel$Qsim; VarSim[!InputsCrit$BoolCrit] <- NA;
##Data_transformation
if("Ind_zeroes" %in% names(InputsCrit) & "epsilon" %in% names(InputsCrit)){ if(length(InputsCrit$Ind_zeroes)>0){
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)
ErrorCrit_KGE <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
## Arguments check
if (!inherits(OutputsModel, "OutputsModel")) {
stop("'OutputsModel' must be of class 'OutputsModel'")
}
##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(sum(!TS_ignore)==1){ OutputsCrit <- list(NA); names(OutputsCrit) <- c("CritValue"); return(OutputsCrit); } ### to avoid a problem in standard deviation computation
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]);
iCrit <- 0
SubCritPrint <- NULL
SubCritNames <- NULL
SubCritValues <- NULL
##SubErrorCrit_____KGE_rPearson__________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," cor(sim, obs, \"pearson\") =", sep = "")
SubCritValues[iCrit] <- NA;
SubCritNames[iCrit] <- "r"
Numer <- sum( (VarObs[!TS_ignore]-meanVarObs)*(VarSim[!TS_ignore]-meanVarSim) );
Deno1 <- sqrt( sum((VarObs[!TS_ignore]-meanVarObs)^2) );
Deno2 <- sqrt( sum((VarSim[!TS_ignore]-meanVarSim)^2) );
if(Numer==0){ if(Deno1==0 & Deno2==0){ Crit <- 1; } else { Crit <- 0; }
} else { Crit <- Numer/(Deno1*Deno2); }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##SubErrorCrit_____KGE_alpha_____________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," sd(sim)/sd(obs) =", sep = "")
SubCritValues[iCrit] <- NA;
SubCritNames[iCrit] <- "alpha"
Numer <- sd(VarSim[!TS_ignore]);
Denom <- sd(VarObs[!TS_ignore]);
if(Numer==0 & Denom==0){ Crit <- 1; } else { Crit <- Numer/Denom ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##SubErrorCrit_____KGE_beta______________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," mean(sim)/mean(obs) =", sep = "")
SubCritValues[iCrit] <- NA;
SubCritNames[iCrit] <- "beta"
if(meanVarSim==0 & meanVarObs==0){ Crit <- 1; } else { Crit <- meanVarSim/meanVarObs ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##ErrorCrit______________________________________
if(sum(is.na(SubCritValues))==0){
CritValue <- ( 1 - sqrt( (SubCritValues[1]-1)^2 + (SubCritValues[2]-1)^2 + (SubCritValues[3]-1)^2 ) );
EC <- .ErrorCrit(InputsCrit = InputsCrit, Crit = "KGE", OutputsModel = OutputsModel, warnings = warnings)
CritValue <- NA
SubCritValues <- rep(NA, 3)
SubCritNames <- c("r", "alpha", "beta")
SubCritPrint <- rep(NA, 3)
if (EC$CritCompute) {
## Other variables preparation
meanVarObs <- mean(EC$VarObs[!EC$TS_ignore])
meanVarSim <- mean(EC$VarSim[!EC$TS_ignore])
## SubErrorCrit KGE rPearson
SubCritPrint[1L] <- paste0(EC$CritName, " cor(sim, obs, \"pearson\") =")
Numer <- sum((EC$VarObs[!EC$TS_ignore] - meanVarObs) * (EC$VarSim[!EC$TS_ignore] - meanVarSim))
Deno1 <- sqrt(sum((EC$VarObs[!EC$TS_ignore] - meanVarObs) ^ 2))
Deno2 <- sqrt(sum((EC$VarSim[!EC$TS_ignore] - meanVarSim) ^ 2))
if (Numer == 0) {
if (Deno1 == 0 & Deno2 == 0) {
Crit <- 1
} else {
Crit <- 0
}
} else {
Crit <- Numer / (Deno1 * Deno2)
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[1L] <- Crit
}
## SubErrorCrit KGE alpha
SubCritPrint[2L] <- paste0(EC$CritName, " sd(sim)/sd(obs) =")
Numer <- sd(EC$VarSim[!EC$TS_ignore])
Denom <- sd(EC$VarObs[!EC$TS_ignore])
if (Numer == 0 & Denom == 0) {
Crit <- 1
} else {
Crit <- Numer / Denom
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[2L] <- Crit
}
## SubErrorCrit KGE beta
SubCritPrint[3L] <- paste0(EC$CritName, " mean(sim)/mean(obs) =")
if (meanVarSim == 0 & meanVarObs == 0) {
Crit <- 1
} else {
Crit <- meanVarSim / meanVarObs
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[3L] <- Crit
}
## ErrorCrit
if (sum(is.na(SubCritValues)) == 0) {
CritValue <- (1 - sqrt((SubCritValues[1L] - 1)^2 + (SubCritValues[2L] - 1)^2 + (SubCritValues[3L] - 1)^2))
}
## Verbose
if (verbose) {
message(sprintf("Crit. %s = %.4f", EC$CritName, CritValue))
message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " "))
}
}
##Verbose______________________________________
if(verbose) {
message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue))
message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " "))
}
## Output
OutputsCrit <- list(CritValue = CritValue,
CritName = EC$CritName,
SubCritValues = SubCritValues,
SubCritNames = SubCritNames,
CritBestValue = EC$CritBestValue,
Multiplier = EC$Multiplier,
Ind_notcomputed = EC$Ind_TS_ignore)
##Output_________________________________________
OutputsCrit <- list(CritValue = CritValue, CritName = CritName,
SubCritValues = SubCritValues, SubCritNames = SubCritNames, CritBestValue = CritBestValue,
Multiplier = Multiplier, Ind_notcomputed = Ind_TS_ignore)
class(OutputsCrit) <- c("KGE", "ErrorCrit")
return(OutputsCrit)
}
class(ErrorCrit_KGE) <- c("FUN_CRIT", class(ErrorCrit_KGE))
ErrorCrit_KGE2 <- function(InputsCrit,OutputsModel, warnings = TRUE, verbose = TRUE){
##Arguments_check________________________________
if(inherits(InputsCrit,"InputsCrit")==FALSE){ stop("InputsCrit must be of class 'InputsCrit' \n"); return(NULL); }
if(inherits(OutputsModel,"OutputsModel")==FALSE){ stop("OutputsModel must be of class 'OutputsModel' \n"); return(NULL); }
##Initialisation_________________________________
CritName <- NA;
if(InputsCrit$transfo=="" ){ CritName <- "KGE'[Q]" ; }
if(InputsCrit$transfo=="sqrt"){ CritName <- "KGE'[sqrt(Q)]"; }
if(InputsCrit$transfo=="log" ){ CritName <- "KGE'[log(Q)]" ; }
if(InputsCrit$transfo=="inv" ){ CritName <- "KGE'[1/Q]" ; }
if(InputsCrit$transfo=="sort"){ CritName <- "KGE'[sort(Q)]"; }
CritValue <- NA;
CritBestValue <- +1;
Multiplier <- -1; ### must be equal to -1 or +1 only
##Data_preparation_______________________________
VarObs <- InputsCrit$Qobs ; VarObs[!InputsCrit$BoolCrit] <- NA;
VarSim <- OutputsModel$Qsim; VarSim[!InputsCrit$BoolCrit] <- NA;
##Data_transformation
if("Ind_zeroes" %in% names(InputsCrit) & "epsilon" %in% names(InputsCrit)){ if(length(InputsCrit$Ind_zeroes)>0){
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(sum(!TS_ignore)==1){ OutputsCrit <- list(NA); names(OutputsCrit) <- c("CritValue"); return(OutputsCrit); } ### to avoid a problem in standard deviation computation
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]);
iCrit <- 0
SubCritPrint <- NULL
SubCritNames <- NULL
SubCritValues <- NULL
##SubErrorCrit_____KGE_rPearson__________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," cor(sim, obs, \"pearson\") =", sep = "")
SubCritNames[iCrit] <- "r"
SubCritValues[iCrit] <- NA;
Numer <- sum( (VarObs[!TS_ignore]-meanVarObs)*(VarSim[!TS_ignore]-meanVarSim) );
Deno1 <- sqrt( sum((VarObs[!TS_ignore]-meanVarObs)^2) );
Deno2 <- sqrt( sum((VarSim[!TS_ignore]-meanVarSim)^2) );
if(Numer==0){ if(Deno1==0 & Deno2==0){ Crit <- 1; } else { Crit <- 0; }
} else { Crit <- Numer/(Deno1*Deno2); }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##SubErrorCrit_____KGE_gamma______________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," cv(sim)/cv(obs) =", sep = "")
SubCritNames[iCrit] <- "gamma"
SubCritValues[iCrit] <- NA;
if(meanVarSim==0){ if(sd(VarSim[!TS_ignore])==0){ CVsim <- 1; } else { CVsim <- 99999; } } else { CVsim <- sd(VarSim[!TS_ignore])/meanVarSim; }
if(meanVarObs==0){ if(sd(VarObs[!TS_ignore])==0){ CVobs <- 1; } else { CVobs <- 99999; } } else { CVobs <- sd(VarObs[!TS_ignore])/meanVarObs; }
if(CVsim==0 & CVobs==0){ Crit <- 1; } else { Crit <- CVsim/CVobs ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##SubErrorCrit_____KGE_beta______________________
iCrit <- iCrit+1;
SubCritPrint[iCrit] <- paste(CritName," mean(sim)/mean(obs) =", sep = "")
SubCritNames[iCrit] <- "beta"
SubCritValues[iCrit] <- NA;
if(meanVarSim==0 & meanVarObs==0){ Crit <- 1; } else { Crit <- meanVarSim/meanVarObs ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
##ErrorCrit______________________________________
if(sum(is.na(SubCritValues))==0){
CritValue <- ( 1 - sqrt( (SubCritValues[1]-1)^2 + (SubCritValues[2]-1)^2 + (SubCritValues[3]-1)^2 ) );
ErrorCrit_KGE2 <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
## Arguments check
if (!inherits(OutputsModel, "OutputsModel")) {
stop("'OutputsModel' must be of class 'OutputsModel'")
}
##Verbose______________________________________
if(verbose) {
message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue))
message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " "))
EC <- .ErrorCrit(InputsCrit = InputsCrit, Crit = "KGE2", OutputsModel = OutputsModel, warnings = warnings)
CritValue <- NA
SubCritValues <- rep(NA, 3)
SubCritNames <- c("r", "gamma", "beta")
SubCritPrint <- rep(NA, 3)
if (EC$CritCompute) {
## Other variables preparation
meanVarObs <- mean(EC$VarObs[!EC$TS_ignore])
meanVarSim <- mean(EC$VarSim[!EC$TS_ignore])
## SubErrorCrit KGE rPearson
SubCritPrint[1L] <- paste0(EC$CritName, " cor(sim, obs, \"pearson\") =")
Numer <- sum((EC$VarObs[!EC$TS_ignore] - meanVarObs) * (EC$VarSim[!EC$TS_ignore] - meanVarSim))
Deno1 <- sqrt(sum((EC$VarObs[!EC$TS_ignore] - meanVarObs)^2))
Deno2 <- sqrt(sum((EC$VarSim[!EC$TS_ignore] - meanVarSim)^2))
if (Numer == 0) {
if (Deno1 == 0 & Deno2 == 0) {
Crit <- 1
} else {
Crit <- 0
}
} else {
Crit <- Numer / (Deno1 * Deno2)
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[1L] <- Crit
}
## SubErrorCrit KGE gamma
SubCritPrint[2L] <- paste0(EC$CritName, " cv(sim)/cv(obs) =")
if (meanVarSim == 0) {
if (sd(EC$VarSim[!EC$TS_ignore]) == 0) {
CVsim <- 1
} else {
CVsim <- 99999
}
} else {
CVsim <- sd(EC$VarSim[!EC$TS_ignore]) / meanVarSim
}
if (meanVarObs == 0) {
if (sd(EC$VarObs[!EC$TS_ignore]) == 0) {
CVobs <- 1
} else {
CVobs <- 99999
}
} else {
CVobs <- sd(EC$VarObs[!EC$TS_ignore]) / meanVarObs
}
if (CVsim == 0 &
CVobs == 0) {
Crit <- 1
} else {
Crit <- CVsim / CVobs
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[2L] <- Crit
}
## SubErrorCrit KGE beta
SubCritPrint[3L] <- paste0(EC$CritName, " mean(sim)/mean(obs) =")
if (meanVarSim == 0 & meanVarObs == 0) {
Crit <- 1
} else {
Crit <- meanVarSim / meanVarObs
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[3L] <- Crit
}
## ErrorCrit
if (sum(is.na(SubCritValues)) == 0) {
CritValue <- (1 - sqrt((SubCritValues[1L] - 1)^2 + (SubCritValues[2L] - 1)^2 + (SubCritValues[3L] - 1)^2))
}
## Verbose
if (verbose) {
message(sprintf("Crit. %s = %.4f", EC$CritName, CritValue))
message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " "))
}
}
##Output_________________________________________
OutputsCrit <- list(CritValue = CritValue, CritName = CritName,
SubCritValues = SubCritValues, SubCritNames = SubCritNames, CritBestValue = CritBestValue,
Multiplier = Multiplier, Ind_notcomputed = Ind_TS_ignore)
## Output
OutputsCrit <- list(CritValue = CritValue,
CritName = EC$CritName,
SubCritValues = SubCritValues,
SubCritNames = SubCritNames,
CritBestValue = EC$CritBestValue,
Multiplier = EC$Multiplier,
Ind_notcomputed = EC$Ind_TS_ignore)
class(OutputsCrit) <- c("KGE2", "ErrorCrit")
return(OutputsCrit)
}
class(ErrorCrit_KGE2) <- c("FUN_CRIT", class(ErrorCrit_KGE2))
ErrorCrit_NSE <- function(InputsCrit,OutputsModel, warnings = TRUE, verbose = TRUE){
ErrorCrit_NSE <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
## Arguments check
if (!inherits(OutputsModel, "OutputsModel")) {
stop("'OutputsModel' must be of class 'OutputsModel'")
}
##Arguments_check________________________________
if(inherits(InputsCrit,"InputsCrit")==FALSE){ stop("InputsCrit must be of class 'InputsCrit' \n"); return(NULL); }
if(inherits(OutputsModel,"OutputsModel")==FALSE){ stop("OutputsModel must be of class 'OutputsModel' \n"); return(NULL); }
EC <- .ErrorCrit(InputsCrit = InputsCrit, Crit = "NSE", OutputsModel = OutputsModel, warnings = warnings)
CritValue <- NA
##Initialisation_________________________________
CritName <- NA;
if(InputsCrit$transfo=="" ){ CritName <- "NSE[Q]" ; }
if(InputsCrit$transfo=="sqrt"){ CritName <- "NSE[sqrt(Q)]"; }
if(InputsCrit$transfo=="log" ){ CritName <- "NSE[log(Q)]" ; }
if(InputsCrit$transfo=="inv" ){ CritName <- "NSE[1/Q]" ; }
if(InputsCrit$transfo=="sort"){ CritName <- "NSE[sort(Q)]"; }
CritValue <- NA;
CritBestValue <- +1;
Multiplier <- -1; ### must be equal to -1 or +1 only
if (EC$CritCompute) {
## ErrorCrit
Emod <- sum((EC$VarSim[!EC$TS_ignore] - EC$VarObs[!EC$TS_ignore])^2)
Eref <- sum((EC$VarObs[!EC$TS_ignore] - mean(EC$VarObs[!EC$TS_ignore]))^2)
if (Emod == 0 & Eref == 0) {
Crit <- 0
} else {
Crit <- (1 - Emod / Eref)
}
if (is.numeric(Crit) & is.finite(Crit)) {
CritValue <- Crit
}
##Data_preparation_______________________________
VarObs <- InputsCrit$Qobs ; VarObs[!InputsCrit$BoolCrit] <- NA;
VarSim <- OutputsModel$Qsim; VarSim[!InputsCrit$BoolCrit] <- NA;
##Data_transformation
if("Ind_zeroes" %in% names(InputsCrit) & "epsilon" %in% names(InputsCrit)){ if(length(InputsCrit$Ind_zeroes)>0){
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)
## Verbose
if (verbose) {
message(sprintf("Crit. %s = %.4f", EC$CritName, CritValue))
}
}
##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))
}
##Output_________________________________________
OutputsCrit <- list(CritValue = CritValue, CritName = CritName,
CritBestValue = CritBestValue,
Multiplier = Multiplier, Ind_notcomputed = Ind_TS_ignore)
## Output
OutputsCrit <- list(CritValue = CritValue,
CritName = EC$CritName,
CritBestValue = EC$CritBestValue,
Multiplier = EC$Multiplier,
Ind_notcomputed = EC$Ind_TS_ignore)
class(OutputsCrit) <- c("NSE", "ErrorCrit")
return(OutputsCrit)
}
class(ErrorCrit_NSE) <- c("FUN_CRIT", class(ErrorCrit_NSE))
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.