Skip to content
Snippets Groups Projects
Commit 5a8caa61 authored by Bock, Andy's avatar Bock, Andy
Browse files

WB fixes

parent 69e2df09
No related branches found
No related tags found
1 merge request!165Fix Duplicate POIs, DA Thresh, Waterbody Fixes
Hide whitespace changes
Inline Side-by-side
workspace/R/NHD_navigate.R
+ 82
− 33
View file @ 5a8caa61
......@@ -409,8 +409,10 @@ movePOI_NA_DA <- function(poi_fix, nhdDF){
# 2/3 - correpsonding segments (both raw and dissolved)
POI_move_down<-function(out_gpkg, POIs, Seg, Seg2, exp, DA_diff){
nexus <- filter(POIs, nexus == TRUE)
POIs <- POIs %>% mutate_if(is.numeric, function(x) ifelse(is.infinite(x), NA, x)) %>%
filter(nexus != TRUE)
filter(nexus != 1)
# downstream POIs and their drainage area ratios and incremental drainage areas
segs_down <- Seg %>%
......@@ -530,7 +532,8 @@ POI_move_down<-function(out_gpkg, POIs, Seg, Seg2, exp, DA_diff){
# Bring new POI data over to new data
fin_POIs <- filter(POIs, !COMID %in% c(merged_POIs$oldPOI, merged_POIs$COMID)) %>%
rbind(merged_POIs %>% select(-c(DirDA, oldPOI)))
rbind(merged_POIs %>% select(-c(DirDA, oldPOI))) %>%
rbind(nexus)
#inner_join(st_drop_geometry(merged_POIs) %>% select(COMID, Type_Gages_A), by = "COMID") %>%
#mutate(Type_Gages = ifelse(!is.na(Type_Gages_A), Type_Gages_A, Type_Gages)) %>% select(-Type_Gages_A)
......@@ -832,6 +835,7 @@ HR_Area_coms <- function(nhd, WBs, data_paths, crs){
# Download NHD HR HUC04 if we dont' have it, other wise load and
# Bind NHDHR waterbodies into one layer
hr_wb <- do.call("rbind", lapply(unique(huc04), function(x){
print(x)
if(!file.exists(file.path(data_paths$nhdplus_dir, vpu,
paste0("NHDPLUS_H_", x, "_HU4_GDB.gdb")))){
......@@ -848,7 +852,8 @@ HR_Area_coms <- function(nhd, WBs, data_paths, crs){
st_zm(.) %>%
st_as_sf() %>%
mutate(source = "NHDHR") %>%
st_transform(crs)
st_transform(crs) #%>%
#st_compatibalize(wb)
# Bind or create new object
if(exists("wb")){
......@@ -897,12 +902,14 @@ WB_event <- function(WBs, nhd_wb, type){
outlet_fl <- filter(nhd_wb, COMID %in% filter(WBs_table, outlet == "outlet")$comid)
# Get the downstream flowline for continuity
ds_fl <- filter(nhd_wb, DnHydroseq %in% outlet_fl$Hydroseq,
ds_fl <- filter(nhd_wb, Hydroseq %in% outlet_fl$DnHydroseq,
LevelPathI %in% outlet_fl$LevelPathI) %>%
rbind(outlet_fl) %>%
group_by(WBAREACOMI) %>%
arrange(desc(Hydroseq)) %>%
group_by(LevelPathI) %>%
# union together
summarize(do_union = T)
summarize(do_union = T) %>%
st_cast("LINESTRING")
WBs_area <- filter(WBs_layer, source == "NHDv2Area")
WBs_HR <- filter(WBs_layer, source == "NHDHR")
......@@ -911,13 +918,29 @@ WB_event <- function(WBs, nhd_wb, type){
if (nrow(WBs_area) > 0){
# Intersect unioned FL with waterbody and get intersecting point
outlet_pnts <- sf::st_intersection(ds_fl, WBs_area) %>%
# Cast to point
outlet_pnts_int <- sf::st_intersection(ds_fl, WBs_area) #%>%
# # Cast to point
# st_cast("POINT")
# group_by(LevelPathI) %>%
# # Get furthest downstream point of intersection
# filter(row_number() == max(row_number(), na.rm = T)) %>%
# ungroup()
outlet_point <- outlet_pnts_int[st_geometry_type(outlet_pnts_int) %in% c("POINT"),]
if(nrow(outlet_point) > 0){
outlet_fl <- outlet_pnts_int[st_geometry_type(outlet_pnts_int) %in% c("LINESTRING", "MULTILINESTRING"),] %>%
st_cast("LINESTRING") %>%
st_cast("POINT")
outlet_pnts_int <- rbind(outlet_point, outlet_fl)
} else {
outlet_pnts_int <- st_cast(outlet_pnts_int, "MULTILINESTRING")
}
outlet_pnts <- outlet_pnts_int %>%
group_by(LevelPathI) %>%
st_cast("POINT") %>%
group_by(WBAREACOMI) %>%
# Get furthest downstream point of intersection
filter(row_number() == max(row_number(), na.rm = T)) %>%
ungroup()
mutate(id = row_number()) %>%# Was comid
filter(id == max(id))
# Derive event from point to use for splitting within hyRefactor
wb_events <- get_flowline_index(nhd_wb, outlet_pnts) %>%
......@@ -1002,29 +1025,50 @@ WB_event <- function(WBs, nhd_wb, type){
inlet_FL <- nhd_wb[lengths(st_intersects(start_pts, WBs_layer)) == 0,] %>%
rbind(filter(nhd_wb, Hydroseq %in% .$DnHydroseq,
LevelPathI %in% .$LevelPathI))
LevelPathI %in% .$LevelPathI)) %>%
arrange(desc(Hydroseq)) %>%
unique()
inlet_ls <- inlet_FL %>%
group_by(LevelPathI) %>%
st_cast("POINT") %>%
summarize(do_union = F) %>%
st_cast("LINESTRING")
st_cast("LINESTRING") %>%
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()
#inlet_pnts <- sf::st_intersection(inlet_ls, st_cast(WBs_layer, "MULTILINESTRING", group_or_split = FALSE)) #%>%
inlet_pnts_int <- sf::st_intersection(inlet_ls, WBs_layer) #%>%
#st_cast("LINESTRING")
int_point <- inlet_pnts_int[st_geometry_type(inlet_pnts_int) %in% c("POINT"),]
if(nrow(int_point) > 0){
inlet_fl <- inlet_pnts_int[st_geometry_type(inlet_pnts_int) %in% c("LINESTRING", "MULTILINESTRING"),] %>%
st_cast("LINESTRING")
inlet_pnts_int <- rbind(int_point, inlet_fl)
}
inlet_pnts <- inlet_pnts_int %>%
group_by(LevelPathI) %>%
st_cast("POINT") %>%
mutate(id = row_number()) %>%# Was comid
filter(id == min(id))
# # the intersection can return linestring features.
# if(sf::st_geometry_type(inlet_pnts, by_geometry = FALSE) %in% c("LINESTRING", "MULTILINESTRING")) {
# inlet_pnts2 <- inlet_pnts %>%
# #st_cast("LINESTRING") %>%
# st_cast("POINT") %>%
# group_by(COMID) %>%
# mutate(id = row_number()) %>%# Was comid
# filter(id == 1)
# #filter(row_number() == 1) %>%
# #filter(row_number() == min(row_number())) %>%
# 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") %>%
......@@ -1062,7 +1106,7 @@ gage_POI_creation <- function(tmp_POIs, gages_info, nhd, combine_meters, reach_m
# Reset gage flag if even is created
gage_POIs_nonevent <- filter(gage_POIs, !Type_Gages %in% events$Type_Gages) %>%
addType(., tmp_POIs, "Gages")
addType(., tmp_POIs, "Gages", nexus = FALSE, bind = TRUE)
tmp_POIs <- data.table::rbindlist(list(gage_POIs_nonevent,
select(events, COMID, Type_Gages, nexus)), fill = TRUE) %>%
......@@ -1206,12 +1250,16 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs,
# Drive the inlet events
nhd_inlet <- filter(nhd_WBs, !COMID %in% wb_outlet_events$COMID)
if(unique(nhd$VPUID) == "12"){
nhd_inlet <- rbind(nhd_inlet, filter(nhd_WBs, COMID == 1304709))
}
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(wbin_POIs) %>%
mutate(nexus = TRUE) %>%
inner_join(select(st_drop_geometry(nhd), COMID, Hydroseq, ToMeas), by = "COMID")
inner_join(select(st_drop_geometry(nhd), COMID, Hydroseq, ToMeas), by = "COMID")
# Determine which events are close enough to the end of an upstream flowline to just use that geometry for the POI
wbin_fl_upstream <- filter(nhd, DnHydroseq %in% filter(nhd, COMID %in% wb_inlet_events$COMID)$Hydroseq) %>%
......@@ -1227,8 +1275,8 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs,
mutate(dsCOMID = COMID, COMID = usCOMID)
if(nrow(wb_inlet_POIs) > 0) {
wbin_POIs <- bind_rows(POI_creation(select(st_drop_geometry(wb_inlet_POIs), dsCOMID, Type_WBIn = WBAREACOMI),
nhd, "WBIn"), wbin_POIs)
wbin_POIs <- bind_rows(wbin_POIs, POI_creation(select(st_drop_geometry(wb_inlet_POIs), dsCOMID, Type_WBIn = WBAREACOMI),
nhd, "WBIn"))
wb_inlet_events <- filter(wb_inlet_events, !COMID %in% wb_inlet_POIs$dsCOMID)
}
......@@ -1237,6 +1285,7 @@ wbin_POIcreation <- function(nhd, WBs_VPU, data_paths, crs,
if(nrow(wb_inlet_events) > 0){
wbin_POIs <- addType(wbin_POIs, tmp_POIs, "WBIn")
wb_inlet_events <- st_compatibalize(wb_inlet_events, wbin_POIs)
wbin_POIs_fin <- data.table::rbindlist(list(wbin_POIs,
select(wb_inlet_events, COMID,
Type_WBIn = WBAREACOMI, nexus)),
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment