Skip to content

Commit

Permalink
progress towards v2.2, general cleanup
Browse files Browse the repository at this point in the history
  • Loading branch information
mikejohnson51 committed May 22, 2024
1 parent 52448a1 commit 1630501
Show file tree
Hide file tree
Showing 19 changed files with 458 additions and 380 deletions.
4 changes: 4 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,7 @@ export(middle_massage)
export(network_metadata)
export(pack_set)
export(pinch_sides)
export(prep_split_events)
export(prepare_network)
export(read_hydrofabric)
export(reconcile_catchment_divides)
Expand All @@ -67,9 +68,11 @@ importFrom(dplyr,"%>%")
importFrom(dplyr,`%>%`)
importFrom(dplyr,add_count)
importFrom(dplyr,arrange)
importFrom(dplyr,bind_cols)
importFrom(dplyr,bind_rows)
importFrom(dplyr,case_when)
importFrom(dplyr,collect)
importFrom(dplyr,contains)
importFrom(dplyr,cur_group_id)
importFrom(dplyr,desc)
importFrom(dplyr,distinct)
Expand Down Expand Up @@ -111,6 +114,7 @@ importFrom(logger,log_level)
importFrom(methods,as)
importFrom(nhdplusTools,add_plus_network_attributes)
importFrom(nhdplusTools,calculate_total_drainage_area)
importFrom(nhdplusTools,get_flowline_index)
importFrom(nhdplusTools,get_node)
importFrom(nhdplusTools,get_sorted)
importFrom(nhdplusTools,get_streamorder)
Expand Down
14 changes: 10 additions & 4 deletions R/add_nonnetwork_divides.R
Original file line number Diff line number Diff line change
Expand Up @@ -15,12 +15,18 @@ add_nonnetwork_divides = function(gpkg = NULL,
divides = NULL,
huc12 = NULL,
reference_gpkg = NULL,
reference_divides = NULL,
verbose = TRUE){

if(is.null(reference_gpkg)){ stop('reference_gpkg cannot be NULL')}
if(is.null(reference_gpkg) & is.null(reference_divides)){ stop('reference_gpkg and reference_divides cannot be NULL')}
if(!is.null(reference_gpkg) & !is.null(reference_divides)){ stop('Either reference_gpkg or reference_divides must be NULL')}

ref_nl = read_hydrofabric(reference_gpkg, verbose = verbose, realization = "catchments")
names(ref_nl$catchments) = tolower(names(ref_nl$catchments))
if(!is.null(reference_gpkg)){
reference_divides = read_hydrofabric(reference_gpkg, verbose = verbose, realization = "catchments")[[1]]
names(reference_divides) = tolower(names(reference_divides))
} else {
names(reference_divides) = tolower(names(reference_divides))
}

catchment_name = grep("divide|catchment", st_layers(gpkg)$name, value = TRUE)
catchment_name = catchment_name[!grepl("network", catchment_name)]
Expand All @@ -37,7 +43,7 @@ add_nonnetwork_divides = function(gpkg = NULL,
u_fl = unique(net$hf_id)

# Reference ND catchments
non_network_divides = filter(ref_nl$catchments, !featureid %in% u_fl) %>%
non_network_divides = filter(reference_divides, !featureid %in% u_fl) %>%
select(id = featureid) %>%
st_transform(st_crs(out_nl$catchments)) %>%
mutate(areasqkm = add_areasqkm(.),
Expand Down
14 changes: 4 additions & 10 deletions R/aggregate_along_mainstems.R
Original file line number Diff line number Diff line change
Expand Up @@ -10,12 +10,6 @@ add_network_type = function(network_list, verbose = TRUE){
mutate(has_flowpath = id %in% network_list$flowpaths$id) %>%
filter(!duplicated(.))

if(verbose){
# message("Has Divide")
# print(table(network_list$flowpaths$has_divide))
# message("\nHas Flowpath")
# print(table(network_list$catchments$has_flowpath))
}

network_list
}
Expand Down Expand Up @@ -70,10 +64,10 @@ aggregate_along_mainstems = function(network_list,
mutate(hl_un = ifelse(hl_un %in% hl_dn, NA, hl_un)) %>%
mutate(
ind = cs_group(
.data$areasqkm,
.data$lengthkm,
.data$hl_dn,
.data$hl_un,
areasqkm,
lengthkm,
hl_dn,
hl_un,
ideal_size_sqkm,
min_area_sqkm,
min_length_km
Expand Down
132 changes: 66 additions & 66 deletions R/aggregate_network_to_outlets.R
Original file line number Diff line number Diff line change
Expand Up @@ -89,32 +89,32 @@ aggregate_network_to_outlets <- function(flowpath, outlets,
# OK to assume since Hydroseq is required in validate flowpath
sets <- unnest_flines(fline_sets, "set") %>%
left_join(select(drop_geometry(flowpath),
toSet = .data$toID, .data$ID, .data$Hydroseq, .data$LevelPathID),
toSet = toID, ID, Hydroseq, LevelPathID),
by = c("set" = "ID"))

# Figure out what the ID of the downstream catchment is.
next_id = sets %>%
group_by(.data$ID) %>%
group_by(ID) %>%
#LEVERAGE PRE SORT!
slice_min(.data$Hydroseq) %>%
slice_min(Hydroseq) %>%
ungroup() %>%
select(.data$ID, .data$toSet, .data$set, .data$LevelPathID) %>%
left_join(select(sets, toID = .data$ID, .data$set), by = c("toSet" = "set")) %>%
select(.data$ID, .data$toID)
select(ID, toSet, set, LevelPathID) %>%
left_join(select(sets, toID = ID, set), by = c("toSet" = "set")) %>%
select(ID, toID)


if(inherits(flowpath, "sf")) {
fline_sets <- select(sets, .data$set, .data$ID) %>%
left_join(select(flowpath, .data$ID), by = c("set" = "ID")) %>%
fline_sets <- select(sets, set, ID) %>%
left_join(select(flowpath, ID), by = c("set" = "ID")) %>%
st_as_sf() %>%
filter(!st_is_empty(.)) %>%
union_linestrings(ID = "ID") %>%
left_join(next_id, by = "ID") %>%
left_join(fline_sets, by = "ID") %>%
left_join(distinct(select(sets, .data$ID, .data$LevelPathID)), by = "ID")
left_join(distinct(select(sets, ID, LevelPathID)), by = "ID")
} else {
fline_sets = left_join(fline_sets, next_id, by = "ID") %>%
left_join(distinct(select(sets, .data$ID, .data$LevelPathID)), by = "ID")
left_join(distinct(select(sets, ID, LevelPathID)), by = "ID")
}

cat_sets <- left_join(cat_sets, next_id, by = "ID")
Expand Down Expand Up @@ -142,40 +142,40 @@ validate_flowpath <- function(flowpath, outlets, post_mortem_file) {

get_lps <- function(flowpath) {
flowpath <- drop_geometry(flowpath) %>%
select(.data$ID, .data$LevelPathID, .data$Hydroseq) %>%
group_by(.data$LevelPathID)
select(ID, LevelPathID, Hydroseq) %>%
group_by(LevelPathID)

headwaters <- filter(flowpath, .data$Hydroseq == max(.data$Hydroseq)) %>%
headwaters <- filter(flowpath, Hydroseq == max(Hydroseq)) %>%
ungroup()

outlets <- filter(flowpath, .data$Hydroseq == min(.data$Hydroseq)) %>%
outlets <- filter(flowpath, Hydroseq == min(Hydroseq)) %>%
ungroup()

left_join(ungroup(flowpath), select(headwaters, head_ID = .data$ID, .data$LevelPathID),
left_join(ungroup(flowpath), select(headwaters, head_ID = ID, LevelPathID),
by = "LevelPathID") %>%
left_join(select(outlets, tail_ID = .data$ID, .data$LevelPathID),
left_join(select(outlets, tail_ID = ID, LevelPathID),
by = "LevelPathID")
}

# Get the levelpath outlet IDs for each of the input outlets.
get_outlets <- function(outlets, lps) {
distinct(left_join(outlets,
select(lps, .data$ID, .data$LevelPathID, .data$tail_ID),
select(lps, ID, LevelPathID, tail_ID),
by = "ID"))
}

# Adds everything that contributes the same recieving catchment as a given tail id.
fix_nexus <- function(flowpath, tail_id, da_thresh = NA, only_larger = FALSE) {
tail <- filter(flowpath, .data$ID == tail_id)
tail <- filter(flowpath, ID == tail_id)

add <- filter(flowpath, .data$toID == tail$toID)
add <- filter(flowpath, toID == tail$toID)

if (only_larger) {
add <- filter(add, .data$LevelPathID <= tail$LevelPathID)
add <- filter(add, LevelPathID <= tail$LevelPathID)
}

if (!is.na(da_thresh)) {
add <- filter(add, .data$TotDASqKM > da_thresh)
add <- filter(add, TotDASqKM > da_thresh)
}
# Can add functionality here to filter which Add IDs to include.
add_ids <- add$ID
Expand Down Expand Up @@ -229,16 +229,16 @@ make_outlets_valid <- function(outlets, flowpath,
outlets <- bind_rows(
outlets,
apply_fix_nexus(filter(get_outlets(outlets, lps),
.data$ID == .data$tail_ID)$ID,
ID == tail_ID)$ID,
flowpath, da_thresh, only_larger)
)

# deduplicate outlets.
outlets <- distinct(outlets) %>%
group_by(.data$ID) %>%
group_by(ID) %>%
filter(!(n() > 1 & # this removes outlets that duplicate terminals.
# they can be added in the above outlets check.
.data$type == "outlet")) %>%
type == "outlet")) %>%
ungroup()

otl <- get_outlets(outlets, lps)
Expand Down Expand Up @@ -269,45 +269,45 @@ make_outlets_valid <- function(outlets, flowpath,
# Need to check that a "next down tributary" in the outlet set has a break along the
# main stem that each outlet contributes to.
otl <- left_join(otl, select(drop_geometry(flowpath),
.data$ID, .data$toID), by = "ID") %>%
left_join(select(lps, .data$ID,
toID_hydroseq = .data$Hydroseq,
toID_tail_ID = .data$tail_ID,
toID_LevelpathID = .data$LevelPathID),
ID, toID), by = "ID") %>%
left_join(select(lps, ID,
toID_hydroseq = Hydroseq,
toID_tail_ID = tail_ID,
toID_LevelpathID = LevelPathID),
by = c("toID" = "ID"))

# this grabs the most downstream if duplicates were generated.
otl <- group_by(otl, .data$ID) %>%
filter(.data$toID_hydroseq == min(.data$toID_hydroseq)) %>%
otl <- group_by(otl, ID) %>%
filter(toID_hydroseq == min(toID_hydroseq)) %>%
ungroup() %>%
# This eliminates groups where only one thing goes to a given tail_ID.
group_by(.data$toID_tail_ID) %>%
group_by(toID_tail_ID) %>%
filter(n() > 1) %>%
ungroup()

# This grabs all the inflows to the nexus that each of these outlets is at.
otl <- left_join(otl, select(drop_geometry(flowpath),
toID_fromID = .data$ID, .data$toID),
toID_fromID = ID, toID),
by = "toID") %>%
mutate(type = "add_outlet")

if (!is.na(da_thresh)) {
otl <- left_join(otl, select(drop_geometry(flowpath),
.data$ID, toID_fromID_TotDASqKM = .data$TotDASqKM),
ID, toID_fromID_TotDASqKM = TotDASqKM),
by = c("toID_fromID" = "ID")) %>%
filter(.data$toID_fromID_TotDASqKM > da_thresh)
filter(toID_fromID_TotDASqKM > da_thresh)
}

if (only_larger) {
otl <- left_join(otl, select(drop_geometry(flowpath),
.data$ID, toID_fromID_lp = .data$LevelPathID),
ID, toID_fromID_lp = LevelPathID),
by = c("toID_fromID" = "ID")) %>%
filter(.data$toID_fromID_lp <= .data$toID_LevelpathID)
filter(toID_fromID_lp <= toID_LevelpathID)
}

otl <- otl %>%
select(ID = .data$toID_fromID, .data$type,
LevelPathID = .data$toID_LevelpathID, .data$tail_ID) %>%
select(ID = toID_fromID, type,
LevelPathID = toID_LevelpathID, tail_ID) %>%
distinct()

# Need to verify that all ID == tail_ID instances are connected.
Expand All @@ -327,8 +327,8 @@ make_outlets_valid <- function(outlets, flowpath,

connected <- FALSE
while (!connected) {
toid <- filter(flowpath, .data$ID == otl[["ID"]][check])$toID
toid_tail_id <- filter(lps, .data$ID == toid)[["tail_ID"]]
toid <- filter(flowpath, ID == otl[["ID"]][check])$toID
toid_tail_id <- filter(lps, ID == toid)[["tail_ID"]]

if (!all(toid_tail_id %in% otl$tail_ID)) {
outlets <- rbind(outlets,
Expand Down Expand Up @@ -359,32 +359,32 @@ sort_outlets <- function(outlets, flowpath) {
# working with a dendritic assumption so each nexus has one and only
# one downstream catchment. Otherwise, we would have to treat things a
# bit differently.
mutate(nexID_stem = .data$ID) %>%
mutate(nexID_stem = ID) %>%
# the nexID_terminal is the to nexus of the terminal flowpath.
# We need to treat this differently because an outlet of type
# "terminal" must be downstream of all other outlets. Without a
# unique nexus ID like this, we can not be gaurunteed a sort
# order that will work correctly.
left_join(select(drop_geometry(flowpath), .data$ID,
nexID_terminal = .data$toNexID), by = "ID") %>%
left_join(select(drop_geometry(flowpath), ID,
nexID_terminal = toNexID), by = "ID") %>%
# Now we can set the nexID (which is used to sort the outlets) appropriately.
mutate(
nexID = case_when(
type == "outlet" ~ .data$nexID_stem,
type == "terminal" ~ .data$nexID_terminal
type == "outlet" ~ nexID_stem,
type == "terminal" ~ nexID_terminal
)
) %>%
select(.data$ID, .data$type, .data$nexID) %>%
select(ID, type, nexID) %>%
distinct()

fp_sort <- nhdplusTools::get_sorted(drop_geometry(select(flowpath, .data$ID, .data$toNexID)))$ID
fp_sort <- nhdplusTools::get_sorted(drop_geometry(select(flowpath, ID, toNexID)))$ID

fp_sort <- c(fp_sort, flowpath$toNexID[flowpath$toNexID < 0])

fp_sort <- fp_sort[fp_sort %in% outlets$nexID]

left_join(data.frame(nexID = as.numeric(fp_sort)), outlets, by = "nexID") %>%
select(-.data$nexID)
select(-nexID)
}

prep_flowpath <- function(flowpath) {
Expand All @@ -401,7 +401,7 @@ prep_flowpath <- function(flowpath) {
}

prep_outlets <- function(outlets, flowpath) {
outlets <- left_join(outlets, select(flowpath, .data$ID, .data$id), by = "ID")
outlets <- left_join(outlets, select(flowpath, ID, id), by = "ID")
outlets$set <- 1:nrow(outlets)
outlets
}
Expand Down Expand Up @@ -530,7 +530,7 @@ get_catchment_sets <- function(flowpath, outlets) {

head_paths <- filter(flowpath, LevelPathID %in% head_outlets$LevelPathID) %>%
left_join(head_outlets, by = "LevelPathID") %>%
filter(Hydroseq >= .data$head_out_Hydroseq) %>%
filter(Hydroseq >= head_out_Hydroseq) %>%
select(ID, set) %>%
group_by(set) %>%
summarise(ID = list(ID))
Expand Down Expand Up @@ -598,27 +598,27 @@ get_minimal_network <- function(flowpath, outlets) {

minimal <- aggregate_network_to_outlets(
flowpath,
outlets = dplyr::filter(outlets, .data$ID %in% flowpath$ID),
outlets = dplyr::filter(outlets, ID %in% flowpath$ID),
da_thresh = NA, only_larger = TRUE)

min_net <- unnest_flines(select(drop_geometry(minimal$fline_sets), -LevelPathID)) %>%
left_join(select(flowpath, .data$ID, .data$LENGTHKM,
.data$AreaSqKM, .data$LevelPathID),
left_join(select(flowpath, ID, LENGTHKM,
AreaSqKM, LevelPathID),
by = c("set" = "ID")) %>%
group_by(ID) %>%
summarise(toID = .data$toID[1],
lengthkm = sum(.data$LENGTHKM),
areasqkm = sum(.data$AreaSqKM),
outlet_levelpath = min(.data$LevelPathID)) %>%
mutate(toID = ifelse(is.na(.data$toID), 0, .data$toID)) %>%
rename(comid = .data$ID,
tocomid = .data$toID,
nameID = .data$outlet_levelpath)%>%
summarise(toID = toID[1],
lengthkm = sum(LENGTHKM),
areasqkm = sum(AreaSqKM),
outlet_levelpath = min(LevelPathID)) %>%
mutate(toID = ifelse(is.na(toID), 0, toID)) %>%
rename(comid = ID,
tocomid = toID,
nameID = outlet_levelpath)%>%
add_plus_network_attributes() %>%
rename(ID = .data$comid, toID = .data$tocomid,
outlet_nhdpv2_levelpath = .data$nameID,
arbolate_sum = .data$weight) %>%
left_join(select(minimal$fline_sets, .data$ID, .data$set),
rename(ID = comid, toID = tocomid,
outlet_nhdpv2_levelpath = nameID,
arbolate_sum = weight) %>%
left_join(select(minimal$fline_sets, ID, set),
by = "ID")

if(inherits(minimal$fline_sets, "sf")) {
Expand Down
Loading

0 comments on commit 1630501

Please sign in to comment.