geointel/scripts/provision_belgium_north_sea_scope.py

"""Provision Belgium and the Belgian North Sea through the canonical API.

This explicit operator downloads a bounded, allowlisted set of authoritative
sources. It never runs at application startup and never writes directly to
PostGIS or vector_features.
"""

from __future__ import annotations

import argparse
import hashlib
import json
import os
import shutil
import sys
import zipfile
from datetime import datetime, timezone
from pathlib import Path
from typing import Any

import requests
from requests.adapters import HTTPAdapter
from shapely.geometry import GeometryCollection, MultiPolygon, Polygon, mapping, shape
from shapely.ops import unary_union
from shapely.validation import make_valid
from urllib3.util.retry import Retry


PROJECT_NAME = "Belgium and North Sea Workbench"
PROJECT_REGION = "Belgium and Belgian North Sea"
DEFAULT_API_URL = "http://127.0.0.1:8000"
DEFAULT_OUTPUT_ROOT = Path("/app/storage/operator-data/geographic-scopes/belgium-north-sea")
NGI_ARCHIVE_URL = (
    "https://ac.ngi.be/remoteclient-open/ngi-standard-open/Vectordata/"
    "TerritorialDivisions/TerritorialDivisions-AdminVector/"
    "fb1e2993-2020-428c-9188-eb5f75e284b9_geopackage+sqlite3_4326.zip"
)
NGI_CATALOG_URL = "https://www.geo.be/catalog/details/fb1e2993-2020-428c-9188-eb5f75e284b9"
NGI_ATTRIBUTION = "National Geographic Institute (NGI), AdminVector"
NGI_LICENSE = "CC BY 4.0"
RBINS_MRU_WFS_URL = "https://spatial.naturalsciences.be/geoserver/od_nature/ows"
RBINS_MRU_METADATA_URL = (
    "https://metadata.naturalsciences.be/geonetwork/srv/api/records/"
    "29f40b0d-2a3e-49a8-870a-e9b4acd4d1e3"
)
RBINS_MRU_LAYER = "od_nature:marine_reporting_units_2024"
RBINS_MSP_WFS_URL = "https://spatial.naturalsciences.be/geoserver/ows"
RBINS_MSP_SOURCE_URL = (
    "https://www.health.belgium.be/en/themes/environment/marine-environment/marine-spatial-plan"
)
MSP_VALID_FROM = "2026-03-20T00:00:00Z"
MAX_ARCHIVE_BYTES = 160 * 1024 * 1024
MAX_EXTRACTED_BYTES = 1024 * 1024 * 1024
MAX_WFS_RESPONSE_BYTES = 64 * 1024 * 1024
WFS_PAGE_SIZE = 5000
GEOJSON_CRS = {"type": "name", "properties": {"name": "EPSG:4326"}}

ADMIN_LAYERS = {
    "belgianterritory": (1, 1),
    "belgianmaritimezone": (1, 1),
    "region": (3, 3),
    "province": (11, 11),
    "municipality": (560, 590),
}

MSP_LAYERS = (
    "imsp26:bmsp_aquaculture_zone",
    "imsp26:bmsp_coastal_protection_experiment_zone",
    "imsp26:bmsp_commercial_industrial_zone",
    "imsp26:bmsp_cultural_heritage",
    "imsp26:bmsp_dredging_zone",
    "imsp26:bmsp_energy_cables_pipelines_zone",
    "imsp26:bmsp_fisheries_zone",
    "imsp26:bmsp_conservation_zone",
    "imsp26:bmsp_measuring_poles",
    "imsp26:bmsp_military_zone",
    "imsp26:bmsp_port_expansion_zone",
    "imsp26:bmsp_radar_towers",
    "imsp26:bmsp_research_recreation_zone",
    "imsp26:bmsp_extraction_zone",
    "imsp26:bmsp_shipping_ports_zone",
)

MARINE_REPORTING_IDS = {
    "belgian_north_sea": "ANS-BE-MS-1",
    "territorial_1_12": "ANS-BE-AA-TEW",
    "coastal_0_1": "ANS-BE-AA-CW",
    "offshore": "ANS-BE-AA-OFFSHORE",
}

AREA_NAMES = {
    "belgium": "Belgium land",
    "flanders": "Flanders",
    "wallonia": "Wallonia",
    "brussels": "Brussels-Capital Region",
    "belgian_north_sea": "Belgian part of the North Sea",
    "territorial_sea": "Belgian territorial sea (0-12 nautical miles)",
    "exclusive_economic_zone": "Belgian exclusive economic zone beyond territorial sea",
    "continental_shelf": "Belgian continental shelf beyond territorial sea",
}

