Skip to content

Commit

Permalink
ready for source
Browse files Browse the repository at this point in the history
  • Loading branch information
@mikejohnson51
mikejohnson51 committed Oct 18, 2023
1 parent 150e5e5 commit eed498f
Showing 1 changed file with 127 additions and 86 deletions.
213 changes: 127 additions & 86 deletions R/subset_network.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,27 +34,43 @@ create_style_row = function (gpkg_path,
)
}

append_style = function (gpkg_path, qml_dir = system.file("qml", package = "hydrofabric"), layer_names) {
append_style = function (gpkg_path,
qml_dir = system.file("qml", package = "hydrofabric"),
layer_names) {
f = list.files(qml_dir, full.names = TRUE)

f = list.files(qml_dir, full.names = TRUE)

good_layers = gsub(".qml", "", basename(f))

layer_names = layer_names[layer_names %in% good_layers]

files = grep(paste(layer_names, collapse = "|"), f, value = TRUE)

styles <- sapply(files, read_qml)
style_names <- sapply(layer_names, paste0, "__hydrofabric_style")
style_rows <- do.call(rbind, mapply(create_style_row, layer_names,
style_names, styles, MoreArgs = list(gpkg_path = gpkg_path),
SIMPLIFY = FALSE))
if ("layer_styles" %in% sf::st_layers(gpkg_path)$name) {
try(st_delete(gpkg_path, "layer_styles"), silent = TRUE)
}
st_write(style_rows, gpkg_path, layer = "layer_styles", append = FALSE,
quiet = TRUE)
return(gpkg_path)
good_layers = gsub(".qml", "", basename(f))

layer_names = layer_names[layer_names %in% good_layers]

files = grep(paste(layer_names, collapse = "|"), f, value = TRUE)

styles <- sapply(files, read_qml)
style_names <-
sapply(layer_names, paste0, "__hydrofabric_style")
style_rows <-
do.call(
rbind,
mapply(
create_style_row,
layer_names,
style_names,
styles,
MoreArgs = list(gpkg_path = gpkg_path),
SIMPLIFY = FALSE
)
)
if ("layer_styles" %in% sf::st_layers(gpkg_path)$name) {
try(sf::st_delete(gpkg_path, "layer_styles"), silent = TRUE)
}
sf::st_write(
style_rows,
gpkg_path,
layer = "layer_styles",
append = FALSE,
quiet = TRUE
)
return(gpkg_path)
}

#' Access Hydrofabric Network
Expand All @@ -68,21 +84,24 @@ append_style = function (gpkg_path, qml_dir = system.file("qml", package = "hydr
get_fabric = function(VPU,
base_s3 = 's3://lynker-spatial/pre-release/',
cache_dir = NULL,
cache_overwrite = FALSE){

cache_overwrite = FALSE) {
Key <- NULL

xx = get_bucket_df(bucket = dirname(base_s3), prefix = basename(base_s3), region = 'us-west-2') %>%
filter(grepl(basename(base_s3), Key) &
grepl(paste0(VPU, ".gpkg$"), Key)) %>%
filter(!grepl("[.]_", Key)) %>%
filter(!grepl("/", dirname(Key)))
xx = aws.s3::get_bucket_df(
bucket = dirname(base_s3),
prefix = basename(base_s3),
region = 'us-west-2'
) |>
dplyr::filter(grepl(basename(base_s3), Key) &
grepl(paste0(VPU, ".gpkg$"), Key)) |>
dplyr::filter(!grepl("[.]_", Key)) |>
dplyr::filter(!grepl("/", dirname(Key)))

if (!is.null(cache_dir)) {
dir.create(cache_dir,
recursive = TRUE,
showWarnings = FALSE)
gpkg = glue("{cache_dir}/{basename(xx$Key)}")
gpkg = glue::glue("{cache_dir}/{basename(xx$Key)}")
if (cache_overwrite) {
unlink(gpkg)
}
Expand All @@ -93,10 +112,12 @@ get_fabric = function(VPU,
}

if (!file.exists(gpkg)) {
save_object(bucket = xx$Bucket,
object = xx$Key,
file = gpkg,
region = 'us-west-2')
aws.s3::save_object(
bucket = xx$Bucket,
object = xx$Key,
file = gpkg,
region = 'us-west-2'
)
}

return(gpkg)
Expand All @@ -105,7 +126,7 @@ get_fabric = function(VPU,

#' Subset Hydrofabric Network
#' @param id hydrofabric id (relevant only to nextgen fabrics)
#' @param comid NHDPlusV2 COMID
#' @param comid NHDPlusV2 COMID
#' @param hl_uri hydrolocation URI (relevant only to nextgen fabrics)
#' @param nldi_feature list with names 'featureSource' and 'featureID' where 'featureSource' is derived from the "source" column of the response of dataRetrieval::get_nldi_sources() and the 'featureID' is a known identifier from the specified 'featureSource'.
#' @param loc Location given as vector of XY in CRS 4326 (long, lat)
Expand Down Expand Up @@ -138,43 +159,40 @@ subset_network = function(id = NULL,
cache_dir = NULL,
qml_dir = system.file("qml", package = "hydrofabric"),
cache_overwrite = FALSE) {
Key <-
hf_hydroseq <-
hf_id <- hydroseq <- member_COMID <- toid <- vpu <- NULL

Key <- hf_hydroseq <- hf_id <- hydroseq <- member_COMID <- toid <- vpu <- NULL

net = tryCatch({
open_dataset(glue(base_s3, "conus_net.parquet")) %>%
select(id, toid, hf_id, hl_uri, hf_hydroseq, hydroseq, vpu) %>%
collect() %>%
distinct()
}, error = function(e) {
NULL
})
net = arrow::open_dataset(glue::glue(base_s3, "conus_net.parquet")) |>
dplyr::select(id, toid, hf_id, hl_uri, hf_hydroseq, hydroseq, vpu) |>
dplyr::collect() |>
dplyr::distinct()

if (!is.null(id) & !is.null(net)) {
comid = filter(net, id == !!id | toid == !!id) %>%
slice_max(hf_hydroseq) %>%
pull(hf_id)
comid = dplyr::filter(net, id == !!id | toid == !!id) |>
dplyr::slice_max(hf_hydroseq) |>
dplyr::pull(hf_id)
}

if (!is.null(nldi_feature)) {
comid = get_nldi_feature(nldi_feature)$comid
comid = nhdplusTools::get_nldi_feature(nldi_feature)$comid
}

if (!is.null(loc)) {
comid = discover_nhdplus_id(point = st_sfc(st_point(c(loc[1], loc[2])), crs = 4326))
comid = nhdplusTools::discover_nhdplus_id(point = sf::st_sfc(sf::st_point(c(loc[1], loc[2])), crs = 4326))
}

if (!is.null(hl_uri) & !is.null(net)) {
origin = filter(net, hl_uri == !!hl_uri) %>%
slice_max(hf_hydroseq) %>%
pull(toid) %>%
origin = dplyr::filter(net, hl_uri == !!hl_uri) |>
dplyr::slice_max(hf_hydroseq) |>
dplyr::pull(toid) |>
unique()
}

if (!is.null(comid) & !is.null(net)) {
origin = filter(net, hf_id == comid) %>%
slice_max(hf_hydroseq) %>%
pull(id) %>%
origin = dplyr::filter(net, hf_id == comid) |>
dplyr::slice_max(hf_hydroseq) |>
dplyr::pull(id) |>
unique()
} else if (is.null(net)) {
origin = comid
Expand All @@ -183,52 +201,75 @@ subset_network = function(id = NULL,
if (is.null(origin)) {
stop("origin not found")
}


if (is.null(net)) {
xx = suppressMessages({
get_nhdplus(comid = comid)
nhdplusTools::get_nhdplus(comid = comid)
})

v = nhdplusTools::vpu_boundaries

vpuid = v$VPUID[which(lengths(st_intersects(st_transform(v, st_crs(xx)), xx)) > 0)]
vpuid = v$VPUID[which(lengths(sf::st_intersects(sf::st_transform(
v, sf::st_crs(xx)
), xx)) > 0)]
} else {
vpuid = unique(pull(filter(net, id == origin |
toid == origin), vpu))
vpuid = unique(dplyr::pull(dplyr::filter(net, id == origin |
toid == origin), vpu))
}


gpkg = get_fabric(VPU = vpuid, base_s3 = base_s3, cache_dir = cache_dir, cache_overwrite = cache_overwrite)
lyrs = lyrs[lyrs %in% st_layers(gpkg)$name]

db <- dbConnect(SQLite(), gpkg)
on.exit(dbDisconnect(db))
gpkg = get_fabric(
VPU = vpuid,
base_s3 = base_s3,
cache_dir = cache_dir,
cache_overwrite = cache_overwrite
)
lyrs = lyrs[lyrs %in% sf::st_layers(gpkg)$name]

db <- DBI::dbConnect(RSQLite::SQLite(), gpkg)
on.exit(DBI::dbDisconnect(db))

if (!is.null(net)) {
sub_net = distinct(select(filter(net, vpu == vpuid), id, toid))
sub_net = dplyr::distinct(dplyr::select(dplyr::filter(net, vpu == vpuid), id, toid))
} else {
lookup <-
c(
id = "ID",
id = "COMID",
toid = "toID",
toid = "toCOMID"
)
dplyr::sub_net = dplyr::tbl(db, lyrs[grepl("flowline|flowpath", lyrs)]) |>
dplyr::select(dplyr::any_of(
c(
"id",
"toid",
'COMID',
'toCOMID',
"ID",
"toID",
"member_COMID"
)
)) |>
dplyr::collect() |>
dplyr::rename(dplyr::any_of(lookup))

lookup <- c(id = "ID", id = "COMID", toid = "toID", toid = "toCOMID")

sub_net = tbl(db, lyrs[grepl("flowline|flowpath", lyrs)]) %>%
select(any_of(c("id", "toid", 'COMID', 'toCOMID', "ID", "toID", "member_COMID"))) %>%
collect() %>%
rename(any_of(lookup))

if("member_COMID" %in% names(sub_net)){
origin = filter(sub_net, grepl(origin, member_COMID)) %>%
pull(id)
if ("member_COMID" %in% names(sub_net)) {
origin = dplyr::filter(sub_net, grepl(origin, member_COMID)) |>
dplyr::pull(id)

sub_net = select(sub_net, -member_COMID)
sub_net = dplyr::select(sub_net, -member_COMID)
}
}

message("Starting from: `", origin, "`")

tmap = suppressWarnings({ get_sorted(distinct(sub_net), outlets = origin) })
tmap = suppressWarnings({
get_sorted(distinct(sub_net), outlets = origin)
})

if(grepl("nex", tail(tmap$id,1))){
if (grepl("nex", tail(tmap$id, 1))) {
tmap = head(tmap, -1)
}

Expand All @@ -241,28 +282,28 @@ subset_network = function(id = NULL,

crs = st_layers(gpkg)$crs

t = tbl(db, lyrs[j]) %>%
filter(if_any(any_of(
t = dplyr::tbl(db, lyrs[j]) |>
dplyr::filter(dplyr::if_any(dplyr::any_of(
c('COMID', 'FEATUREID', 'divide_id', 'id', 'ds_id', "ID")
), ~ . %in% !!ids)) %>%
collect()
), ~ . %in% !!ids)) |>
dplyr::collect()

if (all(!any(is.na(as.character(crs[[j]]))), nrow(t) > 0)) {
if (any(c("geometry", "geom") %in% names(t))) {
t = st_as_sf(t, crs = crs[[j]])
t = sf::st_as_sf(t, crs = crs[[j]])
} else {
t = t
}
}

if (!is.null(outfile)) {
write_sf(t, outfile, lyrs[j])
sf::write_sf(t, outfile, lyrs[j])
} else {
hydrofabric[[lyrs[j]]] = t
}
}

if(is.null(cache_dir)){
if (is.null(cache_dir)) {
unlink(gpkg)
}

Expand Down

0 comments on commit eed498f

Please sign in to comment.