Files
geointel/backend/tests/test_small_selection_raster_analysis.py
T
JensandClaude Opus 5 dd87a62e8f report what an area selection actually measured
Four ways a selection produced a confident number about a different area than
the operator drew:

Flood hazard divided the inundated cells by every cell in the drawn rectangle,
including cells the VMM raster does not model at all. A selection reaching
past the modelled extent therefore reported a diluted risk share, turning
missing data into an implied absence of risk. Terrain, bathymetry and thematic
raster already divided by valid cells; flood hazard was the outlier. It now
reports the three populations separately, states model coverage next to the
drawn area, and returns a null fraction rather than a zero when nothing was
modelled.

geometry_mask selects a cell when its centre falls inside the geometry, so a
rectangle smaller than one cell — or one landing between four centres —
selected nothing and the analysis returned zeros indistinguishable on screen
from "we looked and there is nothing here". On a 100 m population raster a
40 m rectangle over a city block reported no inhabitants. Selection now falls
back to the touched cells and says that it did, since the answer then covers
more ground than was requested. rasterio.mask applies the same centre rule
when cropping, so that call is widened too; the cells that count are still
decided by the centre rule wherever it selects anything.

The object count treated any feature touching the selection as whole, while
intersection_area clipped it — two headline numbers on one panel describing
different populations. The count stays whole-feature, which is what "objecten"
means to an operator, but now reports how many the edge cuts and is marked an
estimate when it does. The area_weighted_sum branch reuses that same count
instead of issuing its own near-identical query.

Partitioned selection de-duplicated the count on source_feature_id but
returned the raw rows, so a building on a municipal boundary was counted once
and drawn twice.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
2026-08-22 14:33:19 +02:00

118 lines
3.9 KiB
Python

"""A selection finer than the source raster must answer, not return zero.
End-to-end counterpart to ``test_raster_cell_selection``: the analysis reads a
real GeoTIFF, so it proves the fallback survives the clip/mask path the service
actually uses rather than only the helper in isolation.
"""
from __future__ import annotations
from pathlib import Path
from uuid import uuid4
import pytest
np = pytest.importorskip("numpy")
rasterio = pytest.importorskip("rasterio")
from pyproj import Transformer
from rasterio.transform import from_origin
from app.core.config import Settings
from app.models import Dataset
from app.schemas.flood_hazard import FloodHazardSelectionRequest
from app.services.flood_hazard_acquisition_service import FloodHazardAcquisitionService
from app.services.flood_hazard_analysis_service import FloodHazardAnalysisService
TO_4326 = Transformer.from_crs("EPSG:31370", "EPSG:4326", always_xy=True)
PRODUCT_KEY = "fluviaal_current_t100"
class FakeSession:
def __init__(self, objects):
self.objects = objects
def get(self, model, item_id):
return self.objects.get((model, item_id))
def _write_raster(path: Path, *, resolution: float, depth: float) -> None:
values = np.full((4, 4), depth, dtype="float32")
with rasterio.open(
path,
"w",
driver="GTiff",
width=4,
height=4,
count=1,
dtype="float32",
crs="EPSG:31370",
transform=from_origin(200_000, 210_000, resolution, resolution),
nodata=-9999.0,
) as output:
output.write(values, 1)
def _dataset(project_id, dataset_id, path: Path) -> Dataset:
return Dataset(
id=dataset_id,
project_id=project_id,
name="vmm-flood.tif",
dataset_type="raster",
source="vmm",
source_name=FloodHazardAcquisitionService.PROVIDER,
status="ready",
storage_path=str(path),
source_metadata={"product_key": PRODUCT_KEY, "normalized_value_unit": "m"},
)
def _bbox_for(min_x: float, min_y: float, max_x: float, max_y: float) -> dict:
left, bottom = TO_4326.transform(min_x, min_y)
right, top = TO_4326.transform(max_x, max_y)
return {"min_x": left, "min_y": bottom, "max_x": right, "max_y": top, "crs": "EPSG:4326"}
def _analyze(tmp_path: Path, bbox: dict, *, resolution: float = 100.0) -> dict:
project_id = uuid4()
dataset_id = uuid4()
path = tmp_path / "flood.tif"
_write_raster(path, resolution=resolution, depth=2.0)
dataset = _dataset(project_id, dataset_id, path)
db = FakeSession({(Dataset, dataset_id): dataset})
return FloodHazardAnalysisService.analyze(
db,
project_id,
dataset_id,
FloodHazardSelectionRequest(bbox=bbox),
settings=Settings(_env_file=None),
)
def test_a_selection_smaller_than_one_cell_reports_the_cell_it_touches(tmp_path: Path) -> None:
# A 40 x 30 m rectangle wholly inside one 100 m cell: no cell centre falls
# inside it, so the centre rule alone would report an empty selection.
result = _analyze(tmp_path, _bbox_for(200_010, 209_960, 200_050, 209_990))
assert result["inundated_cell_count"] == 1
assert result["inundated_fraction"] == pytest.approx(1.0)
assert "kleiner dan één rastercel" in result["coverage_warning"]
def test_a_normal_selection_is_unaffected(tmp_path: Path) -> None:
result = _analyze(tmp_path, _bbox_for(200_000, 209_700, 200_300, 210_000))
assert result["inundated_cell_count"] >= 9
assert result["coverage_warning"] is None
def test_the_reported_area_matches_the_cells_that_were_analysed(tmp_path: Path) -> None:
result = _analyze(tmp_path, _bbox_for(200_010, 209_960, 200_050, 209_990))
metrics = {item["metric_key"]: item["metric_value"] for item in result["summary"]["metrics"]}
# One 100 x 100 m cell, not the 0.12 ha that was drawn.
assert metrics["modelled_inundated_area_ha"] == pytest.approx(1.0)
assert metrics["selection_area_ha"] == pytest.approx(1.0)