Commit dd7c0dbb authored by Thirel Guillaume's avatar Thirel Guillaume Committed by Dorchies David

feat: Modification of Michel's calibration for SD models

Refs #34
parent 52b3cb67
Pipeline #12674 passed with stages
in 11 minutes and 47 seconds
......@@ -56,6 +56,7 @@ export(TransfoParam_GR5H)
export(TransfoParam_GR4J)
export(TransfoParam_GR5J)
export(TransfoParam_GR6J)
export(TransfoParam_LAG)
export(plot.OutputsModel)
exportPattern(".FortranOutputs")
exportPattern(".ErrorCrit")
......
......@@ -57,27 +57,31 @@ Calibration_Michel <- function(InputsModel,
##_check_FUN_TRANSFO
if (is.null(FUN_TRANSFO)) {
if (identical(FUN_MOD, RunModel_GR4H )) {
FUN_TRANSFO <- TransfoParam_GR4H
if (identical(FUN_MOD, RunModel_GR4H) | identical(FUN_MOD, RunModel_CemaNeigeGR4H)) {
FUN1 <- TransfoParam_GR4H
}
if (identical(FUN_MOD, RunModel_GR5H )) {
FUN_TRANSFO <- TransfoParam_GR5H
if (identical(FUN_MOD, RunModel_GR5H) | identical(FUN_MOD, RunModel_CemaNeigeGR5H)) {
FUN1 <- TransfoParam_GR5H
}
if (identical(FUN_MOD, RunModel_GR4J )) {
FUN_TRANSFO <- TransfoParam_GR4J
if (identical(FUN_MOD, RunModel_GR4J) | identical(FUN_MOD, RunModel_CemaNeigeGR4J)) {
FUN1 <- TransfoParam_GR4J
}
if (identical(FUN_MOD, RunModel_GR5J )) {
FUN_TRANSFO <- TransfoParam_GR5J
if (identical(FUN_MOD, RunModel_GR5J) | identical(FUN_MOD, RunModel_CemaNeigeGR5J)) {
FUN1 <- TransfoParam_GR5J
}
if (identical(FUN_MOD, RunModel_GR6J )) {
FUN_TRANSFO <- TransfoParam_GR6J
if (identical(FUN_MOD, RunModel_GR6J) | identical(FUN_MOD, RunModel_CemaNeigeGR6J)) {
FUN1 <- TransfoParam_GR6J
}
if (identical(FUN_MOD, RunModel_GR2M )) {
FUN_TRANSFO <- TransfoParam_GR2M
if (identical(FUN_MOD, RunModel_GR2M )) {
FUN1 <- TransfoParam_GR2M
}
if (identical(FUN_MOD, RunModel_GR1A )) {
FUN_TRANSFO <- TransfoParam_GR1A
if (identical(FUN_MOD, RunModel_GR1A )) {
FUN1 <- TransfoParam_GR1A
}
##_set_FUN_LAG
if (inherits(CalibOptions, "SD")) {
FUN_LAG <- TransfoParam_LAG
}
if (identical(FUN_MOD, RunModel_CemaNeige )) {
if (inherits(CalibOptions, "hysteresis")) {
FUN_TRANSFO <- TransfoParam_CemaNeigeHyst
......@@ -85,29 +89,35 @@ Calibration_Michel <- function(InputsModel,
FUN_TRANSFO <- TransfoParam_CemaNeige
}
}
if (identical(FUN_MOD, RunModel_GR4H) | identical(FUN_MOD, RunModel_GR5H) |
identical(FUN_MOD, RunModel_GR4J) | identical(FUN_MOD, RunModel_GR5J) | identical(FUN_MOD, RunModel_GR6J)) {
if (!inherits(CalibOptions, "SD")) {
FUN_TRANSFO <- FUN1
} else {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam-1)] <- FUN1(ParamIn[, 1:(NParam-1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam ], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
}
if (identical(FUN_MOD, RunModel_CemaNeigeGR4H) | identical(FUN_MOD, RunModel_CemaNeigeGR5H) |
identical(FUN_MOD, RunModel_CemaNeigeGR4J) | identical(FUN_MOD, RunModel_CemaNeigeGR5J) | identical(FUN_MOD, RunModel_CemaNeigeGR6J)) {
if (identical(FUN_MOD, RunModel_CemaNeigeGR4H)) {
FUN1 <- TransfoParam_GR4H
}
if (identical(FUN_MOD, RunModel_CemaNeigeGR5H)) {
FUN1 <- TransfoParam_GR5H
}
if (identical(FUN_MOD, RunModel_CemaNeigeGR4J)) {
FUN1 <- TransfoParam_GR4J
}
if (identical(FUN_MOD, RunModel_CemaNeigeGR5J)) {
FUN1 <- TransfoParam_GR5J
}
if (identical(FUN_MOD, RunModel_CemaNeigeGR6J)) {
FUN1 <- TransfoParam_GR6J
}
if (inherits(CalibOptions, "hysteresis")) {
FUN2 <- TransfoParam_CemaNeigeHyst
} else {
FUN2 <- TransfoParam_CemaNeige
}
if (inherits(CalibOptions, "hysteresis")) {
if (inherits(CalibOptions, "hysteresis") & !inherits(CalibOptions, "SD")) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
......@@ -122,7 +132,8 @@ Calibration_Michel <- function(InputsModel,
}
return(ParamOut)
}
} else {
}
if (!inherits(CalibOptions, "hysteresis") & !inherits(CalibOptions, "SD")) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
......@@ -138,6 +149,40 @@ Calibration_Michel <- function(InputsModel,
return(ParamOut)
}
}
if (inherits(CalibOptions, "hysteresis") & inherits(CalibOptions, "SD")) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam-5)] <- FUN1( ParamIn[, 1:(NParam-5)], Direction)
ParamOut[, (NParam-4):(NParam-1)] <- FUN2( ParamIn[, (NParam-4):(NParam-1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam ], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
if (!inherits(CalibOptions, "hysteresis") & inherits(CalibOptions, "SD")) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam-3)] <- FUN1( ParamIn[, 1:(NParam-3)], Direction)
ParamOut[, (NParam-2):(NParam-1)] <- FUN2( ParamIn[, (NParam-2):(NParam-1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam ], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
}
if (is.null(FUN_TRANSFO)) {
......
......@@ -2,6 +2,7 @@ CreateCalibOptions <- function(FUN_MOD,
FUN_CALIB = Calibration_Michel,
FUN_TRANSFO = NULL,
IsHyst = FALSE,
IsSD = FALSE,
FixedParam = NULL,
SearchRanges = NULL,
StartParamList = NULL,
......@@ -17,6 +18,9 @@ CreateCalibOptions <- function(FUN_MOD,
if (!is.logical(IsHyst) | length(IsHyst) != 1L) {
stop("'IsHyst' must be a logical of length 1")
}
if (!is.logical(IsSD) | length(IsSD) != 1L) {
stop("'IsSD' must be a logical of length 1")
}
##check_FUN_MOD
BOOL <- FALSE
......@@ -75,6 +79,9 @@ CreateCalibOptions <- function(FUN_MOD,
if (IsHyst) {
ObjectClass <- c(ObjectClass, "hysteresis")
}
if (IsSD) {
ObjectClass <- c(ObjectClass, "SD")
}
if (!BOOL) {
stop("incorrect 'FUN_MOD' for use in 'CreateCalibOptions'")
return(NULL)
......@@ -139,11 +146,32 @@ CreateCalibOptions <- function(FUN_MOD,
} else {
FUN2 <- TransfoParam_CemaNeige
}
##_set_FUN_LAG
if (IsSD) {
FUN_LAG <- TransfoParam_LAG
}
##_set_FUN_TRANSFO
if (sum(ObjectClass %in% c("GR4H", "GR5H", "GR4J", "GR5J", "GR6J", "GR2M", "GR1A", "CemaNeige")) > 0) {
FUN_TRANSFO <- FUN1
if (!IsSD) {
FUN_TRANSFO <- FUN1
} else {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam-1)] <- FUN1(ParamIn[, 1:(NParam-1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam ], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
} else {
if (IsHyst) {
if (IsHyst & !IsSD) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
......@@ -151,14 +179,15 @@ CreateCalibOptions <- function(FUN_MOD,
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam - 4)] <- FUN1(ParamIn[, 1:(NParam - 4)], Direction)
ParamOut[, 1:(NParam - 4) ] <- FUN1(ParamIn[, 1:(NParam - 4) ], Direction)
ParamOut[, (NParam - 3):NParam] <- FUN2(ParamIn[, (NParam - 3):NParam], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
} else {
}
if (!IsHyst & !IsSD) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
......@@ -169,7 +198,7 @@ CreateCalibOptions <- function(FUN_MOD,
if (NParam <= 3) {
ParamOut[, 1:(NParam - 2)] <- FUN1(cbind(ParamIn[, 1:(NParam - 2)]), Direction)
} else {
ParamOut[, 1:(NParam - 2)] <- FUN1(ParamIn[, 1:(NParam - 2)], Direction)
ParamOut[, 1:(NParam - 2)] <- FUN1( ParamIn[, 1:(NParam - 2)], Direction)
}
ParamOut[, (NParam - 1):NParam] <- FUN2(ParamIn[, (NParam - 1):NParam], Direction)
if (!Bool) {
......@@ -178,6 +207,44 @@ CreateCalibOptions <- function(FUN_MOD,
return(ParamOut)
}
}
if (IsHyst & IsSD) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
ParamOut[, 1:(NParam - 5) ] <- FUN1( ParamIn[, 1:(NParam - 5) ], Direction)
ParamOut[, (NParam - 4):(NParam - 1)] <- FUN2( ParamIn[, (NParam - 4):(NParam - 1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam ], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
if (!IsHyst & IsSD) {
FUN_TRANSFO <- function(ParamIn, Direction) {
Bool <- is.matrix(ParamIn)
if (!Bool) {
ParamIn <- rbind(ParamIn)
}
ParamOut <- NA * ParamIn
NParam <- ncol(ParamIn)
if (NParam <= 3) {
ParamOut[, 1:(NParam - 3)] <- FUN1(cbind(ParamIn[, 1:(NParam - 3)]), Direction)
} else {
ParamOut[, 1:(NParam - 3)] <- FUN1( ParamIn[, 1:(NParam - 3)], Direction)
}
ParamOut[, (NParam - 2):(NParam - 1)] <- FUN2( ParamIn[, (NParam - 2):(NParam - 1)], Direction)
ParamOut[, NParam ] <- FUN_LAG(ParamIn[, NParam], Direction)
if (!Bool) {
ParamOut <- ParamOut[1, ]
}
return(ParamOut)
}
}
}
}
if (is.null(FUN_TRANSFO)) {
......@@ -228,6 +295,9 @@ CreateCalibOptions <- function(FUN_MOD,
if (IsHyst) {
NParam <- NParam + 2
}
if (IsSD) {
NParam <- NParam + 1
}
##check_FixedParam
if (is.null(FixedParam)) {
......@@ -276,22 +346,22 @@ CreateCalibOptions <- function(FUN_MOD,
if ("GR4H" %in% ObjectClass) {
ParamT <- matrix(c(+5.12, -1.18, +4.34, -9.69,
+5.58, -0.85, +4.74, -9.47,
+6.01, -0.50, +5.14, -8.87), ncol = 4, byrow = TRUE)
+6.01, -0.50, +5.14, -8.87), ncol = 4, byrow = TRUE)
}
if (("GR5H" %in% ObjectClass) & ("interception" %in% ObjectClass)) {
ParamT <- matrix(c(+3.46, -1.25, +4.04, -9.53, -9.34,
+3.74, -0.41, +4.78, -8.94, -3.33,
+4.29, +0.16, +5.39, -7.39, +3.33),ncol=5,byrow=TRUE);
+4.29, +0.16, +5.39, -7.39, +3.33), ncol=5, byrow = TRUE);
}
if (("GR5H" %in% ObjectClass) & !("interception" %in% ObjectClass)) {
ParamT <- matrix(c(+3.28, -0.39, +4.14, -9.54, -7.49,
+3.62, -0.19, +4.80, -9.00, -6.31,
+4.01, -0.04, +5.43, -7.53, -5.33), ncol=5,byrow=TRUE);
+4.01, -0.04, +5.43, -7.53, -5.33), ncol=5, byrow = TRUE);
}
if ("GR4J" %in% ObjectClass) {
ParamT <- matrix(c(+5.13, -1.60, +3.03, -9.05,
+5.51, -0.61, +3.74, -8.51,
+6.07, -0.02, +4.42, -8.06), ncol = 4, byrow = TRUE)
+6.07, -0.02, +4.42, -8.06), ncol = 4, byrow = TRUE)
}
if ("GR5J" %in% ObjectClass) {
ParamT <- matrix(c(+5.17, -1.13, +3.08, -9.37, -7.45,
......@@ -329,12 +399,12 @@ CreateCalibOptions <- function(FUN_MOD,
if (("CemaNeigeGR5H" %in% ObjectClass) & ("interception" %in% ObjectClass)) {
ParamT <- matrix(c(+3.46, -1.25, +4.04, -9.53, -9.34, -9.96, +6.63,
+3.74, -0.41, +4.78, -8.94, -3.33, -9.14, +6.90,
+4.29, +0.16, +5.39, -7.39, +3.33, +4.10, +7.21),ncol=7,byrow=TRUE);
+4.29, +0.16, +5.39, -7.39, +3.33, +4.10, +7.21), ncol = 7, byrow = TRUE);
}
if (("CemaNeigeGR5H" %in% ObjectClass) & !("interception" %in% ObjectClass)) {
ParamT <- matrix(c(+3.28, -0.39, +4.14, -9.54, -7.49, -9.96, +6.63,
+3.62, -0.19, +4.80, -9.00, -6.31, -9.14, +6.90,
+4.01, -0.04, +5.43, -7.53, -5.33, +4.10, +7.21), ncol=7,byrow=TRUE);
+4.01, -0.04, +5.43, -7.53, -5.33, +4.10, +7.21), ncol = 7, byrow = TRUE);
}
if ("CemaNeigeGR4J" %in% ObjectClass) {
ParamT <- matrix(c(+5.13, -1.60, +3.03, -9.05, -9.96, +6.63,
......@@ -358,6 +428,12 @@ CreateCalibOptions <- function(FUN_MOD,
+7.00, +7.00), ncol = 2, byrow = TRUE)
ParamT <- cbind(ParamT, ParamTHyst)
}
if (IsSD) {
ParamTSD <- matrix(c(+1.25,
+2.50,
+5.00), ncol = 1, byrow = TRUE)
ParamT <- cbind(ParamT, ParamTSD)
}
StartParamList <- NULL
StartParamDistrib <- TransfoParam(ParamIn = ParamT, Direction = "TR", FUN_TRANSFO = FUN_TRANSFO)
......
TransfoParam_LAG <- function(ParamIn, Direction) {
## number of model parameters
NParam <- 1L
## check arguments
isVecParamIn <- is.vector(ParamIn)
if (isVecParamIn) {
ParamIn <- matrix(ParamIn, nrow = 1)
}
if (!inherits(ParamIn, "matrix")) {
stop("'ParamIn' must be of class 'matrix'")
}
if (!inherits(Direction, "character") | length(Direction) != 1 | any(!Direction %in% c("RT", "TR"))) {
stop("'Direction' must be a character vector of length 1 equal to 'RT' or 'TR'")
}
if (ncol(ParamIn) != NParam) {
stop(sprintf("the LAG model requires %i parameters", NParam))
}
## transformation
if (Direction == "TR") {
ParamOut <- 20 * (ParamIn + 10) / 20.0
}
if (Direction == "RT") {
ParamOut <- ParamIn * 20.0 / 20 - 10
}
## check output
if (isVecParamIn) {
ParamOut <- as.vector(ParamOut)
}
return(ParamOut)
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment