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    """representative instance = most observed points (best-segmented, validated 8/8 on direction)."""326    return max(insts, key=lambda i: 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    # Robust box-to-box gap. World geometry is gravity aligned by the scene adapter.527    lo_a, hi_a = np.percentile(points_a, [2, 98], axis=0)528    lo_b, hi_b = np.percentile(points_b, [2, 98], axis=0)529    gap = np.maximum(np.maximum(lo_a - hi_b, lo_b - hi_a), 0.0)530    return float(np.linalg.norm(gap))531 532 533def answer_rel_distance(anchor_insts, option_insts):534    """'which option is closest to the anchor?' -> index of the option with min closest-point distance."""535    dists = [536        answer_closest_distance(anchor_insts, oi) if oi is not None else float("inf")537        for oi in option_insts538    ]539    return int(np.argmin(dists)), dists540 541 542# ==========================================================================================543# MASK/DEPTH CLEANUP + BACK-PROJECTION -- per-frame refinement before points enter an544# instance's point cloud. Merged in from perceptual.py, verbatim.545# ==========================================================================================546 547 548def depth_edges(depth_f, valid_f):549    """Per-frame depth-DISCONTINUITY map: gradient magnitude above median + 3*MAD over the valid pixels.550    The threshold is data-derived (robust, universal statistical cut -- NOT benchmark-tuned). These edges are551    where mask-bleed crosses onto an adjacent object at a different depth."""552    if valid_f.sum() < 100:553        return np.zeros_like(depth_f, bool)554    d = np.where(valid_f, depth_f, np.median(depth_f[valid_f]))555    gy, gx = np.gradient(d)556    grad = np.hypot(gx, gy)557    g = grad[valid_f]558    med = np.median(g)559    mad = np.median(np.abs(g - med)) + 1e-9560    return (grad > med + 3.0 * 1.4826 * mad) & valid_f561 562 563try:564    _GUIDED = (565        cv2.ximgproc.guidedFilter566    )  # opencv-contrib; appearance-guided boundary snap567except AttributeError:568    _GUIDED = None569 570 571def refine_mask(mask, rgb):572    """Snap a coarse SAM3 mask boundary to the RGB color edge (appearance-guided). This clips the same-depth573    bleed that depth CANNOT see: where the halo crosses onto a differently-colored neighbor, the color edge574    cuts it. Params are DERIVED, not tuned -- radius from frame scale, smoothness from the image's own color575    variance (same MAD-style principle as depth_edges). No model, no GPU. Falls back to grabCut, then erode.576    """577    if mask.shape[:2] != rgb.shape[:2]:578        mask = cv2.resize(579            mask.astype(np.uint8),580            (rgb.shape[1], rgb.shape[0]),581            interpolation=cv2.INTER_NEAREST,582        ).astype(bool)583    a = int(mask.sum())584    if a < 60:  # too small to refine meaningfully -> 1px erode as before585        me = cv2.erode(mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1).astype(bool)586        return me if me.any() else mask587    if _GUIDED is not None:588        r = max(589            4, int(round(0.02 * float(np.hypot(*mask.shape[:2]))))590        )  # radius ~2% of frame diagonal591        eps = (592            float(np.var(rgb.astype(np.float32) / 255.0)) * 0.01 + 1e-6593        )  # smoothness ~ image color variance594        soft = _GUIDED(rgb, mask.astype(np.float32), r, eps)595        out = soft > 0.5596        return out if out.sum() >= 0.4 * a else mask  # majority guard: don't over-carve597    # fallback (base opencv, no ximgproc): grabCut seeded FG=mask, PR_FG=dilated ring, BG=far exterior598    try:599        gc = np.full(mask.shape[:2], cv2.GC_PR_BGD, np.uint8)600        dil = cv2.dilate(mask.astype(np.uint8), np.ones((15, 15), np.uint8), 1).astype(601            bool602        )603        er = cv2.erode(mask.astype(np.uint8), np.ones((5, 5), np.uint8), 1).astype(bool)604        gc[dil] = cv2.GC_PR_FGD605        gc[er] = cv2.GC_FGD606        bgm = np.zeros((1, 65), np.float64)607        fgm = np.zeros((1, 65), np.float64)608        cv2.grabCut(rgb, gc, None, bgm, fgm, 3, cv2.GC_INIT_WITH_MASK)609        out = (gc == cv2.GC_FGD) | (gc == cv2.GC_PR_FGD)610        return out if 0.4 * a <= out.sum() <= 1.5 * a else mask611    except Exception:612        me = cv2.erode(mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1).astype(bool)613        return me if me.any() else mask614 615 616def backproject_frame(617    depth_f,618    intrinsics,619    c2w_f,620    mask_f,621    conf_f=None,622    conf_thr=0.0,623    valid_f=None,624    return_conf=False,625    edges_f=None,626):627    """Return (M,3) world points for the masked pixels of one frame (c2w_f = cam->world 4x4).628    valid_f: optional precomputed (isfinite & >0) depth mask, reused across all masks of a frame.629    return_conf: also return the (M,) DA3 confidence of each kept point (for per-point noise-aware cleaning).630    edges_f: optional per-frame depth-edge map; if given, keep the mask's largest depth-coherent component631             (cuts mask-bleed onto adjacent objects at the depth boundary, per frame, before back-projection).632    """633    height, width = depth_f.shape634    empty = (635        (np.empty((0, 3), np.float32), np.empty((0,), np.float32))636        if return_conf637        else np.empty((0, 3), np.float32)638    )639    if mask_f.shape != (height, width):640        mask_f = cv2.resize(641            mask_f.astype(np.uint8), (width, height), interpolation=cv2.INTER_NEAREST642        ).astype(bool)643    if valid_f is None:644        valid_f = np.isfinite(depth_f) & (depth_f > 0)645    m = mask_f & valid_f646    if conf_f is not None and conf_thr > 0:647        m &= conf_f >= conf_thr648    if (649        edges_f is not None and m.sum() >= 30650    ):  # keep the largest depth-coherent piece of the mask651        from scipy import ndimage652 653        lab, nlab = ndimage.label(m & ~edges_f)654        if nlab >= 1:655            sizes = np.bincount(lab.ravel())656            sizes[0] = 0657            big = int(sizes.argmax())658            if (659                sizes[big] >= 0.5 * m.sum()660            ):  # object is the MAJORITY piece (bleed is a minority)661                m = lab == big662    if not m.any():663        return empty664    ys, xs = np.nonzero(m)665    z = depth_f[ys, xs]666    if DEPTH_COHERENCE and len(z) >= 8:667        # an object is a coherent depth surface; floor/background BLEED pixels are depth outliers.668        # Drop them via the standard Tukey fence (1.5*IQR) on the masked region's depths -- parameter-free.669        q1, q3 = np.percentile(z, [25, 75])670        iqr = q3 - q1671        keep = (z >= q1 - 1.5 * iqr) & (z <= q3 + 1.5 * iqr)672        if keep.sum() >= 1:673            ys, xs, z = ys[keep], xs[keep], z[keep]674    fx, fy, cx, cy = (675        intrinsics[0, 0],676        intrinsics[1, 1],677        intrinsics[0, 2],678        intrinsics[1, 2],679    )680    camera_points = np.stack(681        [(xs - cx) * z / fx, (ys - cy) * z / fy, z], axis=1682    )  # camera coords683    world_points = (c2w_f[:3, :3] @ camera_points.T).T + c2w_f[:3, 3]  # -> world684    if return_conf:685        cw = (686            conf_f[ys, xs].astype(np.float32)687            if conf_f is not None688            else np.ones(len(ys), np.float32)689        )690        return world_points.astype(np.float32), cw691    return world_points.astype(np.float32)692 693 694# ==========================================================================================695# INSTANCE BUILDING -- oriented extent, 3D box-overlap re-identification, and the main696# build_instances() driver that turns per-frame masks into per-class 3D instances. Merged in697# from perceptual.py, verbatim.698# ==========================================================================================699 700 701def robust_centroid_extent(pts, up_axis=None):702    """median centroid + ORIENTED robust extent.703    If up_axis is given: YAW-ONLY oriented extent -- rotation is found by 2D PCA on the704    floor-projected points only, with the up axis left untouched. This matches VSI-Bench's705    own annotation convention for indoor scans (ScanNet/ARKitScenes OrientedBoundingBox706    objects follow a Manhattan-world assumption: rotated only around the vertical axis,707    never tilted). Unconstrained 3D PCA (the previous behavior) can chase noise on the708    vertical axis for flat/elongated objects and drift away from the true yaw.709    If up_axis is None (unknown at the call site): falls back to unconstrained 3D PCA.710    Either way: parameter-free, rotation-invariant in-plane, p2..p98 robust extent."""711    c = np.median(pts, axis=0)712    centered = pts - c713    if len(centered) > 5000:  # PCA on a sample (deterministic)714        centered = centered[np.random.RandomState(0).choice(len(centered), 5000, False)]715    if up_axis is not None:716        floor_axes = [i for i in range(3) if i != up_axis]717        floor_points = centered[:, floor_axes]718        try:719            _, _, floor_rotation = np.linalg.svd(720                floor_points - floor_points.mean(0), full_matrices=False721            )722            proj_floor = (723                floor_points @ floor_rotation.T724            )  # (N,2) along the object's own floor-plane axes725        except np.linalg.LinAlgError:726            proj_floor = floor_points727        up_col = centered[:, up_axis : up_axis + 1]  # up axis untouched (yaw-only)728        proj = np.concatenate([proj_floor, up_col], axis=1)729    else:730        try:731            _, _, rotation = np.linalg.svd(732                centered - centered.mean(0), full_matrices=False733            )734            proj = centered @ rotation.T  # coordinates along principal axes735        except np.linalg.LinAlgError:736            proj = centered737    lo = np.percentile(proj, 2, axis=0)738    hi = np.percentile(proj, 98, axis=0)739    ext = np.maximum(hi - lo, 0.0)740    dims = np.sort(ext)[::-1]  # the object's 3 oriented side lengths, longest first741    return c.astype(np.float32), float(dims[0]), dims742 743 744def _aabb(pts):745    """robust (p2..p98) axis-aligned 3D box of an instance's world points."""746    return np.percentile(pts, 2, axis=0), np.percentile(pts, 98, axis=0)747 748 749def merge_by_box_overlap(insts, up_axis=None):750    """Parameter-free 3D re-identification, with a temporal-exclusion gate.751 752    SAM3 (a 2D tracker) emits a NEW masklet each time the camera revisits an object, so one753    physical object -> several masklets at the same 3D location. We fuse two same-class masklets754    iff BOTH:755      (a) their measured 3D boxes overlap (or one centroid is inside the other) -- same volume, and756      (b) they NEVER appear in the same frame -- temporal exclusion.757    (b) is the key parameter-free invariant: two masklets co-visible in one frame were tracked by758    SAM3 as distinct objects in that frame, so they ARE distinct (e.g. two chairs around a table);759    we must never merge them, even if depth noise makes their boxes overlap. A revisit-duplicate,760    by contrast, lives in DISJOINT frames. Both tests are exact/measured -- NO tuned threshold.761    """762    n = len(insts)763    if n <= 1:764        return insts765    # Precompute the FULL pairwise box-match as one vectorized boolean matrix. With the per-frame766    # batched detector a class can have 100+ raw instances; the old O(n^3) loop called np.all-based767    # box_match millions of times. Here every pairwise test is one broadcast -> O(1) lookups below.768    lo = np.stack([_aabb(i["pts"])[0] for i in insts]).astype(np.float32)  # (n,3)769    hi = np.stack([_aabb(i["pts"])[1] for i in insts]).astype(np.float32)  # (n,3)770    cents = np.stack([i["centroid"] for i in insts]).astype(np.float32)  # (n,3)771    overlap = (hi[:, None, :] >= lo[None, :, :]).all(-1) & (772        hi[None, :, :] >= lo[:, None, :]773    ).all(-1)774    ins = (cents[:, None, :] >= lo[None, :, :]).all(-1) & (775        cents[:, None, :] <= hi[None, :, :]776    ).all(-1)777    match = (778        overlap | ins | ins.T779    )  # box_match[i,j]: same 3D volume (identical semantics to old code)780    # group-level agglomeration: merge two groups only if their COMBINED frame sets are disjoint781    # (so no two co-visible masklets ever land in one object) AND some cross-pair shares a 3D volume.782    groups = [783        {"members": [i], "frames": set(insts[i].get("frames", set()))} for i in range(n)784    ]785    changed = True786    while changed:787        changed = False788        for a in range(len(groups)):789            for b in range(a + 1, len(groups)):790                if (791                    groups[a]["frames"] & groups[b]["frames"]792                ):  # co-visible -> distinct objects793                    continue794                if match[np.ix_(groups[a]["members"], groups[b]["members"])].any():795                    groups[a]["members"] += groups[b]["members"]796                    groups[a]["frames"] |= groups[b]["frames"]797                    groups.pop(b)798                    changed = True799                    break800            if changed:801                break802    merged = []803    for g in groups:804        idxs = g["members"]805        pts = np.concatenate([insts[k]["pts"] for k in idxs], 0)806        cpts = np.concatenate(807            [808                insts[k].get("conf", np.ones(len(insts[k]["pts"]), np.float32))809                for k in idxs810            ],811            0,812        )813        best = max(814            (insts[k]["best_pts"] for k in idxs), key=len815        )  # largest single obs in the group816        c, longest, dims = robust_centroid_extent(817            best, up_axis818        )  # size+pos+3 oriented dims from best view (#2/#4)819        merged.append(820            {821                "centroid": c,822                "size": longest,823                "dims": dims,824                "first_time": min(insts[k]["first_time"] for k in idxs),825                "pts": pts,826                "conf": cpts,827                "best_pts": best,828                "n": sum(insts[k]["n"] for k in idxs),829                "nframes": len(g["frames"]),830            }831        )  # track persistence (evidence strength)832    return merged833 834 835def build_instances(836    per_class, depth, intr, c2w, conf, frame_times, frame_paths=None, up_axis=None837):838    """-> {class: [ {centroid(3), size(longest dim), first_time, npts} ]}839    up_axis: if known (from room_gravity, computed BEFORE this call), threads through to840    robust_centroid_extent for yaw-only oriented sizing. If None, size falls back to841    unconstrained 3D PCA."""842    out = {}843    stats = {}844    nframes = len(depth)845    # precompute each frame's valid-depth mask ONCE (was recomputed per mask -> per class).846    valid = {}847    edges = {}848    rgb = {}849    used_frames = {fi for frames in per_class.values() for fi in frames if fi < nframes}850    for fi in used_frames:851        valid[fi] = np.isfinite(depth[fi]) & (depth[fi] > 0)852        edges[fi] = depth_edges(depth[fi], valid[fi]) if DEPTH_EDGE_REFINE else None853        if MASK_REFINE and frame_paths and fi < len(frame_paths):854            im = cv2.imread(frame_paths[fi])  # BGR; guidedFilter/grabCut want 3ch uint8855            rgb[fi] = (856                cv2.resize(im, (depth[fi].shape[1], depth[fi].shape[0]))857                if im is not None858                else None859            )860    for cls, frames in per_class.items():861        # peak co-visibility: max distinct masklets SAM3 tracks SIMULTANEOUSLY in any one frame.862        # geometry-free, parameter-free, immune to revisit over-count; provable lower bound on count.863        peak = max((len(objs) for objs in frames.values()), default=0)864        # gather world points + first-seen time + frame set per obj_id (SAM3 track id = masklet)865        pts_by_id, conf_by_id, first_t, frames_by_id = {}, {}, {}, {}866        # px_by_id: (fidx, frame_time, mask_pixel_count) per obj_id, for appearance-order timing.867        # VSI-Bench's own GT defines "first appearance" as the timestamp where an object's pixel868        # count crosses a threshold (paper appendix B.1) -- NOT the first frame with any pixel at869        # all. A single stray mask-bleed/false-positive pixel would otherwise register as "first870        # seen" far too early. The threshold is self-calibrated per instance (that instance's OWN871        # median observed pixel count across its frames), same convention as _clean's median cut.872        px_by_id = {}873        for fidx, objs in frames.items():874            if fidx >= nframes:  # guard: SAM3 frame idx vs DA3 frames875                continue876            for oid, mask in objs.items():877                if MASK_REFINE and rgb.get(fidx) is not None:878                    me = (879                        mask880                        if os.environ.get("VSI_NO_REFINE") == "1"881                        else refine_mask(mask, rgb[fidx])882                    )  # appearance-guided boundary snap (clips same-depth bleed)883                else:884                    # erode 1px to drop mask-edge / background depth bleed885                    me = cv2.erode(886                        mask.astype(np.uint8), np.ones((3, 3), np.uint8), 1887                    ).astype(bool)888                    if not me.any():889                        me = mask890                conf_f = conf[fidx] if conf is not None else None891                conf_thr = (892                    np.percentile(conf_f, CONF_PCT)893                    if (conf_f is not None and CONF_PCT > 0)894                    else 0.0895                )896                world_points, cw = backproject_frame(897                    depth[fidx],898                    intr[fidx],899                    c2w[fidx],900                    me,901                    conf_f,902                    conf_thr=conf_thr,903                    valid_f=valid.get(fidx),904                    return_conf=True,905                    edges_f=edges.get(fidx),906                )907                if len(world_points):908                    pts_by_id.setdefault(oid, []).append(world_points)909                    conf_by_id.setdefault(oid, []).append(910                        cw911                    )  # per-point DA3 confidence (for _clean)912                    frames_by_id.setdefault(oid, set()).add(913                        fidx914                    )  # for co-occurrence gate915                    t = frame_times[fidx]916                    px_by_id.setdefault(oid, []).append((fidx, t, int(me.sum())))917        # appearance_order timing: per instance, first frame at/above its OWN median pixel count918        for oid, obs in px_by_id.items():919            counts = [c for _, _, c in obs]920            thresh = float(np.median(counts))921            crossing = [t for _, t, c in obs if c >= thresh]922            first_t[oid] = min(crossing) if crossing else min(t for _, t, c in obs)923        insts = []924        for oid, plist in pts_by_id.items():925            pts = np.concatenate(plist, 0)926            cpts = np.concatenate(conf_by_id[oid], 0)927            if len(pts) < MIN_INSTANCE_PTS:928                continue929            best = max(930                plist, key=len931            )  # #2/#4: largest single-frame observation (closest/most pixels)932            bc, bsize, bdims = robust_centroid_extent(933                best, up_axis934            )  # size + position + 3 oriented dims from best view935            insts.append(936                {937                    "centroid": bc,938                    "size": bsize,939                    "dims": bdims,940                    "first_time": first_t[oid],941                    "pts": pts,942                    "conf": cpts,943                    "best_pts": best,944                    "n": len(pts),945                    "frames": frames_by_id[oid],946                }947            )948        raw = len(insts)949        # 3D re-ID: fuse same-class masklets that occupy the same measured 3D volume (parameter-free).950        insts = merge_by_box_overlap(insts, up_axis)951        if insts:952            out[cls] = insts953            stats[cls] = {"raw": raw, "merged": len(insts), "peak": peak}954    return out, stats955 956 957# ==========================================================================================958# PER-CLASS SUMMARY + FLOOR AREA -- legacy array-schema row builder, and the room-scale floor959# area calculation. Merged in from perceptual.py, verbatim.960# ==========================================================================================961 962 963# ---- Per-class spatial code (legacy array-schema row) -------------------------------------964def class_spatial_code(insts, peak=0):965    cents = np.stack([i["centroid"] for i in insts], 0)  # (n,3)966    sizes = np.array([i["size"] for i in insts], np.float32)967    n = len(insts)  # merged centroids (for spatial stats)968    count = peak if peak else n  # reported count = peak co-visibility969    x, y, z = cents.mean(0)970    if n == 1:971        e1 = e2 = e3 = px = pz = 0.0972        size_iqr = 0.0973    else:974        cov = np.cov(cents.T)  # 3x3975        vals, vecs = np.linalg.eigh(cov)  # ascending976        order = np.argsort(vals)[::-1]977        vals = np.clip(vals[order], 0, None)978        vecs = vecs[:, order]979        e1, e2, e3 = vals.tolist()980        pv = vecs[:, 0]  # principal eigenvector981        px, pz = float(pv[0]), float(pv[2])982        q1, q3 = np.percentile(sizes, [25, 75])983        size_iqr = float(q3 - q1)984    size_median = float(np.median(sizes))985    first_time = float(min(i["first_time"] for i in insts))986    row = [x, y, z, e1, e2, e3, px, pz, size_median, size_iqr, first_time, count]987    row = [988        (lambda r: 0.0 if r == 0 else r)(round(float(v), 1)) for v in row989    ]  # kill -0.0990    row[-1] = int(count)991    return row992 993 994# ---- floor_area (full-scene min-Y points -> XZ convex hull) -------------------------------995def compute_floor_area(depth, intr, c2w, conf, sky, stride=8, up_vec=None):996    pts = []997    for f in range(depth.shape[0]):998        height, width = depth[f].shape999        ys, xs = np.mgrid[0:height:stride, 0:width:stride]1000        ys = ys.ravel()1001        xs = xs.ravel()1002        z = depth[f][ys, xs]1003        ok = np.isfinite(z) & (z > 0)1004        if sky is not None:1005            ok &= ~sky[f][ys, xs].astype(bool)1006        if conf is not None:1007            ok &= conf[f][ys, xs] >= np.percentile(conf[f], 40)1008        ys, xs, z = ys[ok], xs[ok], z[ok]1009        if not len(z):1010            continue1011        intrinsics = intr[f]1012        fx, fy, cx, cy = (1013            intrinsics[0, 0],1014            intrinsics[1, 1],1015            intrinsics[0, 2],1016            intrinsics[1, 2],1017        )1018        camera_points = np.stack([(xs - cx) * z / fx, (ys - cy) * z / fy, z], 1)1019        world_points = (c2w[f][:3, :3] @ camera_points.T).T + c2w[f][:3, 3]1020        pts.append(world_points.astype(np.float32))1021    if not pts:1022        return 0.01023    points = np.concatenate(pts, 0)1024    if up_vec is not None:1025        # VSI-faithful: area in the plane orthogonal to GRAVITY (RANSAC floor normal), like the1026        # benchmark's gravity-aligned GT meshes. Build an orthonormal in-plane basis (u, v).1027        g = np.asarray(up_vec, np.float64)1028        g /= np.linalg.norm(g) + 1e-121029        a = np.array([1.0, 0.0, 0.0]) if abs(g[0]) < 0.9 else np.array([0.0, 1.0, 0.0])1030        u = np.cross(g, a)1031        u /= np.linalg.norm(u)1032        v = np.cross(g, u)1033        all_floor_points = np.stack([points @ u, points @ v], 1)1034    else:1035        up = int(1036            np.argmin(points.max(0) - points.min(0))1037        )  # legacy: vertical = smallest-extent axis1038        floor_axes = [i for i in range(3) if i != up]1039        all_floor_points = points[:, floor_axes]1040    # VSI-Bench room-size definition = alpha-shape of the floor-plane point cloud (confirmed in their1041    # paper appendix). VSI does not publish the alpha value they use for their own GT mesh, so alpha=21042    # here is NOT a matched/verified constant -- it was chosen empirically for this pipeline's own1043    # (sparser) reconstructed point density. This is the one disclosed benchmark-adjacent tuned constant1044    # in the whole file; everything else is exact/derived or a generic, non-tuned statistical convention.1045    # (Falls back to enclosed-fill below if the alphashape package isn't available.)1046    floor_points = all_floor_points1047    lo = np.percentile(floor_points, 0.5, 0)1048    hi = np.percentile(floor_points, 99.5, 0)  # gentle clip (preserve room extent)1049    floor_points = floor_points[1050        (floor_points[:, 0] >= lo[0])1051        & (floor_points[:, 0] <= hi[0])1052        & (floor_points[:, 1] >= lo[1])1053        & (floor_points[:, 1] <= hi[1])1054    ]1055    if len(floor_points) < 10:1056        return 0.01057    try:1058        import alphashape1059 1060        idx = np.random.RandomState(0).choice(1061            len(floor_points), min(10000, len(floor_points))1062        )1063        return round(1064            float(alphashape.alphashape(floor_points[idx], alpha=2).area), 11065        )  # alpha=2 tuned for recon density1066    except Exception:1067        from scipy import ndimage1068 1069        res = 0.101070        ai = ((floor_points[:, 0] - floor_points[:, 0].min()) / res).astype(int)1071        bi = ((floor_points[:, 1] - floor_points[:, 1].min()) / res).astype(int)1072        grid = np.zeros((ai.max() + 3, bi.max() + 3), np.uint8)1073        grid[ai + 1, bi + 1] = 11074        grid = cv2.morphologyEx(grid, cv2.MORPH_CLOSE, np.ones((7, 7), np.uint8))1075        grid = ndimage.binary_fill_holes(grid).astype(np.uint8)1076        return round(float(grid.sum()) * res * res, 1)1077 1078 1079# ---------------------------------------------------------------------------1080 1081# ==========================================================================================1082# SPATIAL CODE ASSEMBLY -- the top-level entry point this whole file exists for:1083# build_spatial_code() calls everything above to turn already-computed depth/pose/masks into1084# the final spatial code dict. Room outline + JSON writer. This section was already in1085# geometric.py before the perceptual.py merge; build_spatial_code() below is updated to call1086# the geometry functions above DIRECTLY (no `pl.` prefix -- they're plain local functions now1087# that everything is in one file), same logic, unchanged otherwise.1088# ==========================================================================================1089 1090 1091def _room_outline(depth, intr, c2w, conf, bu, bv):1092    """(Currently unemitted -- the one spatial code shape has no room outline field; this1093    math is kept intact for reuse.) Room floor-boundary polygon from the SAME grid as1094    compute_floor_area: floor points1095    projected onto the shared gravity plane (bu, bv), 10cm grid, close 7x7, fill holes,1096    largest contour, 0.2m polygon simplification. Same (bu, bv) as objects/camera, so the1097    outline, object positions, and floor_area all live in one consistent frame."""1098    from scipy import ndimage1099 1100    pts, stride = [], 81101    for f in range(0, depth.shape[0], 3):1102        height, width = depth[f].shape1103        ys, xs = np.mgrid[0:height:stride, 0:width:stride]1104        ys = ys.ravel()1105        xs = xs.ravel()1106        z = depth[f][ys, xs]1107        ok = np.isfinite(z) & (z > 0)1108        if conf is not None:1109            ok &= conf[f][ys, xs] >= np.percentile(conf[f], 40)1110        ys, xs, z = ys[ok], xs[ok], z[ok]1111        if not len(z):1112            continue1113        intrinsics = intr[f]1114        camera_points = np.stack(1115            [1116                (xs - intrinsics[0, 2]) * z / intrinsics[0, 0],1117                (ys - intrinsics[1, 2]) * z / intrinsics[1, 1],1118                z,1119            ],1120            1,1121        )1122        pts.append(1123            ((c2w[f][:3, :3] @ camera_points.T).T + c2w[f][:3, 3]).astype(np.float32)1124        )1125    if not pts:1126        return []1127    world_points = np.concatenate(pts, 0)1128    points = np.stack(1129        [world_points @ bu, world_points @ bv], 11130    )  # gravity-plane projection (same bu,bv as objects/area)1131    lo = np.percentile(points, 0.5, 0)1132    hi = np.percentile(points, 99.5, 0)1133    points = points[1134        (points[:, 0] >= lo[0])1135        & (points[:, 0] <= hi[0])1136        & (points[:, 1] >= lo[1])1137        & (points[:, 1] <= hi[1])1138    ]1139    if len(points) < 10:1140        return []1141    res = 0.101142    x0, y0 = points[:, 0].min(), points[:, 1].min()1143    ai = ((points[:, 0] - x0) / res).astype(int)1144    bi = ((points[:, 1] - y0) / res).astype(int)1145    grid = np.zeros((ai.max() + 3, bi.max() + 3), np.uint8)1146    grid[ai + 1, bi + 1] = 11147    grid = cv2.morphologyEx(grid, cv2.MORPH_CLOSE, np.ones((7, 7), np.uint8))1148    grid = ndimage.binary_fill_holes(grid).astype(np.uint8)1149    cs, _ = cv2.findContours(grid, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)1150    if not cs:1151        return []1152    poly = cv2.approxPolyDP(max(cs, key=cv2.contourArea), 0.2 / res, True)[:, 0, :]1153    # cv2 contour points are (col=bi, row=ai) -> (floor_y, floor_x)1154    return [1155        {1156            "floor_x_meters": round(float((r - 1) * res + x0), 1),1157            "floor_y_meters": round(float((c - 1) * res + y0), 1),1158        }1159        for c, r in poly1160    ]1161 1162 1163def build_spatial_code_raw(depth, intr, c2w, conf, ftimes, per):1164    """Builds THE spatial code -- the one and only shape a spatial code has, everywhere1165    (on disk, in prompts, in this pipeline): unit-strings ("1.4 meters"), spaced keys1166    ("x coordinate"), per-instance position + longest dimension only, room "floor area",1167    "closest classes distance meters from" (rooted per class, distance + closeness rank),1168    and a flat earliest-first "appearance order" list of class names. There is no separate1169    raw/rendered split and no schema flag -- the old v1/v2 branching (VSI_CODE_V2) and the1170    raw intermediate form (floor_x_meters keys, bounding_box, dimensions_meters,1171    seen_in_video_frames, camera_trajectory, room.outline) are gone; every underlying VALUE1172    that survives is computed by exactly the same math as before, only the emitted fields1173    and their formatting changed."""1174    # emission-time class rename: VSI's questions say 'coat rack' while their annotations1175    # (and hence the SAM3 prompt + caches) say 'coat hanger' -- same object, their naming1176    # seam. The model sees questions, so emitted codes follow the question vocabulary.1177    class_aliases = {"coat hanger": "coat rack"}1178    per = {class_aliases.get(k, k): v for k, v in per.items()}1179    inst, stats = build_instances(per, depth, intr, c2w, conf, ftimes)1180    up_vec, up_ax = room_gravity(1181        depth, intr, c2w, conf1182    )  # gravity = RANSAC floor normal (VSI-faithful)1183    bu, bv, bg = _floor_basis(1184        up_vec1185    )  # shared gravity floor frame (bu,bv horizontal, bg up)1186    points = np.concatenate([i["pts"] for cl in inst.values() for i in cl], 0)1187    floor_level = _floor_level(points, bg)1188    fa = compute_floor_area(depth, intr, c2w, conf, None, up_vec=up_vec)1189    code = to_spatial_code(inst, stats, fa, up_ax, up_vec, floor_level)1190    cls = list(inst.keys())1191    class_first = {c: min(i["first_time"] for i in v) for c, v in inst.items()}1192 1193    # Keyed dict + integer ranks (not a sorted list): each question option becomes ONE1194    # direct key access, and "which is closest" = min over small integers -- the filtered1195    # list-scan and decimal comparison were the observed failure modes even on GT data.1196    # 2-decimal distances: 0.1m rounding costs up to ~17% relative error on sub-meter1197    # answers, which fails the strictest MRA thresholds even with perfect values.1198    ccf = {}1199    for a in cls:1200        ds = sorted(

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