bound change detection to the operator's selection

Change detection was the one analysis that ignored the selection entirely. It
compared two datasets in full, loaded every feature of both into Python with no
spatial predicate, and — with include_unchanged defaulting to true — returned a
FeatureCollection holding both datasets. For a regional building layer that is
the wrong answer to "what changed here" and a response no browser should be
asked to hold.

It now accepts bbox and area_id, resolved the way every other analysis resolves
them, and loads through an indexed ST_Intersects predicate.

Features are deliberately not clipped to the selection. A change class
describes a whole object: comparing a clipped earlier footprint against an
unclipped later one would report the selection edge itself as a change. Objects
the edge crosses are compared in full and counted in a warning.

The returned geometry is capped by preview_limit, spending that budget on
modified, added and removed before unchanged, while every count still describes
the whole selection.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
This commit is contained in:
Jens
2026-08-22 14:56:15 +02:00
co-authored by Claude Opus 5
parent d29d572e8c
commit 12aaf1bb4d
4 changed files with 295 additions and 7 deletions
@@ -4,11 +4,12 @@ from datetime import datetime, timezone
from typing import Any
from uuid import UUID
from geoalchemy2.shape import to_shape
from geoalchemy2.shape import from_shape, to_shape
from shapely.geometry import mapping
from shapely.geometry.base import BaseGeometry
from shapely.strtree import STRtree
from shapely.validation import make_valid
from sqlalchemy import func
from sqlalchemy.orm import Session
from app.core.errors import AppError
@@ -30,6 +31,9 @@ class ChangeDetectionService:
iou_threshold: float = 0.8,
include_unchanged: bool = True,
modified_threshold: float = 0.3,
bbox: dict[str, Any] | None = None,
area_id: UUID | None = None,
preview_limit: int = 2_000,
) -> ChangeDetectionSummary:
if source_dataset_id == target_dataset_id:
raise AppError(code="INVALID_PARAMETERS", message="Source and target datasets must differ", status_code=400)
@@ -45,14 +49,19 @@ class ChangeDetectionService:
source_dataset = ChangeDetectionService._get_project_vector_dataset(db, source_dataset_id, project_id, "Source")
target_dataset = ChangeDetectionService._get_project_vector_dataset(db, target_dataset_id, project_id, "Target")
source_features, source_warnings = ChangeDetectionService._load_features(db, source_dataset)
target_features, target_warnings = ChangeDetectionService._load_features(db, target_dataset)
selection_geometry = ChangeDetectionService._selection_geometry(db, project_id, bbox=bbox, area_id=area_id)
source_features, source_warnings = ChangeDetectionService._load_features(db, source_dataset, selection_geometry)
target_features, target_warnings = ChangeDetectionService._load_features(db, target_dataset, selection_geometry)
if not source_features:
raise AppError(code="EMPTY_VECTOR_DATASET", message="Source dataset has no comparable vector features", status_code=422)
if not target_features:
raise AppError(code="EMPTY_VECTOR_DATASET", message="Target dataset has no comparable vector features", status_code=422)
source_features = ChangeDetectionService.restrict_to_selection(source_features, selection_geometry, label="Source")
target_features = ChangeDetectionService.restrict_to_selection(target_features, selection_geometry, label="Target")
classified = ChangeDetectionService._classify_features(
source_features,
target_features,
@@ -79,7 +88,24 @@ class ChangeDetectionService:
if not include_unchanged:
buckets["unchanged"] = []
geojson_features = buckets["added"] + buckets["removed"] + buckets["modified"] + buckets["unchanged"]
geojson_features, preview_truncated = ChangeDetectionService.limit_preview(
buckets["added"] + buckets["removed"] + buckets["modified"] + buckets["unchanged"],
limit=preview_limit,
)
warnings = source_warnings + target_warnings
edge_count = sum(
1 for feature in source_features + target_features if feature.get("partially_covered")
)
if edge_count:
warnings.append(
f"{edge_count} objecten liggen deels buiten de selectie. Ze zijn volledig vergeleken, zodat de "
"selectierand zelf geen wijziging veroorzaakt."
)
if preview_truncated:
warnings.append(
f"De tellingen gelden voor de volledige selectie; de kaart toont maximaal {preview_limit} objecten, "
"wijzigingen eerst."
)
return ChangeDetectionSummary(
source_dataset_id=source_dataset_id,
target_dataset_id=target_dataset_id,
@@ -91,11 +117,109 @@ class ChangeDetectionService:
unchanged_count=unchanged_count,
iou_threshold=iou_threshold,
modified_iou_threshold=modified_threshold,
warnings=source_warnings + target_warnings,
selection_area_id=area_id,
preview_limit=preview_limit,
preview_truncated=preview_truncated,
warnings=warnings,
generated_at=datetime.now(timezone.utc),
geojson={"type": "FeatureCollection", "features": geojson_features},
)
@staticmethod
def _selection_geometry(
db: Session,
project_id: UUID,
*,
bbox: dict[str, Any] | None,
area_id: UUID | None,
) -> BaseGeometry | None:
"""Resolve the drawn rectangle against the named work area, if any."""
from app.models import Area
from shapely.geometry import box as shapely_box
selection = None
if bbox:
selection = shapely_box(
float(bbox["min_x"]), float(bbox["min_y"]), float(bbox["max_x"]), float(bbox["max_y"])
)
if area_id is None:
return selection
area = db.get(Area, 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)
area_geometry = to_shape(area.geometry)
if selection is None:
return area_geometry
intersection = selection.intersection(area_geometry)
if intersection.is_empty or intersection.area <= 0:
raise AppError(
code="CHANGE_DETECTION_SELECTION_OUTSIDE_AREA",
message="Selection does not overlap the selected work area",
status_code=422,
)
return intersection
# Order the preview spends its budget in. An operator asking what changed
# is not helped by a cap filled with unchanged footprints.
PREVIEW_PRIORITY = {"modified": 0, "added": 1, "removed": 2, "unchanged": 3}
@staticmethod
def restrict_to_selection(
features: list[dict[str, Any]],
selection_geometry: BaseGeometry | None,
*,
label: str = "Dataset",
) -> list[dict[str, Any]]:
"""Keep the features a drawn selection reaches, and say which it cuts.
Geometry is deliberately *not* clipped. A change class describes a whole
object: comparing a clipped 2020 footprint against an unclipped 2024 one
would manufacture "modified" along the selection edge. Clipping is right
for an area metric and wrong for an identity comparison.
"""
if selection_geometry is None:
return features
kept: list[dict[str, Any]] = []
for feature in features:
geometry = feature["geometry"]
if not geometry.intersects(selection_geometry):
continue
kept.append({**feature, "partially_covered": not selection_geometry.covers(geometry)})
if not kept:
raise AppError(
code="CHANGE_DETECTION_SELECTION_EMPTY",
message=f"{label} dataset has no features inside this selection",
status_code=422,
)
return kept
@staticmethod
def limit_preview(
features: list[dict[str, Any]],
*,
limit: int,
) -> tuple[list[dict[str, Any]], bool]:
"""Cap the returned geometry without capping the counts.
``include_unchanged`` defaulted to true and nothing bounded the result,
so a regional comparison returned a FeatureCollection holding both
datasets in full. The counts describe the whole selection; the preview
describes what a map can usefully draw.
"""
if limit <= 0 or len(features) <= limit:
return features, False
ordered = sorted(
features,
key=lambda item: ChangeDetectionService.PREVIEW_PRIORITY.get(item["change_type"], 9),
)
return ordered[:limit], True
@staticmethod
def _classify_features(
source_features: list[dict[str, Any]],
@@ -182,8 +306,25 @@ class ChangeDetectionService:
return dataset
@staticmethod
def _load_features(db: Session, dataset: Dataset) -> tuple[list[dict[str, Any]], list[str]]:
rows = db.query(VectorFeature).filter(VectorFeature.dataset_id == dataset.id).all()
def _load_features(
db: Session,
dataset: Dataset,
selection_geometry: BaseGeometry | None = None,
) -> tuple[list[dict[str, Any]], list[str]]:
query = db.query(VectorFeature).filter(VectorFeature.dataset_id == dataset.id)
if selection_geometry is not None and hasattr(query, "filter"):
# Bound the load in the database. Pulling a regional building layer
# into Python to then discard most of it costs memory and time for
# nothing, and the fallback below has no such option.
try:
query = query.filter(
func.ST_Intersects(VectorFeature.geometry, from_shape(selection_geometry, srid=4326))
)
except Exception:
# Lightweight unit-test sessions do not implement every spatial
# predicate; restrict_to_selection still bounds the population.
pass
rows = query.all()
warnings: list[str] = []
if rows:
return [ChangeDetectionService._row_to_feature(row) for row in rows], warnings