"""Build the regions file from one source: Overture Maps divisions — mostly OpenStreetMap (ODbL), with parts from
geoBoundaries, Toitū Te Whenua LINZ and Esri Community Maps contributors (CC BY 4.0), as Overture credits them.

This script is also the "method of making the alterations" the ODbL asks a public user of a derived database to
offer: the site serves it at /map-data/build.py.

    docker run --rm -v "$PWD:/w" -w /w python:3.12-slim sh -c "pip install -q duckdb unidecode && python data/regions/build.py"

Writes data/regions/overture/<CC>.geojsonl.gz — one file per country, one GeoJSON feature per line — and
data/regions/overture/manifest.json (release, licence, counts per country):
  level 1: every country (properties: level, country, name)
  level 2: ONE level of local regions below it, worldwide, as fine as the data allows (level, country, name, population)

The rule for level 2, the same in every country: start from the country and replace any area by its children
(regions, counties, local admin areas, localities whose parent it is) wherever children of one kind cover it —
80 to 125 % of its area (simplified coastlines overlap a little), so they neither leave big holes nor overlap; of several kinds that do, the one with the
most members; if no level of children fits, the descendants of one kind further down may (a broken level in
between is skipped) — then do the same inside each child. An area
stays as it is only where nothing finer covers it. So a state without counties (Berlin) becomes its districts, a
county its municipalities, a municipality its complete subdivisions. Children covering at least half of an area
stay with the rest of it as a region under its name. What is left uncovered (forest, mountain, water, slivers) joins
the neighbour with the longest shared border; a hole bigger than that neighbour becomes a region of its own, named
after the county or state it lies in, a far island after a place on it — so the local level has no holes.

Geometries are simplified to ~30 m (countries ~100 m) and coordinates rounded to 4 decimals (~11 m), to keep the
files and the map tiles small.
OVERTURE_LOCAL=<dir> reads area.parquet / division.parquet (and the hoods_* and dep_* parts) from a local copy instead of the release on S3.
ONLY=CC,CC builds just those countries (several such runs side by side make a full build fast); every finished
country is kept in data/regions/.parts/, and ASSEMBLE=1 then writes manifest.json from them.
"""
from __future__ import annotations

import collections
import gzip
import json
import math
import os
import re
import sys
import time
import unicodedata

import duckdb
from unidecode import unidecode

RELEASE = "2026-09-23.1"
LICENCE = "ODbL 1.0 — © OpenStreetMap contributors (ODbL), Overture Maps Foundation; parts from geoBoundaries, Toitū Te Whenua LINZ and Esri Community Maps contributors (CC BY 4.0)"
S3 = f"s3://overturemaps-us-west-2/release/{RELEASE}/theme=divisions"
OUT = "data/regions/overture"
ONLY = os.environ.get("ONLY", "")   # e.g. ONLY=DE,FR: build these countries only (the rest stays as it is)
ASSEMBLE = os.environ.get("ASSEMBLE") == "1"   # write manifest.json from the per-country parts, build nothing
PARTS = os.path.join("data", "regions", ".parts")   # one JSON per finished country: a stop never loses one
CHILD_COVER = (0.80, 1.25)   # above 100 %: simplified coastlines overlap a little; far above means overlapping kinds
KIDS = ("region", "county", "localadmin", "locality")
# quarters, coarsest first: only for a state or county with no towns mapped in it (the emirate of Dubai is mapped as
# its communities), never to split a town — the local level is the municipality, as a German Gemeinde
HOODS = ("macrohood", "neighborhood", "microhood")
# sizes against the German Gemeinde: quarters only for an area too big to be one local region, and only of a kind
# big enough to gather moments (Ireland's townlands, Japan's blocks are not)
BIG_KM2, MIN_QUARTER_KM2 = 1000, 3   # a German Gemeinde: median ~11 km², the smallest tenth under ~3 km²
# names in the source that are vandalism, not places (checked by hand; the next release may have dropped them)
NOT_PLACES = {"Skibid Jonkler Town"}
# outlines in the source that take in the sea: Macao's is 116 km² for 33 km² of land — its divisions stand alone
BAD_OUTLINE = {"MO"}
# Overture keeps disputed land apart under X codes: Kosovo is a country of its own; the West Bank and Gaza are Palestine
COUNTRY_CODE = {"XK": "XK", "XW": "PS", "XG": "PS", "XH": "SY"}
PARTS_OF = {"PS": {"XW": "West Bank", "XG": "Gaza Strip"}}   # a country made of parts: each part refined on its own
# disputed land joined to the country it belongs to under international law: the Golan Heights (Overture: XH) are
# Syria's — its area becomes part of Syria's, and the land hangs under it as one more region, divided by its own data
JOINED = {"XH": ("SY", "Golan Heights")}
COUNTRY_NAME = {"XK": "Kosovo", "PS": "Palestine"}
# territories Overture maps under their own code but without a country area of their own (Greenland, Hong Kong,
# Puerto Rico …): each is a country here, as in ISO 3166, its area its outline and its divisions together
TERRITORY = {
    "AS": "American Samoa", "AW": "Aruba", "AX": "Åland Islands", "BL": "Saint Barthélemy", "BM": "Bermuda",
    "CC": "Cocos (Keeling) Islands", "CK": "Cook Islands", "CP": "Clipperton Island", "CW": "Curaçao",
    "CX": "Christmas Island", "FK": "Falkland Islands", "FO": "Faroe Islands", "GF": "French Guiana", "GG": "Guernsey",
    "GL": "Greenland", "GP": "Guadeloupe", "GS": "South Georgia and the South Sandwich Islands", "GU": "Guam",
    "HK": "Hong Kong", "IM": "Isle of Man", "JE": "Jersey", "KY": "Cayman Islands", "MO": "Macao",
    "MP": "Northern Mariana Islands", "MQ": "Martinique", "MS": "Montserrat", "NC": "New Caledonia", "NU": "Niue",
    "PF": "French Polynesia", "PM": "Saint Pierre and Miquelon", "PR": "Puerto Rico", "RE": "Réunion",
    "SH": "Saint Helena, Ascension and Tristan da Cunha", "SJ": "Svalbard and Jan Mayen", "TC": "Turks and Caicos Islands",
    "TF": "French Southern Territories", "TK": "Tokelau", "UM": "United States Minor Outlying Islands",
    "VG": "British Virgin Islands", "VI": "United States Virgin Islands", "WF": "Wallis and Futuna", "YT": "Mayotte"}
COUNTRY_NAME.update(TERRITORY)
# the languages whose name stands in, in this order, for one written only in a script that transliterates to nothing
LATIN_LANGS = ("fr", "es", "pt", "de", "it", "nl", "tr", "pl", "id", "ms", "vi", "ro", "cs", "sk", "hu", "hr", "sl", "sq", "lt",
               "lv", "et", "fi", "sv", "da", "no", "az", "uz", "tk", "ku", "so", "sw", "ha", "tl", "mt", "ga", "cy", "eu", "ca",
               "gl", "af")
BIG_TOWN = 100_000   # a town this big, mapped as an area of its own, is never lost inside a bigger region (carve_towns)
PLACEHOLDERS = ("null", "none", "unknown", "unnamed", "n/a", "-")   # names in the data that name nothing
MIN_REGION_KM2 = 0.01   # a local region under a hectare joins its neighbour (drop_tiny)
REMNANT_MIN = 0.10      # a region left with less than this share of itself after the overlaps are cut is dropped   # a local region under a hectare (a wharf, a farmyard filed as a village) joins its neighbour
OUTLINE_REST, OUTLINE_PIECE = 0.05, 0.01   # a territory's land no division claims: from 5 % of it, its pieces from 1 % are regions
PART_MIN = 0.50    # children covering under half of an area are scattered places (village blobs), not a division of it
LATIN = re.compile(r"[A-Za-z]")
def _other_script() -> re.Pattern:
    """Any letter of a script other than Latin (Arabic, Cyrillic, Myanmar, Javanese, Tibetan, syllabics …);
    Latin letters with marks (ə, ł, ő) and modifier letters (the ʻokina, ᵉ) count as Latin."""
    ranges, start = [], None
    for cp in range(0x30000 + 1):
        c = chr(cp) if cp <= 0x2FFFF else "a"
        n = unicodedata.name(c, "")
        other = c.isalpha() and not n.startswith(("LATIN", "MODIFIER LETTER", "MASCULINE ORDINAL", "FEMININE ORDINAL"))
        if other and start is None:
            start = cp
        elif not other and start is not None:
            ranges.append(f"\\U{start:08x}-\\U{cp - 1:08x}")
            start = None
    return re.compile("[" + "".join(ranges) + "]")


