---
title: "Build the zone-set registry (spatial units as vintages, not per-version copies)"
msens:
target_name: build_zone_sets
workflow_type: schema
dependency: [build_cell_grid]
output: data/manifests/build_zone_sets.json
editor_options:
chunk_output_type: console
---
Every release so far baked its own copy of the spatial units — `ply_programareas_2026_v3.gpkg`
through `…_v8.gpkg`, each with its own `zone_cell` re-extracted by `score_zones.qmd`. That makes
the question users actually ask — *how did **this** Program Area's score change from v3 to v8?* —
unanswerable, because nothing guarantees two releases meant the same polygon.
This notebook identifies a spatial unit by its **geometry**, not by the release that used it:
```
zone_set_key = {zone_type}_{YYYY-MM} e.g. programarea_2026-03
```
The vintage is the year-month the geometry **first appeared on disk**, so the label is a checkable
fact rather than an editorial guess, and keys sort chronologically.
> **The per-version files are mostly the same layer copied forward.** Measured across every
> published gpkg, `programarea` is 8 files but only **2** distinct geometries (v2 alone; v3–v8 are
> byte-identical), and `ecoregion` is 10 files but **2**. So `zone_cell` — which depends only on
> (geometry × grid) — is computed **once** per `(zone_set_key, grid_id)` and reused by every release
> on that grid, instead of being re-extracted per version.
## Design
```{mermaid}
%%| label: fig-design
%%| fig-cap: "Fingerprint every zone layer, collapse identical ones into vintages, emit the registry"
flowchart LR
gp["ply_*.gpkg<br/>(top level + every version dir)"] --> gh["msens::zone_geom_hash()<br/>order-invariant WKB fingerprint"]
gh --> gr["group by (zone_type, geom_hash)<br/>vintage = earliest mtime YYYY-MM"]
gr --> vz["msens::validate_zone_sets()<br/>one key ↔ one geometry"]
vz --> csv["data/zone_sets.csv"]
```
## Setup
```{r}
#| label: setup
#| message: false
librarian::shelf(dplyr, glue, knitr, readr, MarineSensitivity/msens, quiet = TRUE)
source(here::here("libs/paths.R"))
# Every zone layer anywhere under derived/ — the top level and every version dir.
# `ply_ecoregion_programarea_*` is a CROSSWALK (ecoregion x programarea slivers),
# not a spatial unit in its own right, and `ply_boem-usa` is the study-area
# outline; neither is a zone set.
#
# Globbed lazily, at the point of USE, because the canonical subregions are
# DERIVED below and must be picked up by the scan. Globbing here at setup ran
# before that file existed, so the first run silently omitted it from the
# registry while still writing the layer.
zone_gpkgs <- function() {
g <- Sys.glob(c(path.expand(glue("{dir_derived}/ply_*.gpkg")),
path.expand(glue("{dir_derived}/*/ply_*.gpkg"))))
grep("ecoregion_programarea|boem-usa", g, value = TRUE, invert = TRUE)
}
message(glue("derived tree: {dir_derived}"))
```
::: {.callout-warning}
This scan is only as complete as the machine it runs on. The **server** holds every version's
`derived/` tree; a laptop holds only the current version's. Render this where the whole tree is,
or the registry will silently omit vintages.
:::
## Derive the canonical subregions
Subregions were never a published layer: each release synthesized them in the database as unions of
**its own** Planning/Program Areas, so they drifted (v1 `AK/AKL48/L48/USA`, v4–v6 `AK/GA/PA/USA`,
v7 adds `FULL`) and v2/v3 carry *two different sets under one `fld`*. Worse, a set derived from the
2026 Program Areas **cannot cover the Atlantic at all** — the 2026 National OCS Program has no
Atlantic areas — so it can never be the frame for comparing releases across all US waters.
The ecoregion layer already solves this. It is **byte-identical across v1–v8**, it spans every US
water, and it carries the rollup as a column: `region_key` / `region_name`. So the canonical
subregions are a dissolve of that column — no hand-written lookup to curate, and nothing to guess:
| region | ecoregions |
|---|---|
| `AK` Alaska | High Arctic, Chukchi & Beaufort Seas, East Bering Sea, **Gulf of Alaska** |
| `AT` Atlantic | NE Continental Shelf, SE Continental Shelf, Puerto Rico & USVI |
| `GA` Gulf of America | **Eastern**, **Western & Central** Gulf of America |
| `PA` Pacific | California Current, Washington/Oregon, Pacific Island Territories |
> Note `GOA` is the Gulf of **Alaska** while `EGOA`/`WCGOA` are the Gulf of **America** — the codes
> invite exactly the wrong inference, which is why this reads `region_key` from the layer rather
> than mapping ecoregion codes by hand.
The vintage is **inherited from the source ecoregions** (not the build date), since the geometry is
a pure function of them — so rebuilding never mints a spurious new vintage.
```{r}
#| label: subregions
#| message: false
librarian::shelf(sf, quiet = TRUE)
eco_src <- path.expand(glue("{dir_derived}/v1/ply_ecoregions_2025.gpkg"))
eco_vin <- format(file.info(eco_src)$mtime, "%Y-%m") # the ecoregion vintage
f_sub <- path.expand(glue("{dir_derived}/ply_subregions_usa_{eco_vin}.gpkg"))
eco <- sf::st_read(eco_src, quiet = TRUE)
stopifnot("ecoregion layer lacks region_key" = all(c("region_key", "region_name") %in% names(eco)))
sub <- eco |>
group_by(subregion_key = region_key, subregion_name = region_name) |>
summarise(n_ecoregions = n(), .groups = "drop") |>
sf::st_make_valid()
sf::st_write(sub, f_sub, delete_dsn = TRUE, quiet = TRUE)
knitr::kable(sf::st_drop_geometry(sub), caption = glue("Canonical subregions -> {basename(f_sub)}"))
```
`USA` (all waters) and `FULL` (the whole map extent) stay **pseudo-subregions**: they are not
polygons but the absence of a filter, so materializing them as geometry would double-count.
## Fingerprint every layer
```{r}
#| label: fingerprint
zone_type_of <- function(f)
dplyr::case_when(
grepl("programarea", f) ~ "programarea",
grepl("planarea", f) ~ "planarea",
grepl("ecoregion", f) ~ "ecoregion",
grepl("subregion", f) ~ "subregion",
.default = NA_character_)
gpkgs <- zone_gpkgs() # re-globbed AFTER the derived subregions exist
message(glue("scanning {length(gpkgs)} zone layer(s)"))
scan <- lapply(gpkgs, function(f) {
# zone_type FIRST: a program-area gpkg carries region_key, planarea_key AND
# programarea_key, so "the first *_key column" picks region_key (3 values for
# 20 features) -- which both mislabels the zones and, being a non-total
# ordering, leaves ties whose order depends on the file's feature order
zt <- zone_type_of(basename(f))
h <- msens::zone_geom_hash(f, zone_type = if (is.na(zt)) NULL else zt)
tibble::tibble(
source = sub(paste0(path.expand(dir_derived), "/"), "", f),
zone_type = zone_type_of(basename(f)),
n_zones = h$n,
key_col = h$key_col,
geom_hash = h$geom_hash,
mtime = file.info(f)$mtime)
}) |> bind_rows() |> filter(!is.na(zone_type))
# which release used it: the version directory it sits in ("" = top level, i.e.
# a staging copy no release scored)
scan <- scan |> mutate(ver = ifelse(grepl("/", source), sub("/.*", "", source), NA_character_))
knitr::kable(scan |> select(source, zone_type, n_zones, ver, geom_hash))
```
## Collapse identical layers into vintages
```{r}
#| label: vintages
vint <- scan |>
group_by(zone_type, geom_hash) |>
summarise(n_zones = first(n_zones),
# First appearance on disk -- EXCEPT the derived canonical subregions,
# whose vintage is inherited from the ecoregions they dissolve (it is
# encoded in their filename). Using their build mtime would mint a
# new vintage on every rebuild for geometry that never changed.
vintage = {
inherited <- sub(".*subregions_usa_([0-9]{4}-[0-9]{2})\\.gpkg$", "\\1",
grep("subregions_usa_", source, value = TRUE)[1])
if (!is.na(inherited) && grepl("^[0-9]{4}-[0-9]{2}$", inherited))
inherited else format(min(mtime), "%Y-%m")
},
n_files = n(),
versions = paste(sort(na.omit(ver)), collapse = " "),
source = source[which.min(mtime)],
.groups = "drop") |>
mutate(zone_set_key = msens::zone_set_key(zone_type, vintage)) |>
arrange(zone_type, vintage)
# Two DISTINCT geometries of one type first appearing in the same month would
# collide on the key. Stop rather than silently merge two different maps.
dup <- vint$zone_set_key[duplicated(vint$zone_set_key)]
if (length(dup))
stop(glue("zone_set_key collision ({paste(dup, collapse=', ')}): two distinct geometries ",
"first appeared in the same month. Disambiguate the vintage by hand (e.g. -a/-b) ",
"before publishing."))
msens::validate_zone_sets(vint)
knitr::kable(vint |> select(zone_set_key, n_zones, n_files, versions, geom_hash))
```
```{r}
#| label: summary
vint |>
count(zone_type, name = "distinct_vintages") |>
left_join(scan |> count(zone_type, name = "files"), by = "zone_type") |>
relocate(files, .after = zone_type) |>
knitr::kable(caption = "How many spatial units really exist, versus how many files")
```
## Write the registry
```{r}
#| label: write
f_csv <- here::here("data/zone_sets.csv")
# generated-but-committed, like a lockfile: deterministic given the files on
# disk, and reviewable in a diff when a new vintage appears
out <- vint |> select(zone_set_key, zone_type, vintage, n_zones, geom_hash, versions, source)
readr::write_csv(out, f_csv)
msens::validate_zone_sets(readr::read_csv(f_csv, show_col_types = FALSE))
message(glue("wrote {f_csv} ({nrow(out)} zone sets)"))
```
## Publish each vintage as PMTiles
Zone outlines are drawn by both apps, and until now came from
`ply_programareas_2026.pmtiles` / `ply_ecoregions_2025.pmtiles` on the file host —
**unversioned**, and inherited artifacts of an archived v7 notebook that nothing
re-runs. So the tiles a v3 map draws were whatever the last hand-run produced.
Built here per vintage, keyed like everything else, and published to S3 beside the
data:
```
zones/{zone_set_key}/zones.pmtiles
```
Vector tiles are versioned by their OWN vintage and shared across MST releases —
six releases share one program-area geometry, so they share one tile set.
```{r}
#| label: pmtiles
do_pmt <- Sys.getenv("ZONE_PMTILES") != "" && nzchar(Sys.which("tippecanoe"))
if (Sys.getenv("ZONE_PMTILES") != "" && !nzchar(Sys.which("tippecanoe")))
message("ZONE_PMTILES set but tippecanoe is not installed - skipping")
pmt <- tibble::tibble(zone_set_key = character(), pmtiles = character())
if (do_pmt) {
dir_pmt <- path.expand(glue("{dir_derived}/zones"))
for (i in seq_len(nrow(vint))) {
z <- vint[i, ]
src <- path.expand(glue("{dir_derived}/{z$source}"))
out <- glue("{dir_pmt}/{z$zone_set_key}/zones.pmtiles")
dir.create(dirname(out), recursive = TRUE, showWarnings = FALSE)
ply <- sf::st_read(src, quiet = TRUE)
kc <- msens::zone_key_col(z$zone_type, names(ply))
# carry the KEY plus any name column, so the app can label without a join
keep <- intersect(c(kc, sub("_key$", "_name", kc)), names(ply))
msens::publish_pmtiles(ply[, keep], out, layer = z$zone_type,
keep_attrs = keep, quiet = TRUE)
pmt <- dplyr::bind_rows(pmt, tibble::tibble(
zone_set_key = z$zone_set_key, pmtiles = as.character(out)))
message(glue(" {z$zone_set_key}: {round(file.info(out)$size/1024^2, 1)} MB"))
}
if (Sys.getenv("ZONE_PMTILES_NO_S3") == "") {
out <- system2("aws", c("s3", "sync", shQuote(dir_pmt),
shQuote(glue("{s3_atlas}/zones")),
"--exclude", shQuote("*"), "--include", shQuote("*.pmtiles"),
"--only-show-errors", "--no-progress"),
stdout = TRUE, stderr = TRUE)
if (!is.null(attr(out, "status")) && attr(out, "status") != 0)
stop("aws s3 sync failed: ", paste(out, collapse = "\n"))
message(glue("published {nrow(pmt)} zone tile set(s) to {s3_atlas}/zones"))
}
knitr::kable(pmt |> dplyr::mutate(mb = round(file.info(pmtiles)$size/1024^2, 1)) |>
dplyr::select(zone_set_key, mb))
} else {
message("ZONE_PMTILES unset - skipping tile build")
}
```
## What each release used, and what replaces it
The per-release subregions this supersedes — note they are not comparable to each other, and v2/v3
are ambiguous even within one database:
| release | subregion keys |
|---|---|
| v1 | `AK`, `AKL48`, `L48`, `USA` |
| v2, v3 | *two different sets under one `fld`* — `AK/AKL48/L48/USA` **and** `AK/GA/PA/USA` |
| v4–v6 | `AK`, `GA`, `PA`, `USA` |
| v7 | `AK`, `FULL`, `GA`, `PA`, `USA` |
| v8 | `AK`, `GA`, `PA`, `USA` |
Every one of them is missing the **Atlantic**, because all were unions of Program Areas and the 2026
National OCS Program has none there. The canonical set derived above does not have that hole, and
being a function of the v1–v8-identical ecoregions, it applies unchanged to every release.
## Target manifest
```{r}
#| label: target-manifest
msens::write_manifest(
here::here("data/manifests/build_zone_sets.json"),
target = "build_zone_sets",
content_hash = digest::digest(out, algo = "xxhash64"),
stats = list(
n_files = nrow(scan),
n_zone_sets = nrow(out),
programarea_n = sum(vint$zone_type == "programarea"),
ecoregion_n = sum(vint$zone_type == "ecoregion"),
planarea_n = sum(vint$zone_type == "planarea"),
subregion_n = sum(vint$zone_type == "subregion")),
force = msens::force_target("build_zone_sets"))
```