From 1b97d68638c2c6e39a51dbd3fe6d31dfa4bac15f Mon Sep 17 00:00:00 2001
From: Andy Bock <abock@usgs.gov>
Date: Mon, 18 Apr 2022 10:17:22 -0600
Subject: [PATCH] Update to elev POIS

---
 workspace/02_NHD_navigate.Rmd | 71 +++++++++++++++++------------------
 1 file changed, 35 insertions(+), 36 deletions(-)

diff --git a/workspace/02_NHD_navigate.Rmd b/workspace/02_NHD_navigate.Rmd
index d6183b8..435bec0 100644
--- a/workspace/02_NHD_navigate.Rmd
+++ b/workspace/02_NHD_navigate.Rmd
@@ -166,7 +166,7 @@ mapview(filter(tmp_POIs, Type_TE != ""), layer.name = "TE Plant POIs", col.regio
 
 ```{r Terminal POIs}
 #  Derive or load Terminal POIs ----------------------
-if(!"Type_Term" %in% colnames(tmp_POIs)) {
+if(all(is.na(tmp_POIs$Type_Term))) {
   
   # Terminal POIs on levelpaths with upstream POIs
   term_paths <- filter(st_drop_geometry(nhd), Hydroseq %in% filter(nhd, COMID %in% tmp_POIs$COMID)$TerminalPa)
@@ -349,7 +349,6 @@ mapview(filter(tmp_POIs, !is.na(Type_WBOut)), layer.name = "Waterbody outlet POI
 
 This chunk breaks up potential aggregated segments where the elevation range of the are contributing to the segment exceeds 500-m
 ```{r Elevation Break POIs}
-
 # derive incremental segments from POIs
 inc_segs <- segment_increment(filter(nhd, minNet == 1), 
                               filter(st_drop_geometry(nhd),
@@ -358,34 +357,34 @@ inc_segs <- segment_increment(filter(nhd, minNet == 1),
   inner_join(select(atts, COMID, DnHydroseq, VA_MA, TOTMA, LENGTHKM, MAXELEVSMO, MINELEVSMO, WBAREACOMI, WBAreaType, FTYPE,
                     AreaSqKM, TotDASqKM), by = "COMID") 
 
-if(all(is.na(tmp_POIs$Type_Elev))) {
-  
-  # Build dataframe for creation of points based on breaks
-  elev_pois_split <- inc_segs %>%
-    group_by(POI_ID) %>%
-    # Get elevation info
-    mutate(MAXELEVSMO = na_if(MAXELEVSMO, -9998), MINELEVSMO = na_if(MINELEVSMO, -9998), 
-           elev_diff_seg = max(MAXELEVSMO) - min(MINELEVSMO)) %>%
-    # Filter out to those incremental segments that need splitting above the defined elevation threshold
-    filter((max(MAXELEVSMO) - min(MINELEVSMO)) > (elev_diff * 100)) %>%
-    # Do not aggregate on fake lfowlines
-    filter(AreaSqKM > 0, TotDASqKM > max_elev_TT_DA) %>% 
-    # Order from upstream to downstream
-    arrange(desc(Hydroseq)) %>%
-    # Get attributes used for splitting
-    # csum_length = cumulative length US to DS along segment, cumsum_elev = cumulative US to DS elev diff along segment
-    # inc_elev = elevation range of each segment
-    # iter = estimated number of times we'd need to split existing agg flowlines, or number of rows in set
-    mutate(csum_length = cumsum(LENGTHKM),
-           inc_elev = cumsum(MAXELEVSMO - MINELEVSMO),
-           #nc_elev_diff = c(inc_elev[1], (inc_elev - lag(inc_elev))[-1]),
-           iter = min(round(max(inc_elev) / (elev_diff * 100)), n()),
-           elev_POI_ID = NA) %>%
-    # Remove segments we can't break
-    filter(iter > 0, n() > 1) %>%
-    # Track original iteration number
-    mutate(orig_iter = iter) %>%
-    ungroup() 
+# Build dataframe for creation of points based on breaks
+elev_pois_split <- inc_segs %>%
+  group_by(POI_ID) %>%
+  # Get elevation info
+  mutate(MAXELEVSMO = na_if(MAXELEVSMO, -9998), MINELEVSMO = na_if(MINELEVSMO, -9998), 
+         elev_diff_seg = max(MAXELEVSMO) - min(MINELEVSMO)) %>%
+  # Filter out to those incremental segments that need splitting above the defined elevation threshold
+  filter((max(MAXELEVSMO) - min(MINELEVSMO)) > (elev_diff * 100)) %>%
+  # Do not aggregate on fake lfowlines
+  filter(AreaSqKM > 0, TotDASqKM > max_elev_TT_DA) %>% 
+  # Order from upstream to downstream
+  arrange(desc(Hydroseq)) %>%
+  # Get attributes used for splitting
+  # csum_length = cumulative length US to DS along segment, cumsum_elev = cumulative US to DS elev diff along segment
+  # inc_elev = elevation range of each segment
+  # iter = estimated number of times we'd need to split existing agg flowlines, or number of rows in set
+  mutate(csum_length = cumsum(LENGTHKM),
+         inc_elev = cumsum(MAXELEVSMO - MINELEVSMO),
+         #nc_elev_diff = c(inc_elev[1], (inc_elev - lag(inc_elev))[-1]),
+         iter = min(round(max(inc_elev) / (elev_diff * 100)), n()),
+         elev_POI_ID = NA) %>%
+  # Remove segments we can't break
+  filter(iter > 0, n() > 1) %>%
+  # Track original iteration number
+  mutate(orig_iter = iter) %>%
+  ungroup() 
+
+if(nrow(elev_pois_split) > 0) {
   
   # For reach iteration
   elev_POIs <- do.call("rbind", lapply(unique(elev_pois_split$POI_ID), split_elev, 
@@ -395,13 +394,14 @@ if(all(is.na(tmp_POIs$Type_Elev))) {
     select(COMID, Type_Elev) %>%
     unique()
   
-  tmp_POIs <- POI_creation(select(elev_POIs, COMID, Type_Elev), nhd, "Elev") %>%
-    mutate(Type_Term = NA) %>%
-    addType(., tmp_POIs, "Elev")
-
+  if(nrow(elev_POIs) > 0){
+    tmp_POIs <- POI_creation(select(elev_POIs, COMID, Type_Elev), nhd, "Elev") %>%
+      addType(., tmp_POIs, "Elev")
+  } 
+  
   write_sf(tmp_POIs, nav_gpkg, nav_poi_layer)
 } else {
-  tmp_POIs <- read_sf(nav_gpkg, nav_poi_layer)
+  tmp_POIs <- read_sf(nav_gpkg, nav_poi_layer) 
   nhd <- read_sf(nav_gpkg, nhd_flowline)
 }
 
@@ -450,7 +450,6 @@ if(all(is.na(tmp_POIs$Type_Travel))) {
     unique()
     
   tmp_POIs <- POI_creation(select(tt_POIs, COMID, Type_Travel), nhd, "Travel") %>%
-    mutate(Type_Term = NA) %>%
     addType(., tmp_POIs, "Travel")
   
   write_sf(tmp_POIs, nav_gpkg, nav_poi_layer)
-- 
GitLab