Newer
Older
Delaigue Olivier
committed
ErrorCrit_KGE2 <- function(InputsCrit,OutputsModel, 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; }
Delaigue Olivier
committed
if(sum(!TS_ignore)<WarningTS & verbose){ warning(paste("\t criterion computed on less than ",WarningTS," time-steps \n",sep="")); }
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
##Other_variables_preparation
meanVarObs <- mean(VarObs[!TS_ignore]);
meanVarSim <- mean(VarSim[!TS_ignore]);
iCrit <- 0;
SubCritNames <- NULL;
SubCritValues <- NULL;
##SubErrorCrit_____KGE_rPearson__________________
iCrit <- iCrit+1;
SubCritNames[iCrit] <- paste(CritName," rPEARSON(sim vs. obs)",sep="");
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_gama______________________
iCrit <- iCrit+1;
SubCritNames[iCrit] <- paste(CritName," CVsim/CVobs",sep="");
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;
SubCritNames[iCrit] <- paste(CritName," MEANsim/MEANobs",sep="");
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 ) );
}
##Output_________________________________________
OutputsCrit <- list(CritValue,CritName,SubCritValues,SubCritNames,CritBestValue,Multiplier,Ind_TS_ignore);
names(OutputsCrit) <- c("CritValue","CritName","SubCritValues","SubCritNames","CritBestValue","Multiplier","Ind_notcomputed");
return(OutputsCrit);
}