Ingest NMFS Core Distribution Areas → global 0.05° cells

Published

2026-07-13 15:00:09

NMFS Core Distribution Areas onto the global 0.05° cell grid. Currently one: Rice’s whale (Balaenoptera ricei), ESA-Endangered, so value = msens::compute_er_score("NMFS:EN") = 100 (extinction-risk score, never hard-coded). mdl_key = "ca_nmfs|{sp_id}". Same vector→cells pattern (msens::cells_from_ranges()) as the other range datasets; add rows to areas (with er_code) for future core areas.

1 Design

Figure 1: NMFS Core Distribution Areas → extinction-risk value on the global 0.05° grid

2 Setup + rasterize

Code 
librarian::shelf(
  arrow, DBI, dplyr, duckdb, fs, glue, here, jsonlite, logger, readr, sf, terra,
  tibble, MarineSensitivity/msens, quiet = T)
source(here("libs/paths.R"))   # dir_raw, dir_big_v, cellid_tif, ver
source(here("libs/vars.R"))    # redo_ingest
options(readr.show_col_types = F)

ds_key   <- "ca_nmfs"
dir_ca   <- glue("{dir_raw}/fisheries.noaa.gov/core-areas")
dir_dist <- glue("{dir_big_v}/marine-atlas/dist/dataset={ds_key}")
manifest <- here("data/manifests/ingest_ca_nmfs.json")
dir_create(c(dir_dist, path_dir(manifest)))
stopifnot("run build_cell_grid.qmd first" = file_exists(cellid_tif))

# core-area shapefile(s) -> species; extend this table as NMFS adds core areas.
# er_code drives the val via msens::compute_er_score() (never hard-code the val):
# Rice's whale is ESA-Endangered (NMFS) -> compute_er_score("NMFS:EN") = 100.
areas <- tribble(
  ~sp_id,               ~sci_name,            ~er_code,   ~shp,
  "Balaenoptera_ricei", "Balaenoptera ricei", "NMFS:EN",
    "shapefile_Rices_whale_core_distribution_area_Jun19_SERO/shapefile_Rices_whale_core_distribution_area_Jun19_SERO.shp") |>
  mutate(mdl_key = mdl_key_raw(ds_key, sp_id))
write_csv(areas |> select(mdl_key, sp_id, sci_name, er_code), glue("{dir_dist}/../model_{ds_key}.csv"))

if (redo_ingest && dir_exists(dir_dist)) { dir_delete(dir_dist); dir_create(dir_dist) }
for (i in seq_len(nrow(areas))) {
  a <- areas[i, ]; out_pq <- fs::path(dir_dist, a$sp_id, ext = "parquet")
  if (file_exists(out_pq)) { log_info("{a$sp_id}: skip"); next }
  ply <- st_read(glue("{dir_ca}/{a$shp}"), quiet = TRUE)
  d   <- cells_from_ranges(ply, cellid_tif, value = compute_er_score(a$er_code))
  if (nrow(d) == 0) { log_warn("{a$sp_id}: no cells"); next }
  msens::write_atlas_parquet(tibble(mdl_key = a$mdl_key, cell_id = d$cell_id, val = d$val), out_pq)
  log_info("{a$sp_id}: {nrow(d)} cells")
}
Spherical geometry (s2) switched off
Spherical geometry (s2) switched on

3 Verify + manifest

Code 
pq  <- dir_ls(dir_dist, glob = "*.parquet")
con <- dbConnect(duckdb())
smry <- dbGetQuery(con, glue(
  "SELECT count(DISTINCT mdl_key) n_models, count(*) n_cells, max(val) v_max
     FROM read_parquet('{dir_dist}/*.parquet')"))
# order-independent content fingerprint of the on-disk surface (not the ingest run)
h <- msens::hash_parquet(glue("{dir_dist}/*.parquet"), con)
dbDisconnect(con, shutdown = TRUE)
msens::report_table(smry, caption = "ca_nmfs model_cell surface")
ca_nmfs model_cell surface
n_models n_cells v_max
1 1969 100
Code 
# content-addressed manifest: deterministic, no wall-clock, no machine paths ->
# downstream targets re-run only when this surface's content actually changes
msens::write_manifest(
  manifest, target = "ingest_ca_nmfs", content_hash = h,
  stats = list(ds_key = ds_key, ver = ver, n_models = nrow(areas),
               n_cells = smry$n_cells, v_max = smry$v_max, n_parquet = length(pq)),
  force = msens::force_target("ingest_ca_nmfs"))