diff --git a/CHANGELOG.md b/CHANGELOG.md
index d5aee0ae..8b95831e 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -9,6 +9,22 @@
## Unreleased - Post-V1 capability completion (2026-07-19)
+- Added the official WALOUS 2018 GeoTIFF as a third live-provisioned Walloon
+ land-cover epoch. Its stable SPW artifact, archive/raster checksums, EPSG:3812
+ identity and published stacked class codes are validated fail-closed. The
+ official view-class crosswalk normalizes stacked codes to the existing
+ 11-class series, while source value `0` is explicitly treated as background
+ nodata and never contributes to area metrics.
+- Added bounded Walloon terrain acquisition from the official SPW MNT
+ 2021-2022 1 m GeoTIFF. Operator provisioning validates archive bounds, safe
+ extraction, CRS, resolution, band count, elevation samples and checksums;
+ selection persistence and terrain metrics retain DNG/EPSG:5710 instead of
+ incorrectly labelling Walloon elevations as TAW.
+- Extended detection-model capabilities with machine-readable training scope,
+ validation scope, validated regions, national-validation status and the
+ operator-review requirement. The configured local model remains bound to
+ Mol/Kempen evidence and cannot become nationally labelled through frontend
+ copy alone.
- Completed live WALOUS 2020/2023 provisioning and fixed signed `int8` source
reads with nodata `-128` across acquisition, analysis and PNG rendering. A
dedicated regression now proves conversion to the persisted `uint8`/`255`
diff --git a/backend/README.md b/backend/README.md
index 6946e8e2..1d3e6660 100644
--- a/backend/README.md
+++ b/backend/README.md
@@ -1603,12 +1603,12 @@ stored as `bounded_selection`.
The Wallonia map flow uses bounded PICC vector products, the queryable legal
SPW flood-hazard polygon layer and provisioned official WALOUS land-cover
-rasters. Provision the 2020 and 2023 source editions once in the persistent
+rasters. Provision the 2018, 2020 and 2023 source editions once in the persistent
storage mount:
```bash
docker exec geointel python /app/scripts/provision_walous_sources.py \
- --years 2020 2023 \
+ --years 2018 2020 2023 \
--destination /app/storage/source-cache/walous
```
@@ -1616,12 +1616,14 @@ The provisioner verifies advertised archive sizes, safe ZIP structure,
EPSG:3812, one band, 1 m cells, the official non-contiguous class codes
`1,2,3,4,5,6,7,8,9,80,90` and SHA-256 checksums. It does not run at
application startup. `GET .../datasets/walous/products` therefore reports
-`source_not_provisioned` until both source files exist.
+`source_not_provisioned` for each edition whose source file is absent.
For a bounded Walloon selection the browser persists the latest edition and
all other configured comparable editions. `POST .../raster/walous/select`
returns cell-area hectares; the temporal API compares the same semantic metric
-keys for 2020 and 2023. WALOUS is land cover, not legal land use, ownership,
+keys for 2018, 2020 and 2023. The 2018 stacked classes use the official visible-
+class crosswalk and retain the earlier-method limitation. WALOUS is land cover,
+not legal land use, ownership,
tree count, timber volume or water volume.
The class semantics follow the official raster codes, not display-list
@@ -1636,10 +1638,13 @@ Settings: `WALOUS_ENABLED`, `WALOUS_SOURCE_DIR`,
`WALOUS_MAX_PIXELS`. The SPW flood polygon adapter uses
`SPW_FLOOD_HAZARD_ENABLED` and `SPW_FLOOD_HAZARD_MAPSERVER_URL`.
-The official Walloon 2021-2022 DTM is currently not an implicit runtime asset:
-the published 1 m whole-region artifact is about 41 GB and the 0.5 m INSPIRE
-artifact about 213 GB. A later operator capacity plan must define storage,
-partitioning and refresh before it can be called operational.
+The official Walloon 2021-2022 1 m MNT is an explicit operator asset. Provision
+it once with `scripts/provision_spw_terrain_source.py`; the runtime then reads
+only bounded windows and persists 5 m analysis derivatives. The full 0.5 m
+artifact remains intentionally excluded because it adds no V1 metric and is
+about 213 GB. Settings: `SPW_TERRAIN_ENABLED`, `SPW_TERRAIN_SOURCE_DIR`,
+`SPW_TERRAIN_ANALYSIS_RESOLUTION_M`, `SPW_TERRAIN_MAX_SIDE_M` and
+`SPW_TERRAIN_MAX_PIXELS`.
Provision the official DOV soil polygons for Mol through the existing vector
upload path:
diff --git a/backend/app/api/routes/datasets.py b/backend/app/api/routes/datasets.py
index b3ea7691..5df2f053 100644
--- a/backend/app/api/routes/datasets.py
+++ b/backend/app/api/routes/datasets.py
@@ -18,6 +18,8 @@ from app.schemas import (
BathymetrySourceRead,
DatasetList,
DhmvProductRead,
+ SpwTerrainAcquireRequest,
+ SpwTerrainProductRead,
Envelope,
FloodHazardProductRead,
FloodHazardSelectionResponse,
@@ -88,6 +90,7 @@ from app.services.grb_acquisition_service import GrbAcquisitionService
from app.services.official_vector_acquisition_service import OfficialVectorAcquisitionService
from app.services.orthophoto_acquisition_service import OrthophotoAcquisitionService
from app.services.dhmv_acquisition_service import DhmvAcquisitionService
+from app.services.spw_terrain_service import SpwTerrainService
from app.services.terrain_analysis_service import TerrainAnalysisService
from app.services.flood_hazard_acquisition_service import FloodHazardAcquisitionService
from app.services.flood_hazard_analysis_service import FloodHazardAnalysisService
@@ -243,6 +246,33 @@ def list_dhmv_products(project_id: UUID, db: Session = Depends(get_db)):
return envelope({"items": items, "total": len(items)})
+@router.post("/datasets/spw-terrain/acquire", response_model=Envelope[JobRead])
+def acquire_bounded_spw_terrain(
+ project_id: UUID,
+ payload: SpwTerrainAcquireRequest,
+ db: Session = Depends(get_db),
+):
+ job = JobService.run_sync_job(
+ db=db,
+ project_id=project_id,
+ job_type="raster.spw-terrain.acquire",
+ parameters=payload.model_dump(mode="json"),
+ operation=lambda: SpwTerrainService.acquire(db, project_id, payload),
+ )
+ return envelope(job)
+
+
+@router.get(
+ "/datasets/spw-terrain/products",
+ response_model=Envelope[ItemList[SpwTerrainProductRead]],
+)
+def list_spw_terrain_products(project_id: UUID, db: Session = Depends(get_db)):
+ if not db.get(Project, project_id):
+ raise AppError(code="PROJECT_NOT_FOUND", message="Project not found", status_code=404)
+ items = SpwTerrainService.list_products()
+ return envelope({"items": items, "total": len(items)})
+
+
@router.post("/datasets/grb/acquire", response_model=Envelope[JobRead])
def acquire_bounded_grb(
project_id: UUID,
diff --git a/backend/app/core/config.py b/backend/app/core/config.py
index ea10c908..6be53d2c 100644
--- a/backend/app/core/config.py
+++ b/backend/app/core/config.py
@@ -294,6 +294,19 @@ class Settings(BaseSettings):
)
walous_max_side_m: float = Field(default=60_000.0, gt=0, validation_alias="WALOUS_MAX_SIDE_M")
walous_max_pixels: int = Field(default=36_000_000, ge=1, validation_alias="WALOUS_MAX_PIXELS")
+ spw_terrain_enabled: bool = Field(default=True, validation_alias="SPW_TERRAIN_ENABLED")
+ spw_terrain_source_dir: str = Field(
+ default="/app/storage/source-cache/spw-terrain",
+ validation_alias="SPW_TERRAIN_SOURCE_DIR",
+ )
+ spw_terrain_analysis_resolution_m: float = Field(
+ default=5.0,
+ ge=1.0,
+ le=10.0,
+ validation_alias="SPW_TERRAIN_ANALYSIS_RESOLUTION_M",
+ )
+ spw_terrain_max_side_m: float = Field(default=20_000.0, gt=0, validation_alias="SPW_TERRAIN_MAX_SIDE_M")
+ spw_terrain_max_pixels: int = Field(default=12_000_000, ge=1, validation_alias="SPW_TERRAIN_MAX_PIXELS")
redis_url: str | None = Field(default=None, validation_alias="REDIS_URL")
log_level: str = Field(default="INFO", validation_alias="GEOINTEL_LOG_LEVEL")
sql_log_level: str = Field(default="WARNING", validation_alias="GEOINTEL_SQL_LOG_LEVEL")
diff --git a/backend/app/schemas/__init__.py b/backend/app/schemas/__init__.py
index f98e568a..60c0870f 100644
--- a/backend/app/schemas/__init__.py
+++ b/backend/app/schemas/__init__.py
@@ -78,6 +78,11 @@ from .dhmv import (
TerrainSelectionResponse,
TerrainSelectionSummary,
)
+from .spw_terrain import (
+ SpwTerrainAcquireRequest,
+ SpwTerrainAcquisitionResult,
+ SpwTerrainProductRead,
+)
from .flood_hazard import (
FloodHazardAcquireRequest,
FloodHazardAcquisitionResult,
@@ -248,6 +253,9 @@ __all__ = [
"DhmvAcquireRequest",
"DhmvAcquisitionResult",
"DhmvProductRead",
+ "SpwTerrainAcquireRequest",
+ "SpwTerrainAcquisitionResult",
+ "SpwTerrainProductRead",
"TerrainMetric",
"TerrainPartitionSelectionRequest",
"TerrainSelectionRequest",
diff --git a/backend/app/schemas/detection.py b/backend/app/schemas/detection.py
index 59ccbe58..feaaad02 100644
--- a/backend/app/schemas/detection.py
+++ b/backend/app/schemas/detection.py
@@ -18,6 +18,11 @@ class DetectionModelCapability(BaseModel):
status: str
limitation_message: str
version: str | None = None
+ training_scope: str | None = None
+ validation_scope: str | None = None
+ validated_regions: list[str] = Field(default_factory=list)
+ nationally_validated: bool = False
+ operator_review_required: bool = True
class DetectionModelsResponse(BaseModel):
diff --git a/backend/app/schemas/spw_terrain.py b/backend/app/schemas/spw_terrain.py
new file mode 100644
index 00000000..b1c4737e
--- /dev/null
+++ b/backend/app/schemas/spw_terrain.py
@@ -0,0 +1,55 @@
+from __future__ import annotations
+
+from uuid import UUID
+
+from pydantic import BaseModel, Field
+
+from .operations import VectorSelectionBBox
+
+
+class SpwTerrainAcquireRequest(BaseModel):
+ bbox: VectorSelectionBBox
+ area_id: UUID | None = None
+ product_key: str = "spw_mnt_1m_2021_2022"
+ resolution_m: float | None = Field(default=None, ge=1.0, le=10.0)
+ force_refresh: bool = False
+
+
+class SpwTerrainProductRead(BaseModel):
+ key: str
+ display_name: str
+ surface_model: str
+ source_filename: str
+ native_resolution_m: float
+ analysis_resolution_m: float
+ source_crs: str
+ vertical_reference: str
+ acquisition_period: str
+ catalog_url: str
+ attribution: str
+ license_note: str
+ limitation_message: str
+ coverage_zones: list[str]
+ configured: bool
+ status: str
+
+
+class SpwTerrainAcquisitionResult(BaseModel):
+ output_dataset_id: UUID
+ reused: bool
+ provider: str
+ product_key: str
+ display_name: str
+ surface_model: str
+ native_resolution_m: float
+ resolution_m: float
+ width: int
+ height: int
+ valid_pixel_count: int
+ nodata_value: float
+ bbox_epsg4326: list[float]
+ bbox_epsg3812: list[float]
+ vertical_reference: str
+ acquisition_period: str
+ attribution: str
+ limitation_message: str
diff --git a/backend/app/services/coverage_registry_service.py b/backend/app/services/coverage_registry_service.py
index abb292ef..84545ca3 100644
--- a/backend/app/services/coverage_registry_service.py
+++ b/backend/app/services/coverage_registry_service.py
@@ -253,10 +253,10 @@ SOURCE_DEFINITIONS = (
license_note="Consult the license of each Geoportail Wallonie product.",
limitation_message=(
"Bounded PICC buildings, road axes and hydrography, the legally current flood-hazard polygons, "
- "and operator-imported SPW bathymetry are operational; other Walloon themes remain separately governed."
+ "operator-imported SPW bathymetry and bounded SPW MNT terrain are operational; other Walloon themes remain separately governed."
),
- materialized_source_names=("spw_picc", "spw_flood_hazard", "spw_walous_land_cover", "spw_bathymetry"),
- operational_themes=("buildings", "roads", "surface_water", "land_cover_use", "flood_climate", "bathymetry"),
+ materialized_source_names=("spw_picc", "spw_flood_hazard", "spw_walous_land_cover", "spw_bathymetry", "spw_terrain"),
+ operational_themes=("buildings", "roads", "surface_water", "land_cover_use", "elevation", "flood_climate", "bathymetry"),
),
_contract(
source_name="urbis",
@@ -376,6 +376,7 @@ REGIONAL_THEME_DATASETS: dict[str, dict[str, dict[str, tuple[str, ...]]]] = {
"roads": {"spw_picc": ("roads",)},
"surface_water": {"spw_picc": ("water",)},
"land_cover_use": {"spw_walous_land_cover": ()},
+ "elevation": {"spw_terrain": ()},
"flood_climate": {"spw_flood_hazard": ("flood_hazard",)},
"bathymetry": {"spw_bathymetry": ()},
},
diff --git a/backend/app/services/model_registry_service.py b/backend/app/services/model_registry_service.py
index 3a980cae..d4b7071b 100644
--- a/backend/app/services/model_registry_service.py
+++ b/backend/app/services/model_registry_service.py
@@ -5,7 +5,10 @@ from typing import Type
from app.core.config import Settings, get_settings
from app.schemas.detection import DetectionModelCapability
-from app.services.segmentation_adapter import SamSegmentationAdapter, YoloSegmentationAdapter
+from app.services.segmentation_adapter import (
+ SamSegmentationAdapter,
+ YoloSegmentationAdapter,
+)
from app.services.yolo_adapter import YoloDetectionAdapter
@@ -39,7 +42,9 @@ class ModelRegistryService:
limitation_message="YOLO/PyTorch inference is not configured in Sprint 8; no model is downloaded or executed.",
version=None,
),
- ModelRegistryService._configured_yolo_capability(resolved_settings, yolo_adapter_class),
+ ModelRegistryService._configured_yolo_capability(
+ resolved_settings, yolo_adapter_class
+ ),
DetectionModelCapability(
model_id="manual-fixture-detector",
display_name="Manual fixture detector",
@@ -104,8 +109,12 @@ class ModelRegistryService:
limitation_message="Fixture segmenter is for explicit tests/demo fixtures only and is not production inference.",
version="fixture-v1",
),
- ModelRegistryService._configured_yolo_seg_capability(resolved_settings, yolo_seg_adapter_class),
- ModelRegistryService._configured_sam_capability(resolved_settings, sam_adapter_class),
+ ModelRegistryService._configured_yolo_seg_capability(
+ resolved_settings, yolo_seg_adapter_class
+ ),
+ ModelRegistryService._configured_sam_capability(
+ resolved_settings, sam_adapter_class
+ ),
]
@staticmethod
@@ -119,7 +128,11 @@ class ModelRegistryService:
"YOLO segmentation is disabled. Set YOLO_SEG_ENABLED=true and YOLO_SEG_MODEL_PATH to a local "
"segmentation model file to enable inference. GeoIntel never downloads model weights automatically."
)
- model_path = Path(settings.yolo_seg_model_path).expanduser() if settings.yolo_seg_model_path else None
+ model_path = (
+ Path(settings.yolo_seg_model_path).expanduser()
+ if settings.yolo_seg_model_path
+ else None
+ )
if settings.yolo_seg_enabled:
if not adapter_class.dependencies_available():
@@ -157,7 +170,11 @@ class ModelRegistryService:
"SAM is disabled. Set SAM_ENABLED=true and SAM_MODEL_PATH to a local SAM-compatible model file to "
"enable class-agnostic segmentation. GeoIntel never downloads model weights automatically."
)
- model_path = Path(settings.sam_model_path).expanduser() if settings.sam_model_path else None
+ model_path = (
+ Path(settings.sam_model_path).expanduser()
+ if settings.sam_model_path
+ else None
+ )
if settings.sam_enabled:
if not adapter_class.dependencies_available():
@@ -192,7 +209,11 @@ class ModelRegistryService:
configured = False
status = "not_configured"
limitation = "YOLO is disabled. Set YOLO_ENABLED=true and YOLO_MODEL_PATH to a local model file to enable inference."
- model_path = Path(settings.yolo_model_path).expanduser() if settings.yolo_model_path else None
+ model_path = (
+ Path(settings.yolo_model_path).expanduser()
+ if settings.yolo_model_path
+ else None
+ )
if settings.yolo_enabled:
if not yolo_adapter_class.dependencies_available():
@@ -217,4 +238,11 @@ class ModelRegistryService:
status=status,
limitation_message=limitation,
version=settings.yolo_model_version,
+ training_scope=(
+ "Operator-managed local weights; the runtime has no nationally governed training-corpus evidence."
+ ),
+ validation_scope="Mol and the Kempen operator evidence; no Belgian national validation matrix is bound.",
+ validated_regions=["flanders_mol_kempen"],
+ nationally_validated=False,
+ operator_review_required=True,
)
diff --git a/backend/app/services/spw_terrain_service.py b/backend/app/services/spw_terrain_service.py
new file mode 100644
index 00000000..bab664e2
--- /dev/null
+++ b/backend/app/services/spw_terrain_service.py
@@ -0,0 +1,467 @@
+from __future__ import annotations
+
+from datetime import UTC, datetime
+import hashlib
+import json
+import math
+from pathlib import Path
+from typing import Any
+from uuid import UUID
+
+from geoalchemy2.shape import to_shape
+from pyproj import Transformer
+from shapely.geometry import box, mapping
+from shapely.ops import transform as shapely_transform
+
+from app.core.config import Settings, get_settings
+from app.core.errors import AppError
+from app.models import Area, Dataset, Project
+from app.schemas.spw_terrain import (
+ SpwTerrainAcquireRequest,
+ SpwTerrainAcquisitionResult,
+ SpwTerrainProductRead,
+)
+from app.services.dataset_service import DatasetService
+
+
+class SpwTerrainService:
+ PROVIDER = "spw_terrain"
+ PRODUCT_KEY = "spw_mnt_1m_2021_2022"
+ DISPLAY_NAME = "SPW terreinmodel (MNT) 2021-2022"
+ SOURCE_FILENAME = "spw_mnt_1m_2021_2022_3812.tif"
+ SOURCE_SHA256_FILENAME = "spw_mnt_1m_2021_2022_3812.sha256"
+ SOURCE_CRS = "EPSG:3812"
+ SOURCE_RESOLUTION_M = 1.0
+ SURFACE_MODEL = "terrain"
+ VERTICAL_REFERENCE = "DNG / Deuxieme Nivellement General (EPSG:5710)"
+ VERTICAL_UNIT_LABEL = "m DNG"
+ ACQUISITION_PERIOD = "2021-02-19/2022-03-05"
+ CATALOG_URL = "https://geoportail.wallonie.be/catalogue/fe13bc84-e371-46ca-9632-8ad4139f1ee5.html"
+ DOWNLOAD_URL = (
+ "https://geoservices.wallonie.be/geotraitement/spwdatadownload/results/"
+ "fe13bc84-e371-46ca-9632-8ad4139f1ee5/RELIEF_WALLONIE_MNT_1M_2021_2022_GEOTIFF_3812.zip"
+ )
+ ATTRIBUTION = (
+ "Service public de Wallonie (SPW) - Relief de la Wallonie MNT 2021-2022"
+ )
+ LICENSE_NOTE = "CC BY 4.0; cite SPW and identify modifications."
+ NODATA = -9999.0
+ LIMITATION = (
+ "GeoIntel leest uitsluitend een begrensd venster uit het checksum-gevalideerde officiele 1 m MNT en "
+ "bewaart een analyse-afgeleide op de gekozen resolutie. Het MNT beschrijft maaiveldhoogte in DNG, niet "
+ "oppervlaktehoogte, afstroming, waterdiepte of watervolume. Kleine bronzones zijn door SPW geinterpoleerd."
+ )
+
+ @staticmethod
+ def _source_path(settings: Settings) -> Path:
+ return Path(settings.spw_terrain_source_dir) / SpwTerrainService.SOURCE_FILENAME
+
+ @staticmethod
+ def list_products(*, settings: Settings | None = None) -> list[dict[str, Any]]:
+ resolved = settings or get_settings()
+ configured = (
+ resolved.spw_terrain_enabled
+ and SpwTerrainService._source_path(resolved).is_file()
+ )
+ product = SpwTerrainProductRead(
+ key=SpwTerrainService.PRODUCT_KEY,
+ display_name=SpwTerrainService.DISPLAY_NAME,
+ surface_model=SpwTerrainService.SURFACE_MODEL,
+ source_filename=SpwTerrainService.SOURCE_FILENAME,
+ native_resolution_m=SpwTerrainService.SOURCE_RESOLUTION_M,
+ analysis_resolution_m=resolved.spw_terrain_analysis_resolution_m,
+ source_crs=SpwTerrainService.SOURCE_CRS,
+ vertical_reference=SpwTerrainService.VERTICAL_REFERENCE,
+ acquisition_period=SpwTerrainService.ACQUISITION_PERIOD,
+ catalog_url=SpwTerrainService.CATALOG_URL,
+ attribution=SpwTerrainService.ATTRIBUTION,
+ license_note=SpwTerrainService.LICENSE_NOTE,
+ limitation_message=SpwTerrainService.LIMITATION,
+ coverage_zones=["wallonia"],
+ configured=configured,
+ status="configured" if configured else "source_not_provisioned",
+ )
+ return [product.model_dump()]
+
+ @staticmethod
+ def _scope_geometry(db, project_id: UUID, payload: SpwTerrainAcquireRequest):
+ if not db.get(Project, project_id):
+ raise AppError(
+ code="PROJECT_NOT_FOUND", message="Project not found", status_code=404
+ )
+ if payload.product_key.strip().lower() != SpwTerrainService.PRODUCT_KEY:
+ raise AppError(
+ code="SPW_TERRAIN_PRODUCT_NOT_SUPPORTED",
+ message="Select the governed SPW MNT 2021-2022 product",
+ details={"product_key": payload.product_key},
+ status_code=422,
+ )
+ if payload.bbox.crs.upper() != "EPSG:4326":
+ raise AppError(
+ code="INVALID_BBOX_CRS",
+ message="SPW terrain acquisition requires EPSG:4326",
+ status_code=400,
+ )
+ values = [
+ payload.bbox.min_x,
+ payload.bbox.min_y,
+ payload.bbox.max_x,
+ payload.bbox.max_y,
+ ]
+ if (
+ not all(math.isfinite(value) for value in values)
+ or values[0] >= values[2]
+ or values[1] >= values[3]
+ ):
+ raise AppError(
+ code="INVALID_BBOX",
+ message="SPW terrain selection must be a finite non-empty rectangle",
+ status_code=400,
+ )
+ selection = box(*values)
+ if payload.area_id is None:
+ return selection, values
+ area = db.get(Area, payload.area_id)
+ if area is None or area.project_id != project_id:
+ raise AppError(
+ code="AREA_NOT_FOUND", message="Area not found", status_code=404
+ )
+ selection = selection.intersection(to_shape(area.geometry))
+ if selection.is_empty or selection.area <= 0:
+ raise AppError(
+ code="SPW_TERRAIN_SELECTION_OUTSIDE_AREA",
+ message="Selection does not overlap the selected work area",
+ status_code=422,
+ )
+ return selection, values
+
+ @staticmethod
+ def _read_source_window(
+ source_path: Path, scope_4326, resolution: float, settings: Settings
+ ) -> tuple[bytes, dict[str, Any]]:
+ try:
+ import numpy as np
+ import rasterio
+ from rasterio.enums import Resampling
+ from rasterio.features import geometry_mask
+ from rasterio.io import MemoryFile
+ from rasterio.transform import from_bounds
+ from rasterio.windows import from_bounds as window_from_bounds
+ except ImportError as exc:
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio and numpy are required for SPW terrain",
+ status_code=503,
+ ) from exc
+
+ scope_metric = shapely_transform(
+ Transformer.from_crs(
+ "EPSG:4326", SpwTerrainService.SOURCE_CRS, always_xy=True
+ ).transform,
+ scope_4326,
+ )
+ try:
+ with rasterio.open(source_path) as source:
+ if (
+ source.crs is None
+ or source.crs.to_epsg() != 3812
+ or source.count != 1
+ ):
+ raise AppError(
+ code="SPW_TERRAIN_SOURCE_INVALID",
+ message="SPW MNT must be a one-band EPSG:3812 raster",
+ status_code=409,
+ )
+ if not all(
+ math.isclose(abs(float(value)), 1.0, abs_tol=0.05)
+ for value in source.res
+ ):
+ raise AppError(
+ code="SPW_TERRAIN_SOURCE_INVALID",
+ message="SPW MNT must retain the official 1 m resolution",
+ status_code=409,
+ )
+ clipped_geometry = scope_metric.intersection(box(*source.bounds))
+ if clipped_geometry.is_empty or clipped_geometry.area <= 0:
+ raise AppError(
+ code="SPW_TERRAIN_SELECTION_OUTSIDE_COVERAGE",
+ message="Selection does not overlap SPW MNT coverage",
+ status_code=422,
+ )
+ min_x, min_y, max_x, max_y = clipped_geometry.bounds
+ bounds = (
+ math.floor(min_x / resolution) * resolution,
+ math.floor(min_y / resolution) * resolution,
+ math.ceil(max_x / resolution) * resolution,
+ math.ceil(max_y / resolution) * resolution,
+ )
+ width_m, height_m = bounds[2] - bounds[0], bounds[3] - bounds[1]
+ if (
+ width_m > settings.spw_terrain_max_side_m
+ or height_m > settings.spw_terrain_max_side_m
+ ):
+ raise AppError(
+ code="SPW_TERRAIN_SELECTION_TOO_LARGE",
+ message="SPW terrain selection exceeds the configured side limit",
+ status_code=422,
+ )
+ width, height = (
+ max(1, round(width_m / resolution)),
+ max(1, round(height_m / resolution)),
+ )
+ if width * height > settings.spw_terrain_max_pixels:
+ raise AppError(
+ code="SPW_TERRAIN_SELECTION_TOO_LARGE",
+ message="SPW terrain selection exceeds the configured cell limit",
+ details={
+ "pixel_count": width * height,
+ "max_pixels": settings.spw_terrain_max_pixels,
+ },
+ status_code=422,
+ )
+ window = window_from_bounds(*bounds, transform=source.transform)
+ band = source.read(
+ 1,
+ window=window,
+ out_shape=(height, width),
+ masked=True,
+ resampling=Resampling.bilinear,
+ )
+ output_transform = from_bounds(*bounds, width, height)
+ outside_scope = geometry_mask(
+ [mapping(clipped_geometry)],
+ out_shape=(height, width),
+ transform=output_transform,
+ invert=False,
+ )
+ values = np.asarray(np.ma.getdata(band), dtype="float32")
+ invalid = (
+ np.ma.getmaskarray(band) | outside_scope | ~np.isfinite(values)
+ )
+ if source.nodata is not None:
+ invalid |= np.isclose(
+ values.astype("float64"), float(source.nodata)
+ )
+ values[invalid] = SpwTerrainService.NODATA
+ valid = values[~invalid]
+ if valid.size == 0:
+ raise AppError(
+ code="SPW_TERRAIN_NO_VALID_DATA",
+ message="SPW MNT contains no valid cells in this selection",
+ status_code=422,
+ )
+ if float(valid.min()) < -100.0 or float(valid.max()) > 1000.0:
+ raise AppError(
+ code="SPW_TERRAIN_SOURCE_INVALID_VALUES",
+ message="SPW MNT contains implausible elevations for Wallonia",
+ details={
+ "minimum": float(valid.min()),
+ "maximum": float(valid.max()),
+ },
+ status_code=409,
+ )
+ profile = {
+ "driver": "GTiff",
+ "width": width,
+ "height": height,
+ "count": 1,
+ "dtype": "float32",
+ "crs": SpwTerrainService.SOURCE_CRS,
+ "transform": output_transform,
+ "nodata": SpwTerrainService.NODATA,
+ "compress": "deflate",
+ "predictor": 3,
+ }
+ with MemoryFile() as memory:
+ with memory.open(**profile) as output:
+ output.write(values, 1)
+ content = memory.read()
+ return content, {
+ "width": width,
+ "height": height,
+ "valid_pixel_count": int(valid.size),
+ "bbox_epsg3812": list(bounds),
+ "source_width": int(source.width),
+ "source_height": int(source.height),
+ "source_nodata": None
+ if source.nodata is None
+ else float(source.nodata),
+ "source_resolution_m": 1.0,
+ "analysis_resolution_m": resolution,
+ "elevation_min_m": float(valid.min()),
+ "elevation_max_m": float(valid.max()),
+ }
+ except AppError:
+ raise
+ except Exception as exc:
+ raise AppError(
+ code="SPW_TERRAIN_SOURCE_READ_FAILED",
+ message="The provisioned SPW MNT could not be read",
+ details={"reason": str(exc)},
+ status_code=500,
+ ) from exc
+
+ @staticmethod
+ def _cached_dataset(db, project_id: UUID, filename: str) -> Dataset | None:
+ candidate = (
+ db.query(Dataset)
+ .filter(
+ Dataset.project_id == project_id,
+ Dataset.name == filename,
+ Dataset.source_name == SpwTerrainService.PROVIDER,
+ Dataset.status == "ready",
+ )
+ .order_by(Dataset.imported_at.desc())
+ .first()
+ )
+ return (
+ candidate
+ if candidate
+ and candidate.storage_path
+ and Path(candidate.storage_path).is_file()
+ else None
+ )
+
+ @staticmethod
+ def acquire(
+ db,
+ project_id: UUID,
+ payload: SpwTerrainAcquireRequest,
+ *,
+ settings: Settings | None = None,
+ ) -> dict[str, Any]:
+ resolved = settings or get_settings()
+ if not resolved.spw_terrain_enabled:
+ raise AppError(
+ code="SPW_TERRAIN_NOT_CONFIGURED",
+ message="SPW terrain bounded analysis is disabled",
+ status_code=503,
+ )
+ source_path = SpwTerrainService._source_path(resolved)
+ if not source_path.is_file():
+ raise AppError(
+ code="SPW_TERRAIN_SOURCE_NOT_PROVISIONED",
+ message="The official SPW MNT source archive has not been provisioned on this runtime",
+ details={
+ "expected_path": str(source_path),
+ "operator_command": "python scripts/provision_spw_terrain_source.py",
+ },
+ status_code=503,
+ )
+ scope, bbox_4326 = SpwTerrainService._scope_geometry(db, project_id, payload)
+ resolution = float(
+ payload.resolution_m or resolved.spw_terrain_analysis_resolution_m
+ )
+ identity = {
+ "product_key": SpwTerrainService.PRODUCT_KEY,
+ "bbox_epsg4326": [round(float(value), 8) for value in bbox_4326],
+ "area_id": str(payload.area_id) if payload.area_id else None,
+ "analysis_resolution_m": resolution,
+ }
+ request_hash = hashlib.sha256(
+ json.dumps(identity, sort_keys=True).encode()
+ ).hexdigest()
+ filename = f"spw_mnt_2021_2022_{request_hash[:12]}_3812.tif"
+ if not payload.force_refresh:
+ cached = SpwTerrainService._cached_dataset(db, project_id, filename)
+ if cached is not None:
+ metadata = cached.source_metadata or {}
+ return SpwTerrainAcquisitionResult(
+ output_dataset_id=cached.id,
+ reused=True,
+ provider=SpwTerrainService.PROVIDER,
+ product_key=SpwTerrainService.PRODUCT_KEY,
+ display_name=SpwTerrainService.DISPLAY_NAME,
+ surface_model=SpwTerrainService.SURFACE_MODEL,
+ native_resolution_m=SpwTerrainService.SOURCE_RESOLUTION_M,
+ resolution_m=resolution,
+ width=int((cached.metadata_json or {}).get("width", 0)),
+ height=int((cached.metadata_json or {}).get("height", 0)),
+ valid_pixel_count=int(metadata.get("valid_pixel_count", 0)),
+ nodata_value=SpwTerrainService.NODATA,
+ bbox_epsg4326=bbox_4326,
+ bbox_epsg3812=list(metadata.get("bbox_epsg3812") or []),
+ vertical_reference=SpwTerrainService.VERTICAL_REFERENCE,
+ acquisition_period=SpwTerrainService.ACQUISITION_PERIOD,
+ attribution=SpwTerrainService.ATTRIBUTION,
+ limitation_message=SpwTerrainService.LIMITATION,
+ ).model_dump(mode="json")
+
+ content, validation = SpwTerrainService._read_source_window(
+ source_path, scope, resolution, resolved
+ )
+ source_sha256_path = source_path.with_name(
+ SpwTerrainService.SOURCE_SHA256_FILENAME
+ )
+ source_sha256 = (
+ source_sha256_path.read_text(encoding="ascii").strip().split()[0]
+ if source_sha256_path.is_file()
+ else None
+ )
+ acquired_at = datetime.now(UTC)
+ dataset = DatasetService.import_raster_bytes(
+ db,
+ project_id=project_id,
+ area_id=payload.area_id,
+ filename=filename,
+ content=content,
+ source="SPW Relief de la Wallonie MNT 2021-2022 operator-provisioned GeoTIFF",
+ source_name=SpwTerrainService.PROVIDER,
+ observed_at=datetime(2022, 3, 5, 23, 59, 59, tzinfo=UTC),
+ valid_from=datetime(2021, 2, 19, tzinfo=UTC),
+ valid_to=datetime(2022, 3, 5, 23, 59, 59, tzinfo=UTC),
+ temporal_granularity="acquisition_period",
+ source_version="RELIEF_WALLONIE_MNT_1M_2021_2022",
+ source_metadata={
+ "provider": SpwTerrainService.PROVIDER,
+ "product_key": SpwTerrainService.PRODUCT_KEY,
+ "product_display_name": SpwTerrainService.DISPLAY_NAME,
+ "surface_model": SpwTerrainService.SURFACE_MODEL,
+ "source_crs": SpwTerrainService.SOURCE_CRS,
+ "source_resolution_m": SpwTerrainService.SOURCE_RESOLUTION_M,
+ "analysis_resolution_m": validation["analysis_resolution_m"],
+ "valid_pixel_count": validation["valid_pixel_count"],
+ "bbox_epsg4326": bbox_4326,
+ "bbox_epsg3812": validation["bbox_epsg3812"],
+ "coverage_zones": ["wallonia"],
+ "vertical_reference": SpwTerrainService.VERTICAL_REFERENCE,
+ "vertical_unit": "m",
+ "vertical_unit_label": SpwTerrainService.VERTICAL_UNIT_LABEL,
+ "acquisition_period": SpwTerrainService.ACQUISITION_PERIOD,
+ "catalog_url": SpwTerrainService.CATALOG_URL,
+ "download_url": SpwTerrainService.DOWNLOAD_URL,
+ "attribution": SpwTerrainService.ATTRIBUTION,
+ "license_note": SpwTerrainService.LICENSE_NOTE,
+ "limitation_message": SpwTerrainService.LIMITATION,
+ },
+ provenance_metadata={
+ "acquisition": "operator_provisioned_official_archive_bounded_window",
+ "acquired_at": acquired_at.isoformat(),
+ "request_hash": request_hash,
+ "source_filename": SpwTerrainService.SOURCE_FILENAME,
+ "source_sha256": source_sha256,
+ "derived_sha256": hashlib.sha256(content).hexdigest(),
+ "resampling": "bilinear",
+ "validation": validation,
+ },
+ )
+ return SpwTerrainAcquisitionResult(
+ output_dataset_id=dataset.id,
+ reused=False,
+ provider=SpwTerrainService.PROVIDER,
+ product_key=SpwTerrainService.PRODUCT_KEY,
+ display_name=SpwTerrainService.DISPLAY_NAME,
+ surface_model=SpwTerrainService.SURFACE_MODEL,
+ native_resolution_m=SpwTerrainService.SOURCE_RESOLUTION_M,
+ resolution_m=validation["analysis_resolution_m"],
+ width=validation["width"],
+ height=validation["height"],
+ valid_pixel_count=validation["valid_pixel_count"],
+ nodata_value=SpwTerrainService.NODATA,
+ bbox_epsg4326=bbox_4326,
+ bbox_epsg3812=validation["bbox_epsg3812"],
+ vertical_reference=SpwTerrainService.VERTICAL_REFERENCE,
+ acquisition_period=SpwTerrainService.ACQUISITION_PERIOD,
+ attribution=SpwTerrainService.ATTRIBUTION,
+ limitation_message=SpwTerrainService.LIMITATION,
+ ).model_dump(mode="json")
diff --git a/backend/app/services/terrain_analysis_service.py b/backend/app/services/terrain_analysis_service.py
index dc5471c1..5cb12f19 100644
--- a/backend/app/services/terrain_analysis_service.py
+++ b/backend/app/services/terrain_analysis_service.py
@@ -22,10 +22,13 @@ from app.schemas.dhmv import (
TerrainSelectionSummary,
)
from app.services.dhmv_acquisition_service import DhmvAcquisitionService
-from app.services.raster_partition_analysis_service import RasterPartitionAnalysisService
+from app.services.raster_partition_analysis_service import (
+ RasterPartitionAnalysisService,
+)
class TerrainAnalysisService:
+ SUPPORTED_PROVIDERS = {DhmvAcquisitionService.PROVIDER, "spw_terrain"}
UNSUPPORTED_METRICS = ["water_depth_m", "water_volume_m3"]
LIMITATION = (
"Hoogte, reliëf en helling zijn afgeleid uit DHMV II. Afstroming vraagt bijkomende hydrologische modellering. "
@@ -36,30 +39,58 @@ class TerrainAnalysisService:
def _load_dataset(db, project_id: UUID, dataset_id: UUID) -> Dataset:
dataset = db.get(Dataset, dataset_id)
if not dataset or dataset.project_id != project_id:
- raise AppError(code="DATASET_NOT_FOUND", message="Dataset not found", status_code=404)
- if dataset.dataset_type != "raster" or dataset.source_name != DhmvAcquisitionService.PROVIDER:
+ raise AppError(
+ code="DATASET_NOT_FOUND", message="Dataset not found", status_code=404
+ )
+ if (
+ dataset.dataset_type != "raster"
+ or dataset.source_name not in TerrainAnalysisService.SUPPORTED_PROVIDERS
+ ):
raise AppError(
code="INVALID_TERRAIN_DATASET",
- message="Terrain analysis requires a governed DHMV raster dataset",
+ message="Terrain analysis requires a governed regional elevation raster",
status_code=400,
)
- if dataset.status != "ready" or not dataset.storage_path or not Path(dataset.storage_path).is_file():
- raise AppError(code="DATASET_FILE_MISSING", message="Persisted DHMV raster file is unavailable", status_code=404)
+ if (
+ dataset.status != "ready"
+ or not dataset.storage_path
+ or not Path(dataset.storage_path).is_file()
+ ):
+ raise AppError(
+ code="DATASET_FILE_MISSING",
+ message="Persisted terrain raster file is unavailable",
+ status_code=404,
+ )
return dataset
@staticmethod
def _selection_geometry(db, project_id: UUID, payload: TerrainSelectionRequest):
- selection = box(payload.bbox.min_x, payload.bbox.min_y, payload.bbox.max_x, payload.bbox.max_y)
+ selection = box(
+ payload.bbox.min_x,
+ payload.bbox.min_y,
+ payload.bbox.max_x,
+ payload.bbox.max_y,
+ )
if payload.area_id is None:
return selection
area = db.get(Area, payload.area_id)
if not area:
- raise AppError(code="AREA_NOT_FOUND", message="Area not found", status_code=404)
+ raise AppError(
+ code="AREA_NOT_FOUND", message="Area not found", status_code=404
+ )
if area.project_id != project_id:
- raise AppError(code="INVALID_DATASET_SCOPE", message="Area does not belong to this project", status_code=400)
+ raise AppError(
+ code="INVALID_DATASET_SCOPE",
+ message="Area does not belong to this project",
+ status_code=400,
+ )
selection = selection.intersection(to_shape(area.geometry))
if selection.is_empty or selection.area <= 0:
- raise AppError(code="TERRAIN_SELECTION_OUTSIDE_AREA", message="Selection does not overlap the selected work area", status_code=422)
+ raise AppError(
+ code="TERRAIN_SELECTION_OUTSIDE_AREA",
+ message="Selection does not overlap the selected work area",
+ status_code=422,
+ )
return selection
@staticmethod
@@ -73,27 +104,51 @@ class TerrainAnalysisService:
) -> dict:
resolved_settings = settings or get_settings()
dataset = TerrainAnalysisService._load_dataset(db, project_id, dataset_id)
- selection_4326 = TerrainAnalysisService._selection_geometry(db, project_id, payload)
+ selection_4326 = TerrainAnalysisService._selection_geometry(
+ db, project_id, payload
+ )
try:
import numpy as np
import rasterio
from rasterio.features import geometry_mask
from rasterio.mask import mask
except ImportError as exc:
- raise AppError(code="RASTER_PROCESSING_UNAVAILABLE", message="Rasterio and numpy are required for terrain analysis", status_code=503) from exc
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio and numpy are required for terrain analysis",
+ status_code=503,
+ ) from exc
source_metadata = dataset.source_metadata or {}
product_key = str(source_metadata.get("product_key") or "")
surface_model = str(source_metadata.get("surface_model") or "")
- if product_key not in DhmvAcquisitionService._products() or surface_model not in {"terrain", "surface"}:
- raise AppError(code="INVALID_TERRAIN_METADATA", message="DHMV product provenance is incomplete", status_code=409)
+ product_is_governed = (
+ product_key in DhmvAcquisitionService._products()
+ if dataset.source_name == DhmvAcquisitionService.PROVIDER
+ else product_key == "spw_mnt_1m_2021_2022"
+ )
+ if not product_is_governed or surface_model not in {"terrain", "surface"}:
+ raise AppError(
+ code="INVALID_TERRAIN_METADATA",
+ message="Regional terrain product provenance is incomplete",
+ status_code=409,
+ )
+ vertical_unit_label = str(source_metadata.get("vertical_unit_label") or "m TAW")
try:
with rasterio.open(dataset.storage_path) as source:
if source.crs is None:
- raise AppError(code="INVALID_DATASET_CRS", message="DHMV raster CRS is missing", status_code=409)
- transformer = Transformer.from_crs("EPSG:4326", source.crs, always_xy=True)
- selection_metric = shapely_transform(transformer.transform, selection_4326)
+ raise AppError(
+ code="INVALID_DATASET_CRS",
+ message="Terrain raster CRS is missing",
+ status_code=409,
+ )
+ transformer = Transformer.from_crs(
+ "EPSG:4326", source.crs, always_xy=True
+ )
+ selection_metric = shapely_transform(
+ transformer.transform, selection_4326
+ )
source_extent = box(*source.bounds)
analysis_geometry = selection_metric.intersection(source_extent)
if analysis_geometry.is_empty or analysis_geometry.area <= 0:
@@ -103,12 +158,17 @@ class TerrainAnalysisService:
status_code=422,
)
min_x, min_y, max_x, max_y = analysis_geometry.bounds
- expected_cells = math.ceil((max_x - min_x) / abs(source.res[0])) * math.ceil((max_y - min_y) / abs(source.res[1]))
+ expected_cells = math.ceil(
+ (max_x - min_x) / abs(source.res[0])
+ ) * math.ceil((max_y - min_y) / abs(source.res[1]))
if expected_cells > resolved_settings.dhmv_max_pixels:
raise AppError(
code="TERRAIN_SELECTION_TOO_LARGE",
message="Terrain analysis exceeds the configured raster cell limit",
- details={"pixel_count": expected_cells, "max_pixels": resolved_settings.dhmv_max_pixels},
+ details={
+ "pixel_count": expected_cells,
+ "max_pixels": resolved_settings.dhmv_max_pixels,
+ },
status_code=422,
)
clipped, clipped_transform = mask(
@@ -133,14 +193,20 @@ class TerrainAnalysisService:
valid_mask = selected_cells & ~np.ma.getmaskarray(elevation) & ~invalid
values = raw[valid_mask]
if values.size == 0:
- raise AppError(code="TERRAIN_NO_VALID_DATA", message="No valid DHMV height cells occur in this selection", status_code=422)
+ raise AppError(
+ code="TERRAIN_NO_VALID_DATA",
+ message="No valid terrain height cells occur in this selection",
+ status_code=422,
+ )
resolution_x = abs(float(source.res[0]))
resolution_y = abs(float(source.res[1]))
slope_values = np.asarray([], dtype="float64")
if raw.shape[0] >= 2 and raw.shape[1] >= 2:
surface = np.where(valid_mask, raw, np.nan)
- gradient_y, gradient_x = np.gradient(surface, resolution_y, resolution_x)
+ gradient_y, gradient_x = np.gradient(
+ surface, resolution_y, resolution_x
+ )
slope = np.degrees(np.arctan(np.hypot(gradient_x, gradient_y)))
slope_values = slope[np.isfinite(slope) & valid_mask]
except AppError:
@@ -148,12 +214,14 @@ class TerrainAnalysisService:
except Exception as exc:
raise AppError(
code="TERRAIN_ANALYSIS_FAILED",
- message="The persisted DHMV raster could not be analysed",
+ message="The persisted terrain raster could not be analysed",
details={"reason": str(exc)},
status_code=500,
) from exc
- def metric(key: str, label: str, value: float, unit: str, method: str) -> TerrainMetric:
+ def metric(
+ key: str, label: str, value: float, unit: str, method: str
+ ) -> TerrainMetric:
return TerrainMetric(
metric_key=key,
metric_label=label,
@@ -163,21 +231,79 @@ class TerrainAnalysisService:
)
prefix = "terrain" if surface_model == "terrain" else "surface"
- elevation_label = "Gemiddelde maaiveldhoogte" if surface_model == "terrain" else "Gemiddelde oppervlaktehoogte"
+ elevation_label = (
+ "Gemiddelde maaiveldhoogte"
+ if surface_model == "terrain"
+ else "Gemiddelde oppervlaktehoogte"
+ )
metrics = [
- metric(f"{prefix}_elevation_mean_m", elevation_label, values.mean(), "m TAW", "mean_valid_cells"),
- metric(f"{prefix}_elevation_min_m", "Laagste hoogte", values.min(), "m TAW", "minimum_valid_cells"),
- metric(f"{prefix}_elevation_max_m", "Hoogste hoogte", values.max(), "m TAW", "maximum_valid_cells"),
- metric(f"{prefix}_elevation_p10_m", "10e percentiel hoogte", np.percentile(values, 10), "m TAW", "percentile_10_valid_cells"),
- metric(f"{prefix}_elevation_p90_m", "90e percentiel hoogte", np.percentile(values, 90), "m TAW", "percentile_90_valid_cells"),
- metric("relief_m", "Reliëfverschil", values.max() - values.min(), "m", "maximum_minus_minimum"),
+ metric(
+ f"{prefix}_elevation_mean_m",
+ elevation_label,
+ values.mean(),
+ vertical_unit_label,
+ "mean_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_min_m",
+ "Laagste hoogte",
+ values.min(),
+ vertical_unit_label,
+ "minimum_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_max_m",
+ "Hoogste hoogte",
+ values.max(),
+ vertical_unit_label,
+ "maximum_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_p10_m",
+ "10e percentiel hoogte",
+ np.percentile(values, 10),
+ vertical_unit_label,
+ "percentile_10_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_p90_m",
+ "90e percentiel hoogte",
+ np.percentile(values, 90),
+ vertical_unit_label,
+ "percentile_90_valid_cells",
+ ),
+ metric(
+ "relief_m",
+ "Reliëfverschil",
+ values.max() - values.min(),
+ "m",
+ "maximum_minus_minimum",
+ ),
]
if slope_values.size:
metrics.extend(
[
- metric("slope_mean_deg", "Gemiddelde helling", slope_values.mean(), "°", "mean_finite_gradient"),
- metric("slope_p90_deg", "90e percentiel helling", np.percentile(slope_values, 90), "°", "percentile_90_finite_gradient"),
- metric("slope_max_deg", "Steilste helling", slope_values.max(), "°", "maximum_finite_gradient"),
+ metric(
+ "slope_mean_deg",
+ "Gemiddelde helling",
+ slope_values.mean(),
+ "°",
+ "mean_finite_gradient",
+ ),
+ metric(
+ "slope_p90_deg",
+ "90e percentiel helling",
+ np.percentile(slope_values, 90),
+ "°",
+ "percentile_90_finite_gradient",
+ ),
+ metric(
+ "slope_max_deg",
+ "Steilste helling",
+ slope_values.max(),
+ "°",
+ "maximum_finite_gradient",
+ ),
]
)
primary = metrics[0]
@@ -194,7 +320,10 @@ class TerrainAnalysisService:
slope_sample_count=int(slope_values.size),
coverage_ratio=round(float(values.size / max(1, selected_cell_count)), 6),
resolution_m=round(max(resolution_x, resolution_y), 4),
- vertical_reference=str(source_metadata.get("vertical_reference") or DhmvAcquisitionService.VERTICAL_REFERENCE),
+ vertical_reference=str(
+ source_metadata.get("vertical_reference")
+ or DhmvAcquisitionService.VERTICAL_REFERENCE
+ ),
summary=TerrainSelectionSummary(
metric_label=primary.metric_label,
metric_value=primary.metric_value,
@@ -204,7 +333,10 @@ class TerrainAnalysisService:
metrics=metrics,
),
unsupported_metrics=TerrainAnalysisService.UNSUPPORTED_METRICS,
- limitation_message=TerrainAnalysisService.LIMITATION,
+ limitation_message=str(
+ source_metadata.get("limitation_message")
+ or TerrainAnalysisService.LIMITATION
+ ),
generated_at=datetime.now(UTC).isoformat(),
)
return response.model_dump(mode="json")
@@ -218,7 +350,9 @@ class TerrainAnalysisService:
settings: Settings | None = None,
) -> dict:
resolved_settings = settings or get_settings()
- product = DhmvAcquisitionService._products().get(payload.product_key.strip().lower())
+ product = DhmvAcquisitionService._products().get(
+ payload.product_key.strip().lower()
+ )
if product is None:
raise AppError(
code="DHMV_PRODUCT_NOT_SUPPORTED",
@@ -226,7 +360,9 @@ class TerrainAnalysisService:
details={"product_key": payload.product_key},
status_code=422,
)
- selection_4326 = TerrainAnalysisService._selection_geometry(db, project_id, payload)
+ selection_4326 = TerrainAnalysisService._selection_geometry(
+ db, project_id, payload
+ )
partition = RasterPartitionAnalysisService.select(
db,
project_id,
@@ -279,7 +415,9 @@ class TerrainAnalysisService:
slope = np.degrees(np.arctan(np.hypot(gradient_x, gradient_y)))
slope_values = slope[np.isfinite(slope) & valid_mask]
- def metric(key: str, label: str, value: float, unit: str, method: str) -> TerrainMetric:
+ def metric(
+ key: str, label: str, value: float, unit: str, method: str
+ ) -> TerrainMetric:
return TerrainMetric(
metric_key=key,
metric_label=label,
@@ -295,19 +433,73 @@ class TerrainAnalysisService:
else "Gemiddelde oppervlaktehoogte"
)
metrics = [
- metric(f"{prefix}_elevation_mean_m", elevation_label, values.mean(), "m TAW", "mean_valid_cells"),
- metric(f"{prefix}_elevation_min_m", "Laagste hoogte", values.min(), "m TAW", "minimum_valid_cells"),
- metric(f"{prefix}_elevation_max_m", "Hoogste hoogte", values.max(), "m TAW", "maximum_valid_cells"),
- metric(f"{prefix}_elevation_p10_m", "10e percentiel hoogte", np.percentile(values, 10), "m TAW", "percentile_10_valid_cells"),
- metric(f"{prefix}_elevation_p90_m", "90e percentiel hoogte", np.percentile(values, 90), "m TAW", "percentile_90_valid_cells"),
- metric("relief_m", "Reliëfverschil", values.max() - values.min(), "m", "maximum_minus_minimum"),
+ metric(
+ f"{prefix}_elevation_mean_m",
+ elevation_label,
+ values.mean(),
+ "m TAW",
+ "mean_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_min_m",
+ "Laagste hoogte",
+ values.min(),
+ "m TAW",
+ "minimum_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_max_m",
+ "Hoogste hoogte",
+ values.max(),
+ "m TAW",
+ "maximum_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_p10_m",
+ "10e percentiel hoogte",
+ np.percentile(values, 10),
+ "m TAW",
+ "percentile_10_valid_cells",
+ ),
+ metric(
+ f"{prefix}_elevation_p90_m",
+ "90e percentiel hoogte",
+ np.percentile(values, 90),
+ "m TAW",
+ "percentile_90_valid_cells",
+ ),
+ metric(
+ "relief_m",
+ "Reliëfverschil",
+ values.max() - values.min(),
+ "m",
+ "maximum_minus_minimum",
+ ),
]
if slope_values.size:
metrics.extend(
[
- metric("slope_mean_deg", "Gemiddelde helling", slope_values.mean(), "°", "mean_finite_gradient"),
- metric("slope_p90_deg", "90e percentiel helling", np.percentile(slope_values, 90), "°", "percentile_90_finite_gradient"),
- metric("slope_max_deg", "Steilste helling", slope_values.max(), "°", "maximum_finite_gradient"),
+ metric(
+ "slope_mean_deg",
+ "Gemiddelde helling",
+ slope_values.mean(),
+ "°",
+ "mean_finite_gradient",
+ ),
+ metric(
+ "slope_p90_deg",
+ "90e percentiel helling",
+ np.percentile(slope_values, 90),
+ "°",
+ "percentile_90_finite_gradient",
+ ),
+ metric(
+ "slope_max_deg",
+ "Steilste helling",
+ slope_values.max(),
+ "°",
+ "maximum_finite_gradient",
+ ),
]
)
primary = metrics[0]
@@ -344,7 +536,9 @@ class TerrainAnalysisService:
return response.model_dump(mode="json")
@staticmethod
- def render_png(db, project_id: UUID, dataset_id: UUID, *, max_dimension: int = 1800) -> bytes:
+ def render_png(
+ db, project_id: UUID, dataset_id: UUID, *, max_dimension: int = 1800
+ ) -> bytes:
dataset = TerrainAnalysisService._load_dataset(db, project_id, dataset_id)
try:
import numpy as np
@@ -352,18 +546,31 @@ class TerrainAnalysisService:
from PIL import Image
from rasterio.enums import Resampling
except ImportError as exc:
- raise AppError(code="RASTER_PROCESSING_UNAVAILABLE", message="Rasterio, numpy and Pillow are required for terrain rendering", status_code=503) from exc
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio, numpy and Pillow are required for terrain rendering",
+ status_code=503,
+ ) from exc
try:
with rasterio.open(dataset.storage_path) as source:
scale = min(1.0, max_dimension / max(source.width, source.height))
width = max(1, round(source.width * scale))
height = max(1, round(source.height * scale))
- data = source.read(1, out_shape=(height, width), masked=True, resampling=Resampling.bilinear)
+ data = source.read(
+ 1,
+ out_shape=(height, width),
+ masked=True,
+ resampling=Resampling.bilinear,
+ )
values = np.asarray(data.filled(np.nan), dtype="float64")
valid = np.isfinite(values) & ~np.ma.getmaskarray(data)
if not valid.any():
- raise AppError(code="TERRAIN_NO_VALID_DATA", message="DHMV raster contains no renderable cells", status_code=422)
+ raise AppError(
+ code="TERRAIN_NO_VALID_DATA",
+ message="Terrain raster contains no renderable cells",
+ status_code=422,
+ )
low, high = np.percentile(values[valid], [2, 98])
if high <= low:
high = low + 1.0
@@ -381,7 +588,9 @@ class TerrainAnalysisService:
)
rgba = np.zeros((height, width, 4), dtype="uint8")
for channel in range(3):
- rgba[:, :, channel] = np.interp(normalized, stops, colors[:, channel]).astype("uint8")
+ rgba[:, :, channel] = np.interp(
+ normalized, stops, colors[:, channel]
+ ).astype("uint8")
rgba[:, :, 3] = np.where(valid, 225, 0).astype("uint8")
output = io.BytesIO()
Image.fromarray(rgba).save(output, format="PNG", optimize=True)
@@ -391,7 +600,7 @@ class TerrainAnalysisService:
except Exception as exc:
raise AppError(
code="TERRAIN_PREVIEW_FAILED",
- message="The persisted DHMV raster could not be rendered",
+ message="The persisted terrain raster could not be rendered",
details={"reason": str(exc)},
status_code=500,
) from exc
diff --git a/backend/app/services/walous_land_cover_service.py b/backend/app/services/walous_land_cover_service.py
index 824207b7..214a2242 100644
--- a/backend/app/services/walous_land_cover_service.py
+++ b/backend/app/services/walous_land_cover_service.py
@@ -40,7 +40,10 @@ class WalousProduct:
catalog_url: str
download_url: str
source_sha256_filename: str
+ attribution: str
accuracy_label: str
+ raw_class_crosswalk: dict[int, int] | None
+ comparability_note: str
observation_start: datetime
observation_end: datetime
@@ -54,7 +57,9 @@ class WalousLandCoverService:
METRIC_KIND = "categorical_area"
NODATA = 255
ATTRIBUTION = "Service public de Wallonie (SPW), Aerospacelab S.A."
- LICENSE_NOTE = "CC BY 4.0; cite the official SPW WALOUS edition and identify modifications."
+ LICENSE_NOTE = (
+ "CC BY 4.0; cite the official SPW WALOUS edition and identify modifications."
+ )
LIMITATION = (
"GeoIntel analyseert een nearest-neighbour afgeleide van het officiele 1 m WALOUS-raster op de "
"geconfigureerde analyseresolutie. Oppervlakten zijn celgebaseerde schattingen; de kaart is landbedekking, "
@@ -88,10 +93,74 @@ class WalousLandCoverService:
80: (78, 125, 70),
90: (107, 164, 87),
}
+ # The original 2018 product retains stacked two-digit codes. The official
+ # "Classe vue" legend resolves those codes to the visible top class. The
+ # only 2018-only visible class, greenhouses (62), is explicitly normalized
+ # to artificial constructions so the stable 11-class series can be used.
+ WALOUS_2018_CLASS_CROSSWALK = {
+ 0: NODATA,
+ 1: 1,
+ 11: 1,
+ 15: 1,
+ 18: 1,
+ 19: 1,
+ 31: 1,
+ 51: 1,
+ 71: 1,
+ 81: 1,
+ 91: 1,
+ 2: 2,
+ 28: 2,
+ 29: 2,
+ 62: 2,
+ 3: 3,
+ 38: 3,
+ 39: 3,
+ 73: 3,
+ 83: 3,
+ 93: 3,
+ 4: 4,
+ 5: 5,
+ 55: 5,
+ 58: 5,
+ 59: 5,
+ 75: 5,
+ 85: 5,
+ 95: 5,
+ 6: 6,
+ 7: 7,
+ 8: 8,
+ 9: 9,
+ 80: 80,
+ 90: 90,
+ }
@staticmethod
def _products() -> dict[str, WalousProduct]:
products = (
+ WalousProduct(
+ key="walous_land_cover_2018",
+ display_name="WALOUS landbedekking 2018",
+ observation_year=2018,
+ source_filename="walous_land_cover_2018_3812.tif",
+ source_version="WALOUS_OCS__2018",
+ catalog_url="https://geoportail.wallonie.be/catalogue/a0ad23a1-1845-4bd5-8c2f-0f62d3f1ec75.html",
+ download_url=(
+ "https://geoservices.wallonie.be/geotraitement/spwdatadownload/results/"
+ "a0ad23a1-1845-4bd5-8c2f-0f62d3f1ec75/WALOUS_OCS__2018_GEOTIFF_3812.zip"
+ ),
+ source_sha256_filename="walous_land_cover_2018_3812.sha256",
+ attribution="Service public de Wallonie (SPW), UCLouvain, ULB, ISSeP",
+ accuracy_label="Officiele globale nauwkeurigheid 91,5%",
+ raw_class_crosswalk=WalousLandCoverService.WALOUS_2018_CLASS_CROSSWALK,
+ comparability_note=(
+ "De 2018-editie gebruikt een eerdere, deels handmatig geconsolideerde methode. GeoIntel past de "
+ "officiele 'Classe vue'-crosswalk toe en groepeert de 2018-only serreklasse bij constructies; "
+ "trends blijven methodologisch begrensde schattingen."
+ ),
+ observation_start=datetime(2018, 1, 1, tzinfo=UTC),
+ observation_end=datetime(2018, 12, 31, 23, 59, 59, tzinfo=UTC),
+ ),
WalousProduct(
key="walous_land_cover_2020",
display_name="WALOUS landbedekking 2020",
@@ -104,7 +173,10 @@ class WalousLandCoverService:
"47b348f1-6e7a-4baa-963c-0232a43c0cff/WAL_OCS_IA__2020_GEOTIFF_3812.zip"
),
source_sha256_filename="walous_land_cover_2020_3812.sha256",
+ attribution=WalousLandCoverService.ATTRIBUTION,
accuracy_label="Officiele globale nauwkeurigheid 83,30%",
+ raw_class_crosswalk=None,
+ comparability_note="",
observation_start=datetime(2020, 4, 1, tzinfo=UTC),
observation_end=datetime(2020, 4, 24, 23, 59, 59, tzinfo=UTC),
),
@@ -120,7 +192,10 @@ class WalousLandCoverService:
"4e780ba1-463c-478e-95df-d2f1963a150d/WAL_OCS_IA__2023_GEOTIFF_3812.zip"
),
source_sha256_filename="walous_land_cover_2023_3812.sha256",
+ attribution=WalousLandCoverService.ATTRIBUTION,
accuracy_label="Officiele globale nauwkeurigheid 87,10%",
+ raw_class_crosswalk=None,
+ comparability_note="",
observation_start=datetime(2023, 5, 27, tzinfo=UTC),
observation_end=datetime(2023, 6, 25, 23, 59, 59, tzinfo=UTC),
),
@@ -136,7 +211,10 @@ class WalousLandCoverService:
resolved = settings or get_settings()
result: list[dict[str, Any]] = []
for product in WalousLandCoverService._products().values():
- configured = resolved.walous_enabled and WalousLandCoverService._source_path(resolved, product).is_file()
+ configured = (
+ resolved.walous_enabled
+ and WalousLandCoverService._source_path(resolved, product).is_file()
+ )
result.append(
ThematicRasterProductRead(
key=product.key,
@@ -151,12 +229,20 @@ class WalousLandCoverService:
observation_year=product.observation_year,
source_version=product.source_version,
catalog_url=product.catalog_url,
- attribution=WalousLandCoverService.ATTRIBUTION,
+ attribution=product.attribution,
license_note=WalousLandCoverService.LICENSE_NOTE,
legend_min_label="WALOUS klasse 1 (kunstmatige bodem)",
legend_max_label="WALOUS klasse 90 (loofbomen tot 3 m)",
included_source_values=list(WalousLandCoverService.CLASS_LABELS),
- limitation_message=f"{WalousLandCoverService.LIMITATION} {product.accuracy_label}.",
+ limitation_message=" ".join(
+ part
+ for part in (
+ WalousLandCoverService.LIMITATION,
+ f"{product.accuracy_label}.",
+ product.comparability_note,
+ )
+ if part
+ ),
coverage_zones=["wallonia"],
configured=configured,
status="configured" if configured else "source_not_provisioned",
@@ -179,21 +265,46 @@ class WalousLandCoverService:
@staticmethod
def _scope_geometry(db, project_id: UUID, payload: ThematicRasterAcquireRequest):
if not db.get(Project, project_id):
- raise AppError(code="PROJECT_NOT_FOUND", message="Project not found", status_code=404)
+ raise AppError(
+ code="PROJECT_NOT_FOUND", message="Project not found", status_code=404
+ )
if payload.bbox.crs.upper() != "EPSG:4326":
- raise AppError(code="INVALID_BBOX_CRS", message="WALOUS acquisition requires EPSG:4326", status_code=400)
- values = [payload.bbox.min_x, payload.bbox.min_y, payload.bbox.max_x, payload.bbox.max_y]
- if not all(math.isfinite(value) for value in values) or values[0] >= values[2] or values[1] >= values[3]:
- raise AppError(code="INVALID_BBOX", message="WALOUS selection must be a finite non-empty rectangle", status_code=400)
+ raise AppError(
+ code="INVALID_BBOX_CRS",
+ message="WALOUS acquisition requires EPSG:4326",
+ status_code=400,
+ )
+ values = [
+ payload.bbox.min_x,
+ payload.bbox.min_y,
+ payload.bbox.max_x,
+ payload.bbox.max_y,
+ ]
+ if (
+ not all(math.isfinite(value) for value in values)
+ or values[0] >= values[2]
+ or values[1] >= values[3]
+ ):
+ raise AppError(
+ code="INVALID_BBOX",
+ message="WALOUS selection must be a finite non-empty rectangle",
+ status_code=400,
+ )
selection = box(*values)
if payload.area_id is None:
return selection, values
area = db.get(Area, payload.area_id)
if area is None or area.project_id != project_id:
- raise AppError(code="AREA_NOT_FOUND", message="Area not found", status_code=404)
+ raise AppError(
+ code="AREA_NOT_FOUND", message="Area not found", status_code=404
+ )
selection = selection.intersection(to_shape(area.geometry))
if selection.is_empty or selection.area <= 0:
- raise AppError(code="WALOUS_SELECTION_OUTSIDE_AREA", message="Selection does not overlap the selected work area", status_code=422)
+ raise AppError(
+ code="WALOUS_SELECTION_OUTSIDE_AREA",
+ message="Selection does not overlap the selected work area",
+ status_code=422,
+ )
return selection, values
@staticmethod
@@ -201,6 +312,7 @@ class WalousLandCoverService:
source_path: Path,
scope_4326,
settings: Settings,
+ product: WalousProduct,
) -> tuple[bytes, dict[str, Any]]:
try:
import numpy as np
@@ -211,20 +323,45 @@ class WalousLandCoverService:
from rasterio.transform import from_bounds
from rasterio.windows import from_bounds as window_from_bounds
except ImportError as exc:
- raise AppError(code="RASTER_PROCESSING_UNAVAILABLE", message="Rasterio and numpy are required for WALOUS", status_code=503) from exc
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio and numpy are required for WALOUS",
+ status_code=503,
+ ) from exc
resolution = float(settings.walous_analysis_resolution_m)
- transformer = Transformer.from_crs("EPSG:4326", WalousLandCoverService.SOURCE_CRS, always_xy=True)
+ transformer = Transformer.from_crs(
+ "EPSG:4326", WalousLandCoverService.SOURCE_CRS, always_xy=True
+ )
scope_metric = shapely_transform(transformer.transform, scope_4326)
try:
with rasterio.open(source_path) as source:
- if source.crs is None or source.crs.to_epsg() != 3812 or source.count != 1:
- raise AppError(code="WALOUS_SOURCE_INVALID", message="WALOUS source must be a one-band EPSG:3812 raster", status_code=409)
- if not all(math.isclose(abs(float(value)), 1.0, abs_tol=0.05) for value in source.res):
- raise AppError(code="WALOUS_SOURCE_INVALID", message="WALOUS source must retain the official 1 m resolution", status_code=409)
+ if (
+ source.crs is None
+ or source.crs.to_epsg() != 3812
+ or source.count != 1
+ ):
+ raise AppError(
+ code="WALOUS_SOURCE_INVALID",
+ message="WALOUS source must be a one-band EPSG:3812 raster",
+ status_code=409,
+ )
+ if not all(
+ math.isclose(abs(float(value)), 1.0, abs_tol=0.05)
+ for value in source.res
+ ):
+ raise AppError(
+ code="WALOUS_SOURCE_INVALID",
+ message="WALOUS source must retain the official 1 m resolution",
+ status_code=409,
+ )
clipped_geometry = scope_metric.intersection(box(*source.bounds))
if clipped_geometry.is_empty or clipped_geometry.area <= 0:
- raise AppError(code="WALOUS_SELECTION_OUTSIDE_COVERAGE", message="Selection does not overlap WALOUS coverage", status_code=422)
+ raise AppError(
+ code="WALOUS_SELECTION_OUTSIDE_COVERAGE",
+ message="Selection does not overlap WALOUS coverage",
+ status_code=422,
+ )
min_x, min_y, max_x, max_y = clipped_geometry.bounds
bounds = (
math.floor(min_x / resolution) * resolution,
@@ -233,20 +370,45 @@ class WalousLandCoverService:
math.ceil(max_y / resolution) * resolution,
)
width_m, height_m = bounds[2] - bounds[0], bounds[3] - bounds[1]
- if width_m > settings.walous_max_side_m or height_m > settings.walous_max_side_m:
+ if (
+ width_m > settings.walous_max_side_m
+ or height_m > settings.walous_max_side_m
+ ):
raise AppError(
code="WALOUS_SELECTION_TOO_LARGE",
message=f"Select no more than {settings.walous_max_side_m:g} by {settings.walous_max_side_m:g} metres",
details={"width_m": width_m, "height_m": height_m},
status_code=422,
)
- width, height = max(1, round(width_m / resolution)), max(1, round(height_m / resolution))
+ width, height = (
+ max(1, round(width_m / resolution)),
+ max(1, round(height_m / resolution)),
+ )
if width * height > settings.walous_max_pixels:
- raise AppError(code="WALOUS_SELECTION_TOO_LARGE", message="WALOUS selection exceeds the configured cell limit", details={"pixel_count": width * height, "max_pixels": settings.walous_max_pixels}, status_code=422)
+ raise AppError(
+ code="WALOUS_SELECTION_TOO_LARGE",
+ message="WALOUS selection exceeds the configured cell limit",
+ details={
+ "pixel_count": width * height,
+ "max_pixels": settings.walous_max_pixels,
+ },
+ status_code=422,
+ )
window = window_from_bounds(*bounds, transform=source.transform)
- band = source.read(1, window=window, out_shape=(height, width), masked=True, resampling=Resampling.nearest)
+ band = source.read(
+ 1,
+ window=window,
+ out_shape=(height, width),
+ masked=True,
+ resampling=Resampling.nearest,
+ )
output_transform = from_bounds(*bounds, width, height)
- outside_scope = geometry_mask([mapping(clipped_geometry)], out_shape=(height, width), transform=output_transform, invert=False)
+ outside_scope = geometry_mask(
+ [mapping(clipped_geometry)],
+ out_shape=(height, width),
+ transform=output_transform,
+ invert=False,
+ )
# The official 2023 GeoTIFF is signed int8 while GDAL exposes
# its nodata sentinel as 255. Filling before widening would
# therefore reject the sentinel as out of range for int8.
@@ -255,11 +417,18 @@ class WalousLandCoverService:
if source.nodata is not None:
invalid |= np.isclose(raw.astype("float64"), float(source.nodata))
raw[invalid] = WalousLandCoverService.NODATA
- valid = raw[raw != WalousLandCoverService.NODATA]
- if valid.size == 0:
- raise AppError(code="WALOUS_NO_VALID_DATA", message="WALOUS contains no valid cells in this selection", status_code=422)
- classes = set(np.unique(valid).astype(int).tolist())
- unexpected = sorted(classes - set(WalousLandCoverService.CLASS_LABELS))
+ source_valid = raw[raw != WalousLandCoverService.NODATA]
+ if source_valid.size == 0:
+ raise AppError(
+ code="WALOUS_NO_VALID_DATA",
+ message="WALOUS contains no valid cells in this selection",
+ status_code=422,
+ )
+ source_classes = set(np.unique(source_valid).astype(int).tolist())
+ governed_source_classes = set(
+ product.raw_class_crosswalk or WalousLandCoverService.CLASS_LABELS
+ )
+ unexpected = sorted(source_classes - governed_source_classes)
if unexpected:
raise AppError(
code="WALOUS_SOURCE_INVALID_VALUES",
@@ -267,6 +436,18 @@ class WalousLandCoverService:
details={"unexpected_classes": unexpected},
status_code=409,
)
+ if product.raw_class_crosswalk:
+ normalized = np.full(
+ raw.shape, WalousLandCoverService.NODATA, dtype="uint8"
+ )
+ for (
+ source_value,
+ normalized_value,
+ ) in product.raw_class_crosswalk.items():
+ normalized[(raw == source_value) & ~invalid] = normalized_value
+ raw = normalized
+ valid = raw[raw != WalousLandCoverService.NODATA]
+ classes = set(np.unique(valid).astype(int).tolist())
profile = {
"driver": "GTiff",
"width": width,
@@ -288,50 +469,87 @@ class WalousLandCoverService:
"height": height,
"valid_pixel_count": int(valid.size),
"classes_present": sorted(classes),
+ "source_classes_present": sorted(source_classes),
+ "class_crosswalk": product.raw_class_crosswalk,
"bbox_epsg3812": list(bounds),
"source_width": int(source.width),
"source_height": int(source.height),
- "source_nodata": None if source.nodata is None else float(source.nodata),
+ "source_nodata": None
+ if source.nodata is None
+ else float(source.nodata),
"source_resolution_m": 1.0,
"analysis_resolution_m": resolution,
}
except AppError:
raise
except Exception as exc:
- raise AppError(code="WALOUS_SOURCE_READ_FAILED", message="The provisioned WALOUS source could not be read", details={"reason": str(exc)}, status_code=500) from exc
+ raise AppError(
+ code="WALOUS_SOURCE_READ_FAILED",
+ message="The provisioned WALOUS source could not be read",
+ details={"reason": str(exc)},
+ status_code=500,
+ ) from exc
@staticmethod
def _cached_dataset(db, project_id: UUID, filename: str) -> Dataset | None:
candidate = (
db.query(Dataset)
- .filter(Dataset.project_id == project_id, Dataset.name == filename, Dataset.source_name == WalousLandCoverService.PROVIDER, Dataset.status == "ready")
+ .filter(
+ Dataset.project_id == project_id,
+ Dataset.name == filename,
+ Dataset.source_name == WalousLandCoverService.PROVIDER,
+ Dataset.status == "ready",
+ )
.order_by(Dataset.imported_at.desc())
.first()
)
- return candidate if candidate and candidate.storage_path and Path(candidate.storage_path).is_file() else None
+ return (
+ candidate
+ if candidate
+ and candidate.storage_path
+ and Path(candidate.storage_path).is_file()
+ else None
+ )
@staticmethod
- def acquire(db, project_id: UUID, payload: ThematicRasterAcquireRequest, *, settings: Settings | None = None) -> dict[str, Any]:
+ def acquire(
+ db,
+ project_id: UUID,
+ payload: ThematicRasterAcquireRequest,
+ *,
+ settings: Settings | None = None,
+ ) -> dict[str, Any]:
resolved = settings or get_settings()
if not resolved.walous_enabled:
- raise AppError(code="WALOUS_NOT_CONFIGURED", message="WALOUS bounded analysis is disabled", status_code=503)
+ raise AppError(
+ code="WALOUS_NOT_CONFIGURED",
+ message="WALOUS bounded analysis is disabled",
+ status_code=503,
+ )
product = WalousLandCoverService._product(payload.product_key)
source_path = WalousLandCoverService._source_path(resolved, product)
if not source_path.is_file():
raise AppError(
code="WALOUS_SOURCE_NOT_PROVISIONED",
message="The official WALOUS source archive has not been provisioned on this runtime",
- details={"expected_path": str(source_path), "operator_command": "python scripts/provision_walous_sources.py --years 2020 2023"},
+ details={
+ "expected_path": str(source_path),
+ "operator_command": "python scripts/provision_walous_sources.py --years 2018 2020 2023",
+ },
status_code=503,
)
- scope, bbox_4326 = WalousLandCoverService._scope_geometry(db, project_id, payload)
+ scope, bbox_4326 = WalousLandCoverService._scope_geometry(
+ db, project_id, payload
+ )
identity = {
"product_key": product.key,
"bbox_epsg4326": [round(float(value), 8) for value in bbox_4326],
"area_id": str(payload.area_id) if payload.area_id else None,
"analysis_resolution_m": resolved.walous_analysis_resolution_m,
}
- request_hash = hashlib.sha256(json.dumps(identity, sort_keys=True).encode()).hexdigest()
+ request_hash = hashlib.sha256(
+ json.dumps(identity, sort_keys=True).encode()
+ ).hexdigest()
filename = f"walous_{product.observation_year}_{request_hash[:12]}_3812.tif"
if not payload.force_refresh:
cached = WalousLandCoverService._cached_dataset(db, project_id, filename)
@@ -345,7 +563,12 @@ class WalousLandCoverService:
display_name=product.display_name,
theme=WalousLandCoverService.THEME,
metric_kind=WalousLandCoverService.METRIC_KIND,
- resolution_m=float(metadata.get("analysis_resolution_m", resolved.walous_analysis_resolution_m)),
+ resolution_m=float(
+ metadata.get(
+ "analysis_resolution_m",
+ resolved.walous_analysis_resolution_m,
+ )
+ ),
width=int((cached.metadata_json or {}).get("width", 0)),
height=int((cached.metadata_json or {}).get("height", 0)),
valid_pixel_count=int(metadata.get("valid_pixel_count", 0)),
@@ -353,16 +576,39 @@ class WalousLandCoverService:
bbox_epsg3812=list(metadata.get("bbox_epsg3812") or []),
observation_year=product.observation_year,
source_value_unit=WalousLandCoverService.SOURCE_VALUE_UNIT,
- attribution=WalousLandCoverService.ATTRIBUTION,
- limitation_message=f"{WalousLandCoverService.LIMITATION} {product.accuracy_label}.",
+ attribution=product.attribution,
+ limitation_message=" ".join(
+ part
+ for part in (
+ WalousLandCoverService.LIMITATION,
+ f"{product.accuracy_label}.",
+ product.comparability_note,
+ )
+ if part
+ ),
).model_dump(mode="json")
- content, validation = WalousLandCoverService._read_source_window(source_path, scope, resolved)
+ content, validation = WalousLandCoverService._read_source_window(
+ source_path, scope, resolved, product
+ )
source_sha256_path = source_path.with_name(product.source_sha256_filename)
- source_sha256 = source_sha256_path.read_text(encoding="ascii").strip().split()[0] if source_sha256_path.is_file() else None
+ source_sha256 = (
+ source_sha256_path.read_text(encoding="ascii").strip().split()[0]
+ if source_sha256_path.is_file()
+ else None
+ )
acquired_at = datetime.now(UTC)
observed_at = product.observation_end
- spatial_series_hash = hashlib.sha256(json.dumps({"bbox": identity["bbox_epsg4326"], "area_id": identity["area_id"], "resolution": identity["analysis_resolution_m"]}, sort_keys=True).encode()).hexdigest()[:24]
+ spatial_series_hash = hashlib.sha256(
+ json.dumps(
+ {
+ "bbox": identity["bbox_epsg4326"],
+ "area_id": identity["area_id"],
+ "resolution": identity["analysis_resolution_m"],
+ },
+ sort_keys=True,
+ ).encode()
+ ).hexdigest()[:24]
dataset = DatasetService.import_raster_bytes(
db,
project_id=project_id,
@@ -394,14 +640,24 @@ class WalousLandCoverService:
"observation_end": product.observation_end.isoformat(),
"valid_pixel_count": validation["valid_pixel_count"],
"classes_present": validation["classes_present"],
+ "source_classes_present": validation["source_classes_present"],
+ "class_crosswalk": validation["class_crosswalk"],
"bbox_epsg4326": bbox_4326,
"bbox_epsg3812": validation["bbox_epsg3812"],
"coverage_zones": ["wallonia"],
"catalog_url": product.catalog_url,
"download_url": product.download_url,
- "attribution": WalousLandCoverService.ATTRIBUTION,
+ "attribution": product.attribution,
"license_note": WalousLandCoverService.LICENSE_NOTE,
- "limitation_message": f"{WalousLandCoverService.LIMITATION} {product.accuracy_label}.",
+ "limitation_message": " ".join(
+ part
+ for part in (
+ WalousLandCoverService.LIMITATION,
+ f"{product.accuracy_label}.",
+ product.comparability_note,
+ )
+ if part
+ ),
},
provenance_metadata={
"acquisition": "operator_provisioned_official_archive_bounded_window",
@@ -430,77 +686,169 @@ class WalousLandCoverService:
bbox_epsg3812=validation["bbox_epsg3812"],
observation_year=product.observation_year,
source_value_unit=WalousLandCoverService.SOURCE_VALUE_UNIT,
- attribution=WalousLandCoverService.ATTRIBUTION,
- limitation_message=f"{WalousLandCoverService.LIMITATION} {product.accuracy_label}.",
+ attribution=product.attribution,
+ limitation_message=" ".join(
+ part
+ for part in (
+ WalousLandCoverService.LIMITATION,
+ f"{product.accuracy_label}.",
+ product.comparability_note,
+ )
+ if part
+ ),
).model_dump(mode="json")
@staticmethod
- def _load_dataset(db, project_id: UUID, dataset_id: UUID) -> tuple[Dataset, WalousProduct]:
+ def _load_dataset(
+ db, project_id: UUID, dataset_id: UUID
+ ) -> tuple[Dataset, WalousProduct]:
dataset = db.get(Dataset, dataset_id)
if not dataset or dataset.project_id != project_id:
- raise AppError(code="DATASET_NOT_FOUND", message="Dataset not found", status_code=404)
- if dataset.dataset_type != "raster" or dataset.source_name != WalousLandCoverService.PROVIDER:
- raise AppError(code="INVALID_WALOUS_DATASET", message="WALOUS analysis requires a governed WALOUS raster", status_code=400)
- if dataset.status != "ready" or not dataset.storage_path or not Path(dataset.storage_path).is_file():
- raise AppError(code="DATASET_FILE_MISSING", message="Persisted WALOUS raster is unavailable", status_code=404)
- product = WalousLandCoverService._product(str((dataset.source_metadata or {}).get("product_key") or ""))
+ raise AppError(
+ code="DATASET_NOT_FOUND", message="Dataset not found", status_code=404
+ )
+ if (
+ dataset.dataset_type != "raster"
+ or dataset.source_name != WalousLandCoverService.PROVIDER
+ ):
+ raise AppError(
+ code="INVALID_WALOUS_DATASET",
+ message="WALOUS analysis requires a governed WALOUS raster",
+ status_code=400,
+ )
+ if (
+ dataset.status != "ready"
+ or not dataset.storage_path
+ or not Path(dataset.storage_path).is_file()
+ ):
+ raise AppError(
+ code="DATASET_FILE_MISSING",
+ message="Persisted WALOUS raster is unavailable",
+ status_code=404,
+ )
+ product = WalousLandCoverService._product(
+ str((dataset.source_metadata or {}).get("product_key") or "")
+ )
return dataset, product
@staticmethod
- def _analysis_geometry(db, project_id: UUID, payload: ThematicRasterSelectionRequest):
- selection = box(payload.bbox.min_x, payload.bbox.min_y, payload.bbox.max_x, payload.bbox.max_y)
+ def _analysis_geometry(
+ db, project_id: UUID, payload: ThematicRasterSelectionRequest
+ ):
+ selection = box(
+ payload.bbox.min_x,
+ payload.bbox.min_y,
+ payload.bbox.max_x,
+ payload.bbox.max_y,
+ )
if payload.area_id is None:
return selection
area = db.get(Area, payload.area_id)
if area is None or area.project_id != project_id:
- raise AppError(code="AREA_NOT_FOUND", message="Area not found", status_code=404)
+ raise AppError(
+ code="AREA_NOT_FOUND", message="Area not found", status_code=404
+ )
selection = selection.intersection(to_shape(area.geometry))
if selection.is_empty or selection.area <= 0:
- raise AppError(code="WALOUS_SELECTION_OUTSIDE_AREA", message="Selection does not overlap the selected work area", status_code=422)
+ raise AppError(
+ code="WALOUS_SELECTION_OUTSIDE_AREA",
+ message="Selection does not overlap the selected work area",
+ status_code=422,
+ )
return selection
@staticmethod
- def analyze(db, project_id: UUID, dataset_id: UUID, payload: ThematicRasterSelectionRequest) -> dict[str, Any]:
- dataset, product = WalousLandCoverService._load_dataset(db, project_id, dataset_id)
- selection_4326 = WalousLandCoverService._analysis_geometry(db, project_id, payload)
+ def analyze(
+ db, project_id: UUID, dataset_id: UUID, payload: ThematicRasterSelectionRequest
+ ) -> dict[str, Any]:
+ dataset, product = WalousLandCoverService._load_dataset(
+ db, project_id, dataset_id
+ )
+ selection_4326 = WalousLandCoverService._analysis_geometry(
+ db, project_id, payload
+ )
try:
import numpy as np
import rasterio
from rasterio.features import geometry_mask
from rasterio.mask import mask
except ImportError as exc:
- raise AppError(code="RASTER_PROCESSING_UNAVAILABLE", message="Rasterio and numpy are required for WALOUS analysis", status_code=503) from exc
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio and numpy are required for WALOUS analysis",
+ status_code=503,
+ ) from exc
try:
with rasterio.open(dataset.storage_path) as source:
- transformer = Transformer.from_crs("EPSG:4326", source.crs, always_xy=True)
- selection_metric = shapely_transform(transformer.transform, selection_4326)
+ transformer = Transformer.from_crs(
+ "EPSG:4326", source.crs, always_xy=True
+ )
+ selection_metric = shapely_transform(
+ transformer.transform, selection_4326
+ )
geometry = selection_metric.intersection(box(*source.bounds))
if geometry.is_empty or geometry.area <= 0:
- raise AppError(code="WALOUS_SELECTION_OUTSIDE_DATASET", message="Selection does not overlap the persisted WALOUS raster", status_code=422)
- clipped, transform = mask(source, [mapping(geometry)], crop=True, filled=False, indexes=[1])
+ raise AppError(
+ code="WALOUS_SELECTION_OUTSIDE_DATASET",
+ message="Selection does not overlap the persisted WALOUS raster",
+ status_code=422,
+ )
+ clipped, transform = mask(
+ source, [mapping(geometry)], crop=True, filled=False, indexes=[1]
+ )
band = np.ma.asarray(clipped[0])
raw = np.asarray(np.ma.getdata(band), dtype="uint8")
- selected = geometry_mask([mapping(geometry)], out_shape=raw.shape, transform=transform, invert=True)
- valid = selected & ~np.ma.getmaskarray(band) & (raw != WalousLandCoverService.NODATA)
+ selected = geometry_mask(
+ [mapping(geometry)],
+ out_shape=raw.shape,
+ transform=transform,
+ invert=True,
+ )
+ valid = (
+ selected
+ & ~np.ma.getmaskarray(band)
+ & (raw != WalousLandCoverService.NODATA)
+ )
values = raw[valid]
selected_count = int(selected.sum())
valid_count = int(values.size)
if not valid_count:
- raise AppError(code="WALOUS_NO_VALID_DATA", message="WALOUS contains no valid cells in this selection", status_code=422)
+ raise AppError(
+ code="WALOUS_NO_VALID_DATA",
+ message="WALOUS contains no valid cells in this selection",
+ status_code=422,
+ )
cell_area_m2 = abs(float(source.res[0]) * float(source.res[1]))
except AppError:
raise
except Exception as exc:
- raise AppError(code="WALOUS_ANALYSIS_FAILED", message="The persisted WALOUS raster could not be analysed", details={"reason": str(exc)}, status_code=500) from exc
+ raise AppError(
+ code="WALOUS_ANALYSIS_FAILED",
+ message="The persisted WALOUS raster could not be analysed",
+ details={"reason": str(exc)},
+ status_code=500,
+ ) from exc
def area_for(classes: set[int]) -> float:
- return float(np.count_nonzero(np.isin(values, list(classes))) * cell_area_m2 / 10_000.0)
+ return float(
+ np.count_nonzero(np.isin(values, list(classes)))
+ * cell_area_m2
+ / 10_000.0
+ )
metric_specs = [
- ("land_cover_observed_area_ha", "Gekarteerde landbedekking", set(WalousLandCoverService.CLASS_LABELS)),
+ (
+ "land_cover_observed_area_ha",
+ "Gekarteerde landbedekking",
+ set(WalousLandCoverService.CLASS_LABELS),
+ ),
("forest_cover_area_ha", "Boom- en bosbedekking", {8, 9, 80, 90}),
("surface_water_area_ha", "Oppervlaktewater", {5}),
- ("artificial_cover_area_ha", "Kunstmatige bedekking en constructies", {1, 2, 3}),
+ (
+ "artificial_cover_area_ha",
+ "Kunstmatige bedekking en constructies",
+ {1, 2, 3},
+ ),
("annual_herbaceous_cover_area_ha", "Jaarlijks wisselende kruidlaag", {6}),
("permanent_herbaceous_cover_area_ha", "Jaarronde kruidlaag", {7}),
("bare_soil_area_ha", "Kale bodem", {4}),
@@ -537,25 +885,52 @@ class WalousLandCoverService:
primary_metric_key=primary.metric_key,
metrics=metrics,
),
- unsupported_metrics=["legal_land_use", "ownership", "tree_count", "timber_volume", "water_volume"],
- limitation_message=f"{WalousLandCoverService.LIMITATION} {product.accuracy_label}.",
+ unsupported_metrics=[
+ "legal_land_use",
+ "ownership",
+ "tree_count",
+ "timber_volume",
+ "water_volume",
+ ],
+ limitation_message=" ".join(
+ part
+ for part in (
+ WalousLandCoverService.LIMITATION,
+ f"{product.accuracy_label}.",
+ product.comparability_note,
+ )
+ if part
+ ),
generated_at=datetime.now(UTC).isoformat(),
).model_dump(mode="json")
@staticmethod
- def render_png(db, project_id: UUID, dataset_id: UUID, *, max_dimension: int = 1800) -> bytes:
- dataset, _product = WalousLandCoverService._load_dataset(db, project_id, dataset_id)
+ def render_png(
+ db, project_id: UUID, dataset_id: UUID, *, max_dimension: int = 1800
+ ) -> bytes:
+ dataset, _product = WalousLandCoverService._load_dataset(
+ db, project_id, dataset_id
+ )
try:
import numpy as np
import rasterio
from PIL import Image
from rasterio.enums import Resampling
except ImportError as exc:
- raise AppError(code="RASTER_PROCESSING_UNAVAILABLE", message="Rasterio, numpy and Pillow are required for WALOUS rendering", status_code=503) from exc
+ raise AppError(
+ code="RASTER_PROCESSING_UNAVAILABLE",
+ message="Rasterio, numpy and Pillow are required for WALOUS rendering",
+ status_code=503,
+ ) from exc
with rasterio.open(dataset.storage_path) as source:
scale = min(1.0, max_dimension / max(source.width, source.height))
- width, height = max(1, round(source.width * scale)), max(1, round(source.height * scale))
- values = source.read(1, out_shape=(height, width), masked=True, resampling=Resampling.nearest)
+ width, height = (
+ max(1, round(source.width * scale)),
+ max(1, round(source.height * scale)),
+ )
+ values = source.read(
+ 1, out_shape=(height, width), masked=True, resampling=Resampling.nearest
+ )
raw = np.asarray(np.ma.getdata(values), dtype="uint8")
rgba = np.zeros((height, width, 4), dtype="uint8")
for value, color in WalousLandCoverService.CLASS_COLORS.items():
diff --git a/backend/tests/test_sprint193_end_user_workbench.py b/backend/tests/test_sprint193_end_user_workbench.py
index 72b05234..f256f399 100644
--- a/backend/tests/test_sprint193_end_user_workbench.py
+++ b/backend/tests/test_sprint193_end_user_workbench.py
@@ -59,6 +59,8 @@ def test_configured_yolo_and_active_asset_are_selected_without_hiding_limitation
assert "getYoloPreflight" in hook
assert 'aria-label="Status gebouwdetectie"' in lab
assert "Nog niet nationaal gevalideerd" in lab
+ assert "selectedDetectionModel?.nationally_validated !== true" in lab
+ assert "selectedDetectionModel?.validation_scope" in lab
assert "vereisen lokale referentiedata en QA" in lab
assert "Modelkalibratie voor beheerders" in lab
diff --git a/backend/tests/test_sprint8b_yolo_foundation.py b/backend/tests/test_sprint8b_yolo_foundation.py
index 49222dd5..bf447f6d 100644
--- a/backend/tests/test_sprint8b_yolo_foundation.py
+++ b/backend/tests/test_sprint8b_yolo_foundation.py
@@ -232,6 +232,10 @@ def test_yolo_configured_model_reports_configured_with_local_model_and_dependenc
assert model.configured is True
assert model.status == "configured"
assert model.version == settings.yolo_model_version
+ assert model.nationally_validated is False
+ assert model.operator_review_required is True
+ assert model.validated_regions == ["flanders_mol_kempen"]
+ assert "Mol and the Kempen" in (model.validation_scope or "")
def test_yolo_dependency_check_uses_real_imports_not_find_spec() -> None:
diff --git a/backend/tests/test_spw_terrain_service.py b/backend/tests/test_spw_terrain_service.py
new file mode 100644
index 00000000..db70da02
--- /dev/null
+++ b/backend/tests/test_spw_terrain_service.py
@@ -0,0 +1,268 @@
+from __future__ import annotations
+
+from pathlib import Path
+from types import SimpleNamespace
+from uuid import uuid4
+
+import numpy as np
+from fastapi.testclient import TestClient
+from pyproj import Transformer
+import rasterio
+from rasterio.transform import from_origin
+
+from app.core.config import Settings
+from app.db.session import get_db
+from app.main import app
+from app.models import Dataset, Job, Project
+from app.schemas.dhmv import TerrainSelectionRequest
+from app.schemas.spw_terrain import SpwTerrainAcquireRequest
+from app.services.dataset_service import DatasetService
+from app.services.spw_terrain_service import SpwTerrainService
+from app.services.terrain_analysis_service import TerrainAnalysisService
+
+
+class FakeQuery:
+ def filter(self, *_args):
+ return self
+
+ def order_by(self, *_args):
+ return self
+
+ def first(self):
+ return None
+
+
+class FakeSession:
+ def __init__(self, project, dataset=None):
+ self.project = project
+ self.dataset = dataset
+ self.added = []
+
+ def get(self, model, row_id):
+ if model is Project and row_id == self.project.id:
+ return self.project
+ if model is Dataset and self.dataset is not None and row_id == self.dataset.id:
+ return self.dataset
+ return next(
+ (
+ item
+ for item in self.added
+ if isinstance(item, model) and item.id == row_id
+ ),
+ None,
+ )
+
+ def query(self, _model):
+ return FakeQuery()
+
+ def add(self, row):
+ self.added.append(row)
+
+ def commit(self):
+ return None
+
+ def rollback(self):
+ return None
+
+ def refresh(self, row):
+ return row
+
+
+def settings(source_dir: Path) -> Settings:
+ return Settings(
+ _env_file=None,
+ SPW_TERRAIN_SOURCE_DIR=str(source_dir),
+ SPW_TERRAIN_ANALYSIS_RESOLUTION_M=5,
+ SPW_TERRAIN_MAX_SIDE_M=20_000,
+ SPW_TERRAIN_MAX_PIXELS=1_000_000,
+ DHMV_MAX_PIXELS=1_000_000,
+ )
+
+
+def make_source(path: Path) -> list[float]:
+ to_3812 = Transformer.from_crs("EPSG:4326", "EPSG:3812", always_xy=True)
+ to_4326 = Transformer.from_crs("EPSG:3812", "EPSG:4326", always_xy=True)
+ x, y = to_3812.transform(4.85, 50.45)
+ values = np.linspace(100.0, 125.0, 40_000, dtype="float32").reshape(200, 200)
+ with rasterio.open(
+ path,
+ "w",
+ driver="GTiff",
+ width=200,
+ height=200,
+ count=1,
+ dtype="float32",
+ crs="EPSG:3812",
+ transform=from_origin(x, y + 200, 1, 1),
+ nodata=-9999.0,
+ ) as target:
+ target.write(values, 1)
+ min_lon, min_lat = to_4326.transform(x, y)
+ max_lon, max_lat = to_4326.transform(x + 200, y + 200)
+ return [min_lon, min_lat, max_lon, max_lat]
+
+
+def test_spw_terrain_registry_reports_real_source_state(tmp_path: Path) -> None:
+ before = SpwTerrainService.list_products(settings=settings(tmp_path))[0]
+ assert before["status"] == "source_not_provisioned"
+
+ make_source(tmp_path / SpwTerrainService.SOURCE_FILENAME)
+ after = SpwTerrainService.list_products(settings=settings(tmp_path))[0]
+
+ assert after["configured"] is True
+ assert after["coverage_zones"] == ["wallonia"]
+ assert after["source_crs"] == "EPSG:3812"
+ assert after["vertical_reference"].endswith("(EPSG:5710)")
+
+
+def test_spw_terrain_acquisition_persists_bounded_dng_raster_and_provenance(
+ tmp_path: Path, monkeypatch
+) -> None:
+ bbox = make_source(tmp_path / SpwTerrainService.SOURCE_FILENAME)
+ project = Project(id=uuid4(), name="Belgium")
+ captured = {}
+
+ def persist(_db, **kwargs):
+ captured.update(kwargs)
+ return SimpleNamespace(id=uuid4())
+
+ monkeypatch.setattr(DatasetService, "import_raster_bytes", persist)
+ result = SpwTerrainService.acquire(
+ FakeSession(project),
+ project.id,
+ SpwTerrainAcquireRequest(
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
+ resolution_m=5,
+ force_refresh=True,
+ ),
+ settings=settings(tmp_path),
+ )
+
+ assert result["provider"] == "spw_terrain"
+ assert result["resolution_m"] == 5
+ assert result["valid_pixel_count"] > 0
+ assert captured["source_metadata"]["vertical_unit_label"] == "m DNG"
+ assert captured["source_metadata"]["coverage_zones"] == ["wallonia"]
+ assert captured["provenance_metadata"]["resampling"] == "bilinear"
+ assert captured["valid_from"].date().isoformat() == "2021-02-19"
+ with rasterio.MemoryFile(captured["content"]) as memory:
+ with memory.open() as derived:
+ assert derived.crs.to_epsg() == 3812
+ assert derived.res == (5.0, 5.0)
+ assert derived.nodata == -9999.0
+
+
+def test_terrain_analysis_preserves_spw_vertical_datum_and_limitations(
+ tmp_path: Path,
+) -> None:
+ bbox = make_source(tmp_path / "derived.tif")
+ project = Project(id=uuid4(), name="Belgium")
+ dataset = Dataset(
+ id=uuid4(),
+ project_id=project.id,
+ name="derived.tif",
+ dataset_type="raster",
+ source="SPW MNT",
+ source_name="spw_terrain",
+ source_metadata={
+ "product_key": SpwTerrainService.PRODUCT_KEY,
+ "surface_model": "terrain",
+ "vertical_reference": SpwTerrainService.VERTICAL_REFERENCE,
+ "vertical_unit_label": "m DNG",
+ "limitation_message": SpwTerrainService.LIMITATION,
+ },
+ storage_path=str(tmp_path / "derived.tif"),
+ status="ready",
+ )
+ result = TerrainAnalysisService.analyze(
+ FakeSession(project, dataset),
+ project.id,
+ dataset.id,
+ TerrainSelectionRequest(
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
+ ),
+ settings=settings(tmp_path),
+ )
+
+ assert result["vertical_reference"] == SpwTerrainService.VERTICAL_REFERENCE
+ assert result["summary"]["metric_unit"] == "m DNG"
+ assert result["summary"]["metrics"][0]["metric_unit"] == "m DNG"
+ assert result["limitation_message"] == SpwTerrainService.LIMITATION
+ assert result["unsupported_metrics"] == ["water_depth_m", "water_volume_m3"]
+
+
+def test_spw_terrain_routes_use_canonical_envelopes(monkeypatch) -> None:
+ project = Project(id=uuid4(), name="Belgium")
+ db = FakeSession(project)
+ monkeypatch.setattr(
+ SpwTerrainService,
+ "acquire",
+ lambda *_args, **_kwargs: {
+ "output_dataset_id": str(uuid4()),
+ "provider": SpwTerrainService.PROVIDER,
+ },
+ )
+ monkeypatch.setattr(
+ SpwTerrainService,
+ "list_products",
+ lambda *_args, **_kwargs: [
+ {
+ "key": SpwTerrainService.PRODUCT_KEY,
+ "display_name": SpwTerrainService.DISPLAY_NAME,
+ "surface_model": "terrain",
+ "source_filename": SpwTerrainService.SOURCE_FILENAME,
+ "native_resolution_m": 1,
+ "analysis_resolution_m": 5,
+ "source_crs": "EPSG:3812",
+ "vertical_reference": SpwTerrainService.VERTICAL_REFERENCE,
+ "acquisition_period": SpwTerrainService.ACQUISITION_PERIOD,
+ "catalog_url": SpwTerrainService.CATALOG_URL,
+ "attribution": SpwTerrainService.ATTRIBUTION,
+ "license_note": SpwTerrainService.LICENSE_NOTE,
+ "limitation_message": SpwTerrainService.LIMITATION,
+ "coverage_zones": ["wallonia"],
+ "configured": True,
+ "status": "configured",
+ }
+ ],
+ )
+ app.dependency_overrides[get_db] = lambda: db
+ try:
+ client = TestClient(app)
+ products = client.get(
+ f"/api/v1/projects/{project.id}/datasets/spw-terrain/products"
+ )
+ acquisition = client.post(
+ f"/api/v1/projects/{project.id}/datasets/spw-terrain/acquire",
+ json={
+ "bbox": {
+ "min_x": 4.8,
+ "min_y": 50.4,
+ "max_x": 4.9,
+ "max_y": 50.5,
+ "crs": "EPSG:4326",
+ },
+ "product_key": SpwTerrainService.PRODUCT_KEY,
+ },
+ )
+ finally:
+ app.dependency_overrides.clear()
+
+ assert products.status_code == 200 and products.json()["data"]["total"] == 1
+ assert (
+ acquisition.status_code == 200
+ and acquisition.json()["data"]["job_type"] == "raster.spw-terrain.acquire"
+ )
+ assert any(isinstance(item, Job) for item in db.added)
diff --git a/backend/tests/test_walous_land_cover_service.py b/backend/tests/test_walous_land_cover_service.py
index a04a06e4..6a9bccbf 100644
--- a/backend/tests/test_walous_land_cover_service.py
+++ b/backend/tests/test_walous_land_cover_service.py
@@ -16,7 +16,10 @@ from app.core.config import Settings
from app.db.session import get_db
from app.main import app
from app.models import Dataset, Job, Project
-from app.schemas.thematic_raster import ThematicRasterAcquireRequest, ThematicRasterSelectionRequest
+from app.schemas.thematic_raster import (
+ ThematicRasterAcquireRequest,
+ ThematicRasterSelectionRequest,
+)
from app.schemas.temporal import TemporalComparisonRequest
from app.services.dataset_service import DatasetService
from app.services.temporal_analysis_service import TemporalAnalysisService
@@ -24,8 +27,12 @@ from app.services.walous_land_cover_service import WalousLandCoverService
def load_provisioner():
- path = Path(__file__).resolve().parents[2] / "scripts" / "provision_walous_sources.py"
- spec = importlib.util.spec_from_file_location("walous_source_provisioner_test", path)
+ path = (
+ Path(__file__).resolve().parents[2] / "scripts" / "provision_walous_sources.py"
+ )
+ spec = importlib.util.spec_from_file_location(
+ "walous_source_provisioner_test", path
+ )
assert spec and spec.loader
module = importlib.util.module_from_spec(spec)
spec.loader.exec_module(module)
@@ -54,7 +61,14 @@ class FakeSession:
return self.project
if model is Dataset and self.dataset is not None and row_id == self.dataset.id:
return self.dataset
- match = next((item for item in self.added if isinstance(item, model) and item.id == row_id), None)
+ match = next(
+ (
+ item
+ for item in self.added
+ if isinstance(item, model) and item.id == row_id
+ ),
+ None,
+ )
if match is not None:
return match
return None
@@ -80,12 +94,13 @@ def make_source(
*,
dtype: str = "uint8",
nodata: int = 255,
+ class_codes: list[int] | None = None,
) -> tuple[list[float], np.ndarray]:
to_3812 = Transformer.from_crs("EPSG:4326", "EPSG:3812", always_xy=True)
to_4326 = Transformer.from_crs("EPSG:3812", "EPSG:4326", always_xy=True)
x, y = to_3812.transform(4.85, 50.45)
transform = from_origin(x, y + 100, 1, 1)
- class_codes = [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
+ class_codes = class_codes or [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
values = np.empty((100, len(class_codes) * 20), dtype=dtype)
for index, class_code in enumerate(class_codes):
values[:, index * 20 : (index + 1) * 20] = class_code
@@ -118,20 +133,40 @@ def settings(source_dir: Path) -> Settings:
def test_walous_registry_reports_real_provisioning_state(tmp_path: Path) -> None:
- before = {item["key"]: item for item in WalousLandCoverService.list_products(settings=settings(tmp_path))}
+ before = {
+ item["key"]: item
+ for item in WalousLandCoverService.list_products(settings=settings(tmp_path))
+ }
assert before["walous_land_cover_2023"]["status"] == "source_not_provisioned"
make_source(tmp_path / "walous_land_cover_2023_3812.tif")
- after = {item["key"]: item for item in WalousLandCoverService.list_products(settings=settings(tmp_path))}
+ after = {
+ item["key"]: item
+ for item in WalousLandCoverService.list_products(settings=settings(tmp_path))
+ }
assert after["walous_land_cover_2023"]["configured"] is True
assert after["walous_land_cover_2023"]["source_crs"] == "EPSG:3812"
assert after["walous_land_cover_2023"]["native_resolution_m"] == 1.0
assert after["walous_land_cover_2023"]["analysis_resolution_m"] == 10.0
assert after["walous_land_cover_2023"]["coverage_zones"] == ["wallonia"]
- assert after["walous_land_cover_2023"]["included_source_values"] == [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
+ assert after["walous_land_cover_2023"]["included_source_values"] == [
+ 1,
+ 2,
+ 3,
+ 4,
+ 5,
+ 6,
+ 7,
+ 8,
+ 9,
+ 80,
+ 90,
+ ]
assert after["walous_land_cover_2023"]["source_value_unit"] == "walous_class_code"
-def test_walous_provisioner_accepts_official_non_contiguous_class_codes(tmp_path: Path) -> None:
+def test_walous_provisioner_accepts_official_non_contiguous_class_codes(
+ tmp_path: Path,
+) -> None:
source_path = tmp_path / "walous_land_cover_2023_3812.tif"
make_source(source_path)
@@ -140,7 +175,86 @@ def test_walous_provisioner_accepts_official_non_contiguous_class_codes(tmp_path
assert validation["sample_classes"] == [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
-def test_walous_acquisition_reads_real_classes_and_persists_provenance(tmp_path: Path, monkeypatch) -> None:
+def test_walous_2018_registry_and_provisioner_accept_official_stacked_classes(
+ tmp_path: Path,
+) -> None:
+ source_codes = sorted(WalousLandCoverService.WALOUS_2018_CLASS_CROSSWALK)
+ source_path = tmp_path / "walous_land_cover_2018_3812.tif"
+ make_source(source_path, class_codes=source_codes)
+
+ validation = load_provisioner().validate_raster(source_path)
+ registry = {
+ item["key"]: item
+ for item in WalousLandCoverService.list_products(settings=settings(tmp_path))
+ }
+
+ assert validation["sample_classes"] == source_codes
+ assert validation["implicit_source_nodata_values"] == [0]
+ assert registry["walous_land_cover_2018"]["configured"] is True
+ assert registry["walous_land_cover_2018"]["observation_year"] == 2018
+ assert "crosswalk" in registry["walous_land_cover_2018"]["limitation_message"]
+
+
+def test_walous_2018_acquisition_normalizes_stacked_classes_with_explicit_provenance(
+ tmp_path: Path, monkeypatch
+) -> None:
+ source_codes = sorted(WalousLandCoverService.WALOUS_2018_CLASS_CROSSWALK)
+ bbox, _values = make_source(
+ tmp_path / "walous_land_cover_2018_3812.tif",
+ class_codes=source_codes,
+ )
+ project = Project(id=uuid4(), name="Belgium")
+ captured = {}
+
+ def persist(_db, **kwargs):
+ captured.update(kwargs)
+ return SimpleNamespace(id=uuid4())
+
+ monkeypatch.setattr(DatasetService, "import_raster_bytes", persist)
+ result = WalousLandCoverService.acquire(
+ FakeSession(project),
+ project.id,
+ ThematicRasterAcquireRequest(
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
+ product_key="walous_land_cover_2018",
+ force_refresh=True,
+ ),
+ settings=settings(tmp_path),
+ )
+
+ assert result["observation_year"] == 2018
+ assert captured["source_metadata"]["classes_present"] == [
+ 1,
+ 2,
+ 3,
+ 4,
+ 5,
+ 6,
+ 7,
+ 8,
+ 9,
+ 80,
+ 90,
+ ]
+ assert captured["source_metadata"]["source_classes_present"] == source_codes
+ assert captured["source_metadata"]["class_crosswalk"][62] == 2
+ assert captured["source_metadata"]["class_crosswalk"][0] == 255
+ assert (
+ captured["source_metadata"]["attribution"]
+ == "Service public de Wallonie (SPW), UCLouvain, ULB, ISSeP"
+ )
+ assert captured["observed_at"].date().isoformat() == "2018-12-31"
+
+
+def test_walous_acquisition_reads_real_classes_and_persists_provenance(
+ tmp_path: Path, monkeypatch
+) -> None:
bbox, _values = make_source(tmp_path / "walous_land_cover_2023_3812.tif")
project = Project(id=uuid4(), name="Belgium")
db = FakeSession(project)
@@ -156,7 +270,13 @@ def test_walous_acquisition_reads_real_classes_and_persists_provenance(tmp_path:
db,
project.id,
ThematicRasterAcquireRequest(
- bbox={"min_x": bbox[0], "min_y": bbox[1], "max_x": bbox[2], "max_y": bbox[3], "crs": "EPSG:4326"},
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
product_key="walous_land_cover_2023",
force_refresh=True,
),
@@ -166,7 +286,19 @@ def test_walous_acquisition_reads_real_classes_and_persists_provenance(tmp_path:
assert result["output_dataset_id"] == str(output_id)
assert result["resolution_m"] == 10
assert captured["source_name"] == "spw_walous_land_cover"
- assert captured["source_metadata"]["classes_present"] == [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
+ assert captured["source_metadata"]["classes_present"] == [
+ 1,
+ 2,
+ 3,
+ 4,
+ 5,
+ 6,
+ 7,
+ 8,
+ 9,
+ 80,
+ 90,
+ ]
assert captured["provenance_metadata"]["resampling"] == "nearest"
assert captured["temporal_series_key"].startswith("spw:walous:land-cover:")
assert captured["observed_at"].date().isoformat() == "2023-06-25"
@@ -174,7 +306,9 @@ def test_walous_acquisition_reads_real_classes_and_persists_provenance(tmp_path:
assert captured["valid_to"] == captured["observed_at"]
-def test_walous_acquisition_accepts_official_signed_int8_nodata(tmp_path: Path, monkeypatch) -> None:
+def test_walous_acquisition_accepts_official_signed_int8_nodata(
+ tmp_path: Path, monkeypatch
+) -> None:
bbox, _values = make_source(
tmp_path / "walous_land_cover_2023_3812.tif",
dtype="int8",
@@ -193,7 +327,13 @@ def test_walous_acquisition_accepts_official_signed_int8_nodata(tmp_path: Path,
FakeSession(project),
project.id,
ThematicRasterAcquireRequest(
- bbox={"min_x": bbox[0], "min_y": bbox[1], "max_x": bbox[2], "max_y": bbox[3], "crs": "EPSG:4326"},
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
product_key="walous_land_cover_2023",
force_refresh=True,
),
@@ -201,14 +341,28 @@ def test_walous_acquisition_accepts_official_signed_int8_nodata(tmp_path: Path,
)
assert result["output_dataset_id"] == str(output_id)
- assert captured["source_metadata"]["classes_present"] == [1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90]
+ assert captured["source_metadata"]["classes_present"] == [
+ 1,
+ 2,
+ 3,
+ 4,
+ 5,
+ 6,
+ 7,
+ 8,
+ 9,
+ 80,
+ 90,
+ ]
with rasterio.MemoryFile(captured["content"]) as memory:
with memory.open() as derived:
assert derived.dtypes == ("uint8",)
assert derived.nodata == 255
-def test_walous_analysis_returns_semantic_area_metrics(tmp_path: Path, monkeypatch) -> None:
+def test_walous_analysis_returns_semantic_area_metrics(
+ tmp_path: Path, monkeypatch
+) -> None:
bbox, _values = make_source(tmp_path / "walous_land_cover_2023_3812.tif")
project = Project(id=uuid4(), name="Belgium")
output_id = uuid4()
@@ -221,7 +375,13 @@ def test_walous_analysis_returns_semantic_area_metrics(tmp_path: Path, monkeypat
monkeypatch.setattr(DatasetService, "import_raster_bytes", persist)
db = FakeSession(project)
payload = ThematicRasterAcquireRequest(
- bbox={"min_x": bbox[0], "min_y": bbox[1], "max_x": bbox[2], "max_y": bbox[3], "crs": "EPSG:4326"},
+ bbox={
+ "min_x": bbox[0],
+ "min_y": bbox[1],
+ "max_x": bbox[2],
+ "max_y": bbox[3],
+ "crs": "EPSG:4326",
+ },
product_key="walous_land_cover_2023",
force_refresh=True,
)
@@ -247,7 +407,10 @@ def test_walous_analysis_returns_semantic_area_metrics(tmp_path: Path, monkeypat
output_id,
ThematicRasterSelectionRequest(bbox=payload.bbox),
)
- metrics = {item["metric_key"]: item["metric_value"] for item in result["summary"]["metrics"]}
+ metrics = {
+ item["metric_key"]: item["metric_value"]
+ for item in result["summary"]["metrics"]
+ }
assert result["metric_kind"] == "categorical_area"
assert metrics["land_cover_observed_area_ha"] > 0
@@ -270,7 +433,10 @@ def test_walous_render_png_uses_governed_class_colours(tmp_path: Path) -> None:
dataset_type="raster",
source="SPW WALOUS",
source_name="spw_walous_land_cover",
- source_metadata={"product_key": "walous_land_cover_2023", "bbox_epsg4326": bbox},
+ source_metadata={
+ "product_key": "walous_land_cover_2023",
+ "bbox_epsg4326": bbox,
+ },
storage_path=str(tmp_path / "walous_land_cover_2023_3812.tif"),
status="ready",
)
@@ -281,7 +447,9 @@ def test_walous_render_png_uses_governed_class_colours(tmp_path: Path) -> None:
assert rendered.startswith(b"\x89PNG\r\n\x1a\n")
-def test_walous_temporal_comparison_reuses_persisted_raster_metrics(monkeypatch) -> None:
+def test_walous_temporal_comparison_reuses_persisted_raster_metrics(
+ monkeypatch,
+) -> None:
project_id = uuid4()
earlier = Dataset(
id=uuid4(),
@@ -320,14 +488,16 @@ def test_walous_temporal_comparison_reuses_persisted_raster_metrics(monkeypatch)
"metric_unit": "ha",
"aggregation_method": "nearest_resampled_cells_times_cell_area",
"primary_metric_key": "land_cover_observed_area_ha",
- "metrics": [{
- "metric_key": "land_cover_observed_area_ha",
- "metric_label": "Gekarteerde landbedekking",
- "metric_value": value,
- "metric_unit": "ha",
- "aggregation_method": "nearest_resampled_cells_times_cell_area",
- "is_estimate": True,
- }],
+ "metrics": [
+ {
+ "metric_key": "land_cover_observed_area_ha",
+ "metric_label": "Gekarteerde landbedekking",
+ "metric_value": value,
+ "metric_unit": "ha",
+ "aggregation_method": "nearest_resampled_cells_times_cell_area",
+ "is_estimate": True,
+ }
+ ],
},
"limitation_message": "Cell-based estimate.",
}
@@ -336,10 +506,18 @@ def test_walous_temporal_comparison_reuses_persisted_raster_metrics(monkeypatch)
payload = TemporalComparisonRequest(
earlier_dataset_id=earlier.id,
later_dataset_id=later.id,
- bbox={"min_x": 4.8, "min_y": 50.4, "max_x": 4.9, "max_y": 50.5, "crs": "EPSG:4326"},
+ bbox={
+ "min_x": 4.8,
+ "min_y": 50.4,
+ "max_x": 4.9,
+ "max_y": 50.5,
+ "crs": "EPSG:4326",
+ },
)
- result = TemporalAnalysisService.compare(TemporalSession(), project_id=project_id, payload=payload)
+ result = TemporalAnalysisService.compare(
+ TemporalSession(), project_id=project_id, payload=payload
+ )
assert result.metric.earlier_value == 4.0
assert result.metric.later_value == 5.5
@@ -354,7 +532,10 @@ def test_walous_api_routes_use_canonical_envelopes(monkeypatch) -> None:
monkeypatch.setattr(
WalousLandCoverService,
"acquire",
- lambda *_args, **_kwargs: {"output_dataset_id": str(dataset_id), "provider": WalousLandCoverService.PROVIDER},
+ lambda *_args, **_kwargs: {
+ "output_dataset_id": str(dataset_id),
+ "provider": WalousLandCoverService.PROVIDER,
+ },
)
monkeypatch.setattr(
WalousLandCoverService,
@@ -364,7 +545,13 @@ def test_walous_api_routes_use_canonical_envelopes(monkeypatch) -> None:
"product_key": "walous_land_cover_2023",
"theme": "land_cover_use",
"metric_kind": "categorical_area",
- "selection_bbox": {"min_x": 4.8, "min_y": 50.4, "max_x": 4.9, "max_y": 50.5, "crs": "EPSG:4326"},
+ "selection_bbox": {
+ "min_x": 4.8,
+ "min_y": 50.4,
+ "max_x": 4.9,
+ "max_y": 50.5,
+ "crs": "EPSG:4326",
+ },
"selected_cell_count": 100,
"valid_cell_count": 100,
"coverage_ratio": 1.0,
@@ -390,20 +577,37 @@ def test_walous_api_routes_use_canonical_envelopes(monkeypatch) -> None:
acquisition = client.post(
f"/api/v1/projects/{project.id}/datasets/walous/acquire",
json={
- "bbox": {"min_x": 4.8, "min_y": 50.4, "max_x": 4.9, "max_y": 50.5, "crs": "EPSG:4326"},
+ "bbox": {
+ "min_x": 4.8,
+ "min_y": 50.4,
+ "max_x": 4.9,
+ "max_y": 50.5,
+ "crs": "EPSG:4326",
+ },
"product_key": "walous_land_cover_2023",
},
)
selection = client.post(
f"/api/v1/projects/{project.id}/datasets/{dataset_id}/raster/walous/select",
- json={"bbox": {"min_x": 4.8, "min_y": 50.4, "max_x": 4.9, "max_y": 50.5, "crs": "EPSG:4326"}},
+ json={
+ "bbox": {
+ "min_x": 4.8,
+ "min_y": 50.4,
+ "max_x": 4.9,
+ "max_y": 50.5,
+ "crs": "EPSG:4326",
+ }
+ },
)
finally:
app.dependency_overrides.clear()
assert products.status_code == 200 and set(products.json()) == {"data"}
- assert products.json()["data"]["total"] == 2
+ assert products.json()["data"]["total"] == 3
assert acquisition.status_code == 200 and set(acquisition.json()) == {"data"}
assert acquisition.json()["data"]["job_type"] == "raster.walous.acquire"
- assert selection.status_code == 200 and selection.json()["data"]["theme"] == "land_cover_use"
+ assert (
+ selection.status_code == 200
+ and selection.json()["data"]["theme"] == "land_cover_use"
+ )
assert any(isinstance(item, Job) for item in db.added)
diff --git a/deploy/unraid/geointel-unraid-template.xml b/deploy/unraid/geointel-unraid-template.xml
index d2afde7f..2d5979c5 100644
--- a/deploy/unraid/geointel-unraid-template.xml
+++ b/deploy/unraid/geointel-unraid-template.xml
@@ -105,11 +105,16 @@
30000000
300
160
- true
+ true
/app/storage/source-cache/walous
10
60000
36000000
+ true
+ /app/storage/source-cache/spw-terrain
+ 5
+ 20000
+ 12000000
false
/app/models
{detectionQualityInterpretation(selectedOperatorProfile?.f1)}
- {selectedOperatorProfile && !selectedOperatorProfile.nationallyValidated ? (
+ {selectedOperatorProfile && selectedDetectionModel?.nationally_validated !== true ? (
Nog niet nationaal gevalideerd
- Dit model is operationeel voor gecontroleerde beeldanalyse, maar de gemeten kwaliteit geldt alleen voor {selectedOperatorProfile.validationScope}.
+ Dit model is operationeel voor gecontroleerde beeldanalyse, maar de gemeten kwaliteit geldt alleen voor {selectedDetectionModel?.validation_scope ?? selectedOperatorProfile.validationScope}.
Resultaten elders in Belgie of op zee vereisen lokale referentiedata en QA voordat ze als betrouwbaar kunnen worden vrijgegeven.
diff --git a/frontend/src/components/detection/detectionProfiles.ts b/frontend/src/components/detection/detectionProfiles.ts
index b3d6b66e..075c863f 100644
--- a/frontend/src/components/detection/detectionProfiles.ts
+++ b/frontend/src/components/detection/detectionProfiles.ts
@@ -11,7 +11,6 @@ export interface DetectionOperatorProfile {
positiveSampleCount: number
maxBackgroundDetections: number
validationScope: string
- nationallyValidated: boolean
description: string
limitationMessage: string
}
@@ -30,7 +29,6 @@ export const DETECTION_OPERATOR_PROFILES: DetectionOperatorProfile[] = [
positiveSampleCount: 7,
maxBackgroundDetections: 0,
validationScope: '7 onafhankelijke testgebieden in Mol en de Kempen',
- nationallyValidated: false,
description:
'Aanbevolen profiel met een evenwicht tussen gevonden en gemiste kleine gebouwen, opnieuw gemeten over zeven onafhankelijke testgebieden in Mol en de Kempen.',
limitationMessage:
@@ -49,7 +47,6 @@ export const DETECTION_OPERATOR_PROFILES: DetectionOperatorProfile[] = [
positiveSampleCount: 7,
maxBackgroundDetections: 0,
validationScope: '7 onafhankelijke testgebieden in Mol en de Kempen',
- nationallyValidated: false,
description: 'Voorgaand profiel voor controles waarbij minder foutieve vondsten belangrijker zijn dan maximale dekking.',
limitationMessage:
'De lege-achtergrondtest is geslaagd. Dit profiel vindt minder onterechte objecten, maar mist meer kleine gebouwen dan het aanbevolen profiel.',
@@ -67,7 +64,6 @@ export const DETECTION_OPERATOR_PROFILES: DetectionOperatorProfile[] = [
positiveSampleCount: 7,
maxBackgroundDetections: 0,
validationScope: '7 onafhankelijke testgebieden in Mol en de Kempen',
- nationallyValidated: false,
description: 'Profiel met hoge precisie voor controles waarbij zo weinig mogelijk foutieve vondsten zwaarder wegen dan volledige dekking.',
limitationMessage:
'Goedgekeurd na de lege-achtergrondtest. Resultaten in dun bebouwde context blijven altijd controlebewijs en geen automatische waarheid.',
diff --git a/frontend/src/components/map/MapWorkspace.tsx b/frontend/src/components/map/MapWorkspace.tsx
index d10adac5..a6b956b1 100644
--- a/frontend/src/components/map/MapWorkspace.tsx
+++ b/frontend/src/components/map/MapWorkspace.tsx
@@ -108,7 +108,7 @@ function productSupportsSelection(product: OnDemandMapProduct, bbox: VectorSelec
const scale = selectionAnalysisScale(bbox)
if (scale === 'overview') return false
const dimensions = selectionDimensions(bbox)
- if (product.kind === 'dhmv' || product.kind === 'flood_hazard') {
+ if (product.kind === 'dhmv' || product.kind === 'spw_terrain' || product.kind === 'flood_hazard') {
return dimensions.areaSquareMetres <= 280_000_000
}
if (product.kind === 'thematic_raster' || product.kind === 'walous') {
@@ -338,7 +338,7 @@ function datasetAvailabilityLabel(
): string {
const partitionCount = partitions.length
const regionalSuffix = partitionCount > 1 ? ` · ${partitionCount} gemeenten` : ''
- if (dataset.dataset_type === 'raster' && dataset.source_name === 'digitaal_vlaanderen_dhmv') {
+ if (dataset.dataset_type === 'raster' && ['digitaal_vlaanderen_dhmv', 'spw_terrain'].includes(dataset.source_name ?? '')) {
const resolution = Number(dataset.source_metadata?.['analysis_resolution_m'])
return `${Number.isFinite(resolution) ? `${resolution.toLocaleString('nl-BE')} m` : 'Raster'} hoogtegrid${regionalSuffix}`
}
@@ -416,7 +416,7 @@ function datasetMatchesTheme(dataset: DatasetCreateResponse, theme: DataTheme):
if (dataset.source_name === 'spw_bathymetry') {
return theme.id === 'bathymetry'
}
- if (dataset.source_name === 'digitaal_vlaanderen_dhmv') {
+ if (['digitaal_vlaanderen_dhmv', 'spw_terrain'].includes(dataset.source_name ?? '')) {
return theme.id === 'elevation'
}
if (dataset.source_name === 'department_omgeving_thematic_raster') {
@@ -444,7 +444,7 @@ function datasetMatchesTheme(dataset: DatasetCreateResponse, theme: DataTheme):
function isPartitionedRaster(dataset: DatasetCreateResponse | null | undefined): boolean {
return Boolean(
dataset?.dataset_type === 'raster'
- && ['digitaal_vlaanderen_dhmv', 'vmm_flood_hazard'].includes(dataset.source_name ?? ''),
+ && ['digitaal_vlaanderen_dhmv', 'spw_terrain', 'vmm_flood_hazard'].includes(dataset.source_name ?? ''),
)
}
@@ -544,6 +544,7 @@ function pickThemeDataset(
(dataset.source_name === 'spw_walous_land_cover' ? 5_000_000 : 0) +
(dataset.source_name === 'digitaal_vlaanderen_buildings_addresses_register' ? 5_000_000 : 0) +
(dataset.source_name === 'digitaal_vlaanderen_dhmv' ? 5_000_000 : 0) +
+ (dataset.source_name === 'spw_terrain' ? 5_000_000 : 0) +
(dataset.source_name === 'vmm_flood_hazard' ? 5_000_000 : 0) +
(dataset.source_name === 'vmm_vha_bathymetry_profiles' ? 5_000_000 : 0) +
(dataset.source_name === 'spw_bathymetry' ? 5_100_000 : 0) +
@@ -873,6 +874,7 @@ export function MapWorkspace({
selectedCoverageZones?.includes('flanders')
|| (!selectedCoverageZones && activeScopeProject?.name === FLANDERS_WORKSPACE_PROJECT_NAME),
)
+ const walloniaScopeSelected = Boolean(selectedCoverageZones?.includes('wallonia'))
const {
themeInsights,
themeInsightsLoading: themeResultsLoading,
@@ -992,6 +994,14 @@ export function MapWorkspace({
&& datasetCoversSelectedArea(dataset, selectedMapAreaId, selectedMapArea?.name, regionalScopeSelected),
) ?? null
}
+ if (walloniaScopeSelected && officialMapProducts.spwTerrain.some((product) => product.configured)) {
+ result.elevation = availableMapDatasets.find(
+ (dataset) =>
+ dataset.source_name === 'spw_terrain'
+ && datasetProductKey(dataset) === 'spw_mnt_1m_2021_2022'
+ && datasetCoversSelectedArea(dataset, selectedMapAreaId, selectedMapArea?.name, regionalScopeSelected),
+ ) ?? null
+ }
if (flandersScopeSelected && officialMapProducts.thematic.length > 0) {
for (const product of officialMapProducts.thematic) {
if (!result[product.theme]) {
@@ -1021,6 +1031,7 @@ export function MapWorkspace({
flandersScopeSelected,
floodHazardDatasets,
officialMapProducts.dhmv.length,
+ officialMapProducts.spwTerrain,
officialMapProducts.floodHazard.length,
officialMapProducts.grb,
officialMapProducts.officialVector,
@@ -1032,6 +1043,7 @@ export function MapWorkspace({
selectedMapArea?.name,
selectedMapAreaId,
selectedCoverageZones,
+ walloniaScopeSelected,
])
const themePartitionMap = useMemo(
() =>
@@ -1066,6 +1078,7 @@ export function MapWorkspace({
effectiveZones?.includes('flanders')
|| (!effectiveZones && activeScopeProject?.name === FLANDERS_WORKSPACE_PROJECT_NAME),
)
+ const includesWallonia = Boolean(effectiveZones?.includes('wallonia'))
if (includesFlanders) {
for (const product of officialMapProducts.thematic) {
result.push({
@@ -1161,6 +1174,21 @@ export function MapWorkspace({
coverageZones: ['flanders'],
})
}
+ const spwTerrainProduct = includesWallonia
+ ? officialMapProducts.spwTerrain.find((product) => product.configured)
+ : null
+ if (spwTerrainProduct) {
+ result.push({
+ kind: 'spw_terrain',
+ productKey: spwTerrainProduct.key,
+ displayName: spwTerrainProduct.display_name,
+ theme: 'elevation',
+ availabilityLabel: `${spwTerrainProduct.analysis_resolution_m} m analyse · ${spwTerrainProduct.acquisition_period} · automatisch bij selectie`,
+ attribution: spwTerrainProduct.attribution,
+ limitationMessage: spwTerrainProduct.limitation_message,
+ coverageZones: spwTerrainProduct.coverage_zones,
+ })
+ }
const floodProduct = includesFlanders
? officialMapProducts.floodHazard.find(
(product) => product.key === selectedFloodHazardProductKey,
@@ -1319,7 +1347,7 @@ export function MapWorkspace({
activeTheme.id === 'elevation' && selectedProjectId
? activeThemePartitions.flatMap((dataset) => {
const bounds = dataset.source_metadata?.['bbox_epsg4326']
- return dataset.source_name === 'digitaal_vlaanderen_dhmv'
+ return ['digitaal_vlaanderen_dhmv', 'spw_terrain'].includes(dataset.source_name ?? '')
&& Array.isArray(bounds)
&& bounds.length === 4
? [{
@@ -1530,7 +1558,7 @@ export function MapWorkspace({
const terrainReliefMetric = activeSupportingMetrics.find((metric) => metric.metric_key === 'relief_m')
const populationDensityMetric = activeSupportingMetrics.find((metric) => metric.metric_key === 'population_density_mean_per_ha')
const scoreMedianMetric = activeSupportingMetrics.find((metric) => metric.metric_key.endsWith('_median'))
- const activeSecondaryMetric = activeMetricUnit === 'm TAW'
+ const activeSecondaryMetric = activeMetricUnit === 'm TAW' || activeMetricUnit === 'm DNG'
? terrainReliefMetric ? `${terrainReliefMetric.metric_value.toLocaleString('nl-BE', { maximumFractionDigits: 2 })} m reliëf` : null
: activeMetricUnit === 'inwoners'
? populationDensityMetric ? selectionMetricLabel(populationDensityMetric) : null
@@ -1543,7 +1571,7 @@ export function MapWorkspace({
: null
const activeSecondaryLabel = activeMetricUnit === 'ha'
? 'Aandeel selectie'
- : activeMetricUnit === 'm TAW'
+ : activeMetricUnit === 'm TAW' || activeMetricUnit === 'm DNG'
? 'Reliëf'
: activeMetricUnit === 'inwoners'
? 'Gemiddelde dichtheid'
@@ -2301,7 +2329,7 @@ export function MapWorkspace({
{officialMapProductsError}
) : null}
- {analysisMode === 'current' && activeTheme.id === 'elevation' && officialMapProducts.dhmv.length > 0 ? (
+ {analysisMode === 'current' && activeTheme.id === 'elevation' && flandersScopeSelected && officialMapProducts.dhmv.length > 0 ? (
Hoogtemodel
DTM meet het maaiveld; DSM bevat ook gebouwen en vegetatie.
) : null}
+ {analysisMode === 'current' && activeTheme.id === 'elevation' && walloniaScopeSelected && officialMapProducts.spwTerrain.some((product) => product.configured) ? (
+ SPW MNT 2021-2022 · 1 m bron · 5 m begrensde analyse · hoogte in m DNG.
+ ) : null}
{analysisMode === 'current' && activeTheme.id === 'flood_hazard' && officialMapProducts.floodHazard.length > 0 ? (
diff --git a/frontend/src/components/map/mapWorkspaceUtils.ts b/frontend/src/components/map/mapWorkspaceUtils.ts
index 6f56af55..5daf91a4 100644
--- a/frontend/src/components/map/mapWorkspaceUtils.ts
+++ b/frontend/src/components/map/mapWorkspaceUtils.ts
@@ -108,7 +108,7 @@ export function persistedDatasetSupportsSelection(
const scale = selectionAnalysisScale(bbox)
if (scale === 'overview') return false
const dimensions = selectionDimensions(bbox)
- if (dataset.source_name === 'digitaal_vlaanderen_dhmv' || dataset.source_name === 'vmm_flood_hazard') {
+ if (dataset.source_name === 'digitaal_vlaanderen_dhmv' || dataset.source_name === 'spw_terrain' || dataset.source_name === 'vmm_flood_hazard') {
return dimensions.areaSquareMetres <= 280_000_000
}
if (dataset.source_name === 'department_omgeving_thematic_raster') {
diff --git a/frontend/src/hooks/useMapSelectionExtract.ts b/frontend/src/hooks/useMapSelectionExtract.ts
index 41f3aec4..9f9d787b 100644
--- a/frontend/src/hooks/useMapSelectionExtract.ts
+++ b/frontend/src/hooks/useMapSelectionExtract.ts
@@ -42,7 +42,7 @@ export function useMapSelectionExtract({
setMapSelectionError('Open a vector dataset before extracting a map area.')
return null
}
- const terrainDataset = selectedDataset.dataset_type === 'raster' && selectedDataset.source_name === 'digitaal_vlaanderen_dhmv'
+ const terrainDataset = selectedDataset.dataset_type === 'raster' && ['digitaal_vlaanderen_dhmv', 'spw_terrain'].includes(selectedDataset.source_name ?? '')
const floodHazardDataset = selectedDataset.dataset_type === 'raster' && selectedDataset.source_name === 'vmm_flood_hazard'
const thematicRasterDataset = selectedDataset.dataset_type === 'raster' && selectedDataset.source_name === 'department_omgeving_thematic_raster'
const bathymetryRasterDataset = selectedDataset.dataset_type === 'raster' && selectedDataset.source_name === 'spw_bathymetry'
diff --git a/frontend/src/hooks/useMapThemeSelectionInsights.ts b/frontend/src/hooks/useMapThemeSelectionInsights.ts
index fa86960f..f79bf2e4 100644
--- a/frontend/src/hooks/useMapThemeSelectionInsights.ts
+++ b/frontend/src/hooks/useMapThemeSelectionInsights.ts
@@ -11,6 +11,7 @@ export type MapThemeAcquisitionKind =
| 'thematic_raster'
| 'walous'
| 'dhmv'
+ | 'spw_terrain'
| 'flood_hazard'
| 'grb'
| 'official_vector'
@@ -136,6 +137,11 @@ export function useMapThemeSelectionInsights(
...commonPayload,
product_key: acquisition.productKey as 'dtm_1m' | 'dsm_1m',
})
+ : acquisition.kind === 'spw_terrain'
+ ? await datasetsApi.acquireSpwTerrain(selectedProjectId, {
+ ...commonPayload,
+ product_key: 'spw_mnt_1m_2021_2022',
+ })
: acquisition.kind === 'flood_hazard'
? await datasetsApi.acquireFloodHazard(selectedProjectId, {
...commonPayload,
@@ -198,9 +204,9 @@ export function useMapThemeSelectionInsights(
const acquiredDatasetIds = acquiredDatasets.map((item) => item.id)
const selectedDatasetIds = acquiredDatasetIds.length > 0 ? acquiredDatasetIds : datasetIds ?? []
const acquiredAsPartitions = acquiredDatasetIds.length > 1
- const result = dataset.dataset_type === 'raster' && dataset.source_name === 'digitaal_vlaanderen_dhmv'
+ const result = dataset.dataset_type === 'raster' && ['digitaal_vlaanderen_dhmv', 'spw_terrain'].includes(dataset.source_name ?? '')
? terrainSelectionToMapSelection(
- partitioned || acquiredAsPartitions
+ dataset.source_name === 'digitaal_vlaanderen_dhmv' && (partitioned || acquiredAsPartitions)
? await datasetsApi.selectTerrainPartitions(selectedProjectId, {
bbox,
area_id: areaId,
diff --git a/frontend/src/hooks/useOfficialMapProducts.ts b/frontend/src/hooks/useOfficialMapProducts.ts
index 2f987f47..833dd289 100644
--- a/frontend/src/hooks/useOfficialMapProducts.ts
+++ b/frontend/src/hooks/useOfficialMapProducts.ts
@@ -5,6 +5,7 @@ import type {
BathymetrySourceRead,
CoverageResolveResponse,
DhmvProductRead,
+ SpwTerrainProductRead,
FloodHazardProductRead,
GrbProductRead,
OfficialVectorProductRead,
@@ -16,6 +17,7 @@ export interface OfficialMapProducts {
thematic: ThematicRasterProductRead[]
walous: ThematicRasterProductRead[]
dhmv: DhmvProductRead[]
+ spwTerrain: SpwTerrainProductRead[]
floodHazard: FloodHazardProductRead[]
grb: GrbProductRead[]
officialVector: OfficialVectorProductRead[]
@@ -26,6 +28,7 @@ const EMPTY_PRODUCTS: OfficialMapProducts = {
thematic: [],
walous: [],
dhmv: [],
+ spwTerrain: [],
floodHazard: [],
grb: [],
officialVector: [],
@@ -54,17 +57,19 @@ export function useOfficialMapProducts(selectedProjectId: string | null) {
datasetsApi.listThematicRasterProducts(selectedProjectId),
datasetsApi.listWalousProducts(selectedProjectId),
datasetsApi.listDhmvProducts(selectedProjectId),
+ datasetsApi.listSpwTerrainProducts(selectedProjectId),
datasetsApi.listFloodHazardProducts(selectedProjectId),
datasetsApi.listGrbProducts(selectedProjectId),
datasetsApi.listOfficialVectorProducts(selectedProjectId),
datasetsApi.listBathymetrySources(selectedProjectId),
])
- .then(([thematic, walous, dhmv, floodHazard, grb, officialVector, bathymetry]) => {
+ .then(([thematic, walous, dhmv, spwTerrain, floodHazard, grb, officialVector, bathymetry]) => {
if (!cancelled) {
setProducts({
thematic: thematic.items,
walous: walous.items,
dhmv: dhmv.items,
+ spwTerrain: spwTerrain.items,
floodHazard: floodHazard.items,
grb: grb.items,
officialVector: officialVector.items,
@@ -75,7 +80,7 @@ export function useOfficialMapProducts(selectedProjectId: string | null) {
.catch((requestError) => {
if (!cancelled) {
setProducts(EMPTY_PRODUCTS)
- setError(formatError(requestError, 'De officiële Vlaamse kaartcatalogi konden niet worden geladen.'))
+ setError(formatError(requestError, 'De officiële regionale kaartcatalogi konden niet worden geladen.'))
}
})
.finally(() => {
diff --git a/frontend/src/lib/datasetCapabilities.ts b/frontend/src/lib/datasetCapabilities.ts
index 8b4aed88..afdc0957 100644
--- a/frontend/src/lib/datasetCapabilities.ts
+++ b/frontend/src/lib/datasetCapabilities.ts
@@ -3,6 +3,7 @@ import type { DatasetCreateResponse } from '../types'
const NON_IMAGERY_RASTER_SOURCES = new Set([
'department_omgeving_thematic_raster',
'digitaal_vlaanderen_dhmv',
+ 'spw_terrain',
'vmm_flood_hazard',
'spw_bathymetry',
'spw_walous_land_cover',
diff --git a/frontend/src/lib/datasetDisplay.ts b/frontend/src/lib/datasetDisplay.ts
index 82a02332..e14dc2f3 100644
--- a/frontend/src/lib/datasetDisplay.ts
+++ b/frontend/src/lib/datasetDisplay.ts
@@ -29,6 +29,7 @@ const DATASET_SOURCE_LABELS: Record = {
digitaal_vlaanderen_orthophoto: 'Digitaal Vlaanderen',
digitaal_vlaanderen_buildings_addresses_register: 'Digitaal Vlaanderen',
digitaal_vlaanderen_dhmv: 'Digitaal Vlaanderen',
+ spw_terrain: 'Service public de Wallonie',
grb: 'GRB',
historical_landuse: 'Digitaal Vlaanderen',
inbo_bwk_natura2000: 'INBO',
@@ -54,6 +55,10 @@ export function getDatasetDisplayName(dataset: DatasetCreateResponse): string {
const productName = dataset.source_metadata?.['product_display_name']
return typeof productName === 'string' && productName.trim() ? productName : 'DHMV II hoogtemodel'
}
+ if (dataset.source_name === 'spw_terrain') {
+ const productName = dataset.source_metadata?.['product_display_name']
+ return typeof productName === 'string' && productName.trim() ? productName : 'SPW terreinmodel 2021-2022'
+ }
if (dataset.source_name === 'vmm_flood_hazard') {
const productName = dataset.source_metadata?.['product_display_name']
return typeof productName === 'string' && productName.trim() ? productName : 'VMM-overstromingsscenario'
diff --git a/frontend/src/lib/sourcePortfolio.ts b/frontend/src/lib/sourcePortfolio.ts
index 38dd92db..7056afd7 100644
--- a/frontend/src/lib/sourcePortfolio.ts
+++ b/frontend/src/lib/sourcePortfolio.ts
@@ -199,10 +199,10 @@ export const OFFICIAL_SOURCE_PORTFOLIO: OfficialSourceDefinition[] = [
owner: 'Digitaal Vlaanderen',
coverage: 'DTM/DSM, opname 2013-2015',
value: 'Hoogte, reliëf en helling als gemeten terreinbasis.',
- metricExamples: 'hoogte m TAW, reliëf en helling in graden',
+ metricExamples: 'hoogte in de regionale verticale referentie, reliëf en helling in graden',
priority: 'next',
url: 'https://www.vlaanderen.be/digitaal-vlaanderen/onze-diensten-en-platformen/earth-observation-data-science-eodas/het-digitaal-hoogtemodel/digitaal-hoogtemodel-vlaanderen-ii',
- matches: (dataset) => sourceNameIs(dataset, 'digitaal_vlaanderen_dhmv'),
+ matches: (dataset) => sourceNameIs(dataset, 'digitaal_vlaanderen_dhmv') || sourceNameIs(dataset, 'spw_terrain'),
},
{
key: 'soil_map',
diff --git a/frontend/src/services/api/datasets.ts b/frontend/src/services/api/datasets.ts
index e3bc3c3d..ad02cb96 100644
--- a/frontend/src/services/api/datasets.ts
+++ b/frontend/src/services/api/datasets.ts
@@ -27,6 +27,8 @@ import type {
BathymetrySourceProbeRead,
BathymetrySourceRead,
DhmvProductRead,
+ SpwTerrainAcquireRequest,
+ SpwTerrainProductRead,
GrbAcquireRequest,
GrbProductRead,
OfficialVectorAcquireRequest,
@@ -151,6 +153,10 @@ export const datasetsApi = {
apiPost(`/api/v1/projects/${projectId}/datasets/dhmv/acquire`, payload),
listDhmvProducts: (projectId: string): Promise<{ items: DhmvProductRead[]; total: number }> =>
apiGet<{ items: DhmvProductRead[]; total: number }>(`/api/v1/projects/${projectId}/datasets/dhmv/products`),
+ acquireSpwTerrain: (projectId: string, payload: SpwTerrainAcquireRequest): Promise =>
+ apiPost(`/api/v1/projects/${projectId}/datasets/spw-terrain/acquire`, payload),
+ listSpwTerrainProducts: (projectId: string): Promise<{ items: SpwTerrainProductRead[]; total: number }> =>
+ apiGet<{ items: SpwTerrainProductRead[]; total: number }>(`/api/v1/projects/${projectId}/datasets/spw-terrain/products`),
acquireGrb: (projectId: string, payload: GrbAcquireRequest): Promise =>
apiPost(`/api/v1/projects/${projectId}/datasets/grb/acquire`, payload),
listGrbProducts: (projectId: string): Promise<{ items: GrbProductRead[]; total: number }> =>
diff --git a/frontend/src/types.ts b/frontend/src/types.ts
index 98d4b229..062622d3 100644
--- a/frontend/src/types.ts
+++ b/frontend/src/types.ts
@@ -362,6 +362,33 @@ export interface DhmvProductRead {
limitation_message: string
}
+export interface SpwTerrainAcquireRequest {
+ bbox: VectorSelectionBBox
+ area_id?: string | null
+ product_key?: 'spw_mnt_1m_2021_2022'
+ resolution_m?: number | null
+ force_refresh?: boolean
+}
+
+export interface SpwTerrainProductRead {
+ key: 'spw_mnt_1m_2021_2022'
+ display_name: string
+ surface_model: 'terrain'
+ source_filename: string
+ native_resolution_m: number
+ analysis_resolution_m: number
+ source_crs: 'EPSG:3812'
+ vertical_reference: string
+ acquisition_period: string
+ catalog_url: string
+ attribution: string
+ license_note: string
+ limitation_message: string
+ coverage_zones: string[]
+ configured: boolean
+ status: string
+}
+
export interface GrbAcquireRequest {
bbox: VectorSelectionBBox
area_id?: string | null
@@ -1156,6 +1183,11 @@ export interface DetectionModelCapability {
status: string
limitation_message: string
version?: string | null
+ training_scope?: string | null
+ validation_scope?: string | null
+ validated_regions: string[]
+ nationally_validated: boolean
+ operator_review_required: boolean
}
export interface DetectionModelsResponse {
diff --git a/scripts/README.md b/scripts/README.md
index ea7a73e6..06640967 100644
--- a/scripts/README.md
+++ b/scripts/README.md
@@ -4,22 +4,39 @@ Setup-, import-, demo- en maintenance-scripts voor GeoIntel.
## WALOUS source provisioning
-Run the networked operator only after checking at least 2 GB of archive space
+Run the networked operator only after checking at least 3 GB of archive space
plus room for the extracted official GeoTIFFs:
```bash
python scripts/provision_walous_sources.py \
- --years 2020 2023 \
+ --years 2018 2020 2023 \
--destination storage/source-cache/walous
```
-The command accepts only the hard-coded official SPW 2020/2023 archives,
-streams with a 1 GB per-archive cap, rejects changed content lengths, extracts
+The command accepts only the hard-coded official SPW 2018/2020/2023 archives,
+streams with a 1.25 GB per-archive cap, rejects changed content lengths, extracts
only the single GeoTIFF by basename, validates the raster contract and writes
checksums plus `provisioning-report.json`. Existing valid sources are reused;
`--force` performs a new download. This is an operator acquisition, not an
application startup task.
+## SPW Wallonia terrain source provisioning
+
+The official 1 m MNT is a large operator asset, never an implicit startup
+download. Reserve at least 90 GB temporarily for archive plus extraction and
+run:
+
+```bash
+python scripts/provision_spw_terrain_source.py \
+ --destination storage/source-cache/spw-terrain
+```
+
+The provisioner accepts only the fixed official SPW artifact, enforces a
+30-60 GB archive range and a single safe GeoTIFF member, validates EPSG:3812,
+one band, native 1 m cells and representative elevation samples, writes
+source/archive SHA-256 evidence and removes the archive after successful
+extraction unless `--keep-archive` is supplied.
+
## Runtime verification
Inspect interrupted runtime state without changing it:
diff --git a/scripts/provision_spw_terrain_source.py b/scripts/provision_spw_terrain_source.py
new file mode 100644
index 00000000..37f3b421
--- /dev/null
+++ b/scripts/provision_spw_terrain_source.py
@@ -0,0 +1,242 @@
+#!/usr/bin/env python3
+"""Provision the official Wallonia 1 m MNT for bounded GeoIntel analysis."""
+
+from __future__ import annotations
+
+import argparse
+import hashlib
+import json
+from pathlib import Path
+import shutil
+import sys
+from urllib.request import Request, urlopen
+from zipfile import BadZipFile, ZipFile
+
+
+SOURCE_URL = (
+ "https://geoservices.wallonie.be/geotraitement/spwdatadownload/results/"
+ "fe13bc84-e371-46ca-9632-8ad4139f1ee5/RELIEF_WALLONIE_MNT_1M_2021_2022_GEOTIFF_3812.zip"
+)
+TARGET_FILENAME = "spw_mnt_1m_2021_2022_3812.tif"
+ARCHIVE_FILENAME = "spw_mnt_1m_2021_2022_3812.zip"
+MIN_ARCHIVE_BYTES = 30_000_000_000
+MAX_ARCHIVE_BYTES = 60_000_000_000
+MAX_EXTRACTED_BYTES = 80_000_000_000
+
+
+def sha256_file(path: Path) -> str:
+ digest = hashlib.sha256()
+ with path.open("rb") as handle:
+ while chunk := handle.read(16 * 1024 * 1024):
+ digest.update(chunk)
+ return digest.hexdigest()
+
+
+def download(destination: Path) -> str:
+ temporary = destination.with_suffix(destination.suffix + ".part")
+ temporary.unlink(missing_ok=True)
+ digest = hashlib.sha256()
+ received = 0
+ request = Request(
+ SOURCE_URL,
+ headers={"User-Agent": "GeoIntel/1.0 SPW-terrain-source-provisioner"},
+ )
+ try:
+ with urlopen(request, timeout=300) as response, temporary.open("wb") as output:
+ content_length = int(response.headers.get("Content-Length") or 0)
+ if (
+ content_length
+ and not MIN_ARCHIVE_BYTES <= content_length <= MAX_ARCHIVE_BYTES
+ ):
+ raise RuntimeError(
+ f"official archive size is outside the governed range: {content_length}"
+ )
+ while chunk := response.read(16 * 1024 * 1024):
+ received += len(chunk)
+ if received > MAX_ARCHIVE_BYTES:
+ raise RuntimeError(
+ "official archive exceeds the governed 60 GB transfer limit"
+ )
+ digest.update(chunk)
+ output.write(chunk)
+ if received % (1024 * 1024 * 1024) < len(chunk):
+ print(
+ f" downloaded {received / 1024 / 1024 / 1024:.1f} GiB",
+ flush=True,
+ )
+ if received < MIN_ARCHIVE_BYTES:
+ raise RuntimeError(
+ f"official archive is unexpectedly small: {received} bytes"
+ )
+ temporary.replace(destination)
+ return digest.hexdigest()
+ except Exception:
+ temporary.unlink(missing_ok=True)
+ raise
+
+
+def extract_single_geotiff(archive: Path, target: Path) -> None:
+ try:
+ with ZipFile(archive) as bundle:
+ candidates = [
+ item
+ for item in bundle.infolist()
+ if not item.is_dir()
+ and item.filename.lower().endswith((".tif", ".tiff"))
+ ]
+ if len(candidates) != 1:
+ raise RuntimeError(
+ f"archive must contain exactly one GeoTIFF, found {len(candidates)}"
+ )
+ member = candidates[0]
+ if member.file_size <= 0 or member.file_size > MAX_EXTRACTED_BYTES:
+ raise RuntimeError(
+ f"GeoTIFF uncompressed size is outside the governed limit: {member.file_size}"
+ )
+ temporary = target.with_suffix(target.suffix + ".part")
+ temporary.unlink(missing_ok=True)
+ with bundle.open(member) as source, temporary.open("wb") as output:
+ shutil.copyfileobj(source, output, length=16 * 1024 * 1024)
+ temporary.replace(target)
+ except BadZipFile as exc:
+ raise RuntimeError("official SPW MNT archive is not a valid ZIP file") from exc
+
+
+def validate_raster(path: Path) -> dict:
+ try:
+ import numpy as np
+ import rasterio
+ from rasterio.windows import Window
+ except ImportError as exc:
+ raise RuntimeError(
+ "rasterio and numpy are required to validate the SPW MNT source"
+ ) from exc
+ with rasterio.open(path) as source:
+ if source.crs is None or source.crs.to_epsg() != 3812:
+ raise RuntimeError(f"SPW MNT must use EPSG:3812, found {source.crs}")
+ if source.count != 1:
+ raise RuntimeError(f"SPW MNT must have one band, found {source.count}")
+ if not all(abs(abs(float(value)) - 1.0) <= 0.05 for value in source.res):
+ raise RuntimeError(f"SPW MNT must retain 1 m cells, found {source.res}")
+ sample_windows = []
+ sample_size = 512
+ for x_fraction, y_fraction in (
+ (0.1, 0.1),
+ (0.5, 0.5),
+ (0.9, 0.9),
+ (0.1, 0.9),
+ (0.9, 0.1),
+ ):
+ col = max(
+ 0,
+ min(
+ source.width - sample_size,
+ round(source.width * x_fraction - sample_size / 2),
+ ),
+ )
+ row = max(
+ 0,
+ min(
+ source.height - sample_size,
+ round(source.height * y_fraction - sample_size / 2),
+ ),
+ )
+ sample_windows.append(
+ Window(
+ col,
+ row,
+ min(sample_size, source.width),
+ min(sample_size, source.height),
+ )
+ )
+ samples = [
+ source.read(1, window=window, masked=True).compressed().astype("float64")
+ for window in sample_windows
+ ]
+ values = np.concatenate([sample for sample in samples if sample.size])
+ values = values[np.isfinite(values)]
+ if not values.size:
+ raise RuntimeError(
+ "SPW MNT validation samples contain no finite elevation values"
+ )
+ if float(values.min()) < -100.0 or float(values.max()) > 1000.0:
+ raise RuntimeError(
+ f"SPW MNT samples contain implausible values: {values.min()}..{values.max()}"
+ )
+ return {
+ "path": str(path),
+ "crs": str(source.crs),
+ "width": int(source.width),
+ "height": int(source.height),
+ "resolution": [float(value) for value in source.res],
+ "bounds": [float(value) for value in source.bounds],
+ "nodata": None if source.nodata is None else float(source.nodata),
+ "dtype": source.dtypes[0],
+ "sample_min_m": float(values.min()),
+ "sample_max_m": float(values.max()),
+ }
+
+
+def provision(destination: Path, force: bool, keep_archive: bool) -> dict:
+ destination.mkdir(parents=True, exist_ok=True)
+ target = destination / TARGET_FILENAME
+ archive = destination / ARCHIVE_FILENAME
+ archive_digest = None
+ if target.is_file() and not force:
+ print(f"SPW MNT: validating existing source {target}")
+ else:
+ if archive.is_file():
+ archive_size = archive.stat().st_size
+ if not MIN_ARCHIVE_BYTES <= archive_size <= MAX_ARCHIVE_BYTES:
+ raise RuntimeError(
+ f"existing archive size is outside the governed range: {archive_size}"
+ )
+ print(f"SPW MNT: using existing archive {archive}")
+ archive_digest = sha256_file(archive)
+ else:
+ print("SPW MNT: downloading official archive")
+ archive_digest = download(archive)
+ print(f"SPW MNT: archive sha256 {archive_digest}")
+ extract_single_geotiff(archive, target)
+ if not keep_archive:
+ archive.unlink(missing_ok=True)
+ validation = validate_raster(target)
+ source_digest = sha256_file(target)
+ target.with_suffix(".sha256").write_text(
+ f"{source_digest} {target.name}\n", encoding="ascii"
+ )
+ validation.update(
+ {
+ "source_sha256": source_digest,
+ "archive_sha256": archive_digest,
+ "download_url": SOURCE_URL,
+ "catalog_url": "https://geoportail.wallonie.be/catalogue/fe13bc84-e371-46ca-9632-8ad4139f1ee5.html",
+ }
+ )
+ report_path = destination / "provisioning-report.json"
+ report_path.write_text(json.dumps(validation, indent=2) + "\n", encoding="utf-8")
+ print(f"SPW MNT: ready ({target.stat().st_size / 1024 / 1024 / 1024:.1f} GiB)")
+ print(f"Provisioning report: {report_path}")
+ return validation
+
+
+def main() -> int:
+ parser = argparse.ArgumentParser(
+ description="Provision the official Wallonia 1 m MNT 2021-2022 GeoTIFF."
+ )
+ parser.add_argument(
+ "--destination", type=Path, default=Path("storage/source-cache/spw-terrain")
+ )
+ parser.add_argument("--force", action="store_true")
+ parser.add_argument("--keep-archive", action="store_true")
+ args = parser.parse_args()
+ provision(args.destination.resolve(), args.force, args.keep_archive)
+ return 0
+
+
+if __name__ == "__main__":
+ try:
+ raise SystemExit(main())
+ except Exception as exc:
+ print(f"SPW_TERRAIN_PROVISIONING_FAILED: {exc}", file=sys.stderr)
+ raise SystemExit(1) from exc
diff --git a/scripts/provision_walous_sources.py b/scripts/provision_walous_sources.py
index c9be618a..0235f4f3 100644
--- a/scripts/provision_walous_sources.py
+++ b/scripts/provision_walous_sources.py
@@ -14,6 +14,14 @@ from zipfile import BadZipFile, ZipFile
SOURCES = {
+ 2018: {
+ "url": (
+ "https://geoservices.wallonie.be/geotraitement/spwdatadownload/results/"
+ "a0ad23a1-1845-4bd5-8c2f-0f62d3f1ec75/WALOUS_OCS__2018_GEOTIFF_3812.zip"
+ ),
+ "expected_archive_bytes": 1_122_785_133,
+ "target": "walous_land_cover_2018_3812.tif",
+ },
2020: {
"url": (
"https://geoservices.wallonie.be/geotraitement/spwdatadownload/results/"
@@ -31,9 +39,46 @@ SOURCES = {
"target": "walous_land_cover_2023_3812.tif",
},
}
-MAX_ARCHIVE_BYTES = 1_000_000_000
+MAX_ARCHIVE_BYTES = 1_250_000_000
MAX_EXTRACTED_BYTES = 50_000_000_000
-WALOUS_CLASS_CODES = {1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90}
+WALOUS_CLASS_CODES = {
+ 0,
+ 1,
+ 2,
+ 3,
+ 4,
+ 5,
+ 6,
+ 7,
+ 8,
+ 9,
+ 11,
+ 15,
+ 18,
+ 19,
+ 28,
+ 29,
+ 31,
+ 38,
+ 39,
+ 51,
+ 55,
+ 58,
+ 59,
+ 62,
+ 71,
+ 73,
+ 75,
+ 80,
+ 81,
+ 83,
+ 85,
+ 90,
+ 91,
+ 93,
+ 95,
+}
+WALOUS_CANONICAL_CLASS_CODES = {1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 90}
def sha256_file(path: Path) -> str:
@@ -49,22 +94,30 @@ def download(url: str, destination: Path, expected_bytes: int) -> str:
temporary.unlink(missing_ok=True)
digest = hashlib.sha256()
received = 0
- request = Request(url, headers={"User-Agent": "GeoIntel/1.0 WALOUS-source-provisioner"})
+ request = Request(
+ url, headers={"User-Agent": "GeoIntel/1.0 WALOUS-source-provisioner"}
+ )
try:
with urlopen(request, timeout=300) as response, temporary.open("wb") as output:
content_length = int(response.headers.get("Content-Length") or 0)
if content_length and content_length != expected_bytes:
- raise RuntimeError(f"official archive size changed: expected {expected_bytes}, advertised {content_length}")
+ raise RuntimeError(
+ f"official archive size changed: expected {expected_bytes}, advertised {content_length}"
+ )
while chunk := response.read(8 * 1024 * 1024):
received += len(chunk)
if received > MAX_ARCHIVE_BYTES:
- raise RuntimeError("official archive exceeds the governed 1 GB transfer limit")
+ raise RuntimeError(
+ "official archive exceeds the governed 1 GB transfer limit"
+ )
digest.update(chunk)
output.write(chunk)
if received % (128 * 1024 * 1024) < len(chunk):
print(f" downloaded {received / 1024 / 1024:.0f} MiB", flush=True)
if received != expected_bytes:
- raise RuntimeError(f"archive is incomplete: expected {expected_bytes} bytes, received {received}")
+ raise RuntimeError(
+ f"archive is incomplete: expected {expected_bytes} bytes, received {received}"
+ )
temporary.replace(destination)
return digest.hexdigest()
except Exception:
@@ -75,13 +128,25 @@ def download(url: str, destination: Path, expected_bytes: int) -> str:
def extract_single_geotiff(archive: Path, target: Path) -> None:
try:
with ZipFile(archive) as bundle:
- candidates = [item for item in bundle.infolist() if not item.is_dir() and item.filename.lower().endswith((".tif", ".tiff"))]
+ candidates = [
+ item
+ for item in bundle.infolist()
+ if not item.is_dir()
+ and item.filename.lower().endswith((".tif", ".tiff"))
+ ]
if len(candidates) != 1:
- raise RuntimeError(f"archive must contain exactly one GeoTIFF, found {len(candidates)}")
+ raise RuntimeError(
+ f"archive must contain exactly one GeoTIFF, found {len(candidates)}"
+ )
member = candidates[0]
if member.file_size <= 0 or member.file_size > MAX_EXTRACTED_BYTES:
- raise RuntimeError(f"GeoTIFF uncompressed size is outside the governed limit: {member.file_size}")
- if Path(member.filename).name != member.filename.replace("\\", "/").split("/")[-1]:
+ raise RuntimeError(
+ f"GeoTIFF uncompressed size is outside the governed limit: {member.file_size}"
+ )
+ if (
+ Path(member.filename).name
+ != member.filename.replace("\\", "/").split("/")[-1]
+ ):
# Nested paths are accepted only by basename; extraction never trusts archive paths.
pass
temporary = target.with_suffix(target.suffix + ".part")
@@ -99,21 +164,37 @@ def validate_raster(path: Path) -> dict:
import rasterio
from rasterio.enums import Resampling
except ImportError as exc:
- raise RuntimeError("rasterio and numpy are required to validate WALOUS sources") from exc
+ raise RuntimeError(
+ "rasterio and numpy are required to validate WALOUS sources"
+ ) from exc
with rasterio.open(path) as source:
if source.crs is None or source.crs.to_epsg() != 3812:
raise RuntimeError(f"WALOUS raster must use EPSG:3812, found {source.crs}")
if source.count != 1:
- raise RuntimeError(f"WALOUS raster must have one band, found {source.count}")
+ raise RuntimeError(
+ f"WALOUS raster must have one band, found {source.count}"
+ )
if not all(abs(abs(float(value)) - 1.0) <= 0.05 for value in source.res):
- raise RuntimeError(f"WALOUS raster must retain 1 m cells, found {source.res}")
+ raise RuntimeError(
+ f"WALOUS raster must retain 1 m cells, found {source.res}"
+ )
sample_height = min(2048, source.height)
sample_width = min(2048, source.width)
- sample = source.read(1, out_shape=(sample_height, sample_width), masked=True, resampling=Resampling.nearest)
+ sample = source.read(
+ 1,
+ out_shape=(sample_height, sample_width),
+ masked=True,
+ resampling=Resampling.nearest,
+ )
values = np.unique(sample.compressed()).astype(int).tolist()
- unexpected = sorted(set(values) - WALOUS_CLASS_CODES)
+ allowed_codes = (
+ WALOUS_CLASS_CODES if "2018" in path.name else WALOUS_CANONICAL_CLASS_CODES
+ )
+ unexpected = sorted(set(values) - allowed_codes)
if unexpected:
- raise RuntimeError(f"WALOUS sample contains classes outside the official 11-class code set: {unexpected}")
+ raise RuntimeError(
+ f"WALOUS sample contains classes outside the official 11-class code set: {unexpected}"
+ )
return {
"path": str(path),
"crs": str(source.crs),
@@ -123,6 +204,9 @@ def validate_raster(path: Path) -> dict:
"bounds": [float(value) for value in source.bounds],
"nodata": None if source.nodata is None else float(source.nodata),
"sample_classes": values,
+ "implicit_source_nodata_values": [0]
+ if "2018" in path.name and 0 in values
+ else [],
}
@@ -135,8 +219,17 @@ def provision(year: int, destination: Path, force: bool) -> dict:
validation = validate_raster(target)
digest = sha256_file(target)
else:
- print(f"WALOUS {year}: downloading official archive")
- archive_digest = download(source["url"], archive, source["expected_archive_bytes"])
+ if (
+ archive.is_file()
+ and archive.stat().st_size == source["expected_archive_bytes"]
+ ):
+ print(f"WALOUS {year}: using existing official archive {archive}")
+ archive_digest = sha256_file(archive)
+ else:
+ print(f"WALOUS {year}: downloading official archive")
+ archive_digest = download(
+ source["url"], archive, source["expected_archive_bytes"]
+ )
print(f"WALOUS {year}: archive sha256 {archive_digest}")
extract_single_geotiff(archive, target)
validation = validate_raster(target)
@@ -150,15 +243,40 @@ def provision(year: int, destination: Path, force: bool) -> dict:
def main() -> int:
- parser = argparse.ArgumentParser(description="Provision official WALOUS 2020/2023 GeoTIFF sources.")
- parser.add_argument("--years", nargs="+", type=int, choices=sorted(SOURCES), default=sorted(SOURCES))
- parser.add_argument("--destination", type=Path, default=Path("storage/source-cache/walous"))
+ parser = argparse.ArgumentParser(
+ description="Provision official WALOUS 2018/2020/2023 GeoTIFF sources."
+ )
+ parser.add_argument(
+ "--years", nargs="+", type=int, choices=sorted(SOURCES), default=sorted(SOURCES)
+ )
+ parser.add_argument(
+ "--destination", type=Path, default=Path("storage/source-cache/walous")
+ )
parser.add_argument("--force", action="store_true")
args = parser.parse_args()
args.destination.mkdir(parents=True, exist_ok=True)
- report = [provision(year, args.destination.resolve(), args.force) for year in args.years]
+ report = [
+ provision(year, args.destination.resolve(), args.force) for year in args.years
+ ]
report_path = args.destination / "provisioning-report.json"
- report_path.write_text(json.dumps({"sources": report}, indent=2) + "\n", encoding="utf-8")
+ retained: dict[int, dict] = {}
+ if report_path.is_file():
+ try:
+ retained = {
+ int(item["year"]): item
+ for item in json.loads(report_path.read_text(encoding="utf-8")).get(
+ "sources", []
+ )
+ if isinstance(item, dict) and item.get("year") in SOURCES
+ }
+ except (OSError, ValueError, TypeError):
+ retained = {}
+ retained.update({int(item["year"]): item for item in report})
+ report_path.write_text(
+ json.dumps({"sources": [retained[year] for year in sorted(retained)]}, indent=2)
+ + "\n",
+ encoding="utf-8",
+ )
print(f"Provisioning report: {report_path}")
return 0