"""Manual Ground-Truth (§7 Hierarchy Level 2) Held-Out Point Evaluation Harness. Computes: - Per-point displacement (in pixels and meters) - RMSE (x, y, total) - Success Rate with tolerance 1x GSD and 2x GSD (<= 1 px and <= 2 px) - Confidence tier breakdown (Tier 1: Certain, Tier 2: Pretty Sure, Tier 3: Best Guess) - Feature type breakdown (crater_rim, boulder, crater_center, ridge) Note: NEVER mix held-out check points into any fit/transformation calculation. """ import argparse import csv import json import math import sys from pathlib import Path from typing import Dict, List, Tuple, Optional def load_points_from_csv(file_path: Path) -> List[Dict]: """Load points from a QGIS-exported CSV file.""" if not file_path.exists(): raise FileNotFoundError(f"CSV file not found: {file_path}") points = [] with open(file_path, "r", encoding="utf-8-sig") as f: reader = csv.DictReader(f) for idx, row in enumerate(reader): # QGIS exports geometry as X, Y or WKT, or custom x, y columns x_val = row.get("X") or row.get("x") or row.get("lon") or row.get("X84") y_val = row.get("Y") or row.get("y") or row.get("lat") or row.get("Y84") if x_val is None or y_val is None: # Check if WKT geometry field exists (e.g. "POINT (x y)") wkt = row.get("WKT") or row.get("wkt") or row.get("geometry") if wkt and "POINT" in wkt.upper(): try: coords = wkt.replace("POINT", "").replace("(", "").replace(")", "").strip().split() x_val, y_val = coords[0], coords[1] except Exception: pass if x_val is None or y_val is None: raise ValueError(f"Could not parse X/Y coordinates at row {idx+1} in {file_path}. Row keys: {list(row.keys())}") feat_type = row.get("feature_type") or row.get("type") or "unknown" conf = row.get("confidence") or row.get("conf") or "1" notes = row.get("notes") or row.get("note") or "" points.append({ "index": idx, "x": float(x_val), "y": float(y_val), "feature_type": str(feat_type).strip().lower(), "confidence": int(float(conf)), "notes": str(notes).strip() }) return points def evaluate_manual_pairs( ref_points: List[Dict], prod_points: List[Dict], gsd_m: float = 0.25 ) -> Dict: """Compute per-point displacement, RMSE, and Success Rates from paired reference/product CSVs.""" if len(ref_points) != len(prod_points): raise ValueError( f"Row count mismatch: reference CSV has {len(ref_points)} rows, " f"product CSV has {len(prod_points)} rows. Both files must pair by row order!" ) if len(ref_points) == 0: raise ValueError("No points found to evaluate!") per_point_results = [] for i, (ref, prod) in enumerate(zip(ref_points, prod_points)): dx = prod["x"] - ref["x"] dy = prod["y"] - ref["y"] disp_px = math.hypot(dx, dy) disp_m = disp_px * gsd_m # Use confidence and feature_type from reference point (or fallback to product) conf = ref.get("confidence", prod.get("confidence", 1)) feat_type = ref.get("feature_type", prod.get("feature_type", "unknown")) per_point_results.append({ "point_id": i + 1, "ref_x": ref["x"], "ref_y": ref["y"], "prod_x": prod["x"], "prod_y": prod["y"], "dx_px": dx, "dy_px": dy, "disp_px": disp_px, "disp_m": disp_m, "confidence": conf, "feature_type": feat_type, "pass_1x_gsd": disp_px <= 1.0, "pass_2x_gsd": disp_px <= 2.0, }) def compute_stats_for_subset(subset: List[Dict]) -> Dict: if not subset: return { "count": 0, "rmse_x_px": None, "rmse_y_px": None, "rmse_total_px": None, "rmse_m": None, "success_rate_1x_gsd": None, "success_rate_2x_gsd": None, } n = len(subset) sum_dx2 = sum(p["dx_px"] ** 2 for p in subset) sum_dy2 = sum(p["dy_px"] ** 2 for p in subset) sum_disp2 = sum(p["disp_px"] ** 2 for p in subset) rmse_x = math.sqrt(sum_dx2 / n) rmse_y = math.sqrt(sum_dy2 / n) rmse_total = math.sqrt(sum_disp2 / n) rmse_m = rmse_total * gsd_m n_pass_1x = sum(1 for p in subset if p["pass_1x_gsd"]) n_pass_2x = sum(1 for p in subset if p["pass_2x_gsd"]) return { "count": n, "rmse_x_px": round(rmse_x, 4), "rmse_y_px": round(rmse_y, 4), "rmse_total_px": round(rmse_total, 4), "rmse_m": round(rmse_m, 4), "success_rate_1x_gsd": round(n_pass_1x / n, 4), "success_rate_2x_gsd": round(n_pass_2x / n, 4), "n_pass_1x_gsd": n_pass_1x, "n_pass_2x_gsd": n_pass_2x, } # Overall stats overall_stats = compute_stats_for_subset(per_point_results) # By Confidence Tier (1, 2, 3) by_confidence = {} for tier in [1, 2, 3]: sub = [p for p in per_point_results if p["confidence"] == tier] tier_names = {1: "1_certain", 2: "2_pretty_sure", 3: "3_best_guess"} by_confidence[tier_names[tier]] = compute_stats_for_subset(sub) # By Feature Type feature_types = set(p["feature_type"] for p in per_point_results) by_feature = {} for ft in sorted(feature_types): sub = [p for p in per_point_results if p["feature_type"] == ft] by_feature[ft] = compute_stats_for_subset(sub) return { "hierarchy_level": 2, "evaluation_type": "Manual Ground-Truth Check Set", "gsd_m": gsd_m, "total_points": len(per_point_results), "overall": overall_stats, "by_confidence_tier": by_confidence, "by_feature_type": by_feature, "per_point_details": per_point_results, } def generate_markdown_report(report: Dict) -> str: """Generate clean Markdown table output for human review.""" lines = [] lines.append("# Manual Ground-Truth Held-Out Point Evaluation (§7 Hierarchy Level 2)") lines.append("") lines.append(f"- **Total Check Points**: `{report['total_points']}`") lines.append(f"- **GSD (Ground Sample Distance)**: `{report['gsd_m']} m/px`") lines.append(f"- **Evaluation Protocol**: Row-paired manual features on Reference vs Registered Product") lines.append("") # Overall summary table lines.append("## Overall Metrics") lines.append("") ov = report["overall"] lines.append("| Metric | Value (px) | Value (meters) | Notes / Criteria |") lines.append("|---|---|---|---|") lines.append(f"| **RMSE (Total)** | `{ov['rmse_total_px']} px` | `{ov['rmse_m']} m` | Radial displacement |") lines.append(f"| **RMSE (X)** | `{ov['rmse_x_px']} px` | `{round(ov['rmse_x_px']*report['gsd_m'], 4)} m` | Easting / Along-row |") lines.append(f"| **RMSE (Y)** | `{ov['rmse_y_px']} px` | `{round(ov['rmse_y_px']*report['gsd_m'], 4)} m` | Northing / Along-col |") lines.append(f"| **Success Rate (1× GSD)** | `{ov['success_rate_1x_gsd']*100:.1f}%` | — | `{ov['n_pass_1x_gsd']}/{ov['count']}` points $\\le 1.0$ px |") lines.append(f"| **Success Rate (2× GSD)** | `{ov['success_rate_2x_gsd']*100:.1f}%` | — | `{ov['n_pass_2x_gsd']}/{ov['count']}` points $\\le 2.0$ px |") lines.append("") # Confidence Tier breakdown table lines.append("## Metrics by Confidence Tier") lines.append("") lines.append("| Tier | Label | Count | RMSE (px) | RMSE (m) | Success (1× GSD) | Success (2× GSD) |") lines.append("|---|---|---|---|---|---|---|") for tier_key, stats in report["by_confidence_tier"].items(): if stats["count"] == 0: lines.append(f"| `{tier_key}` | — | 0 | N/A | N/A | N/A | N/A |") else: lines.append( f"| `{tier_key}` | `{tier_key.split('_', 1)[1].replace('_', ' ').title()}` | {stats['count']} | " f"`{stats['rmse_total_px']} px` | `{stats['rmse_m']} m` | " f"`{stats['success_rate_1x_gsd']*100:.1f}%` | `{stats['success_rate_2x_gsd']*100:.1f}%` |" ) lines.append("") # Feature Type breakdown table lines.append("## Metrics by Feature Type") lines.append("") lines.append("| Feature Type | Count | RMSE (px) | RMSE (m) | Success (1× GSD) | Success (2× GSD) |") lines.append("|---|---|---|---|---|---|") for ft_key, stats in report["by_feature_type"].items(): lines.append( f"| `{ft_key}` | {stats['count']} | " f"`{stats['rmse_total_px']} px` | `{stats['rmse_m']} m` | " f"`{stats['success_rate_1x_gsd']*100:.1f}%` | `{stats['success_rate_2x_gsd']*100:.1f}%` |" ) lines.append("") lines.append("> [!IMPORTANT]") lines.append("> **Circularity Guard**: These ground-truth points were picked manually and evaluated exclusively on the final registered product. They were NEVER included in the registration fitting set.") lines.append("") return "\n".join(lines) def main(): parser = argparse.ArgumentParser(description="Manual Ground-Truth Held-Out Point Evaluation") parser.add_argument("--ref", type=str, help="Path to heldout_ref.csv (reference image picks)") parser.add_argument("--prod", type=str, help="Path to heldout_prod.csv (registered product picks)") parser.add_argument("--single", type=str, help="Path to single paired CSV (e.g. ref_x, ref_y, prod_x, prod_y)") parser.add_argument("--gsd", type=float, default=0.25, help="Ground Sample Distance in meters per pixel (default: 0.25)") parser.add_argument("--out-json", type=str, default="outputs/heldout_metrics.json", help="Path to output JSON") parser.add_argument("--out-md", type=str, default="outputs/heldout_report.md", help="Path to output Markdown") args = parser.parse_args() if args.single: # Load single paired file single_path = Path(args.single) if not single_path.exists(): print(f"Error: Single CSV file '{single_path}' does not exist.") sys.exit(1) ref_pts, prod_pts = [], [] with open(single_path, "r", encoding="utf-8-sig") as f: reader = csv.DictReader(f) for idx, r in enumerate(reader): ref_pts.append({ "x": float(r["ref_x"]), "y": float(r["ref_y"]), "confidence": int(float(r.get("confidence", 1))), "feature_type": r.get("feature_type", "unknown") }) prod_pts.append({ "x": float(r["prod_x"]), "y": float(r["prod_y"]), "confidence": int(float(r.get("confidence", 1))), "feature_type": r.get("feature_type", "unknown") }) elif args.ref and args.prod: ref_pts = load_points_from_csv(Path(args.ref)) prod_pts = load_points_from_csv(Path(args.prod)) else: print("Error: Must specify either (--ref AND --prod) or --single") sys.exit(1) res = evaluate_manual_pairs(ref_pts, prod_pts, gsd_m=args.gsd) md = generate_markdown_report(res) Path(args.out_json).parent.mkdir(parents=True, exist_ok=True) with open(args.out_json, "w", encoding="utf-8") as f: json.dump(res, f, indent=2) with open(args.out_md, "w", encoding="utf-8") as f: f.write(md) print(md) print(f"\nSaved metrics to {args.out_json} and report to {args.out_md}") if __name__ == "__main__": main()