diff --git a/R/Calibration.GRiwrmInputsModel.R b/R/Calibration.GRiwrmInputsModel.R index 5eb509c86f0633055e50999e8c54c26a8b74449a..88715a5aea436db88fef884a112c8e25d09a8c73 100644 --- a/R/Calibration.GRiwrmInputsModel.R +++ b/R/Calibration.GRiwrmInputsModel.R @@ -48,6 +48,7 @@ Calibration.GRiwrmInputsModel <- function(InputsModel, l <- updateParameters4Ungauged(id, InputsModel, RunOptions, + CalibOptions, OutputsModel, useUpstreamQsim) IM <- l$InputsModel @@ -84,14 +85,24 @@ Calibration.GRiwrmInputsModel <- function(InputsModel, subBasinAreas <- calcSubBasinAreas(IM) } for (uId in Ids) { - # Add OutputsCalib for ungauged nodes - OutputsCalib[[uId]] <- OutputsCalib[[id]] - # Copy parameters and transform X4 relatively to the sub-basin area - OutputsCalib[[uId]]$ParamFinalR <- - OutputsCalib[[uId]]$ParamFinalR[IM[[uId]]$model$indexParamUngauged] - if(IM[[id]]$model$hasX4) { - OutputsCalib[[uId]]$ParamFinalR[IM[[uId]]$model$iX4] <- - X4 * (subBasinAreas[uId] / sum(subBasinAreas)) ^ 0.3 + if(!IM[[uId]]$isReservoir) { + # Add OutputsCalib for ungauged nodes + OutputsCalib[[uId]] <- OutputsCalib[[id]] + # Copy parameters and transform X4 relatively to the sub-basin area + OutputsCalib[[uId]]$ParamFinalR <- + OutputsCalib[[uId]]$ParamFinalR[IM[[uId]]$model$indexParamUngauged] + if(IM[[id]]$model$hasX4) { + OutputsCalib[[uId]]$ParamFinalR[IM[[uId]]$model$iX4] <- + X4 * (subBasinAreas[uId] / sum(subBasinAreas, na.rm = TRUE)) ^ 0.3 + } + } else { + OutputsCalib[[uId]] <- Calibration( + InputsModel = IM[[uId]], + RunOptions = RunOptions[[uId]], + InputsCrit = IC, + CalibOptions = CalibOptions[[uId]], + ... + ) } } IM <- IM[[id]] @@ -171,6 +182,7 @@ reduceGRiwrmObj4Ungauged <- function(griwrm, obj) { updateParameters4Ungauged <- function(GaugedId, InputsModel, RunOptions, + CalibOptions, OutputsModel, useUpstreamQsim) { @@ -188,9 +200,15 @@ updateParameters4Ungauged <- function(GaugedId, ### Modify InputsModel for the reduced network ### # Remove nodes outside of reduced network InputsModel <- reduceGRiwrmObj4Ungauged(g, InputsModel) + # Copy fixed parameters for Reservoirs + for (id in names(InputsModel)) { + if (InputsModel[[id]]$isReservoir) { + InputsModel[[id]]$FixedParam <- CalibOptions[[id]]$FixedParam + } + } # Update griwrm attr(InputsModel, "GRiwrm") <- g - # Update Qupstream already modelled in the reduced network upstream nodes + # Update Qupstream already modeled in the reduced network upstream nodes idIM <- unique(g$down[g$id %in% upIds]) for (id in idIM) { if(useUpstreamQsim && any(InputsModel[[id]]$UpstreamIsModeled)) { @@ -258,9 +276,12 @@ calcSubBasinAreas <- function(IM) { #' @noRd RunModel_Ungauged <- function(InputsModel, RunOptions, Param) { InputsModel$FUN_MOD <- NULL - SBVI <- sum(calcSubBasinAreas(InputsModel)) + SBVI <- sum(calcSubBasinAreas(InputsModel), na.rm = TRUE) # Compute Param for each sub-basin P <- lapply(InputsModel, function(IM) { + if (IM$isReservoir) { + return(IM$FixedParam) + } p <- Param[IM$model$indexParamUngauged] if(IM$model$hasX4) { p[IM$model$iX4] <- Param[IM$model$iX4] * (IM$BasinAreas[length(IM$BasinAreas)] / SBVI) ^ 0.3 diff --git a/tests/testthat/test-RunModel_Reservoir.R b/tests/testthat/test-RunModel_Reservoir.R index f4b781eb934ed033b2ff1782dfc396127455116f..f1c36289b3b79795939b7df7268d2d6ffdbd8ec9 100644 --- a/tests/testthat/test-RunModel_Reservoir.R +++ b/tests/testthat/test-RunModel_Reservoir.R @@ -79,4 +79,22 @@ test_that("Calibration with ungauged node and reservoir in the middle works", { expect_equal(g$donor[g$id == "54095"], "54001") + e <- setupRunModel(griwrm = g, runRunModel = FALSE, Qobs2 = Qobs2) + for(x in ls(e)) assign(x, get(x, e)) + + InputsCrit <- CreateInputsCrit(InputsModel, + ErrorCrit_KGE2, + RunOptions = RunOptions, + Obs = Qobs[IndPeriod_Run, ]) + CalibOptions <- CreateCalibOptions(InputsModel) + CalibOptions[["Dam"]]$FixedParam <- c(650E6, 1) + OC <- Calibration( + InputsModel = InputsModel, + RunOptions = RunOptions, + InputsCrit = InputsCrit, + CalibOptions = CalibOptions + ) + # X1, X2, X3 are identical + expect_equal(OC$`54001`$ParamFinalR[2:4], OC$`54095`$ParamFinalR[1:3]) + expect_equal(OC$Dam$ParamFinalR, CalibOptions[["Dam"]]$FixedParam) })