OTHER_SCRIPT = _other_script()
# scripts whose transliteration reads as nothing (Arabic and Persian lose their vowels, Thai its word breaks):
# a region named only in one of them takes the name of a town inside it instead
POOR_SCRIPT = re.compile("[" + "".join(
    f"\\U{cp:08x}" for cp in range(0x2FFFF) if chr(cp).isalpha()
    and unicodedata.name(chr(cp), "").startswith(("ARABIC", "HEBREW", "SYRIAC", "THAANA", "NKO", "THAI", "LAO", "KHMER",
                                                 "MYANMAR", "TIFINAGH"))) + "]")
WARD = re.compile(r"(?i)\bward\s*\d+")   # an electoral ward named by its number is no place
TOL_LOCAL, TOL_COUNTRY = 0.0003, 0.001


def tidy(name: str) -> str:
    """The marks a source or a transliteration leaves in a name: underscores, stray symbols, SHOUTED or s p a c e d
    letters, "de Nako Department"."""
    name = re.sub(r"[`@~\[\]{}|^]", "", name.replace("_", " "))
    name = re.sub(r"\s+", " ", name).strip(" -–,")
    if re.fullmatch(r"(?:[A-Z] )+[A-Z]", name):
        name = name.replace(" ", "")
    letters = [c for c in name if c.isalpha()]
    if len(letters) >= 4 and all(c.isupper() for c in letters):
        name = " ".join(w.capitalize() if w.isalpha() or "-" in w or "'" in w else w for w in name.split())
    return re.sub(r"^(?:de |d')(.+ Department)$", r"\1", name)


def connect() -> duckdb.DuckDBPyConnection:
    con = duckdb.connect()
    con.execute(f"SET memory_limit='{os.environ.get('BUILD_MEMORY', '8GB')}'; SET threads=4; SET preserve_insertion_order=false; INSTALL spatial; LOAD spatial;")
    local = os.environ.get("OVERTURE_LOCAL")
    if local:
        con.execute(f"CREATE VIEW area AS SELECT * FROM read_parquet(['{local}/area.parquet', '{local}/hoods_area.parquet', '{local}/dep_area.parquet'])")
        con.execute(f"CREATE TABLE division AS SELECT * FROM read_parquet(['{local}/division.parquet', '{local}/hoods_division.parquet', '{local}/dep_division.parquet'])")
        con.execute(f"CREATE TABLE country_en AS SELECT * FROM read_parquet('{local}/country_en.parquet')")
        con.execute(f"CREATE TABLE division_en AS SELECT * FROM read_parquet(['{local}/division_en.parquet', '{local}/hoods_en.parquet', '{local}/dep_en.parquet'])")
        con.execute(f"CREATE TABLE division_latin AS SELECT * FROM read_parquet('{local}/division_latin.parquet')")
        return con
    con.execute("INSTALL httpfs; LOAD httpfs; SET s3_region='us-west-2';")
    con.execute(f"""CREATE TABLE division AS SELECT id, parent_division_id AS parent, subtype, country, names.primary AS name, population
        FROM read_parquet('{S3}/type=division/*', hive_partitioning=1) WHERE subtype IN ('country','dependency','region','county','localadmin','locality','macrohood','neighborhood','microhood')""")
    con.execute(f"""CREATE TABLE country_en AS SELECT country, names.common['en'] AS name_en
        FROM read_parquet('{S3}/type=division/*', hive_partitioning=1) WHERE subtype = 'country'""")
    con.execute(f"""CREATE TABLE division_en AS SELECT id, names.common['en'] AS name_en
        FROM read_parquet('{S3}/type=division/*', hive_partitioning=1)
        WHERE subtype IN ('dependency','region','county','localadmin','locality','macrohood','neighborhood','microhood') AND names.common['en'] IS NOT NULL""")
    # no English name: the first one in Latin letters in another language (Tunisia's sectors in French), for a name
    # whose own script transliterates to nothing
    latin = "list_filter([" + ", ".join(f"names.common['{l}']" for l in LATIN_LANGS) + "], x -> x IS NOT NULL AND " \
            "regexp_matches(x, '^[\\p{Latin}\\p{M}0-9 .,''’()/&-]+$') AND regexp_matches(x, '\\p{L}'))[1]"
    con.execute(f"""CREATE TABLE division_latin AS SELECT id, {latin} AS name_latin
        FROM read_parquet('{S3}/type=division/*', hive_partitioning=1)
        WHERE subtype IN ('country','dependency','region','county','localadmin','locality','macrohood','neighborhood','microhood')
          AND names.common['en'] IS NULL AND names.common IS NOT NULL AND {latin} IS NOT NULL""")
    con.execute(f"""CREATE TABLE area_t AS SELECT division_id, country, subtype,
        ST_AsWKB(ST_SimplifyPreserveTopology(geometry, CASE WHEN subtype IN ('country','dependency') THEN {TOL_COUNTRY} ELSE {TOL_LOCAL} END)) AS wkb
        FROM read_parquet('{S3}/type=division_area/*', hive_partitioning=1)
        WHERE class='land' AND subtype IN ('country','dependency','region','county','localadmin','locality','macrohood','neighborhood','microhood')""")
    con.execute("CREATE VIEW area AS SELECT * FROM area_t")
    return con


def fill_gaps(con: duckdb.DuckDBPyConnection, country_wkb: bytes, leaves: list, namers: list, island_name=None) -> tuple[dict[str, bytes], list]:
    """fill_gaps_once, and when the geometry library trips over an invalid source geometry (a TopologyException),
    once more with every input snapped to a 1 m grid, which repairs it; failing that, the country keeps its
    local regions as they are and the build goes on."""
    for grid in (0.0, 0.00001):
        try:
            return fill_gaps_once(con, country_wkb, leaves, namers, grid, island_name)
        except duckdb.Error as exc:
            print(f"  gap fill {'failed' if grid else 'retrying on a 1 m grid'}: {str(exc)[:90]}", file=sys.stderr)
    return {}, []


