get minimal network function

Author: dblodgett-usgs

NAMESPACE

@@ -9,6 +9,7 @@ export(collapse_flowlines)
 export(download_elev)
 export(download_fdr_fac)
 export(flowpaths_to_linestrings)
+export(get_minimal_network)
 export(map_outlet_ids)
 export(reconcile_catchment_divides)
 export(reconcile_collapsed_flowlines)
@@ -43,6 +44,7 @@ importFrom(httr,write_disk)
 importFrom(magrittr,"%>%")
 importFrom(methods,is)
 importFrom(methods,slot)
+importFrom(nhdplusTools,add_plus_network_attributes)
 importFrom(nhdplusTools,get_node)
 importFrom(nhdplusTools,get_vaa)
 importFrom(nhdplusTools,prepare_nhdplus)

R/aggregate_network.R

@@ -496,3 +496,88 @@ get_catchment_sets <- function(flowpath, outlets) {
   list(fline_sets, cat_sets)
 }
 
+#' Get Minimal Network
+#' @description Given a set of outlets, will generate a minimal network by 
+#' calling \code{\link{aggregate_network}} and adding nhdplus attributes to the result.
+#' 
+#' If geometry is included with the network, it will be merged and returned.
+#' 
+#' @inheritParams aggregate_network
+#' @return a data.frame (potentially including an sfc list column) with
+#' attributes generated by \code{\link[nhdplusTools]{add_plus_network_attributes}} 
+#' and a list column "set" containing members of each output flowpath.
+#' @importFrom nhdplusTools add_plus_network_attributes
+#' @export
+#' @examples
+#' source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
+#' fline <- walker_flowline
+#' 
+#' outlets <- data.frame(ID = c(5329357, 5329317, 5329365, 5329435, 5329817),
+#'                       type = c("outlet", "outlet", "outlet", "outlet", "outlet"))
+#' 
+#' #' Add toCOMID
+#' fline[["toCOMID"]] <- nhdplusTools::get_tocomid(fline)
+#' 
+#' # get attributes set
+#' fline <- dplyr::select(fline, ID = COMID, toID = toCOMID,
+#'                        levelpathid = LevelPathI, hydroseq = Hydroseq,
+#'                        areasqkm = AreaSqKM, lengthkm = LENGTHKM)
+#' 
+#' min_net <- get_minimal_network(fline, outlets) 
+#' 
+#' plot(sf::st_geometry(fline), col = "blue")
+#' plot(sf::st_geometry(min_net), lwd = 2, add = TRUE)
+#' plot(sf::st_geometry(nhdplusTools::get_node(min_net)), add = TRUE)
+#' 
+get_minimal_network <- function(flowpath, outlets) {
+  
+  outlets <- add_terminals(flowpath, outlets)
+
+  minimal <- hyRefactor::aggregate_network(
+    flowpath, dplyr::filter(outlets, .data$ID %in% flowpath$ID), 
+    da_thresh = NA, only_larger = TRUE)
+
+  min_net <- tidyr::unnest_longer(drop_geometry(minimal$fline_sets), 
+                                  col = .data$set) %>%
+    left_join(select(flowpath, .data$ID, .data$lengthkm, 
+                     .data$areasqkm, .data$levelpathid), 
+              by = c("set" = "ID")) %>%
+    group_by(ID) %>%
+    summarise(toID = .data$toID[1],
+              lengthkm = sum(.data$lengthkm),
+              areasqkm = sum(.data$areasqkm),
+              outlet_levelpath = min(.data$levelpathid)) %>%
+    mutate(toID = ifelse(is.na(.data$toID), 0, .data$toID)) %>%
+    rename(comid = .data$ID, 
+           tocomid = .data$toID, 
+           nameID = .data$outlet_levelpath)%>%
+    add_plus_network_attributes() %>%
+    rename(ID = .data$comid, toID = .data$tocomid, 
+           outlet_nhdpv2_levelpath = .data$nameID, 
+           arbolate_sum = .data$weight) %>%
+    left_join(select(minimal$fline_sets, .data$ID, .data$set),
+                       by = "ID")
+  
+  if(inherits(minimal$fline_sets, "sf")) {
+    sf::st_sf(min_net)
+  } else {
+    min_net
+  }
+}
+
+add_terminals <- function(flowpath, outlets) {
+  flowpath_sort <- left_join(
+    data.frame(ID = flowpath$ID),
+    nhdplusTools::get_sorted(flowpath[, c("ID", "toID"), drop = TRUE], 
+                             split = TRUE), by = "ID")
+  
+  # Grab terminal paths that matter and combine with outlets..
+  terminal_paths <- unique(flowpath_sort$terminalID[flowpath_sort$ID %in% outlets$ID])
+  flowpath <- flowpath[flowpath_sort$terminalID %in% terminal_paths, ]
+  terminal_paths <- flowpath$ID[flowpath_sort$ID %in% terminal_paths]
+  
+  rbind(outlets,
+        data.frame(ID = terminal_paths,
+                   type = "terminal"))
+}
+

man/get_minimal_network.Rd

@@ -28,3 +28,30 @@ test_that("example runs", {
 
   expect_error(aggregate_network(fline, outlets), "Terminal paths must have an NA or 0 toID")
 })
+
+test_that("minimal network", {
+  source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
+  fline <- walker_flowline
+
+  outlets <- data.frame(ID = c(5329357, 5329317, 5329365, 5329435, 5329817),
+                        type = c("outlet", "outlet", "outlet", "outlet", "outlet"))
+
+  #' Add toCOMID
+  fline[["toCOMID"]] <- nhdplusTools::get_tocomid(fline)
+
+  fline <- dplyr::select(fline, ID = COMID, toID = toCOMID,
+                         levelpathid = LevelPathI, hydroseq = Hydroseq,
+                         areasqkm = AreaSqKM, lengthkm = LENGTHKM)
+
+  min_net <- get_minimal_network(fline, outlets)
+  
+  expect_equal(nrow(min_net), 8)
+  
+  expect_s3_class(min_net, "sf")  
+  
+  min_net <- get_minimal_network(sf::st_drop_geometry(fline), outlets)
+  
+  expect_s3_class(min_net, c("tbl_df","tbl","data.frame"), exact = TRUE) 
+  
+  expect_true(all(outlets$ID %in% min_net$ID))
+})