diff --git a/workspace/00_get_data.Rmd b/workspace/00_get_data.Rmd index 91bfcd70011d9af2fb72eff91e33f127d093503d..cdeef4f85971abec57863c3b5022835c7196a3dd 100644 --- a/workspace/00_get_data.Rmd +++ b/workspace/00_get_data.Rmd @@ -300,8 +300,6 @@ waterbodies_path <- file.path(nhdplus_dir, "nhdplus_waterbodies.rds") if(!file.exists(waterbodies_path)) { message("formatting NHDPlus watebodies...") - nhdplus_path(nhdplus_gdb) - # Read the feature class wbSF <- read_sf(nhdplus_gdb, "NHDWaterbody") %>% st_transform(5070) diff --git a/workspace/01_Gage_Selection.Rmd b/workspace/01_Gage_Selection.Rmd index 8a9a6ecc6c84c10d8129b9b466c295ab4fb0e7f9..30546dfa191eb278ad7b1ed6e6727da40e4e99b2 100644 --- a/workspace/01_Gage_Selection.Rmd +++ b/workspace/01_Gage_Selection.Rmd @@ -36,7 +36,7 @@ source("R/gage_mod.R") # Load user-defined data paths data_paths <- jsonlite::read_json(file.path("cache", "data_paths.json")) # Load NHDPlus Data if precprocessed to RDS format -nhdplus_path(data_paths$nhdplus_gdb) + staged_nhd <- stage_national_data(nhdplus_data = data_paths$nhdplus_gdb, output_path = data_paths$nhdplus_dir) diff --git a/workspace/R/utils.R b/workspace/R/utils.R index cf022aa6366d61897b2864726a5a639f2b45d10a..0c8535912efc3143b8d2b42b84f695c7a056f41f 100644 --- a/workspace/R/utils.R +++ b/workspace/R/utils.R @@ -10,117 +10,6 @@ if(!(require(hyfabric, quietly = TRUE) && packageVersion("hyfabric") == hyfabric repos = NULL, type = "source") } -fix_headwaters <- function(nhd_rds, out_gpkg, new_atts = NULL, - nhdpath = nhdplusTools::nhdplus_path(), - par = 6) { - - if(!file.exists(out_gpkg)) { - nhd <- readRDS(nhd_rds) - - # Need to replace headwater flowlines with the trimmed burn line event. - # Some headwater flowlines are not fully within the catchment they go with. - ble <- sf::read_sf(nhdpath, "BurnLineEvent") - - ble <- dplyr::left_join(dplyr::select(sf::st_drop_geometry(nhd), COMID), - ble, by = "COMID") - - ble <- sf::st_zm(sf::st_as_sf(ble)) - nhd <- sf::st_zm(nhd) - - sf::st_geometry(nhd)[!sf::st_is_empty(sf::st_geometry(ble)) & (nhd$StartFlag == 1)] <- - sf::st_geometry(ble)[!sf::st_is_empty(sf::st_geometry(ble)) & (nhd$StartFlag == 1)] - - if(!is.null(new_atts)) { - new_atts <- data.table::fread(new_atts, - integer64 = "character", - showProgress = FALSE) - - if("weight" %in% names(new_atts)) new_atts <- rename(new_atts, - ArbolateSu = weight) - - nhd <- select(nhd, COMID, FromNode, ToNode, - StartFlag, StreamCalc, Divergence, DnMinorHyd) %>% - left_join(new_atts, by = c("COMID" = "comid")) %>% - nhdplusTools::align_nhdplus_names() - } - - nhd <- sf::st_sf(nhd) - - m_to_km <- (1/1000) - - nhd$LENGTHKM <- as.numeric(sf::st_length(sf::st_transform(nhd, 5070))) * m_to_km - - custom_net <- select(sf::st_drop_geometry(nhd), - COMID, FromNode, ToNode, Divergence) - - custom_net <- nhdplusTools::get_tocomid(custom_net, remove_coastal = FALSE) - - custom_net <- select(custom_net, comid, override_tocomid = tocomid) - - nhd <- left_join(nhd, custom_net, by = c("COMID" = "comid")) - - nhd <- mutate(nhd, override_tocomid = ifelse(toCOMID == 0, override_tocomid, toCOMID)) - - # is a headwater and for sure flows to something. - check <- !nhd$COMID %in% nhd$override_tocomid & - !(nhd$override_tocomid == 0 | is.na(nhd$override_tocomid) | - !nhd$override_tocomid %in% nhd$COMID) - check_direction <- dplyr::filter(nhd, check) - - if(!all(check_direction$override_tocomid[check_direction$override_tocomid != 0] %in% nhd$COMID)) - stop("this won't work") - - matcher <- match(check_direction$override_tocomid, nhd$COMID) - - matcher <- nhd[matcher, ] - - fn_list <- Map(list, - split(check_direction, seq_len(nrow(check_direction))), - split(matcher, seq_len(nrow(matcher)))) - - cl <- parallel::makeCluster(16) - - # check and fix these. - new_geom <- pbapply::pblapply(cl = cl, - X = fn_list, - FUN = function(x) { - nhdplusTools::fix_flowdir(x[[1]]$COMID, - fn_list = list(flowline = x[[1]], - network = x[[2]], - check_end = "end")) - }) - - parallel::stopCluster(cl) - - # One is empty. Deal with it and remove from the replacement. - error_index <- sapply(new_geom, inherits, what = "try-error") - errors <- filter(nhd, COMID %in% check_direction$COMID[error_index]) - check[which(nhd$COMID %in% errors$COMID)] <- FALSE - - if(!all(sapply(sf::st_geometry(errors), sf::st_is_empty))) { - stop("Errors other than empty geometry from fix flowdir") - } - - new_geom <- new_geom[!error_index] - new_geom <- do.call(c, new_geom) - st_geometry(nhd)[check] <- new_geom - - nhd <- dplyr::filter(nhd, COMID %in% new_atts$comid) - - nhd <- select(nhd, -override_tocomid) - - # remove remnant attributes - nhd <- select(nhd, -FromNode, -ToNode, -StartFlag -StreamCalc, -Divergence, -DnMinorHyd) - - warning("naive duplication!") - nhd <- nhd[!duplicated(nhd$COMID), ] - - sf::write_sf(nhd, out_gpkg, "reference_flowlines") - - } - return(out_gpkg) -} - ## Is available in sbtools now, but will leave till on cran. #' Retrieves files from SB in facet file_structure #' @param sbitem (character) ScienceBase item ID