def fill_gaps_once(con: duckdb.DuckDBPyConnection, country_wkb: bytes, leaves: list, namers: list, grid: float, island_name=None) -> tuple[dict[str, bytes], list]:
    """No hole in the local level: every part of the country no local region covers (a lake, a forest, land no
    municipality claims, the slivers simplification leaves between neighbours) is joined to the local region it
    shares the longest border with — or, for an island of its own, the nearest one. Returns id → new geometry
    for the regions that grew."""
    con.execute("CREATE OR REPLACE TEMP TABLE lv (id VARCHAR, g GEOMETRY)")
    snap = (lambda g: f"ST_MakeValid(ST_ReducePrecision({g}, {grid}))") if grid else (lambda g: f"ST_MakeValid({g})")
    con.executemany(f"INSERT INTO lv VALUES (?, {snap('ST_GeomFromWKB(?)')})", [(l[0], l[6]) for l in leaves])
    con.execute("""CREATE OR REPLACE TEMP TABLE gp AS
        SELECT row_number() OVER () AS gid, d.geom AS g FROM (
            SELECT unnest(ST_Dump(ST_Difference(""" + snap("ST_GeomFromWKB(?)") + """, (SELECT ST_Union_Agg(g) FROM lv)))) AS d)
        WHERE ST_Area(d.geom) > least(1e-6, 0.001 * ST_Area(ST_GeomFromWKB(?)))""", [country_wkb, country_wkb])   # about a hectare (less in a
    # country of a few km², Monaco): smaller slivers are invisible on the map
    if not con.execute("SELECT count(*) FROM gp").fetchone()[0]:
        return {}, []
    # the neighbour: the region the gap's thin halo overlaps most (≈ the longest shared border); failing any, the
    # nearest region within half a degree
    con.execute("""CREATE OR REPLACE TEMP TABLE pick AS
        WITH halo AS (SELECT gid, ST_Buffer(g, 0.002) AS h FROM gp),   -- buffer the small gap, never the big neighbour
        touch AS (
            SELECT halo.gid, lv.id, ST_Area(ST_Intersection(halo.h, lv.g)) AS shared
            FROM halo JOIN lv ON ST_Intersects(halo.h, lv.g)),
        best AS (SELECT gid, arg_max(id, shared) AS id FROM touch GROUP BY gid),
        lone AS (
            SELECT gp.gid, arg_min(lv.id, ST_Distance(gp.g, lv.g)) AS id
            FROM gp JOIN lv ON ST_DWithin(gp.g, lv.g, 0.5) WHERE gp.gid NOT IN (SELECT gid FROM best) GROUP BY gp.gid)
        SELECT * FROM best UNION ALL SELECT * FROM lone""")
    # a hole bigger than the region that would take it is not a gap but land of its own: it becomes a region named
    # after the state or county it lies in (never the country), instead of swelling its neighbour many times over
    extra = []
    big = con.execute("""SELECT gp.gid, ST_AsWKB(gp.g), ST_Area(gp.g) FROM pick JOIN gp USING (gid) JOIN lv USING (id)
                         WHERE ST_Area(gp.g) > ST_Area(lv.g)""").fetchall()
    by_namer: dict[int, list] = {}
    if big and namers:
        con.execute("CREATE OR REPLACE TEMP TABLE nm (i INTEGER, g GEOMETRY)")
        con.executemany("INSERT INTO nm VALUES (?, ST_MakeValid(ST_GeomFromWKB(?)))", [(i, n[6]) for i, n in enumerate(namers)])
        for gid, wkb, area in big:
            # the smallest division holding at least half of the hole (the county, not the state it lies in);
            # failing that, the one holding most of it
            hit = con.execute("""SELECT i FROM nm WHERE ST_Intersects(g, ST_GeomFromWKB(?))
                                 ORDER BY (ST_Area(ST_Intersection(g, ST_GeomFromWKB(?))) >= 0.5 * ?) DESC, ST_Area(g) ASC LIMIT 1""",
                              [wkb, wkb, area or 1e-12]).fetchone()
            if hit:
                by_namer.setdefault(hit[0], []).append(bytes(wkb))
                con.execute("DELETE FROM pick WHERE gid = ?", [gid])
    for i, pieces in by_namer.items():   # all the big pieces of one area become ONE region under its name
        n = namers[i]
        w, a = con.execute("SELECT ST_AsWKB(g), ST_Area(g) FROM (SELECT ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(x))) AS g "
                           "FROM (SELECT unnest(?::BLOB[]) AS x))", [pieces]).fetchone()
        extra.append((f"{n[0]}#rest", n[1], n[2], n[3], None, a, bytes(w)))
    # land with no region anywhere near (a far island): a region of its own, named after a place on it
    if island_name:
        islands: dict[str, list] = {}
        for gid, wkb, area in con.execute("""SELECT gid, ST_AsWKB(g), ST_Area(g) FROM gp
                                             WHERE gid NOT IN (SELECT gid FROM pick)""").fetchall():
            if any(gid == g for g, _, _ in big):
                continue
            name = island_name(bytes(wkb))
            if name:
                islands.setdefault(name, []).append(bytes(wkb))
        for k, (name, pieces) in enumerate(islands.items()):   # the islets around one place are one region
            w, a = con.execute("SELECT ST_AsWKB(g), ST_Area(g) FROM (SELECT ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(x))) AS g "
                               "FROM (SELECT unnest(?::BLOB[]) AS x))", [pieces]).fetchone()
            extra.append((f"island#{k}", "locality", None, name, None, a, bytes(w)))
    rows = con.execute("""
        SELECT lv.id, ST_AsWKB(ST_MakeValid(ST_Union_Agg(x.g))) FROM (
            SELECT id, g FROM lv WHERE id IN (SELECT id FROM pick)
            UNION ALL SELECT pick.id, gp.g FROM pick JOIN gp USING (gid)) x JOIN lv USING (id) GROUP BY lv.id""").fetchall()
    return {i: bytes(w) for i, w in rows}, extra


def cut_overlaps(con: duckdb.DuckDBPyConnection, leaves: list) -> list:
    """No two local regions on top of each other: where the data lets two overlap (quarters filed twice, a town
    under one area and the municipality around it under another), the smaller keeps the overlap and is cut out of
    the bigger — the same one a point there belongs to (regions_at). A region left with nothing is dropped."""
    if len(leaves) < 2:
        return leaves
    con.execute("CREATE OR REPLACE TEMP TABLE lf (id VARCHAR, a DOUBLE, g GEOMETRY)")
    con.executemany("INSERT INTO lf VALUES (?, ?, ST_MakeValid(ST_GeomFromWKB(?)))", [(l[0], l[5] or 0, l[6]) for l in leaves])
    rows = con.execute("""
        WITH pairs AS (
            SELECT b.id, s.g FROM lf b JOIN lf s ON ST_Intersects(b.g, s.g)
            WHERE s.id <> b.id AND (s.a < b.a OR (s.a = b.a AND s.id < b.id))
              AND ST_Area(ST_Intersection(b.g, s.g)) > 1e-7)
        SELECT lf.id, ST_AsWKB(ST_MakeValid(ST_Difference(lf.g, ST_Union_Agg(pairs.g)))) AS w
        FROM lf JOIN pairs USING (id) GROUP BY lf.id, lf.g""").fetchall()
    if not rows:
        return leaves
    cut = {i: bytes(w) for i, w in rows}
    out = []
    for l in leaves:
        if l[0] not in cut:
            out.append(l)
            continue
        a = con.execute("SELECT ST_Area(ST_CollectionExtract(ST_GeomFromWKB(?), 3))", [cut[l[0]]]).fetchone()[0] or 0
        # a region the others cover almost whole (a community filed over its neighbours) is no region of its own:
        # its remnant joins a neighbour in the gap fill
        if a > 1e-7 and a >= REMNANT_MIN * (l[5] or 0):
            w = con.execute("SELECT ST_AsWKB(ST_CollectionExtract(ST_GeomFromWKB(?), 3))", [cut[l[0]]]).fetchone()[0]
            out.append((*l[:5], a, bytes(w)))
    return out


def carve_towns(con: duckdb.DuckDBPyConnection, leaves: list, rows: list) -> list:
    """A town of BIG_TOWN people or more that the data maps as an area of its own, but that lies inside a bigger local
    region of another name (Bouaké in the Gbêkê Region, Espoo in the Helsinki sub-region), becomes a local region of
    its own; cut_overlaps then takes it out of the region around it."""
    have = {l[0] for l in leaves}
    towns = [r for r in rows if r[1] == "locality" and (r[4] or 0) >= BIG_TOWN and r[0] not in have and r[6]]
    if not towns or not leaves:
        return leaves
    con.execute("CREATE OR REPLACE TEMP TABLE ct (id VARCHAR, name VARCHAR, a DOUBLE, g GEOMETRY)")
    con.executemany("INSERT INTO ct VALUES (?, ?, ?, ST_MakeValid(ST_GeomFromWKB(?)))", [(t[0], t[3], t[5] or 0, t[6]) for t in towns])
    con.execute("CREATE OR REPLACE TEMP TABLE cl (id VARCHAR, name VARCHAR, a DOUBLE, g GEOMETRY)")
    con.executemany("INSERT INTO cl VALUES (?, ?, ?, ST_MakeValid(ST_GeomFromWKB(?)))", [(l[0], l[3], l[5] or 0, l[6]) for l in leaves])
    # the region holding most of the town: the town lies at least 80 % in it and is under half its size
    hits = con.execute("""SELECT ct.id, arg_max(cl.name, ST_Area(ST_Intersection(ct.g, cl.g))), max(ST_Area(ST_Intersection(ct.g, cl.g))) / max(ct.a),
                                 arg_max(cl.a, ST_Area(ST_Intersection(ct.g, cl.g))), max(ct.a)
                          FROM ct JOIN cl ON ST_Intersects(ct.g, cl.g) GROUP BY ct.id""").fetchall()
    by_id = {t[0]: t for t in towns}
    # also inside a region of its own name (Manaus in the municipality of Manaus): the town is the town, the rest its
    # surroundings
    add = [by_id[i] for i, around, share, around_a, a in hits if share >= 0.8 and a < 0.5 * (around_a or 0)]
    return leaves + add


