test_that("airGR::Calibration should work", { ## loading catchment data data(L0123001) ## preparation of InputsModel object InputsModel <- CreateInputsModel(RunModel_GR4J, DatesR = BasinObs$DatesR, Precip = BasinObs$P, PotEvap = BasinObs$E) ## calibration period selection Ind_Run <- seq(which(format(BasinObs$DatesR, format = "%Y-%m-%d")=="1990-01-01"), which(format(BasinObs$DatesR, format = "%Y-%m-%d")=="1999-12-31")) Ind_WarmUp <- seq(which(format(BasinObs$DatesR, format = "%Y-%m-%d")=="1989-01-01"), which(format(BasinObs$DatesR, format = "%Y-%m-%d")=="1989-12-31")) ## preparation of RunOptions object RunOptions <- CreateRunOptions(RunModel_GR4J, InputsModel = InputsModel, IndPeriod_Run = Ind_Run, IndPeriod_WarmUp = Ind_WarmUp) ## calibration criterion: preparation of the InputsCrit object InputsCrit <- CreateInputsCrit(ErrorCrit_NSE, InputsModel = InputsModel, RunOptions = RunOptions, Obs = BasinObs$Qmm[Ind_Run]) ## preparation of CalibOptions object CalibOptions <- CreateCalibOptions(RunModel_GR4J, FUN_CALIB = Calibration_Michel) ## calibration OutputsCalib <- Calibration(InputsModel = InputsModel, RunOptions = RunOptions, InputsCrit = InputsCrit, CalibOptions = CalibOptions, FUN_MOD = RunModel_GR4J, FUN_CALIB = Calibration_Michel) expect_length(OutputsCalib$ParamFinalR, 4) }) # data set up e <- runCalibration(runRunModel = TRUE, FUN_CRIT = ErrorCrit_NSE) for(x in ls(e)) assign(x, get(x, e)) test_that("Calibrated parameters remains unchanged", { skip_on_cran() InputsCrit <- CreateInputsCrit( InputsModel = InputsModel, RunOptions = RunOptions, Obs = Qobs[IndPeriod_Run,] ) OC <- Calibration( InputsModel = InputsModel, RunOptions = RunOptions, InputsCrit = InputsCrit, CalibOptions = CalibOptions ) ParamFinalR <- extractParam(OutputsCalib) lapply(names(ParamFinalR), function(id) expect_equal(ParamFinalR[[!!id]], ParamMichel[[id]])) }) skip_on_cran() test_that("Calibration with regularization is OK", { InputsCrit <- CreateInputsCrit( InputsModel = InputsModel, RunOptions = RunOptions, Obs = Qobs[IndPeriod_Run,], AprioriIds = c( "54057" = "54032", "54032" = "54001", "54001" = "54095" ), transfo = "sqrt" ) OC <- Calibration( InputsModel = InputsModel, RunOptions = RunOptions, InputsCrit = InputsCrit, CalibOptions = CalibOptions ) ParamLavenne <- extractParam(OC) expect_equal(OC[["54095"]]$CritFinal, ErrorCrit( InputsCrit[["54095"]], RunModel(InputsModel, RunOptions, ParamLavenne)[["54095"]] )$CritValue) OM <- RunModel(InputsModel, RunOptions, ParamLavenne) lapply(names(OC), function(id) { expect_gt( ErrorCrit( InputsCrit[[id]], OM[[id]] )$CritValue, 0.89 ) }) }) test_that("Calibration with Diversion works", { n_div <- rbind(nodes, data.frame(id = "54029", down = "54002", length = 50, area = NA, model = "Diversion")) g_div <- CreateGRiwrm(n_div) Qmin = matrix(1E5, nrow = length(DatesR), ncol = 1) colnames(Qmin) = "54029" Qdiv <- -Qmin IM_div <- CreateInputsModel(g_div, DatesR, Precip, PotEvap, Qinf = Qdiv, Qmin = Qmin) RO_div <- setupRunOptions(IM_div)$RunOptions P_div <- ParamMichel P_div$`54002` <- c(1, ParamMichel$`54002`) IC_div <- CreateInputsCrit( InputsModel = IM_div, RunOptions = RO_div, Obs = Qobs[IndPeriod_Run,], ) CO_div <- CreateCalibOptions(IM_div) OC <- Calibration( InputsModel = IM_div, RunOptions = RO_div, InputsCrit = IC_div, CalibOptions = CO_div ) expect_length(OC$`54002`$ParamFinalR, 5) }) test_that("Derivation and normal connection should return same calibration", { n_2ol <- nodes[nodes$id %in% c("54095", "54001"), ] n_2ol[n_2ol$id %in% c("54095", "54001"), c("down", "length")] <- c(NA, NA) meteoIds <- n_2ol$id n_2ol$area[n_2ol$id == "54001"] <- n_2ol$area[n_2ol$id == "54001"] - n_2ol$area[n_2ol$id == "54095"] n_2ol <- rbind(n_2ol, data.frame(id = "54095", down = "54001", length = 42, area = NA, model = "Diversion"), data.frame(id = "upstream", down = "54095", length = 0, area = NA, model = NA)) g_2ol <- CreateGRiwrm(n_2ol) # Add upstream flow on 54095 that is removed by the Diversion # and derive previously simulated flow in order to get the same Qsim as before Qinf = matrix(0, nrow = length(DatesR), ncol = 2) Qinf[IndPeriod_Run, 1] <- OM_GriwrmInputs[["54095"]]$Qsim_m3 Qinf[IndPeriod_Run, 2] <- - OM_GriwrmInputs[["54095"]]$Qsim_m3 Qinf[IndPeriod_WarmUp, 1] <- OM_GriwrmInputs[["54095"]]$RunOptions$WarmUpQsim_m3 Qinf[IndPeriod_WarmUp, 2] <- - OM_GriwrmInputs[["54095"]]$RunOptions$WarmUpQsim_m3 colnames(Qinf) <- c("upstream", "54095") Qmin = matrix(0, nrow = length(DatesR), ncol = 1) colnames(Qmin) <- "54095" IM_2ol <- CreateInputsModel(g_2ol, DatesR, Precip[, meteoIds], PotEvap[, meteoIds], Qinf = Qinf, Qmin = Qmin) # Copy area of upstream node to downstream node in order to get # correct conversion of Qsim in mm IM_2ol[["54001"]]$BasinAreas[1] <- tail(IM_2ol[["54095"]]$BasinAreas, 1) RO_2ol <- setupRunOptions(IM_2ol)$RunOptions IC_2ol <- CreateInputsCrit( InputsModel = IM_2ol, RunOptions = RO_2ol, Obs = Qobs[IndPeriod_Run,], ) CO_2ol <- CreateCalibOptions(IM_2ol) CO_2ol[["54095"]]$FixedParam[1] <- 1 OC_2ol <- Calibration( InputsModel = IM_2ol, RunOptions = RO_2ol, InputsCrit = IC_2ol, CalibOptions = CO_2ol ) ParamRef <- ParamMichel[names(IM_2ol)] ParamRef[["54095"]] <- c(1, ParamRef[["54095"]]) ParamFinalR <- extractParam(OC_2ol) lapply(names(ParamFinalR), function(id) expect_equal(OC_2ol[[id]]$CritFinal, OutputsCalib[[id]]$CritFinal, tolerance = 1E-5)) #Excepted parameter #2 of GR4J all others are equal (precision 3/1000) lapply(names(ParamFinalR), function(id) expect_equal(ParamFinalR[[!!id]][-3] / ParamRef[[!!id]][-3], rep(1, 4), tolerance = 3E-3)) })