Merge branch '130-ungauged-node-diversion-to-reservoir-crashes-calibration-2' into 'dev'

Resolve "Ungauged node: Diversion to Reservoir crashes Calibration"

Closes #130

See merge request !77

NAMESPACE

@@ -23,6 +23,8 @@ S3method(RunModel,GRiwrmInputsModel)
 S3method(RunModel,InputsModel)
 S3method(RunModel,SD)
 S3method(RunModel,Supervisor)
+S3method(isNodeDownstream,GRiwrm)
+S3method(isNodeDownstream,GRiwrmInputsModel)
 S3method(plot,GRiwrm)
 S3method(plot,GRiwrmOutputsModel)
 S3method(plot,OutputsModelReservoir)
@@ -44,13 +46,16 @@ export(getNodeProperties)
 export(getNodeRanking)
 export(getSD_Ids)
 export(isNodeDownstream)
+export(isNodeUpstream)
 export(plot.Qm3s)
 import(airGR)
+import(dplyr)
 importFrom(dplyr,"%>%")
 importFrom(grDevices,rainbow)
 importFrom(graphics,matplot)
 importFrom(graphics,par)
 importFrom(graphics,plot)
 importFrom(graphics,title)
+importFrom(stats,setNames)
 importFrom(utils,read.table)
 importFrom(utils,tail)

R/CreateCalibOptions.GRiwrmInputsModel.R

 #' @rdname CreateCalibOptions
 #' @export
