Jens/geointel
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

correct the tiled inference chain and move runs off the request thread

Browse Source 
Tile handling produced results that were wrong before any model quality
question arose:

- orthophoto tiles reached the model through PIL convert("RGB"), which
  truncates the high byte of a 16-bit product and treats a 4-band RGB+NIR
  tile's infrared channel as colour. Tiles are now read with rasterio, the
  visible bands are chosen explicitly, and values are percentile-stretched
  across all three bands together so hue is preserved;
- an object wider than the tile overlap was truncated by both tiles into two
  boxes that barely intersect, so IoU suppression kept both: two false
  positives and one missed footprint per seam building. Suppression now also
  compares overlap against the smaller box, and boxes cut by an interior tile
  edge are dropped in favour of the neighbouring tile's complete view;
- georeferencing fell back to an assumed EPSG:4326 when a manifest carried no
  CRS, producing geometry that renders plausibly in the wrong place. QA
  already refused such a tile; inference now fails closed too.

Segmentation QA scored candidates against every reference feature in the
dataset, so every building outside the inferred tiles counted as a false
negative. It now applies the same persisted tile coverage that detection QA
has always used, including the indexed ST_Intersects prefilter.

Duplicate suppression uses an STRtree instead of the O(n^2) scan, tiles are
predicted in batches of YOLO_BATCH_SIZE (a setting that existed but was never
read), and detection/segmentation runs can be queued through /run-async for a
polling background worker rather than holding an HTTP worker thread for
minutes of GPU work.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
This commit is contained in:
Jens
2026-08-22 14:32:44 +02:00
co-authored by Claude Opus 5
parent 2b968b74cf
commit 08188005bd
19 changed files with 1727 additions and 112 deletions
Download Patch File Download Diff File Expand all files Collapse all files
backend/app/services/detection_service.py
+317 -26
View File
@@ -9,7 +9,11 @@ from typing import Any
from typing import Type
from geoalchemy2.shape import from_shape, to_shape
from pyproj import Transformer
from shapely.geometry import box as shapely_box
from shapely.geometry import mapping, shape
from shapely.ops import transform as shapely_transform
from shapely.strtree import STRtree
from sqlalchemy import func
from app.core.config import Settings, get_settings
@@ -18,6 +22,7 @@ from app.core.request_context import get_request_id
from app.models import AnalysisRun, Area, Dataset, Detection, Job, Project, VectorFeature
from app.schemas.detection import DetectionListResponse, DetectionRead, DetectionRunListResponse, DetectionRunRead, DetectionRunResponse
from app.services.detection_georeferencing import pixel_bbox_to_epsg4326_polygon
from app.services.detection_metrics_service import DetectionMetricsService
from app.services.detection_qa_service import DetectionQaService
from app.services.dataset_consumption_gate_service import DatasetConsumptionGate
from app.services.model_asset_catalog_service import ModelAssetCatalogService
@@ -51,22 +56,11 @@ class DetectionService:
parameters_json: dict[str, Any] | None = None,
settings: Settings | None = None,
yolo_adapter_class: Type[YoloDetectionAdapter] = YoloDetectionAdapter,
existing_job: Job | None = None,
) -> DetectionRunResponse:
parameters = dict(parameters_json or {})
resolved_settings = settings or get_settings()
project = db.get(Project, project_id)
if not project:
raise AppError(code="PROJECT_NOT_FOUND", message="Project not found", status_code=404)
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":
raise AppError(
code="INVALID_DATASET_TYPE",
message="Detection requires a raster dataset",
details={"dataset_type": dataset.dataset_type},
status_code=400,
)
dataset = DetectionService._validate_run_request(db, project_id=project_id, dataset_id=dataset_id)
TemporalCompatibilityService.ensure_detection_source_supported(dataset)
selected_model_asset = None
@@ -122,7 +116,7 @@ class DetectionService:
"tile_manifest_path": tile_manifest_path,
"parameters_json": parameters,
}
job = DetectionService._create_job(db, project_id, dataset_id, run_parameters)
job = DetectionService._create_job(db, project_id, dataset_id, run_parameters, existing_job=existing_job)
analysis_run = DetectionService._create_analysis_run(db, project_id, dataset_id, job.id, model, run_parameters)
logger.info(
"detection_started request_id=%s project_id=%s dataset_id=%s job_id=%s analysis_run_id=%s model_id=%s",
@@ -419,7 +413,19 @@ class DetectionService:
class_name=class_name,
min_confidence=min_confidence,
)
raw_candidate_geometries = [({"id": str(row.id), "class_name": row.class_name}, to_shape(row.geometry)) for row in detections]
raw_candidate_geometries = [
(
{
"id": str(row.id),
"class_name": row.class_name,
# Confidence lets the matcher rank candidates the way
# detection benchmarks do instead of by row order.
"confidence": row.confidence,
},
to_shape(row.geometry),
)
for row in detections
]
candidate_geometries = raw_candidate_geometries
coverage = None
@@ -510,10 +516,31 @@ class DetectionService:
reference_envelopes,
iou_threshold,
)
candidate_geometry_mode = DetectionQaService.candidate_geometry_mode(candidate_geometries)
box_to_footprint_diagnostics = DetectionQaService.box_to_footprint_diagnostics(
evidence,
envelope_evidence,
iou_threshold=iou_threshold,
candidate_geometry_mode=candidate_geometry_mode,
)
box_to_footprint_diagnostics["envelope_precision_recall_curve"] = (
DetectionMetricsService.precision_recall_curve(
candidate_geometries,
reference_envelopes,
iou_threshold=iou_threshold,
)
)
if candidate_geometry_mode == "axis_aligned_boxes":
coverage_warnings.append(
"Candidates are axis-aligned detector boxes; strict footprint IoU cannot reach 1 for "
"rotated or non-rectangular buildings. See box_to_footprint_diagnostics."
)
# Threshold-independent view of the same populations, so the run can be
# compared with another model instead of only with itself.
precision_recall_curve = DetectionMetricsService.precision_recall_curve(
candidate_geometries,
reference_geometries,
iou_threshold=iou_threshold,
)
mean_iou = None if not evidence.match_iou_values else sum(evidence.match_iou_values) / len(evidence.match_iou_values)
precision = evidence.matches / (evidence.matches + evidence.false_positives) if evidence.matches + evidence.false_positives > 0 else None
@@ -549,6 +576,7 @@ class DetectionService:
"coverage": coverage_summary,
"temporal_compatibility": temporal_compatibility,
"box_to_footprint_diagnostics": box_to_footprint_diagnostics,
"precision_recall_curve": precision_recall_curve,
"match_evidence": evidence.match_evidence,
"false_positive_evidence": evidence.false_positive_evidence,
"false_negative_evidence": evidence.false_negative_evidence,
@@ -560,6 +588,9 @@ class DetectionService:
"mean_iou": mean_iou,
"false_positive_count": evidence.false_positives,
"false_negative_count": evidence.false_negatives,
"average_precision": precision_recall_curve["average_precision"],
"best_f1": precision_recall_curve["best_f1"],
"best_f1_threshold": precision_recall_curve["best_f1_threshold"],
},
)
logger.info(
@@ -594,13 +625,32 @@ class DetectionService:
"coverage": coverage_summary,
"temporal_compatibility": temporal_compatibility,
"box_to_footprint_diagnostics": box_to_footprint_diagnostics,
"precision_recall_curve": precision_recall_curve,
"match_evidence": evidence.match_evidence,
"false_positive_evidence": evidence.false_positive_evidence,
"false_negative_evidence": evidence.false_negative_evidence,
}
@staticmethod
def _create_job(db, project_id: uuid.UUID, dataset_id: uuid.UUID, parameters: dict[str, Any]) -> Job:
def _create_job(
db,
project_id: uuid.UUID,
dataset_id: uuid.UUID,
parameters: dict[str, Any],
existing_job: Job | None = None,
) -> Job:
if existing_job is not None:
# A queued job already represents this run; reuse it so the client
# keeps polling one identifier from request to result.
existing_job.status = "running"
existing_job.dataset_id = dataset_id
existing_job.input_dataset_id = dataset_id
existing_job.parameters_json = {**(existing_job.parameters_json or {}), **parameters}
existing_job.started_at = DetectionService._now()
db.add(existing_job)
db.commit()
db.refresh(existing_job)
return existing_job
job = Job(
id=uuid.uuid4(),
job_type="detection.run",
@@ -616,6 +666,78 @@ class DetectionService:
db.refresh(job)
return job
@staticmethod
def enqueue_detection(
db,
project_id: uuid.UUID,
dataset_id: uuid.UUID,
model_id: str,
confidence_threshold: float,
model_asset_id: str | None = None,
class_filter: list[str] | None = None,
tile_manifest_path: str | None = None,
parameters_json: dict[str, Any] | None = None,
) -> Job:
"""Accept a detection run for background execution.
Everything cheap enough to answer inside the request is checked here,
so an operator learns about a missing dataset or an unvalidated class
immediately rather than from a job that fails minutes later.
"""
DetectionService._validate_run_request(
db,
project_id=project_id,
dataset_id=dataset_id,
)
job = Job(
id=uuid.uuid4(),
job_type="detection.run",
status="queued",
project_id=project_id,
dataset_id=dataset_id,
input_dataset_id=dataset_id,
parameters_json={
"project_id": str(project_id),
"dataset_id": str(dataset_id),
"model_id": model_id,
"model_asset_id": model_asset_id,
"confidence_threshold": confidence_threshold,
"class_filter": class_filter or [],
"tile_manifest_path": tile_manifest_path,
"parameters_json": dict(parameters_json or {}),
},
)
db.add(job)
db.commit()
db.refresh(job)
logger.info(
"detection_queued request_id=%s project_id=%s dataset_id=%s job_id=%s model_id=%s",
get_request_id(),
project_id,
dataset_id,
job.id,
model_id,
)
return job
@staticmethod
def _validate_run_request(db, *, project_id: uuid.UUID, dataset_id: uuid.UUID) -> Dataset:
project = db.get(Project, project_id)
if not project:
raise AppError(code="PROJECT_NOT_FOUND", message="Project not found", status_code=404)
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":
raise AppError(
code="INVALID_DATASET_TYPE",
message="Detection requires a raster dataset",
details={"dataset_type": dataset.dataset_type},
status_code=400,
)
return dataset
@staticmethod
def _query_detection_rows(
db,
@@ -634,7 +756,14 @@ class DetectionService:
query = query.filter(Detection.class_name == class_name)
if min_confidence is not None:
query = query.filter(Detection.confidence >= min_confidence)
return query.order_by(Detection.created_at.desc()).all()
# ``created_at`` defaults to the transaction timestamp, so every
# detection in a run shares one value and ordering by it alone leaves
# the row order undefined. Confidence first, id as a stable tiebreak.
return query.order_by(
Detection.confidence.desc(),
Detection.created_at.desc(),
Detection.id.asc(),
).all()
@staticmethod
def _detection_properties(detection: Detection) -> dict[str, Any]:
@@ -792,10 +921,20 @@ class DetectionService:
model = adapter.load_model(model_path)
allowed_classes = {DetectionService._canonical_class_name(value) for value in class_filter if DetectionService._canonical_class_name(value)}
candidates: list[dict[str, Any]] = []
manifest_crs = manifest.get("crs") or manifest.get("source_crs") or manifest.get("dataset_crs") or "EPSG:4326"
for tile in manifest["tiles"]:
tile_path = DetectionService._resolve_tile_path(tile, Path(tile_manifest_path or "").expanduser())
for raw in adapter.predict_tile(model, tile_path, confidence_threshold):
manifest_crs = DetectionService._require_manifest_crs(manifest)
raster_bounds = DetectionService._bounds_to_epsg4326(manifest.get("bounds"), manifest_crs)
tiles = list(manifest["tiles"])
tile_paths = [
DetectionService._resolve_tile_path(tile, Path(tile_manifest_path or "").expanduser()) for tile in tiles
]
# Batched so the GPU is not idle between tiles; each tile keeps its own
# transform for georeferencing, so results stay per tile and in order.
detections_per_tile = adapter.predict_tiles(model, tile_paths, confidence_threshold)
for tile, tile_path, raw_detections in zip(tiles, tile_paths, detections_per_tile):
tile_crs = tile.get("crs") or manifest_crs
tile_bounds_4326 = DetectionService._bounds_to_epsg4326(tile.get("bounds"), tile_crs)
tile_edge_tolerance = DetectionService._tile_edge_tolerance(tile, tile_bounds_4326)
for raw in raw_detections:
model_class_name = str(raw.get("class_name") or "").strip()
class_name = DetectionService._canonical_class_name(model_class_name)
confidence = float(raw.get("confidence", 0.0))
@@ -806,7 +945,7 @@ class DetectionService:
bbox = raw.get("bbox")
if not isinstance(bbox, list):
raise AppError(code="DETECTION_INVALID_BBOX", message="YOLO adapter returned a detection without bbox", status_code=422)
geometry = pixel_bbox_to_epsg4326_polygon(bbox=bbox, tile=tile, crs=tile.get("crs") or manifest_crs)
geometry = pixel_bbox_to_epsg4326_polygon(bbox=bbox, tile=tile, crs=tile_crs)
properties = dict(raw.get("properties") or {})
if model_class_name and model_class_name != class_name:
properties.setdefault("model_class_name", model_class_name)
@@ -818,10 +957,19 @@ class DetectionService:
"bbox": bbox,
"source_tile_path": str(tile_path),
"properties": {**properties, "tile_index": tile.get("index")},
"tile_bounds": tile_bounds_4326,
"tile_edge_tolerance": tile_edge_tolerance,
}
)
edge_filtered_candidates = candidates
if settings.yolo_suppress_tile_edge_detections:
edge_filtered_candidates = DetectionService._drop_tile_edge_truncations(
candidates,
raster_bounds=raster_bounds,
tolerance=0.0,
)
filtered_candidates = DetectionService._suppress_duplicate_candidates(
candidates,
edge_filtered_candidates,
iou_threshold=float(settings.yolo_duplicate_iou_threshold),
)
persisted: list[Detection] = []
@@ -858,7 +1006,9 @@ class DetectionService:
return persisted, {
"raw_detection_count": len(candidates),
"suppressed_detection_count": len(candidates) - len(filtered_candidates),
"tile_edge_truncated_count": len(candidates) - len(edge_filtered_candidates),
"duplicate_iou_threshold": float(settings.yolo_duplicate_iou_threshold),
"containment_suppression_threshold": DetectionService.CONTAINMENT_SUPPRESSION_THRESHOLD,
"runtime_model_provenance": runtime_model_provenance.as_dict(),
}
@@ -887,21 +1037,90 @@ class DetectionService:
def _canonical_class_name(value: Any) -> str:
return str(value or "").strip().casefold()
# An object wider than the tile overlap is truncated by both tiles, so the
# two halves barely intersect and IoU alone never suppresses them. Overlap
# measured against the smaller box catches that case; the threshold is
# deliberately strict so that terraced houses stay separate detections.
CONTAINMENT_SUPPRESSION_THRESHOLD = 0.85
@staticmethod
def _suppress_duplicate_candidates(candidates: list[dict[str, Any]], iou_threshold: float) -> list[dict[str, Any]]:
if iou_threshold <= 0 or len(candidates) < 2:
return candidates
ordered = sorted(
candidates,
key=lambda item: (-float(item["confidence"]), str(item.get("source_tile_path") or "")),
)
kept: list[dict[str, Any]] = []
for candidate in sorted(candidates, key=lambda item: float(item["confidence"]), reverse=True):
kept_geometries: list[Any] = []
tree = None
for candidate in ordered:
geometry = candidate["geometry"]
duplicate = False
for kept_candidate in kept:
# Only geometries that actually touch this candidate can suppress
# it, so an index keeps a dense AOI from turning into an O(n^2) scan.
neighbour_indexes = range(len(kept)) if tree is None else (int(index) for index in tree.query(geometry))
for index in neighbour_indexes:
kept_candidate = kept[index]
if candidate["class_name"] != kept_candidate["class_name"]:
continue
if DetectionService._geometry_iou(candidate["geometry"], kept_candidate["geometry"]) >= iou_threshold:
other = kept_geometries[index]
if DetectionService._geometry_iou(geometry, other) >= iou_threshold:
duplicate = True
break
if (
DetectionService._geometry_containment(geometry, other)
>= DetectionService.CONTAINMENT_SUPPRESSION_THRESHOLD
):
duplicate = True
break
if not duplicate:
kept.append(candidate)
kept_geometries.append(geometry)
tree = STRtree(kept_geometries)
return kept
@staticmethod
def _drop_tile_edge_truncations(
candidates: list[dict[str, Any]],
*,
raster_bounds: tuple[float, float, float, float] | None,
tolerance: float,
) -> list[dict[str, Any]]:
"""Discard boxes cut off by an interior tile edge.
Such a box describes only the part of the object that fell inside its
tile. Because tiles overlap, the neighbouring tile saw the object whole
and contributed the box worth keeping. A box against the outer raster
edge has no such neighbour and is kept.
"""
if raster_bounds is None or tolerance <= 0:
return candidates
raster_left, raster_bottom, raster_right, raster_top = raster_bounds
kept: list[dict[str, Any]] = []
for candidate in candidates:
tile_bounds = candidate.get("tile_bounds")
if not tile_bounds or len(tuple(tile_bounds)) != 4:
kept.append(candidate)
continue
tile_left, tile_bottom, tile_right, tile_top = (float(value) for value in tile_bounds)
left, bottom, right, top = candidate["geometry"].bounds
# A pixel-sized tolerance per tile: a fixed degree value would be
# wrong for both a 10 cm orthophoto and a coarse thematic raster.
tolerance = float(candidate.get("tile_edge_tolerance") or 0.0) or tolerance
touches_interior_edge = (
(abs(left - tile_left) <= tolerance and abs(tile_left - raster_left) > tolerance)
or (abs(right - tile_right) <= tolerance and abs(tile_right - raster_right) > tolerance)
or (abs(bottom - tile_bottom) <= tolerance and abs(tile_bottom - raster_bottom) > tolerance)
or (abs(top - tile_top) <= tolerance and abs(tile_top - raster_top) > tolerance)
)
if not touches_interior_edge:
kept.append(candidate)
return kept
@staticmethod
@@ -916,6 +1135,78 @@ class DetectionService:
return 0.0
return intersection_area / union_area
@staticmethod
def _geometry_containment(left, right) -> float:
"""Intersection over the smaller of the two areas."""
if left.is_empty or right.is_empty:
return 0.0
smaller_area = min(left.area, right.area)
if smaller_area <= 0:
return 0.0
intersection_area = left.intersection(right).area
if intersection_area <= 0:
return 0.0
return intersection_area / smaller_area
@staticmethod
def _require_manifest_crs(manifest: dict[str, Any]) -> str:
"""Refuse to georeference inference output against a guessed CRS.
Detection QA already rejects a tile without explicit CRS metadata.
Silently assuming EPSG:4326 on the inference side produced geometry
that looks plausible on a map but sits in the wrong place.
"""
raw_crs = manifest.get("crs") or manifest.get("source_crs") or manifest.get("dataset_crs")
if not isinstance(raw_crs, str) or not raw_crs.strip():
raise AppError(
code="DETECTION_TILE_MANIFEST_INVALID",
message="Raster tile manifest requires explicit CRS metadata for georeferencing",
status_code=422,
)
return raw_crs.strip()
@staticmethod
def _bounds_to_epsg4326(bounds: Any, crs: str | None) -> tuple[float, float, float, float] | None:
if not isinstance(bounds, (list, tuple)) or len(bounds) != 4:
return None
try:
left, bottom, right, top = (float(value) for value in bounds)
except (TypeError, ValueError):
return None
if left >= right or bottom >= top:
return None
if not crs or str(crs).strip().upper() in {"EPSG:4326", "4326"}:
return (left, bottom, right, top)
try:
transformer = Transformer.from_crs(crs, "EPSG:4326", always_xy=True)
# Transform the whole rectangle, not just two corners: a projected
# box does not stay axis-aligned after reprojection.
projected = shapely_transform(transformer.transform, shapely_box(left, bottom, right, top))
return projected.bounds
except Exception:
return None
@staticmethod
def _tile_edge_tolerance(tile: dict[str, Any], tile_bounds_4326: tuple[float, float, float, float] | None) -> float:
"""One and a half pixels, expressed in the degrees the boxes live in."""
if tile_bounds_4326 is None:
return 0.0
pixel_window = tile.get("pixel_window")
if not (isinstance(pixel_window, (list, tuple)) and len(pixel_window) == 4):
return 0.0
try:
width = float(pixel_window[2])
height = float(pixel_window[3])
except (TypeError, ValueError):
return 0.0
if width <= 0 or height <= 0:
return 0.0
left, bottom, right, top = tile_bounds_4326
return 1.5 * max((right - left) / width, (top - bottom) / height)
@staticmethod
def _load_tile_manifest(tile_manifest_path: str | None, max_tiles: int) -> dict[str, Any]:
if not tile_manifest_path: