diff --git a/NAMESPACE b/NAMESPACE index b5558de..8e74001 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -6,6 +6,7 @@ S3method(dbplyr_edition,OGRSQLConnection) S3method(st_as_sf,tbl_OGRSQLConnection) export(OGRSQL) export(accumulate_downstream) +export(add_hydroseq) export(as_ogr) export(collect) export(st_as_sf) @@ -58,6 +59,8 @@ importFrom(httr,GET) importFrom(httr,progress) importFrom(httr,write_disk) importFrom(igraph,as_ids) +importFrom(igraph,dfs) +importFrom(igraph,graph_from_data_frame) importFrom(igraph,graph_from_edgelist) importFrom(igraph,is_dag) importFrom(igraph,topo_sort) diff --git a/R/network_properties.R b/R/network_properties.R index 614e1cb..015b9cc 100644 --- a/R/network_properties.R +++ b/R/network_properties.R @@ -92,3 +92,60 @@ accumulate_downstream <- function(x, id = "flowpath_id", toid = "flowpath_toid # Return totals aligned to input rows as.numeric(total[idx]) } + + +#' Compute and add the hydrosequence to a directed acyclic network. +#' +#' @param topology A data frame (or tibble) containing at least the identifier column +#' given by `id` and the downstream pointer column given by `toid`. +#' @param id Character scalar. Column name in `topology` with unique node identifiers. +#' Defaults to `"flowpath_id"`. +#' @param toid Character scalar. Column name in `topology` with the *downstream* node +#' identifier for each row. Use `NA` or `0` for outlets/terminals. +#' Defaults to `"flowpath_toid"`. +#' @param colname Character scalar. Column name to use in result. +#' Defaults to `"hydroseq"` +#' +#' @returns The data frame `topology` with an additional column, named `colname`, representing the hydrosequence. +#' @importFrom igraph dfs graph_from_data_frame +#' @export +add_hydroseq <- function(topology, id = "flowpath_id", toid = "flowpath_toid", colname = "hydroseq") { + # Create a _transposed_ network, where traversing the network + # is equivalent to traversing the hydrological network upstream. + # + # This assumes the outlets of this network all connect to an + # ephemeral "0" node (forming a rooted tree network). + edgelist <- topology[, c(toid, id)] + names(edgelist) <- c("id", "toid") + edgelist$id[is.na(edgelist$id)] <- "0" + + # TODO: Check if multiple components exist. If they do, then + # we need to add "0" edges for each component not rooted on "0". + + # Perform DFS from each terminal upstream to get a + # distinct topological sort for the hydrosequence. + sorted <- data.frame( + node = as.integer( + names( + igraph::dfs( + igraph::graph_from_data_frame(edgelist), + root = "0", + mode = "out" + )$order + ) + ) + ) + + sorted$hydroseq <- c(0, seq_len(nrow(sorted) - 1)) + + # Merge the initial hydrosequence to the edgelist and handle ties in the hydrosequence. + result <- merge(edgelist, sorted, by.x = "id", by.y = "node", all.x = TRUE) + result <- result[!is.na(result$hydroseq), ] + result <- result[order(result$hydroseq, result$id, result$toid), c("toid", "id")] + result$hydroseq <- seq_len(nrow(result)) + names(result) <- c(id, toid, "hydroseq") + + # Arrange into input order + topology[[colname]] <- result$hydroseq[match(topology[[id]], result[[id]])] + topology +} \ No newline at end of file diff --git a/man/add_hydroseq.Rd b/man/add_hydroseq.Rd new file mode 100644 index 0000000..e58ff16 --- /dev/null +++ b/man/add_hydroseq.Rd @@ -0,0 +1,33 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/network_properties.R +\name{add_hydroseq} +\alias{add_hydroseq} +\title{Compute and add the hydrosequence to a directed acyclic network.} +\usage{ +add_hydroseq( + topology, + id = "flowpath_id", + toid = "flowpath_toid", + colname = "hydroseq" +) +} +\arguments{ +\item{topology}{A data frame (or tibble) containing at least the identifier column +given by `id` and the downstream pointer column given by `toid`.} + +\item{id}{Character scalar. Column name in `topology` with unique node identifiers. +Defaults to `"flowpath_id"`.} + +\item{toid}{Character scalar. Column name in `topology` with the *downstream* node +identifier for each row. Use `NA` or `0` for outlets/terminals. +Defaults to `"flowpath_toid"`.} + +\item{colname}{Character scalar. Column name to use in result. +Defaults to `"hydroseq"`} +} +\value{ +The data frame `topology` with an additional column, named `colname`, representing the hydrosequence. +} +\description{ +Compute and add the hydrosequence to a directed acyclic network. +}