CHANGELOG.md

@@ -2,6 +2,13 @@
 
 All notable changes to **History** are listed here, newest first.
 
+## v0.13.0 — 2026-10-12 08:30
+- **Replaced the geometric sub-region grid with real named sub-national administrative regions (#16, PR 2).** The choropleth source `public/regions.geo.json` now contains real named administrative regions — states, provinces, governorates, oblasts, … (Natural Earth admin-1, 1:10m) — instead of the abstract ~100–200 km geometric grid. Each region keeps a stable `id` and a `name`. The map's data contract is **unchanged** (`properties.parent` → the parent country's ISO_A3, `properties.name` → the region name), so no change to the map code or the LLM contract: `WorldMap.tsx`, `App.tsx`, `shared/protocol.ts`, and `server/map.ts` all keep working as-is.
+- **Generator rewritten** — `scripts/generate-regions.mjs` (same path, dependency-free, `npm run generate:regions`) now groups the admin-1 features by their `adm0_a3` (the parent's ISO_A3), keeps only parents in the 177-country set, and id's each region with its ISO 3166-2 code when present (e.g. `US-CA`, `FR-75`, `RU-MOW`) else `<ISO_A3>-<nn>`; within a country, ids are disambiguated with a `-2/-3/…` suffix on collision (ISO 3166-2 codes are not always unique within a country — e.g. all nine of Bosnia's entities report `BA-BIH`). A country with no admin-1 data would collapse to one whole-country region `<ISO_A3>-01` (0 are needed at 1:10m, which covers all 177). Coordinates are rounded to 4 decimal places, and the file is Douglas–Peucker-simplified only if it exceeds the ~8 MB guard. Output is sorted by (parent, id) and fully deterministic.
+- **The dataset is 4,138 regions covering all 177 countries**, 2,831 KB (2.76 MB) after Douglas–Peucker (the 4-dp base is ~23 MB). Per-country region count: min 1 / median 15 / max 232 (United Kingdom). There are **0 whole-country fallbacks** — 1:10m admin-1 covers every one of the 177 (the five single-region countries are SOL, SAH, CYN, PRI, FLK). Samples: `US-CA` California and `US-NY` New York (US); `FR-75` Paris and `FR-13` Bouches-du-Rhône (France); `RU-MOW` Moskovskaya (Russia).
+- **Validation test updated** — `scripts/regions.test.ts` was rewritten for the new dataset: well-formed features, globally unique ids, all 177 countries covered (≥1 region each), sane region areas (median ≥ 100 km², no non-Antarctica region ≥ 5,000,000 km² — Antarctica is exempt because its full-360° wrap ring breaks the planar area model), and each country's region area summing to within [40%, 175%] of its own area (median within [70%, 130%] — the band is wider than a naive [40%, 110%] because Natural Earth's admin-1 polygons are digitized independently of the admin-0 country polygons and are not clipped to them). Runs in O(total vertices).
+- Part of #16.
+
 ## v0.12.0 — 2026-10-12 01:03
 - **Map now renders the sub-national regional dataset (#16, PR 2).** The choropleth source switches from `public/countries.geo.json` (177 countries) to `public/regions.geo.json` (**2,597 sub-regions**): each country is subdivided into ~100–200 km regions, shown as a finer choropleth with thin internal border lines. The MapLibre source/layer ids are renamed `countries`/`countries-fill`/`countries-line` → `regions`/`regions-fill`/`regions-line`, and the border line is thinned to `0.75` (from `1.25`) so the subdivision reads within each country.
 - **Ownership stays at country granularity.** A nation still owns whole countries, so every sub-region of a country shares that country's owning-nation color (the fill is "same colors as the country map, but subdivided"). The owner of a sub-region is resolved from the RegionMap via the feature's `properties.parent` (the parent country's ISO_A3) — never the sub-region id (`FRA-03`) nor the repacked top-level `feature.id`. Clicking a sub-region switches the dossier to its parent country's owning nation (unchanged behavior); the hover tooltip now also shows the sub-region name + its parent country; a small "Selected …" label appears in the dossier after a map click.

package.json

@@ -1,6 +1,6 @@
 {
   "name": "history",
-  "version": "0.12.0",
+  "version": "0.13.0",
   "private": true,
   "type": "module",
   "description": "AI-driven alternate-history simulator — the LLM is the world engine (Pax Historia style).",

scripts/generate-regions.mjs

@@ -2,49 +2,45 @@
 /**
  * scripts/generate-regions.mjs — dependency-free generator for `public/regions.geo.json`.
  *
- * Subdivides every country in `public/countries.geo.json` into sub-national
- * regions sized ~100–200 km on a side (adaptive: a regular lat/lon grid
- * clipped to the country, so vast/sparse countries — Siberia, Sahara,
- * Antarctica — end up with bigger regions while Europe gets smaller ones).
- * Each country gets a few to tens of regions (capped at 36); tiny countries
- * (and tiny antimeridian islet groups) stay a single region.
+ * Builds the sub-national region choropleth from REAL named administrative
+ * regions (Natural Earth admin-1: states, provinces, governorates, oblasts, …)
+ * instead of the old geometric grid the user rejected ("some stupid geometric
+ * grid" — it looked abstract, not like recognizable territories).
  *
- * Model: a country's lat/lon bbox is tiled with a cols × rows grid of
- * ~square-in-km cells. Each grid cell yields at most ONE region feature:
- * every sub-polygon of the country is clipped against the cell with
- * Sutherland–Hodgman (edges in order: left, bottom, right, top), clipped
- * pieces smaller than 300 km² are dropped, and the kept pieces become the
- * rings of that cell's region (a MultiPolygon when there is more than one).
- * Region ids are `<ISO_A3>-<nn>`, nn = 1..k in row-major order (rows
- * top→bottom by descending lat, left→right by ascending lon); names are
- * `<Country> R<nn>`. Coordinates are rounded to 4 decimal places (≈11 m).
- * Because a cell produces at most one region, the per-country count can
- * never exceed the (capped) cell count ≤ 36, even for archipelagos.
+ * Sources (both GeoJSON):
+ *   (a) the Natural Earth admin-1 file — real, named sub-national boundaries.
+ *       It is read from a NON-committed temp path (default /tmp/ne10m_admin1.geojson,
+ *       override with the NE_ADMIN1 env var) so the repo only stores the derived
+ *       output, not the 40 MB source. Natural Earth is public domain (CC0).
+ *   (b) the committed `public/countries.geo.json` (177 countries) — the
+ *       authoritative ISO_A3 set, the country names, and the whole-country
+ *       geometries used as the fallback for any country with no admin-1 data.
  *
- * Area: local planar approximation, x = lon·111.32·cos(midLat), y = lat·
- * 110.574 (km), shoelace per sub-polygon ring, summed in abs. Target region
- * area T = 12,000 km² (≈110 km side); target count n = clamp(round(area /
- * T), 1, 36). S = sqrt(bboxAreaKm2 / n); cols = max(1, round(wKm/S)),
- * rows = max(1, round(hKm/S)); the grid is then reduced if cols·rows > 36.
+ * Output contract (UNCHANGED from the old generator, so `WorldMap.tsx`,
+ * `App.tsx`, `shared/protocol.ts` and `server/map.ts` all keep working —
+ * ownership is resolved via `properties.parent`, the parent country's ISO_A3):
+ *   - a GeoJSON FeatureCollection
+ *   - each feature: a unique top-level `id`, `properties = { name, parent }`
+ *     (`parent` = the country ISO_A3 = the RegionMap key; `name` = the region's
+ *     human name), and a Polygon/MultiPolygon geometry.
  *
- * Antimeridian handling (three cases):
- *  - Small antimeridian island groups (area ≤ 50,000 km² — Fiji, Kiribati,
- *    Tuvalu, …): one region = the whole country (`<ISO_A3>-01`, plain name).
- *  - Large antimeridian countries made of ordinary sub-polygons (Russia):
- *    the sub-polygons are shifted into longitude space [0, 360) (lons < 0 →
- *    +360) so the country sits in one contiguous band, tiled there, and
- *    output rings that straddle ±180° are split at the antimeridian into a
- *    MultiPolygon with each ring fully on one side — standard −180..180
- *    GeoJSON that renders without cross-map lines.
- *  - Full-360° "polar" wrap rings (Antarctica's main ring, whose only
- *    longitude jump sits at the pole and which crosses the prime meridian
- *    once): such a ring self-intersects in any single planar representation,
- *    so it is replaced by two simple "wedge" rings — the eastern and western
- *    hemispheric halves, each closed along the prime meridian, the pole
- *    (lat −90), and the antimeridian. The wedges tile −180..180 space
- *    without overlap and are tiled like any other ring (no shifting needed).
+ * For each country in `countries.geo.json` (in file order, deterministic):
+ *   - if it has admin-1 sub-regions -> emit each as a feature, in the source
+ *     file's order. Region `id` = its ISO 3166-2 code (e.g. `US-CA`, `FR-75`,
+ *     `RU-MOW`) when present, else `<ISO_A3>-<nn>` (1-based index in file
+ *     order). Because ISO 3166-2 codes are not always unique WITHIN a country
+ *     (e.g. all nine of Bosnia's entities report `BA-BIH`), ids are
+ *     disambiguated with a `-2/-3/…` suffix on collision (BA-BIH, BA-BIH-2, …).
+ *   - if it has none (tiny / not covered by admin-1) -> emit ONE whole-country
+ *     region: `id = <ISO_A3>-01`, `name` = the country's name, geometry = the
+ *     country's geometry.
+ *
+ * Every coordinate is rounded to 4 decimal places (0.0001° ≈ 11 m). If the
+ * resulting file exceeds ~8 MB, a light Douglas–Peucker pass (0.05°, escalating
+ * through a small table if still too big) trims it under the guard while keeping
+ * shapes recognizable, and the before/after sizes are reported. Features are
+ * sorted by (parent, id). Fully deterministic: same inputs -> same bytes.
  *
- * Deterministic: same input → byte-identical output.
  * Run: `npm run generate:regions`
  */
 
@@ -55,466 +51,202 @@ import { fileURLToPath } from 'node:url';
 const here = path.dirname(fileURLToPath(import.meta.url));
 const root = path.resolve(here, '..');
 
-// --- parameters ------------------------------------------------------------------
-const T = 12_000; // target region area, km² (≈110 km side)
-const MAX_REGIONS = 36; // per-country cap ("tens")
-const SLIVER_KM2 = 300; // drop clipped pieces smaller than this
-const BIG_KM2 = 50_000; // antimeridian countries at/above this are tiled, not 1 region
-const KM_LON = 111.32; // km per degree of longitude (at the equator)
-const KM_LAT = 110.574; // km per degree of latitude
+// --- parameters ----------------------------------------------------------------
+const SIZE_GUARD_KB = 8192; // ~8 MB file-size guard (task spec)
+const SIZE_TARGET_KB = 6144; // land under ~6 MB once simplified (margin under the guard)
+const DP_STEPS = [0.05, 0.1, 0.2, 0.3, 0.5, 0.75, 1.0]; // Douglas–Peucker tolerances (degrees)
+const R4 = 4; // base coordinate precision (decimal places, ≈11 m)
 
-// --- geo helpers (same planar approximation as scripts/regions.test.ts) -----------
-const rad = (deg) => (deg * Math.PI) / 180;
-const r4 = (v) => Math.round(v * 1e4) / 1e4;
+// Natural Earth admin-1 source (NON-committed temp file; override with NE_ADMIN1).
+const INPUT = process.env.NE_ADMIN1 || '/tmp/ne10m_admin1.geojson';
 
-function polygons(geom) {
-  return geom.type === 'Polygon' ? [geom.coordinates] : geom.coordinates;
-}
+// --- small geo helpers ---------------------------------------------------------
+const rN = (v, N) => Math.round(v * Math.pow(10, N)) / Math.pow(10, N);
 
-function allCoords(geom) {
+/** Round a ring to N decimals, drop consecutive duplicates, and close it. */
+function roundRing(ring, N) {
   const out = [];
-  for (const poly of polygons(geom)) for (const ring of poly) for (const c of ring) out.push(c);
-  return out;
-}
-
-function meanLat(geom) {
-  const cs = allCoords(geom);
-  let s = 0;
-  for (const c of cs) s += c[1];
-  return cs.length ? s / cs.length : 0;
-}
-
-function bbox(rings) {
-  let minLon = Infinity,
-    maxLon = -Infinity,
-    minLat = Infinity,
-    maxLat = -Infinity;
-  for (const ring of rings)
-    for (const [lo, la] of ring) {
-      if (lo < minLon) minLon = lo;
-      if (lo > maxLon) maxLon = lo;
-      if (la < minLat) minLat = la;
-      if (la > maxLat) maxLat = la;
-    }
-  return { minLon, maxLon, minLat, maxLat };
-}
-
-/** |shoelace| of a (possibly open) ring in km², using one mid-latitude for x. */
-function ringAreaKm2(ring, midLat) {
-  const kx = KM_LON * Math.cos(rad(midLat));
-  let a = 0;
-  const n = ring.length;
-  for (let i = 0; i < n; i++) {
-    const [lon1, lat1] = ring[i];
-    const [lon2, lat2] = ring[(i + 1) % n];
-    a += lon1 * kx * (lat2 * KM_LAT) - lon2 * kx * (lat1 * KM_LAT);
-  }
-  return Math.abs(a / 2);
-}
-
-/** Remove consecutive duplicate points and a trailing closing point equal to the first. */
-function dedupeRing(ring) {
-  const out = [];
-  for (const p of ring) {
+  for (const [x, y] of ring) {
+    const p = [rN(x, N), rN(y, N)];
     const last = out[out.length - 1];
     if (!last || last[0] !== p[0] || last[1] !== p[1]) out.push(p);
   }
+  // drop a trailing point equal to the first (a redundant close), then re-close
   if (
     out.length > 1 &&
     out[0][0] === out[out.length - 1][0] &&
     out[0][1] === out[out.length - 1][1]
   )
     out.pop();
-  return out;
+  if (out.length < 3) return out; // degenerate; caller drops it
+  return [...out, [...out[0]]];
 }
 
-/** Close an open ring by appending a copy of the first point. */
-function closeRing(ring) {
-  return ring.length ? [...ring, [...ring[0]]] : ring;
+/** Perpendicular distance (degrees) from point p to the segment a→b. */
+function perpDist(p, a, b) {
+  const dx = b[0] - a[0];
+  const dy = b[1] - a[1];
+  const l2 = dx * dx + dy * dy;
+  if (l2 === 0) return Math.hypot(p[0] - a[0], p[1] - a[1]);
+  let t = ((p[0] - a[0]) * dx + (p[1] - a[1]) * dy) / l2;
+  t = Math.max(0, Math.min(1, t));
+  return Math.hypot(p[0] - (a[0] + t * dx), p[1] - (a[1] + t * dy));
 }
 
-/**
- * The tiling/area subject rings of a feature, in standard −180..180 space:
- * one open ring per sub-polygon (outer ring only — holes are ignored by
- * design: Sutherland–Hodgman clips simple rings and any enclosed water is
- * at most a lake, which the area-sanity tolerance in the test absorbs).
- * A full-360° wrap ring is replaced by its two hemispheric wedges (see
- * makeWedges); if that is not possible the ring is returned as-is (the
- * caller may still shift it).
- */
-function subjectRings(geom) {
-  const out = [];
-  for (const poly of polygons(geom)) {
-    const ring = dedupeRing(poly[0]);
-    if (hasAMJump(ring)) {
-      const wedges = makeWedges(ring);
-      if (wedges) out.push(...wedges);
-      else out.push(ring);
-    } else {
-      out.push(ring);
+/** Douglas–Peucker simplification of an OPEN point sequence (keeps endpoints). */
+function dpOpen(pts, tol) {
+  if (pts.length < 3) return pts;
+  let maxD = 0;
+  let idx = 0;
+  const a = pts[0];
+  const b = pts[pts.length - 1];
+  for (let i = 1; i < pts.length - 1; i++) {
+    const d = perpDist(pts[i], a, b);
+    if (d > maxD) {
+      maxD = d;
+      idx = i;
     }
   }
-  return out;
-}
-
-// --- Sutherland–Hodgman: clip an open subject ring against one half-plane --------
-function clipEdge(poly, inside, intersect) {
-  const out = [];
-  for (let i = 0; i < poly.length; i++) {
-    const a = poly[i];
-    const b = poly[(i + 1) % poly.length];
-    const aIn = inside(a);
-    const bIn = inside(b);
-    if (aIn && bIn) out.push(b);
-    else if (aIn && !bIn) out.push(intersect(a, b));
-    else if (!aIn && bIn) {
-      out.push(intersect(a, b));
-      out.push(b);
-    }
+  if (maxD > tol) {
+    const left = dpOpen(pts.slice(0, idx + 1), tol);
+    const right = dpOpen(pts.slice(idx), tol);
+    return left.concat(right.slice(1));
   }
-  return out;
+  return [a, b];
 }
 
-/**
- * Clip an open subject ring against an axis-aligned (lon/lat) rectangle.
- * Edges are processed in order: left, bottom, right, top. Returns the clipped
- * open ring, or [] if empty.
- */
-function clipRect(ring, [x0, y0, x1, y1]) {
-  let poly = clipEdge(
-    ring,
-    (p) => p[0] >= x0,
-    (a, b) => {
-      const t = (x0 - a[0]) / (b[0] - a[0]);
-      return [x0, a[1] + (b[1] - a[1]) * t];
-    },
-  );
-  if (!poly.length) return [];
-  poly = clipEdge(
-    poly,
-    (p) => p[1] >= y0,
-    (a, b) => {
-      const t = (y0 - a[1]) / (b[1] - a[1]);
-      return [a[0] + (b[0] - a[0]) * t, y0];
-    },
-  );
-  if (!poly.length) return [];
-  poly = clipEdge(
-    poly,
-    (p) => p[0] <= x1,
-    (a, b) => {
-      const t = (x1 - a[0]) / (b[0] - a[0]);
-      return [x1, a[1] + (b[1] - a[1]) * t];
-    },
-  );
-  if (!poly.length) return [];
-  poly = clipEdge(
-    poly,
-    (p) => p[1] <= y1,
-    (a, b) => {
-      const t = (y1 - a[1]) / (b[1] - a[1]);
-      return [a[0] + (b[0] - a[0]) * t, y1];
-    },
-  );
-  return poly;
+/** Douglas–Peucker a CLOSED ring (first == last) at tolerance `tol` (0 = skip). */
+function dpRing(ring, tol) {
+  if (tol <= 0 || ring.length < 3) return ring;
+  let pts = ring;
+  if (
+    pts[0][0] === pts[pts.length - 1][0] &&
+    pts[0][1] === pts[pts.length - 1][1]
+  )
+    pts = pts.slice(0, -1);
+  const simplified = dpOpen(pts, tol);
+  return [...simplified, [...simplified[0]]];
 }
 
-// --- antimeridian handling ---------------------------------------------------------
-
-/** True if the ring has an edge with |Δlon| > 180 (an antimeridian jump). */
-function hasAMJump(ring) {
-  const n = ring.length;
-  for (let i = 0; i < n; i++)
-    if (Math.abs(ring[i][0] - ring[(i + 1) % n][0]) > 180) return true;
-  return false;
+/** Round (and optionally simplify) every ring of a Polygon/MultiPolygon geometry. */
+function transformGeom(geom, tol) {
+  const polys = geom.type === 'Polygon' ? [geom.coordinates] : geom.coordinates;
+  const out = polys
+    .map((poly) => poly.map((ring) => roundRing(ring, R4)).map((ring) => (tol > 0 ? dpRing(ring, tol) : ring)))
+    .filter((poly) => poly.length >= 1 && poly.every((ring) => ring.length >= 3));
+  if (!out.length) return geom; // nothing survived simplification — keep the rounded original
+  if (geom.type === 'Polygon') return { type: 'Polygon', coordinates: out[0] };
+  return { type: 'MultiPolygon', coordinates: out };
 }
 
-function crossesAntimeridian(geom) {
-  let minLon = Infinity,
-    maxLon = -Infinity;
-  for (const [lo] of allCoords(geom)) {
-    if (lo < minLon) minLon = lo;
-    if (lo > maxLon) maxLon = lo;
-  }
-  return maxLon > 0 && minLon < 0 && maxLon - minLon > 180;
-}
-
-/**
- * Replace a full-360° wrap ring (one antimeridian jump, one prime-meridian
- * crossing) by two simple hemispheric "wedge" rings, each closed along the
- * prime meridian (x = 0), the pole (lat = −90), and the antimeridian
- * (x = ±180). Both wedges lie in standard −180..180 space, are disjoint, and
- * together cover exactly the original ring's enclosed region. Returns
- * [westRing, eastRing] (open rings), or null if the ring has an unexpected
- * shape (no wrap handling possible).
- */
-function makeWedges(ring) {
-  const n = ring.length;
-  if (n < 3) return null;
-  const jAM = [];
-  const jPM = [];
-  for (let i = 0; i < n; i++) {
-    const a = ring[i];
-    const b = ring[(i + 1) % n];
-    if (Math.abs(a[0] - b[0]) > 180) jAM.push(i);
-    else if (a[0] !== 0 && b[0] !== 0 && (a[0] < 0) !== (b[0] < 0)) jPM.push(i);
-  }
-  if (jAM.length !== 1 || jPM.length !== 1) return null;
-
-  const iAM = jAM[0];
-  const a = ring[iAM]; // the +180 end of the jump
-  const b = ring[(iAM + 1) % n]; // the −180 end of the jump
-  if (Math.abs(a[0] - 180) > 1e-6 || Math.abs(b[0] + 180) > 1e-6) return null;
-
-  const iPM = jPM[0];
-  const p = ring[iPM];
-  const q = ring[(iPM + 1) % n];
-  const t = (0 - p[0]) / (q[0] - p[0]);
-  const P = [0, p[1] + (q[1] - p[1]) * t]; // the exact prime-meridian crossing
-
-  // Chain from b (−180) forward to P: everything in between must have lon ≤ 0.
-  const west = [];
-  for (let i = (iAM + 1) % n; ; i = (i + 1) % n) {
-    const pt = ring[i];
-    if (pt[0] > 0) return null;
-    west.push(pt);
-    if (i === iPM) break;
-  }
-  west.push(P);
-  // Chain from P forward to a (+180): everything in between must have lon ≥ 0.
-  const east = [P];
-  for (let i = (iPM + 1) % n; ; i = (i + 1) % n) {
-    const pt = ring[i];
-    if (pt[0] < 0) return null;
-    east.push(pt);
-    if (i === iAM) break;
-  }
-
-  const pole = (x) => [x, -90];
-  // West wedge: P → pole along x=0 → pole along the pole → up the −180 side →
-  // the western coast back to P. East wedge: the mirror image on x ≥ 0.
-  const westWedge = dedupeRing([P, pole(0), pole(-180), b, ...west.slice(1)]);
-  // The east chain q→…→a runs in *increasing* lon, so to trace the coast back
-  // from a to q (keeping the wedge simple) we walk it in reverse.
-  const eastWedge = dedupeRing([P, pole(0), pole(180), a, ...east.slice(1, -1).reverse()]);
-  if (westWedge.length < 3 || eastWedge.length < 3) return null;
-  return [westWedge, eastWedge];
-}
-
-/**
- * Split an open ring (in shifted [0,360) space) at the antimeridian (x = 180).
- *
- * Instead of walking the ring between crossings (which breaks on jagged
- * coastlines that cross the antimeridian many times), clip the ring against
- * the two half-planes x ≤ 180 and x ≥ 180 with Sutherland–Hodgman. Each piece
- * is a proper closed-on-the-split-line polygon, the two pieces share only the
- * x = 180 boundary, and the split is area-preserving (the boundary has zero
- * area). Both pieces are returned in shifted space; the caller remaps any
- * x > 180 piece back into −180..0.
- */
-function splitAtAM(ring) {
-  if (ring.length < 3) return [];
-  const cross = (a, b) => {
-    const t = (180 - a[0]) / (b[0] - a[0]);
-    return [180, a[1] + (b[1] - a[1]) * t];
-  };
-  const left = clipEdge(ring, (p) => p[0] <= 180, cross); // x in [0, 180]
-  // Right piece: x in [180, 360). Remap EVERY vertex by -360 (including the
-  // shared x = 180 boundary vertex → -180) so the piece sits fully on the west
-  // side of the antimeridian with no ±180 jump inside the ring.
-  const right = clipEdge(ring, (p) => p[0] >= 180, cross).map(
-    ([lo, la]) => [lo - 360, la],
-  );
-  const pieces = [];
-  if (left.length >= 3) pieces.push(left);
-  if (right.length >= 3) pieces.push(right);
-  return pieces;
-}
-
-// --- feature processing ------------------------------------------------------------
-
-function makeFeature(id, name, parent, rings) {
-  // Each ring is an INDEPENDENT piece (a separate island or a separate clipped
-  // cell fragment), never a hole. GeoJSON: a MultiPolygon is an array of
-  // polygons, each polygon an array of rings — so wrap every ring in its own
-  // (single-ring) polygon. ([rings] would instead make ring 0 the outer and
-  // every other ring a hole, swallowing e.g. Fiji's other islands.)
-  return {
-    type: 'Feature',
-    id,
-    properties: { name, parent },
-    geometry:
-      rings.length === 1
-        ? { type: 'Polygon', coordinates: [rings[0]] }
-        : { type: 'MultiPolygon', coordinates: rings.map((r) => [r]) },
-  };
-}
-
-/**
- * Whole-country single region (no tiling, no sliver drops): all sub-polygons
- * in one feature, coordinates rounded. Used for tiny countries (n = 1) and
- * small antimeridian-crossing island groups.
- */
-function wholeCountryFeature(feature) {
-  const id = feature.id;
-  const name = feature.properties.name;
-  const rings = polygons(feature.geometry).map((poly) =>
-    closeRing(dedupeRing(poly[0]).map(([lo, la]) => [r4(lo), r4(la)])),
-  );
-  return [makeFeature(`${id}-01`, name, id, rings)];
-}
-
-/**
- * Tile the country's bbox with a cols × rows grid; each non-empty cell yields
- * one region feature (its kept clipped pieces as a Polygon/MultiPolygon).
- * `shifted` = true for large non-wrapping antimeridian countries (Russia):
- * work in [0, 360) longitude space, then split output rings at ±180°.
- */
-function tileFeature(feature, rings, shifted) {
-  const id = feature.id;
-  const name = feature.properties.name;
-  const midLat = meanLat(feature.geometry);
-  const area = rings.reduce((s, r) => s + ringAreaKm2(r, midLat), 0);
-
-  const n = Math.min(MAX_REGIONS, Math.max(1, Math.round(area / T)));
-  if (n === 1) return wholeCountryFeature(feature);
-
-  const shiftLon = (lo) => (shifted && lo < 0 ? lo + 360 : lo);
-  const work = rings.map((r) => r.map(([lo, la]) => [shiftLon(lo), la]));
-  const bb = bbox(work);
-
-  const kx = KM_LON * Math.cos(rad(midLat));
-  const wKm = (bb.maxLon - bb.minLon) * kx;
-  const hKm = (bb.maxLat - bb.minLat) * KM_LAT;
-  const S = Math.sqrt((wKm * hKm) / n);
-  let cols = Math.max(1, Math.round(wKm / S));
-  let rows = Math.max(1, Math.round(hKm / S));
-  while (cols * rows > MAX_REGIONS) {
-    if (cols >= rows) cols = Math.max(1, Math.floor(MAX_REGIONS / rows));
-    else rows = Math.max(1, Math.floor(MAX_REGIONS / cols));
-  }
-
-  const dLon = (bb.maxLon - bb.minLon) / cols;
-  const dLat = (bb.maxLat - bb.minLat) / rows;
-
-  const out = [];
-  let nn = 0;
-  for (let r = 0; r < rows; r++) {
-    const yTop = bb.maxLat - r * dLat;
-    const yBot = bb.maxLat - (r + 1) * dLat;
-    for (let c = 0; c < cols; c++) {
-      const xLeft = bb.minLon + c * dLon;
-      const xRight = bb.minLon + (c + 1) * dLon;
-      // One region per cell: clip every sub-ring, keep the pieces that clear
-      // the sliver threshold, and merge them into a single feature.
-      const cellPieces = [];
-      for (const ring of work) {
-        const clipped = clipRect(ring, [xLeft, yBot, xRight, yTop]);
-        if (clipped.length < 3) continue;
-        if (ringAreaKm2(clipped, midLat) < SLIVER_KM2) continue;
-        // splitAtAM returns pieces already in standard −180..180 space (the
-        // x > 180 half remapped), so no further per-vertex remap is needed.
-        const pieces = shifted
-          ? splitAtAM(dedupeRing(clipped)).map((pc) => dedupeRing(pc))
-          : [dedupeRing(clipped)];
-        for (const pc of pieces)
-          if (pc.length >= 3) cellPieces.push(pc.map(([lo, la]) => [r4(lo), r4(la)]));
-      }
-      if (!cellPieces.length) continue;
-      const closed = cellPieces
-        .map(dedupeRing)
-        .filter((pc) => pc.length >= 3)
-        .map(closeRing);
-      if (!closed.length) continue;
-      nn += 1;
-      const nnStr = String(nn).padStart(2, '0');
-      out.push(makeFeature(`${id}-${nnStr}`, `${name} R${nnStr}`, id, closed));
-    }
-  }
-  return out;
-}
-
-function processFeature(feature) {
-  const geom = feature.geometry;
-  const am = crossesAntimeridian(geom);
-  if (am && featureAreaKm2(feature) <= BIG_KM2) {
-    // Small antimeridian island group (Fiji, Kiribati, Tuvalu, …): one region.
-    return wholeCountryFeature(feature);
-  }
-  // Wrapping rings (Antarctica) are already in standard space as wedges — no
-  // shifting; non-wrapping antimeridian countries (Russia) are shifted.
-  const rings = subjectRings(geom);
-  const hasWrap = polygons(geom).some((p) => hasAMJump(dedupeRing(p[0])));
-  const shifted = am && !hasWrap;
-  return tileFeature(feature, rings, shifted);
-}
-
-/** Total planar area of a feature's tiling rings (same convention as tileFeature). */
-function featureAreaKm2(feature) {
-  const midLat = meanLat(feature.geometry);
-  return subjectRings(feature.geometry).reduce((s, r) => s + ringAreaKm2(r, midLat), 0);
-}
-
-// --- main ---------------------------------------------------------------------------
-
+// --- main ----------------------------------------------------------------------
 function main() {
-  const src = readFileSync(path.join(root, 'public', 'countries.geo.json'), 'utf8');
-  const countries = JSON.parse(src);
+  const countries = JSON.parse(readFileSync(path.join(root, 'public', 'countries.geo.json'), 'utf8'));
+  const admin1 = JSON.parse(readFileSync(INPUT, 'utf8'));
+  const cset = new Set(countries.features.map((c) => c.id));
 
-  const features = [];
-  const amSingletons = [];
-  const amTiled = [];
-  for (const f of countries.features) {
-    if (crossesAntimeridian(f.geometry)) {
-      (featureAreaKm2(f) <= BIG_KM2 ? amSingletons : amTiled).push(f.id);
+  // Group admin-1 features by parent ISO_A3 (their `adm0_a3`), in file order,
+  // keeping only parents that are one of our 177 countries.
+  const byParent = new Map();
+  for (const f of admin1.features) {
+    const p = f.properties.adm0_a3;
+    if (typeof p === 'string' && cset.has(p)) {
+      if (!byParent.has(p)) byParent.set(p, []);
+      byParent.get(p).push(f);
+    }
+  }
+
+  // Build the base specs (raw geometry + unique id + name + parent) in
+  // countries.geo.json file order, then per-region source order — deterministic.
+  const specs = [];
+  const seenIds = new Set();
+  const nextId = (base) => {
+    let id = base;
+    let k = 1;
+    while (seenIds.has(id)) {
+      k += 1;
+      id = `${base}-${k}`;
+    }
+    seenIds.add(id);
+    return id;
+  };
+
+  let fallbackCount = 0;
+  const fallbackCountries = [];
+  for (const c of countries.features) {
+    const regs = byParent.get(c.id) ?? [];
+    if (regs.length) {
+      let nn = 0;
+      for (const r of regs) {
+        nn += 1;
+        const iso = r.properties.iso_3166_2;
+        const base = typeof iso === 'string' && iso ? iso : `${c.id}-${String(nn).padStart(2, '0')}`;
+        const name = typeof r.properties.name === 'string' && r.properties.name ? r.properties.name : c.properties.name;
+        specs.push({ id: nextId(base), name, parent: c.id, geom: r.geometry });
+      }
+    } else {
+      fallbackCount += 1;
+      fallbackCountries.push(c.id);
+      specs.push({ id: nextId(`${c.id}-01`), name: c.properties.name, parent: c.id, geom: c.geometry });
+    }
+  }
+
+  const sortById = (a, b) => (a.properties.parent !== b.properties.parent ? (a.properties.parent < b.properties.parent ? -1 : 1) : a.id.localeCompare(b.id, 'en', { numeric: true }));
+
+  const build = (tol) => {
+    const fc = {
+      type: 'FeatureCollection',
+      features: specs.map((s) => ({
+        type: 'Feature',
+        id: s.id,
+        properties: { name: s.name, parent: s.parent },
+        geometry: transformGeom(s.geom, tol),
+      })),
+    };
+    fc.features.sort(sortById);
+    return fc;
+  };
+
+  // Base: 4 decimal places, no simplification.
+  let fc = build(0);
+  let json = JSON.stringify(fc);
+  let sizeKB = Buffer.byteLength(json) / 1024;
+  console.log(`base (4dp, no DP): ${sizeKB.toFixed(1)} KB  (${fc.features.length} features)`);
+
+  // File-size guard: if the base is over the guard, simplify until it is under
+  // the target (or we exhaust the tolerance table).
+  let usedTol = 0;
+  if (sizeKB > SIZE_GUARD_KB) {
+    for (const t of DP_STEPS) {
+      fc = build(t);
+      json = JSON.stringify(fc);
+      sizeKB = Buffer.byteLength(json) / 1024;
+      usedTol = t;
+      console.log(`DP ${t}: ${sizeKB.toFixed(1)} KB`);
+      if (sizeKB < SIZE_TARGET_KB) break;
     }
-    features.push(...processFeature(f));
   }
-  // Sorted by (parent, nn) — the processing order already guarantees this, but
-  // sort explicitly so the contract holds regardless of traversal changes.
-  features.sort((a, b) => {
-    if (a.properties.parent !== b.properties.parent)
-      return a.properties.parent < b.properties.parent ? -1 : 1;
-    return a.id.localeCompare(b.id, 'en', { numeric: true });
-  });
 
-  const out = { type: 'FeatureCollection', features };
-  const json = JSON.stringify(out);
   writeFileSync(path.join(root, 'public', 'regions.geo.json'), json + '\n');
 
-  // --- summary ----------------------------------------------------------------------
+  // --- summary ----------------------------------------------------------------
   const counts = new Map();
-  for (const f of features)
-    counts.set(f.properties.parent, (counts.get(f.properties.parent) ?? 0) + 1);
+  for (const f of fc.features) counts.set(f.properties.parent, (counts.get(f.properties.parent) ?? 0) + 1);
   const sortedCounts = [...counts.values()].sort((a, b) => a - b);
   const median =
     sortedCounts.length % 2
       ? sortedCounts[Math.floor(sortedCounts.length / 2)]
       : (sortedCounts[sortedCounts.length / 2 - 1] + sortedCounts[sortedCounts.length / 2]) / 2;
+  const maxCountry = [...counts.entries()].sort((a, b) => b[1] - a[1])[0];
+  const maxCountryName = countries.features.find((c) => c.id === maxCountry[0])?.properties.name ?? '';
 
-  // Max region area (planar approximation, same mid-latitude convention as the
-  // tiling: the parent country's mean latitude).
-  const midLatByParent = new Map(countries.features.map((f) => [f.id, meanLat(f.geometry)]));
-  let maxArea = 0;
-  let maxId = '';
-  for (const f of features) {
-    const midLat = midLatByParent.get(f.properties.parent);
-    let a = 0;
-    for (const poly of polygons(f.geometry)) for (const ring of poly) a += ringAreaKm2(ring, midLat);
-    if (a > maxArea) {
-      maxArea = a;
-      maxId = f.id;
-    }
-  }
-
-  console.log(`regions: ${features.length} total (${countries.features.length} countries)`);
-  console.log(`file:    public/regions.geo.json — ${(Buffer.byteLength(json + '\n') / 1024).toFixed(1)} KB`);
-  console.log(
-    `count:   min ${sortedCounts[0]} / median ${median} / max ${sortedCounts[sortedCounts.length - 1]} per country`,
-  );
-  console.log(`max region: ${maxId} ≈ ${Math.round(maxArea).toLocaleString('en-US')} km²`);
-  if (amSingletons.length) console.log(`AM → 1 region: ${amSingletons.join(', ')}`);
-  if (amTiled.length) console.log(`AM → tiled (shifted/wedges): ${amTiled.join(', ')}`);
+  const finalKB = Buffer.byteLength(json) / 1024;
+  console.log(`\nwrote public/regions.geo.json — ${finalKB.toFixed(1)} KB (${(finalKB / 1024).toFixed(2)} MB)`);
+  console.log(`features: ${fc.features.length} total across ${counts.size} of ${countries.features.length} countries`);
+  console.log(`per-country count: min ${sortedCounts[0]} / median ${median} / max ${sortedCounts[sortedCounts.length - 1]}`);
+  console.log(`max-count country: ${maxCountry[0]} (${maxCountryName}) with ${maxCountry[1]}`);
+  if (usedTol) console.log(`simplification: Douglas–Peucker @ ${usedTol}° (base was ${SIZE_GUARD_KB}+ KB)`);
+  console.log(`whole-country fallback: ${fallbackCount} ${fallbackCount ? `(${fallbackCountries.join(', ')})` : '(none)'}`);
 }
 
 main();

scripts/regions.test.ts

@@ -1,15 +1,22 @@
 /**
  * scripts/regions.test.ts — validation for the committed `public/regions.geo.json`.
  *
- * Self-contained (no import of the generator): it re-implements the generator's
- * exact planar area model so the assertions measure the same quantity the
- * generator used to size and clip regions:
+ * The dataset is now REAL named sub-national administrative regions (Natural
+ * Earth admin-1: states, provinces, governorates, oblasts, …) rather than the
+ * old geometric grid — one feature per named region, each carrying
+ * `properties = { name, parent }` where `parent` is the country's ISO_A3 (the
+ * RegionMap key the map resolves ownership through). Tiny or uncovered
+ * countries collapse to a single whole-country region (`<ISO_A3>-01`).
+ *
+ * This test is self-contained (no import of the generator): it re-implements
+ * the same planar area model the generator/analysis use, so the area-based
+ * assertions measure the same quantity:
  *
  *   x = lon · 111.32 · cos(midLat),  y = lat · 110.574  (km)
  *
- * with a shoelace per outer ring, where each ring is measured at its own mean
- * latitude. Outer rings only — the generator never emits holes (every ring is an
- * independent piece: a separate island or a separate clipped cell fragment).
+ * with a shoelace per OUTER ring, each ring measured at its own mean latitude.
+ * Outer rings only — a region's inner rings (holes) are lakes and are ignored
+ * by design, the same way the coverage band already absorbed them.
  *
  * It reads each file exactly once and runs in O(total vertices).
  */
@@ -21,7 +28,7 @@ import { describe, expect, it } from 'vitest';
 const here = path.dirname(fileURLToPath(import.meta.url));
 const root = path.resolve(here, '..');
 
-// --- area model (must match scripts/generate-regions.mjs exactly) -----------------
+// --- area model (planar approximation; matches scripts/generate-regions.mjs) --
 const KM_LON = 111.32;
 const KM_LAT = 110.574;
 const rad = (deg: number) => (deg * Math.PI) / 180;
@@ -59,7 +66,14 @@ function featureAreaKm2(geom: Geometry): number {
   }
   return total;
 }
-// -------------------------------------------------------------------------------
+
+function median(nums: number[]): number {
+  if (!nums.length) return 0;
+  const s = [...nums].sort((a, b) => a - b);
+  const m = Math.floor(s.length / 2);
+  return s.length % 2 ? s[m] : (s[m - 1] + s[m]) / 2;
+}
+// ------------------------------------------------------------------------------
 
 interface RegionFeature {
   type: 'Feature';
@@ -73,96 +87,116 @@ interface CountryFeature {
   properties: { name: string };
   geometry: Geometry;
 }
-interface FeatureCollection<T> {
-  type: 'FeatureCollection';
-  features: T[];
-}
 
 // Read each file exactly once.
 const regions = JSON.parse(
   readFileSync(path.join(root, 'public', 'regions.geo.json'), 'utf8'),
-) as FeatureCollection<RegionFeature>;
+) as { type: 'FeatureCollection'; features: RegionFeature[] };
 const countries = JSON.parse(
   readFileSync(path.join(root, 'public', 'countries.geo.json'), 'utf8'),
-) as FeatureCollection<CountryFeature>;
+) as { type: 'FeatureCollection'; features: CountryFeature[] };
 
-// Precompute areas once (O(total vertices)).
+// Precompute once (O(total vertices)).
+const countryIds = new Set(countries.features.map((c) => c.id));
 const countryArea = new Map<string, number>();
 for (const c of countries.features) countryArea.set(c.id, featureAreaKm2(c.geometry));
 const regionArea = new Map<string, number>();
 for (const f of regions.features) regionArea.set(f.id, featureAreaKm2(f.geometry));
-const regionSumByParent = new Map<string, number>();
-for (const f of regions.features) {
-  regionSumByParent.set(
-    f.properties.parent,
-    (regionSumByParent.get(f.properties.parent) ?? 0) + (regionArea.get(f.id) ?? 0),
-  );
-}
+const sumByParent = new Map<string, number>();
+for (const f of regions.features)
+  sumByParent.set(f.properties.parent, (sumByParent.get(f.properties.parent) ?? 0) + (regionArea.get(f.id) ?? 0));
 const countByParent = new Map<string, number>();
-for (const f of regions.features) {
+for (const f of regions.features)
   countByParent.set(f.properties.parent, (countByParent.get(f.properties.parent) ?? 0) + 1);
-}
 
 describe('regions.geo.json', () => {
-  it('is a valid GeoJSON FeatureCollection', () => {
+  it('is a valid GeoJSON FeatureCollection of well-formed features', () => {
     expect(regions.type).toBe('FeatureCollection');
     expect(Array.isArray(regions.features)).toBe(true);
     expect(regions.features.length).toBeGreaterThan(0);
-  });
-
-  it('every region is a well-formed feature with a valid id and a valid parent country', () => {
-    const countryIds = new Set(countries.features.map((c) => c.id));
-    const idRe = /^[A-Z]{3}-\d+$/;
     for (const f of regions.features) {
       expect(f.type, `${f.id ?? '?'}: not a Feature`).toBe('Feature');
-      expect(idRe.test(f.id), `${f.id ?? '?'}: bad id`).toBe(true);
+      expect(typeof f.id, `${f.id ?? '?'}: id must be a string`).toBe('string');
+      expect(f.id.length, `${f.id ?? '?'}: empty id`).toBeGreaterThan(0);
       expect(typeof f.properties?.name, `${f.id}: missing name`).toBe('string');
       expect(f.properties.name.length, `${f.id}: empty name`).toBeGreaterThan(0);
       expect(typeof f.properties.parent, `${f.id}: missing parent`).toBe('string');
       expect(countryIds.has(f.properties.parent), `${f.id}: unknown parent ${f.properties.parent}`).toBe(true);
-      expect(f.properties.parent, `${f.id}: parent != ISO_A3 prefix`).toBe(f.id.slice(0, 3));
       const g = f.geometry;
       expect(
         !!g && (g.type === 'Polygon' || g.type === 'MultiPolygon'),
         `${f.id}: bad geometry type ${g?.type}`,
       ).toBe(true);
-    }
-  });
-
-  it('gives every country 1..36 regions, and >=2 for any country > 50,000 km²', () => {
-    for (const c of countries.features) {
-      const n = countByParent.get(c.id) ?? 0;
-      expect(n, `country ${c.id} has no regions`).toBeGreaterThanOrEqual(1);
-      expect(n, `country ${c.id} has ${n} regions (> 36)`).toBeLessThanOrEqual(36);
-      if ((countryArea.get(c.id) ?? 0) > 50_000) {
-        expect(n, `large country ${c.id} should have >= 2 regions`).toBeGreaterThanOrEqual(2);
+      if (g) {
+        const polys = g.type === 'Polygon' ? [g.coordinates] : g.coordinates;
+        for (const poly of polys)
+          for (const ring of poly)
+            expect(ring.length, `${f.id}: ring with ${ring.length} pts`).toBeGreaterThanOrEqual(3);
       }
     }
   });
 
-  it('every region has a sane area: >= 200 km² (microstates exempt) and < 5,000,000 km²', () => {
-    let maxArea = 0;
-    let maxId = '';
+  it('has globally unique feature ids (ISO 3166-2 disambiguated on within-country collision)', () => {
+    const seen = new Set<string>();
+    const dups = new Set<string>();
     for (const f of regions.features) {
-      const a = regionArea.get(f.id) ?? 0;
-      const parentArea = countryArea.get(f.properties.parent) ?? 0;
-      // A whole-country singleton of a sub-200 km² microstate is exempt.
-      expect(a >= 200 || parentArea < 200, `${f.id}: area ${a.toFixed(1)} km² is too small`).toBe(true);
-      expect(a < 5_000_000, `${f.id}: area ${Math.round(a).toLocaleString('en-US')} km² is absurdly large`).toBe(true);
-      if (a > maxArea) {
-        maxArea = a;
-        maxId = f.id;
-      }
+      if (seen.has(f.id)) dups.add(f.id);
+      seen.add(f.id);
     }
-    console.log(`\n  max region area: ${maxId} ≈ ${Math.round(maxArea).toLocaleString('en-US')} km²`);
+    expect(
+      dups.size,
+      `duplicate ids: ${[...dups].slice(0, 10).join(', ')}`,
+    ).toBe(0);
+    expect(seen.size).toBe(regions.features.length);
   });
 
-  it('covers each country to within [45%, 108%] of its own area', () => {
-    // The band is a planar-approximation sanity check. The measured range is
-    // ~[64%, 105%]: Antarctica's hemispheric wedges under-fill to ~64%, and
-    // Norway's fjord-heavy, high-latitude coast over-fills to ~105.1% under the
-    // single-mid-latitude approximation. [45%, 108%] brackets that with margin
-    // while still catching a genuinely missing or doubled-up region.
+  it('covers every one of the 177 countries with at least one region', () => {
+    const counts = countries.features.map((c) => countByParent.get(c.id) ?? 0);
+    for (const c of countries.features) {
+      expect(countByParent.get(c.id) ?? 0, `country ${c.id} has no regions`).toBeGreaterThanOrEqual(1);
+    }
+    const min = Math.min(...counts);
+    const max = Math.max(...counts);
+    const med = median(counts);
+    const maxCountry = [...countByParent.entries()].sort((a, b) => b[1] - a[1])[0];
+    const maxCountryName = countries.features.find((c) => c.id === maxCountry[0])?.properties.name ?? '';
+    console.log(
+      `\n  per-country region count: min ${min} / median ${med} / max ${max}` +
+        ` (max: ${maxCountry[0]} ${maxCountryName} with ${maxCountry[1]})`,
+    );
+  });
+
+  it('has sane region areas: median >= 100 km2, and no non-Antarctica region >= 5,000,000 km2', () => {
+    const areas = regions.features.map((f) => regionArea.get(f.id) ?? 0);
+    const med = median(areas);
+    const maxRegion = regions.features.reduce((m, f) => ((regionArea.get(f.id) ?? 0) > (regionArea.get(m.id) ?? 0) ? f : m));
+    const minRegion = regions.features.reduce((m, f) => ((regionArea.get(f.id) ?? 0) < (regionArea.get(m.id) ?? 0) ? f : m));
+    const tiny = areas.filter((a) => a < 1).length;
+    // Antarctica (ATA) is a full-360° wrap ring, so the planar model can't
+    // measure it (it reports ~21.8M km2 for a ~14M km2 landmass) — exempt it.
+    for (const f of regions.features) {
+      if (f.properties.parent === 'ATA') continue;
+      const a = regionArea.get(f.id) ?? 0;
+      expect(a < 5_000_000, `${f.id}: area ${Math.round(a).toLocaleString('en-US')} km2 is absurdly large`).toBe(true);
+    }
+    expect(med, `median region area ${med.toFixed(1)} km2 is too small`).toBeGreaterThanOrEqual(100);
+    console.log(
+      `\n  region area: min ${minRegion.id} (${minRegion.properties.name}) ≈ ${(regionArea.get(minRegion.id) ?? 0).toExponential(2)} km2` +
+        ` / median ${med.toLocaleString('en-US')} km2` +
+        ` / max ${maxRegion.id} (${maxRegion.properties.name}) ≈ ${Math.round(regionArea.get(maxRegion.id) ?? 0).toLocaleString('en-US')} km2` +
+        ` — ${tiny} regions < 1 km2 (legit micro-regions: tiny islands, small urban districts, degenerate disputed enclaves)`,
+    );
+  });
+
+  it('sums to within [40%, 175%] of each country area (median in [70%, 130%])', () => {
+    // The band is a planar-approximation sanity check. It is wider than a naive
+    // [40%, 110%] because Natural Earth's admin-1 polygons are digitized
+    // independently of the admin-0 country polygons (they are NOT clipped to
+    // them), so for a handful of small/island nations the region union
+    // over/under-fills the country: measured range is ~[52%, 151%]
+    // (Brunei under-fills, Vanuatu over-fills). The band brackets that with
+    // margin while still catching a genuinely missing or doubled-up region.
+    const ratios: number[] = [];
     let minRatio = Infinity;
     let maxRatio = -Infinity;
     let minId = '';
@@ -170,22 +204,24 @@ describe('regions.geo.json', () => {
     for (const c of countries.features) {
       const ca = countryArea.get(c.id) ?? 0;
       if (ca <= 0) continue;
-      const ratio = (regionSumByParent.get(c.id) ?? 0) / ca;
-      if (ratio < minRatio) {
-        minRatio = ratio;
-        minId = c.id;
-      }
-      if (ratio > maxRatio) {
-        maxRatio = ratio;
-        maxId = c.id;
-      }
+      const ratio = (sumByParent.get(c.id) ?? 0) / ca;
+      ratios.push(ratio);
+      if (ratio < minRatio) { minRatio = ratio; minId = c.id; }
+      if (ratio > maxRatio) { maxRatio = ratio; maxId = c.id; }
       expect(
-        ratio >= 0.45 && ratio <= 1.08,
-        `${c.id}: coverage ${(100 * ratio).toFixed(1)}% is outside [45%, 108%]`,
+        ratio >= 0.4 && ratio <= 1.75,
+        `${c.id}: coverage ${(100 * ratio).toFixed(1)}% is outside [40%, 175%]`,
       ).toBe(true);
     }
+    const medRatio = median(ratios);
+    expect(
+      medRatio >= 0.7 && medRatio <= 1.3,
+      `median coverage ${(100 * medRatio).toFixed(1)}% is outside [70%, 130%]`,
+    ).toBe(true);
     console.log(
-      `\n  coverage range: min ${minId} ${(100 * minRatio).toFixed(1)}% / max ${maxId} ${(100 * maxRatio).toFixed(1)}% (band 45–108%)`,
+      `\n  coverage range: min ${minId} ${(100 * minRatio).toFixed(1)}%` +
+        ` / median ${(100 * medRatio).toFixed(1)}%` +
+        ` / max ${maxId} ${(100 * maxRatio).toFixed(1)}% (band 40–175%)`,
     );
   });
 });