diff --git a/workspace/07-2_NonDend.Rmd b/workspace/07-2_NonDend.Rmd
index fe535cc41ff00ff2b8dd3fa52a202305e5e6f1c4..04288603a3450654931d4468d9ef5b216a63a71b 100644
--- a/workspace/07-2_NonDend.Rmd
+++ b/workspace/07-2_NonDend.Rmd
@@ -38,6 +38,8 @@ if (vpu == "20"){
full_nhd <- readRDS(file.path(data_paths$nhdplus_dir, "nhdplus_flowline.rds"))
}
+elev <- data_paths$elev_cm[grepl(paste0("Ned", substr(vpu, 1, 2)), data_paths$elev_cm, ignore.case = TRUE)]
+
# HUC extraction for specific NHD+ vpus
if(vpu == "02"){
grep_exp <-"^02|^04"
@@ -45,6 +47,7 @@ if(vpu == "02"){
grep_exp <- "^03|^08"
} else {
grep_exp <- paste0("^", substr(vpu, start = 1, stop = 2))
+ elev <- append(elev, list(rpu_03g = data_paths$elev_cm$rpu_03g))
}
cat_rpu_table <- readRDS(data_paths$fullcats_table)
@@ -57,8 +60,6 @@ full_nhd <- full_nhd %>%
vpu_WBD <- readRDS(file.path(data_paths$nhdplus_dir, "HUC12.rds")) %>%
filter(grepl(paste0("^", grep_exp, ".*"), .data$HUC_12))
-elev <- data_paths$elev_cm[grepl(paste0("Ned", vpu), data_paths$elev_cm, ignore.case = TRUE)]
-
nhd <- st_transform(read_sf(gf_gpkg, nhd_flowline), crs)
cats <- st_transform(read_sf(gf_gpkg, nhd_catchment), crs)
divides <- st_transform(read_sf(gf_gpkg, divide_layer), crs)
@@ -98,7 +99,7 @@ if(needs_layer(ND_gpkg, xwalk_layer)){
# Bring over divides/HUC12 intersection information into divides layer
xwalk_nhd_wbd <- st_drop_geometry(nhd_wbd_int$cats_HUC12) %>%
- select(-c(ACRES, SOURCEFC, Shape_Length, Shape_Area, HU_12_MOD))
+ select(-c(ACRES, HU_12_MOD))
rm(nhd_wbd_int)
@@ -350,7 +351,7 @@ if(needs_layer(ND_gpkg, "missing_cats")){
st_make_valid()
# Prob HRU - filter(all_hrus, POI_ID == 140402000209)
- all_hrus <- divides_lu %>%
+ all_hrus <- filter(divides_lu, !is.na(POI_ID)) %>%
group_by(POI_ID) %>%
summarize(do_union = T) %>%
sfheaders::sf_remove_holes(.) %>%
diff --git a/workspace/R/non_dend.R b/workspace/R/non_dend.R
index fc3a2565de1e0f402314ce1c48c58a23dc13c089..2d2943ea2e9a1cd62192ec01e8fff2736bceb81a 100644
--- a/workspace/R/non_dend.R
+++ b/workspace/R/non_dend.R
@@ -82,7 +82,7 @@ assign_HUC12 <- function(divides_poi, HUC12_xwalk, nhd, CAC_num){
ungroup() %>%
group_by(FEATUREID, HUC_12_int) %>%
summarize(n_NHD = n(),
- NHD_area = mean(AreaSqKM),
+ NHD_area = mean(AREASQKM),
HUC_12_area = sum(AreaHUC12),
intarea = sum(intArea),
CAC = intarea / ((NHD_area - intarea) + (HUC_12_area - intarea) + intarea)) %>%
@@ -118,7 +118,7 @@ assign_HUC12 <- function(divides_poi, HUC12_xwalk, nhd, CAC_num){
ungroup() %>%
group_by(HUC_12_int) %>%
summarize(n_HUC12 = n(),
- NHD_area = sum(AreaSqKM),
+ NHD_area = sum(AREASQKM),
HUC_12_area = mean(AreaHUC12),
intarea = sum(intArea),
CAC = intarea / ((NHD_area - intarea) + (HUC_12_area - intarea) + intarea)) %>%
@@ -233,7 +233,7 @@ assign_HUC10 <- function(divides, HUC12_xwalk, nhd, CAC_num){
ungroup() %>%
group_by(HUC_10) %>%
mutate(n_NHD = n(),
- NHD_area = mean(AreaSqKM),
+ NHD_area = mean(AREASQKM),
HUC_10_area = mean(AreaHUC10),
intarea = sum(intArea),
CAC = intarea / ((NHD_area - intarea) + (HUC_10_area - intarea) + intarea)) %>%
@@ -849,35 +849,46 @@ miss_term_assign <- function(term_outlets, divides_poi, nhd, elev){
st_intersection(divides_lu_poi) %>%
filter(!member_COMID == incomid)
- # convert cat buf to format for extract input
- buf <- st_as_sf(cats_buff)
-
- #**************************************************************
- # Then do a zonal statitics type operation of the DEM within each
- # buffered polygon to get the min cell value, this is a hard coded path to where the dems reside
-
- elev <- as.character(elev)
- elev <- elev[grepl(unique(missing_cats$rpu), elev)]
-
- dem <- terra::rast(elev)
- ex <- terra::extract(dem, terra::vect(buf), min)
- buf_df <- as.data.frame(ex)[,2, drop = FALSE]
-
- # Attribute with values
- buf_df$outlet_COMID <- incomid
- buf_df$neighbor_COMID <- as.character(cats_buff$member_COMID)
- buf_df$POI_ID = as.character(cats_buff$POI_ID)
-
- colnames(buf_df) <- c("Elev","outlet_COMID", "neighbor_COMID", "POI_ID")
- buf_df2 <- buf_df[which.min(buf_df$Elev), ]
-
- return(buf_df2)
+ if(nrow(cats_buff) == 0){
+ miss_cat <- dplyr::select(st_drop_geometry(missing_cats),
+ outlet_COMID = member_COMID) %>%
+ dplyr::mutate(Elev = NA, neighbor_COMID = NA, POI_ID = NA) %>%
+ dplyr::select(Elev, outlet_COMID, neighbor_COMID, POI_ID)
+ return(miss_cat)
+ } else {
+ # convert cat buf to format for extract input
+ buf <- st_as_sf(cats_buff)
+
+ #**************************************************************
+ # Then do a zonal statitics type operation of the DEM within each
+ # buffered polygon to get the min cell value, this is a hard coded path to where the dems reside
+
+ elev <- as.character(elev)
+ elev <- elev[grepl(unique(missing_cats$rpu), elev)]
+
+ dem <- terra::rast(elev)
+ ex <- terra::extract(dem, terra::vect(buf), min)
+ buf_df <- as.data.frame(ex)[,2, drop = FALSE]
+
+ # Attribute with values
+ buf_df$outlet_COMID <- incomid
+ buf_df$neighbor_COMID <- as.character(cats_buff$member_COMID)
+ buf_df$POI_ID = as.character(cats_buff$POI_ID)
+
+ colnames(buf_df) <- c("Elev","outlet_COMID", "neighbor_COMID", "POI_ID")
+ buf_df2 <- buf_df[which.min(buf_df$Elev), ]
+
+ return(buf_df2)
+ }
}
- out_df <- do.call(rbind,
- pbapply::pblapply(unique(term_outlets$outlet_COMID),
+ out_df <- do.call(rbind,
+ pbapply::pblapply(unique(term_outlets$outlet_COMID),
assign_func, divides_poi, nhd, elev, cl = NULL))
-
+
+ # out_df <- lapply(unique(term_outlets$outlet_COMID),
+ # assign_func, divides_poi, nhd, elev)
+
out_df <- out_df %>%
left_join(dplyr::select(st_drop_geometry(divides_poi), HUC12_neighbor = HUC_12_int, member_COMID),
by = c("neighbor_COMID" = "member_COMID")) %>%