THEME_BY_LAYER = {
    "belgium_land_boundary": "administrative",
    "belgium_regions": "administrative",
    "belgium_provinces": "administrative",
    "belgium_municipalities": "administrative",
    "marine_legal_scopes": "marine_environment",
    "marine_spatial_plan_2026": "maritime_planning",
}


def parse_args() -> argparse.Namespace:
    parser = argparse.ArgumentParser(description="Provision the Belgium and Belgian North Sea scope.")
    parser.add_argument("--base-url", default=os.environ.get("GEOINTEL_INTERNAL_API_URL", DEFAULT_API_URL))
    parser.add_argument(
        "--output-root",
        type=Path,
        default=Path(os.environ.get("GEOINTEL_NATIONAL_SCOPE_OUTPUT_ROOT", DEFAULT_OUTPUT_ROOT)),
    )
    parser.add_argument("--request-timeout", type=int, default=300)
    parser.add_argument("--import-timeout", type=int, default=3600)
    parser.add_argument("--force", action="store_true")
    parser.add_argument("--fetch-only", action="store_true")
    return parser.parse_args()


def utc_now() -> str:
    return datetime.now(timezone.utc).isoformat()


def sha256_file(path: Path) -> str:
    digest = hashlib.sha256()
    with path.open("rb") as handle:
        for chunk in iter(lambda: handle.read(1024 * 1024), b""):
            digest.update(chunk)
    return digest.hexdigest()


def write_json_atomic(path: Path, payload: dict[str, Any], *, pretty: bool = False) -> None:
    path.parent.mkdir(parents=True, exist_ok=True)
    temporary = path.with_suffix(f"{path.suffix}.partial")
    temporary.write_text(
        json.dumps(
            payload,
            ensure_ascii=False,
            indent=2 if pretty else None,
            separators=None if pretty else (",", ":"),
            sort_keys=pretty,
        ),
        encoding="utf-8",
    )
    temporary.replace(path)


def build_session() -> requests.Session:
    retry = Retry(
        total=5,
        connect=5,
        read=5,
        status=5,
        backoff_factor=1.0,
        status_forcelist=(429, 500, 502, 503, 504),
        allowed_methods=frozenset({"GET"}),
        raise_on_status=True,
    )
    session = requests.Session()
    adapter = HTTPAdapter(max_retries=retry)
    session.mount("https://", adapter)
    session.mount("http://", adapter)
    session.headers.update({"User-Agent": "GeoIntel-Belgium-North-Sea-Operator/1.0"})
    return session


def download_limited(session: requests.Session, url: str, target: Path, timeout: int, max_bytes: int) -> Path:
    if url != NGI_ARCHIVE_URL:
        raise RuntimeError("Archive URL is not in the fixed operator allowlist")
    temporary = target.with_suffix(f"{target.suffix}.partial")
    target.parent.mkdir(parents=True, exist_ok=True)
    size = 0
    with session.get(url, stream=True, timeout=timeout) as response:
        response.raise_for_status()
        content_length = int(response.headers.get("content-length") or 0)
        if content_length > max_bytes:
            raise RuntimeError(f"NGI archive declares {content_length} bytes, above the {max_bytes} byte limit")
        with temporary.open("wb") as handle:
            for chunk in response.iter_content(chunk_size=1024 * 1024):
                if not chunk:
                    continue
                size += len(chunk)
                if size > max_bytes:
                    raise RuntimeError(f"NGI archive exceeded the {max_bytes} byte limit")
                handle.write(chunk)
    if size == 0:
        raise RuntimeError("NGI archive download was empty")
    temporary.replace(target)
    return target


def extract_single_geopackage(archive_path: Path, output_dir: Path) -> Path:
    extraction_root = output_dir / "adminvector"
    temporary_root = output_dir / "adminvector.partial"
    if temporary_root.exists():
        shutil.rmtree(temporary_root)
    temporary_root.mkdir(parents=True)
    with zipfile.ZipFile(archive_path) as archive:
        members = [member for member in archive.infolist() if not member.is_dir()]
        gpkg_members = [member for member in members if Path(member.filename).suffix.lower() == ".gpkg"]
        if len(gpkg_members) != 1:
            raise RuntimeError(f"Expected one GeoPackage in NGI archive, received {len(gpkg_members)}")
        member = gpkg_members[0]
        member_path = Path(member.filename)
        if member_path.is_absolute() or ".." in member_path.parts:
            raise RuntimeError("NGI archive contains an unsafe GeoPackage path")
        if member.file_size <= 0 or member.file_size > MAX_EXTRACTED_BYTES:
            raise RuntimeError("NGI GeoPackage size is empty or above the extraction limit")
        target = temporary_root / "adminvector_4326.gpkg"
        with archive.open(member) as source, target.open("wb") as destination:
            shutil.copyfileobj(source, destination, length=1024 * 1024)
    if extraction_root.exists():
        shutil.rmtree(extraction_root)
    temporary_root.replace(extraction_root)
    return extraction_root / "adminvector_4326.gpkg"


