Build global 0.05° cell grid (v8 sampling unit)

Global ocean 0.05° raster cells in [-180,180] with env covariates + US study-area / Program-Area membership

Published

2026-07-13 14:50:44

The v8 spatial sampling unit is the global 0.05° raster cell in [-180,180] (rolled back from H3 hex). The grid is anchored to the shared 0.05°, 3600×7200, EPSG:4326 topology of Bio-Oracle v3 AND the AquaX SDM TIFs — so cell_id is the raster cell index (1:ncell), and any SDM raster on this grid (Bio-Oracle, AquaX, resampled AquaMaps) maps to cell_id by position with no projection.

This notebook builds the canonical cell table (ocean cells only) — cell_id, the 12 Bio-Oracle env covariates, area_km2, and boolean in_usa / in_pra membership — plus the global cell-id COG that titiler paints (and that ingest_aquamaps.qmd uses to build its bilinear weight table).

The global build (631 MB Bio-Oracle → 25.9M cells → ~17M-row cell table + 47 MB COG) is skipped when the cell-id COG already exists (resume); force it with REDO_CELL_GRID=1 (libs/vars.R). A bbox prototype run (CELL_BBOX=...) always rebuilds into its own CELL_DB.

1 Design

Figure 1: Bio-Oracle raster → global 0.05° cell grid: cell table + cell-id COG, content-addressed manifest
Code 
librarian::shelf(DBI, dplyr, duckdb, glue, here, jsonlite, sf, terra, quiet = TRUE)
source(here("libs/paths.R"))   # ver, sdm_db, bio_oracle_tif, cellid_tif, ply_usa_gpkg, ply_pra_gpkg
source(here("libs/vars.R"))    # redo_cell_grid
stopifnot(file.exists(bio_oracle_tif), file.exists(ply_usa_gpkg), file.exists(ply_pra_gpkg))

# optional bbox restriction for prototype runs: CELL_BBOX="xmin,xmax,ymin,ymax"
cell_bbox <- Sys.getenv("CELL_BBOX", "")
cell_db   <- Sys.getenv("CELL_DB", sdm_db)
dir.create(dirname(cell_db), recursive = TRUE, showWarnings = FALSE)
dir.create(here("data/manifests"), recursive = TRUE, showWarnings = FALSE)

# skip the expensive global build if the cell-id COG already exists (resume),
# unless forced; a bbox prototype always rebuilds into its own CELL_DB
need_build <- redo_cell_grid || nzchar(cell_bbox) || !file.exists(cellid_tif)

2 Read Bio-Oracle, assign global cell_id, membership + area

Code 
if (need_build) {
  r <- rast(bio_oracle_tif)   # global [-180,180], 3600 x 7200, 12 env layers
  r$cell_id  <- 1:ncell(r)    # global raster-position id (shared with AquaX/AquaMaps grids)
  r$area_km2 <- cellSize(r$depth_mean, unit = "km")

  # membership: rasterize the study-area + Program-Area polygons onto the same grid
  r$in_usa <- rasterize(vect(st_read(ply_usa_gpkg, quiet = TRUE)), r$depth_mean, field = 1, background = 0)
  r$in_pra <- rasterize(vect(st_read(ply_pra_gpkg, quiet = TRUE)), r$depth_mean, field = 1, background = 0)

  if (nzchar(cell_bbox)) {
    bb <- as.numeric(strsplit(cell_bbox, ",")[[1]])
    r  <- crop(r, ext(bb[1], bb[2], bb[3], bb[4]))
    message(glue("restricted to bbox [{cell_bbox}]"))
  }
}

3 Write cell table (ocean cells) + global cell-id COG

Code 
if (need_build) {
  # ocean cells only (depth_mean non-NA); cell_id stays the global raster index
  d_cell <- as.data.frame(r, xy = TRUE, na.rm = FALSE) |>
    filter(!is.na(depth_mean)) |>
    transmute(
      cell_id = as.integer(cell_id),
      lon = x, lat = y,
      across(c(depth_mean, depth_min, depth_max, oxy_b_mean, oxy_mean, prim_prod_mean,
               ice_con_ann, salinity_b_mean, salinity_mean, sbt_an_mean, sst_an_mean, fao_area_m)),
      area_km2,
      in_usa = in_usa == 1,
      in_pra = in_pra == 1)

  con <- dbConnect(duckdb(cell_db))
  dbWriteTable(con, "cell", d_cell, overwrite = TRUE)
  dbExecute(con, "ALTER TABLE cell ADD PRIMARY KEY (cell_id)")
  dbDisconnect(con, shutdown = TRUE)

  # global cell-id COG in [-180,180] for titiler (uint32; no overviews to keep ids intact)
  writeRaster(r$cell_id, cellid_tif, filetype = "COG", datatype = "INT4U", overwrite = TRUE,
              gdal = c("COMPRESS=LZW", "OVERVIEWS=NONE"))
} else {
  message(glue("cell grid up to date (cellid_tif exists); set REDO_CELL_GRID=1 to rebuild"))
}
cell grid up to date (cellid_tif exists); set REDO_CELL_GRID=1 to rebuild

4 Summary + manifest

Code 
con  <- dbConnect(duckdb(cell_db, read_only = TRUE))
smry <- dbGetQuery(con, "
  SELECT count(*) AS n_cell,
         count(*) FILTER (WHERE in_usa) AS n_usa,
         count(*) FILTER (WHERE in_pra) AS n_pra,
         round(min(area_km2), 3) AS area_min_km2,
         round(max(area_km2), 3) AS area_max_km2,
         round(min(lon),1) AS lon_min, round(max(lon),1) AS lon_max,
         count(*) FILTER (WHERE depth_mean IS NULL) AS n_null_depth
  FROM cell")
# content fingerprint of the cell table (before disconnect)
h <- msens::hash_query(con, "cell")
dbDisconnect(con, shutdown = TRUE)
msens::report_table(smry, caption = "global 0.05° cell grid summary")
global 0.05° cell grid summary
n_cell n_usa n_pra area_min_km2 area_max_km2 lon_min lon_max n_null_depth
17072105 634208 343576 0.014 30.773 -180 180 0
Code 
# content-addressed manifest: deterministic (no wall-clock, no machine paths) ->
# downstream targets re-run only when the cell table's content actually changes
msens::write_manifest(
  here("data/manifests/build_cell_grid.json"), target = "build_cell_grid", content_hash = h,
  stats = c(as.list(smry), list(version = ver, bbox = cell_bbox,
            source = basename(bio_oracle_tif))),
  force = msens::force_target("build_cell_grid"))
message(glue("cell grid: {smry$n_cell} cells ({smry$n_usa} in US, {smry$n_pra} in PRAs) -> {cell_db}"))
cell grid: 17072105 cells (634208 in US, 343576 in PRAs) -> ~/_big/msens/derived/v8/sdm.duckdb