diff --git a/workspace/02_NHD_navigate.Rmd b/workspace/02_NHD_navigate.Rmd index d6183b8ac13a39f3a9ba99a05fb3e155dbfb233c..435bec078d1a0e24e62034e26c87daf55f78b68e 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)