+#' @importFrom stats setNames
 CreateCalibOptions.GRiwrmInputsModel <- function(x, FixedParam = NULL, ...) {
   dots <- list(...)
   if ("IsHyst" %in% names(dots)) {

R/CreateGRiwrm.R

@@ -125,34 +125,42 @@ CheckColumnTypes <- function(df, coltypes, keep_all) {
 }
 
 
-#' Sorting of the nodes from upstream to downstream
+#' Sorting of the nodes from upstream to downstream for RunModel and Calibration
 #'
 #' @param griwrm \[object of class `GRiwrm`\] see [CreateGRiwrm] for details
 #'
-#' @return [numeric] ordered node names
+#' @return [numeric] ordered node ids
 #' @export
+#' @import dplyr
 getNodeRanking <- function(griwrm) {
   if (!inherits(griwrm, "GRiwrm")) {
     stop("getNodeRanking: griwrm argument should be of class GRiwrm")
   }
   # Remove upstream nodes without model (direct flow connections)
-  griwrm <- griwrm[!is.na(griwrm$model), ]
-  # Rank 1
-  rank <- setdiff(griwrm$id, griwrm$down)
-  ranking <- rank
-  # Next ranks
-  while (any(griwrm$id %in% rank)) {
-    rank <- griwrm$down[griwrm$id %in% rank]
-    ranking <- c(ranking, rank)
+  g <- griwrm[!is.na(griwrm$model), ]
+  r <- c()
+  while (nrow(g) > 0) {
+    # Search for gauged ids or ungauged with upstream donor
+    upIds <- unique(g$id[!g$id %in% g$down & (g$id == g$donor | !g$donor %in% g$id)])
+    r <- c(r, upIds)
+    g <- g[!g$id %in% upIds, ]
+    #Search for ungauged ids
+    upIds <- unique(g$id[!g$id %in% g$down & g$id != g$donor])
+    while(length(upIds) > 0) {
+      upId <- upIds[1]
+      #Browse the ungauged sub-network until the donor
+      upDonor <- g$donor[g$id == upId]
+      g2 <- g %>% filter(donor == upDonor)
+      g2$donor <- g2$id
+      ungaugedIds <- getNodeRanking(g2)
+      upIds <- upIds[!upIds %in% ungaugedIds]
+      r <- c(r, ungaugedIds)
+      g <- g[!g$id %in% ungaugedIds, ]
+    }
   }
-  ranking <- unique(ranking, fromLast = TRUE)
-  ranking <- ranking[-length(ranking)]
-  # Remove intermediate nodes without model (direct flow connections)
-  ranking <- ranking[ranking %in% griwrm$id]
-  return(ranking)
+  return(r)
 }
 
-
 checkNetworkConsistency <- function(db) {
   db2 <- db[getDiversionRows(db, TRUE), ]
   if (any(duplicated(db2$id))) {
@@ -268,10 +276,7 @@ setDonor <- function(griwrm) {
   sapply(seq(nrow(griwrm)), function(i) {
     id <- griwrm$id[i]
     model <- griwrm$model[i]
-    if (is.na(model) || model == "Diversion") {
-      # Diversion and Direct injection are "Non Applicable"
-      return(NA)
-    } else if(model == "RunModel_Reservoir" && is.na(griwrm$down[i])){
+    if(model == "RunModel_Reservoir" && is.na(griwrm$down[i])){
       # RunModel_Reservoir needs to be its own "donor" only if at downstream
       # Otherwise we search the first gauged station downstream to allow
       # calibration with ungauged upstream nodes

R/CreateInputsModel.GRiwrm.R

@@ -431,24 +431,14 @@ getInputBV <- function(x, id, unset = NULL) {
 #'
 #' @noRd
 hasUngaugedNodes <- function(id, griwrm) {
-  nps <- getAllNodesProperties(griwrm)
-  upNodes <- griwrm[!is.na(griwrm$down) & griwrm$down == id, ]
-  upIds <- upNodes$id[upNodes$model != "Diversion"]
-  # No upstream nodes
-  if(length(upIds) == 0) return(FALSE)
-  # At least one upstream node is ungauged
-  UngNodes <- griwrm$model[griwrm$id %in% upIds] == "Ungauged"
-  UngNodes <- UngNodes[!is.na(UngNodes)]
-  if(length(UngNodes) > 0 && any(UngNodes)) return(TRUE)
-
-  upNps <- nps[nps$id %in% upIds, ]
-  if(any(upNps$DirectInjection) || any(upNps$Reservoir)) {
-    # At least one node's model is NA or Reservoir, we need to investigate next level
-    out <- sapply(upNps$id[upNps$DirectInjection | upNps$Reservoir],
-                  hasUngaugedNodes,
-                  griwrm = griwrm)
-    return(any(out))
-  }
+  g <- griwrm[!is.na(griwrm$model), ]
+  idsWithCurrentAsDonor <- g$id[g$id != id & g$donor == id]
+  if (length(idsWithCurrentAsDonor) == 0) return(FALSE)
+  areNodesUpstream <- sapply(idsWithCurrentAsDonor,
+                             function(x) isNodeUpstream(g, id, x))
+  if (!any(areNodesUpstream)) return(FALSE)
+  g_red <- g[g$id %in% idsWithCurrentAsDonor[areNodesUpstream], ]
+  if (any(g_red$model == "Ungauged")) return(TRUE)
   return(FALSE)
 }
 

R/globals.R

+#' getNodeRanking: no visible binding for global variable 'donor'
+#' updateParameters4Ungauged: no visible binding for global variable
+#' down'
+#' updateParameters4Ungauged: no visible binding for global variable
+#' model'
+utils::globalVariables(c("donor", "down ", "model "))

R/utils.GRiwrm.R

@@ -33,18 +33,40 @@ getNoSD_Ids <- function(InputsModel, include_diversion = TRUE) {
 }
 
 
-#' Check if a node is downstream another one
+#' Check if a node is downstream or upstream another one
 #'
-#' @param InputsModel \[`GRiwrmInputsModel` object\] see [CreateInputsModel.GRiwrm] for details
+#' @param x \[`GRiwrmInputsModel` object\] (see [CreateInputsModel.GRiwrm]) or
+#'        \[`GRiwrm` object\] (See [CreateGRiwrm])
 #' @param current_node [character] with the id of the current node
-#' @param down_node [character] with the id of the node for which we want to know if it is downstream `current_node`
+#' @param candidate_node [character] with the id of the node for which we want
+#'        to know if it is downstream or upstream `current_node`
 #'
-#' @return [logical] `TRUE` if the node with the id `down_node` is downstream the node with the id `current_node`
+#' @return [logical] `TRUE` if the node with the id `down_candidate` is downstream
+#'         or upstream the node with the id `current_node`
 #' @export
+#' @rdname isNodeDownstream
 #'
-isNodeDownstream <- function(InputsModel, current_node, down_node) {
-  current_down_node <- InputsModel[[current_node]]$down
-  if (is.na(current_down_node)) return(FALSE)
-  if (current_down_node == down_node) return(TRUE)
-  return(isNodeDownstream(InputsModel, current_down_node, down_node))
+isNodeDownstream <- function(x, current_node, candidate_node) {
+  UseMethod("isNodeDownstream", x)
+}
+
+#' @export
+#' @rdname isNodeDownstream
+isNodeDownstream.GRiwrmInputsModel <- function(x, current_node, candidate_node) {
+  isNodeDownstream(attr(x, "GRiwrm"), current_node, candidate_node)
+}
+
+#' @export
+#' @rdname isNodeDownstream
+isNodeDownstream.GRiwrm <- function(x, current_node, candidate_node) {
+  current_down_node <- x$down[x$id %in% current_node]
+  if (all(is.na(current_down_node))) return(FALSE)
+  if (any(current_down_node == candidate_node)) return(TRUE)
+  return(isNodeDownstream(x, current_down_node, candidate_node))
+}
+
+#' @export
+#' @rdname isNodeDownstream
+isNodeUpstream <- function(x, current_node, candidate_node) {
+  !isNodeDownstream(x, current_node, candidate_node)
 }

man/getNodeRanking.Rd

@@ -2,7 +2,7 @@
 % Please edit documentation in R/CreateGRiwrm.R
 \name{getNodeRanking}
 \alias{getNodeRanking}
-\title{Sorting of the nodes from upstream to downstream}
+\title{Sorting of the nodes from upstream to downstream for RunModel and Calibration}
 \usage{
 getNodeRanking(griwrm)
 }
@@ -10,8 +10,8 @@ getNodeRanking(griwrm)
 \item{griwrm}{[object of class \code{GRiwrm}] see \link{CreateGRiwrm} for details}
 }
 \value{
-\link{numeric} ordered node names
+\link{numeric} ordered node ids
 }
 \description{
-Sorting of the nodes from upstream to downstream
+Sorting of the nodes from upstream to downstream for RunModel and Calibration
 }

man/isNodeDownstream.Rd

@@ -2,20 +2,32 @@
 % Please edit documentation in R/utils.GRiwrm.R
 \name{isNodeDownstream}
 \alias{isNodeDownstream}
-\title{Check if a node is downstream another one}
+\alias{isNodeDownstream.GRiwrmInputsModel}
+\alias{isNodeDownstream.GRiwrm}
+\alias{isNodeUpstream}
+\title{Check if a node is downstream or upstream another one}
 \usage{
-isNodeDownstream(InputsModel, current_node, down_node)
+isNodeDownstream(x, current_node, candidate_node)
+
+\method{isNodeDownstream}{GRiwrmInputsModel}(x, current_node, candidate_node)
+
+\method{isNodeDownstream}{GRiwrm}(x, current_node, candidate_node)
+
+isNodeUpstream(x, current_node, candidate_node)
 }
 \arguments{
-\item{InputsModel}{[\code{GRiwrmInputsModel} object] see \link{CreateInputsModel.GRiwrm} for details}
+\item{x}{[\code{GRiwrmInputsModel} object] (see \link{CreateInputsModel.GRiwrm}) or
+[\code{GRiwrm} object] (See \link{CreateGRiwrm})}
 
 \item{current_node}{\link{character} with the id of the current node}
 
-\item{down_node}{\link{character} with the id of the node for which we want to know if it is downstream \code{current_node}}
+\item{candidate_node}{\link{character} with the id of the node for which we want
+to know if it is downstream or upstream \code{current_node}}
 }
 \value{
-\link{logical} \code{TRUE} if the node with the id \code{down_node} is downstream the node with the id \code{current_node}
+\link{logical} \code{TRUE} if the node with the id \code{down_candidate} is downstream
+or upstream the node with the id \code{current_node}
 }
 \description{
-Check if a node is downstream another one
+Check if a node is downstream or upstream another one
 }

tests/testthat/helper_RunModel_Reservoir.R

@@ -56,3 +56,20 @@ getGriwrmDerivedReservoirUngauged <- function(donorByDerivation) {
     g$donor[g$id == "54095"] <- "54029"
   return(g)
 }
+
+testDerivedUngauged <- function(donorByDerivation) {
+  g <- getGriwrmDerivedReservoirUngauged(donorByDerivation)
+  Qobs2 <- matrix(-1E9, ncol = 2, nrow = 11536)
+  colnames(Qobs2) <- c("54095", "Dam")
+  Qobs2[, "54095"] <- -1E9
+  Qobs2[, "Dam"] <- 1E9
+  e <- setupRunModel(griwrm = g, runRunModel = FALSE, Qobs2 = Qobs2)
+  for (x in ls(e)) assign(x, get(x, e))
+
+  CalibOptions <- CreateCalibOptions(InputsModel,
+                                     FixedParam = list(Dam = c(650E6, 1)))
+  e <- runCalibration(g, Qobs2 = Qobs2, CalibOptions = CalibOptions)
+  for(x in ls(e)) assign(x, get(x, e))
+  expect_equal(Param[["54095"]][1:3],
+               Param[[ifelse(donorByDerivation, "54029", "54001")]][2:4])
+}

tests/testthat/test-RunModel_Ungauged.R

@@ -290,27 +290,10 @@ test_that("Cemaneige with hysteresis works",  {
                c("54057" = 9, "54032" = 9, "54001" = 8))
 })
 
-testDerivdedUngauged <- function(donorByDerivation) {
-  g <- getGriwrmDerivedReservoirUngauged(donorByDerivation)
-  Qobs2 <- matrix(-1E9, ncol = 2, nrow = 11536)
-  colnames(Qobs2) <- c("54095", "Dam")
-  Qobs2[, "54095"] <- -1E9
-  Qobs2[, "Dam"] <- 1E9
-  e <- setupRunModel(griwrm = g, runRunModel = FALSE, Qobs2 = Qobs2)
-  for (x in ls(e)) assign(x, get(x, e))
-
-  CalibOptions <- CreateCalibOptions(InputsModel,
-                                     FixedParam = list(Dam = c(650E6, 1)))
-  e <- runCalibration(g, Qobs2 = Qobs2, CalibOptions = CalibOptions)
-  for(x in ls(e)) assign(x, get(x, e))
-  expect_equal(Param[["54095"]][1:3],
-               Param[[ifelse(donorByDerivation, "54029", "54001")]][2:4])
-}
-
 test_that("Ungauged node with derivation to reservoir should work", {
-  testDerivdedUngauged(FALSE)
+  testDerivedUngauged(FALSE)
 })
 
-# test_that("Ungauged node with donor by derivation through reservoir should work", {
-#   testDerivdedUngauged(TRUE)
-# })
+test_that("Ungauged node with donor by derivation through reservoir should work", {
+  testDerivedUngauged(TRUE)
+})

tests/testthat/test-createGRiwrm.R

@@ -68,7 +68,7 @@ test_that("Derivated ungauged node without downstream node should have derivated
   nodes <- rbind(nodes,
                  data.frame(id = "54001", down = "54032", length = 45, area = NA, model = "Diversion"))
   g <- CreateGRiwrm(nodes)
-  expect_equal(g$donor, c(rep("54032", 3), NA))
+  expect_equal(g$donor, rep("54032", 4))
 })
 
 test_that("Reservoir between ungauged and gauged node should have the first downstream node as donor", {

tests/testthat/test-getNodeRanking.R

@@ -24,7 +24,12 @@ test_that("Check ranking with Diversion", {
   expect_lt(which(r == "54029"), which(r == "54002"))
 })
 
-# test_that("Check ranking with Ungauged node, reservoir, and Diversion #130", {
-#   g <- getGriwrmDerivedReservoirUngauged(TRUE)
-#   expect_equal(getNodeRanking(g), c("54029", "Dam", "54029", "54001", "54032"))
-# })
+test_that("Check ranking with Ungauged node, reservoir, and Diversion #130", {
+  g <- getGriwrmDerivedReservoirUngauged(FALSE)
+  expect_equal(getNodeRanking(g), c("54095", "54001", "Dam", "54029", "54032"))
+})
+
+test_that("Check ranking with Ungauged node, reservoir, and Diversion #130", {
+  g <- getGriwrmDerivedReservoirUngauged(TRUE)
+  expect_equal(getNodeRanking(g), c("54095", "Dam", "54029", "54001", "54032"))
+})