From 5d184fe0e3853e52d8f39bf416537e2e87de2f68 Mon Sep 17 00:00:00 2001
From: David <david.dorchies@inrae.fr>
Date: Tue, 20 Jun 2023 15:02:45 +0200
Subject: [PATCH] test: calibration with ungauged nodes and hysteresis version
of cemaneige
Refs #134
---
tests/testthat/helper_1_RunModel.R | 13 ++++++-
tests/testthat/helper_Calibration.R | 12 ++++--
tests/testthat/test-RunModel_Ungauged.R | 49 +++++++++++++++++++++++++
3 files changed, 69 insertions(+), 5 deletions(-)
diff --git a/tests/testthat/helper_1_RunModel.R b/tests/testthat/helper_1_RunModel.R
index 565d840..3e9acbe 100644
--- a/tests/testthat/helper_1_RunModel.R
+++ b/tests/testthat/helper_1_RunModel.R
@@ -15,7 +15,8 @@ setupRunModel <-
runRunOptions = TRUE,
runRunModel = TRUE,
griwrm = NULL,
- Qobs2 = NULL) {
+ Qobs2 = NULL,
+ IsHyst = FALSE) {
data(Severn)
@@ -43,12 +44,20 @@ setupRunModel <-
# Convert meteo data to SD (remove upstream areas)
Precip <- ConvertMeteoSD(griwrm, PrecipTot)
PotEvap <- ConvertMeteoSD(griwrm, PotEvapTot)
+ if (IsHyst) {
+ TempMean <- PotEvap+5 # Fake temperatures
+ } else {
+ TempMean <- NULL
+ }
# set up inputs
if (!runInputsModel)
return(environment())
InputsModel <-
- suppressWarnings(CreateInputsModel(griwrm, DatesR, Precip, PotEvap, Qobs = Qobs2))
+ suppressWarnings(CreateInputsModel(griwrm, DatesR, Precip, PotEvap,
+ TempMean = TempMean,
+ Qobs = Qobs2,
+ IsHyst = IsHyst))
# RunOptions
if (!runRunOptions)
diff --git a/tests/testthat/helper_Calibration.R b/tests/testthat/helper_Calibration.R
index c292811..08dc636 100644
--- a/tests/testthat/helper_Calibration.R
+++ b/tests/testthat/helper_Calibration.R
@@ -1,14 +1,19 @@
runCalibration <- function(nodes = NULL,
Qobs2 = NULL,
InputsCrit = NULL,
+ CalibOptions = NULL,
FUN_CRIT = ErrorCrit_KGE2,
- runRunModel = FALSE) {
+ runRunModel = FALSE,
+ IsHyst = FALSE) {
if (is.null(nodes)) {
griwrm <- NULL
} else {
griwrm <- CreateGRiwrm(nodes)
}
- e <- setupRunModel(griwrm = griwrm, runRunModel = runRunModel, Qobs2 = Qobs2)
+ e <- setupRunModel(griwrm = griwrm,
+ runRunModel = runRunModel,
+ Qobs2 = Qobs2,
+ IsHyst = IsHyst)
for(x in ls(e)) assign(x, get(x, e))
rm(e)
np <- getAllNodesProperties(griwrm)
@@ -22,7 +27,8 @@ runCalibration <- function(nodes = NULL,
)
}
- CalibOptions <- CreateCalibOptions(InputsModel)
+ if (is.null(CalibOptions)) CalibOptions <- CreateCalibOptions(InputsModel)
+
OutputsCalib <- Calibration(InputsModel, RunOptions, InputsCrit, CalibOptions)
Param <- sapply(OutputsCalib, "[[", "ParamFinalR")
return(environment())
diff --git a/tests/testthat/test-RunModel_Ungauged.R b/tests/testthat/test-RunModel_Ungauged.R
index 1a03a4f..1030fd1 100644
--- a/tests/testthat/test-RunModel_Ungauged.R
+++ b/tests/testthat/test-RunModel_Ungauged.R
@@ -233,3 +233,52 @@ test_that("Donor node with diversion should work", {
for(x in ls(e)) assign(x, get(x, e))
expect_equal(OC_ref$`54032`$CritFinal, OutputsCalib$`54032`$CritFinal, tolerance = 1E-3)
})
+
+test_that("Cemaneige with hysteresis works", {
+ nodes <- loadSevernNodes()
+ nodes <- nodes[nodes$id %in% c("54057", "54032", "54001"), ]
+ nodes$model <- "RunModel_CemaNeigeGR4J"
+ nodes$model[nodes$id != 54057] <- "Ungauged"
+ griwrm <- CreateGRiwrm(nodes)
+
+ # # The custom ErrorCrit function !!!
+ ErrorCrit_KGE3 <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
+ OutputsCritQ <- suppressMessages(
+ ErrorCrit_KGE2(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE)
+ )
+ InputsCrit$Obs <- InputsCrit$SCA #a adapter
+ OutputsModel$Qsim <- OutputsModel$CemaNeigeLayers[[1]]$Gratio #a adapter
+ OutputsCritSCA <- suppressMessages(
+ ErrorCrit_KGE2(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE)
+ )
+ OutputsCritQ$CritValue <-
+ (OutputsCritQ$CritValue + OutputsCritSCA$CritValue) / 2
+ OutputsCritQ$CritName <- "(0.5 * KGE2[Q] + 0.5 * KGE2[SCA])"
+ return(OutputsCritQ)
+ }
+ class(ErrorCrit_KGE3) <- c("FUN_CRIT", class(ErrorCrit_KGE3))
+
+ e <- setupRunModel(griwrm = griwrm, runRunModel = FALSE, IsHyst = TRUE)
+ for(x in ls(e)) assign(x, get(x, e))
+ np <- getAllNodesProperties(griwrm)
+ InputsCrit <- CreateInputsCrit(
+ InputsModel,
+ FUN_CRIT = ErrorCrit_KGE3,
+ RunOptions = RunOptions,
+ Obs = Qobs[IndPeriod_Run, np$id[np$RunOff & np$calibration == "Gauged"], drop = FALSE],
+ )
+ InputsCrit$`54057`$SCA <- runif(length(IndPeriod_Run)) # Fake SCA
+ CalibOptions <- CreateCalibOptions(InputsModel)
+ CO <- lapply(CalibOptions, function(x) {
+ x$StartParamList <- matrix(
+ c(0.605, 320.596, -0.042, 37.991, 2.221, 0.705, 6.764, 85.000, 0.850),
+ nrow = 1)
+ x$StartParamDistrib <- NULL
+ x
+ })
+ class(CO) <- class(CalibOptions)
+ e <- runCalibration(nodes, InputsCrit = InputsCrit, CalibOptions = CO, IsHyst = TRUE)
+ for(x in ls(e)) assign(x, get(x, e))
+ expect_equal(sapply(Param, length),
+ c("54057" = 9, "54032" = 9, "54001" = 8))
+})
--
GitLab