Skip to content

Instantly share code, notes, and snippets.

@mdsumner
Created November 8, 2024 01:17
Show Gist options
  • Save mdsumner/95f3edf1909e10e82d69a9bdcee51ca2 to your computer and use it in GitHub Desktop.
Save mdsumner/95f3edf1909e10e82d69a9bdcee51ca2 to your computer and use it in GitHub Desktop.
library(terra)

options(parallelly.fork.enable = TRUE, future.rng.onMisuse = "ignore")
library(furrr)



topo <- project(rast(sds::gebco()), rast(), by_util = TRUE)
plot(topo)
cells <- sample(which(values(topo)[,1] > 0), 1000)

p2ex <- function(x) {
  x <- floor(x)
  cbind(x[,1], x[,1] + 1, x[,2], x[,2] + 1)
  #rep(floor(x), 2) + c(0, 0, 1, 1)
}
extents <- xyFromCell(topo, cells) |> p2ex() #|> vaster::plot_extent(add = T)

l <- vector("list")

get_stac <- function(ex) {
  h <- try(tacky:::hrefs(sds::stacit(ex, "2024-11")), silent = TRUE)
  if (!inherits(h, "try-error") && nrow(h) > 0) {
    h$cell <- cells[i]
    return(h)
  }
  return(NULL)
}
plan(multicore)
l <- future_map_dfr(split(t(extents), rep(1:nrow(extents), each = 4)), 
               get_stac)



bb2llb <- function(x) {
  exy <- do.call(rbind, x)
  wch <- exy[,3] < exy[,1]
  if (any(wch)) {
    exy[wch, 3] <-  exy[wch, 3] + 360
  }
  wk::rct(exy[,1], exy[,2], exy[,3], exy[,4])
}
x <- bb2llb(l$bbox)
#pl```ot(x)


library(wk)
rc <- as_rct(grd(rct(-180, -90, 180, 90), dx = 1, dy = 1, type = "polygons"))


tree <- geos_basic_strtree(rc)
idx <- geos_basic_strtree_query(tree, x)

ic <- sample(1:nrow(idx), 1)
plot(c(x[idx[ic, 1]], rc[idx[ic, 2]]), col = c("grey", "transparent"))
abline(v = 180, lty = 2)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment