diff --git a/R/ConvertMeteoSD.R b/R/ConvertMeteoSD.R index 3266d0a0924a6d53295a4040b4493e1a704195d4..7c1ba3fc729865a70cf36dfef42ca1e77ff26486 100644 --- a/R/ConvertMeteoSD.R +++ b/R/ConvertMeteoSD.R @@ -32,14 +32,14 @@ ConvertMeteoSD.GRiwrm <- function(x, meteo, ...) { #' @export #' @rdname ConvertMeteoSD ConvertMeteoSD.character <- function(x, griwrm, meteo, ...) { - upperBasins <- !is.na(griwrm$down) & griwrm$down == x & !is.na(griwrm$area) - if(all(!upperBasins)) { + griwrm <- griwrm[getDiversionRows(griwrm, inverse = TRUE), ] + upperIDs <- getUpstreamRunOffIds(x, griwrm) + if(length(upperIDs) == 1) { return(meteo[,x]) } - upperIDs <- griwrm$id[upperBasins] - areas <- griwrm$area[match(c(x, upperIDs), griwrm$id)] + areas <- griwrm$area[match(upperIDs, griwrm$id)] output <- ConvertMeteoSD( - meteo[,c(x, upperIDs), drop = FALSE], + meteo[, upperIDs, drop = FALSE], areas = areas ) return(output) @@ -79,3 +79,22 @@ ConvertMeteoSD.matrix <- function(x, areas, temperature = FALSE, ...) { return(as.matrix(meteoDown, ncol = 1)) } +getUpstreamRunOffIds <- function(id, griwrm) { + griwrm <- griwrm[getDiversionRows(griwrm, inverse = TRUE), ] + upstreamNodeIds <- griwrm$id[griwrm$down == id & !is.na(griwrm$down)] + upstreamRunOffIds <- griwrm$id[griwrm$id %in% upstreamNodeIds & !is.na(griwrm$area)] + upstreamNaAreaIds <- upstreamNodeIds[!upstreamNodeIds %in% upstreamRunOffIds] + if(length(upstreamNaAreaIds) > 0) { + upstreamRunOffIds <- c( + upstreamRunOffIds, + unlist(sapply(upstreamNaAreaIds, getUpstreamRunOffIds, griwrm = griwrm)) + ) + upstreamRunOffIds <- upstreamRunOffIds[!is.na(griwrm$area[griwrm$id %in% upstreamRunOffIds])] + } + + if (is.na(griwrm$area[griwrm$id == id])) { + return(upstreamRunOffIds) + } + + return(c(id, upstreamRunOffIds)) +} diff --git a/R/CreateGRiwrm.R b/R/CreateGRiwrm.R index fe5e93fc366b98624b6216368ddf899e1ec56d29..8d888e312fcf6a251c85113b281fdd4fb17a859a 100644 --- a/R/CreateGRiwrm.R +++ b/R/CreateGRiwrm.R @@ -229,9 +229,8 @@ nodeError <- function(node, s) { #' #' @noRd getGaugedId <- function(id, griwrm) { - np <- getNodeProperties(id, griwrm) - if (np$RunOff && np$calibration == "Gauged") { - # Match with au gauged station! + if (isNodeGauged(id, griwrm, skip_reservoirs = TRUE)) { + # Match with a gauged station! return(id) } else { # Otherwise we need to search downstream on the natural network diff --git a/R/getNodeProperties.R b/R/getNodeProperties.R index 8a6e6509a32ef14f261bf7ef033580b2dec45c3e..2f144b28b3cb98b203ac2dd1c7cebfb997c4c1e4 100644 --- a/R/getNodeProperties.R +++ b/R/getNodeProperties.R @@ -1,39 +1,56 @@ #' Get list of properties of a GRiwrm network node #' +#' #' @param id [character] Id of the node in the GRiwrm object #' @param griwrm \[GRiwrm object\] describing the network of the semi-distributed model (See [CreateGRiwrm]) #' #' @return A [list] with the following items: #' - "position" ([character]): Position of the node in the network ("Upstream" or "Intermediate") -#' - "calibration" ([character]): describe if the node is a "Gauged", or an "Ungauged" station, -#' modelled with an hydrological model, or "NA" otherwise -#' - "Upstream" ([logical]): is the node an upstream node? #' - "DirectInjection" ([logical]): is the node a Direct Injection node? #' - "Diversion" ([logical]): is the node a Diversion node? +#' - "Reservoir" ([logical]): is the node a Reservoir? +#' - "airGR" ([logical]): is the node contains an airGR model? +#' - "calibration" ([character]): describe if the node is a "Gauged", or an "Ungauged" station, +#' (see details), or "NA" otherwise +#' - "Upstream" ([logical]): is the node an upstream node? +#' - "RunOff" ([logical]): is the node contains an hydrological model? +#' +#' @details +#' A "Gauged" node is either a node containing a model that is already +#' calibrated (parameters are already fixed) or a node containing a model where +#' observations are available for calibration. +#' +#' A "Ungauged" node is a node containing a model which derives its parameters from +#' another "donor" node. #' #' @export #' #' @example man-examples/getNodeProperties.R getNodeProperties <- function(id, griwrm) { - stopifnot(inherits(griwrm, "GRiwrm")) + stopifnot(inherits(griwrm, "GRiwrm"), + "donor" %in% names(griwrm)) g2 <- griwrm[getDiversionRows(griwrm, TRUE), , drop = FALSE] upstreamIds <- griwrm$id[!griwrm$id %in% griwrm$down] model <- g2$model[g2$id == id] + if (is.na(model)) { + donor_model <- "DirectInjection" + } else { + donor_model <- g2$model[g2$id == g2$donor[g2$id == id]] + } p <- list( position = ifelse(id %in% upstreamIds, "Upstream", "Intermediate"), DirectInjection = is.na(model), Diversion = "Diversion" %in% griwrm$model[griwrm$id == id], - Reservoir = !is.na(model) && model == "RunModel_Reservoir" + Reservoir = !is.na(model) && model == "RunModel_Reservoir", + airGR = grepl("RunModel_", donor_model) ) if (p$DirectInjection) { p$calibration <- "NA" - } else if (model == "Ungauged") { - p$calibration <- "Ungauged" } else { - p$calibration <- "Gauged" + p$calibration <- ifelse(isNodeGauged(id, griwrm), "Gauged", "Ungauged") } p$Upstream <- p$position == "Upstream" - p$RunOff <- !p$DirectInjection && !p$Reservoir + p$RunOff <- !p$DirectInjection && !p$Reservoir && donor_model != "RunModel_Lag" return(p) } @@ -58,3 +75,19 @@ getAllNodesProperties <- function(griwrm) { rownames(df) <- uids return(df) } + + +#' Does this node contains a gauged model? +#' +#' @inheritParams getNodeProperties +#' +#' @return `TRUE` if the node contains a gauged model, `FALSE` if not. +#' @noRd +#' +isNodeGauged <- function(id, griwrm, skip_reservoirs = FALSE) { + g2 <- griwrm[getDiversionRows(griwrm, inverse = TRUE), , drop = FALSE] + model <- g2$model[g2$id == id] + ungaugedModels <- c("Ungauged") + if (skip_reservoirs) ungaugedModels <- c(ungaugedModels, "RunModel_Reservoir") + return(!is.na(model) && !model %in% ungaugedModels) +} diff --git a/man/getNodeProperties.Rd b/man/getNodeProperties.Rd index b8f1126ec6cf0424ba5b100148899e5fd92bea88..d9f7fd38ac3a7d9618832c753ccb12511f17185f 100644 --- a/man/getNodeProperties.Rd +++ b/man/getNodeProperties.Rd @@ -15,16 +15,27 @@ getNodeProperties(id, griwrm) A \link{list} with the following items: \itemize{ \item "position" (\link{character}): Position of the node in the network ("Upstream" or "Intermediate") -\item "calibration" (\link{character}): describe if the node is a "Gauged", or an "Ungauged" station, -modelled with an hydrological model, or "NA" otherwise -\item "Upstream" (\link{logical}): is the node an upstream node? \item "DirectInjection" (\link{logical}): is the node a Direct Injection node? \item "Diversion" (\link{logical}): is the node a Diversion node? +\item "Reservoir" (\link{logical}): is the node a Reservoir? +\item "airGR" (\link{logical}): is the node contains an airGR model? +\item "calibration" (\link{character}): describe if the node is a "Gauged", or an "Ungauged" station, +(see details), or "NA" otherwise +\item "Upstream" (\link{logical}): is the node an upstream node? +\item "RunOff" (\link{logical}): is the node contains an hydrological model? } } \description{ Get list of properties of a GRiwrm network node } +\details{ +A "Gauged" node is either a node containing a model that is already +calibrated (parameters are already fixed) or a node containing a model where +observations are available for calibration. + +A "Ungauged" node is a node containing a model which derives its parameters from +another "donor" node. +} \examples{ ############################################################################### # Severn network with : # diff --git a/tests/testthat/test-ConvertMeteoSD.R b/tests/testthat/test-ConvertMeteoSD.R index 4054650db15eca35f8f56a420a225b408101eafa..18d0b6995e8534937f506d21d5be7dfba1931331 100644 --- a/tests/testthat/test-ConvertMeteoSD.R +++ b/tests/testthat/test-ConvertMeteoSD.R @@ -71,8 +71,26 @@ nodes <- ) griwrm <- CreateGRiwrm(nodes) +test_that("getUpstreamRunoffIds works", { + expect_equal(getUpstreamRunOffIds("Up", griwrm), "Up") + expect_equal(getUpstreamRunOffIds("Down", griwrm), c("Down", "Up")) + + nodes <- + data.frame( + id = c("Up1", "Up2", "Inter", "Up3", "Down"), + down = c("Inter", "Inter", "Down", "Down", NA), + area = c(1.1, 1.2, NA, 2.1, 5), + length = c(0, 0, 0, 0, NA), + model = "RunModel_GR4J", + stringsAsFactors = FALSE + ) + nodes$model[3] <- NA + griwrm <- CreateGRiwrm(nodes) + + expect_equal(sort(getUpstreamRunOffIds("Down", griwrm)), sort(nodes$id[nodes$id != "Inter"])) +}) + test_that("Downstream data should return 2", { - resultMatrix[,] <- 2 meteo <- matrix(rep(c(1, 1.5), dataNumRows), ncol = 2, byrow = TRUE) colnames(meteo) <- c("Up", "Down") expect_equivalent( diff --git a/tests/testthat/test-RunModel.R b/tests/testthat/test-RunModel.R index 611caec457960c86c2967db8ffeaeac093da49bf..d3f39504b41aa91046e89fda5fd958e149742467 100644 --- a/tests/testthat/test-RunModel.R +++ b/tests/testthat/test-RunModel.R @@ -223,3 +223,26 @@ test_that("Upstream node - equal Diversion should return same results", { expect_equal(OM_airGR$Qsimdown, OM_GriwrmInputs[[!!id]]$Qsimdown) expect_equal(OM_airGR$Qsim_m3, OM_GriwrmInputs[[!!id]]$Qsim_m3) }) + +test_that("A transparent upstream node with area=NA should return same result #124", { + nodes <- loadSevernNodes() + nodes <- nodes[nodes$id %in% c("54095", "54001"), ] + nodes[nodes$id == "54001", c("down", "length")] <- c(NA, NA) + nodes$down[nodes$id == "54095"] <- "P" + nodes$length[nodes$id == "54095"] <- 0 + nodes <- rbind( + nodes, + data.frame(id = "P", down = "54001", length = 42, area = NA, model = "RunModel_Lag") + ) + g <- CreateGRiwrm(nodes) + e <- setupRunModel(griwrm = g, runRunModel = FALSE) + for(x in ls(e)) assign(x, get(x, e)) + Param <- ParamMichel[c("54095", "54001")] + Param[["P"]] <- 1 + OM <- RunModel(InputsModel, + RunOptions = RunOptions, + Param = Param) + expect_equal(OM[["P"]]$Qsim_m3, OM[["54095"]]$Qsim_m3) + names(OM[["54001"]]$StateEnd$SD[[1]]) <- "54095" # For matching upstream IDs with ref + expect_equal(OM[["54001"]], OM_GriwrmInputs[["54001"]]) +}) diff --git a/tests/testthat/test-RunModel_Ungauged.R b/tests/testthat/test-RunModel_Ungauged.R index c85301e9fabf2fc2d2bbd9c57a4bba5fba5e7cac..dee9d632cc6a680cd8daaada13c5ca85d5aefd05 100644 --- a/tests/testthat/test-RunModel_Ungauged.R +++ b/tests/testthat/test-RunModel_Ungauged.R @@ -225,3 +225,67 @@ test_that("Ungauged node with diversion outside the sub-network should work", { expect_equal(OC1[[id]]$CritFinal, CritValue) }) }) + +test_that("Ungauged node with upstream node with diversion should work", { + nodes <- loadSevernNodes() + nodes <- nodes[nodes$id %in% c("54095", "54001", "54032"), ] + nodes[nodes$id == "54032", c("down", "length")] <- c(NA, NA) + nodes$model[nodes$id == "54001"] <- "Ungauged" + nodes$down[nodes$id == "54095"] <- "P" + nodes$length[nodes$id == "54095"] <- 0 + nodes <- rbind(nodes, + data.frame(id = "P", down = "54001", length = 42, area = NA, model = "RunModel_Lag"), + data.frame(id = c("54095", "P", "54001"), + down = NA, + length = NA, + area = NA, + model = "Diversion")) + g <- CreateGRiwrm(nodes) + Qobs2 <- matrix(0, ncol = length(g$id[g$model == "Diversion"]), nrow = 11536) + colnames(Qobs2) <- g$id[g$model == "Diversion"] + e <- setupRunModel(griwrm = g, runRunModel = FALSE, Qobs2 = Qobs2) + for(x in ls(e)) assign(x, get(x, e)) + np <- getAllNodesProperties(g) + + IC <- CreateInputsCrit( + InputsModel, + FUN_CRIT = ErrorCrit_KGE2, + RunOptions = RunOptions, + Obs = Qobs[IndPeriod_Run, np$id[np$RunOff & np$calibration == "Gauged"], drop = FALSE], + ) + + CO <- CreateCalibOptions(InputsModel) + CO[["P"]]$FixedParam = 1 + OC_Lag <- Calibration(InputsModel, RunOptions, IC, CO) + Param_Lag <- sapply(OC_Lag, "[[", "ParamFinalR") + Param_Lag[["P"]] <- NULL + expect_equal(Param_Lag, Param) +}) + +test_that("Donor node with diversion should work", { + nodes <- rbind(nodes, + data.frame(id = "54032", + down = NA, + length = NA, + area = NA, + model = "Diversion")) + g <- CreateGRiwrm(nodes) + Qobs2 <- matrix(0, ncol = length(g$id[g$model == "Diversion"]), nrow = 11536) + colnames(Qobs2) <- g$id[g$model == "Diversion"] + e <- setupRunModel(griwrm = g, runRunModel = FALSE, Qobs2 = Qobs2) + for(x in ls(e)) assign(x, get(x, e)) + np <- getAllNodesProperties(g) + + IC <- CreateInputsCrit( + InputsModel, + FUN_CRIT = ErrorCrit_KGE2, + RunOptions = RunOptions, + Obs = Qobs[IndPeriod_Run, np$id[np$RunOff & np$calibration == "Gauged"], drop = FALSE], + ) + + CO <- CreateCalibOptions(InputsModel) + OC_Div <- Calibration(InputsModel, RunOptions, IC, CO) + Param_Div <- sapply(OC_Div, "[[", "ParamFinalR") + expect_equal(Param_Div, Param) + expect_equal(OC$`54032`$CritFinal, OC_Div$`54032`$CritFinal) +})