def _polygonal(geometry):
    if geometry is None or geometry.is_empty:
        return None
    if not geometry.is_valid:
        geometry = make_valid(geometry)
    if isinstance(geometry, (Polygon, MultiPolygon)):
        return geometry
    if isinstance(geometry, GeometryCollection):
        parts = [part for part in geometry.geoms if isinstance(part, (Polygon, MultiPolygon))]
        merged = unary_union(parts) if parts else None
        if isinstance(merged, (Polygon, MultiPolygon)) and not merged.is_empty:
            return merged
    return None


def read_adminvector_layers(gpkg_path: Path) -> dict[str, dict[str, Any]]:
    try:
        import geopandas
        import pandas
        import pyogrio
    except ImportError as exc:
        raise RuntimeError("GeoPandas and Pyogrio are required for the national scope operator") from exc

    available_layers = {str(row[0]) for row in pyogrio.list_layers(gpkg_path)}
    missing = sorted(set(ADMIN_LAYERS) - available_layers)
    if missing:
        raise RuntimeError(f"NGI AdminVector is missing expected layers: {', '.join(missing)}")

    payloads: dict[str, dict[str, Any]] = {}
    for layer_name, (minimum, maximum) in ADMIN_LAYERS.items():
        frame = geopandas.read_file(gpkg_path, layer=layer_name)
        if frame.crs is None:
            raise RuntimeError(f"NGI layer {layer_name} has no CRS")
        frame = frame.to_crs(4326)
        if not minimum <= len(frame) <= maximum:
            raise RuntimeError(
                f"NGI layer {layer_name} has {len(frame)} records; expected between {minimum} and {maximum}"
            )
        payload = json.loads(frame.to_json(drop_id=False, default=str))
        source_max_modification = None
        if "modifdate" in frame.columns:
            parsed = pandas.to_datetime(frame["modifdate"], errors="coerce", utc=True).dropna()
            if not parsed.empty:
                source_max_modification = parsed.max().isoformat()
        payload.update(
            {
                "name": f"NGI AdminVector {layer_name}",
                "crs": GEOJSON_CRS,
                "source_url": NGI_CATALOG_URL,
                "attribution": NGI_ATTRIBUTION,
                "license": NGI_LICENSE,
                "source_max_modification": source_max_modification,
            }
        )
        for feature in payload["features"]:
            properties = feature.setdefault("properties", {})
            properties.update(
                {
                    "source_name": "ngi_adminvector",
                    "source_layer": layer_name,
                    "source_feature_id": str(feature.get("id") or properties.get("tgid") or ""),
                    "authority_level": "authoritative",
                    "attribution": NGI_ATTRIBUTION,
                    "source_url": NGI_CATALOG_URL,
                }
            )
        payloads[layer_name] = payload
    return payloads


def _response_json_limited(response: requests.Response, max_bytes: int) -> dict[str, Any]:
    response.raise_for_status()
    content = response.content
    if len(content) > max_bytes:
        raise RuntimeError(f"WFS response exceeded the {max_bytes} byte limit")
    try:
        payload = response.json()
    except ValueError as exc:
        raise RuntimeError("WFS returned non-JSON content") from exc
    if payload.get("type") != "FeatureCollection" or not isinstance(payload.get("features"), list):
        raise RuntimeError("WFS response is not a GeoJSON FeatureCollection")
    return payload