def drop_tiny(con: duckdb.DuckDBPyConnection, leaves: list) -> list:
    """No local region under MIN_REGION_KM2 (a wharf, a farmyard filed as a village): it joins its neighbour in the
    gap fill. A country made only of such keeps its biggest."""
    if len(leaves) < 2:
        return leaves
    con.execute("CREATE OR REPLACE TEMP TABLE td (id VARCHAR, g GEOMETRY)")
    con.executemany("INSERT INTO td VALUES (?, ST_GeomFromWKB(?))", [(l[0], l[6]) for l in leaves])
    km2 = dict(con.execute("SELECT id, ST_Area(g) * 12392 * greatest(cos(radians(ST_Y(ST_Centroid(g)))), 0.01) FROM td").fetchall())
    kept = [l for l in leaves if (km2.get(l[0]) or 0) >= MIN_REGION_KM2]
    return kept or [max(leaves, key=lambda l: km2.get(l[0]) or 0)]


def assemble() -> None:
    """manifest.json from the parts the builders left (several may run side by side on different countries)."""
    per_country = {f[:-5]: json.load(open(os.path.join(PARTS, f), encoding="utf-8"))
                   for f in sorted(os.listdir(PARTS)) if f.endswith(".json")}
    with open(os.path.join(OUT, "manifest.json"), "w", encoding="utf-8") as fh:
        json.dump({"source": "Overture Maps Foundation, divisions theme", "release": RELEASE, "licence": LICENCE,
                   "countries": len(per_country), "local_regions": sum(v["local_regions"] for v in per_country.values()),
                   "per_country": per_country}, fh, ensure_ascii=False, indent=1, sort_keys=True)
    print(f"manifest: {len(per_country)} countries, {sum(v['local_regions'] for v in per_country.values()):,} local regions", file=sys.stderr)


