From bbc65718128329c89999656c4dd0d982a644ddc9 Mon Sep 17 00:00:00 2001
From: Andy Bock <abock@usgs.gov>
Date: Tue, 16 Mar 2021 18:12:17 -0600
Subject: [PATCH] stye changes hydReg -> VPU
---
workspace/02_NHD_navigate.Rmd | 58 +++++++++++++++++++----------------
1 file changed, 31 insertions(+), 27 deletions(-)
diff --git a/workspace/02_NHD_navigate.Rmd b/workspace/02_NHD_navigate.Rmd
index 4c6dfb6..faed329 100644
--- a/workspace/02_NHD_navigate.Rmd
+++ b/workspace/02_NHD_navigate.Rmd
@@ -15,19 +15,17 @@ knitr::opts_chunk$set(error = TRUE)
```{r setup}
# Load data and libraries ----------------------
-ptm<-proc.time()
-
library(nhdplusTools)
library(dplyr)
library(sf)
library(hyRefactor)
library(tidyr)
-#library(hyfabric)
+library(hyfabric)
# Load custom functions
source("R/utils.R")
source("R/NHD_navigate.R")
-source("R/hyRefactor_config.R")
+source("R/hyrefactor_config.R")
# Load user-defined data paths
data_paths <- jsonlite::read_json(file.path("cache", "data_paths.json"))
@@ -35,26 +33,32 @@ data_paths <- jsonlite::read_json(file.path("cache", "data_paths.json"))
nhdplus_path(data_paths$nhdplus_gdb)
staged_nhd <- stage_national_data(nhdplus_data = data_paths$nhdplus_gdb,
output_path = data_paths$nhdplus_dir)
-
```
```{r huc12 POIs}
# Derive or load HUC12 POIs ----------------------
# Still need to deal with some duplicate HUC12s
-if(needs_layer(nav_gpkg, temp_POIs) | needs_layer(nav_gpkg, nhdflowline)) {
-
- # Subset NHD by VPU
- nhd <- VPU_Subset(staged_nhd$flowline, VPU)
-
- # Filter and write dendritic/non-coastal subset to gpkg
- # This will be iterated over to produce the final network after POIs identified
- non_dend <- unique(unlist(lapply(filter(nhd, TerminalFl == 1 & TotDASqKM < min_da_km)
- %>% pull(COMID), NetworkNav, st_drop_geometry(nhd))))
+if(needs_layer(nav_gpkg, temp_POIs)) {
- # Add fields to note dendritic and POI flowlines
- nhd <- nhd %>%
- mutate(dend = ifelse(!COMID %in% non_dend, 1, 0),
- poi = ifelse(!COMID %in% non_dend & TotDASqKM >= min_da_km, 1, 0))
+ if (needs_layer(nav_gpkg, nhdflowline)){
+ # Subset NHD by VPU
+ nhd <- VPU_Subset(staged_nhd$flowline, VPU)
+
+ # Filter and write dendritic/non-coastal subset to gpkg
+ # This will be iterated over to produce the final network after POIs identified
+ non_dend <- unique(unlist(lapply(filter(nhd, TerminalFl == 1 & TotDASqKM < min_da_km)
+ %>% pull(COMID), NetworkNav, st_drop_geometry(nhd))))
+
+ # Add fields to note dendritic and POI flowlines
+ nhd <- nhd %>%
+ mutate(dend = ifelse(!COMID %in% non_dend, 1, 0),
+ poi = ifelse(!COMID %in% non_dend & TotDASqKM >= min_da_km, 1, 0))
+
+ write_sf(nhd, nav_gpkg, nhdflowline)
+ } else {
+ nhd <- read_sf(nav_gpkg, nhdflowline)
+ nhd <- select(nhd, -c(minNet, WB, struct_POI, struct_Net, POI_ID))
+ }
# Join Dave's HUC12 outlets with NHD
HUC12_COMIDs <- read_sf(data_paths$hu12_points_path, "hu_outlets") %>%
@@ -69,9 +73,9 @@ if(needs_layer(nav_gpkg, temp_POIs) | needs_layer(nav_gpkg, nhdflowline)) {
huc12_POIs <- POI_creation(st_drop_geometry(HUC12_COMIDs), filter(nhd, poi == 1), "HUC12")
# Write out geopackage layer representing POIs for given theme
- write_sf(nhd, nav_gpkg, nhdflowline)
write_sf(huc12_POIs, nav_gpkg, temp_POIs)
tmp_POIs <- huc12_POIs
+
} else {
# Load HUC12 POIs as the tmpPOIs if they already exist
tmp_POIs <- read_sf(nav_gpkg, temp_POIs)
@@ -86,7 +90,7 @@ if(needs_layer(nav_gpkg, temp_POIs) | needs_layer(nav_gpkg, nhdflowline)) {
if(all(is.na(tmp_POIs$Type_Gages))) {
# Previously identified streamgages within Gage_Selection.Rmd
- streamgages <- read_sf("cache/Gage_info.gpkg", "Potential_Gages") %>%
+ streamgages <- read_sf("cache/Gages_info.gpkg", "Gages") %>%
select(COMID, Gage_no) %>%
filter(COMID %in% nhd$COMID) %>%
st_drop_geometry() %>%
@@ -213,7 +217,7 @@ if(all(is.na(tmp_POIs$Type_WBOut))) {
switchDiv(., nhd) %>%
select(COMID, WBAREACOMI)
- wbin_COMIDs <- filter(nhd, WB == 1 & DnHydroseq %in% filter(nhd, WB==1)$Hydroseq) %>%
+ wbin_COMIDs <- filter(nhd, WB == 0, DnHydroseq %in% filter(nhd, WB == 1)$Hydroseq) %>%
select(-WBAREACOMI) %>%
switchDiv(., nhd) %>%
inner_join(st_drop_geometry(filter(nhd, minNet == 1)) %>%
@@ -223,10 +227,10 @@ if(all(is.na(tmp_POIs$Type_WBOut))) {
slice(n = 1)
tmp_POIs <- POI_creation(filter(wbout_COMIDs, !COMID %in% wbin_COMIDs$COMID), filter(nhd, poi == 1), "WBOut") %>%
- addType(., tmp_POIs, "WBOut", bind = FALSE)
+ addType(., tmp_POIs, "WBOut")
tmp_POIs <- POI_creation(filter(wbin_COMIDs, !COMID %in% wbout_COMIDs$COMID), filter(nhd, poi == 1), "WBIn") %>%
- addType(., tmp_POIs, "WBIn", bind = FALSE)
+ addType(., tmp_POIs, "WBIn")
write_sf(nhd, nav_gpkg, nhdflowline)
write_sf(tmp_POIs, nav_gpkg, temp_POIs)
@@ -280,11 +284,11 @@ if(needs_layer(nav_gpkg, n_segments)) {
# Sort POIs by Levelpath and Hydrosequence in upstream to downstream order
seg_POIs <- filter(st_drop_geometry(nhd), COMID %in% tmp_POIs$COMID, COMID %in% filter(nhd, minNet == 1)$COMID)
-
- # Attribute reaches that contribute locally to each POI
+ # derive incremental segments from POIs
inc_segs <- segment_increment(filter(nhd, minNet == 1), seg_POIs)
+
nhd_Final <- nhd %>%
- left_join(select(st_drop_geometry(inc_segs), COMID, POI_ID), by = "COMID")
+ left_join(select(inc_segs, COMID, POI_ID), by = "COMID")
# create and write out final dissolved segments
nsegments_fin <- segment_creation(nhd_Final, xWalk)
@@ -326,7 +330,7 @@ if(needs_layer(pois_merge)) {
# Write out geopackage layer representing POIs for given theme
write_sf(out_HUC12$allPOIs, nav_gpkg, pois_merge)
- write_sf(out_HUC12$FL, nav_gpkg, nhdflowline)
+ write_sf(nhd_Final, nav_gpkg, nhdflowline)
write_sf(out_HUC12$segs, nav_gpkg, n_segments)
} else {
final_POIs <- read_sf(col_gpkg, pois_merge)
--
GitLab