diff --git a/backend/festival4d/frames.py b/backend/festival4d/frames.py index 5bfabe9..3e98edb 100644 --- a/backend/festival4d/frames.py +++ b/backend/festival4d/frames.py @@ -11,18 +11,78 @@ from __future__ import annotations import logging from pathlib import Path +import numpy as np + log = logging.getLogger("festival4d.frames") def sharpness(image) -> float: - """Variance of the Laplacian of a grayscale image (OpenCV) — higher is sharper.""" - raise NotImplementedError("lane B (M2): implement sharpness") + """Variance of the Laplacian of an image (OpenCV) — higher is sharper. + + Accepts a BGR or grayscale ``ndarray`` (as returned by ``cv2.VideoCapture.read``). + The variance of the Laplacian is the standard focus/blur measure: a sharp frame has + strong high-frequency edges (high variance), a blurred one is smooth (low variance). + """ + import cv2 + + img = np.asarray(image) + if img.ndim == 3: + img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY) + return float(cv2.Laplacian(img, cv2.CV_64F).var()) def sample_frames(video_path: Path, video_id: int, out_dir: Path, target_fps: float = 2.0, window_s: float = 0.5) -> list[Path]: """Sample the sharpest frame per ``window_s`` window from ``video_path`` (spec M2). - Returns the written JPEG paths (named ``{video_id}_{frame_idx}.jpg``). + Decodes every frame, buckets them into non-overlapping ``window_s``-second windows, and + keeps the single sharpest frame (variance of Laplacian) in each window. This yields + ~``1 / window_s`` frames per second (≈ ``target_fps`` at the defaults), biased toward the + in-focus frames COLMAP needs. Writes JPEGs named ``{video_id}_{frame_idx}.jpg`` (the + original decoded ``frame_idx``, so ``t_video = frame_idx / fps``) into ``out_dir``. + + Returns the written JPEG paths, ordered by frame index. """ - raise NotImplementedError("lane B (M2): implement sample_frames") + import cv2 + + out_dir = Path(out_dir) + out_dir.mkdir(parents=True, exist_ok=True) + + cap = cv2.VideoCapture(str(video_path)) + if not cap.isOpened(): + raise RuntimeError(f"could not open video for frame sampling: {video_path}") + try: + fps = float(cap.get(cv2.CAP_PROP_FPS)) + if not np.isfinite(fps) or fps <= 0.0: + log.warning("frames: %s reported fps=%s; falling back to 30.0", video_path, fps) + fps = 30.0 + bucket_s = float(window_s) if window_s and window_s > 0 else 0.5 + + # window_idx -> (sharpness, frame_idx, frame_bgr). Only the current best per window + # is retained, so memory stays ~O(number of windows), not O(frames). + best: dict[int, tuple[float, int, np.ndarray]] = {} + frame_idx = 0 + while True: + ok, frame = cap.read() + if not ok: + break + t = frame_idx / fps + win = int(t / bucket_s) + s = sharpness(frame) + cur = best.get(win) + if cur is None or s > cur[0]: + best[win] = (s, frame_idx, frame) + frame_idx += 1 + finally: + cap.release() + + written: list[Path] = [] + for win in sorted(best): + _, fidx, frame = best[win] + path = out_dir / f"{video_id}_{fidx}.jpg" + if not cv2.imwrite(str(path), frame): + raise RuntimeError(f"failed to write frame JPEG: {path}") + written.append(path) + log.info("frames: %s -> %d frames sampled into %s", Path(video_path).name, + len(written), out_dir) + return written diff --git a/backend/festival4d/geometry.py b/backend/festival4d/geometry.py index 2d5825c..0f670ad 100644 --- a/backend/festival4d/geometry.py +++ b/backend/festival4d/geometry.py @@ -121,7 +121,8 @@ def colmap_to_threejs(q: ArrayLike, t: ArrayLike) -> tuple[NDArray[np.float64], # --------------------------------------------------------------------------- -# Lane B stubs (spec M2 + M8). Signatures FROZEN; bodies raise NotImplementedError. +# Lane B (spec M2 + M8): pose interpolation, pixel ray casting, two-view triangulation, +# and the single-view nearest-point fallback. Signatures FROZEN; bodies implemented below. # --------------------------------------------------------------------------- def slerp_pose( q0: ArrayLike, @@ -149,7 +150,39 @@ def slerp_pose( ------- (q, t) : the interpolated quaternion ``[w, x, y, z]`` and translation ``[x, y, z]``. """ - raise NotImplementedError("lane B (M2): implement slerp_pose") + a = float(alpha) + q0 = np.asarray(q0, dtype=np.float64).reshape(4) + q1 = np.asarray(q1, dtype=np.float64).reshape(4) + n0 = np.linalg.norm(q0) + n1 = np.linalg.norm(q1) + if n0 < 1e-12 or n1 < 1e-12: + raise ValueError("slerp endpoint quaternion has near-zero norm") + q0 = q0 / n0 + q1 = q1 / n1 + + # Double cover: pick the sign of q1 that lies on the same hemisphere as q0, so slerp + # takes the shorter arc (a rotation and its negation are the same orientation). + dot = float(np.dot(q0, q1)) + if dot < 0.0: + q1 = -q1 + dot = -dot + dot = min(1.0, max(-1.0, dot)) + + if dot > 0.9995: + # Endpoints almost coincide: nlerp is numerically safe and visually identical. + q_interp = q0 + a * (q1 - q0) + else: + theta_0 = np.arccos(dot) + sin_0 = np.sin(theta_0) + s0 = np.sin((1.0 - a) * theta_0) / sin_0 + s1 = np.sin(a * theta_0) / sin_0 + q_interp = s0 * q0 + s1 * q1 + q_interp = q_interp / np.linalg.norm(q_interp) + + t0 = np.asarray(t0, dtype=np.float64).reshape(3) + t1 = np.asarray(t1, dtype=np.float64).reshape(3) + t_interp = (1.0 - a) * t0 + a * t1 + return q_interp, t_interp def ray_from_pixel( @@ -174,7 +207,20 @@ def ray_from_pixel( direction : ndarray, shape (3,) Unit ray direction in world coords, pointing into the scene. """ - raise NotImplementedError("lane B (M8): implement ray_from_pixel") + R = quat_to_mat(q) # world -> cam + t_vec = np.asarray(t, dtype=np.float64).reshape(3) + origin = -R.T @ t_vec # camera center in world coords + + # Pinhole back-projection. A world point X projects with x_cam = R X + t and + # px = fx * x_cam.x / x_cam.z + cx, py = fy * x_cam.y / x_cam.z + cy + # (COLMAP camera axes: +x right, +y down, +z forward). So the camera-space direction + # through pixel (px, py) is [(px-cx)/fx, (py-cy)/fy, 1], pointing forward into the scene. + d_cam = np.array([(px - cx) / fx, (py - cy) / fy, 1.0], dtype=np.float64) + d_world = R.T @ d_cam # rotate direction cam -> world + norm = np.linalg.norm(d_world) + if norm < 1e-12: + raise ValueError("degenerate ray direction") + return origin, d_world / norm def triangulate_rays( @@ -194,7 +240,37 @@ def triangulate_rays( reject the triangulation when the rays are near-parallel or ``gap`` exceeds the spec threshold (0.5 scene units). """ - raise NotImplementedError("lane B (M8): implement triangulate_rays") + oa = np.asarray(origin_a, dtype=np.float64).reshape(3) + ob = np.asarray(origin_b, dtype=np.float64).reshape(3) + da = np.asarray(dir_a, dtype=np.float64).reshape(3) + db = np.asarray(dir_b, dtype=np.float64).reshape(3) + na, nb = np.linalg.norm(da), np.linalg.norm(db) + if na < 1e-12 or nb < 1e-12: + raise ValueError("triangulate_rays: zero-length direction") + da = da / na + db = db / nb + + # Shortest segment between two lines P(s)=oa+s*da, Q(u)=ob+u*db. Minimize |P-Q|^2. + # With unit directions: b = da.db, denom = 1 - b^2 (0 when parallel). + w0 = oa - ob + b = float(np.dot(da, db)) + d = float(np.dot(da, w0)) + e = float(np.dot(db, w0)) + denom = 1.0 - b * b + if denom < 1e-9: + # Near-parallel: no unique closest pair. Anchor on oa, take the closest point on + # line b to it; gap is the line-to-line perpendicular distance. Callers reject + # near-parallel rays up front, so this branch just stays numerically safe. + s = 0.0 + u = e + else: + s = (b * e - d) / denom + u = (e - b * d) / denom + pa = oa + s * da + pb = ob + u * db + point = 0.5 * (pa + pb) + gap = float(np.linalg.norm(pa - pb)) + return point, gap def nearest_point_on_ray( @@ -214,4 +290,26 @@ def nearest_point_on_ray( point : ndarray shape (3,) or None The selected point-cloud point, or ``None`` if none lie within ``radius``. """ - raise NotImplementedError("lane B (M8): implement nearest_point_on_ray") + o = np.asarray(origin, dtype=np.float64).reshape(3) + d = np.asarray(direction, dtype=np.float64).reshape(3) + nd = np.linalg.norm(d) + if nd < 1e-12: + raise ValueError("nearest_point_on_ray: zero-length direction") + d = d / nd + + pts = np.asarray(points, dtype=np.float64).reshape(-1, 3) + if len(pts) == 0: + return None + + v = pts - o # origin -> each point + proj = v @ d # signed distance along the ray + perp = v - np.outer(proj, d) # component perpendicular to the ray + perp_dist = np.linalg.norm(perp, axis=1) + + # In front of the origin and inside the cylinder of the given radius. + mask = (proj > 0.0) & (perp_dist <= radius) + if not np.any(mask): + return None + idx_in = np.where(mask)[0] + best = idx_in[np.argmin(perp_dist[idx_in])] + return pts[best].copy() diff --git a/backend/festival4d/sfm.py b/backend/festival4d/sfm.py index dc5f902..da9a2a1 100644 --- a/backend/festival4d/sfm.py +++ b/backend/festival4d/sfm.py @@ -1,33 +1,61 @@ """COLMAP orchestration, model parsing, scene normalization, pose export (spec M2, lane B). -STUB — lane B fills the bodies; the public signatures below are frozen (``cli.py`` calls -``run_reconstruct``). - -Pipeline (spec M2, for lane B's reference): +Pipeline (spec M2): - Sample frames (``frames.sample_frames``), then drive the COLMAP CLI via subprocess: - ``feature_extractor`` (OPENCV model, single camera per folder) -> ``sequential_matcher`` - within each video + ``exhaustive_matcher`` across -> ``mapper`` -> ``model_converter`` to - a TXT model. Parse ``images.txt`` / ``cameras.txt`` / ``points3D.txt`` by hand (no - pycolmap). + ``feature_extractor`` (OPENCV model, single camera per folder) -> ``exhaustive_matcher`` + (small frame counts; exhaustive covers within- and cross-video pairs) -> ``mapper`` -> + ``image_undistorter`` + ``model_converter`` to a PINHOLE TXT model. Parse ``images.txt`` / + ``cameras.txt`` / ``points3D.txt`` by hand (no pycolmap). - **Normalize** the scene: centroid -> origin, camera bounding sphere radius -> 10, world-up - -> average camera up. Apply the same similarity transform to poses and points. + -> average camera up. The same similarity transform is applied to poses and points. - **Interpolate** poses for unregistered frames (slerp rotation, lerp translation, mark - ``registered=False``); never extrapolate past the first/last registered frame. -- Export ``config.POINTS_PLY`` (binary little-endian, matching ``synthetic.write_ply``'s - format) and fill ``camera_poses`` via ``db.set_poses``. -- **Failure handling** (spec M2): COLMAP optional at runtime (``shutil.which("colmap")``); - if it registers < 60% of frames or < 2 videos into one model, print a clear diagnostic and - leave existing poses untouched — never corrupt the DB. The app degrades to "synced videos, - no 3D". + ``registered=False``); never extrapolate past the first/last registered frame of a video. +- Export ``config.POINTS_PLY`` (binary little-endian, matching ``synthetic.write_ply``) and + fill ``camera_poses`` via ``db.set_poses``. +- **Failure handling** (spec M2): COLMAP optional at runtime (``shutil.which("colmap")``); if + it registers < 60% of sampled frames or < 2 videos into one model, print a clear diagnostic + and leave existing poses untouched — never corrupt the DB. The app degrades to "synced + videos, no 3D". """ from __future__ import annotations import logging +import functools +import struct +import subprocess from pathlib import Path +import numpy as np + +from festival4d import config, db, frames +from festival4d.geometry import mat_to_quat, quat_to_mat, slerp_pose + log = logging.getLogger("festival4d.sfm") +# Minimum share of sampled frames COLMAP must register, and minimum videos in one model, +# for the reconstruction to be trusted (spec M2 failure handling). +MIN_REGISTERED_FRACTION = 0.60 +MIN_VIDEOS_IN_MODEL = 2 +NORMALIZED_SPHERE_RADIUS = 10.0 + +# How COLMAP camera models lay out PARAMS -> (fx, fy, cx, cy). Distortion terms are ignored: +# the stored intrinsics are the pinhole part (exact after image_undistorter's PINHOLE output, +# approximate if a distorted model is parsed directly). Covers the models COLMAP emits here. +_INTRINSICS_FROM_PARAMS = { + "SIMPLE_PINHOLE": lambda p: (p[0], p[0], p[1], p[2]), + "PINHOLE": lambda p: (p[0], p[1], p[2], p[3]), + "SIMPLE_RADIAL": lambda p: (p[0], p[0], p[1], p[2]), + "RADIAL": lambda p: (p[0], p[0], p[1], p[2]), + "SIMPLE_RADIAL_FISHEYE": lambda p: (p[0], p[0], p[1], p[2]), + "RADIAL_FISHEYE": lambda p: (p[0], p[0], p[1], p[2]), + "OPENCV": lambda p: (p[0], p[1], p[2], p[3]), + "OPENCV_FISHEYE": lambda p: (p[0], p[1], p[2], p[3]), + "FULL_OPENCV": lambda p: (p[0], p[1], p[2], p[3]), + "FOV": lambda p: (p[0], p[1], p[2], p[3]), + "THIN_PRISM_FISHEYE": lambda p: (p[0], p[1], p[2], p[3]), +} + def colmap_available() -> bool: """Whether the COLMAP binary is on PATH (spec M2 optional-at-runtime rule).""" @@ -36,29 +64,416 @@ def colmap_available() -> bool: return shutil.which("colmap") is not None +# --------------------------------------------------------------------------- +# COLMAP TXT-model parsers (hand-written; no pycolmap dependency) +# --------------------------------------------------------------------------- +def _iter_data_lines(path: Path): + """Yield non-empty, non-comment lines from a COLMAP TXT file.""" + with open(path, "r", encoding="utf-8") as f: + for line in f: + line = line.strip() + if line and not line.startswith("#"): + yield line + + def parse_images_txt(path: Path) -> list[dict]: - """Parse a COLMAP ``images.txt`` into per-image pose records (spec M2, lane B).""" - raise NotImplementedError("lane B (M2): implement parse_images_txt") + """Parse a COLMAP ``images.txt`` into per-image pose records (spec M2, lane B). + + Format: two lines per image. The first is + ``IMAGE_ID QW QX QY QZ TX TY TZ CAMERA_ID NAME``; the second is the 2D keypoint list + (ignored here). Returns one dict per image with the world->camera quaternion + ``[w,x,y,z]``, translation, camera id, and image name. + """ + records: list[dict] = [] + take_pose_line = True # image data alternates: pose line, then points2D line + for line in _iter_data_lines(path): + if take_pose_line: + parts = line.split() + if len(parts) < 10: + raise ValueError(f"malformed images.txt pose line: {line!r}") + records.append({ + "image_id": int(parts[0]), + "qw": float(parts[1]), "qx": float(parts[2]), + "qy": float(parts[3]), "qz": float(parts[4]), + "tx": float(parts[5]), "ty": float(parts[6]), "tz": float(parts[7]), + "camera_id": int(parts[8]), + "name": " ".join(parts[9:]), # NAME may itself contain spaces + }) + take_pose_line = not take_pose_line + return records def parse_cameras_txt(path: Path) -> dict[int, dict]: - """Parse a COLMAP ``cameras.txt`` into ``camera_id -> intrinsics`` (spec M2, lane B).""" - raise NotImplementedError("lane B (M2): implement parse_cameras_txt") + """Parse a COLMAP ``cameras.txt`` into ``camera_id -> intrinsics`` (spec M2, lane B). + + Format per line: ``CAMERA_ID MODEL WIDTH HEIGHT PARAMS[]``. Returns, per camera id, a + dict with ``model``, ``width``, ``height``, pinhole ``fx, fy, cx, cy`` and the raw + ``params`` list. + """ + cameras: dict[int, dict] = {} + for line in _iter_data_lines(path): + parts = line.split() + if len(parts) < 4: + raise ValueError(f"malformed cameras.txt line: {line!r}") + cam_id = int(parts[0]) + model = parts[1] + width, height = int(parts[2]), int(parts[3]) + params = [float(x) for x in parts[4:]] + if model not in _INTRINSICS_FROM_PARAMS: + raise ValueError(f"unsupported COLMAP camera model {model!r}") + try: + fx, fy, cx, cy = _INTRINSICS_FROM_PARAMS[model](params) + except IndexError: + raise ValueError( + f"cameras.txt: model {model} has too few PARAMS: {line!r}" + ) from None + cameras[cam_id] = { + "model": model, "width": width, "height": height, + "fx": float(fx), "fy": float(fy), "cx": float(cx), "cy": float(cy), + "params": params, + } + return cameras -def parse_points3d_txt(path: Path): - """Parse a COLMAP ``points3D.txt`` into ``(points Nx3, colors Nx3)`` (spec M2, lane B).""" - raise NotImplementedError("lane B (M2): implement parse_points3d_txt") +def parse_points3d_txt(path: Path) -> tuple[np.ndarray, np.ndarray]: + """Parse a COLMAP ``points3D.txt`` into ``(points Nx3 float32, colors Nx3 uint8)``. + + Format per line: ``POINT3D_ID X Y Z R G B ERROR TRACK[]`` (track ignored). + """ + xyz: list[tuple[float, float, float]] = [] + rgb: list[tuple[int, int, int]] = [] + for line in _iter_data_lines(path): + parts = line.split() + if len(parts) < 7: + raise ValueError(f"malformed points3D.txt line: {line!r}") + xyz.append((float(parts[1]), float(parts[2]), float(parts[3]))) + rgb.append((int(parts[4]), int(parts[5]), int(parts[6]))) + if not xyz: + return (np.zeros((0, 3), dtype=np.float32), np.zeros((0, 3), dtype=np.uint8)) + return (np.asarray(xyz, dtype=np.float32), np.asarray(rgb, dtype=np.uint8)) + + +# --------------------------------------------------------------------------- +# Scene normalization +# --------------------------------------------------------------------------- +def _rotation_aligning(a: np.ndarray, b: np.ndarray) -> np.ndarray: + """Shortest-arc rotation matrix that maps unit vector ``a`` onto unit vector ``b``.""" + a = a / np.linalg.norm(a) + b = b / np.linalg.norm(b) + v = np.cross(a, b) + c = float(np.dot(a, b)) + s = float(np.linalg.norm(v)) + if s < 1e-12: + if c > 0: # already aligned + return np.eye(3) + # antiparallel: rotate 180 deg about any axis perpendicular to a + perp = np.array([1.0, 0.0, 0.0]) + if abs(a[0]) > 0.9: + perp = np.array([0.0, 1.0, 0.0]) + axis = np.cross(a, perp) + axis /= np.linalg.norm(axis) + return 2.0 * np.outer(axis, axis) - np.eye(3) + vx = np.array([[0.0, -v[2], v[1]], [v[2], 0.0, -v[0]], [-v[1], v[0], 0.0]]) + return np.eye(3) + vx + vx @ vx * ((1.0 - c) / (s * s)) + + +def _camera_center(pose: dict) -> np.ndarray: + """World-space camera center ``C = -R^T t`` for a world->camera pose dict.""" + R = quat_to_mat([pose["qw"], pose["qx"], pose["qy"], pose["qz"]]) + t = np.array([pose["tx"], pose["ty"], pose["tz"]], dtype=np.float64) + return -R.T @ t + + +def _camera_up(pose: dict) -> np.ndarray: + """World-space camera up vector. COLMAP camera +y is *down*, so up = R^T @ (0,-1,0).""" + R = quat_to_mat([pose["qw"], pose["qx"], pose["qy"], pose["qz"]]) + return R.T @ np.array([0.0, -1.0, 0.0]) def normalize_scene(points, poses): - """Similarity transform: centroid->origin, camera sphere radius->10, up->+Y (spec M2).""" - raise NotImplementedError("lane B (M2): implement normalize_scene") + """Similarity transform: centroid->origin, camera sphere radius->10, up->+Y (spec M2). + + ``points`` is an ``(N,3)`` array; ``poses`` a list of world->camera pose dicts (keys + ``qw,qx,qy,qz,tx,ty,tz`` plus any extras, which are preserved). The single similarity + transform ``X' = s * Rn @ (X - c)`` (translate to the point centroid, rotate the average + camera-up onto +Y, scale so the camera bounding sphere has radius 10) is applied to both + the points and the poses. + + Returns ``(points_new float32 (N,3), poses_new list[dict])``. + """ + pts = np.asarray(points, dtype=np.float64).reshape(-1, 3) + poses = list(poses) + + # 1. Translate to the point-cloud centroid. + centroid = pts.mean(axis=0) if len(pts) else np.zeros(3) + + # 2. Rotate so the average camera-up aligns with +Y (phones held roughly upright). + if poses: + ups = np.array([_camera_up(p) for p in poses]) + up_avg = ups.mean(axis=0) + if np.linalg.norm(up_avg) > 1e-9: + Rn = _rotation_aligning(up_avg, np.array([0.0, 1.0, 0.0])) + else: + Rn = np.eye(3) + else: + Rn = np.eye(3) + + # 3. Scale so the farthest camera center from the centroid sits at radius 10. + if poses: + centers = np.array([_camera_center(p) for p in poses]) + max_dist = float(np.max(np.linalg.norm(centers - centroid, axis=1))) + scale = NORMALIZED_SPHERE_RADIUS / max_dist if max_dist > 1e-9 else 1.0 + else: + scale = 1.0 + + # Apply X' = s * Rn @ (X - c) to points. + points_new = (scale * (pts - centroid) @ Rn.T).astype(np.float32) + + # Apply the matching transform to each world->camera pose: + # R' = R @ Rn^T , t' = s * (t + R @ c) (derived so projections are unchanged). + poses_new: list[dict] = [] + for p in poses: + R = quat_to_mat([p["qw"], p["qx"], p["qy"], p["qz"]]) + t = np.array([p["tx"], p["ty"], p["tz"]], dtype=np.float64) + R_new = R @ Rn.T + t_new = scale * (t + R @ centroid) + q_new = mat_to_quat(R_new) + out = dict(p) + out["qw"], out["qx"], out["qy"], out["qz"] = ( + float(q_new[0]), float(q_new[1]), float(q_new[2]), float(q_new[3])) + out["tx"], out["ty"], out["tz"] = float(t_new[0]), float(t_new[1]), float(t_new[2]) + poses_new.append(out) + + return points_new, poses_new -def run_colmap(frames_dir: Path, workspace: Path) -> Path: - """Drive the COLMAP CLI and return the exported TXT-model directory (spec M2, lane B).""" - raise NotImplementedError("lane B (M2): implement run_colmap") +# --------------------------------------------------------------------------- +# Pose interpolation +# --------------------------------------------------------------------------- +def interpolate_poses(registered: list[dict], sampled: list[dict]) -> list[dict]: + """Fill unregistered sampled frames of one video between registered ones (spec M2). + + ``registered`` are solved poses (with ``registered=True``); ``sampled`` are the frames + that were fed to COLMAP as ``{frame_idx, t_video_s, ...}``. For each sampled frame that + lies strictly between the first and last registered frame but wasn't solved, slerp the + rotation and lerp the translation of the two bracketing registered poses (intrinsics + inherited from the earlier neighbor), marking it ``registered=False``. Frames outside the + registered span are dropped (no extrapolation). Returns all poses sorted by ``t_video_s``. + """ + reg = sorted(registered, key=lambda p: p["t_video_s"]) + if not reg: + return [] + reg_frames = {int(p["frame_idx"]) for p in reg} + t_first, t_last = reg[0]["t_video_s"], reg[-1]["t_video_s"] + + out: list[dict] = list(reg) + for sf in sampled: + fidx = int(sf["frame_idx"]) + tv = float(sf["t_video_s"]) + if fidx in reg_frames: + continue + if tv <= t_first or tv >= t_last: + continue # no extrapolation beyond the registered span + # bracketing registered neighbors + lo = max((p for p in reg if p["t_video_s"] <= tv), key=lambda p: p["t_video_s"]) + hi = min((p for p in reg if p["t_video_s"] >= tv), key=lambda p: p["t_video_s"]) + span = hi["t_video_s"] - lo["t_video_s"] + alpha = 0.0 if span <= 1e-12 else (tv - lo["t_video_s"]) / span + q, t = slerp_pose( + [lo["qw"], lo["qx"], lo["qy"], lo["qz"]], [lo["tx"], lo["ty"], lo["tz"]], + [hi["qw"], hi["qx"], hi["qy"], hi["qz"]], [hi["tx"], hi["ty"], hi["tz"]], + alpha, + ) + out.append({ + "frame_idx": fidx, "t_video_s": tv, + "qw": float(q[0]), "qx": float(q[1]), "qy": float(q[2]), "qz": float(q[3]), + "tx": float(t[0]), "ty": float(t[1]), "tz": float(t[2]), + "fx": lo["fx"], "fy": lo["fy"], "cx": lo["cx"], "cy": lo["cy"], + "registered": False, + }) + out.sort(key=lambda p: p["t_video_s"]) + return out + + +# --------------------------------------------------------------------------- +# COLMAP CLI orchestration +# --------------------------------------------------------------------------- +@functools.lru_cache(maxsize=None) +def _colmap_help(command: str) -> str: + """Cached ``colmap --help`` text, for build-aware option detection. + + COLMAP prints its option list to stderr, so both streams are captured. + """ + try: + proc = subprocess.run( + ["colmap", command, "--help"], capture_output=True, text=True + ) + return (proc.stdout or "") + "\n" + (proc.stderr or "") + except Exception: + return "" + + +def _cpu_flag(command: str) -> list[str]: + """Force CPU SIFT for ``command`` if this COLMAP build exposes a ``use_gpu`` option. + + The option prefix differs by version: COLMAP 3.x uses ``SiftExtraction`` / + ``SiftMatching``, 4.x uses ``FeatureExtraction`` / ``FeatureMatching``. We read the option + name from ``--help`` so the pipeline runs headless on either (and on CPU-only builds that + omit the option entirely, we pass nothing). Returns e.g. ``["--FeatureExtraction.use_gpu", + "0"]`` or ``[]``. + """ + for line in _colmap_help(command).splitlines(): + line = line.strip() + if line.startswith("--") and ".use_gpu" in line: + return [line.split()[0], "0"] + return [] + + +def _run_colmap_step(args: list[str]) -> bool: + """Run one ``colmap `` step; return True on success, False (logged) on failure.""" + cmd = ["colmap", *args] + log.info("colmap: %s", " ".join(args[:2])) + try: + proc = subprocess.run(cmd, capture_output=True, text=True) + except FileNotFoundError: + log.error("colmap binary vanished from PATH") + return False + if proc.returncode != 0: + log.error("colmap %s failed (rc=%d): %s", args[0], proc.returncode, + (proc.stderr or proc.stdout or "").strip()[-500:]) + return False + return True + + +def _registered_image_count(model_dir: Path) -> int | None: + """Number of registered images in a COLMAP model dir, or None if it has no image list. + + The mapper writes binary models by default; ``images.bin`` begins with a little-endian + ``uint64`` giving the registered-image count, which is the exact metric we want (file + size would instead track total keypoint observations and mis-rank multi-component runs). + Falls back to counting the TXT model's two-lines-per-image list. + """ + images_txt = model_dir / "images.txt" + images_bin = model_dir / "images.bin" + if images_txt.exists(): + return sum(1 for _ in _iter_data_lines(images_txt)) // 2 + if images_bin.exists(): + with open(images_bin, "rb") as f: + header = f.read(8) + return struct.unpack(" Path | None: + """Return the reconstruction subdir (``0``, ``1``, ...) with the most registered images.""" + candidates = [d for d in sparse_dir.iterdir() if d.is_dir()] if sparse_dir.exists() else [] + best, best_n = None, -1 + for d in candidates: + n = _registered_image_count(d) + if n is None: + continue + if n > best_n: + best, best_n = d, n + return best + + +def run_colmap(frames_dir: Path, workspace: Path) -> Path | None: + """Drive the COLMAP CLI and return the exported TXT-model directory (spec M2, lane B). + + ``frames_dir`` holds one subfolder of JPEGs per video (single camera per folder). Runs + feature extraction (CPU SIFT, OPENCV model), exhaustive matching, incremental mapping, + then undistortion + TXT conversion. Returns the directory containing ``images.txt`` / + ``cameras.txt`` / ``points3D.txt``, or ``None`` if COLMAP produced no usable model + (caller degrades gracefully). + """ + workspace = Path(workspace) + # Start each run from a clean workspace: a stale database.db makes feature_extractor + # fail ("images already exist"), and a leftover sparse/ model would poison the pick. + if workspace.exists(): + import shutil + shutil.rmtree(workspace) + workspace.mkdir(parents=True, exist_ok=True) + database = workspace / "database.db" + sparse = workspace / "sparse" + sparse.mkdir(exist_ok=True) + + ok = _run_colmap_step([ + "feature_extractor", + "--database_path", str(database), + "--image_path", str(frames_dir), + "--ImageReader.single_camera_per_folder", "1", + "--ImageReader.camera_model", "OPENCV", + *_cpu_flag("feature_extractor"), + ]) + if not ok: + return None + + # Exhaustive matching covers within- and cross-video pairs; frame counts are small. + if not _run_colmap_step([ + "exhaustive_matcher", + "--database_path", str(database), + *_cpu_flag("exhaustive_matcher"), + ]): + return None + + if not _run_colmap_step([ + "mapper", + "--database_path", str(database), + "--image_path", str(frames_dir), + "--output_path", str(sparse), + ]): + return None + + model = _largest_model_dir(sparse) + if model is None: + log.error("colmap mapper produced no reconstruction") + return None + + # Undistort -> PINHOLE model, then convert to TXT. Fall back to converting the raw + # (possibly distorted) sparse model directly if undistortion fails. + txt_dir = workspace / "model_txt" + txt_dir.mkdir(exist_ok=True) + dense = workspace / "dense" + undistorted = _run_colmap_step([ + "image_undistorter", + "--image_path", str(frames_dir), + "--input_path", str(model), + "--output_path", str(dense), + "--output_type", "COLMAP", + ]) + convert_input = (dense / "sparse") if undistorted and (dense / "sparse").exists() else model + if not _run_colmap_step([ + "model_converter", + "--input_path", str(convert_input), + "--output_path", str(txt_dir), + "--output_type", "TXT", + ]): + return None + if not (txt_dir / "images.txt").exists(): + log.error("colmap model_converter produced no images.txt") + return None + return txt_dir + + +# --------------------------------------------------------------------------- +# Full reconstruction entrypoint +# --------------------------------------------------------------------------- +def _video_id_and_frame(name: str) -> tuple[int, int] | None: + """Recover ``(video_id, frame_idx)`` from a ``{video_id}_{frame_idx}.jpg`` image name.""" + stem = Path(name).stem + if "_" not in stem: + return None + vid_s, frame_s = stem.rsplit("_", 1) + try: + return int(vid_s), int(frame_s) + except ValueError: + return None + + +def _diagnostic(msg: str) -> None: + """Emit a clear operator-facing diagnostic (spec M2 graceful degradation).""" + log.warning("reconstruct: %s", msg) + print(f"[reconstruct] {msg}") def run_reconstruct() -> dict: @@ -66,6 +481,128 @@ def run_reconstruct() -> dict: Entrypoint for ``python -m festival4d reconstruct``. Degrades gracefully when COLMAP is absent or the reconstruction is too weak (spec M2 failure handling): existing poses are - left untouched. Returns a summary dict. + left untouched and the DB is never corrupted. Returns a summary dict with a ``status``. """ - raise NotImplementedError("lane B (M2): implement run_reconstruct") + db.init_engine() + db.init_db() + videos = db.get_videos() + if not videos: + _diagnostic("no videos in the project — run `synthetic` or `ingest` first.") + return {"status": "no_videos"} + + if not colmap_available(): + _diagnostic( + "COLMAP not found on PATH — skipping 3D reconstruction. Install it " + "(`brew install colmap` on macOS) to enable pose/point-cloud solving. " + "Existing (synthetic/previous) poses are left untouched; the app still runs as " + "synced videos without 3D." + ) + return {"status": "skipped_no_colmap", "videos": len(videos)} + + # 1. Sample sharp frames per video into one subfolder each (single camera per folder). + frames_root = config.FRAMES_DIR + per_video_sampled: dict[int, list[dict]] = {} + for v in videos: + video_path = config.RAW_DIR / v.filename + if not video_path.exists(): + _diagnostic(f"video file missing, skipping: {video_path}") + continue + out_dir = frames_root / str(v.id) + paths = frames.sample_frames(video_path, v.id, out_dir) + fps = float(v.fps) if v.fps and v.fps > 0 else 30.0 + sampled = [] + for p in paths: + parsed = _video_id_and_frame(p.name) + if parsed is None: + continue + _, fidx = parsed + sampled.append({"frame_idx": fidx, "t_video_s": fidx / fps}) + per_video_sampled[v.id] = sorted(sampled, key=lambda s: s["frame_idx"]) + + total_sampled = sum(len(s) for s in per_video_sampled.values()) + if total_sampled == 0: + _diagnostic("no frames could be sampled from any video — nothing to reconstruct.") + return {"status": "no_frames"} + + # 2. Run COLMAP. + model_dir = run_colmap(frames_root, config.COLMAP_DIR) + if model_dir is None: + _diagnostic( + "COLMAP produced no usable reconstruction. Likely causes: too-dark footage, " + "motion blur, or insufficient view overlap between cameras. Existing poses left " + "untouched — the app degrades to synced videos, no 3D." + ) + return {"status": "failed_no_model", "sampled": total_sampled} + + # 3. Parse the model. + images = parse_images_txt(model_dir / "images.txt") + cameras = parse_cameras_txt(model_dir / "cameras.txt") + points, colors = parse_points3d_txt(model_dir / "points3D.txt") + + # 4. Attach each registered image to its video + intrinsics. + registered_by_video: dict[int, list[dict]] = {} + for img in images: + parsed = _video_id_and_frame(img["name"]) + if parsed is None: + continue + video_id, frame_idx = parsed + cam = cameras.get(img["camera_id"]) + if cam is None: + continue + fps = next((float(v.fps) for v in videos if v.id == video_id), 30.0) or 30.0 + registered_by_video.setdefault(video_id, []).append({ + "frame_idx": frame_idx, "t_video_s": frame_idx / fps, + "qw": img["qw"], "qx": img["qx"], "qy": img["qy"], "qz": img["qz"], + "tx": img["tx"], "ty": img["ty"], "tz": img["tz"], + "fx": cam["fx"], "fy": cam["fy"], "cx": cam["cx"], "cy": cam["cy"], + "registered": True, "video_id": video_id, + }) + + registered_count = sum(len(v) for v in registered_by_video.values()) + videos_in_model = len(registered_by_video) + frac = registered_count / total_sampled if total_sampled else 0.0 + if frac < MIN_REGISTERED_FRACTION or videos_in_model < MIN_VIDEOS_IN_MODEL: + _diagnostic( + f"weak reconstruction: {registered_count}/{total_sampled} frames " + f"({frac:.0%}) registered across {videos_in_model} video(s); need " + f">= {MIN_REGISTERED_FRACTION:.0%} of frames and >= {MIN_VIDEOS_IN_MODEL} " + "videos. Likely causes: too-dark footage, motion blur, or too little view " + "overlap. Existing poses left untouched — degrading to synced videos, no 3D." + ) + return { + "status": "failed_weak", "sampled": total_sampled, + "registered": registered_count, "videos_in_model": videos_in_model, + } + + # 5. Normalize the whole scene (points + all registered poses) with one transform. + all_registered = [p for v in registered_by_video.values() for p in v] + norm_points, norm_registered = normalize_scene(points, all_registered) + norm_by_video: dict[int, list[dict]] = {} + for p in norm_registered: + norm_by_video.setdefault(p["video_id"], []).append(p) + + # 6. Interpolate unregistered sampled frames, then write poses per video (atomic replace). + interpolated_total = 0 + for video_id, reg in norm_by_video.items(): + sampled = per_video_sampled.get(video_id, []) + final = interpolate_poses(reg, sampled) + interpolated_total += sum(1 for p in final if not p.get("registered", True)) + db.set_poses(video_id, final) + + # 7. Export the normalized point cloud in the frozen PLY format. + from festival4d.synthetic import write_ply + write_ply(config.POINTS_PLY, norm_points, colors) + + summary = { + "status": "ok", + "videos_in_model": videos_in_model, + "registered": registered_count, + "interpolated": interpolated_total, + "points": int(len(norm_points)), + "points_ply": str(config.POINTS_PLY), + } + log.info("reconstruct: done — %s", summary) + print(f"[reconstruct] reconstructed {videos_in_model} videos, " + f"{registered_count} registered + {interpolated_total} interpolated poses, " + f"{len(norm_points)} points -> {config.POINTS_PLY}") + return summary diff --git a/backend/tests/test_geometry.py b/backend/tests/test_geometry.py index a5aaadc..0ca0f30 100644 --- a/backend/tests/test_geometry.py +++ b/backend/tests/test_geometry.py @@ -108,3 +108,226 @@ def test_quat_to_mat_matches_scipy_oracle(): np.testing.assert_allclose( quat_to_mat(q_wxyz), Rotation.from_quat(q_xyzw).as_matrix(), atol=1e-9 ) + + +# =========================================================================== +# Lane B (spec M2 + M8) tests — below the frozen block. These exercise the +# stubs foundation left: slerp_pose, ray_from_pixel, triangulate_rays, +# nearest_point_on_ray. They do NOT touch the frozen colmap_to_threejs cases above. +# =========================================================================== +from scipy.spatial.transform import Slerp # noqa: E402 + +from festival4d.geometry import ( # noqa: E402 + colmap_to_threejs, + nearest_point_on_ray, + ray_from_pixel, + slerp_pose, + triangulate_rays, +) + + +def _wxyz_to_xyzw(q): + q = np.asarray(q, dtype=float) + return np.array([q[1], q[2], q[3], q[0]]) + + +def _rand_unit_quat(rng): + q = rng.normal(size=4) + return q / np.linalg.norm(q) + + +# --- slerp_pose ------------------------------------------------------------ +def test_slerp_pose_endpoints(): + q0 = _rand_unit_quat(np.random.default_rng(1)) + q1 = _rand_unit_quat(np.random.default_rng(2)) + t0 = np.array([1.0, -2.0, 3.0]) + t1 = np.array([-4.0, 5.0, 6.0]) + + q_a, t_a = slerp_pose(q0, t0, q1, t1, 0.0) + q_b, t_b = slerp_pose(q0, t0, q1, t1, 1.0) + # endpoints recover the endpoint *rotations* (quaternion up to sign) and translations + np.testing.assert_allclose(quat_to_mat(q_a), quat_to_mat(q0), atol=1e-12) + np.testing.assert_allclose(quat_to_mat(q_b), quat_to_mat(q1), atol=1e-12) + np.testing.assert_allclose(t_a, t0, atol=1e-12) + np.testing.assert_allclose(t_b, t1, atol=1e-12) + + +def test_slerp_pose_matches_scipy_oracle(): + rng = np.random.default_rng(20260716) + for _ in range(30): + q0 = _rand_unit_quat(rng) + q1 = _rand_unit_quat(rng) + oracle = Slerp([0.0, 1.0], Rotation.from_quat( + [_wxyz_to_xyzw(q0), _wxyz_to_xyzw(q1)])) + for alpha in (0.1, 0.25, 0.5, 0.73, 0.9): + q, _ = slerp_pose(q0, [0, 0, 0], q1, [0, 0, 0], alpha) + np.testing.assert_allclose( + quat_to_mat(q), oracle(alpha).as_matrix(), atol=1e-9) + + +def test_slerp_pose_translation_is_linear(): + q = np.array([1.0, 0.0, 0.0, 0.0]) + t0 = np.array([0.0, 0.0, 0.0]) + t1 = np.array([10.0, -4.0, 2.0]) + for alpha in (0.0, 0.3, 0.5, 1.0): + _, t = slerp_pose(q, t0, q, t1, alpha) + np.testing.assert_allclose(t, (1 - alpha) * t0 + alpha * t1, atol=1e-12) + + +def test_slerp_pose_double_cover_takes_short_arc(): + """q1 and -q1 are the same rotation; slerp must yield the same (short-arc) result.""" + rng = np.random.default_rng(99) + q0 = _rand_unit_quat(rng) + q1 = _rand_unit_quat(rng) + for alpha in (0.2, 0.5, 0.8): + qa, _ = slerp_pose(q0, [0, 0, 0], q1, [0, 0, 0], alpha) + qb, _ = slerp_pose(q0, [0, 0, 0], -q1, [0, 0, 0], alpha) + np.testing.assert_allclose(quat_to_mat(qa), quat_to_mat(qb), atol=1e-12) + + +def test_slerp_pose_midpoint_is_half_angle(): + """A 180-degree-ish pair: the midpoint rotation angle is half the endpoint angle.""" + q0 = np.array([1.0, 0.0, 0.0, 0.0]) # identity + ang = np.radians(100.0) + q1 = np.array([np.cos(ang / 2), 0.0, np.sin(ang / 2), 0.0]) # yaw 100 deg about Y + qm, _ = slerp_pose(q0, [0, 0, 0], q1, [0, 0, 0], 0.5) + # rotation angle of qm relative to identity should be ~50 deg + angle = 2.0 * np.arccos(min(1.0, abs(qm[0]))) + assert abs(np.degrees(angle) - 50.0) < 1e-6 + + +# --- ray_from_pixel -------------------------------------------------------- +def test_ray_from_pixel_origin_is_camera_center(): + """The ray origin equals the Three.js camera center from the frozen contract.""" + rng = np.random.default_rng(5) + for _ in range(20): + q = _rand_unit_quat(rng) + t = rng.uniform(-5, 5, size=3) + position, _ = colmap_to_threejs(q, t) + origin, _ = ray_from_pixel(q, t, 600, 600, 320, 180, 320, 180) + np.testing.assert_allclose(origin, position, atol=1e-9) + + +def test_ray_from_pixel_center_is_forward_axis(): + """The center pixel unprojects along the camera forward (+z) axis, in world coords.""" + q = np.array([1.0, 0.0, 0.0, 0.0]) # identity world->cam + t = np.array([0.0, 0.0, -10.0]) # camera center at (0,0,10) + origin, direction = ray_from_pixel(q, t, 500, 500, 320, 180, 320, 180) + np.testing.assert_allclose(origin, [0.0, 0.0, 10.0], atol=1e-9) + np.testing.assert_allclose(direction, [0.0, 0.0, 1.0], atol=1e-9) + + +def _project_colmap(q, t, fx, fy, cx, cy, X): + """Forward COLMAP pinhole projection of a world point X -> pixel (px, py).""" + R = quat_to_mat(q) + xc = R @ np.asarray(X, float) + np.asarray(t, float) + return fx * xc[0] / xc[2] + cx, fy * xc[1] / xc[2] + cy, xc[2] + + +def test_ray_from_pixel_is_projection_inverse(): + """A world point in front of the camera lies exactly on the ray through its pixel.""" + rng = np.random.default_rng(7) + for _ in range(50): + q = _rand_unit_quat(rng) + t = rng.uniform(-3, 3, size=3) + fx = fy = rng.uniform(400, 800) + cx, cy = 320.0, 180.0 + C = -quat_to_mat(q).T @ t # camera center + forward = quat_to_mat(q).T @ np.array([0, 0, 1.0]) + X = C + rng.uniform(2, 8) * forward + rng.uniform(-1, 1, size=3) # in front + px, py, zc = _project_colmap(q, t, fx, fy, cx, cy, X) + if zc <= 0.1: + continue + origin, direction = ray_from_pixel(q, t, fx, fy, cx, cy, px, py) + to_X = X - origin + # X - origin must be parallel to direction and in front (positive projection) + cross = np.cross(to_X, direction) + assert np.linalg.norm(cross) < 1e-6 * (1 + np.linalg.norm(to_X)) + assert np.dot(to_X, direction) > 0 + + +# --- triangulate_rays ------------------------------------------------------ +def test_triangulate_rays_intersecting(): + P = np.array([1.0, 2.0, 3.0]) + oa = np.array([0.0, 0.0, 0.0]) + ob = np.array([4.0, 0.0, 0.0]) + point, gap = triangulate_rays(oa, P - oa, ob, P - ob) + np.testing.assert_allclose(point, P, atol=1e-9) + assert gap < 1e-9 + + +def test_triangulate_rays_skew_known_geometry(): + """Ray A along +x at z=0; ray B along +y at z=1. Closest points (0,0,0),(0,0,1).""" + oa = np.array([0.0, 0.0, 0.0]); da = np.array([1.0, 0.0, 0.0]) + ob = np.array([0.0, 0.0, 1.0]); db = np.array([0.0, 1.0, 0.0]) + point, gap = triangulate_rays(oa, da, ob, db) + np.testing.assert_allclose(point, [0.0, 0.0, 0.5], atol=1e-9) + assert abs(gap - 1.0) < 1e-9 + + +def test_triangulate_rays_parallel_is_safe(): + """Parallel rays must not divide-by-zero; gap ~ their separation.""" + oa = np.array([0.0, 0.0, 0.0]); da = np.array([1.0, 0.0, 0.0]) + ob = np.array([0.0, 2.0, 0.0]); db = np.array([1.0, 0.0, 0.0]) + point, gap = triangulate_rays(oa, da, ob, db) + assert np.all(np.isfinite(point)) + assert abs(gap - 2.0) < 1e-6 + + +def test_triangulate_rays_from_two_cameras_recovers_point(): + """Two cameras looking at a world point; triangulating their pixel rays recovers it.""" + rng = np.random.default_rng(123) + P = np.array([0.5, 1.0, -0.5]) + fx = fy = 600.0; cx, cy = 320.0, 180.0 + rays = [] + for center in ([6.0, 2.0, 6.0], [-6.0, 2.5, 6.0]): + C = np.array(center, float) + z = P - C; z /= np.linalg.norm(z) + up = np.array([0.0, 1.0, 0.0]) + up = up - np.dot(up, z) * z; up /= np.linalg.norm(up) + x = np.cross(z, up); y = -up + R_c2w = np.column_stack([x, y, z]); R = R_c2w.T + t = -R @ C + q = mat_to_quat(R) + px, py, _ = _project_colmap(q, t, fx, fy, cx, cy, P) + rays.append(ray_from_pixel(q, t, fx, fy, cx, cy, px, py)) + point, gap = triangulate_rays(rays[0][0], rays[0][1], rays[1][0], rays[1][1]) + assert np.linalg.norm(point - P) < 0.2 # spec M8 acceptance tolerance + assert gap < 1e-6 + + +# --- nearest_point_on_ray -------------------------------------------------- +def test_nearest_point_on_ray_hits_within_radius(): + o = np.array([0.0, 0.0, 0.0]); d = np.array([0.0, 0.0, 1.0]) + points = np.array([ + [0.1, 0.0, 5.0], # perp 0.1, in front -> candidate + [0.05, 0.0, 3.0], # perp 0.05, in front -> closest + [2.0, 0.0, 5.0], # perp 2.0 -> outside radius + ]) + got = nearest_point_on_ray(o, d, points, radius=0.3) + np.testing.assert_allclose(got, [0.05, 0.0, 3.0], atol=1e-12) + + +def test_nearest_point_on_ray_excludes_outside_and_behind(): + o = np.array([0.0, 0.0, 0.0]); d = np.array([0.0, 0.0, 1.0]) + behind = np.array([[0.05, 0.0, -3.0]]) # in cylinder but behind origin + outside = np.array([[1.0, 0.0, 3.0]]) # in front but outside radius + assert nearest_point_on_ray(o, d, behind, radius=0.3) is None + assert nearest_point_on_ray(o, d, outside, radius=0.3) is None + assert nearest_point_on_ray(o, d, np.zeros((0, 3)), radius=0.3) is None + + +def test_nearest_point_on_ray_synthetic_stage_corner(): + """M8 single-view fallback: a ray toward a stage corner selects that corner point.""" + from festival4d import synthetic + + corner = np.array(synthetic.STAGE_CORNERS["Stage FL"], float) + points, _ = synthetic.generate_point_cloud() + # a camera looking straight at the corner; ray through the image center hits it + q, t = synthetic.look_at_colmap((5.0, 3.0, 9.0), corner) + intr = synthetic.intrinsics(640, 360) + origin, direction = ray_from_pixel( + q, t, intr["fx"], intr["fy"], intr["cx"], intr["cy"], intr["cx"], intr["cy"]) + got = nearest_point_on_ray(origin, direction, points.astype(float), radius=0.3) + assert got is not None + np.testing.assert_allclose(got, corner, atol=1e-4) diff --git a/backend/tests/test_sfm.py b/backend/tests/test_sfm.py new file mode 100644 index 0000000..29b3c4f --- /dev/null +++ b/backend/tests/test_sfm.py @@ -0,0 +1,494 @@ +"""Tests for lane B: frames.py (sampling) + sfm.py (parsers, normalization, +interpolation, and the graceful-degradation reconstruct pipeline). + +Everything here is verifiable without COLMAP: parsers run on hand-written TXT snippets, the +export path runs on a hand-built (synthetic-pose) COLMAP model injected in place of a real +COLMAP run, and the failure paths are exercised by monkeypatching. Spec M2 acceptance: +"on the synthetic fixture (which skips COLMAP and injects known poses) the whole export path +runs; if COLMAP is installed, the pipeline runs end-to-end without crashing." +""" + +from __future__ import annotations + +import struct +from pathlib import Path + +import numpy as np +import pytest + +from festival4d import config, db, frames, sfm, synthetic +from festival4d.geometry import quat_to_mat + + +# =========================================================================== +# COLMAP TXT parsers (hand-written snippets) +# =========================================================================== +IMAGES_TXT = """\ +# Image list with two lines of data per image: +# IMAGE_ID, QW, QX, QY, QZ, TX, TY, TZ, CAMERA_ID, NAME +# POINTS2D[] as (X, Y, POINT3D_ID) +# Number of images: 2, mean observations per image: 2 +1 0.9998 0.0 0.02 0.0 -1.5 0.3 8.0 1 1/1_0.jpg +480.0 270.0 5 500.1 260.2 -1 +2 0.7071 0.0 0.7071 0.0 2.0 0.1 7.5 2 2/2_30.jpg +100.0 120.0 -1 300.0 200.0 5 +""" + +CAMERAS_TXT = """\ +# Camera list with one line of data per camera: +# CAMERA_ID, MODEL, WIDTH, HEIGHT, PARAMS[] +# Number of cameras: 2 +1 PINHOLE 640 360 500.0 500.0 320.0 180.0 +2 OPENCV 640 360 510.0 511.0 321.0 181.0 0.01 -0.02 0.0 0.0 +""" + +POINTS3D_TXT = """\ +# 3D point list with one line of data per point: +# POINT3D_ID, X, Y, Z, R, G, B, ERROR, TRACK[] as (IMAGE_ID, POINT2D_IDX) +# Number of points: 2, mean track length: 2 +5 1.0 2.0 3.0 200 100 50 0.5 1 0 2 1 +9 -1.0 0.0 4.0 10 20 30 0.8 1 1 +""" + + +def test_parse_images_txt(tmp_path): + p = tmp_path / "images.txt" + p.write_text(IMAGES_TXT) + recs = sfm.parse_images_txt(p) + assert len(recs) == 2 + assert recs[0]["image_id"] == 1 + assert recs[0]["camera_id"] == 1 + assert recs[0]["name"] == "1/1_0.jpg" + np.testing.assert_allclose( + [recs[0]["qw"], recs[0]["qx"], recs[0]["qy"], recs[0]["qz"]], + [0.9998, 0.0, 0.02, 0.0]) + np.testing.assert_allclose( + [recs[1]["tx"], recs[1]["ty"], recs[1]["tz"]], [2.0, 0.1, 7.5]) + + +def test_parse_cameras_txt(tmp_path): + p = tmp_path / "cameras.txt" + p.write_text(CAMERAS_TXT) + cams = sfm.parse_cameras_txt(p) + assert set(cams) == {1, 2} + assert cams[1]["model"] == "PINHOLE" + assert (cams[1]["fx"], cams[1]["fy"], cams[1]["cx"], cams[1]["cy"]) == \ + (500.0, 500.0, 320.0, 180.0) + # OPENCV: first 4 params are fx, fy, cx, cy; distortion ignored + assert (cams[2]["fx"], cams[2]["fy"], cams[2]["cx"], cams[2]["cy"]) == \ + (510.0, 511.0, 321.0, 181.0) + + +def test_parse_points3d_txt(tmp_path): + p = tmp_path / "points3D.txt" + p.write_text(POINTS3D_TXT) + points, colors = sfm.parse_points3d_txt(p) + assert points.shape == (2, 3) and colors.shape == (2, 3) + np.testing.assert_allclose(points[0], [1.0, 2.0, 3.0]) + assert tuple(colors[0]) == (200, 100, 50) + assert points.dtype == np.float32 and colors.dtype == np.uint8 + + +def test_parse_points3d_empty(tmp_path): + p = tmp_path / "points3D.txt" + p.write_text("# header only\n") + points, colors = sfm.parse_points3d_txt(p) + assert points.shape == (0, 3) and colors.shape == (0, 3) + + +# =========================================================================== +# Scene normalization +# =========================================================================== +def _all_synth_poses(): + poses = [] + for i in range(3): + poses.extend(synthetic.camera_track(i, 20.0, 30, 640, 360)) + return poses + + +def test_normalize_scene_invariants(): + points, _ = synthetic.generate_point_cloud() + poses = _all_synth_poses() + npoints, nposes = sfm.normalize_scene(points.astype(float), poses) + + # centroid at origin + assert np.linalg.norm(npoints.mean(axis=0)) < 1e-3 + # camera bounding sphere (from origin = point centroid) has radius 10 + centers = np.array([sfm._camera_center(p) for p in nposes]) + assert abs(float(np.max(np.linalg.norm(centers, axis=1))) - 10.0) < 1e-6 + # average camera-up aligns with +Y + ups = np.array([sfm._camera_up(p) for p in nposes]) + up_avg = ups.mean(axis=0) + up_avg /= np.linalg.norm(up_avg) + np.testing.assert_allclose(up_avg, [0.0, 1.0, 0.0], atol=1e-6) + + +def test_normalize_scene_preserves_projection(): + """A world point projects to the same pixel before and after normalization + (a similarity transform of the world must not change any image).""" + points, _ = synthetic.generate_point_cloud() + poses = _all_synth_poses() + intr = synthetic.intrinsics(640, 360) + for p in poses: + p.update(intr) + X = np.array([0.7, 0.9, 0.3]) # arbitrary world point + + def project(pose, Xw): + R = quat_to_mat([pose["qw"], pose["qx"], pose["qy"], pose["qz"]]) + t = np.array([pose["tx"], pose["ty"], pose["tz"]]) + xc = R @ Xw + t + return np.array([pose["fx"] * xc[0] / xc[2] + pose["cx"], + pose["fy"] * xc[1] / xc[2] + pose["cy"]]) + + npoints, nposes = sfm.normalize_scene(points.astype(float), poses) + # transform the same world point by X' = s*Rn@(X-c): recover s,Rn,c from the point map. + # Easier: the invariant is projection equality, so map X through the same transform used + # on the cloud by fitting it — instead, just check pixel equality using each pose pair. + # Reconstruct transform from three non-collinear cloud points is overkill; use the fact + # that projection is invariant, so compare original point's pixel to the transformed + # point's pixel where the transform is inferred from the point cloud centroid+scale. + c = points.astype(float).mean(axis=0) + centers = np.array([sfm._camera_center(p) for p in poses]) + scale = 10.0 / float(np.max(np.linalg.norm(centers - c, axis=1))) + ups = np.array([sfm._camera_up(p) for p in poses]) + up_avg = ups.mean(axis=0) + Rn = sfm._rotation_aligning(up_avg, np.array([0.0, 1.0, 0.0])) + Xn = scale * Rn @ (X - c) + + for p_old, p_new in zip(poses, nposes): + np.testing.assert_allclose(project(p_old, X), project(p_new, Xn), atol=1e-6) + + +# =========================================================================== +# Pose interpolation +# =========================================================================== +def _pose(frame_idx, t, q, t3, registered=True): + return {"frame_idx": frame_idx, "t_video_s": t, + "qw": q[0], "qx": q[1], "qy": q[2], "qz": q[3], + "tx": t3[0], "ty": t3[1], "tz": t3[2], + "fx": 500.0, "fy": 500.0, "cx": 320.0, "cy": 180.0, + "registered": registered, "video_id": 1} + + +def test_interpolate_poses_fills_gaps_no_extrapolation(): + q0 = [1.0, 0.0, 0.0, 0.0] + ang = np.radians(60.0) + q1 = [np.cos(ang / 2), 0.0, np.sin(ang / 2), 0.0] # 60 deg yaw about Y + registered = [_pose(0, 0.0, q0, [0, 0, 0]), _pose(60, 2.0, q1, [6, 0, 0])] + sampled = [{"frame_idx": f, "t_video_s": f / 30.0} for f in (0, 30, 60, 90)] + + out = sfm.interpolate_poses(registered, sampled) + by_frame = {p["frame_idx"]: p for p in out} + assert set(by_frame) == {0, 30, 60} # frame 90 dropped (beyond last registered) + assert by_frame[30]["registered"] is False # interpolated + assert by_frame[0]["registered"] is True and by_frame[60]["registered"] is True + # midpoint translation is the lerp; rotation is ~30 deg + np.testing.assert_allclose( + [by_frame[30]["tx"], by_frame[30]["ty"], by_frame[30]["tz"]], [3, 0, 0], atol=1e-9) + mid_angle = 2.0 * np.degrees(np.arccos(min(1.0, abs(by_frame[30]["qw"])))) + assert abs(mid_angle - 30.0) < 1e-6 + + +def test_interpolate_poses_empty_registered(): + assert sfm.interpolate_poses([], [{"frame_idx": 0, "t_video_s": 0.0}]) == [] + + +# =========================================================================== +# frames.py — sharpness + windowed sampling +# =========================================================================== +def test_sharpness_sharp_beats_blurred(): + import cv2 + + rng = np.random.default_rng(0) + noise = (rng.integers(0, 256, size=(120, 120), dtype=np.uint8)) # high-frequency + flat = np.full((120, 120), 128, dtype=np.uint8) # no edges + blurred = cv2.GaussianBlur(noise, (0, 0), sigmaX=5) + assert frames.sharpness(noise) > frames.sharpness(blurred) > frames.sharpness(flat) + # accepts BGR too + bgr = cv2.cvtColor(noise, cv2.COLOR_GRAY2BGR) + assert frames.sharpness(bgr) > 0 + + +def test_sample_frames_picks_sharpest_per_window(tmp_path): + import cv2 + + w = h = 128 + fps = 30.0 + n_frames = 30 # 1.0 s -> two 0.5 s windows + vid_path = tmp_path / "clip.mp4" + writer = cv2.VideoWriter(str(vid_path), cv2.VideoWriter_fourcc(*"mp4v"), fps, (w, h)) + if not writer.isOpened(): + pytest.skip("no usable VideoWriter codec in this environment") + checker = np.indices((h, w)).sum(axis=0) % 2 # fine checkerboard -> sharp + sharp = (checker * 255).astype(np.uint8) + flat = np.full((h, w), 128, dtype=np.uint8) + sharp_frames = {7, 22} + for i in range(n_frames): + g = sharp if i in sharp_frames else flat + writer.write(cv2.cvtColor(g, cv2.COLOR_GRAY2BGR)) + writer.release() + + out = frames.sample_frames(vid_path, video_id=5, out_dir=tmp_path / "frames") + got_frames = sorted(int(Path(p).stem.split("_")[1]) for p in out) + + # ground-truth: decode the written stream and take the per-window argmax ourselves + cap = cv2.VideoCapture(str(vid_path)) + sh, idx = [], 0 + while True: + ok, fr = cap.read() + if not ok: + break + sh.append(frames.sharpness(fr)); idx += 1 + cap.release() + bucket = {} + for i, s in enumerate(sh): + win = int((i / fps) / 0.5) + if win not in bucket or s > bucket[win][0]: + bucket[win] = (s, i) + expected = sorted(v[1] for v in bucket.values()) + + assert got_frames == expected + assert all(Path(p).name.startswith("5_") for p in out) + + +# =========================================================================== +# run_reconstruct — graceful degradation + full export path (no real COLMAP) +# =========================================================================== +def _snapshot_poses(): + return { + v.id: [(p.frame_idx, round(p.qw, 9), round(p.tx, 9), p.registered) + for p in db.get_poses(v.id)] + for v in db.get_videos() + } + + +def _write_colmap_model(model_dir: Path, records: list[dict], cameras: dict, + points, colors) -> None: + """Write a COLMAP TXT model (images/cameras/points3D) from pose records.""" + model_dir.mkdir(parents=True, exist_ok=True) + with open(model_dir / "images.txt", "w") as f: + f.write("# image list\n") + for i, r in enumerate(records, start=1): + f.write(f"{i} {r['qw']} {r['qx']} {r['qy']} {r['qz']} " + f"{r['tx']} {r['ty']} {r['tz']} {r['camera_id']} {r['name']}\n") + f.write("100.0 100.0 -1\n") # dummy points2D line + with open(model_dir / "cameras.txt", "w") as f: + f.write("# camera list\n") + for cid, cam in cameras.items(): + f.write(f"{cid} PINHOLE {cam['width']} {cam['height']} " + f"{cam['fx']} {cam['fy']} {cam['cx']} {cam['cy']}\n") + with open(model_dir / "points3D.txt", "w") as f: + f.write("# point list\n") + for j, (p, c) in enumerate(zip(points, colors), start=1): + f.write(f"{j} {p[0]} {p[1]} {p[2]} {int(c[0])} {int(c[1])} {int(c[2])} 0.5 1 0\n") + + +def test_reconstruct_no_videos(): + db.init_engine() + db.reset_db() + result = sfm.run_reconstruct() + assert result["status"] == "no_videos" + + +def test_reconstruct_skips_without_colmap(monkeypatch): + synthetic.build(run_ffmpeg=False) + monkeypatch.setattr(sfm, "colmap_available", lambda: False) + before = _snapshot_poses() + result = sfm.run_reconstruct() + assert result["status"] == "skipped_no_colmap" + assert _snapshot_poses() == before # poses untouched + + +def _fake_sampling(monkeypatch, frames_per_video): + """Make sampling hermetic: placeholder raw files + fake sample_frames output.""" + for v in db.get_videos(): + (config.RAW_DIR / v.filename).parent.mkdir(parents=True, exist_ok=True) + (config.RAW_DIR / v.filename).write_bytes(b"stub") + + def fake_sample(video_path, video_id, out_dir, **kw): + return [Path(f"{video_id}_{f}.jpg") for f in frames_per_video] + + monkeypatch.setattr(sfm.frames, "sample_frames", fake_sample) + monkeypatch.setattr(sfm, "colmap_available", lambda: True) + + +def test_reconstruct_failed_no_model(monkeypatch): + synthetic.build(run_ffmpeg=False) + _fake_sampling(monkeypatch, [0, 30, 60]) + monkeypatch.setattr(sfm, "run_colmap", lambda frames_dir, workspace: None) + before = _snapshot_poses() + result = sfm.run_reconstruct() + assert result["status"] == "failed_no_model" + assert _snapshot_poses() == before # poses untouched + + +def test_reconstruct_failed_weak(monkeypatch, tmp_path): + synthetic.build(run_ffmpeg=False) + _fake_sampling(monkeypatch, [0, 30, 60]) # 3 videos x 3 = 9 sampled + # model registers only one image of one video -> ~11%, < 60% and < 2 videos + track = synthetic.camera_track(0, 20.0, 30, 640, 360) + rec = dict(track[0]); rec["camera_id"] = 1; rec["name"] = "1_0.jpg" + cams = {1: {"width": 640, "height": 360, **synthetic.intrinsics(640, 360)}} + pts, cols = synthetic.generate_point_cloud() + model_dir = tmp_path / "weak" + _write_colmap_model(model_dir, [rec], cams, pts[:50], cols[:50]) + monkeypatch.setattr(sfm, "run_colmap", lambda frames_dir, workspace: model_dir) + before = _snapshot_poses() + result = sfm.run_reconstruct() + assert result["status"] == "failed_weak" + assert _snapshot_poses() == before # poses untouched + + +def test_reconstruct_full_export_path(monkeypatch, tmp_path): + """The whole export path runs on injected known poses (spec M2 acceptance).""" + synthetic.build(run_ffmpeg=False) + videos = db.get_videos() + id_by_cam = {i: videos[i].id for i in range(len(videos))} # cam index -> db id + sampled_frames = [0, 15, 30, 45, 60] + _fake_sampling(monkeypatch, sampled_frames) + + # Build a model: videos for cam0/cam1 fully registered; cam2 registers a subset so + # interpolation must fill its gaps. All frame indices come from the synthetic tracks. + records, cameras = [], {} + reg_plan = {0: sampled_frames, 1: sampled_frames, 2: [0, 30, 60]} + for cam_idx, frame_list in reg_plan.items(): + vid = id_by_cam[cam_idx] + cameras[vid] = {"width": 640, "height": 360, **synthetic.intrinsics(640, 360)} + track = {p["frame_idx"]: p for p in synthetic.camera_track(cam_idx, 20.0, 30, 640, 360)} + for f in frame_list: + p = track[f] + records.append({**p, "camera_id": vid, "name": f"{vid}_{f}.jpg"}) + pts, cols = synthetic.generate_point_cloud() + model_dir = tmp_path / "good" + _write_colmap_model(model_dir, records, cameras, pts, cols) + monkeypatch.setattr(sfm, "run_colmap", lambda frames_dir, workspace: model_dir) + + result = sfm.run_reconstruct() + assert result["status"] == "ok" + assert result["registered"] == 13 # 5 + 5 + 3 + assert result["videos_in_model"] == 3 + assert result["interpolated"] == 2 # cam2 frames 15 and 45 + + # cam2 poses: 5 total (3 registered + 2 interpolated), sorted, with interp flags + cam2_poses = db.get_poses(id_by_cam[2]) + assert len(cam2_poses) == 5 + assert sum(1 for p in cam2_poses if not p.registered) == 2 + assert [p.frame_idx for p in cam2_poses] == [0, 15, 30, 45, 60] + + # normalized point cloud written in the frozen PLY format; camera sphere radius ~ 10 + assert config.POINTS_PLY.exists() + npoints, _ = synthetic.read_ply(config.POINTS_PLY) + assert len(npoints) == len(pts) + all_centers = [] + for v in videos: + for p in db.get_poses(v.id): + R = quat_to_mat([p.qw, p.qx, p.qy, p.qz]) + all_centers.append(-R.T @ np.array([p.tx, p.ty, p.tz])) + max_r = float(np.max(np.linalg.norm(np.array(all_centers), axis=1))) + assert abs(max_r - 10.0) < 0.5 # registered cams normalized to radius 10 + + +def test_colmap_model_parser_roundtrip(tmp_path): + """Writing then parsing a model recovers the poses/intrinsics/points.""" + track = synthetic.camera_track(0, 20.0, 30, 640, 360)[:3] + records = [{**p, "camera_id": 1, "name": f"1_{p['frame_idx']}.jpg"} for p in track] + cams = {1: {"width": 640, "height": 360, **synthetic.intrinsics(640, 360)}} + pts, cols = synthetic.generate_point_cloud() + model_dir = tmp_path / "rt" + _write_colmap_model(model_dir, records, cams, pts[:20], cols[:20]) + + images = sfm.parse_images_txt(model_dir / "images.txt") + cameras = sfm.parse_cameras_txt(model_dir / "cameras.txt") + ppoints, pcolors = sfm.parse_points3d_txt(model_dir / "points3D.txt") + assert len(images) == 3 + np.testing.assert_allclose(images[0]["qw"], track[0]["qw"], atol=1e-6) + assert cameras[1]["fx"] == pytest.approx(synthetic.intrinsics(640, 360)["fx"]) + assert len(ppoints) == 20 + + +def test_parse_cameras_txt_truncated_params(tmp_path): + """A camera line with too few PARAMS gives a clean ValueError, not a raw IndexError.""" + p = tmp_path / "cameras.txt" + p.write_text("# header\n7 PINHOLE 640 360 500.0 500.0 320.0\n") # PINHOLE needs 4 params + with pytest.raises(ValueError, match="too few PARAMS"): + sfm.parse_cameras_txt(p) + + +def test_largest_model_dir_picks_by_registered_count(tmp_path): + """Binary models are ranked by the images.bin header count, not file size.""" + sparse = tmp_path / "sparse" + # model 0: 10 registered images but a LARGE images.bin (many keypoint observations) + (sparse / "0").mkdir(parents=True) + (sparse / "0" / "images.bin").write_bytes(struct.pack(" dropped, no extrapolation + assert frames_out == {30, 45, 60} + by = {p["frame_idx"]: p for p in out} + assert by[45]["registered"] is False # interior gap still interpolated + + +def test_reconstruct_subset_leaves_others_untouched(monkeypatch, tmp_path): + """When only a subset of videos registers, the rest keep their existing poses (DB safety).""" + synthetic.build(run_ffmpeg=False) + videos = db.get_videos() + id_by_cam = {i: videos[i].id for i in range(len(videos))} + sampled_frames = [0, 15, 30, 45, 60] + _fake_sampling(monkeypatch, sampled_frames) + + # Register cam0 + cam1 fully (10/15 = 67% >= 60%, 2 videos -> ok); cam2 NOT in the model. + records, cameras = [], {} + for cam_idx in (0, 1): + vid = id_by_cam[cam_idx] + cameras[vid] = {"width": 640, "height": 360, **synthetic.intrinsics(640, 360)} + track = {p["frame_idx"]: p for p in synthetic.camera_track(cam_idx, 20.0, 30, 640, 360)} + for f in sampled_frames: + records.append({**track[f], "camera_id": vid, "name": f"{vid}_{f}.jpg"}) + pts, cols = synthetic.generate_point_cloud() + model_dir = tmp_path / "subset" + _write_colmap_model(model_dir, records, cameras, pts, cols) + monkeypatch.setattr(sfm, "run_colmap", lambda frames_dir, workspace: model_dir) + + cam2_before = [(p.frame_idx, p.qw, p.tx, p.registered) for p in db.get_poses(id_by_cam[2])] + result = sfm.run_reconstruct() + assert result["status"] == "ok" + assert result["videos_in_model"] == 2 + # cam2 absent from the model: its rows must be byte-identical afterward + cam2_after = [(p.frame_idx, p.qw, p.tx, p.registered) for p in db.get_poses(id_by_cam[2])] + assert cam2_after == cam2_before + assert len(cam2_before) == 41 # original synthetic poses intact + assert len(db.get_poses(id_by_cam[0])) == 5 # cam0 replaced with sampled frames + + +def test_reconstruct_no_frames_when_videos_missing(monkeypatch): + """All raw files absent -> videos skipped -> no_frames, poses untouched.""" + synthetic.build(run_ffmpeg=False) + for f in config.RAW_DIR.glob("*.mp4"): # clear stubs left by other tests + f.unlink() + monkeypatch.setattr(sfm, "colmap_available", lambda: True) + before = _snapshot_poses() + result = sfm.run_reconstruct() + assert result["status"] == "no_frames" + assert _snapshot_poses() == before + + +def test_reconstruct_no_frames_empty_sampling(monkeypatch): + """Raw files exist but sampling yields nothing -> no_frames, poses untouched.""" + synthetic.build(run_ffmpeg=False) + _fake_sampling(monkeypatch, []) # sample_frames returns [] + before = _snapshot_poses() + result = sfm.run_reconstruct() + assert result["status"] == "no_frames" + assert _snapshot_poses() == before diff --git a/plan/status/lane-B.md b/plan/status/lane-B.md new file mode 100644 index 0000000..3c3b380 --- /dev/null +++ b/plan/status/lane-B.md @@ -0,0 +1,59 @@ +# Status — lane-B + +## Round 0 — 2026-07-16 — STATUS: ready_to_merge + +**Directives acknowledged:** round 0 of plan/DIRECTIVES.md (run order C→A→D→B; merge when +acceptance passes, don't wait on siblings; maintain this status file; evidence discipline). + +**Environment:** COLMAP 4.1.0 (CPU-only) + ffmpeg both installed, so the pipeline was +exercised against **real COLMAP end-to-end**, not just the synthetic export path. + +**Isolation note (see blocker):** worked in a dedicated git worktree +(`/Users/m3ultra/Documents/festifun-laneB` on `lane/b-recon`) because all lane agents were +sharing one working directory. Only the 6 lane-B-owned files were touched. + +**Acceptance checklist** (spec M2 + M8 geometry + lane brief) — all pass: +- [x] `geometry.py` M8 stubs implemented: `slerp_pose` (shortest-arc, double-cover), `ray_from_pixel` + (COLMAP +y-down back-projection), `triangulate_rays` (closest-point + parallel guard), + `nearest_point_on_ray` (in-front radius cylinder). Frozen `colmap_to_threejs`/`quat_to_mat`/ + `mat_to_quat` and the 3 frozen POSE_TEST_VECTORS untouched (`git diff 5fa7301` shows 0 removed + lines in test_geometry.py; the frozen funcs are unchanged). +- [x] `frames.py`: `sharpness` (variance of Laplacian) + windowed sharpest-frame `sample_frames`. +- [x] COLMAP TXT parsers (images/cameras/points3D) with hand-written-snippet tests + a truncated-PARAMS + guard (clean ValueError). +- [x] `normalize_scene`: centroid≈0, camera-sphere radius≈10, up≈+Y — verified by `test_normalize_scene_invariants` + and a projection-invariance test (`test_normalize_scene_preserves_projection`). +- [x] pose interpolation (slerp+lerp, `registered=False`, NO extrapolation past first/last registered) — + both guard directions tested. +- [x] PLY export via frozen `synthetic.write_ply`; poses written via atomic `db.set_poses`. +- [x] graceful degradation: COLMAP absent / <60% frames / <2 videos / no frames / missing files → + diagnostic, existing poses UNTOUCHED, DB uncorrupted. Mixed-subset case (some videos registered, + others kept) tested for the DB-safety invariant. +- [x] real end-to-end: `synthetic` → `reconstruct` ran the full COLMAP CLI (feature_extractor → + exhaustive_matcher → mapper → image_undistorter → model_converter → my parsers), produced 6 + components (2/15/11/16/17/37 imgs), correctly selected the largest by `images.bin` header count, + judged 37/120 (31%) weak, degraded gracefully, **left all 41 synthetic poses/video untouched, exit 0.** +- [x] `pytest` green — **61 passed** (24 foundation + 37 lane-B; test_geometry.py + test_sfm.py). +- [x] no edits outside owned files (`git diff 5fa7301 --name-only` = frames/geometry/sfm/test_geometry only; + untracked: test_sfm.py, this status file); no new deps. + +**Bugs found & fixed while validating on real COLMAP:** +- COLMAP 4.x renamed `SiftExtraction`/`SiftMatching` → `FeatureExtraction`/`FeatureMatching`; option + detected from `--help` (which COLMAP prints to **stderr**) so the pipeline runs on 3.x and 4.x. +- `_largest_model_dir` ranked binary models by file size (tracks keypoints, not image count) → now reads + the exact registered-image count from the `images.bin` uint64 header (validated on the real 6-component run). +- workspace cleared before each run so a stale `database.db` can't fail reruns. + +**Adversarial review:** 5-lens review + verify pass (23 agents) → 18 raw findings, 5 CONFIRMED (2 code, 3 +test-coverage), all applied; 13 refuted (incl. triangulate/nearest-point "behind origin" — correct-as-written). + +**Blockers / questions for coordinator:** +- **Shared working tree collision.** All four lane agents are operating in the *same* working directory + (`/Users/m3ultra/Documents/festifun`), not separate clones/worktrees as `plan/README.md` requires. + Mid-session the shared branch was switched to `lane/d-events` and the tree accumulated uncommitted edits + from lanes A/B/C/D at once. I isolated lane B in its own worktree so nothing is clobbered. **Recommend the + coordinator have each lane use `git worktree`/separate clones, and adjudicate merges centrally**, since a + naive `git branch -f main` from any lane could silently drop a sibling's merge. I did NOT auto-merge to + `main` for this reason — lane B is ready and awaiting a safe merge. + +**Next:** commit lane B to `lane/b-recon`; coordinator to merge to `main` (safely, given the shared-tree issue).