"""Sweep ChangeStar probability and component-area thresholds without rerunning inference."""
|
|
from __future__ import annotations
|
|
import argparse
|
import csv
|
import json
|
from pathlib import Path
|
|
import cv2
|
import numpy as np
|
import rasterio
|
from PIL import Image, ImageDraw, ImageFont
|
|
import importlib.util
|
|
|
def load_demo():
|
path = Path(__file__).with_name("run_change_detection.py")
|
spec = importlib.util.spec_from_file_location("change_detection_demo", path)
|
if spec is None or spec.loader is None:
|
raise RuntimeError("Unable to load change-detection helpers.")
|
module = importlib.util.module_from_spec(spec)
|
spec.loader.exec_module(module)
|
return module
|
|
|
def clean_mask(raw: np.ndarray, minimum_area: int) -> tuple[np.ndarray, int]:
|
cleaned = raw.copy()
|
count, labels, stats, _ = cv2.connectedComponentsWithStats((cleaned > 0).astype(np.uint8), connectivity=8)
|
for component in range(1, count):
|
if int(stats[component, cv2.CC_STAT_AREA]) < minimum_area:
|
cleaned[labels == component] = 0
|
return cleaned, count - 1
|
|
|
def preview(image: np.ndarray, mask: np.ndarray, max_dimension: int = 720) -> Image.Image:
|
height, width = image.shape[:2]
|
scale = min(1.0, max_dimension / max(height, width))
|
size = (max(1, int(round(width * scale))), max(1, int(round(height * scale))))
|
base = Image.fromarray(image).resize(size, Image.Resampling.LANCZOS)
|
small_mask = Image.fromarray((mask > 0).astype(np.uint8) * 255).resize(size, Image.Resampling.NEAREST)
|
base_array = np.asarray(base).copy()
|
mask_array = np.asarray(small_mask) > 0
|
red = np.zeros_like(base_array)
|
red[..., 0] = 255
|
base_array[mask_array] = (base_array[mask_array].astype(np.float32) * 0.45 + red[mask_array].astype(np.float32) * 0.55).astype(np.uint8)
|
return Image.fromarray(base_array)
|
|
|
def run_scan(run_dir: Path, output_dir: Path, thresholds: list[float], areas: list[int]) -> dict:
|
demo = load_demo()
|
probability_path = run_dir / "change_probability.tif"
|
input_path = run_dir.parents[3] / "data" / "processed" / "00-change-detection" / run_dir.name / "after_registered.tif"
|
if not probability_path.is_file() or not input_path.is_file():
|
raise FileNotFoundError("The source run must contain change_probability.tif and processed after_registered.tif.")
|
with rasterio.open(probability_path) as probability_dataset:
|
probability = probability_dataset.read(1).astype(np.float32)
|
transform = probability_dataset.transform
|
crs = probability_dataset.crs.to_string() if probability_dataset.crs else None
|
with rasterio.open(input_path) as input_dataset:
|
image = np.transpose(input_dataset.read([1, 2, 3]), (1, 2, 0))
|
output_dir.mkdir(parents=True, exist_ok=False)
|
rows: list[dict] = []
|
contact_items: list[tuple[str, Image.Image]] = []
|
for threshold in thresholds:
|
for area in areas:
|
label = f"threshold-{threshold:.2f}_area-{area}"
|
item_dir = output_dir / label
|
item_dir.mkdir()
|
raw = (probability >= threshold).astype(np.uint8) * 255
|
cleaned, raw_components = clean_mask(raw, area)
|
cleaned_components = int(cv2.connectedComponents((cleaned > 0).astype(np.uint8), connectivity=8)[0] - 1)
|
mask_path = item_dir / "change_mask.tif"
|
with rasterio.open(mask_path, "w", driver="GTiff", height=cleaned.shape[0], width=cleaned.shape[1], count=1, dtype="uint8", transform=transform, crs=crs, compress="lzw") as destination:
|
destination.write(cleaned, 1)
|
# Vectorization is intentionally deferred: full-resolution GeoTIFF
|
# polygons are expensive and are not needed to compare thresholds.
|
features = demo._feature_summary(cleaned, probability, transform, crs is not None)
|
(item_dir / "regions.json").write_text(json.dumps({"features": features}, ensure_ascii=False, indent=2), encoding="utf-8")
|
overlay = preview(image, cleaned)
|
overlay.save(item_dir / "overlay_preview.jpg", quality=90)
|
row = {
|
"threshold": threshold,
|
"minimum_area_pixels": area,
|
"raw_components": raw_components,
|
"cleaned_components": cleaned_components,
|
"raw_changed_pixels": int((raw > 0).sum()),
|
"changed_pixels": int((cleaned > 0).sum()),
|
"changed_pixel_ratio": round(float((cleaned > 0).mean()), 6),
|
"vector_feature_count": len(features),
|
"directory": label,
|
}
|
(item_dir / "summary.json").write_text(json.dumps(row, ensure_ascii=False, indent=2), encoding="utf-8")
|
rows.append(row)
|
contact_items.append((f"T={threshold:.2f} / A={area} | {len(features)} regions", overlay))
|
(output_dir / "scan_summary.json").write_text(json.dumps({"source_run": run_dir.name, "results": rows}, ensure_ascii=False, indent=2), encoding="utf-8")
|
with (output_dir / "scan_summary.csv").open("w", newline="", encoding="utf-8-sig") as stream:
|
writer = csv.DictWriter(stream, fieldnames=list(rows[0]))
|
writer.writeheader()
|
writer.writerows(rows)
|
tile_width = max(image.width for _, image in contact_items)
|
tile_height = max(image.height for _, image in contact_items) + 30
|
sheet = Image.new("RGB", (tile_width * 3, tile_height * 4), "#202b25")
|
draw = ImageDraw.Draw(sheet)
|
for index, (label, item) in enumerate(contact_items):
|
x = (index % 3) * tile_width
|
y = (index // 3) * tile_height
|
sheet.paste(item, (x, y + 28))
|
draw.text((x + 6, y + 6), label, fill="white")
|
sheet.save(output_dir / "parameter_scan_contact_sheet.jpg", quality=92)
|
return {"source_run": run_dir.name, "output_dir": str(output_dir), "results": rows}
|
|
|
def main() -> None:
|
parser = argparse.ArgumentParser()
|
parser.add_argument("--run-dir", type=Path, required=True)
|
parser.add_argument("--output", type=Path, required=True)
|
parser.add_argument("--threshold", type=float, action="append", dest="thresholds")
|
parser.add_argument("--minimum-area", type=int, action="append", dest="areas")
|
args = parser.parse_args()
|
thresholds = args.thresholds or [0.3, 0.4, 0.5]
|
areas = args.areas or [16, 64, 256, 686]
|
result = run_scan(args.run_dir, args.output, thresholds, areas)
|
print(json.dumps(result, ensure_ascii=False, indent=2))
|
|
|
if __name__ == "__main__":
|
main()
|