def fetch_wfs_layer(
    session: requests.Session,
    *,
    service_url: str,
    layer_name: str,
    timeout: int,
    page_size: int = WFS_PAGE_SIZE,
) -> dict[str, Any]:
    if (service_url, layer_name) not in {
        (RBINS_MRU_WFS_URL, RBINS_MRU_LAYER),
        *((RBINS_MSP_WFS_URL, layer) for layer in MSP_LAYERS),
    }:
        raise RuntimeError(f"WFS layer is not in the fixed operator allowlist: {layer_name}")
    features: list[dict[str, Any]] = []
    source_ids: set[str] = set()
    start_index = 0
    matched: int | None = None
    while True:
        response = session.get(
            service_url,
            params={
                "service": "WFS",
                "version": "2.0.0",
                "request": "GetFeature",
                "typeNames": layer_name,
                "outputFormat": "application/json",
                "srsName": "EPSG:4326",
                "count": page_size,
                "startIndex": start_index,
            },
            timeout=timeout,
        )
        payload = _response_json_limited(response, MAX_WFS_RESPONSE_BYTES)
        page = payload["features"]
        if matched is None and str(payload.get("numberMatched", "")).isdigit():
            matched = int(payload["numberMatched"])
        for feature in page:
            source_id = str(feature.get("id") or "")
            if source_id and source_id in source_ids:
                raise RuntimeError(f"WFS layer {layer_name} returned duplicate feature id {source_id}")
            if source_id:
                source_ids.add(source_id)
            features.append(feature)
        if not page or len(page) < page_size or (matched is not None and len(features) >= matched):
            break
        start_index += len(page)
    if matched is not None and len(features) != matched:
        raise RuntimeError(f"WFS layer {layer_name} returned {len(features)} of {matched} matched features")
    return {
        "type": "FeatureCollection",
        "name": layer_name,
        "crs": GEOJSON_CRS,
        "features": features,
    }


def _reporting_id(feature: dict[str, Any]) -> str:
    properties = feature.get("properties") or {}
    return str(
        properties.get("MarineReportingUnitId")
        or properties.get("marineReportingUnitId")
        or properties.get("localId")
        or ""
    )


def derive_marine_scope_payload(reporting_units: dict[str, Any]) -> dict[str, Any]:
    by_id = {_reporting_id(feature): feature for feature in reporting_units.get("features") or []}
    missing = [identifier for identifier in MARINE_REPORTING_IDS.values() if identifier not in by_id]
    if missing:
        raise RuntimeError(f"Marine reporting units are missing required identifiers: {', '.join(missing)}")

    def geometry(identifier: str):
        result = _polygonal(shape(by_id[identifier].get("geometry")))
        if result is None:
            raise RuntimeError(f"Marine reporting unit {identifier} has no valid polygon geometry")
        return result

    bpns = geometry(MARINE_REPORTING_IDS["belgian_north_sea"])
    territorial = _polygonal(
        unary_union(
            (
                geometry(MARINE_REPORTING_IDS["territorial_1_12"]),
                geometry(MARINE_REPORTING_IDS["coastal_0_1"]),
            )
        )
    )
    offshore = geometry(MARINE_REPORTING_IDS["offshore"])
    if territorial is None or not bpns.covers(territorial.representative_point()) or not bpns.covers(offshore.representative_point()):
        raise RuntimeError("Derived marine legal scopes do not align with the official BPNS reporting unit")

    definitions = (
        (
            "belgian_north_sea",
            AREA_NAMES["belgian_north_sea"],
            bpns,
            "marine_water_and_seabed_scope",
            [MARINE_REPORTING_IDS["belgian_north_sea"]],
        ),
        (
            "territorial_sea",
            AREA_NAMES["territorial_sea"],
            territorial,
            "territorial_water_column_and_seabed",
            [MARINE_REPORTING_IDS["coastal_0_1"], MARINE_REPORTING_IDS["territorial_1_12"]],
        ),
        (
            "exclusive_economic_zone",
            AREA_NAMES["exclusive_economic_zone"],
            offshore,
            "water_column_rights_beyond_territorial_sea",
            [MARINE_REPORTING_IDS["offshore"]],
        ),
        (
            "continental_shelf",
            AREA_NAMES["continental_shelf"],
            offshore,
            "seabed_and_subsoil_rights_beyond_territorial_sea",
            [MARINE_REPORTING_IDS["offshore"]],
        ),
    )
    features = []
    for zone, name, legal_geometry, legal_domain, source_ids in definitions:
        features.append(
            {
                "type": "Feature",
                "id": f"belgian-marine-scope:{zone}",
                "geometry": mapping(legal_geometry),
                "properties": {
                    "name": name,
                    "coverage_zone": zone,
                    "legal_domain": legal_domain,
                    "derived_from_reporting_unit_ids": source_ids,
                    "source_name": "rbins_marine_reporting_units",
                    "source_layer": RBINS_MRU_LAYER,
                    "authority_level": "authoritative",
                    "attribution": "Royal Belgian Institute of Natural Sciences (RBINS), BMDC",
                    "source_url": RBINS_MRU_METADATA_URL,
                    "derivation": (
                        "Official reporting-unit geometry reused with explicit legal semantics"
                        if len(source_ids) == 1
                        else "Union of official 0-1 nm coastal waters and 1-12 nm territorial waters"
                    ),
                },
            }
        )
    return {
        "type": "FeatureCollection",
        "name": "Belgian marine legal scopes",
        "crs": GEOJSON_CRS,
        "features": features,
        "source_url": RBINS_MRU_METADATA_URL,
    }