def build() -> None:
    if ASSEMBLE:
        return assemble()
    con = connect()
    t0 = time.time()
    # one row per division with an area: its parent, name, population, area (deg²) and geometry
    con.execute("""CREATE TABLE a AS
        SELECT ar.division_id AS id, CASE ar.country {xcase} ELSE ar.country END AS country, ar.subtype, d.parent, d.name,
               d.population, ST_Area(ST_GeomFromWKB(ar.wkb)) AS area, ar.wkb
        FROM area ar JOIN division d ON d.id = ar.division_id
        WHERE ar.country IS NOT NULL AND length(ar.country) = 2 AND (NOT ar.country LIKE 'X%' OR ar.country IN ({xin}))""".format(
        xcase=" ".join(f"WHEN '{x}' THEN '{c}'" for x, c in COUNTRY_CODE.items()), xin=", ".join(f"'{x}'" for x in COUNTRY_CODE)))
    # a territory without a country area: one is made of its top divisions (those whose parent lies in another
    # country), which then hang under it, and of its outline where Overture has one (a dependency: Greenland's
    # national park belongs to no municipality)
    for cc in TERRITORY:
        if con.execute("SELECT count(*) FROM a WHERE country = ? AND subtype = 'country'", [cc]).fetchone()[0]:
            continue
        top = [i for (i,) in con.execute("SELECT id FROM a WHERE country = ? AND subtype <> 'dependency' AND (parent IS NULL "
                                         "OR parent NOT IN (SELECT id FROM a WHERE country = ? AND subtype <> 'dependency'))",
                                         [cc, cc]).fetchall()]
        outline = [i for (i,) in con.execute("SELECT id FROM a WHERE country = ? AND subtype = 'dependency'", [cc]).fetchall()]
        if not top:
            continue
        cid = f"country:{cc}"
        marks = ", ".join("?" * len(top))
        con.execute(f"INSERT INTO a SELECT ?, ?, 'country', NULL, ?, NULL, ST_Area(g), ST_AsWKB(g) FROM "
                    f"(SELECT ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(wkb))) AS g FROM a WHERE id IN ({marks}{', ?' * len(outline)}))",
                    [cid, cc, TERRITORY[cc], *top, *outline])
        # land of the outline that no division claims (Greenland's national park, Nordaustlandet, most of Aruba): where
        # it is at least OUTLINE_REST of the territory, its pieces of at least OUTLINE_PIECE are regions of their own
        # among the top divisions, one per side it lies on — "Northeast Greenland"; smaller slivers go to the gap fill
        pieces = con.execute(f"""WITH t AS (SELECT ST_GeomFromWKB(wkb) AS t FROM a WHERE id = ?),
                                     p AS (SELECT unnest(ST_Dump(ST_MakeValid(ST_Difference((SELECT t FROM t),
                                           (SELECT ST_Union_Agg(ST_GeomFromWKB(wkb)) FROM a WHERE id IN ({marks})))))).geom AS g)
                                 SELECT ST_AsWKB(g), ST_X(ST_Centroid(g)) - ST_X(ST_Centroid(t)), ST_Y(ST_Centroid(g)) - ST_Y(ST_Centroid(t))
                                 FROM p, t WHERE ST_Area(g) >= ? * ST_Area(t)""", [cid, *top, OUTLINE_PIECE]).fetchall() \
            if outline and cc not in BAD_OUTLINE else []
        if pieces and con.execute("SELECT sum(ST_Area(ST_GeomFromWKB(x))) FROM (SELECT unnest(?::BLOB[]) AS x)",
                                  [[bytes(p[0]) for p in pieces]]).fetchone()[0] < OUTLINE_REST * con.execute(
                                  "SELECT area FROM a WHERE id = ?", [cid]).fetchone()[0]:
            pieces = []   # a few slivers, no land of its own
        if outline and cc in BAD_OUTLINE:   # its divisions alone, as before
            con.execute(f"UPDATE a SET wkb = (SELECT ST_AsWKB(ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(wkb)))) FROM a WHERE id IN ({marks})) "
                        "WHERE id = ?", [*top, cid])
            con.execute("UPDATE a SET area = ST_Area(ST_GeomFromWKB(wkb)) WHERE id = ?", [cid])
        own = con.execute("SELECT name_en FROM division_en WHERE id = ?", outline[:1]).fetchone() if outline else None
        by_side: dict[str, list] = {}
        for wkb, dx, dy in pieces:
            side = ("East", "Northeast", "North", "Northwest", "West", "Southwest", "South", "Southeast")[
                round(math.degrees(math.atan2(dy, dx)) / 45) % 8]
            by_side.setdefault(side, []).append(bytes(wkb))
        for side, parts in by_side.items():
            rid, rname = f"outline:{cc}:{side}", f"{side} {own[0] if own and own[0] else TERRITORY[cc]}"
            con.execute("INSERT INTO a SELECT ?, ?, 'region', ?, ?, NULL, ST_Area(g), ST_AsWKB(g) FROM (SELECT ST_MakeValid("
                        "ST_Union_Agg(ST_GeomFromWKB(x))) AS g FROM (SELECT unnest(?::BLOB[]) AS x))", [rid, cc, cid, rname, parts])
            con.execute("INSERT INTO division (id, parent, subtype, country, name, population) VALUES (?, ?, 'region', ?, ?, NULL)",
                        [rid, cid, cc, rname])
        con.execute(f"UPDATE a SET parent = ? WHERE id IN ({marks})", [cid, *top])
        con.execute("INSERT INTO division (id, parent, subtype, country, name, population) VALUES (?, NULL, 'country', ?, ?, NULL)", [cid, cc, TERRITORY[cc]])
        con.execute(f"UPDATE division SET parent = ? WHERE id IN ({marks})", [cid, *top])
    con.execute("DELETE FROM a WHERE subtype = 'dependency'")   # an outline is only a territory's shape, never a region
    # a country Overture keeps in two parts (the West Bank and Gaza): one country made of their union, the parts
    # under it as its first areas, so each is divided as far as its own data goes
    for cc, named_parts in PARTS_OF.items():
        cid = f"country:{cc}"
        parts = {x: i for x, i in con.execute("SELECT d.country, a.id FROM a JOIN division d USING (id) "
                                              "WHERE a.country = ? AND a.subtype = 'country'", [cc]).fetchall()}
        if not parts:
            continue
        marks = ", ".join("?" * len(parts))
        con.execute(f"INSERT INTO a SELECT ?, ?, 'country', NULL, ?, NULL, ST_Area(g), ST_AsWKB(g) FROM "
                    f"(SELECT ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(wkb))) AS g FROM a WHERE id IN ({marks}))",
                    [cid, cc, COUNTRY_NAME[cc], *parts.values()])
        con.execute("INSERT INTO division (id, parent, subtype, country, name, population) VALUES (?, NULL, 'country', ?, ?, NULL)",
                    [cid, cc, COUNTRY_NAME[cc]])
        for x, i in parts.items():
            con.execute("UPDATE a SET subtype = 'region', parent = ?, name = ? WHERE id = ?", [cid, named_parts.get(x, x), i])
            con.execute("UPDATE division SET subtype = 'region', parent = ?, name = ? WHERE id = ?", [cid, named_parts.get(x, x), i])
    for x, (cc, part_name) in JOINED.items():
        part, whole = (con.execute("SELECT a.id FROM a JOIN division d USING (id) WHERE d.country = ? AND a.subtype = 'country' "
                                   "ORDER BY a.area DESC LIMIT 1", [c]).fetchone() for c in (x, cc))
        if not part or not whole:
            continue
        wkb, area = con.execute("SELECT ST_AsWKB(g), ST_Area(g) FROM (SELECT ST_MakeValid(ST_Union_Agg(ST_GeomFromWKB(wkb))) AS g "
                                "FROM a WHERE id IN (?, ?))", [whole[0], part[0]]).fetchone()
        con.execute("UPDATE a SET wkb = ?, area = ? WHERE id = ?", [bytes(wkb), area, whole[0]])
        con.execute("UPDATE a SET subtype = 'region', parent = ?, name = ? WHERE id = ?", [whole[0], part_name, part[0]])
        con.execute("UPDATE division SET subtype = 'region', parent = ?, name = ? WHERE id = ?", [whole[0], part_name, part[0]])
    countries =[r[0] for r in con.execute("SELECT DISTINCT country FROM a WHERE subtype = 'country' ORDER BY 1").fetchall()]
    names = dict(con.execute("SELECT id, name FROM division").fetchall())
    parent_of = dict(con.execute("SELECT id, parent FROM division WHERE parent IS NOT NULL").fetchall())
    country_ids = {i for (i,) in con.execute("SELECT id FROM division WHERE subtype = 'country'").fetchall()}
    # countries by their English name, for a board read worldwide; local regions keep their own official name
    english = {c: n for c, n in con.execute("SELECT country, name_en FROM country_en WHERE name_en IS NOT NULL").fetchall()}
    english.update(COUNTRY_NAME)
    english_local = dict(con.execute("SELECT id, name_en FROM division_en").fetchall())
    latin_local = dict(con.execute("SELECT id, name_latin FROM division_latin").fetchall())
    # a division without an English name borrows the one of a namesake in the same country that has it
    # (the county "München" has none, the city "München" is "Munich")
    english_alike = {(c, n): e for c, n, e in con.execute(
        "SELECT d.country, d.name, arg_max(e.name_en, coalesce(d.population, 0)) FROM division d JOIN division_en e USING (id) "
        "WHERE d.name IS NOT NULL GROUP BY 1, 2").fetchall()}
    country_of = dict(con.execute("SELECT id, country FROM division").fetchall())

    def source_name(division_id: str | None, primary: str | None) -> str | None:
        """The English name where the data has one (Munich, Al Khor), else the name itself — or, for a name only in a
        script that transliterates to nothing, its name in Latin letters in another language (Tunisia: French)."""
        division_id = division_id.split("#")[0] if division_id else division_id   # a rest or a hole of r is named like r
        name = ((english_local.get(division_id) or english_alike.get((country_of.get(division_id), primary)))
                if division_id else None) or primary
        if division_id and name and POOR_SCRIPT.search(name) and latin_local.get(division_id):
            return latin_local[division_id]
        return name

    def unreadable(division_id: str | None, primary: str | None) -> bool:
        """A name only in a script whose transliteration reads as nothing ("dhstn ql`h qfh")."""
        name = source_name(division_id, primary)
        if not name or not POOR_SCRIPT.search(name):
            return False
        return not any(LATIN.search(t) or (OTHER_SCRIPT.search(t) and not POOR_SCRIPT.search(t)) for t in name.split())

    def local_name(division_id: str | None, primary: str | None) -> str | None:
        """Every name in English, for a board read worldwide: the English name where the data has one (Munich,
        Al Khor); else the name itself if it is in Latin letters; else its Latin part (a name that strings several
        scripts together); else the name transliterated into Latin letters — leaving out a part in a script that
        transliterates to nothing when another part reads."""
        name = source_name(division_id, primary)
        if not name or not OTHER_SCRIPT.search(name):
            return tidy(name) if name else name
        latin = re.sub(r"\s*\([^)]*\)", lambda m: "" if OTHER_SCRIPT.search(m.group(0)) else m.group(0), name)
        latin = " ".join(t for t in latin.split() if not OTHER_SCRIPT.search(t)).strip(" -–,")
        if LATIN.search(latin):
            return tidy(latin)
        readable = " ".join(t for t in name.split() if not POOR_SCRIPT.search(t))
        readable = readable if re.search(r"[^\W\d_]", readable) else name   # a number left alone is no name
        t = tidy(unidecode(readable).strip()) or name
        return t[:1].upper() + t[1:]

    place_cache: dict[str, list] = {}

    def places(cc: str) -> list:
        """Named places of a country with their position — only fetched for a country with numbered regions."""
        if cc not in place_cache:
            local = os.environ.get("OVERTURE_LOCAL")
            codes = ", ".join(f"'{c}'" for c in [cc] + [x for x, c in COUNTRY_CODE.items() if c == cc and x != cc])
            src = f"(SELECT * FROM read_parquet('{local}/places_all.parquet') WHERE country IN ({codes}))" if local else f"""(SELECT country, subtype,
                coalesce(names.common['en'], names.primary) AS name, population, cartography.prominence AS prominence,
                ST_AsWKB(geometry) AS wkb FROM read_parquet('{S3}/type=division/*', hive_partitioning=1)
                WHERE country IN ({codes}) AND subtype IN ('locality','neighborhood','macrohood','microhood','localadmin','county','region'))"""
            place_cache[cc] = con.execute(f"SELECT subtype, name, coalesce(prominence, 0), coalesce(population, 0), wkb FROM {src} "
                                          f"WHERE name IS NOT NULL AND lower(name) NOT IN {PLACEHOLDERS}").fetchall()   # a placeholder is no name
            place_cache[cc] = [p for p in place_cache[cc] if p[1] not in NOT_PLACES]
        return place_cache[cc]

    SPECIFIC = {"macrohood": 0, "neighborhood": 1, "microhood": 2, "locality": 3, "localadmin": 4, "county": 5, "region": 6}

    def place_inside(cc: str, leaf) -> str | None:
        """The most specific named place inside a region, the most prominent first."""
        inside = [p for p in places(cc)
                  if con.execute("SELECT ST_Contains(ST_GeomFromWKB(?), ST_GeomFromWKB(?))", [leaf[6], p[4]]).fetchone()[0]]
        if inside:
            return local_name(None, min(inside, key=lambda p: (SPECIFIC.get(p[0], 9), -p[2], -p[3]))[1])
        return None

    def towns_inside(cc: str, leaves: list) -> dict[str, str]:
        """For regions named only in a script that transliterates to nothing (Iran's rural districts, Tunisia's
        sectors): the town inside each with a name in Latin letters, the most prominent first — all in one query."""
        if not leaves:
            return {}
        latin = [(p[0], p[1], p[2], p[3], p[4]) for p in places(cc) if not OTHER_SCRIPT.search(p[1]) and LATIN.search(p[1])]
        if not latin:
            return {}
        con.execute("CREATE OR REPLACE TEMP TABLE tw (subtype VARCHAR, name VARCHAR, prominence INTEGER, population BIGINT, g GEOMETRY)")
        con.executemany("INSERT INTO tw VALUES (?, ?, ?, ?, ST_GeomFromWKB(?))", latin)
        con.execute("CREATE OR REPLACE TEMP TABLE tl (id VARCHAR, g GEOMETRY)")
        con.executemany("INSERT INTO tl VALUES (?, ST_MakeValid(ST_GeomFromWKB(?)))", [(l[0], l[6]) for l in leaves])
        rows = con.execute("""SELECT tl.id, arg_max(tw.name, (CASE WHEN tw.subtype IN ('locality', 'localadmin') THEN 1e12 ELSE 0 END)
                                                         + tw.prominence * 1e9 + least(tw.population, 999999999))
                              FROM tl JOIN tw ON ST_Contains(tl.g, tw.g) GROUP BY tl.id""").fetchall()
        # none inside: the nearest town within ~30 km
        rows += con.execute("""SELECT tl.id, arg_min(tw.name, ST_Distance(tl.g, tw.g)) FROM tl JOIN tw ON ST_DWithin(tl.g, tw.g, 0.3)
                               WHERE tw.subtype IN ('locality', 'localadmin') AND tl.id NOT IN (SELECT id FROM tl JOIN tw
                                   ON ST_Contains(tl.g, tw.g)) GROUP BY tl.id""").fetchall()
        return {i: tidy(n) for i, n in rows}

    def places_in(cc: str, leaves: list) -> dict[str, list[str]]:
        """For twins: the named places inside each, readable ones only, the most specific and prominent first — one
        query for all."""
        if not leaves:
            return {}
        ps = [(p[0], p[1], p[2], p[3], p[4]) for p in places(cc)
              if not POOR_SCRIPT.search(p[1]) and LATIN.search(local_name(None, p[1]) or "")]
        if not ps:
            return {}
        con.execute("CREATE OR REPLACE TEMP TABLE pw (subtype VARCHAR, name VARCHAR, prominence INTEGER, population BIGINT, g GEOMETRY)")
        con.executemany("INSERT INTO pw VALUES (?, ?, ?, ?, ST_GeomFromWKB(?))", ps)
        con.execute("CREATE OR REPLACE TEMP TABLE pl (id VARCHAR, g GEOMETRY)")
        con.executemany("INSERT INTO pl VALUES (?, ST_MakeValid(ST_GeomFromWKB(?)))", [(l[0], l[6]) for l in leaves])
        spec = "CASE pw.subtype " + " ".join(f"WHEN '{k}' THEN {v}" for k, v in SPECIFIC.items()) + " ELSE 9 END"
        out: dict[str, list[str]] = collections.defaultdict(list)
        for i, n in con.execute(f"""SELECT pl.id, pw.name FROM pl JOIN pw ON ST_Contains(pl.g, pw.g)
                                    ORDER BY pl.id, {spec}, pw.prominence DESC, pw.population DESC""").fetchall():
            t = tidy(local_name(None, n) or "")
            if t and t not in out[i]:
                out[i].append(t)
        return out

    def numbered_name(cc: str, leaf, parent: str | None) -> str:
        """A region the data names only by a number (Qatar's zones): the most specific named place inside it —
        "Al Sadd (Zone 23)" — or, with none inside, "Zone 23, Doha"."""
        best = place_inside(cc, leaf)
        if best:
            return f"{best} (Zone {leaf[3]})"
        return f"Zone {leaf[3]}, {parent}" if parent else f"Zone {leaf[3]}"
    stats: collections.Counter = collections.Counter()
    per_country: dict[str, dict] = {}
    os.makedirs(OUT, exist_ok=True)
    for old in ([] if ONLY else os.listdir(OUT)):   # a trial run of some countries leaves the rest alone
        os.remove(os.path.join(OUT, old))

    def rounded(c):
        return [rounded(x) for x in c] if isinstance(c[0], list) else [round(c[0], 4), round(c[1], 4)]

    def feature(props: dict, wkb: bytes) -> None:
        # snapped to the 4 decimals the file keeps, and valid at that precision (rounding alone can make a ring touch itself)
        g = json.loads(con.execute("SELECT ST_AsGeoJSON(ST_CollectionExtract(ST_MakeValid(ST_ReducePrecision(ST_GeomFromWKB(?), 0.0001)), 3))",
                                   [wkb]).fetchone()[0])
        g["coordinates"] = rounded(g["coordinates"])
        out.write(json.dumps({"type": "Feature", "properties": props, "geometry": g}, ensure_ascii=False, separators=(",", ":")) + "\n")

    for cc in countries:
        if ONLY and cc not in ONLY.split(","):
            continue
        out = gzip.open(os.path.join(OUT, f"{cc}.geojsonl.gz"), "wt", encoding="utf-8")
        rows = con.execute("SELECT id, subtype, parent, name, population, area, wkb FROM a WHERE country = ? "
                           f"AND (subtype = 'country' OR lower(coalesce(name, '')) NOT IN {PLACEHOLDERS})", [cc]).fetchall()   # a placeholder name is no region
        rows = [r for r in rows if r[3] not in NOT_PLACES]
        country = max((r for r in rows if r[1] == "country"), key=lambda r: r[5])
        c_area = country[5] or 1e-12
        feature({"level": 1, "country": cc, "name": english.get(cc) or country[3] or cc}, country[6])

        by_parent: dict[str, list] = collections.defaultdict(list)
        for r in rows:
            if r[1] != "country" and r[2]:
                by_parent[r[2]].append(r)
        def below(r, kind, depth=0):
            """Every descendant of one kind under r, at whatever depth (a level in between may be broken)."""
            out = []
            for k in by_parent.get(r[0], []):
                if k[1] == kind:
                    out.append(k)
                elif depth < 6:
                    out.extend(below(k, kind, depth + 1))
            return out

        def km2(deg2, wkb):
            """An area in deg² as km², at the latitude of the shape it lies in."""
            lat = con.execute("SELECT ST_Y(ST_Centroid(ST_GeomFromWKB(?)))", [wkb]).fetchone()[0] or 0
            return (deg2 or 0) * 12392 * max(math.cos(math.radians(lat)), 0.01)

        def quarters(r, kind, depth=0):
            """The quarters of one kind under r, reached through states, counties and coarser quarters only."""
            out = []
            for k in by_parent.get(r[0], []):
                if k[1] == kind:
                    out.append(k)
                elif k[1] not in ("locality", "localadmin") and depth < 6:
                    out.extend(quarters(k, kind, depth + 1))
            return out

        clipped: dict[tuple, float] = {}

        def share(group, r):
            """How much of r the group covers. Above CHILD_COVER the children are measured inside r only: Finland's
            municipalities reach out over the sea past their sub-region's outline, which is not an overlap."""
            raw = sum(k[5] for k in group) / (r[5] or 1e-12)
            if raw <= CHILD_COVER[1] or len(group) > 5000:
                return raw
            key = (r[0], tuple(sorted(k[0] for k in group)))
            if key not in clipped:
                con.execute("CREATE OR REPLACE TEMP TABLE sk (g GEOMETRY)")
                con.executemany("INSERT INTO sk VALUES (ST_MakeValid(ST_GeomFromWKB(?)))", [(k[6],) for k in group])
                clipped[key] = con.execute("SELECT sum(ST_Area(ST_Intersection(sk.g, ST_MakeValid(ST_GeomFromWKB(?))))) FROM sk",
                                           [r[6]]).fetchone()[0] / (r[5] or 1e-12)
            return clipped[key]

        def clip(group, r):
            """A group that covers r only when measured inside r (clipped in share): each member cut to r, so what
            reaches out over the sea (an atoll's municipality and its ocean) never becomes part of a local region."""
            if sum(k[5] for k in group) / (r[5] or 1e-12) <= CHILD_COVER[1]:
                return group
            con.execute("CREATE OR REPLACE TEMP TABLE ck (id VARCHAR, g GEOMETRY)")
            con.executemany("INSERT INTO ck VALUES (?, ST_MakeValid(ST_GeomFromWKB(?)))", [(k[0], k[6]) for k in group])
            cut = {i: (a, bytes(w)) for i, a, w in con.execute(
                "SELECT id, ST_Area(x), ST_AsWKB(x) FROM (SELECT id, ST_CollectionExtract(ST_MakeValid(ST_Intersection(ck.g, "
                "ST_MakeValid(ST_GeomFromWKB(?)))), 3) AS x FROM ck) WHERE NOT ST_IsEmpty(x)", [r[6]]).fetchall()}
            return [(*k[:5], *cut[k[0]]) for k in group if k[0] in cut]

        def unnest(group, r):
            """Children of one kind that overlap because the data files two levels as one (the Federation of Bosnia
            and Herzegovina: its cantons and its municipalities, all counties under it): the bigger ones are the
            level, and each finer one moves under the one it lies in, to divide it there."""
            # never among towns: a town filed twice with its villages under it (Drei Gleichen) stays one town
            if len(group) < 2 or group[0][1] == "locality" or share(group, r) <= CHILD_COVER[1] or len(group) > 5000:
                return group
            con.execute("CREATE OR REPLACE TEMP TABLE ng (id VARCHAR, a DOUBLE, g GEOMETRY)")
            con.executemany("INSERT INTO ng VALUES (?, ?, ST_MakeValid(ST_GeomFromWKB(?)))", [(k[0], k[5] or 0, k[6]) for k in group])
            inside = con.execute("""SELECT o.id, arg_min(m.id, m.a) FROM ng m JOIN ng o ON o.id <> m.id AND o.a < m.a
                                        AND ST_Contains(m.g, ST_PointOnSurface(o.g)) GROUP BY o.id""").fetchall()
            if not inside:
                return group
            moved = dict(inside)
            for k in group:
                if k[0] in moved:
                    by_parent[moved[k[0]]].append(k)
            by_parent[r[0]] = [k for k in by_parent[r[0]] if k[0] not in moved]
            return [k for k in group if k[0] not in moved]

        def nested(group) -> dict[str, str]:
            """Members of a group that lie inside a bigger member (a town and the municipality around it, both children
            of one state; a village filed as a town inside its town): id → the smallest member holding it. Two such
            must never both be local regions — they would lie on top of each other."""
            if len(group) < 2 or len(group) > 20000:
                return {}
            con.execute("CREATE OR REPLACE TEMP TABLE ng (id VARCHAR, a DOUBLE, g GEOMETRY)")
            con.executemany("INSERT INTO ng VALUES (?, ?, ST_MakeValid(ST_GeomFromWKB(?)))", [(k[0], k[5] or 0, k[6]) for k in group])
            return dict(con.execute("""SELECT o.id, arg_min(m.id, m.a) FROM ng m JOIN ng o ON o.id <> m.id AND o.a < m.a
                                           AND ST_Contains(m.g, ST_PointOnSurface(o.g))
                                           AND ST_Area(ST_Intersection(m.g, o.g)) > 0.5 * o.a GROUP BY o.id""").fetchall())

        def denest(group):
            """The group without its nested members, and where each of them goes."""
            inside = nested(group)
            return [k for k in group if k[0] not in inside], inside

        def take(r, group, inside):
            """The chosen group's nested members move under the member holding them, to divide it there if they
            can — but a village inside a town stays part of the town (Drei Gleichen)."""
            kind = {k[0]: k[1] for k in group}
            for k in group:
                if k[0] in inside and not (k[1] == "locality" and kind.get(inside[k[0]]) == "locality"):
                    by_parent[inside[k[0]]].append(k)
            if inside:
                by_parent[r[0]] = [k for k in by_parent.get(r[0], []) if k[0] not in inside]
            return [k for k in group if k[0] not in inside]

        def settle(r, group, leaves):
            """What r's children leave of r: each piece smaller than the child beside it joins that child; the pieces
            bigger than their neighbour become ONE region under r's name (England's unparished towns, the land between
            Lebanon's municipalities). At country level nothing here — the country-wide fill does it, never under the
            country's name."""
            if r[1] == "country" or share(group, r) >= 0.99:
                return leaves
            grown, rest = fill_gaps(con, r[6], leaves, [r])
            leaves = [(*l[:5], con.execute("SELECT ST_Area(ST_GeomFromWKB(?))", [grown[l[0]]]).fetchone()[0], grown[l[0]])
                      if l[0] in grown else l for l in leaves]
            return leaves + rest

        def refine(r, depth=0):
            kids = [k for k in by_parent.get(r[0], []) if k[1] in KIDS]
            if kids and depth < 8:
                # 1. children of one kind (or all together) that cover r — checked for r itself, never pooled over a
                #    bigger area; of several, the one with the most members
                # (every group without its nested members: a member inside another is never a region beside it)
                full = [unnest([k for k in kids if k[1] == st], r) for st in KIDS] + [kids]
                cut = [denest(g) for g in full]
                direct = [kept for kept, _ in cut]
                fits = [i for i, g in enumerate(direct) if g and CHILD_COVER[0] <= share(g, r) <= CHILD_COVER[1]]
                if fits:
                    i = max(fits, key=lambda i: len(direct[i]))
                    group = clip(take(r, full[i], cut[i][1]), r)
                    return settle(r, group, [leaf for k in group for leaf in refine(k, depth + 1)])
                # 2. only if no level of children fits: descendants of one kind further down (a broken level in between
                #    is skipped — Morocco's regions reach into Western Sahara and over-cover the country)
                deep = [denest(g)[0] for g in (below(r, st) for st in KIDS) if g]
                fits2 = [g for g in deep if g and CHILD_COVER[0] <= share(g, r) <= CHILD_COVER[1]]
                if fits2:
                    group = clip(max(fits2, key=len), r)
                    return settle(r, group, [leaf for k in group for leaf in refine(k, depth + 1)])
                # 3. children of one kind that cover at least half of r (a division with holes, not scattered places):
                #    they stay at their true size, and the rest of r becomes a region of its own under r's name. Not at
                #    country level: what the regions leave there is filled into neighbours, so no country appears as a
                #    local region.
                part = [i for i, g in enumerate(direct[:-1]) if g and PART_MIN <= share(g, r) < CHILD_COVER[0]]
                if part:
                    i = max(part, key=lambda i: share(direct[i], r))
                    group = clip(take(r, full[i], cut[i][1]), r)
                    return settle(r, group, [leaf for k in group for leaf in refine(k, depth + 1)])
            # 4. a state or county with no towns mapped in it: its quarters, the coarsest kind that covers at least
            #    half of it (Dubai's communities) — found under it without passing through a town, so no town is split
            if r[1] in ("country", "region", "county") and km2(r[5], r[6]) > BIG_KM2:
                for kind in HOODS:
                    group = quarters(r, kind)
                    if sum(1 for k in group if WARD.search(k[3] or "")) * 2 >= len(group):   # electoral wards, not places
                        continue
                    if group and PART_MIN <= share(group, r) <= CHILD_COVER[1] and \
                            km2(sorted(k[5] for k in group)[len(group) // 2], r[6]) >= MIN_QUARTER_KM2:
                        return settle(r, group, group)
            return [r]

        leaves = refine(country)
        leaves = [l for l in leaves if l[1] != "country"] or [country]   # a country with nothing finer is its own one region
        if len(leaves) == 1 and leaves[0] is not country:   # one region for the whole country (Monaco): all of the country
            leaves = [(*leaves[0][:5], country[5], country[6])]
        leaves = drop_tiny(con, cut_overlaps(con, carve_towns(con, leaves, rows)))
        before = {l[0]: l[5] or 0 for l in leaves}   # each region's own area, to see what the gap fill added
        if len(leaves) > 1:
            def island_name(wkb):
                """The most prominent named place on a far island, or the nearest one within ~30 km."""
                inside = [p for p in places(cc)
                          if con.execute("SELECT ST_Contains(ST_GeomFromWKB(?), ST_GeomFromWKB(?))", [wkb, p[4]]).fetchone()[0]]
                if not inside:
                    near = [(con.execute("SELECT ST_Distance(ST_GeomFromWKB(?), ST_GeomFromWKB(?))", [wkb, p[4]]).fetchone()[0], p)
                            for p in places(cc)]
                    inside = [p for d, p in sorted(near, key=lambda x: x[0])[:1] if d < 0.3]
                if not inside:
                    return None
                return local_name(None, min(inside, key=lambda p: (-p[2], -p[3]))[1])

            grown, holes = fill_gaps(con, country[6], leaves, [r for r in rows if r[1] in ("region", "county")], island_name)
            # a region that grew carries its new area too, so the checks in the manifest measure the filled level
            leaves = [(*l[:5], con.execute("SELECT ST_Area(ST_GeomFromWKB(?))", [grown[l[0]]]).fetchone()[0], grown[l[0]])
                      if l[0] in grown else l for l in leaves]
            by_id: dict[str, tuple] = {}
            for l in leaves + holes:
                if l[0] in by_id:   # an area's rest from settle() and its big holes from here are one region
                    o = by_id[l[0]]
                    w, a = con.execute("SELECT ST_AsWKB(g), ST_Area(g) FROM (SELECT ST_MakeValid(ST_Union("
                                       "ST_GeomFromWKB(?), ST_GeomFromWKB(?))) AS g)", [o[6], l[6]]).fetchone()
                    by_id[l[0]] = (*o[:5], a, bytes(w))
                    before[l[0]] = before.get(l[0], 0) + (l[5] or 0)   # the hole is its own land, not growth
                else:
                    by_id[l[0]] = l
            leaves = list(by_id.values())
            stats["gaps_filled_regions"] += len(grown)
            stats["holes_as_regions"] += len(holes)
        def named(up):   # an ancestor's name in English; never the country (the board names it); none for one that reads as nothing
            if up in country_ids:
                return None
            return None if unreadable(up, names.get(up)) else local_name(up, names.get(up))

        def ancestors(l, base):
            """The names above a region, nearest first, skipping any that repeats the region's own name or one already
            given ("Armagh City, Banbridge and Craigavon" after "… District")."""
            out, up = [], l[2]
            while up and len(out) < 3:
                n = named(up)
                if n and n != base and not any(re.search(rf"(?<!\w){re.escape(a)}(?!\w)", b)
                                               for o in out for a, b in ((n, o), (o, n))):
                    out.append(n)
                up = parent_of.get(up)
            return out

        town_of = towns_inside(cc, [l for l in leaves if unreadable(l[0], l[3])])
        base_of = {}
        for l in leaves:
            name = town_of.get(l[0]) or local_name(l[0], l[3]) or names.get(l[2]) or cc
            if not re.search(r"[^\W\d_]", name):   # a number, not a name
                name = numbered_name(cc, l, local_name(l[2], names.get(l[2])))
            elif name.lower() in PLACEHOLDERS:   # a placeholder in the data (El Salvador's "Null", Italy's "None"): a place inside, else the area above
                name = place_inside(cc, l) or local_name(l[2], names.get(l[2])) or english.get(cc) or cc
            base_of[l[0]] = name
        # twins get as many ancestors as it takes to tell them apart: "Mitte (Berlin)",
        # "Jackson (Jackson County, Ohio)"; true namesakes in one and the same place stay alike
        label = dict(base_of)
        for depth in (1, 2):
            seen = collections.Counter(label.values())
            for l in leaves:
                if seen[label[l[0]]] > 1:
                    up = ancestors(l, base_of[l[0]])[:depth]
                    if up:
                        label[l[0]] = f"{base_of[l[0]]} ({', '.join(up)})"
        # twins in one and the same district (townlands, a borough and its township) are told apart by where
        # they lie towards each other: "Largy (…, north)", "Largy (…, south)"
        def add(name, word):
            return f"{name[:-1]}, {word})" if name.endswith(")") else f"{name} ({word})"

        # an area's land outside its towns, when it shares the name of its town: "Larnaca (surroundings)"
        seen = collections.Counter(label.values())
        towns = {label[l[0]] for l in leaves if "#rest" not in l[0]}
        for l in leaves:
            if "#rest" in l[0] and seen[label[l[0]]] > 1 and label[l[0]] in towns:
                label[l[0]] = add(label[l[0]], "surroundings")
        # a municipality around a town of the same name (Poland's rural gminas, a US township around its borough):
        # the ring is the town's surroundings, "Radomsko (…, surroundings)"
        def twin_groups():
            seen = collections.Counter(label.values())
            groups = collections.defaultdict(list)
            for l in leaves:
                if seen[label[l[0]]] > 1:
                    groups[label[l[0]]].append(l)
            return groups

        for name, group in twin_groups().items():
            outer = [o for o in group if any(
                (i[5] or 0) < (o[5] or 0) and con.execute(
                    "SELECT ST_Contains(ST_ConvexHull(ST_GeomFromWKB(?)), ST_PointOnSurface(ST_GeomFromWKB(?)))",
                    [o[6], i[6]]).fetchone()[0] for i in group if i is not o)]
            if len(outer) == 1:
                label[outer[0][0]] = add(name, "surroundings")
        # twins still alike are told apart by a named place inside each, where every one has its own: "Kairouan
        # (Medina)", "Kairouan (Sahbi)"; only then by where they lie
        twins = twin_groups()
        inside = places_in(cc, [l for g in twins.values() for l in g])
        for name, group in twins.items():
            used: set[str] = set()
            for l in sorted(group, key=lambda l: -(l[5] or 0)):
                pick = next((p for p in inside.get(l[0], []) if p not in used and p.lower() not in name.lower()), None)
                if pick:
                    used.add(pick)
                    label[l[0]] = add(name, pick)
        twins = twin_groups()
        compass = ("east", "north-east", "north", "north-west", "west", "south-west", "south", "south-east")
        for name, group in twins.items():
            at = [con.execute("SELECT ST_X(c), ST_Y(c) FROM (SELECT ST_Centroid(ST_GeomFromWKB(?)) c)", [l[6]]).fetchone()
                  for l in group]
            mx, my = sum(x for x, _ in at) / len(at), sum(y for _, y in at) / len(at)
            if len(group) == 2:   # two: the axis they are furthest apart on
                (x0, y0), (x1, y1) = at
                words = (["north", "south"] if y0 > y1 else ["south", "north"]) if abs(y0 - y1) >= abs(x0 - x1) else \
                        (["east", "west"] if x0 > x1 else ["west", "east"])
            else:
                words = [compass[round(math.degrees(math.atan2(y - my, x - mx)) / 45) % 8] for x, y in at]
            taken = collections.Counter(words)
            for l, w, k in zip(group, words, range(len(group))):
                w = w if taken[w] == 1 else f"{w} {sum(1 for o in words[:k + 1] if o == w)}"
                label[l[0]] = add(name, w)
        name_of: dict[str, str] = {}
        for l in leaves:
            name = label[l[0]]
            name_of[l[0]] = name
            feature({"level": 2, "country": cc, "name": name, "population": l[4]}, l[6])
            stats[l[1]] += 1
        out.close()
        stats["countries"] += 1
        areas = sorted((l[5] or 0) for l in leaves)
        growth = [((l[5] or 0) / before[l[0]]) if before.get(l[0]) else 1.0 for l in leaves]
        written = collections.Counter(name_of[l[0]] for l in leaves)
        rests = [l for l in leaves if "#" in l[0]]
        per_country[cc] = {"name": english.get(cc) or country[3] or cc, "local_regions": len(leaves),
                           "kinds": dict(collections.Counter(l[1] for l in leaves)),
                           # the checks, read per country to find a division that makes no sense:
                           "cover": round(sum(areas) / c_area, 3),                         # share of the country covered
                           "largest_share": round(areas[-1] / c_area, 3) if areas else 0,  # the biggest region's share
                           "largest_to_median": round(areas[-1] / (areas[len(areas) // 2] or 1e-12), 1) if areas else 0,
                           "inflated": sum(1 for g in growth if g > 2),                    # regions the gap fill more than doubled
                           "max_growth": round(max(growth), 1) if growth else 1,
                           "rests": len(rests), "rest_share": round(sum(l[5] or 0 for l in rests) / c_area, 3),
                           "duplicate_names": sum(n - 1 for n in written.values() if n > 1)}
        os.makedirs(PARTS, exist_ok=True)
        with open(os.path.join(PARTS, f"{cc}.json"), "w", encoding="utf-8") as fh:
            json.dump(per_country[cc], fh, ensure_ascii=False)
        print(f"{cc}: {len(leaves)} local regions", file=sys.stderr)

    if ONLY:   # a rebuild of some countries: their entries in the manifest are replaced, the rest stays
        path = os.path.join(OUT, "manifest.json")
        if not os.path.exists(path):
            print(f"trial run {ONLY}: {dict(stats)} — no manifest to update", file=sys.stderr)
            return
        old = json.load(open(path, encoding="utf-8"))
        old["per_country"].update(per_country)
        old["local_regions"] = sum(v["local_regions"] for v in old["per_country"].values())
        with open(path, "w", encoding="utf-8") as fh:
            json.dump(old, fh, ensure_ascii=False, indent=1, sort_keys=True)
        print(f"rebuilt {ONLY}: {dict(stats)} — manifest updated, {old['local_regions']:,} local regions", file=sys.stderr)
        return
    with open(os.path.join(OUT, "manifest.json"), "w", encoding="utf-8") as fh:
        json.dump({"source": "Overture Maps Foundation, divisions theme", "release": RELEASE,
                   "licence": LICENCE,
                   "countries": stats["countries"], "local_regions": sum(v["local_regions"] for v in per_country.values()),
                   "per_country": per_country}, fh, ensure_ascii=False, indent=1, sort_keys=True)
    size = sum(os.path.getsize(os.path.join(OUT, f)) for f in os.listdir(OUT))
    print(f"done in {round(time.time() - t0)} s — {dict(stats)} — {OUT} {size / 1e6:.0f} MB", file=sys.stderr)


if __name__ == "__main__":
    build()