diff --git a/workspace/02_NHD_navigate.Rmd b/workspace/02_NHD_navigate.Rmd index 22769bcfebe605cc6e87d2a432b413e14a0f785c..cb9ef9b386be04bedbf9d0556a187d24454227a2 100644 --- a/workspace/02_NHD_navigate.Rmd +++ b/workspace/02_NHD_navigate.Rmd @@ -344,7 +344,7 @@ mapview(filter(tmp_POIs, Type_NID != ""), layer.name = "NID POIs", col.regions = if(all(is.na(tmp_POIs$Type_WBIn))) { wb_layers <- wbin_POIcreation(nhd, WBs_VPU, data_paths, crs, split_layer) - tmp_POIs <- wb_layers$POIs + tmp_POIs <- addType(wb_layers$POIs, tmp_POIs, "WBIn") if(!all(is.na(wb_layers$events))) { wb_inlet_events <- select(wb_layers$events, -c(id, offset, Hydroseq, ToMeas)) %>% diff --git a/workspace/R/NHD_navigate.R b/workspace/R/NHD_navigate.R index dda22721e7498d39839542dbbdc44b44a8e5c332..49715be4658fedb74a9ec8676f1fbfb1e5af3b53 100644 --- a/workspace/R/NHD_navigate.R +++ b/workspace/R/NHD_navigate.R @@ -713,6 +713,10 @@ split_tt <- function(in_POI_ID, split_DF, tt_diff, max_DA){ filter(AreaSqKM > 0) %>% mutate(inc_seg_area = c(TotDASqKM[1], (TotDASqKM - lag(TotDASqKM))[-1])) + if(nrow(split_tt_DF) == 0) { + return(NULL) + } + first_iter <- unique(split_tt_DF$iter) # Iterate through segs that need splitting for (i in c(first_iter:max(split_DF$iter, na.rm = T))){ @@ -1005,15 +1009,22 @@ WB_event <- function(WBs, nhd_wb, type){ summarize(do_union = F) %>% st_cast("LINESTRING") - inlet_pnts <- sf::st_intersection(inlet_ls, WBs_layer) %>% - #group_by("LINESTRING") %>% - st_cast("POINT") %>% - st_as_sf() - #mutate(row_id = row_number()) %>% - #group_by(LevelPathI) %>% - #filter(row_id == min(row_id)) %>% - #ungroup() - + inlet_pnts <- sf::st_intersection(inlet_ls, WBs_layer) + + # the intersection can return linestring features. + if(sf::st_geometry_type(inlet_pnts, by_geometry = FALSE) == "LINESTRING") { + inlet_pnts <- inlet_pnts %>% + st_cast("POINT") %>% + st_as_sf() %>% + group_by(COMID) %>% + filter(row_number() == 1) %>% + ungroup() + } else { + inlet_pnts <- inlet_pnts %>% + st_cast("POINT") %>% + st_as_sf() + } + wb_events <- get_flowline_index(nhd_wb, inlet_pnts) %>% inner_join(select(st_drop_geometry(nhd_wb), COMID, WBAREACOMI, LevelPathI), by = "COMID") %>% mutate(event_type = type) %>% @@ -1167,9 +1178,8 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs, split_layer){ } wbout_COMIDs <- filter(tmp_POIs, !is.na(Type_WBOut)) - tmp_POIs <- POI_creation(filter(st_drop_geometry(wbin_COMIDs), !COMID %in% wbout_COMIDs$COMID), - nhd, "WBIn") %>% - addType(., tmp_POIs, "WBIn") + wbin_POIs <- POI_creation(filter(st_drop_geometry(wbin_COMIDs), !COMID %in% wbout_COMIDs$COMID), + nhd, "WBIn") # Get NHDArea and HR waterbodies WBs_VPU_areaHR <- HR_Area_coms(nhd, WBs_VPU, data_paths, crs) @@ -1193,7 +1203,7 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs, split_layer){ if (nrow(nhd_inlet) > 0){ # Get inlet events and bind with outlets wb_inlet_events <- WB_event(WBs_VPU_areaHR, nhd_inlet, "inlet") %>% - st_compatibalize(tmp_POIs) %>% + st_compatibalize(wbin_POIs) %>% mutate(nexus = TRUE) %>% inner_join(select(st_drop_geometry(nhd), COMID, Hydroseq, ToMeas), by = "COMID") @@ -1210,27 +1220,28 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs, split_layer){ inner_join(select(st_drop_geometry(wbin_fl_upstream), usCOMID = COMID, toCOMID), by = c("COMID" = "toCOMID")) %>% mutate(dsCOMID = COMID, COMID = usCOMID) - tmp_POIs <- POI_creation(select(st_drop_geometry(wb_inlet_POIs), COMID, Type_WBIn = WBAREACOMI), - nhd, "WBIn") %>% - addType(., tmp_POIs, "WBIn") + wbin_POIs <- bind_rows(POI_creation(select(st_drop_geometry(wb_inlet_POIs), COMID, Type_WBIn = WBAREACOMI), + nhd, "WBIn"), wbin_POIs) wb_inlet_events <- filter(wb_inlet_events, !COMID %in% wb_inlet_POIs$dsCOMID) } if(nrow(wb_inlet_events) > 0){ - tmp_POIs <- data.table::rbindlist(list(tmp_POIs, - select(wb_inlet_events, COMID, Type_WBIn = WBAREACOMI, nexus)), fill = TRUE) %>% + wbin_POIs <- data.table::rbindlist(list(wbin_POIs, + select(wb_inlet_events, COMID, + Type_WBIn = WBAREACOMI, nexus)), + fill = TRUE) %>% st_as_sf() } - return(list(POIs = tmp_POIs, events = wb_inlet_events)) + return(list(POIs = wbin_POIs, events = wb_inlet_events)) } else { print("no waterbody inlets events") - return(list(POIs = tmp_POIs, events = NA)) + return(list(POIs = NA, events = NA)) } } else { - return(list(POIs = tmp_POIs, events = NA)) + return(list(POIs = NA, events = NA)) } }