def build_msp_payload(layer_payloads: dict[str, dict[str, Any]]) -> dict[str, Any]:
    if set(layer_payloads) != set(MSP_LAYERS):
        missing = sorted(set(MSP_LAYERS) - set(layer_payloads))
        raise RuntimeError(f"Marine Spatial Plan payload is incomplete: {', '.join(missing)}")
    features: list[dict[str, Any]] = []
    for layer_name in MSP_LAYERS:
        short_name = layer_name.split(":", 1)[1]
        for index, source_feature in enumerate(layer_payloads[layer_name]["features"]):
            feature = dict(source_feature)
            properties = dict(feature.get("properties") or {})
            source_id = str(feature.get("id") or f"{short_name}.{index}")
            properties.update(
                {
                    "source_name": "rbins_msp_2026",
                    "source_layer": layer_name,
                    "source_feature_id": source_id,
                    "authority_level": "authoritative",
                    "valid_from": MSP_VALID_FROM,
                    "valid_to": "2034-12-31T23:59:59Z",
                    "attribution": "Belgian federal Marine Environment service and RBINS",
                    "source_url": RBINS_MSP_SOURCE_URL,
                }
            )
            feature["id"] = f"{short_name}:{source_id}"
            feature["properties"] = properties
            features.append(feature)
    return {
        "type": "FeatureCollection",
        "name": "Belgian Marine Spatial Plan 2026-2034",
        "crs": GEOJSON_CRS,
        "features": features,
        "source_url": RBINS_MSP_SOURCE_URL,
        "valid_from": MSP_VALID_FROM,
        "valid_to": "2034-12-31T23:59:59Z",
    }


def _max_admin_modification(payloads: dict[str, dict[str, Any]], fallback: str) -> str:
    values = [
        str(payload["source_max_modification"])
        for payload in payloads.values()
        if payload.get("source_max_modification")
    ]
    value = max(values) if values else fallback
    return f"{value[:10]}T00:00:00Z"


def prepare_artifacts(args: argparse.Namespace) -> tuple[dict[str, Path], dict[str, Any]]:
    output_root = args.output_root
    output_root.mkdir(parents=True, exist_ok=True)
    archive_path = output_root / "adminvector_4326.zip"
    generated_at = utc_now()
    with build_session() as session:
        if args.force or not archive_path.is_file():
            download_limited(session, NGI_ARCHIVE_URL, archive_path, args.request_timeout, MAX_ARCHIVE_BYTES)
        gpkg_path = extract_single_geopackage(archive_path, output_root)
        admin_payloads = read_adminvector_layers(gpkg_path)
        reporting_units = fetch_wfs_layer(
            session,
            service_url=RBINS_MRU_WFS_URL,
            layer_name=RBINS_MRU_LAYER,
            timeout=args.request_timeout,
        )
        marine_scopes = derive_marine_scope_payload(reporting_units)
        msp_payloads = {
            layer: fetch_wfs_layer(
                session,
                service_url=RBINS_MSP_WFS_URL,
                layer_name=layer,
                timeout=args.request_timeout,
            )
            for layer in MSP_LAYERS
        }
        msp = build_msp_payload(msp_payloads)

    artifact_payloads = {
        "belgium_land_boundary": admin_payloads["belgianterritory"],
        "belgium_regions": admin_payloads["region"],
        "belgium_provinces": admin_payloads["province"],
        "belgium_municipalities": admin_payloads["municipality"],
        "ngi_maritime_zone": admin_payloads["belgianmaritimezone"],
        "marine_legal_scopes": marine_scopes,
        "marine_spatial_plan_2026": msp,
    }
    artifacts: dict[str, Path] = {}
    for key, payload in artifact_payloads.items():
        path = output_root / f"{key}.geojson"
        write_json_atomic(path, payload)
        artifacts[key] = path

    admin_observed_at = _max_admin_modification(admin_payloads, generated_at)
    mru_versions = [
        str((feature.get("properties") or {}).get("beginLife") or "")
        for feature in reporting_units["features"]
        if (feature.get("properties") or {}).get("beginLife")
    ]
    marine_observed_at = f"{max(mru_versions)[:10]}T00:00:00Z" if mru_versions else generated_at
    manifest = {
        "schema_version": 1,
        "status": "complete",
        "generated_at": generated_at,
        "project_name": PROJECT_NAME,
        "scope": "belgium-and-belgian-north-sea",
        "ngi_archive_url": NGI_ARCHIVE_URL,
        "ngi_archive_sha256": sha256_file(archive_path),
        "ngi_geopackage_sha256": sha256_file(gpkg_path),
        "admin_observed_at": admin_observed_at,
        "marine_reporting_units_url": RBINS_MRU_WFS_URL,
        "marine_reporting_units_feature_count": len(reporting_units["features"]),
        "marine_observed_at": marine_observed_at,
        "msp_wfs_url": RBINS_MSP_WFS_URL,
        "msp_layer_count": len(MSP_LAYERS),
        "msp_feature_count": len(msp["features"]),
        "msp_valid_from": MSP_VALID_FROM,
        "artifacts": {
            key: {
                "filename": path.name,
                "sha256": sha256_file(path),
                "feature_count": len(artifact_payloads[key]["features"]),
            }
            for key, path in artifacts.items()
        },
    }
    write_json_atomic(output_root / "manifest.json", manifest, pretty=True)
    return artifacts, manifest


