From 16be7564f377640785525c40344c7baafda8e9fc Mon Sep 17 00:00:00 2001 From: Codex Date: Wed, 15 Jul 2026 13:20:13 +0200 Subject: [PATCH] feat: add governed Mol nature value layer --- CHANGELOG.md | 14 + backend/README.md | 15 + .../app/services/vector_feature_service.py | 42 +- .../tests/test_sprint204_bwk_natura2000.py | 270 ++++++ deploy/unraid/Dockerfile.all-in-one | 1 + docs/API_CONTRACTS.md | 8 + docs/DATA_SOURCES.md | 21 + docs/DATA_SPECIFICATION.md | 15 + docs/STORAGE_ARCHITECTURE.md | 7 + frontend/README.md | 6 + .../datasets/SourceCatalogPanel.tsx | 12 +- frontend/src/components/map/MapWorkspace.tsx | 11 +- frontend/src/styles/app.css | 1 + scripts/README.md | 19 + scripts/provision_mol_bwk_natura2000.py | 815 ++++++++++++++++++ scripts/run_readiness_check.sh | 1 + 16 files changed, 1251 insertions(+), 7 deletions(-) create mode 100644 backend/tests/test_sprint204_bwk_natura2000.py create mode 100644 scripts/provision_mol_bwk_natura2000.py diff --git a/CHANGELOG.md b/CHANGELOG.md index ed10bb6b..7ccb3211 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,20 @@ # Changelog +## Sprint 204 Governed BWK and Natura 2000 for Mol (2026-07-15) + +- Added an explicit operator for the official INBO BWK/Natura 2000 state 2025 + WFS with complete link/start-index pagination, retained raw page checksums and + exact Mol clipping in EPSG:31370. +- Preserved original BWK evaluation, mapping-unit, origin, habitat, percentage + and habitat-origin fields while adding readable, source-faithful properties. +- Added filtered PostGIS selection metrics for each official BWK value class, + Natura 2000 habitat shares, regional biotopes and uncertain habitat gaps. +- Added Nature Value as a seventh map theme and moved BWK/Natura 2000 from the + follow-up catalogue to loaded evidence only after a real Dataset exists. +- Kept the 2025 edition honest: it is a map state rather than one uniform 2025 + survey, and PHAB-based partial-selection hectares remain labelled estimates. + ## Sprint 203 Governed hydrology and historical imagery (2026-07-15) - Added an explicit Waterinfo KiWIS operator for annual water-level and diff --git a/backend/README.md b/backend/README.md index 567ef067..834b1fc2 100644 --- a/backend/README.md +++ b/backend/README.md @@ -1158,3 +1158,18 @@ observations through the canonical dataset upload API. Every station has its own temporal-series key. Water levels and discharges remain Point measurements; they are never averaged across stations or presented as municipal water volume. Use `--fetch-only` to prepare and audit artifacts without persistence. + +## BWK/Natura 2000 state 2025 + +Run the governed Mol operator after the regional workspace and Mol Area exist: + +```bash +docker exec geointel python /app/scripts/provision_mol_bwk_natura2000.py +``` + +The command fetches the official INBO WFS, retains raw checksummed pages, +clips in EPSG:31370 and imports through DatasetService. `--fetch-only` builds +evidence without persistence. A conflicting checksum for an already persisted +state-2025 Mol Dataset fails closed instead of creating a silent replacement. +PostGIS selection summaries keep BWK value classes separate and label +PHAB-derived habitat hectares as estimates. diff --git a/backend/app/services/vector_feature_service.py b/backend/app/services/vector_feature_service.py index 9b6e2fad..43d1596d 100644 --- a/backend/app/services/vector_feature_service.py +++ b/backend/app/services/vector_feature_service.py @@ -24,6 +24,7 @@ FULL_AREA_CLIPPED_OPERATOR_TOOLS = { "provision_regional_grb_buildings.py", "provision_regional_grb_context.py", "provision_waterinfo_station_history.py", + "provision_mol_bwk_natura2000.py", } @@ -84,6 +85,7 @@ SEMANTIC_SELECTION_METRICS: dict[str, tuple[dict[str, Any], ...]] = { "warning": "GRB-percelen zijn een grafische referentie en vormen geen juridische grensopmeting.", }, ), + "nature_value": (), } SEMANTIC_COUNT_LABELS = { @@ -93,6 +95,7 @@ SEMANTIC_COUNT_LABELS = { "water": "Waterobjecten", "roads": "Wegsegmenten", "parcels": "Percelen", + "nature_value": "BWK-kaartvlakken", } @@ -114,6 +117,10 @@ class VectorFeatureService: "waterways": "water", "road": "roads", "parcel": "parcels", + "nature": "nature_value", + "biodiversity": "nature_value", + "bwk": "nature_value", + "natura2000": "nature_value", } for candidate in candidates: if not isinstance(candidate, str) or not candidate.strip(): @@ -356,6 +363,17 @@ class VectorFeatureService: primary_config = semantic_metrics[0] metric_configs = [primary_config] + configured_metrics = source_metadata.get("selection_metrics") + if isinstance(configured_metrics, list): + existing_metric_keys = {str(primary_config.get("metric_key") or "")} + for configured_item in configured_metrics: + if not isinstance(configured_item, dict): + continue + metric_key = str(configured_item.get("metric_key") or "").strip() + if not metric_key or metric_key in existing_metric_keys: + continue + metric_configs.append(dict(configured_item)) + existing_metric_keys.add(metric_key) for semantic_metric in semantic_metrics: signature = (semantic_metric["method"], semantic_metric["unit"]) existing = { @@ -419,6 +437,20 @@ class VectorFeatureService: metric_filter = selection_filter if dimension in {1, 2}: metric_filter += (func.ST_Dimension(VectorFeature.geometry) == int(dimension),) + filter_property = str(config.get("filter_property") or "").strip() + filter_values = config.get("filter_values") + if filter_property: + if not isinstance(filter_values, list) or not filter_values: + raise AppError( + code="INVALID_SELECTION_AGGREGATION", + message="Dataset selection metric filter requires one or more values", + details={"dataset_id": str(dataset.id), "filter_property": filter_property}, + status_code=500, + ) + normalized_filter_values = [str(value) for value in filter_values] + metric_filter += ( + VectorFeature.properties_json.op("->>")(filter_property).in_(normalized_filter_values), + ) if method == "intersection_area": measured_geometry = ( @@ -461,7 +493,7 @@ class VectorFeatureService: aggregate_function = func.avg if method == "mean" else func.sum aggregate_value = ( db.query(func.coalesce(aggregate_function(value_expression), 0.0)) - .filter(*selection_filter) + .filter(*metric_filter) .filter(VectorFeature.properties_json.op("->>")(property_name).isnot(None)) .scalar() ) @@ -469,16 +501,16 @@ class VectorFeatureService: if method == "area_weighted_sum" and not full_dataset_area: partial_feature_count = ( db.query(func.count(VectorFeature.id)) - .filter(*selection_filter) + .filter(*metric_filter) .filter(coverage_ratio < 0.999999) .scalar() ) - is_estimate = bool(partial_feature_count) + is_estimate = bool(config.get("is_estimate", False)) or bool(partial_feature_count) if not is_estimate and config.get("warning_only_when_estimate", True): warning = None elif method == "area_weighted_sum": - is_estimate = False - if config.get("warning_only_when_estimate", True): + is_estimate = bool(config.get("is_estimate", False)) + if not is_estimate and config.get("warning_only_when_estimate", True): warning = None elif method != "feature_count": raise AppError( diff --git a/backend/tests/test_sprint204_bwk_natura2000.py b/backend/tests/test_sprint204_bwk_natura2000.py new file mode 100644 index 00000000..04065954 --- /dev/null +++ b/backend/tests/test_sprint204_bwk_natura2000.py @@ -0,0 +1,270 @@ +from __future__ import annotations + +import importlib.util +import json +import sys +from pathlib import Path +from uuid import uuid4 + +import pytest +from shapely.geometry import box, shape +from shapely.ops import transform as transform_geometry + +from app.models import Dataset +from app.schemas.operations import VectorSelectionSummary +from app.services.vector_feature_service import VectorFeatureService + + +ROOT = Path(__file__).resolve().parents[2] +BBOX = {"min_x": 5.0, "min_y": 51.1, "max_x": 5.2, "max_y": 51.3, "crs": "EPSG:4326"} + + +def load_operator(): + script_path = ROOT / "scripts" / "provision_mol_bwk_natura2000.py" + spec = importlib.util.spec_from_file_location("bwk_natura2000_operator", script_path) + assert spec is not None + assert spec.loader is not None + module = importlib.util.module_from_spec(spec) + sys.modules[spec.name] = module + spec.loader.exec_module(module) + return module + + +class FakeResponse: + ok = True + status_code = 200 + text = "" + + def __init__(self, payload: dict, url: str): + self.payload = payload + self.url = url + self.content = json.dumps(payload).encode("utf-8") + + def raise_for_status(self): + return None + + def json(self): + return self.payload + + +class FakeSession: + def __init__(self, responses: list[FakeResponse]): + self.responses = iter(responses) + self.calls: list[tuple[str, dict | None]] = [] + + def get(self, url, *, params=None, timeout): # noqa: ANN001, ARG002 + self.calls.append((url, params)) + return next(self.responses) + + +class ScalarQuery: + def __init__(self, value: float): + self.value = value + + def filter(self, *args): # noqa: ANN002, ARG002 + return self + + def scalar(self): + return self.value + + +class SequenceScalarSession: + def __init__(self, values: list[float]): + self.values = iter(values) + + def query(self, *args): # noqa: ANN002, ARG002 + return ScalarQuery(next(self.values)) + + +def test_official_bwk_contract_and_class_labels_are_fixed() -> None: + module = load_operator() + + assert module.WFS_URL == "https://geo.api.vlaanderen.be/BWK/wfs" + assert module.TYPE_NAME == "BWK:Bwkhab" + assert module.SOURCE_VERSION == "2025" + assert module.ATTRIBUTION == "Bron: INBO" + assert module.EVALUATION_LABELS == { + "z": "Biologisch zeer waardevol", + "w": "Biologisch waardevol", + "m": "Biologisch minder waardevol", + "wz": "Complex van waardevolle en zeer waardevolle elementen", + "mwz": "Complex van minder waardevolle, waardevolle en zeer waardevolle elementen", + "mz": "Complex van minder waardevolle en zeer waardevolle elementen", + "mw": "Complex van minder waardevolle en waardevolle elementen", + } + + +def test_wfs_pagination_follows_server_next_links() -> None: + module = load_operator() + page_one = { + "type": "FeatureCollection", + "features": [{"id": "one"}], + "numberReturned": 1, + "links": [{"rel": "next", "href": "https://geo.api.vlaanderen.be/BWK/wfs?STARTINDEX=1"}], + } + page_two = {"type": "FeatureCollection", "features": [], "numberReturned": 0, "links": []} + session = FakeSession( + [ + FakeResponse(page_one, "https://geo.api.vlaanderen.be/BWK/wfs?first"), + FakeResponse(page_two, "https://geo.api.vlaanderen.be/BWK/wfs?STARTINDEX=1"), + ] + ) + + pages = list(module.iter_wfs_pages(session, (5.0, 51.1, 5.2, 51.3), page_limit=1000, timeout=30)) + + assert len(pages) == 2 + assert session.calls[0][1]["sortBy"] == "UIDN" + assert session.calls[0][1]["srsName"] == "EPSG:4326" + assert session.calls[1] == ("https://geo.api.vlaanderen.be/BWK/wfs?STARTINDEX=1", None) + + +def test_wfs_page_limit_without_next_link_uses_controlled_start_index_fallback() -> None: + module = load_operator() + session = FakeSession( + [ + FakeResponse( + {"type": "FeatureCollection", "features": [{"id": "one"}], "numberReturned": 1}, + "https://geo.api.vlaanderen.be/BWK/wfs", + ), + FakeResponse( + {"type": "FeatureCollection", "features": [], "numberReturned": 0}, + "https://geo.api.vlaanderen.be/BWK/wfs?startIndex=1", + ), + ] + ) + + pages = list(module.iter_wfs_pages(session, (5.0, 51.1, 5.2, 51.3), page_limit=1, timeout=30)) + + assert len(pages) == 2 + assert session.calls[1][1]["startIndex"] == "1" + + +def test_feature_is_clipped_in_lambert72_and_keeps_bwk_habitat_provenance() -> None: + module = load_operator() + boundary_wgs84 = box(5.10, 51.20, 5.11, 51.21) + boundary_lambert72 = transform_geometry(module.TO_LAMBERT72.transform, boundary_wgs84) + feature = { + "type": "Feature", + "id": "Bwkhab.42", + "geometry": mapping_box(5.095, 51.195, 5.105, 51.205), + "properties": { + "UIDN": 42, + "EVAL": "wz", + "BWKLABEL": "qb + qs", + "EENH1": "qb", + "EENH2": "qs", + "HERK": "225", + "HAB1": "9190", + "PHAB1": 60, + "HAB2": "rbbppm", + "PHAB2": 30, + "HAB3": "gh", + "PHAB3": 10, + "HABLEGENDE": "phab", + "HERKHAB": "225", + "HERKPHAB": "a", + }, + } + + normalized, was_clipped = module.normalize_feature(feature, boundary_lambert72) + + assert normalized is not None + assert was_clipped is True + normalized_geometry = shape(normalized["geometry"]) + assert normalized_geometry.difference(boundary_wgs84.buffer(1e-7)).area < 1e-12 + properties = normalized["properties"] + assert properties["bwk_evaluation_code"] == "wz" + assert properties["bwk_evaluation_label"].startswith("Complex van waardevolle") + assert properties["natura2000_codes"] == "9190" + assert properties["regional_biotope_codes"] == "rbbppm" + assert properties["natura2000_area_ha"] == pytest.approx(properties["clipped_area_ha"] * 0.6) + assert properties["regional_biotope_area_ha"] == pytest.approx(properties["clipped_area_ha"] * 0.3) + assert properties["habitat_share_origin_code"] == "a" + + +def mapping_box(min_x: float, min_y: float, max_x: float, max_y: float) -> dict: + return { + "type": "Polygon", + "coordinates": [[ + [min_x, min_y], + [max_x, min_y], + [max_x, max_y], + [min_x, max_y], + [min_x, min_y], + ]], + } + + +def test_uncertain_habitat_status_is_not_presented_as_confirmed_habitat() -> None: + module = load_operator() + + entries, natura_share, regional_share, uncertain_share = module.habitat_breakdown( + {"HAB1": "gh", "PHAB1": 100, "HABLEGENDE": "ohab"} + ) + + assert entries == [{"code": "gh", "share_percent": 100.0}] + assert natura_share == 0 + assert regional_share == 0 + assert uncertain_share == 100 + + +def test_nature_value_summary_returns_separate_official_classes_and_habitat_metrics() -> None: + module = load_operator() + dataset = Dataset( + id=uuid4(), + project_id=uuid4(), + name="bwk_natura2000_2025_mol.geojson", + dataset_type="vector", + dataset_role="reference", + source_name="inbo_bwk_natura2000", + reference_layer_name="nature_value", + source_metadata={ + "theme": "nature_value", + "semantic_metrics": False, + "selection_aggregation": { + "metric_key": "bwk_mapped_area", + "method": "intersection_area", + "label": "BWK-gekarteerde oppervlakte", + "unit": "ha", + "geometry_dimension": 2, + }, + "selection_metrics": module.selection_metrics(), + }, + ) + + result = VectorFeatureService.summarize_features_by_bbox( + SequenceScalarSession([100_000.0, 10_000.0, 20_000.0, 30_000.0, 40_000.0, 5.5, 2.5, 1.5]), + dataset=dataset, + bbox=BBOX, + total_feature_count=125, + full_dataset_area=True, + ) + + assert result["primary_metric_key"] == "bwk_mapped_area" + assert result["metric_value"] == 10.0 + metrics = {item["metric_key"]: item for item in result["metrics"]} + assert metrics["bwk_very_valuable_area"]["metric_value"] == 1.0 + assert metrics["bwk_valuable_area"]["metric_value"] == 2.0 + assert metrics["bwk_less_valuable_area"]["metric_value"] == 3.0 + assert metrics["bwk_mixed_value_area"]["metric_value"] == 4.0 + assert metrics["natura2000_area"]["metric_value"] == 5.5 + assert metrics["natura2000_area"]["is_estimate"] is True + assert metrics["regional_biotope_area"]["metric_value"] == 2.5 + assert metrics["uncertain_habitat_area"]["metric_value"] == 1.5 + assert metrics["feature_count"]["metric_value"] == 125 + VectorSelectionSummary(**result) + + +def test_operator_is_packaged_readiness_checked_and_wired_to_map() -> None: + dockerfile = (ROOT / "deploy" / "unraid" / "Dockerfile.all-in-one").read_text(encoding="utf-8") + readiness = (ROOT / "scripts" / "run_readiness_check.sh").read_text(encoding="utf-8") + service = (ROOT / "backend" / "app" / "services" / "vector_feature_service.py").read_text(encoding="utf-8") + map_workspace = (ROOT / "frontend" / "src" / "components" / "map" / "MapWorkspace.tsx").read_text(encoding="utf-8") + source_catalog = (ROOT / "frontend" / "src" / "components" / "datasets" / "SourceCatalogPanel.tsx").read_text(encoding="utf-8") + + assert "COPY scripts/provision_mol_bwk_natura2000.py" in dockerfile + assert "py_compile scripts/provision_mol_bwk_natura2000.py" in readiness + assert '"provision_mol_bwk_natura2000.py"' in service + assert "id: 'nature_value'" in map_workspace + assert "Natuurwaarde" in map_workspace + assert "source.key === 'bwk'" in source_catalog diff --git a/deploy/unraid/Dockerfile.all-in-one b/deploy/unraid/Dockerfile.all-in-one index ad1ca8f2..a12a8100 100644 --- a/deploy/unraid/Dockerfile.all-in-one +++ b/deploy/unraid/Dockerfile.all-in-one @@ -78,6 +78,7 @@ COPY scripts/provision_mol_population_history.py /app/scripts/provision_mol_popu COPY scripts/provision_mol_historical_landuse.py /app/scripts/provision_mol_historical_landuse.py COPY scripts/provision_official_landuse_timeseries.py /app/scripts/provision_official_landuse_timeseries.py COPY scripts/provision_waterinfo_station_history.py /app/scripts/provision_waterinfo_station_history.py +COPY scripts/provision_mol_bwk_natura2000.py /app/scripts/provision_mol_bwk_natura2000.py COPY scripts/provision_regional_timeseries.py /app/scripts/provision_regional_timeseries.py COPY scripts/geographic_scopes.py /app/scripts/geographic_scopes.py COPY scripts/provision_geographic_scope.py /app/scripts/provision_geographic_scope.py diff --git a/docs/API_CONTRACTS.md b/docs/API_CONTRACTS.md index 6891ee37..75b379e5 100644 --- a/docs/API_CONTRACTS.md +++ b/docs/API_CONTRACTS.md @@ -496,6 +496,11 @@ Rules: - Results are generated from persisted PostGIS `vector_features`, not from client-side map data. - `feature_count` is the number of GeoJSON features returned in the bounded preview. `total_feature_count` is the exact number of persisted rows intersecting the requested bbox or persisted Area geometry. - `summary` keeps one backwards-compatible primary metric and exposes all relevant measurements in `metrics`. Known themes use metric PostGIS calculations: building/forest/water/parcel surfaces in hectares, road and watercourse lengths in kilometres, population in inhabitants and intersecting feature counts as supporting evidence. +- Governed datasets may declare additional `source_metadata.selection_metrics`. + Each metric may constrain one persisted feature-property to an explicit value + allowlist before the same PostGIS aggregation runs. The BWK/Natura 2000 + dataset uses this only for official `EVAL` classes; it does not collapse mixed + classes into a made-up score. - Area and length calculations transform geometry to Belgian Lambert 72 (`EPSG:31370`); they are never calculated in geographic degrees. - Water volume is not inferred from 2D GRB geometry. It remains unavailable until a source provides reliable depth or bathymetry with compatible spatial coverage and provenance. - The response is capped by `limit` and returns `truncated=true` when `total_feature_count` exceeds the returned preview. @@ -1582,6 +1587,9 @@ aggregation method, feature count, estimate status and an optional warning. Supported PostGIS aggregations are feature count, intersection area, intersection length, numeric sum and area-weighted numeric sum. Area and length are measured after transformation to EPSG:31370. +Configured supplemental metrics can use `filter_property` plus +`filter_values`. Filters are server-owned dataset metadata, not arbitrary +client SQL or request expressions. ### PATCH `/api/v1/projects/{project_id}/datasets/{dataset_id}/temporal` diff --git a/docs/DATA_SOURCES.md b/docs/DATA_SOURCES.md index 3303fc84..307682fa 100644 --- a/docs/DATA_SOURCES.md +++ b/docs/DATA_SOURCES.md @@ -227,6 +227,27 @@ until a governed operator import, provenance record and validation pass exist. LiDAR, suitable for elevation, slope and drainage. It does not provide water depth. +## Governed BWK and Natura 2000 state 2025 + +`scripts/provision_mol_bwk_natura2000.py` uses the official production WFS +`https://geo.api.vlaanderen.be/BWK/wfs`, feature type `BWK:Bwkhab`. It follows +server `next` links and uses WFS `startIndex` only when a full page omits that +link. A short final page is required for completeness. Raw response pages, +request URLs and SHA256 checksums are retained. + +The source bbox is reduced to the exact persisted Mol Area in Belgian Lambert +72. Only valid polygonal intersections are emitted, transformed back to +EPSG:4326 and imported through DatasetService with `source_name` equal to +`inbo_bwk_natura2000` and `reference_layer_name` equal to `nature_value`. +Required attribution is `Bron: INBO`. + +`EVAL` classes are reported separately. Natura 2000 and regionally important +biotope hectares use official `PHAB` shares; uncertain `ohab` polygons remain +knowledge gaps. The official report warns that parts of the 2025 edition, +especially outside Habitats Directive areas, can still be based on older +mapping. GeoIntel therefore preserves origin fields and never presents this +single edition as a historical trend. + ## Governed Waterinfo station history `scripts/provision_waterinfo_station_history.py` uses the public Waterinfo diff --git a/docs/DATA_SPECIFICATION.md b/docs/DATA_SPECIFICATION.md index 24eae66d..46856531 100644 --- a/docs/DATA_SPECIFICATION.md +++ b/docs/DATA_SPECIFICATION.md @@ -213,6 +213,21 @@ series may not be merged into one area-wide value. Point water level and discharge may not be converted to water volume without governed compatible depth/profile data. +### BWK and Natura 2000 polygons + +The governed state-2025 import uses `BWK:Bwkhab` polygons. Geometry is fetched +in EPSG:4326, validated and clipped against the persisted Mol boundary in +EPSG:31370, then persisted in EPSG:4326. Required retained fields include +`EVAL`, `EENH1..8`, `V1..3`, `HERK`, `BWKLABEL`, `HAB1..5`, `PHAB1..5`, +`HERKHAB`, `HERKPHAB` and `HABLEGENDE`. + +BWK value classes remain independent. Natura 2000 codes, regionally important +biotope codes and uncertain `ohab` knowledge gaps are also independent. PHAB +shares can be theoretical source allocations; derived hectares therefore keep +an estimate warning, especially for partial rectangle intersections. State +2025 is a product edition and must not be presented as a uniform 2025 field +survey or an annual time series. + ### Vector Ondersteund: diff --git a/docs/STORAGE_ARCHITECTURE.md b/docs/STORAGE_ARCHITECTURE.md index 42ff547b..b4a4976b 100644 --- a/docs/STORAGE_ARCHITECTURE.md +++ b/docs/STORAGE_ARCHITECTURE.md @@ -109,6 +109,13 @@ product/layer, request/spatial hash, temporal validity and limitations are held in source/provenance metadata. Browser PNG rendering is derived on request and does not replace the stored GeoTIFF. +BWK/Natura 2000 evidence lives under +`storage/operator-evidence/bwk-natura2000-2025/mol/`. The `raw/` directory +contains immutable WFS pages; the adjacent manifest records their URLs, +checksums, boundary checksum, page/feature counts, clipping diagnostics and +official class totals. The final GeoJSON is the retained source artifact while +`datasets` and `vector_features` remain the canonical queryable PostGIS state. + Offline demo export artifacts can be inspected and cleaned with: ```bash diff --git a/frontend/README.md b/frontend/README.md index 2679e7a6..206a6cc6 100644 --- a/frontend/README.md +++ b/frontend/README.md @@ -470,6 +470,12 @@ The workbench shell clamps page-level horizontal overflow on mobile while keepin The QA/QC result list includes client-side search, status and check-type filters plus a latest-results cap so long-lived demo projects remain scan-friendly without changing the API response shape. +The current-state Map explorer includes `Natuurwaarde` when the governed INBO +BWK/Natura 2000 state-2025 Dataset is loaded. It shows exact mapped BWK area, +separate official value classes, PHAB-derived Natura 2000 and regional-biotope +hectares, uncertain habitat knowledge gaps and source provenance. This single +edition is not exposed as a fabricated annual series in Evolution. + ## Useful repository scripts - `bash scripts/frontend_install.sh` diff --git a/frontend/src/components/datasets/SourceCatalogPanel.tsx b/frontend/src/components/datasets/SourceCatalogPanel.tsx index fcb8cfe1..96adf0c8 100644 --- a/frontend/src/components/datasets/SourceCatalogPanel.tsx +++ b/frontend/src/components/datasets/SourceCatalogPanel.tsx @@ -8,6 +8,7 @@ const THEME_LABELS: Record = { buildings: 'Bebouwing', population: 'Bevolking', forest: 'Bos', + nature_value: 'Natuurwaarde', water: 'Water', roads: 'Wegen en transport', parcels: 'Percelen', @@ -103,9 +104,11 @@ export function SourceCatalogPanel({ datasets }: SourceCatalogPanelProps): JSX.E dataset.source_name === 'digitaal_vlaanderen_orthophoto' && String(dataset.source_metadata?.['product_key'] ?? 'most_recent') !== 'most_recent', ) + const bwkDatasets = ready.filter((dataset) => dataset.source_name === 'inbo_bwk_natura2000') const pendingSources = AVAILABLE_SOURCES.filter((source) => { if (source.key === 'waterinfo') return waterinfoDatasets.length === 0 if (source.key === 'historical_orthophoto') return historicalOrthophotos.length === 0 + if (source.key === 'bwk') return bwkDatasets.length === 0 return true }) const themes = Object.keys(THEME_LABELS).map((theme) => { @@ -165,7 +168,7 @@ export function SourceCatalogPanel({ datasets }: SourceCatalogPanelProps): JSX.E ))} - {waterinfoDatasets.length > 0 || historicalOrthophotos.length > 0 ? ( + {waterinfoDatasets.length > 0 || historicalOrthophotos.length > 0 || bwkDatasets.length > 0 ? (
{waterinfoDatasets.length > 0 ? (
@@ -181,6 +184,13 @@ export function SourceCatalogPanel({ datasets }: SourceCatalogPanelProps): JSX.E

Officiële jaargangen en periodes zijn via de kaart beschikbaar zonder actuele GRB-validatie.

) : null} + {bwkDatasets.length > 0 ? ( +
+ BWK / Natura 2000 + Toestand 2025 · {bwkDatasets.reduce((total, dataset) => total + (dataset.feature_count ?? 0), 0).toLocaleString('nl-BE')} kaartvlakken +

Biologische waardering en habitataandelen blijven afzonderlijke, brongetrouwe metingen met herkomstinformatie.

+
+ ) : null}
) : null} diff --git a/frontend/src/components/map/MapWorkspace.tsx b/frontend/src/components/map/MapWorkspace.tsx index bcd4e44f..ffa400d2 100644 --- a/frontend/src/components/map/MapWorkspace.tsx +++ b/frontend/src/components/map/MapWorkspace.tsx @@ -13,7 +13,7 @@ const EMPTY_TEMPORAL_SERIES: DatasetCreateResponse[] = [] const MOL_PROJECT_NAME = 'Mol Municipality Workbench' const KEMPEN_PROJECT_NAME = 'Kempen Regional Workbench' -type DataThemeId = 'buildings' | 'population' | 'forest' | 'water' | 'roads' | 'parcels' +type DataThemeId = 'buildings' | 'population' | 'forest' | 'nature_value' | 'water' | 'roads' | 'parcels' interface DataTheme { id: DataThemeId @@ -51,6 +51,13 @@ const DATA_THEMES: DataTheme[] = [ description: 'Bos, natuur en groenbedekking uit een ingeladen vectorbron.', tokens: ['forest', 'forestry', 'woodland', 'bos', 'groen', 'vegetation'], }, + { + id: 'nature_value', + label: 'Natuurwaarde', + shortLabel: 'BWK-oppervlakte', + description: 'Biologische waardering, Natura 2000-habitat en regionaal belangrijke biotopen uit de BWK.', + tokens: ['nature_value', 'nature value', 'natuurwaarde', 'bwk', 'natura2000', 'natura 2000', 'biodiversity'], + }, { id: 'water', label: 'Water', @@ -78,6 +85,7 @@ const DATA_THEME_MAP_STYLES: Record buildings: { fill: '#d45f3d', line: '#9f3e24' }, population: { fill: '#7559a6', line: '#5b3f88' }, forest: { fill: '#347950', line: '#225f3b' }, + nature_value: { fill: '#9a4f64', line: '#74364a' }, water: { fill: '#2676a8', line: '#155b85' }, roads: { fill: '#6b7280', line: '#4b5563' }, parcels: { fill: '#a7792f', line: '#7d571f' }, @@ -111,6 +119,7 @@ function pickThemeDataset(datasets: DatasetCreateResponse[], theme: DataTheme): (dataset.reference_layer_name && theme.tokens.includes(dataset.reference_layer_name.toLowerCase()) ? 1_000_000 : 0) + (dataset.source_name === 'grb' ? 100_000 : 0) + (dataset.source_name === 'department_omgeving_land_use' ? 90_000 : 0) + + (dataset.source_name === 'inbo_bwk_natura2000' ? 95_000 : 0) + (dataset.dataset_role === 'reference' ? 10_000 : 0) + (dataset.observed_at ? new Date(dataset.observed_at).getTime() / 100_000_000 : 0) + (dataset.feature_count ?? dataset.vector_summary?.feature_count ?? 0) diff --git a/frontend/src/styles/app.css b/frontend/src/styles/app.css index 6db34373..080c46ce 100644 --- a/frontend/src/styles/app.css +++ b/frontend/src/styles/app.css @@ -5753,6 +5753,7 @@ section { .geo-theme-symbol-buildings { background: #d45f3d; } .geo-theme-symbol-population { background: #7559a6; } .geo-theme-symbol-forest { background: #347950; } +.geo-theme-symbol-nature_value { background: #9a4f64; } .geo-theme-symbol-water { background: #2676a8; } .geo-theme-symbol-roads { background: #6b7280; } .geo-theme-symbol-parcels { background: #a7792f; } diff --git a/scripts/README.md b/scripts/README.md index 4058d225..3d22bc61 100644 --- a/scripts/README.md +++ b/scripts/README.md @@ -1438,6 +1438,25 @@ missing discharge series is reported without synthesizing values. Different stations remain separate temporal series and may not be treated as area-wide water level or volume. +## Governed BWK/Natura 2000 state 2025 for Mol + +Build and validate raw WFS evidence without importing: + +```bash +docker exec geointel python /app/scripts/provision_mol_bwk_natura2000.py --fetch-only +``` + +Import after reviewing the manifest: + +```bash +docker exec geointel python /app/scripts/provision_mol_bwk_natura2000.py +``` + +The operator uses the official `BWK:Bwkhab` layer, follows complete WFS +pagination, retains raw JSON/checksums, clips polygon geometry in EPSG:31370 +and writes only through the Dataset upload flow. `--force` refetches source +evidence but never silently replaces a conflicting persisted 2025 artifact. + ## Tower deployment Push the local branch to Gitea, then rebuild the Unraid/Tower Docker runtime: diff --git a/scripts/provision_mol_bwk_natura2000.py b/scripts/provision_mol_bwk_natura2000.py new file mode 100644 index 00000000..7f277a92 --- /dev/null +++ b/scripts/provision_mol_bwk_natura2000.py @@ -0,0 +1,815 @@ +"""Provision the official 2025 BWK/Natura 2000 dataset for Mol. + +The operator follows every WFS page, retains raw checksummed responses, clips +polygon geometry to the persisted Mol Area in EPSG:31370 and imports the final +GeoJSON through the canonical DatasetService upload route. It never writes +directly to vector_features. +""" + +from __future__ import annotations + +import argparse +import hashlib +import json +import os +import re +import sys +from collections import defaultdict +from datetime import datetime, timezone +from pathlib import Path +from typing import Any, Iterable + +import requests +from pyproj import Transformer +from requests.adapters import HTTPAdapter +from shapely.geometry import MultiPolygon, Polygon, mapping, shape +from shapely.ops import transform as transform_geometry +from shapely.ops import unary_union +from shapely.validation import make_valid +from urllib3.util.retry import Retry + + +WFS_URL = "https://geo.api.vlaanderen.be/BWK/wfs" +CATALOG_URL = ( + "https://www.vlaanderen.be/datavindplaats/catalogus/" + "biologische-waarderingskaart-en-natura-2000-habitatkaart-toestand-2025" +) +REPORT_URL = "https://doi.org/10.21436/inbor.129502912" +TYPE_NAME = "BWK:Bwkhab" +SOURCE_VERSION = "2025" +PUBLICATION_DATE = "2025-12-10T00:00:00Z" +ATTRIBUTION = "Bron: INBO" +DEFAULT_API_URL = "http://127.0.0.1:8000" +DEFAULT_OUTPUT_DIR = "/app/storage/operator-evidence/bwk-natura2000-2025/mol" +DEFAULT_PROJECT_NAME = "Kempen Regional Workbench" +DEFAULT_AREA_NAME = "Gemeente Mol - officiele grens" +DATASET_FILENAME = "bwk_natura2000_2025_mol.geojson" +MANIFEST_FILENAME = "bwk_natura2000_2025_mol.manifest.json" +SCHEMA_VERSION = 1 + +EVALUATION_LABELS = { + "z": "Biologisch zeer waardevol", + "w": "Biologisch waardevol", + "m": "Biologisch minder waardevol", + "wz": "Complex van waardevolle en zeer waardevolle elementen", + "mwz": "Complex van minder waardevolle, waardevolle en zeer waardevolle elementen", + "mz": "Complex van minder waardevolle en zeer waardevolle elementen", + "mw": "Complex van minder waardevolle en waardevolle elementen", +} + +HABITAT_STATUS_LABELS = { + "gh": "Geen Natura 2000-habitat aangeduid", + "hab": "Volledig habitatwaardig kaartvlak", + "phab": "Gedeeltelijk habitatwaardig kaartvlak", + "ohab": "Onzeker habitat of kennislacune", +} + +TO_LAMBERT72 = Transformer.from_crs("EPSG:4326", "EPSG:31370", always_xy=True) +TO_WGS84 = Transformer.from_crs("EPSG:31370", "EPSG:4326", always_xy=True) + + +def parse_args() -> argparse.Namespace: + parser = argparse.ArgumentParser(description="Provision official BWK/Natura 2000 state 2025 for Mol.") + parser.add_argument("--base-url", default=os.environ.get("GEOINTEL_INTERNAL_API_URL", DEFAULT_API_URL)) + parser.add_argument("--project-name", default=DEFAULT_PROJECT_NAME) + parser.add_argument("--area-name", default=DEFAULT_AREA_NAME) + parser.add_argument( + "--output-dir", + type=Path, + default=Path(os.environ.get("GEOINTEL_BWK_OUTPUT_DIR", DEFAULT_OUTPUT_DIR)), + ) + parser.add_argument("--page-limit", type=int, default=1000) + parser.add_argument("--max-features", type=int, default=30_000) + parser.add_argument("--request-timeout", type=int, default=180) + parser.add_argument("--import-timeout", type=int, default=1800) + parser.add_argument("--force", action="store_true", help="Refetch source pages and rebuild the artifact.") + parser.add_argument("--fetch-only", action="store_true", help="Build evidence without importing a Dataset.") + return parser.parse_args() + + +def utc_now() -> str: + return datetime.now(timezone.utc).isoformat() + + +def sha256_bytes(value: bytes) -> str: + return hashlib.sha256(value).hexdigest() + + +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_bytes_atomic(path: Path, value: bytes) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + temporary = path.with_suffix(path.suffix + ".tmp") + temporary.write_bytes(value) + temporary.replace(path) + + +def write_json_atomic(path: Path, value: Any, *, pretty: bool = False) -> None: + encoded = json.dumps( + value, + ensure_ascii=False, + indent=2 if pretty else None, + separators=None if pretty else (",", ":"), + ).encode("utf-8") + write_bytes_atomic(path, encoded) + + +def source_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() + session.headers.update({"User-Agent": "GeoIntel-BWK-Natura2000-Mol-Operator/1.0"}) + adapter = HTTPAdapter(max_retries=retry) + session.mount("https://", adapter) + session.mount("http://", adapter) + return session + + +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[:300]}") from exc + if not response.ok: + raise RuntimeError( + f"GeoIntel API failed ({response.status_code}): {json.dumps(payload, ensure_ascii=False)[:800]}" + ) + 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 paginated_items(session: requests.Session, url: str, *, timeout: int) -> list[dict[str, Any]]: + items: list[dict[str, Any]] = [] + offset = 0 + total: int | None = None + while total is None or offset < total: + page = response_data(session.get(url, params={"limit": 200, "offset": offset}, timeout=timeout)) + page_items = list(page.get("items") or []) + page_total = int(page.get("total") or 0) + if total is None: + total = page_total + elif page_total != total: + raise RuntimeError("GeoIntel pagination total changed while reading the workspace") + items.extend(page_items) + if not page_items: + break + offset += len(page_items) + if total is not None and len(items) != total: + raise RuntimeError(f"GeoIntel pagination returned {len(items)} of {total} items") + return items + + +def polygonal_geometry(geometry): + if geometry is None or geometry.is_empty: + return None + if not geometry.is_valid: + geometry = make_valid(geometry) + polygons: list[Polygon] = [] + + def collect(candidate) -> None: + if candidate is None or candidate.is_empty: + return + if isinstance(candidate, Polygon): + polygons.append(candidate) + return + if isinstance(candidate, MultiPolygon): + polygons.extend(part for part in candidate.geoms if not part.is_empty) + return + if hasattr(candidate, "geoms"): + for part in candidate.geoms: + collect(part) + + collect(geometry) + if not polygons: + return None + result = unary_union(polygons) + if not result.is_valid: + result = make_valid(result) + return result if not result.is_empty and result.is_valid else None + + +def locate_workspace( + session: requests.Session, + base_url: str, + project_name: str, + area_name: str, + timeout: int, +) -> tuple[str, str, Any, list[dict[str, Any]]]: + projects = paginated_items(session, f"{base_url}/api/v1/projects", timeout=timeout) + project = next((item for item in projects if item.get("name") == project_name), None) + if not project: + raise RuntimeError(f"Project {project_name!r} is missing") + project_id = str(project["id"]) + areas = paginated_items(session, f"{base_url}/api/v1/projects/{project_id}/areas", timeout=timeout) + normalized_target = area_name.casefold().replace("officiele", "officiële") + area = next( + ( + item + for item in areas + if str(item.get("name") or "").casefold().replace("officiele", "officiële") == normalized_target + ), + None, + ) + if not area: + area = next((item for item in areas if "gemeente mol" in str(item.get("name") or "").casefold()), None) + if not area or not area.get("geometry"): + raise RuntimeError(f"Persisted Mol Area {area_name!r} with geometry is missing") + boundary_wgs84 = polygonal_geometry(shape(area["geometry"])) + if boundary_wgs84 is None: + raise RuntimeError("Persisted Mol Area geometry is invalid or not polygonal") + datasets = paginated_items(session, f"{base_url}/api/v1/projects/{project_id}/datasets", timeout=timeout) + return project_id, str(area["id"]), boundary_wgs84, datasets + + +def next_page_url(payload: dict[str, Any]) -> str | None: + for link in payload.get("links") or []: + if str(link.get("rel") or "").lower() == "next" and link.get("href"): + return str(link["href"]) + return None + + +def iter_wfs_pages( + session: requests.Session, + bbox: tuple[float, float, float, float], + *, + page_limit: int, + timeout: int, +) -> Iterable[tuple[dict[str, Any], str, bytes]]: + url: str | None = WFS_URL + base_params: dict[str, str] = { + "service": "WFS", + "version": "2.0.0", + "request": "GetFeature", + "typeNames": TYPE_NAME, + "srsName": "EPSG:4326", + "bbox": ",".join(f"{value:.8f}" for value in bbox) + ",EPSG:4326", + "count": str(page_limit), + "sortBy": "UIDN", + "outputFormat": "application/json", + } + params: dict[str, str] | None = dict(base_params) + seen: set[str] = set() + next_start_index = 0 + while url: + request_key = requests.Request("GET", url, params=params).prepare().url or url + if request_key in seen: + raise RuntimeError("BWK WFS pagination loop detected") + seen.add(request_key) + response = session.get(url, params=params, timeout=timeout) + params = None + response.raise_for_status() + payload = response.json() + if not isinstance(payload, dict) or payload.get("type") != "FeatureCollection": + raise RuntimeError("BWK WFS returned an invalid FeatureCollection") + yield payload, response.url, response.content + next_url = next_page_url(payload) + number_returned = int(payload.get("numberReturned") or len(payload.get("features") or [])) + next_start_index += number_returned + if next_url: + url = next_url + params = None + elif number_returned >= page_limit: + url = WFS_URL + params = {**base_params, "startIndex": str(next_start_index)} + else: + url = None + params = None + + +def habitat_breakdown(properties: dict[str, Any]) -> tuple[list[dict[str, Any]], float, float, float]: + entries: list[dict[str, Any]] = [] + natura_share = 0.0 + regional_share = 0.0 + uncertain_share = 0.0 + for index in range(1, 6): + code = str(properties.get(f"HAB{index}") or "").strip() + try: + share = max(0.0, min(100.0, float(properties.get(f"PHAB{index}") or 0.0))) + except (TypeError, ValueError): + share = 0.0 + if not code or share <= 0: + continue + normalized = code.lower() + entries.append({"code": code, "share_percent": share}) + if re.match(r"^\d", normalized): + natura_share += share + elif normalized.startswith("rbb"): + regional_share += share + elif normalized == "ohab": + uncertain_share += share + status = str(properties.get("HABLEGENDE") or "").strip().lower() + if status == "ohab" and uncertain_share <= 0: + uncertain_share = 100.0 + return entries, min(natura_share, 100.0), min(regional_share, 100.0), min(uncertain_share, 100.0) + + +def normalize_feature(feature: dict[str, Any], boundary_lambert72) -> tuple[dict[str, Any] | None, bool]: + geometry_payload = feature.get("geometry") + if not geometry_payload: + return None, False + source_wgs84 = polygonal_geometry(shape(geometry_payload)) + if source_wgs84 is None: + return None, False + source_lambert72 = polygonal_geometry(transform_geometry(TO_LAMBERT72.transform, source_wgs84)) + if source_lambert72 is None or not source_lambert72.intersects(boundary_lambert72): + return None, False + was_clipped = not source_lambert72.within(boundary_lambert72) + clipped_lambert72 = polygonal_geometry(source_lambert72.intersection(boundary_lambert72)) + if clipped_lambert72 is None or clipped_lambert72.area <= 0: + return None, was_clipped + clipped_wgs84 = polygonal_geometry(transform_geometry(TO_WGS84.transform, clipped_lambert72)) + if clipped_wgs84 is None: + return None, was_clipped + + raw_properties = dict(feature.get("properties") or {}) + source_id = str(feature.get("id") or raw_properties.get("UIDN") or raw_properties.get("OIDN") or "").strip() + if not source_id: + source_id = sha256_bytes(json.dumps(geometry_payload, sort_keys=True).encode("utf-8")) + stable_id = f"BWK:Bwkhab:{raw_properties.get('UIDN') or source_id}" + evaluation_code = str(raw_properties.get("EVAL") or "").strip().lower() + habitat_status = str(raw_properties.get("HABLEGENDE") or "").strip().lower() + habitats, natura_share, regional_share, uncertain_share = habitat_breakdown(raw_properties) + clipped_area_ha = float(clipped_lambert72.area) / 10_000.0 + natura_codes = [entry["code"] for entry in habitats if re.match(r"^\d", str(entry["code"]))] + regional_codes = [entry["code"] for entry in habitats if str(entry["code"]).lower().startswith("rbb")] + + properties = { + **raw_properties, + "source_name": "inbo_bwk_natura2000", + "source_collection": TYPE_NAME, + "source_feature_id": stable_id, + "reference_layer_name": "nature_value", + "theme": "nature_value", + "authority_level": "authoritative", + "coverage_scope": "municipality", + "municipality": "Mol", + "nis_code": "13025", + "source_version": SOURCE_VERSION, + "attribution": ATTRIBUTION, + "bwk_evaluation_code": evaluation_code or "unknown", + "bwk_evaluation_label": EVALUATION_LABELS.get(evaluation_code, "Onbekende of ontbrekende BWK-waardering"), + "bwk_label": str(raw_properties.get("BWKLABEL") or "").strip(), + "bwk_units": ", ".join( + str(raw_properties.get(f"EENH{index}") or "").strip() + for index in range(1, 9) + if str(raw_properties.get(f"EENH{index}") or "").strip() + ), + "bwk_origin_code": str(raw_properties.get("HERK") or "").strip(), + "habitat_status_code": habitat_status or "unknown", + "habitat_status_label": HABITAT_STATUS_LABELS.get(habitat_status, "Onbekende of ontbrekende habitatstatus"), + "natura2000_codes": ", ".join(natura_codes), + "regional_biotope_codes": ", ".join(regional_codes), + "habitat_entries": habitats, + "habitat_origin_code": str(raw_properties.get("HERKHAB") or "").strip(), + "habitat_share_origin_code": str(raw_properties.get("HERKPHAB") or "").strip(), + "clipped_area_ha": round(clipped_area_ha, 8), + "natura2000_share_percent": natura_share, + "regional_biotope_share_percent": regional_share, + "uncertain_habitat_share_percent": uncertain_share, + "natura2000_area_ha": round(clipped_area_ha * natura_share / 100.0, 8), + "regional_biotope_area_ha": round(clipped_area_ha * regional_share / 100.0, 8), + "uncertain_habitat_area_ha": round(clipped_area_ha * uncertain_share / 100.0, 8), + "habitat_share_spatial_limitation": ( + "PHAB percentages apply to the full source polygon; a partial map selection scales them by intersected area." + ), + } + return { + "type": "Feature", + "id": stable_id, + "geometry": mapping(clipped_wgs84), + "properties": properties, + }, was_clipped + + +def reusable_artifact(output_dir: Path) -> tuple[Path, Path, dict[str, Any]] | None: + artifact_path = output_dir / DATASET_FILENAME + manifest_path = output_dir / MANIFEST_FILENAME + if not artifact_path.is_file() or not manifest_path.is_file(): + return None + try: + manifest = json.loads(manifest_path.read_text(encoding="utf-8")) + except (OSError, ValueError): + return None + if manifest.get("schema_version") != SCHEMA_VERSION or manifest.get("source_version") != SOURCE_VERSION: + return None + if manifest.get("artifact_sha256") != sha256_file(artifact_path): + return None + for page in manifest.get("raw_pages") or []: + page_path = output_dir / str(page.get("path") or "") + if not page_path.is_file() or page.get("sha256") != sha256_file(page_path): + return None + return artifact_path, manifest_path, manifest + + +def prepare_artifact( + session: requests.Session, + boundary_wgs84, + output_dir: Path, + *, + page_limit: int, + max_features: int, + timeout: int, + force: bool, +) -> tuple[Path, Path, dict[str, Any]]: + if not force: + reusable = reusable_artifact(output_dir) + if reusable: + return reusable + + output_dir.mkdir(parents=True, exist_ok=True) + raw_dir = output_dir / "raw" + raw_dir.mkdir(parents=True, exist_ok=True) + boundary_lambert72 = polygonal_geometry(transform_geometry(TO_LAMBERT72.transform, boundary_wgs84)) + if boundary_lambert72 is None: + raise RuntimeError("Mol boundary could not be transformed to EPSG:31370") + + retained: list[dict[str, Any]] = [] + raw_pages: list[dict[str, Any]] = [] + source_urls: list[str] = [] + seen_ids: set[str] = set() + raw_feature_count = 0 + duplicate_count = 0 + rejected_count = 0 + clipped_count = 0 + + for page_number, (payload, source_url, raw_bytes) in enumerate( + iter_wfs_pages(session, boundary_wgs84.bounds, page_limit=page_limit, timeout=timeout), + start=1, + ): + page_path = raw_dir / f"bwk_bwkhab_page_{page_number:05d}.json" + write_bytes_atomic(page_path, raw_bytes) + page_features = list(payload.get("features") or []) + raw_feature_count += len(page_features) + if raw_feature_count > max_features: + raise RuntimeError( + f"BWK WFS returned more than the {max_features} feature safety limit; refusing a truncated import" + ) + raw_pages.append( + { + "path": str(page_path.relative_to(output_dir)), + "sha256": sha256_bytes(raw_bytes), + "size_bytes": len(raw_bytes), + "feature_count": len(page_features), + "source_url": source_url, + } + ) + source_urls.append(source_url) + for source_feature in page_features: + source_identity = str( + source_feature.get("id") + or (source_feature.get("properties") or {}).get("UIDN") + or "" + ) + if source_identity and source_identity in seen_ids: + duplicate_count += 1 + continue + if source_identity: + seen_ids.add(source_identity) + normalized, was_clipped = normalize_feature(source_feature, boundary_lambert72) + if normalized is None: + rejected_count += 1 + continue + if was_clipped: + clipped_count += 1 + retained.append(normalized) + + if not retained: + raise RuntimeError("BWK WFS returned no valid polygon features inside the persisted Mol Area") + + evaluation_area: dict[str, float] = defaultdict(float) + habitat_status_area: dict[str, float] = defaultdict(float) + natura_area = 0.0 + regional_area = 0.0 + uncertain_area = 0.0 + for feature in retained: + properties = feature["properties"] + area = float(properties["clipped_area_ha"]) + evaluation_area[str(properties["bwk_evaluation_code"])] += area + habitat_status_area[str(properties["habitat_status_code"])] += area + natura_area += float(properties["natura2000_area_ha"]) + regional_area += float(properties["regional_biotope_area_ha"]) + uncertain_area += float(properties["uncertain_habitat_area_ha"]) + + generated_at = utc_now() + artifact = { + "type": "FeatureCollection", + "name": "BWK en Natura 2000 - toestand 2025 - Gemeente Mol", + "features": retained, + "source": "INBO BWK/Natura 2000 WFS", + "source_version": SOURCE_VERSION, + "attribution": ATTRIBUTION, + "catalog_url": CATALOG_URL, + "report_url": REPORT_URL, + "coverage_scope": "municipality", + "municipality": "Mol", + "nis_code": "13025", + "reference_truncated": False, + "generated_at": generated_at, + } + artifact_path = output_dir / DATASET_FILENAME + write_json_atomic(artifact_path, artifact) + manifest = { + "schema_version": SCHEMA_VERSION, + "source_version": SOURCE_VERSION, + "source_type_name": TYPE_NAME, + "wfs_url": WFS_URL, + "catalog_url": CATALOG_URL, + "report_url": REPORT_URL, + "attribution": ATTRIBUTION, + "generated_at": generated_at, + "crs_source": "EPSG:4326", + "crs_clip": "EPSG:31370", + "crs_persisted": "EPSG:4326", + "boundary_sha256": sha256_bytes(json.dumps(mapping(boundary_wgs84), sort_keys=True).encode("utf-8")), + "boundary_bbox_wgs84": list(boundary_wgs84.bounds), + "page_limit": page_limit, + "page_count": len(raw_pages), + "raw_source_feature_count": raw_feature_count, + "feature_count": len(retained), + "duplicate_count": duplicate_count, + "rejected_or_outside_count": rejected_count, + "clipped_feature_count": clipped_count, + "reference_truncated": False, + "raw_pages": raw_pages, + "source_urls": source_urls, + "evaluation_area_ha": {key: round(value, 6) for key, value in sorted(evaluation_area.items())}, + "habitat_status_area_ha": {key: round(value, 6) for key, value in sorted(habitat_status_area.items())}, + "natura2000_area_ha": round(natura_area, 6), + "regional_biotope_area_ha": round(regional_area, 6), + "uncertain_habitat_area_ha": round(uncertain_area, 6), + "artifact_path": str(artifact_path), + "artifact_sha256": sha256_file(artifact_path), + "artifact_size_bytes": artifact_path.stat().st_size, + "limitations": [ + "The 2025 edition is the best available map state, not one uniform 2025 field survey.", + "PHAB percentages may be theoretical shares; partial selections scale shares by intersected polygon area.", + "The real field situation remains authoritative for policy and legal use.", + ], + } + manifest_path = output_dir / MANIFEST_FILENAME + write_json_atomic(manifest_path, manifest, pretty=True) + return artifact_path, manifest_path, manifest + + +def selection_metrics() -> list[dict[str, Any]]: + share_warning = ( + "Oppervlakte op basis van de officiele PHAB-aandelen. Automatisch verdeelde aandelen kunnen lokaal " + "afwijken van de terreinsituatie." + ) + return [ + { + "metric_key": "bwk_very_valuable_area", + "method": "intersection_area", + "label": "Biologisch zeer waardevol", + "unit": "ha", + "geometry_dimension": 2, + "filter_property": "bwk_evaluation_code", + "filter_values": ["z"], + }, + { + "metric_key": "bwk_valuable_area", + "method": "intersection_area", + "label": "Biologisch waardevol", + "unit": "ha", + "geometry_dimension": 2, + "filter_property": "bwk_evaluation_code", + "filter_values": ["w"], + }, + { + "metric_key": "bwk_less_valuable_area", + "method": "intersection_area", + "label": "Biologisch minder waardevol", + "unit": "ha", + "geometry_dimension": 2, + "filter_property": "bwk_evaluation_code", + "filter_values": ["m"], + }, + { + "metric_key": "bwk_mixed_value_area", + "method": "intersection_area", + "label": "Gemengde BWK-waardering", + "unit": "ha", + "geometry_dimension": 2, + "filter_property": "bwk_evaluation_code", + "filter_values": ["wz", "mwz", "mz", "mw"], + }, + { + "metric_key": "natura2000_area", + "method": "area_weighted_sum", + "property": "natura2000_area_ha", + "label": "Natura 2000-habitat", + "unit": "ha", + "is_estimate": True, + "warning": share_warning, + "warning_only_when_estimate": False, + }, + { + "metric_key": "regional_biotope_area", + "method": "area_weighted_sum", + "property": "regional_biotope_area_ha", + "label": "Regionaal belangrijk biotoop", + "unit": "ha", + "is_estimate": True, + "warning": share_warning, + "warning_only_when_estimate": False, + }, + { + "metric_key": "uncertain_habitat_area", + "method": "area_weighted_sum", + "property": "uncertain_habitat_area_ha", + "label": "Onzeker habitat / kennislacune", + "unit": "ha", + "is_estimate": True, + "warning": "Dit is een maximale potentiele oppervlakte van kaartvlakken met de status ohab, geen bevestigd habitat.", + "warning_only_when_estimate": False, + }, + ] + + +def upload_artifact( + session: requests.Session, + *, + base_url: str, + project_id: str, + area_id: str, + artifact_path: Path, + manifest_path: Path, + manifest: dict[str, Any], + timeout: int, +) -> dict[str, Any]: + source_metadata = { + "provider": "Instituut voor Natuur- en Bosonderzoek", + "theme": "nature_value", + "layer_name": "BWK/Natura 2000 toestand 2025", + "authority_level": "authoritative", + "coverage_scope": "municipality", + "municipality": "Mol", + "nis_code": "13025", + "feature_count": manifest["feature_count"], + "geometry_clipped_to_area": True, + "semantic_metrics": False, + "attribution": ATTRIBUTION, + "catalog_url": CATALOG_URL, + "report_url": REPORT_URL, + "selection_aggregation": { + "metric_key": "bwk_mapped_area", + "method": "intersection_area", + "label": "BWK-gekarteerde oppervlakte", + "unit": "ha", + "geometry_dimension": 2, + "warning": ( + "Uitgave 2025 is de best beschikbare kaarttoestand, maar niet elk kaartvlak is in 2025 op terrein gekarteerd." + ), + }, + "selection_metrics": selection_metrics(), + } + provenance_metadata = { + "operator_tool": "provision_mol_bwk_natura2000.py", + "operator_explicit_fetch": True, + "geometry_clipped_to_area": True, + "source_type_name": TYPE_NAME, + "wfs_url": WFS_URL, + "catalog_url": CATALOG_URL, + "report_url": REPORT_URL, + "manifest_path": str(manifest_path), + "artifact_sha256": manifest["artifact_sha256"], + "raw_page_checksums": {page["path"]: page["sha256"] for page in manifest["raw_pages"]}, + "source_urls": manifest["source_urls"], + "reference_truncated": False, + "generated_at": manifest["generated_at"], + "limitations": manifest["limitations"], + } + with artifact_path.open("rb") as handle: + response = session.post( + f"{base_url}/api/v1/projects/{project_id}/datasets/upload", + data={ + "dataset_type": "vector", + "source": "operator_official_import", + "dataset_role": "reference", + "source_name": "inbo_bwk_natura2000", + "reference_layer_name": "nature_value", + "source_metadata_json": json.dumps(source_metadata, ensure_ascii=False), + "provenance_metadata_json": json.dumps(provenance_metadata, ensure_ascii=False), + "area_id": area_id, + "observed_at": PUBLICATION_DATE, + "temporal_granularity": "snapshot", + "source_version": SOURCE_VERSION, + }, + files={"file": (artifact_path.name, handle, "application/geo+json")}, + timeout=timeout, + ) + return response_data(response) + + +def main() -> int: + args = parse_args() + if args.page_limit < 1 or args.page_limit > 5000 or args.max_features < args.page_limit: + print(json.dumps({"status": "error", "message": "Invalid page or feature safety limits"}), file=sys.stderr) + return 2 + base_url = args.base_url.rstrip("/") + try: + with requests.Session() as api_session: + project_id, area_id, boundary, existing = locate_workspace( + api_session, + base_url, + args.project_name, + args.area_name, + args.import_timeout, + ) + with source_session() as official_session: + artifact_path, manifest_path, manifest = prepare_artifact( + official_session, + boundary, + args.output_dir, + page_limit=args.page_limit, + max_features=args.max_features, + timeout=args.request_timeout, + force=args.force, + ) + existing_dataset = next( + ( + item + for item in existing + if item.get("source_name") == "inbo_bwk_natura2000" + and item.get("source_version") == SOURCE_VERSION + and str(item.get("area_id") or "") == area_id + ), + None, + ) + if existing_dataset: + persisted_checksum = str(existing_dataset.get("checksum_sha256") or "") + if persisted_checksum and persisted_checksum != manifest["artifact_sha256"]: + raise RuntimeError( + "A different BWK 2025 Mol artifact is already persisted; refusing a silent duplicate or replacement" + ) + result = { + "status": "existing", + "dataset_id": existing_dataset["id"], + "feature_count": existing_dataset.get("feature_count") or manifest["feature_count"], + } + elif args.fetch_only: + result = { + "status": "prepared", + "artifact_path": str(artifact_path), + "feature_count": manifest["feature_count"], + } + else: + dataset = upload_artifact( + api_session, + base_url=base_url, + project_id=project_id, + area_id=area_id, + artifact_path=artifact_path, + manifest_path=manifest_path, + manifest=manifest, + timeout=args.import_timeout, + ) + result = { + "status": "imported", + "dataset_id": dataset["id"], + "feature_count": dataset.get("feature_count"), + } + except (KeyError, OSError, RuntimeError, ValueError, requests.RequestException) as exc: + print(json.dumps({"status": "error", "message": str(exc)}, ensure_ascii=False), file=sys.stderr) + return 1 + + print( + json.dumps( + { + "status": "ok", + "project": args.project_name, + "area": args.area_name, + "source_version": SOURCE_VERSION, + "manifest_path": str(manifest_path), + "metrics": { + "evaluation_area_ha": manifest["evaluation_area_ha"], + "natura2000_area_ha": manifest["natura2000_area_ha"], + "regional_biotope_area_ha": manifest["regional_biotope_area_ha"], + "uncertain_habitat_area_ha": manifest["uncertain_habitat_area_ha"], + }, + "result": result, + }, + ensure_ascii=False, + indent=2, + ) + ) + return 0 + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/scripts/run_readiness_check.sh b/scripts/run_readiness_check.sh index 2ae93f38..5fac8864 100755 --- a/scripts/run_readiness_check.sh +++ b/scripts/run_readiness_check.sh @@ -48,6 +48,7 @@ ${PYTHON_BIN} -m py_compile scripts/provision_mol_population_history.py ${PYTHON_BIN} -m py_compile scripts/provision_mol_historical_landuse.py ${PYTHON_BIN} -m py_compile scripts/provision_official_landuse_timeseries.py ${PYTHON_BIN} -m py_compile scripts/provision_waterinfo_station_history.py +${PYTHON_BIN} -m py_compile scripts/provision_mol_bwk_natura2000.py ${PYTHON_BIN} -m py_compile scripts/provision_regional_timeseries.py ${PYTHON_BIN} -m py_compile scripts/geographic_scopes.py ${PYTHON_BIN} -m py_compile scripts/provision_geographic_scope.py