Ingest USFWS ESA Critical Habitat → global 0.05° cells

USFWS ESA critical-habitat polygons onto the global 0.05° cell grid. Like FWS ranges, critical habitat is designated for many terrestrial species (840 total), so we match v7’s marine-relevant subset (29 species) for comparability (data/v7_ch_fws_spcodes.csv). Value = msens::compute_er_score("FWS:{code}") from the ESA listing_st (Endangered→EN=100, Threatened→TN=50), never hard-coded. Whole range via cells_from_ranges. mdl_key = "ch_fws|{spcode}".

1 Design

Figure 1: USFWS critical habitat (v7-matched subset) → extinction-risk value on the global 0.05° grid

2 Ingest

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)
sf_use_s2(FALSE)

Spherical geometry (s2) switched off

ds_key   <- "ch_fws"
poly_shp <- glue("{dir_raw}/fws.gov/crithab_all_layers/crithab_poly.shp")
dir_dist <- glue("{dir_big_v}/marine-atlas/dist/dataset={ds_key}")
manifest <- here("data/manifests/ingest_ch_fws.json")
dir_create(c(dir_dist, path_dir(manifest)))
stopifnot("run build_cell_grid.qmd first" = file_exists(cellid_tif), file_exists(poly_shp))
v7_spcodes <- read_csv(here("data/v7_ch_fws_spcodes.csv"))$spcode   # 29 v7 ch_fws species

ply <- st_read(poly_shp, quiet = TRUE)
sp <- ply |> st_drop_geometry() |>
  filter(spcode %in% v7_spcodes) |>
  group_by(spcode, sciname, comname) |>
  summarise(status = if_else(any(grepl("Endangered", listing_st)), "EN", "TN"), .groups = "drop") |>
  mutate(er_code  = paste0("FWS:", status),
         er_score = compute_er_score(er_code),
         mdl_key  = mdl_key_raw(ds_key, spcode)) |>
  arrange(sciname)
write_csv(sp, glue("{dir_dist}/../model_{ds_key}.csv"))

if (redo_ingest && dir_exists(dir_dist)) { dir_delete(dir_dist); dir_create(dir_dist) }
done <- path_ext_remove(path_file(dir_ls(dir_dist, glob = "*.parquet")))
todo <- sp |> filter(!spcode %in% done)
log_info("{ds_key}: {length(done)} done; rasterizing {nrow(todo)} of {nrow(sp)} (v7-matched) ...")
for (i in seq_len(nrow(todo))) {
  r      <- todo[i, ]
  ply_sp <- ply |> filter(spcode == r$spcode) |> summarise(.groups = "drop")
  d <- cells_from_ranges(ply_sp, cellid_tif, value = r$er_score)
  if (nrow(d) == 0) { log_warn("  {r$spcode}: no cells"); next }
  msens::write_atlas_parquet(tibble(mdl_key = r$mdl_key, cell_id = d$cell_id, val = d$val),
                fs::path(dir_dist, r$spcode, ext = "parquet"))
}

although coordinates are longitude/latitude, st_union assumes that they are
planar

although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar
although coordinates are longitude/latitude, st_union assumes that they are
planar

log_info("{ds_key}: done")

3 Verify + manifest

pq  <- dir_ls(dir_dist, glob = "*.parquet")
con <- dbConnect(duckdb())
smry <- dbGetQuery(con, glue(
  "SELECT count(DISTINCT mdl_key) n_models, count(*) n_cells, min(val) v_min, 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 = "ch_fws model_cell surface")

ch_fws model_cell surface

n_models	n_cells	v_min	v_max
29	70190	50	100

# 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_ch_fws", content_hash = h,
  stats = list(ds_key = ds_key, ver = ver, n_models = nrow(sp),
               n_cells = smry$n_cells, v_min = smry$v_min, v_max = smry$v_max,
               n_parquet = length(pq)),
  force = msens::force_target("ingest_ch_fws"))