def response_data(response: requests.Response) -> Any:
    try:
        payload = response.json()
    except ValueError as exc:
        raise RuntimeError(f"GeoIntel API returned non-JSON ({response.status_code}): {response.text[:500]}") from exc
    if not response.ok:
        raise RuntimeError(f"GeoIntel API failed ({response.status_code}): {json.dumps(payload, ensure_ascii=False)[:1000]}")
    if not isinstance(payload, dict) or "data" not in payload:
        raise RuntimeError("GeoIntel API response does not use the canonical data envelope")
    return payload["data"]


def list_paginated(session: requests.Session, url: str, timeout: int) -> list[dict[str, Any]]:
    items: list[dict[str, Any]] = []
    offset = 0
    while True:
        page = response_data(session.get(url, params={"limit": 200, "offset": offset}, timeout=timeout))
        page_items = list(page.get("items") or [])
        items.extend(page_items)
        total = int(page.get("total") or len(items))
        if not page_items or len(items) >= total:
            break
        offset += len(page_items)
    return items


def _feature_geometry(payload: dict[str, Any], *, feature_id: str | None = None) -> dict[str, Any]:
    features = payload["features"]
    feature = next((item for item in features if str(item.get("id")) == feature_id), None) if feature_id else features[0]
    if not feature or not feature.get("geometry"):
        raise RuntimeError(f"Expected geometry is missing for {feature_id or payload.get('name')}")
    return feature["geometry"]


def _canonical_region_name(feature: dict[str, Any]) -> str:
    properties = feature.get("properties") or {}
    nis_code = str(properties.get("niscode") or "").lstrip("0")
    if nis_code == "2000":
        return "flanders"
    if nis_code == "3000":
        return "wallonia"
    if nis_code == "4000":
        return "brussels"
    names = " ".join(str(value) for key, value in properties.items() if "name" in key.lower() or "naam" in key.lower())
    normalized = names.casefold()
    if "vlaams" in normalized or "flam" in normalized:
        return "flanders"
    if "wallon" in normalized:
        return "wallonia"
    if "brux" in normalized or "brussel" in normalized:
        return "brussels"
    raise RuntimeError(f"Cannot map NGI region feature {feature.get('id')} to a canonical region")


def _upload_dataset(
    session: requests.Session,
    *,
    base_url: str,
    project_id: str,
    area_id: str,
    path: Path,
    source_name: str,
    reference_layer_name: str,
    coverage_zones: list[str],
    observed_at: str,
    valid_from: str,
    valid_to: str | None,
    source_version: str,
    artifact_sha256: str,
    source_url: str,
    attribution: str,
    license_note: str,
    timeout: int,
) -> dict[str, Any]:
    source_metadata = {
        "provider": source_name,
        "authority_level": "authoritative",
        "theme": THEME_BY_LAYER[reference_layer_name],
        "coverage_zones": coverage_zones,
        "reference_layer_name": reference_layer_name,
        "source_url": source_url,
        "attribution": attribution,
        "license_note": license_note,
    }
    provenance = {
        "operator_tool": "provision_belgium_north_sea_scope.py",
        "operator_explicit_fetch": True,
        "artifact_sha256": artifact_sha256,
        "source_url": source_url,
        "direct_vector_feature_write": False,
    }
    data = {
        "dataset_type": "vector",
        "source": "operator_official_import",
        "dataset_role": "reference",
        "source_name": source_name,
        "reference_layer_name": reference_layer_name,
        "source_metadata_json": json.dumps(source_metadata, ensure_ascii=False),
        "provenance_metadata_json": json.dumps(provenance, ensure_ascii=False),
        "area_id": area_id,
        "temporal_series_key": f"{source_name}:{reference_layer_name}",
        "observed_at": observed_at,
        "valid_from": valid_from,
        "temporal_granularity": "period" if valid_to else "snapshot",
        "source_version": source_version,
    }
    if valid_to:
        data["valid_to"] = valid_to
    with path.open("rb") as handle:
        response = session.post(
            f"{base_url}/api/v1/projects/{project_id}/datasets/upload",
            data=data,
            files={"file": (path.name, handle, "application/geo+json")},
            timeout=timeout,
        )
    return response_data(response)


