AntonioJun/workspace

1"""Geometric assembler (runs in the pipeline venv). Imported by encoder/render.py.2 3Pure math/assembly -- builds the per-instance spatial code (object positions/sizes, pairwise4distances, closeness ranks, room outline, camera trajectory, appearance order) FROM5already-computed depth/pose/masks. Does NOT run DA3 or SAM3, and does not call cache_or_load()6-- that's run.py's job entirely (the only file that calls the actual model-inference functions).7Callers may provide raw depth/intr/c2w/conf/ftimes/per inputs or canonical world-space8geometry. Both paths emit the same compact spatial-code schema.9 10Formerly this file called into perceptual.py (as a dynamically-loaded `pl` module) for its own11geometry helpers -- build_instances, room_gravity, compute_floor_area, answer_closest_distance,12to_spatial_code, and everything else in this file below _room_outline(). Those functions are13now merged in directly, verbatim, since they were never DA3/SAM3 calls -- they're geometric14computations over already-extracted depth/masks, which is exactly this file's job.15perceptual.py's OTHER half (the actual model-calling functions) moved to run.py instead;16perceptual.py itself no longer exists.17 18The spatial code is the sole spatial representation the downstream VLM sees; no answer engine19is computed here (the harness runs the model).20"""21 22import os23import json24import numpy as np25import cv226 27# ==========================================================================================28# CONSTANTS -- geometry-cleanup knobs, merged in from perceptual.py. Env-tunable levers that29# affect only the math below (build_instances/backproject_frame/etc.), never model inference.30# ==========================================================================================31 32SCHEMA = [33    "x",34    "y",35    "z",36    "e1",37    "e2",38    "e3",39    "px",40    "pz",41    "size_median",42    "size_IQR",43    "time",44    "n",45]  # legacy row schema, see to_labeled()46KD_WORKERS = max(1, int(os.environ.get("VSI_KD_WORKERS", "4")))47# Batch workers set VSI_KD_WORKERS from their per-process CPU budget. Keeping a48# bounded default avoids pathological all-core thread spawning for small queries.49MIN_INSTANCE_PTS = (50    1  # bare geometry floor only: need >=1 valid depth pixel to place a 3D point.51)52# NO quality filtering / dedup -- count = exactly SAM3's tracked masklets (honest)53FLOOR_BAND = (54    0.15  # m above the floor to count as floor-level (LEGACY -- currently unreferenced;55)56# compute_floor_area now uses the RANSAC gravity plane directly instead)57# ---- geometry-cleanup levers (improve abs_distance etc.; env-tunable) ----58# DEPTH_COHERENCE (Tukey-fence bleed removal): SAFE + helpful everywhere -> default ON.59#   ablation: ARKit abs_distance 0.757->0.729 (no harm), ScanNet++ d755 0.443->0.486 (helps).60# CONF_PCT (per-frame confidence percentile): noisy-dataset ONLY -> default OFF.61#   helps ScanNet++ distance more (~0.56) but DESTROYS clean ARKit (0.757->0.429). Opt in via VSI_CONF_PCT=55.62CONF_PCT = float(os.environ.get("VSI_CONF_PCT", "0"))63DEPTH_COHERENCE = os.environ.get("VSI_DEPTH_COHERENCE", "1") == "1"64# cut mask-bleed at per-frame depth edges. Default OFF: it's a NO-OP on the dominant failure (same-depth bleed65# -- adjacent objects at similar range have no depth edge), and only helps depth-SEPARATED bleed. Enable per-need.66DEPTH_EDGE_REFINE = os.environ.get("VSI_DEPTH_EDGE_REFINE", "0") == "1"67# MASK_REFINE (appearance-guided boundary snap): uses the RGB color edge to clip same-depth mask bleed that68#   depth can't see. Principled + cheap (CPU, no model). Default OFF until validated; enable via VSI_MASK_REFINE=1.69MASK_REFINE = os.environ.get("VSI_MASK_REFINE", "0") == "1"70 71 72# ==========================================================================================73# ROOM/OBJECT GEOMETRY -- gravity, floor basis, per-instance spatial-code records. Merged in74# from perceptual.py, verbatim.75# ==========================================================================================76 77 78def to_labeled(objects, floor_area):79    """LEGACY per-class summary (kept for --format array-compat). Superseded by to_spatial_code."""80    out = {"objects": {}}81    for cls, row in objects.items():82        d = dict(zip(SCHEMA, row))83        out["objects"][cls] = {84            "count": int(d["n"]),85            "centroid_meters": {"x": d["x"], "y": d["y"], "z": d["z"]},86            "longest_dimension_meters": {87                "median": d["size_median"],88                "iqr": d["size_IQR"],89            },90            "centroid_spread": {91                "eigenvalues": [d["e1"], d["e2"], d["e3"]],92                "principal_axis_xz": [d["px"], d["pz"]],93            },94            "first_seen_seconds": d["time"],95        }96    out["room"] = {"floor_area_square_meters": floor_area}97    return out98 99 100def room_up_axis(instances, c2w):101    """up axis = smallest-extent axis of all object points; sign from gravity (floor->camera).102    Floor = densest horizontal slab; cameras are always above it, which fixes the sign.103    Returns (axis_index, signed_unit_vector)."""104    points = np.concatenate([i["pts"] for v in instances.values() for i in v], 0)105    ext = np.percentile(points, 98, 0) - np.percentile(points, 2, 0)106    up = int(np.argmin(ext))107    h, edges = np.histogram(points[:, up], bins=80)108    floor = 0.5 * (edges[h.argmax()] + edges[h.argmax() + 1])  # densest slab = floor109    cam_up = c2w[:, :3, 3][:, up].mean()110    e = np.zeros(3, np.float32)111    e[up] = 1.0 if cam_up > floor else -1.0112    return up, e113 114 115def room_gravity(116    depth, intr, c2w, conf, conf_pct=40, stride=12, fstride=15, iters=300, thr=0.05117):118    """Robust UP vector = normal of the RANSAC floor plane (a real physical plane), oriented toward119    the cameras. Works when the reconstruction is tilted/drifted (ScanNet++) where argmin-extent fails.120    Floor = large planar support with most non-inlier mass on ONE side. Returns (gravity_unit_vec, axis).121    """122    points = []123    for f in range(0, len(depth), fstride):124        height, width = depth[f].shape125        ys, xs = np.mgrid[0:height:stride, 0:width:stride]126        ys = ys.ravel()127        xs = xs.ravel()128        z = depth[f][ys, xs]129        ok = (z > 0) & np.isfinite(z)130        if conf is not None and conf_pct > 0:131            ok &= conf[f][ys, xs] >= np.percentile(conf[f], conf_pct)132        ys, xs, z = ys[ok], xs[ok], z[ok]133        if not len(z):134            continue135        intrinsics = intr[f]136        camera_points = np.stack(137            [138                (xs - intrinsics[0, 2]) * z / intrinsics[0, 0],139                (ys - intrinsics[1, 2]) * z / intrinsics[1, 1],140                z,141            ],142            1,143        )144        points.append((c2w[f][:3, :3] @ camera_points.T).T + c2w[f][:3, 3])145    points = np.concatenate(points).astype(np.float64) if points else np.zeros((0, 3))146    cam = c2w[:, :3, 3].mean(0)147    if len(points) < 100:  # fallback to axis-extent148        ext = (149            np.percentile(points, 98, 0) - np.percentile(points, 2, 0)150            if len(points)151            else np.ones(3)152        )153        ax = int(np.argmin(ext))154        g = np.zeros(3)155        g[ax] = 1.0156        return g, ax157    rng = np.random.RandomState(0)158    best = None159    best_score = -1160    for _ in range(iters):161        a, b, c = points[rng.choice(len(points), 3, False)]162        nrm = np.cross(b - a, c - a)163        ln = np.linalg.norm(nrm)164        if ln < 1e-6:165            continue166        nrm /= ln167        d = -nrm @ a168        side = points @ nrm + d169        ninl = int((np.abs(side) < thr).sum())170        if ninl < 50:171            continue172        score = ninl * max(173            np.mean(side > thr), np.mean(side < -thr)174        )  # big + one-sided = floor175        if score > best_score:176            best_score = score177            best = (nrm, d)178    if best is None:179        ext = np.percentile(points, 98, 0) - np.percentile(points, 2, 0)180        ax = int(np.argmin(ext))181        g = np.zeros(3)182        g[ax] = 1.0183        return g, ax184    nrm, d = best185    if (cam @ nrm + d) < 0:186        nrm = -nrm  # orient toward cameras (up)187    return nrm.astype(np.float32), int(np.argmax(np.abs(nrm)))188 189 190def pos3(rec):191    """Read either legacy numeric or current unit-string position formatting."""192    p = rec.get("position") or {}193 194    def meters(current, legacy):195        value = p.get(current, p.get(legacy, 0.0))196        if isinstance(value, str):197            if value.endswith(" meters"):198                value = value[: -len(" meters")]199            value = value.strip()200        return float(value)201 202    return [203        meters("x coordinate", "floor_x_meters"),204        meters("y coordinate", "floor_y_meters"),205        meters("height above floor", "height_above_floor_meters"),206    ]207 208 209def _floor_basis(up_vec):210    """Orthonormal floor basis (u, v horizontal; g = up) from the gravity vector. u is the OLD211    floor_x world axis projected onto the gravity plane (v the old floor_y axis), so this frame212    differs from the old axis-drop frame ONLY by the tilt correction -- NO arbitrary in-plane213    rotation (aligned scenes stay put; only tilt gets corrected). Objects / camera / room_outline214    all share this one gravity-plane frame; floor_area is rotation-invariant so it matches too.215    """216    g = np.asarray(up_vec, np.float64)217    g = g / (np.linalg.norm(g) + 1e-12)218    up_ax = int(np.argmax(np.abs(g)))219    floor_axes = [220        a for a in range(3) if a != up_ax221    ]  # the two world axes the old frame used222    e0 = np.zeros(3)223    e0[floor_axes[0]] = 1.0  # old floor_x world axis224    u = e0 - (e0 @ g) * g225    u = u / (np.linalg.norm(u) + 1e-12)  # project it into the gravity plane226    v = np.cross(g, u)227    if v[floor_axes[1]] < 0:228        v = -v  # keep floor_y sign aligned with the old axis229    return u, v, g230 231 232def _floor_level(points, gravity, v2=None):233    """Reference floor estimator shared by every model representation.234 235    The supplied geometry uses the densest gravity-height slab for its reproducible v1236    behavior and the robust second percentile for v2.  Keep that switch here, after model237    adapters have produced world points, so it cannot become model-specific.238    """239    heights = np.asarray(points, np.float64) @ np.asarray(gravity, np.float64)240    heights = heights[np.isfinite(heights)]241    if not len(heights):242        return 0.0243    if v2 is None:244        v2 = os.environ.get("VSI_CODE_V2") == "1"245    if v2:246        return float(np.percentile(heights, 2))247    counts, edges = np.histogram(heights, bins=80)248    index = int(counts.argmax())249    return float(0.5 * (edges[index] + edges[index + 1]))250 251 252def _object_records(insts, count, u, v, g, floor_level):253    """Up to `count` instances, strongest-evidence first (most observed points = best-segmented,254    closest, most geometry). `count` (peak co-visibility) decides HOW MANY; total observed points255    decide WHICH -- no threshold. Positions are projected onto the shared gravity floor basis256    (u, v horizontal; g up; height 0 = floor_level) and emitted directly in THE final spatial257    code shape: unit-strings ("1.4 meters"), spaced keys ("x coordinate"), and exactly two258    fields per instance (position + longest dimension) -- there is no separate raw form.259 260    Deliberately does NOT report a per-instance first_seen_seconds: the reported instances are261    chosen by STRONGEST evidence (most points/frames), but the class's true first appearance can262    come from a weaker, earlier masklet that never makes this cut (confirmed empirically -- e.g. a263    brief early detection with few points, superseded here by a longer later observation of264    presumably the same object). A per-instance timestamp here would silently describe a DIFFERENT265    detection than the class-level "first seen" a reader would assume it means. appearance_order266    (built below in build_spatial_code() from min(first_time) over ALL detected masklets, not just267    the reported ones) is the sole reliable source for first-appearance timing."""268    ranked = sorted(insts, key=lambda i: (i["n"], i.get("nframes", 0)), reverse=True)[269        : max(count, 1)270    ]271    recs = []272    for it in ranked:273        c = np.asarray(it["centroid"], np.float64)274        recs.append(275            {276                "position": {277                    "x coordinate": f"{round(float(c @ u), 2)} meters",278                    "y coordinate": f"{round(float(c @ v), 2)} meters",279                    "height above floor": f"{round(float(c @ g - floor_level), 2)} meters",280                },281                "longest dimension": f"{round(float(it['size']), 2)} meters",282            }283        )284    return recs285 286 287def to_spatial_code(instances, stats, floor_area, up_axis, up_vec, floor_level):288    """Per-instance spatial code, emitted directly in THE final shape: class -> {count (peak289    co-visibility), instances:[{position, longest dimension}]} plus room -> {"floor area"}.290    Positions are projected onto the GRAVITY floor plane (u, v horizontal via _floor_basis -- the291    SAME frame as compute_floor_area); "height above floor" is along gravity with 0 = on the292    floor. floor_level is the gravity-height of the floor. (up_axis is retained for signature293    compatibility; the frame now derives from up_vec.)"""294    u, v, g = _floor_basis(up_vec)295    out = {"objects": {}}296    for cls, insts in instances.items():297        cnt = int(stats[cls]["peak"])298        out["objects"][cls] = {299            "count": cnt,300            "instances": _object_records(insts, cnt, u, v, g, floor_level),301        }302    out["room"] = {"floor area": f"{floor_area} square meters"}303    return out304 305 306# ==========================================================================================307# DETERMINISTIC ANSWER LAYER (parameter-free; validated on VSI GT). Reads the in-memory308# instances (pos for direction/route, point clouds for distance). Merged in from309# perceptual.py, verbatim.310# ==========================================================================================311 312 313def _find_cls(name, classes):314    name = name.strip().lower()315    for c in classes:316        if c == name or c.replace(" ", "") == name.replace(" ", ""):317            return c318    for c in classes:319        if name in c or c in name:320            return c321    return None322 323 324def _rep(insts):325    """Prefer a persistent track, using point support as the tie-breaker."""326    return max(insts, key=lambda i: (i.get("nframes", 0), i["n"]))327 328 329def answer_rel_direction(point_a, point_b, point_c, up_vec, up_ax, mode="hard"):330    """Standing at A facing B, where is C? front/back=dot(C-A,fwd); left/right=dot(C-A, up x fwd).331    Right-handed world (OpenCV cam frame + det+1 c2w) makes up x fwd = left a fixed identity.332    Projection uses the gravity VECTOR (v-(v.g)g), so a tilted floor (ScanNet++) is handled; for an333    axis-aligned up this reduces to zeroing that axis (ARKit unchanged)."""334    g = np.asarray(up_vec, np.float64)335    g = g / (np.linalg.norm(g) + 1e-12)336 337    def fl(v):338        w = v.astype(np.float64)339        return w - (w @ g) * g340 341    fwd = fl(point_b - point_a)342    n = np.linalg.norm(fwd)343    if n < 1e-6:344        return None345    fwd /= n346    left = np.cross(g, fwd)347    left /= np.linalg.norm(left) + 1e-9348    d = fl(point_c - point_a)349    f = float(d @ fwd)350    lateral = float(d @ left)351    if mode == "medium":352        if abs(np.degrees(np.arctan2(lateral, f))) >= 135:353            return "back"354        return "left" if lateral > 0 else "right"355    return f"{'front' if f > 0 else 'back'}-{'left' if lateral > 0 else 'right'}"356 357 358def _classify_turn(h_in, h_out, up_vec, up_ax):359    """rotation h_in->h_out in floor plane -> turn left/right/back (135deg = VSI's own 'back' cutoff)."""360    g = np.asarray(up_vec, np.float64)361    g = g / (np.linalg.norm(g) + 1e-12)362 363    def fl(v):364        w = v.astype(np.float64)365        return w - (w @ g) * g366 367    a, b = fl(h_in), fl(h_out)368    na, nb = np.linalg.norm(a), np.linalg.norm(b)369    if na < 1e-6 or nb < 1e-6:370        return None371    a /= na372    b /= nb373    ang = np.degrees(np.arctan2(float(up_vec @ np.cross(a, b)), float(a @ b)))374    if abs(ang) >= 135:375        return "turn back"376    return "turn left" if ang > 0 else "turn right"377 378 379def answer_route(ql, cents, up_vec, up_ax):380    """Chain the turn primitive over the waypoint sequence -> ordered ['turn left/right/back', ...]."""381    import re as _re382 383    m = _re.search(r"beginning at the (.+?) (?:and )?facing the (.+?)\.", ql)384    if not m:385        return None386 387    def position(name):388        class_name = _find_cls(name, cents)389        return cents[class_name] if class_name else None390 391    steps_txt = ql.split(":", 1)[1] if ":" in ql else ql392    steps = _re.findall(393        r"\d+\.\s*(\[please fill in\]|Go forward until the [^0-9\[]+?)(?=\s*\d+\.|$)",394        steps_txt,395    )396    cur_pos = position(m.group(1).strip())397    if cur_pos is None:398        return None399    fac = position(m.group(2).strip())400    cur_head = (fac - cur_pos) if fac is not None else None401    turns = []402    i = 0403    while i < len(steps):404        s = steps[i].strip()405        if s.startswith("Go forward"):406            tp = position(_re.sub(r"^Go forward until the ", "", s).strip().rstrip("."))407            if tp is not None:408                cur_head = tp - cur_pos409                cur_pos = tp410        else:411            nxt = next(412                (413                    _re.sub(r"^Go forward until the ", "", steps[j].strip())414                    .strip()415                    .rstrip(".")416                    for j in range(i + 1, len(steps))417                    if steps[j].strip().startswith("Go forward")418                ),419                None,420            )421            tp = position(nxt) if nxt else None422            if tp is None or cur_head is None:423                turns.append(None)424            else:425                turns.append(_classify_turn(cur_head, tp - cur_pos, up_vec, up_ax))426                cur_head = tp - cur_pos427        i += 1428    return turns429 430 431# ==========================================================================================432# POINT-CLOUD CLEANING + DISTANCE ANSWERS -- outlier removal, closest-distance queries. Merged433# in from perceptual.py, verbatim.434# ==========================================================================================435 436 437def _sor(pts, k=16, std=2.0, cap=4000):438    """Statistical outlier removal: drop points whose mean distance to their k nearest neighbors exceeds439    mean + std*sigma. Removes mask-bleed / depth-speckle points that corrupt a literal closest-point min.440    k=16 and std=2 are universal robust-statistics defaults -- NOT tuned to VSI (no benchmark-fit knob).441    """442    from scipy.spatial import cKDTree443 444    if len(pts) < k + 2:445        return pts446    rs = np.random.RandomState(0)447    points = pts if len(pts) <= cap else pts[rs.choice(len(pts), cap, False)]448    d, _ = cKDTree(points).query(points, k=k + 1, workers=KD_WORKERS)449    md = d[:, 1:].mean(1)450    return points[md <= md.mean() + std * md.std()]451 452 453def _main_cluster(pts):454    """Keep the object's dominant spatial cluster. Coherent mask-bleed onto a SPATIALLY-SEPARATED adjacent455    object forms a disconnected component (a gap separates two objects); the true object is the largest one.456    Linkage scale self-calibrates from the cloud's own nearest-neighbor spacing -- no fixed distance.457    """458    from scipy.spatial import cKDTree459 460    if len(pts) < 30:461        return pts462    tree = cKDTree(pts)463    nn, _ = tree.query(pts, k=2)464    eps = 3.0 * float(np.median(nn[:, 1]))  # 3x median NN gap (data-derived)465    pairs = tree.query_pairs(eps, output_type="ndarray")466    if len(pairs) == 0:467        return pts468    parent = np.arange(len(pts))469 470    def find(x):471        r = x472        while parent[r] != r:473            r = parent[r]474        while parent[x] != r:475            parent[x], x = r, parent[x]476        return r477 478    for a, b in pairs:479        ra, rb = find(int(a)), find(int(b))480        if ra != rb:481            parent[ra] = rb482    roots = np.array([find(i) for i in range(len(pts))])483    v, cnt = np.unique(roots, return_counts=True)484    return pts[roots == v[cnt.argmax()]]485 486 487def _clean(inst, cap=4000):488    """Outlier removal self-calibrated from the pipeline's OWN per-object signals -- no universal constant,489    no benchmark knob. Cached on the instance so the O(classes^2) distance table cleans each object once.490      (1) per-point DA3 confidence: boundary mask-bleed = depth discontinuity = LOW conf -> drop below the491          object's OWN median confidence (data-derived cut);492      (2) statistical density outlier removal on the survivors;493      (3) dominant spatial cluster: coherent bleed onto a separated adjacent object is a disconnected cluster.494    """495    if inst.get("_cleanpts") is not None:496        return inst["_cleanpts"]497    pts = inst["pts"]498    conf = inst.get("conf")499    rs = np.random.RandomState(0)500    if len(pts) > cap:501        idx = rs.choice(len(pts), cap, False)502        pts = pts[idx]503        conf = conf[idx] if conf is not None else None504    if conf is not None and len(conf) > 20:  # self-calibrating: object's OWN median505        keep = conf >= np.median(conf)506        if keep.sum() >= 10:507            pts = pts[keep]508    pts = _sor(pts, cap=cap)  # density outlier removal on the survivors509    # NOTE: spatial de-bleeding (_main_cluster) was tried and REVERTED -- "largest cluster = the object" is510    # not guaranteed; when an object's near edge splits off it drops the true closest part (overshoot). You511    # cannot post-hoc recover the true object from a bleeding mask in 3D -- that needs better SOURCE masks.512    inst["_cleanpts"] = pts513    return pts514 515 516def answer_closest_distance(instances_a, instances_b, k=4000):517    """Closest distance between the two objects' point clouds ('closest point of each object'). Points are518    cleaned by _clean (self-calibrated per-object confidence + density) then exact nearest-neighbor min via519    KD-tree, so mask-bleed/speckle can't collapse the answer. NOTE: VSI's GT is box-to-box on CLEAN annotation520    boxes; the residual gap is object-segmentation quality (bleed), not this formula -- see _clean docstring.521    """522    points_a = _clean(_rep(instances_a), cap=k)523    points_b = _clean(_rep(instances_b), cap=k)524    if len(points_a) == 0 or len(points_b) == 0:525        return float("inf")526    from scipy.spatial import cKDTree527 528    if len(points_a) <= len(points_b):529        d, _ = cKDTree(points_a).query(points_b, k=1, workers=KD_WORKERS)530    else:531        d, _ = cKDTree(points_b).query(points_a, k=1, workers=KD_WORKERS)532    return float(d.min())533 534 535def answer_rel_distance(anchor_insts, option_insts):536    """'which option is closest to the anchor?' -> index of the option with min closest-point distance."""537    dists = [538        answer_closest_distance(anchor_insts, oi) if oi is not None else float("inf")539        for oi in option_insts540    ]541    return int(np.argmin(dists)), dists542 543 544# ==========================================================================================545# MASK/DEPTH CLEANUP + BACK-PROJECTION -- per-frame refinement before points enter an546# instance's point cloud. Merged in from perceptual.py, verbatim.547# ==========================================================================================548 549 550def depth_edges(depth_f, valid_f):551    """Per-frame depth-DISCONTINUITY map: gradient magnitude above median + 3*MAD over the valid pixels.552    The threshold is data-derived (robust, universal statistical cut -- NOT benchmark-tuned). These edges are553    where mask-bleed crosses onto an adjacent object at a different depth."""554    if valid_f.sum() < 100:555        return np.zeros_like(depth_f, bool)556    d = np.where(valid_f, depth_f, np.median(depth_f[valid_f]))557    gy, gx = np.gradient(d)558    grad = np.hypot(gx, gy)559    g = grad[valid_f]560    med = np.median(g)561    mad = np.median(np.abs(g - med)) + 1e-9562    return (grad > med + 3.0 * 1.4826 * mad) & valid_f563 564 565try:566    _GUIDED = (567        cv2.ximgproc.guidedFilter568    )  # opencv-contrib; appearance-guided boundary snap569except AttributeError:570    _GUIDED = None571 572 573def refine_mask(mask, rgb):574    """Snap a coarse SAM3 mask boundary to the RGB color edge (appearance-guided). This clips the same-depth575    bleed that depth CANNOT see: where the halo crosses onto a differently-colored neighbor, the color edge576    cuts it. Params are DERIVED, not tuned -- radius from frame scale, smoothness from the image's own color577    variance (same MAD-style principle as depth_edges). No model, no GPU. Falls back to grabCut, then erode.578    """579    if mask.shape[:2] != rgb.shape[:2]:580        mask = cv2.resize(581            mask.astype(np.uint8),582            (rgb.shape[1], rgb.shape[0]),583            interpolation=cv2.INTER_NEAREST,584        ).astype(bool)585    a = int(mask.sum())586    if a < 60:  # too small to refine meaningfully -> 1px erode as before587        me = cv2.erode(mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1).astype(bool)588        return me if me.any() else mask589    if _GUIDED is not None:590        r = max(591            4, int(round(0.02 * float(np.hypot(*mask.shape[:2]))))592        )  # radius ~2% of frame diagonal593        eps = (594            float(np.var(rgb.astype(np.float32) / 255.0)) * 0.01 + 1e-6595        )  # smoothness ~ image color variance596        soft = _GUIDED(rgb, mask.astype(np.float32), r, eps)597        out = soft > 0.5598        return out if out.sum() >= 0.4 * a else mask  # majority guard: don't over-carve599    # fallback (base opencv, no ximgproc): grabCut seeded FG=mask, PR_FG=dilated ring, BG=far exterior600    try:601        gc = np.full(mask.shape[:2], cv2.GC_PR_BGD, np.uint8)602        dil = cv2.dilate(mask.astype(np.uint8), np.ones((15, 15), np.uint8), 1).astype(603            bool604        )605        er = cv2.erode(mask.astype(np.uint8), np.ones((5, 5), np.uint8), 1).astype(bool)606        gc[dil] = cv2.GC_PR_FGD607        gc[er] = cv2.GC_FGD608        bgm = np.zeros((1, 65), np.float64)609        fgm = np.zeros((1, 65), np.float64)610        cv2.grabCut(rgb, gc, None, bgm, fgm, 3, cv2.GC_INIT_WITH_MASK)611        out = (gc == cv2.GC_FGD) | (gc == cv2.GC_PR_FGD)612        return out if 0.4 * a <= out.sum() <= 1.5 * a else mask613    except Exception:614        me = cv2.erode(mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1).astype(bool)615        return me if me.any() else mask616 617 618def backproject_frame(619    depth_f,620    intrinsics,621    c2w_f,622    mask_f,623    conf_f=None,624    conf_thr=0.0,625    valid_f=None,626    return_conf=False,627    edges_f=None,628):629    """Return (M,3) world points for the masked pixels of one frame (c2w_f = cam->world 4x4).630    valid_f: optional precomputed (isfinite & >0) depth mask, reused across all masks of a frame.631    return_conf: also return the (M,) DA3 confidence of each kept point (for per-point noise-aware cleaning).632    edges_f: optional per-frame depth-edge map; if given, keep the mask's largest depth-coherent component633             (cuts mask-bleed onto adjacent objects at the depth boundary, per frame, before back-projection).634    """635    height, width = depth_f.shape636    empty = (637        (np.empty((0, 3), np.float32), np.empty((0,), np.float32))638        if return_conf639        else np.empty((0, 3), np.float32)640    )641    if mask_f.shape != (height, width):642        mask_f = cv2.resize(643            mask_f.astype(np.uint8), (width, height), interpolation=cv2.INTER_NEAREST644        ).astype(bool)645    if valid_f is None:646        valid_f = np.isfinite(depth_f) & (depth_f > 0)647    m = mask_f & valid_f648    if conf_f is not None and conf_thr > 0:649        m &= conf_f >= conf_thr650    if (651        edges_f is not None and m.sum() >= 30652    ):  # keep the largest depth-coherent piece of the mask653        from scipy import ndimage654 655        lab, nlab = ndimage.label(m & ~edges_f)656        if nlab >= 1:657            sizes = np.bincount(lab.ravel())658            sizes[0] = 0659            big = int(sizes.argmax())660            if (661                sizes[big] >= 0.5 * m.sum()662            ):  # object is the MAJORITY piece (bleed is a minority)663                m = lab == big664    if not m.any():665        return empty666    ys, xs = np.nonzero(m)667    z = depth_f[ys, xs]668    if DEPTH_COHERENCE and len(z) >= 8:669        # an object is a coherent depth surface; floor/background BLEED pixels are depth outliers.670        # Drop them via the standard Tukey fence (1.5*IQR) on the masked region's depths -- parameter-free.671        q1, q3 = np.percentile(z, [25, 75])672        iqr = q3 - q1673        keep = (z >= q1 - 1.5 * iqr) & (z <= q3 + 1.5 * iqr)674        if keep.sum() >= 1:675            ys, xs, z = ys[keep], xs[keep], z[keep]676    fx, fy, cx, cy = (677        intrinsics[0, 0],678        intrinsics[1, 1],679        intrinsics[0, 2],680        intrinsics[1, 2],681    )682    camera_points = np.stack(683        [(xs - cx) * z / fx, (ys - cy) * z / fy, z], axis=1684    )  # camera coords685    world_points = (c2w_f[:3, :3] @ camera_points.T).T + c2w_f[:3, 3]  # -> world686    if return_conf:687        cw = (688            conf_f[ys, xs].astype(np.float32)689            if conf_f is not None690            else np.ones(len(ys), np.float32)691        )692        return world_points.astype(np.float32), cw693    return world_points.astype(np.float32)694 695 696# ==========================================================================================697# INSTANCE BUILDING -- oriented extent, 3D box-overlap re-identification, and the main698# build_instances() driver that turns per-frame masks into per-class 3D instances. Merged in699# from perceptual.py, verbatim.700# ==========================================================================================701 702 703def robust_centroid_extent(pts, up_axis=None):704    """median centroid + ORIENTED robust extent.705    If up_axis is given: YAW-ONLY oriented extent -- rotation is found by 2D PCA on the706    floor-projected points only, with the up axis left untouched. This matches VSI-Bench's707    own annotation convention for indoor scans (ScanNet/ARKitScenes OrientedBoundingBox708    objects follow a Manhattan-world assumption: rotated only around the vertical axis,709    never tilted). Unconstrained 3D PCA (the previous behavior) can chase noise on the710    vertical axis for flat/elongated objects and drift away from the true yaw.711    If up_axis is None (unknown at the call site): falls back to unconstrained 3D PCA.712    Either way: parameter-free, rotation-invariant in-plane, p2..p98 robust extent."""713    c = np.median(pts, axis=0)714    centered = pts - c715    if len(centered) > 5000:  # PCA on a sample (deterministic)716        centered = centered[np.random.RandomState(0).choice(len(centered), 5000, False)]717    if up_axis is not None:718        floor_axes = [i for i in range(3) if i != up_axis]719        floor_points = centered[:, floor_axes]720        try:721            _, _, floor_rotation = np.linalg.svd(722                floor_points - floor_points.mean(0), full_matrices=False723            )724            proj_floor = (725                floor_points @ floor_rotation.T726            )  # (N,2) along the object's own floor-plane axes727        except np.linalg.LinAlgError:728            proj_floor = floor_points729        up_col = centered[:, up_axis : up_axis + 1]  # up axis untouched (yaw-only)730        proj = np.concatenate([proj_floor, up_col], axis=1)731    else:732        try:733            _, _, rotation = np.linalg.svd(734                centered - centered.mean(0), full_matrices=False735            )736            proj = centered @ rotation.T  # coordinates along principal axes737        except np.linalg.LinAlgError:738            proj = centered739    lo = np.percentile(proj, 2, axis=0)740    hi = np.percentile(proj, 98, axis=0)741    ext = np.maximum(hi - lo, 0.0)742    dims = np.sort(ext)[::-1]  # the object's 3 oriented side lengths, longest first743    return c.astype(np.float32), float(dims[0]), dims744 745 746def _aabb(pts):747    """robust (p2..p98) axis-aligned 3D box of an instance's world points."""748    return np.percentile(pts, 2, axis=0), np.percentile(pts, 98, axis=0)749 750 751def merge_by_box_overlap(insts, up_axis=None):752    """Parameter-free 3D re-identification, with a temporal-exclusion gate.753 754    SAM3 (a 2D tracker) emits a NEW masklet each time the camera revisits an object, so one755    physical object -> several masklets at the same 3D location. We fuse two same-class masklets756    iff BOTH:757      (a) their measured 3D boxes overlap (or one centroid is inside the other) -- same volume, and758      (b) they NEVER appear in the same frame -- temporal exclusion.759    (b) is the key parameter-free invariant: two masklets co-visible in one frame were tracked by760    SAM3 as distinct objects in that frame, so they ARE distinct (e.g. two chairs around a table);761    we must never merge them, even if depth noise makes their boxes overlap. A revisit-duplicate,762    by contrast, lives in DISJOINT frames. Both tests are exact/measured -- NO tuned threshold.763    """764    n = len(insts)765    if n <= 1:766        return insts767    # Precompute the FULL pairwise box-match as one vectorized boolean matrix. With the per-frame768    # batched detector a class can have 100+ raw instances; the old O(n^3) loop called np.all-based769    # box_match millions of times. Here every pairwise test is one broadcast -> O(1) lookups below.770    lo = np.stack([_aabb(i["pts"])[0] for i in insts]).astype(np.float32)  # (n,3)771    hi = np.stack([_aabb(i["pts"])[1] for i in insts]).astype(np.float32)  # (n,3)772    cents = np.stack([i["centroid"] for i in insts]).astype(np.float32)  # (n,3)773    overlap = (hi[:, None, :] >= lo[None, :, :]).all(-1) & (774        hi[None, :, :] >= lo[:, None, :]775    ).all(-1)776    ins = (cents[:, None, :] >= lo[None, :, :]).all(-1) & (777        cents[:, None, :] <= hi[None, :, :]778    ).all(-1)779    match = (780        overlap | ins | ins.T781    )  # box_match[i,j]: same 3D volume (identical semantics to old code)782    # group-level agglomeration: merge two groups only if their COMBINED frame sets are disjoint783    # (so no two co-visible masklets ever land in one object) AND some cross-pair shares a 3D volume.784    groups = [785        {"members": [i], "frames": set(insts[i].get("frames", set()))} for i in range(n)786    ]787    changed = True788    while changed:789        changed = False790        for a in range(len(groups)):791            for b in range(a + 1, len(groups)):792                if (793                    groups[a]["frames"] & groups[b]["frames"]794                ):  # co-visible -> distinct objects795                    continue796                if match[np.ix_(groups[a]["members"], groups[b]["members"])].any():797                    groups[a]["members"] += groups[b]["members"]798                    groups[a]["frames"] |= groups[b]["frames"]799                    groups.pop(b)800                    changed = True801                    break802            if changed:803                break804    merged = []805    for g in groups:806        idxs = g["members"]807        pts = np.concatenate([insts[k]["pts"] for k in idxs], 0)808        cpts = np.concatenate(809            [810                insts[k].get("conf", np.ones(len(insts[k]["pts"]), np.float32))811                for k in idxs812            ],813            0,814        )815        best = max(816            (insts[k]["best_pts"] for k in idxs), key=len817        )  # largest single obs in the group818        c, longest, dims = robust_centroid_extent(819            best, up_axis820        )  # size+pos+3 oriented dims from best view (#2/#4)821        merged.append(822            {823                "centroid": c,824                "size": longest,825                "dims": dims,826                "first_time": min(insts[k]["first_time"] for k in idxs),827                "pts": pts,828                "conf": cpts,829                "best_pts": best,830                "n": sum(insts[k]["n"] for k in idxs),831                "nframes": len(g["frames"]),832            }833        )  # track persistence (evidence strength)834    return merged835 836 837def build_instances(838    per_class, depth, intr, c2w, conf, frame_times, frame_paths=None, up_axis=None839):840    """-> {class: [ {centroid(3), size(longest dim), first_time, npts} ]}841    up_axis: if known (from room_gravity, computed BEFORE this call), threads through to842    robust_centroid_extent for yaw-only oriented sizing. If None, size falls back to843    unconstrained 3D PCA."""844    out = {}845    stats = {}846    nframes = len(depth)847    # precompute each frame's valid-depth mask ONCE (was recomputed per mask -> per class).848    valid = {}849    edges = {}850    rgb = {}851    used_frames = {fi for frames in per_class.values() for fi in frames if fi < nframes}852    for fi in used_frames:853        valid[fi] = np.isfinite(depth[fi]) & (depth[fi] > 0)854        edges[fi] = depth_edges(depth[fi], valid[fi]) if DEPTH_EDGE_REFINE else None855        if MASK_REFINE and frame_paths and fi < len(frame_paths):856            im = cv2.imread(frame_paths[fi])  # BGR; guidedFilter/grabCut want 3ch uint8857            rgb[fi] = (858                cv2.resize(im, (depth[fi].shape[1], depth[fi].shape[0]))859                if im is not None860                else None861            )862    for cls, frames in per_class.items():863        # peak co-visibility: max distinct masklets SAM3 tracks SIMULTANEOUSLY in any one frame.864        # geometry-free, parameter-free, immune to revisit over-count; provable lower bound on count.865        peak = max((len(objs) for objs in frames.values()), default=0)866        # gather world points + first-seen time + frame set per obj_id (SAM3 track id = masklet)867        pts_by_id, conf_by_id, first_t, frames_by_id = {}, {}, {}, {}868        # px_by_id: (fidx, frame_time, mask_pixel_count) per obj_id, for appearance-order timing.869        # VSI-Bench's own GT defines "first appearance" as the timestamp where an object's pixel870        # count crosses a threshold (paper appendix B.1) -- NOT the first frame with any pixel at871        # all. A single stray mask-bleed/false-positive pixel would otherwise register as "first872        # seen" far too early. The threshold is self-calibrated per instance (that instance's OWN873        # median observed pixel count across its frames), same convention as _clean's median cut.874        px_by_id = {}875        for fidx, objs in frames.items():876            if fidx >= nframes:  # guard: SAM3 frame idx vs DA3 frames877                continue878            for oid, mask in objs.items():879                if MASK_REFINE and rgb.get(fidx) is not None:880                    me = (881                        mask882                        if os.environ.get("VSI_NO_REFINE") == "1"883                        else refine_mask(mask, rgb[fidx])884                    )  # appearance-guided boundary snap (clips same-depth bleed)885                else:886                    # erode 1px to drop mask-edge / background depth bleed887                    me = cv2.erode(888                        mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1889                    ).astype(bool)890                    if not me.any():891                        me = mask892                conf_f = conf[fidx] if conf is not None else None893                conf_thr = (894                    np.percentile(conf_f, CONF_PCT)895                    if (conf_f is not None and CONF_PCT > 0)896                    else 0.0897                )898                world_points, cw = backproject_frame(899                    depth[fidx],900                    intr[fidx],901                    c2w[fidx],902                    me,903                    conf_f,904                    conf_thr=conf_thr,905                    valid_f=valid.get(fidx),906                    return_conf=True,907                    edges_f=edges.get(fidx),908                )909                if len(world_points):910                    pts_by_id.setdefault(oid, []).append(world_points)911                    conf_by_id.setdefault(oid, []).append(912                        cw913                    )  # per-point DA3 confidence (for _clean)914                    frames_by_id.setdefault(oid, set()).add(915                        fidx916                    )  # for co-occurrence gate917                    t = frame_times[fidx]918                    px_by_id.setdefault(oid, []).append((fidx, t, int(me.sum())))919        # appearance_order timing: per instance, first frame at/above its OWN median pixel count920        for oid, obs in px_by_id.items():921            counts = [c for _, _, c in obs]922            thresh = float(np.median(counts))923            crossing = [t for _, t, c in obs if c >= thresh]924            first_t[oid] = min(crossing) if crossing else min(t for _, t, c in obs)925        insts = []926        for oid, plist in pts_by_id.items():927            pts = np.concatenate(plist, 0)928            cpts = np.concatenate(conf_by_id[oid], 0)929            if len(pts) < MIN_INSTANCE_PTS:930                continue931            best = max(932                plist, key=len933            )  # #2/#4: largest single-frame observation (closest/most pixels)934            bc, bsize, bdims = robust_centroid_extent(935                best, up_axis936            )  # size + position + 3 oriented dims from best view937            insts.append(938                {939                    "centroid": bc,940                    "size": bsize,941                    "dims": bdims,942                    "first_time": first_t[oid],943                    "pts": pts,944                    "conf": cpts,945                    "best_pts": best,946                    "n": len(pts),947                    "frames": frames_by_id[oid],948                }949            )950        raw = len(insts)951        # 3D re-ID: fuse same-class masklets that occupy the same measured 3D volume (parameter-free).952        insts = merge_by_box_overlap(insts, up_axis)953        if insts:954            out[cls] = insts955            stats[cls] = {"raw": raw, "merged": len(insts), "peak": peak}956    return out, stats957 958 959# ==========================================================================================960# PER-CLASS SUMMARY + FLOOR AREA -- legacy array-schema row builder, and the room-scale floor961# area calculation. Merged in from perceptual.py, verbatim.962# ==========================================================================================963 964 965# ---- Per-class spatial code (legacy array-schema row) -------------------------------------966def class_spatial_code(insts, peak=0):967    cents = np.stack([i["centroid"] for i in insts], 0)  # (n,3)968    sizes = np.array([i["size"] for i in insts], np.float32)969    n = len(insts)  # merged centroids (for spatial stats)970    count = peak if peak else n  # reported count = peak co-visibility971    x, y, z = cents.mean(0)972    if n == 1:973        e1 = e2 = e3 = px = pz = 0.0974        size_iqr = 0.0975    else:976        cov = np.cov(cents.T)  # 3x3977        vals, vecs = np.linalg.eigh(cov)  # ascending978        order = np.argsort(vals)[::-1]979        vals = np.clip(vals[order], 0, None)980        vecs = vecs[:, order]981        e1, e2, e3 = vals.tolist()982        pv = vecs[:, 0]  # principal eigenvector983        px, pz = float(pv[0]), float(pv[2])984        q1, q3 = np.percentile(sizes, [25, 75])985        size_iqr = float(q3 - q1)986    size_median = float(np.median(sizes))987    first_time = float(min(i["first_time"] for i in insts))988    row = [x, y, z, e1, e2, e3, px, pz, size_median, size_iqr, first_time, count]989    row = [990        (lambda r: 0.0 if r == 0 else r)(round(float(v), 1)) for v in row991    ]  # kill -0.0992    row[-1] = int(count)993    return row994 995 996# ---- floor_area (full-scene min-Y points -> XZ convex hull) -------------------------------997def compute_floor_area(depth, intr, c2w, conf, sky, stride=8, up_vec=None):998    pts = []999    for f in range(depth.shape[0]):1000        height, width = depth[f].shape1001        ys, xs = np.mgrid[0:height:stride, 0:width:stride]1002        ys = ys.ravel()1003        xs = xs.ravel()1004        z = depth[f][ys, xs]1005        ok = np.isfinite(z) & (z > 0)1006        if sky is not None:1007            ok &= ~sky[f][ys, xs].astype(bool)1008        if conf is not None:1009            ok &= conf[f][ys, xs] >= np.percentile(conf[f], 40)1010        ys, xs, z = ys[ok], xs[ok], z[ok]1011        if not len(z):1012            continue1013        intrinsics = intr[f]1014        fx, fy, cx, cy = (1015            intrinsics[0, 0],1016            intrinsics[1, 1],1017            intrinsics[0, 2],1018            intrinsics[1, 2],1019        )1020        camera_points = np.stack([(xs - cx) * z / fx, (ys - cy) * z / fy, z], 1)1021        world_points = (c2w[f][:3, :3] @ camera_points.T).T + c2w[f][:3, 3]1022        pts.append(world_points.astype(np.float32))1023    if not pts:1024        return 0.01025    points = np.concatenate(pts, 0)1026    if up_vec is not None:1027        # VSI-faithful: area in the plane orthogonal to GRAVITY (RANSAC floor normal), like the1028        # benchmark's gravity-aligned GT meshes. Build an orthonormal in-plane basis (u, v).1029        g = np.asarray(up_vec, np.float64)1030        g /= np.linalg.norm(g) + 1e-121031        a = np.array([1.0, 0.0, 0.0]) if abs(g[0]) < 0.9 else np.array([0.0, 1.0, 0.0])1032        u = np.cross(g, a)1033        u /= np.linalg.norm(u)1034        v = np.cross(g, u)1035        all_floor_points = np.stack([points @ u, points @ v], 1)1036    else:1037        up = int(1038            np.argmin(points.max(0) - points.min(0))1039        )  # legacy: vertical = smallest-extent axis1040        floor_axes = [i for i in range(3) if i != up]1041        all_floor_points = points[:, floor_axes]1042    # VSI-Bench room-size definition = alpha-shape of the floor-plane point cloud (confirmed in their1043    # paper appendix). VSI does not publish the alpha value they use for their own GT mesh, so alpha=21044    # here is NOT a matched/verified constant -- it was chosen empirically for this pipeline's own1045    # (sparser) reconstructed point density. This is the one disclosed benchmark-adjacent tuned constant1046    # in the whole file; everything else is exact/derived or a generic, non-tuned statistical convention.1047    # (Falls back to enclosed-fill below if the alphashape package isn't available.)1048    floor_points = all_floor_points1049    lo = np.percentile(floor_points, 0.5, 0)1050    hi = np.percentile(floor_points, 99.5, 0)  # gentle clip (preserve room extent)1051    floor_points = floor_points[1052        (floor_points[:, 0] >= lo[0])1053        & (floor_points[:, 0] <= hi[0])1054        & (floor_points[:, 1] >= lo[1])1055        & (floor_points[:, 1] <= hi[1])1056    ]1057    if len(floor_points) < 10:1058        return 0.01059    try:1060        import alphashape1061 1062        idx = np.random.RandomState(0).choice(1063            len(floor_points), min(10000, len(floor_points))1064        )1065        return round(1066            float(alphashape.alphashape(floor_points[idx], alpha=2).area), 11067        )  # alpha=2 tuned for recon density1068    except Exception:1069        from scipy import ndimage1070 1071        res = 0.101072        ai = ((floor_points[:, 0] - floor_points[:, 0].min()) / res).astype(int)1073        bi = ((floor_points[:, 1] - floor_points[:, 1].min()) / res).astype(int)1074        grid = np.zeros((ai.max() + 3, bi.max() + 3), np.uint8)1075        grid[ai + 1, bi + 1] = 11076        grid = cv2.morphologyEx(grid, cv2.MORPH_CLOSE, np.ones((7, 7), np.uint8))1077        grid = ndimage.binary_fill_holes(grid).astype(np.uint8)1078        return round(float(grid.sum()) * res * res, 1)1079 1080 1081# ---------------------------------------------------------------------------1082 1083# ==========================================================================================1084# SPATIAL CODE ASSEMBLY -- the top-level entry point this whole file exists for:1085# build_spatial_code() calls everything above to turn already-computed depth/pose/masks into1086# the final spatial code dict. Room outline + JSON writer. This section was already in1087# geometric.py before the perceptual.py merge; build_spatial_code() below is updated to call1088# the geometry functions above DIRECTLY (no `pl.` prefix -- they're plain local functions now1089# that everything is in one file), same logic, unchanged otherwise.1090# ==========================================================================================1091 1092 1093def _room_outline(depth, intr, c2w, conf, bu, bv):1094    """(Currently unemitted -- the one spatial code shape has no room outline field; this1095    math is kept intact for reuse.) Room floor-boundary polygon from the SAME grid as1096    compute_floor_area: floor points1097    projected onto the shared gravity plane (bu, bv), 10cm grid, close 7x7, fill holes,1098    largest contour, 0.2m polygon simplification. Same (bu, bv) as objects/camera, so the1099    outline, object positions, and floor_area all live in one consistent frame."""1100    from scipy import ndimage1101 1102    pts, stride = [], 81103    for f in range(0, depth.shape[0], 3):1104        height, width = depth[f].shape1105        ys, xs = np.mgrid[0:height:stride, 0:width:stride]1106        ys = ys.ravel()1107        xs = xs.ravel()1108        z = depth[f][ys, xs]1109        ok = np.isfinite(z) & (z > 0)1110        if conf is not None:1111            ok &= conf[f][ys, xs] >= np.percentile(conf[f], 40)1112        ys, xs, z = ys[ok], xs[ok], z[ok]1113        if not len(z):1114            continue1115        intrinsics = intr[f]1116        camera_points = np.stack(1117            [1118                (xs - intrinsics[0, 2]) * z / intrinsics[0, 0],1119                (ys - intrinsics[1, 2]) * z / intrinsics[1, 1],1120                z,1121            ],1122            1,1123        )1124        pts.append(1125            ((c2w[f][:3, :3] @ camera_points.T).T + c2w[f][:3, 3]).astype(np.float32)1126        )1127    if not pts:1128        return []1129    world_points = np.concatenate(pts, 0)1130    points = np.stack(1131        [world_points @ bu, world_points @ bv], 11132    )  # gravity-plane projection (same bu,bv as objects/area)1133    lo = np.percentile(points, 0.5, 0)1134    hi = np.percentile(points, 99.5, 0)1135    points = points[1136        (points[:, 0] >= lo[0])1137        & (points[:, 0] <= hi[0])1138        & (points[:, 1] >= lo[1])1139        & (points[:, 1] <= hi[1])1140    ]1141    if len(points) < 10:1142        return []1143    res = 0.101144    x0, y0 = points[:, 0].min(), points[:, 1].min()1145    ai = ((points[:, 0] - x0) / res).astype(int)1146    bi = ((points[:, 1] - y0) / res).astype(int)1147    grid = np.zeros((ai.max() + 3, bi.max() + 3), np.uint8)1148    grid[ai + 1, bi + 1] = 11149    grid = cv2.morphologyEx(grid, cv2.MORPH_CLOSE, np.ones((7, 7), np.uint8))1150    grid = ndimage.binary_fill_holes(grid).astype(np.uint8)1151    cs, _ = cv2.findContours(grid, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)1152    if not cs:1153        return []1154    poly = cv2.approxPolyDP(max(cs, key=cv2.contourArea), 0.2 / res, True)[:, 0, :]1155    # cv2 contour points are (col=bi, row=ai) -> (floor_y, floor_x)1156    return [1157        {1158            "floor_x_meters": round(float((r - 1) * res + x0), 1),1159            "floor_y_meters": round(float((c - 1) * res + y0), 1),1160        }1161        for c, r in poly1162    ]1163 1164 1165def build_spatial_code_raw(depth, intr, c2w, conf, ftimes, per):1166    """Builds THE spatial code -- the one and only shape a spatial code has, everywhere1167    (on disk, in prompts, in this pipeline): unit-strings ("1.4 meters"), spaced keys1168    ("x coordinate"), per-instance position + longest dimension only, room "floor area",1169    "closest classes distance meters from" (rooted per class, distance + closeness rank),1170    and a flat earliest-first "appearance order" list of class names. There is no separate1171    raw/rendered split and no schema flag -- the old v1/v2 branching (VSI_CODE_V2) and the1172    raw intermediate form (floor_x_meters keys, bounding_box, dimensions_meters,1173    seen_in_video_frames, camera_trajectory, room.outline) are gone; every underlying VALUE1174    that survives is computed by exactly the same math as before, only the emitted fields1175    and their formatting changed."""1176    # emission-time class rename: VSI's questions say 'coat rack' while their annotations1177    # (and hence the SAM3 prompt + caches) say 'coat hanger' -- same object, their naming1178    # seam. The model sees questions, so emitted codes follow the question vocabulary.1179    class_aliases = {"coat hanger": "coat rack"}1180    per = {class_aliases.get(k, k): v for k, v in per.items()}1181    inst, stats = build_instances(per, depth, intr, c2w, conf, ftimes)1182    up_vec, up_ax = room_gravity(1183        depth, intr, c2w, conf1184    )  # gravity = RANSAC floor normal (VSI-faithful)1185    bu, bv, bg = _floor_basis(1186        up_vec1187    )  # shared gravity floor frame (bu,bv horizontal, bg up)1188    points = np.concatenate([i["pts"] for cl in inst.values() for i in cl], 0)1189    floor_level = _floor_level(points, bg)1190    fa = compute_floor_area(depth, intr, c2w, conf, None, up_vec=up_vec)1191    code = to_spatial_code(inst, stats, fa, up_ax, up_vec, floor_level)1192    cls = list(inst.keys())1193    class_first = {c: min(i["first_time"] for i in v) for c, v in inst.items()}1194 1195    # Keyed dict + integer ranks (not a sorted list): each question option becomes ONE1196    # direct key access, and "which is closest" = min over small integers -- the filtered1197    # list-scan and decimal comparison were the observed failure modes even on GT data.1198    # 2-decimal distances: 0.1m rounding costs up to ~17% relative error on sub-meter1199    # answers, which fails the strictest MRA thresholds even with perfect values.1200    ccf = {}

Showing the first 1,200 of 1590 lines. Download the file for the rest.