# define output files and see if they have already been generated
out_file_mm <- file.path(samp_dir, "mm", "samp_pts_new.shp")
out_file_nk <- file.path(samp_dir, "nk", "samp_pts_new.shp")
if (!all(file.exists(out_file_mm), file.exists(out_file_nk))){
# loop over plots and find most similar pixels -- mm
n_mm <- nrow(unburned_mm_pc)
for (priority in c("B", "C")){
for (i in 1:n_mm) {
plot_id <- unburned_mm_pc$plot_id[i]
dist_rast <- dist_stack_mm[[paste0("dist_", plot_id)]]
min_dist <- global(dist_rast, "min", na.rm = T)[1,1]
min_rast <- ifel(abs(dist_rast - min_dist) < 1e-9, 1, NA)
min_pt <- as.points(min_rast)[1,]
min_pt$dist <- extract(dist_rast, min_pt, ID = F)[1,1]
min_pt$priority <- priority
min_pt$plot_id <- plot_id
min_pt <- min_pt[,c("plot_id", "priority", "dist")]
min_pt$type <- "new"
min_pts_mm <- rbind(min_pts_mm, min_pt)
eraser <- buffer(min_pt, 100)
dist_stack_mm <- mask(dist_stack_mm, eraser, inverse = T)
}
}
# loop over plots and find most similar pixels -- nk
n_nk <- nrow(unburned_nk_pc)
for (priority in c("B", "C")){
for (i in 1:n_nk) {
plot_id <- unburned_nk_pc$plot_id[i]
dist_rast <- dist_stack_nk[[paste0("dist_", plot_id)]]
min_dist <- global(dist_rast, "min", na.rm = T)[1,1]
min_rast <- ifel(abs(dist_rast - min_dist) < 1e-9, 1, NA)
min_pt <- as.points(min_rast)[1,]
min_pt$dist <- extract(dist_rast, min_pt, ID = F)[1,1]
min_pt$priority <- priority
min_pt$plot_id <- plot_id
min_pt <- min_pt[,c("plot_id", "priority", "dist")]
min_pt$type <- "new"
min_pts_nk <- rbind(min_pts_nk, min_pt)
eraser <- buffer(min_pt, 100)
dist_stack_nk <- mask(dist_stack_nk, eraser, inverse = T)
}
}
# write to files
writeVector(min_pts_mm, out_file_mm, overwrite = T)
writeVector(min_pts_nk, out_file_nk, overwrite = T)
} else {
# read the existing files in
min_pts_mm <- vect(out_file_mm)
min_pts_nk <- vect(out_file_nk)
}
# plot the data out -- mm
par(mfrow = c(1,2))
plot(als_mm, main = "Monroe Mountain")
plot(unburned_mm, col = 1, add = T)
plot(burned_mm, col = 2, add = T)
cols <- NA
cols[min_pts_mm$priority == "A"] <- 3
cols[min_pts_mm$priority == "B"] <- 4
cols[min_pts_mm$priority == "C"] <- 5
plot(min_pts_mm, col = cols, add = T)
add_legend(
"bottomright",
legend = c(
"2023 unburned",
"2023 burned",
"2026 sample A",
"2026 sample B",
"2026 sample C"
),
pch = 16,
col = c(1,2,3,4,5),
cex = 2/3
)
# plot the data out -- nk
plot(als_nk, main = "North Kaibab")
plot(unburned_nk, col = 1, add = T)
plot(burned_nk, col = 2, add = T)
cols <- NA
cols[min_pts_nk$priority == "A"] <- 3
cols[min_pts_nk$priority == "B"] <- 4
cols[min_pts_nk$priority == "C"] <- 5
plot(min_pts_nk, col = cols, add = T)
add_legend(
"bottomleft",
legend = c(
"2023 unburned",
"2023 burned",
"2026 sample A",
"2026 sample B",
"2026 sample C"
),
pch = 16,
col = c(1,2,3,4,5),
cex = 2/3
)