def provision(args: argparse.Namespace, artifacts: dict[str, Path], manifest: dict[str, Any]) -> dict[str, Any]:
    base_url = args.base_url.rstrip("/")
    payloads = {key: json.loads(path.read_text(encoding="utf-8")) for key, path in artifacts.items()}
    with requests.Session() as session:
        projects = list_paginated(session, f"{base_url}/api/v1/projects", args.import_timeout)
        project = next((item for item in projects if item.get("name") == PROJECT_NAME), None)
        if project is None:
            project = response_data(
                session.post(
                    f"{base_url}/api/v1/projects",
                    json={
                        "name": PROJECT_NAME,
                        "description": (
                            "Authoritative national and maritime scope for Belgium, its territorial sea, "
                            "exclusive economic zone and continental shelf."
                        ),
                        "region": PROJECT_REGION,
                    },
                    timeout=args.import_timeout,
                )
            )
        project_id = str(project["id"])
        existing_areas = list_paginated(
            session,
            f"{base_url}/api/v1/projects/{project_id}/areas",
            args.import_timeout,
        )
        areas_by_name = {str(area.get("name")): area for area in existing_areas}

        def ensure_area(name: str, geometry: dict[str, Any]) -> dict[str, Any]:
            if name in areas_by_name:
                existing = areas_by_name[name]
                existing_geometry = existing.get("geometry")
                if not existing_geometry or not shape(existing_geometry).equals(shape(geometry)):
                    raise RuntimeError(
                        f"Persisted Area {name!r} differs from the current authoritative geometry; "
                        "provision a fresh versioned national project instead of silently mutating it"
                    )
                return existing
            area = response_data(
                session.post(
                    f"{base_url}/api/v1/projects/{project_id}/areas",
                    json={"name": name, "crs": "EPSG:4326", "geometry": geometry},
                    timeout=args.import_timeout,
                )
            )
            areas_by_name[name] = area
            return area

        land_area = ensure_area(AREA_NAMES["belgium"], _feature_geometry(payloads["belgium_land_boundary"]))
        region_features = {_canonical_region_name(feature): feature for feature in payloads["belgium_regions"]["features"]}
        for zone in ("flanders", "wallonia", "brussels"):
            ensure_area(AREA_NAMES[zone], region_features[zone]["geometry"])
        marine_features = {
            str((feature.get("properties") or {}).get("coverage_zone")): feature
            for feature in payloads["marine_legal_scopes"]["features"]
        }
        for zone in ("belgian_north_sea", "territorial_sea", "exclusive_economic_zone", "continental_shelf"):
            ensure_area(AREA_NAMES[zone], marine_features[zone]["geometry"])

        datasets = list_paginated(
            session,
            f"{base_url}/api/v1/projects/{project_id}/datasets",
            args.import_timeout,
        )
        specs = (
            (
                "belgium_land_boundary",
                "ngi_adminvector",
                "belgium_land_boundary",
                ["belgium"],
                str(land_area["id"]),
                manifest["admin_observed_at"],
                manifest["admin_observed_at"],
                None,
                manifest["admin_observed_at"][:10],
                NGI_CATALOG_URL,
                NGI_ATTRIBUTION,
                NGI_LICENSE,
            ),
            (
                "belgium_regions",
                "ngi_adminvector",
                "belgium_regions",
                ["belgium", "flanders", "wallonia", "brussels"],
                str(land_area["id"]),
                manifest["admin_observed_at"],
                manifest["admin_observed_at"],
                None,
                manifest["admin_observed_at"][:10],
                NGI_CATALOG_URL,
                NGI_ATTRIBUTION,
                NGI_LICENSE,
            ),
            (
                "belgium_provinces",
                "ngi_adminvector",
                "belgium_provinces",
                ["belgium", "flanders", "wallonia"],
                str(land_area["id"]),
                manifest["admin_observed_at"],
                manifest["admin_observed_at"],
                None,
                manifest["admin_observed_at"][:10],
                NGI_CATALOG_URL,
                NGI_ATTRIBUTION,
                NGI_LICENSE,
            ),
            (
                "belgium_municipalities",
                "ngi_adminvector",
                "belgium_municipalities",
                ["belgium", "flanders", "wallonia", "brussels"],
                str(land_area["id"]),
                manifest["admin_observed_at"],
                manifest["admin_observed_at"],
                None,
                manifest["admin_observed_at"][:10],
                NGI_CATALOG_URL,
                NGI_ATTRIBUTION,
                NGI_LICENSE,
            ),
            (
                "marine_legal_scopes",
                "rbins_marine_reporting_units",
                "marine_legal_scopes",
                [
                    "belgian_north_sea",
                    "territorial_sea",
                    "exclusive_economic_zone",
                    "continental_shelf",
                ],
                str(areas_by_name[AREA_NAMES["belgian_north_sea"]]["id"]),
                manifest["marine_observed_at"],
                manifest["marine_observed_at"],
                None,
                manifest["marine_observed_at"][:10],
                RBINS_MRU_METADATA_URL,
                "Royal Belgian Institute of Natural Sciences (RBINS), BMDC",
                "See source metadata",
            ),
            (
                "marine_spatial_plan_2026",
                "rbins_msp_2026",
                "marine_spatial_plan_2026",
                [
                    "belgian_north_sea",
                    "territorial_sea",
                    "exclusive_economic_zone",
                    "continental_shelf",
                ],
                str(areas_by_name[AREA_NAMES["belgian_north_sea"]]["id"]),
                MSP_VALID_FROM,
                MSP_VALID_FROM,
                "2034-12-31T23:59:59Z",
                "2026-2034",
                RBINS_MSP_SOURCE_URL,
                "Belgian federal Marine Environment service and RBINS",
                "See source metadata",
            ),
        )
        persisted = []
        for (
            artifact_key,
            source_name,
            layer_name,
            coverage_zones,
            area_id,
            observed_at,
            valid_from,
            valid_to,
            source_version,
            source_url,
            attribution,
            license_note,
        ) in specs:
            checksum = manifest["artifacts"][artifact_key]["sha256"]
            existing = next(
                (
                    dataset
                    for dataset in datasets
                    if dataset.get("source_name") == source_name
                    and dataset.get("reference_layer_name") == layer_name
                    and dataset.get("source_version") == source_version
                ),
                None,
            )
            if existing:
                persisted_checksum = (existing.get("provenance_metadata") or {}).get("artifact_sha256")
                if persisted_checksum != checksum:
                    raise RuntimeError(
                        f"Existing immutable dataset {layer_name} has checksum {persisted_checksum}, expected {checksum}"
                    )
                persisted.append(existing)
                continue
            created = _upload_dataset(
                session,
                base_url=base_url,
                project_id=project_id,
                area_id=area_id,
                path=artifacts[artifact_key],
                source_name=source_name,
                reference_layer_name=layer_name,
                coverage_zones=coverage_zones,
                observed_at=observed_at,
                valid_from=valid_from,
                valid_to=valid_to,
                source_version=source_version,
                artifact_sha256=checksum,
                source_url=source_url,
                attribution=attribution,
                license_note=license_note,
                timeout=args.import_timeout,
            )
            persisted.append(created)

    return {
        "project_id": project_id,
        "project_name": PROJECT_NAME,
        "area_count": len(AREA_NAMES),
        "dataset_ids": [str(dataset["id"]) for dataset in persisted],
    }


def main() -> int:
    args = parse_args()
    try:
        artifacts, manifest = prepare_artifacts(args)
        workspace = None if args.fetch_only else provision(args, artifacts, manifest)
    except (OSError, RuntimeError, ValueError, KeyError, requests.RequestException, zipfile.BadZipFile) as exc:
        print(json.dumps({"status": "error", "message": str(exc)}, ensure_ascii=False), file=sys.stderr)
        return 1
    print(
        json.dumps(
            {
                "status": "ok",
                "mode": "fetch_only" if args.fetch_only else "provisioned",
                "scope": manifest["scope"],
                "manifest_path": str(args.output_root / "manifest.json"),
                "artifact_count": len(artifacts),
                "msp_feature_count": manifest["msp_feature_count"],
                "workspace": workspace,
            },
            ensure_ascii=False,
            indent=2,
        )
    )
    return 0


if __name__ == "__main__":
    sys.exit(main())