From 99408de04a60a21ad25d79e4a729e25fd959eae8 Mon Sep 17 00:00:00 2001 From: Yogesh Chawla Date: Wed, 22 Jul 2026 15:01:50 -0700 Subject: [PATCH 01/42] Add apple tree collision-aware build script and camera capture rig apple_tree.py builds apple trees with fruit-aware soft collision avoidance enabled and prints each tree's bounding-box extent. apple_tree_cameras.py reuses build_apple_tree() to render a same-distance, fixed-elevation-arc camera rig (above/level/below) per tree, intended as input views for Gaussian splatting. Co-Authored-By: Claude Sonnet 5 --- apple_tree.py | 146 ++++++++++++++++++++++++++++++++++++++++++ apple_tree_cameras.py | 76 ++++++++++++++++++++++ 2 files changed, 222 insertions(+) create mode 100644 apple_tree.py create mode 100644 apple_tree_cameras.py diff --git a/apple_tree.py b/apple_tree.py new file mode 100644 index 0000000..96343e9 --- /dev/null +++ b/apple_tree.py @@ -0,0 +1,146 @@ +import sys + +from pyhelios import Context, PlantArchitecture, Visualizer +from pyhelios.types import RGBcolor, vec3 + + +def build_apple_tree(plantarch, position=None, age_days=365.0, build_parameters=None): + """Load the apple library model and build one tree instance. + + Args: + plantarch: PlantArchitecture instance bound to a Context. + position: Base of the trunk (default: origin). + age_days: Plant age in days. Older = larger / more developed. + build_parameters: Optional training overrides, e.g. + {'trunk_height': 0.8, 'num_scaffolds': 4, 'scaffold_angle': 40} + + Returns: + plant_id for the created tree. + """ + if position is None: + position = vec3(0, 0, 0) + + if build_parameters and build_parameters.get("trunk_height", 0.8) > 0.8: + raise ValueError( + "The built-in apple model requires trunk_height <= 0.80 m " + "(20 nodes at 0.04 m per internode)." + ) + + plantarch.loadPlantModelFromLibrary("apple") + plant_id = plantarch.buildPlantInstanceFromLibrary( + base_position=position, + age=age_days, + build_parameters=build_parameters, + ) + return plant_id + + +def visualize_trees(context, plantarch, plant_ids): + """Open one interactive window showing all built trees together.""" + all_uuids = [] + for plant_id in plant_ids: + uuids = plantarch.getAllPlantUUIDs(plant_id) + if not uuids: + raise RuntimeError(f"Plant {plant_id} has no geometry to visualize") + all_uuids.extend(uuids) + + if not all_uuids: + raise RuntimeError("No tree geometry to visualize") + + # Fit the camera to the combined scene bounding box + x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(all_uuids) + center = vec3( + 0.5 * (x_bounds.x + x_bounds.y), + 0.5 * (y_bounds.x + y_bounds.y), + 0.5 * (z_bounds.x + z_bounds.y), + ) + extent = max( + x_bounds.y - x_bounds.x, + y_bounds.y - y_bounds.x, + z_bounds.y - z_bounds.x, + 1.0, + ) + + with Visualizer(width=1000, height=800) as visualizer: + visualizer.buildContextGeometry(context, uuids=all_uuids) + visualizer.setCameraPosition( + position=vec3( + center.x + 1.8 * extent, + center.y + 1.8 * extent, + center.z + 0.6 * extent, + ), + lookAt=center, + ) + visualizer.setBackgroundColor(RGBcolor(0.70, 0.85, 1.0)) + visualizer.setLightingModel("phong_shadowed") + visualizer.plotInteractive() + + +if __name__ == "__main__": + # Age is the main "size" control for library plants (days). + # Try ~90 for a young tree, ~365 for a year-old tree. + AGE_DAYS = [720.0] * 10 + + POSITIONS = [ + vec3(0, 0, 0), + vec3(1.5, 0, 0), + vec3(3, 0, 0), + vec3(9, 0, 0), + vec3(12, 0, 0), + vec3(15, 0, 0), + vec3(18, 0, 0), + vec3(21, 0, 0), + vec3(24, 0, 0), + vec3(27, 0, 0), + ] + + # Optional apple training-system overrides (omit to use defaults): + BUILD_PARAMS_LIST = [ + {"trunk_height": 0.72, "num_scaffolds": 4, "scaffold_angle": 33}, + {"trunk_height": 0.78, "num_scaffolds": 6, "scaffold_angle": 47}, + {"trunk_height": 0.79, "num_scaffolds": 5, "scaffold_angle": 39}, + {"trunk_height": 0.74, "num_scaffolds": 4, "scaffold_angle": 44}, + {"trunk_height": 0.70, "num_scaffolds": 6, "scaffold_angle": 31}, + {"trunk_height": 0.80, "num_scaffolds": 5, "scaffold_angle": 50}, + {"trunk_height": 0.76, "num_scaffolds": 4, "scaffold_angle": 36}, + {"trunk_height": 0.77, "num_scaffolds": 5, "scaffold_angle": 42}, + {"trunk_height": 0.71, "num_scaffolds": 6, "scaffold_angle": 35}, + {"trunk_height": 0.74, "num_scaffolds": 5, "scaffold_angle": 48}, + ] + + with Context() as context: + with PlantArchitecture(context) as plantarch: + plantarch.loadPlantModelFromLibrary("apple") + + # Include fruit in collision detection and steer growth away from + # overlaps (fruit is excluded by default) so apples on neighboring + # trees/branches don't grow into each other. + plantarch.setCollisionRelevantOrgans( + include_internodes=True, + include_leaves=True, + include_fruit=True, + ) + plantarch.enableSoftCollisionAvoidance(enable_fruit_collision=True) + + plant_ids = [] + + for i in range(3): + plant_id_i = build_apple_tree( + plantarch, + position=POSITIONS[i], + age_days=AGE_DAYS[i], + build_parameters=BUILD_PARAMS_LIST[i], + ) + plant_ids.append(plant_id_i) + + object_ids = plantarch.getAllPlantObjectIDs(plant_id_i) + uuids = plantarch.getAllPlantUUIDs(plant_id_i) + x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(uuids) + print(f"Apple tree plant_id={plant_id_i}") + print(f" objects: {len(object_ids)}") + print(f" primitives: {len(uuids)}") + print(f" age: {AGE_DAYS[i]} days") + print(f" extent: x={x_bounds.y - x_bounds.x:.3f}m, y={y_bounds.y - y_bounds.x:.3f}m, z={z_bounds.y - z_bounds.x:.3f}m") + + if "--no-visualization" not in sys.argv: + visualize_trees(context, plantarch, plant_ids) diff --git a/apple_tree_cameras.py b/apple_tree_cameras.py new file mode 100644 index 0000000..c07410c --- /dev/null +++ b/apple_tree_cameras.py @@ -0,0 +1,76 @@ +import os +import sys + +from apple_tree import build_apple_tree +from pyhelios import Context, PlantArchitecture, Visualizer +from pyhelios.types import RGBcolor, vec3 + +# Gaussian-splat capture rig: all 3 cameras sit at the same horizontal +# distance in front of the tree (same x, same y offset) and differ only in +# elevation (z) -- one above the tree looking down, one level with the +# tree's center, one low near the ground looking up. The same rig is +# reused for every tree, just re-centered on that tree's position/height. +CAMERA_DISTANCE = 2.5 # front offset from the tree, in meters, same for all rigs + +CAMERA_RIGS = [ + {"name": "above", "height_frac": 1.15}, + {"name": "level", "height_frac": 0.5}, + {"name": "below", "height_frac": 0.05}, +] + +OUTPUT_DIR = "renders" + + +def camera_position_for_tree(context, plantarch, plant_id, base_position, rig): + """Compute a camera position/lookAt in front of one tree for a given rig.""" + uuids = plantarch.getAllPlantUUIDs(plant_id) + x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(uuids) + tree_height = z_bounds.y - z_bounds.x + + look_at = vec3(base_position.x, base_position.y, z_bounds.x + 0.5 * tree_height) + position = vec3( + base_position.x, + base_position.y - CAMERA_DISTANCE, + z_bounds.x + rig["height_frac"] * tree_height, + ) + return position, look_at + + +if __name__ == "__main__": + AGE_DAYS = 365.0 + POSITIONS = [vec3(0, 0, 0), vec3(3, 0, 0), vec3(6, 0, 0)] + + os.makedirs(OUTPUT_DIR, exist_ok=True) + + with Context() as context: + with PlantArchitecture(context) as plantarch: + plant_ids = [] + for i, position in enumerate(POSITIONS): + plant_id = build_apple_tree(plantarch, position=position, age_days=AGE_DAYS) + plant_ids.append(plant_id) + print(f"Built apple tree plant_id={plant_id} at {position}") + + all_uuids = [] + for plant_id in plant_ids: + all_uuids.extend(plantarch.getAllPlantUUIDs(plant_id)) + + with Visualizer(width=1000, height=800, headless=True) as visualizer: + visualizer.buildContextGeometry(context, uuids=all_uuids) + visualizer.setBackgroundColor(RGBcolor(0.70, 0.85, 1.0)) + visualizer.setLightingModel("phong_shadowed") + + for plant_id, position in zip(plant_ids, POSITIONS): + for rig in CAMERA_RIGS: + camera_pos, look_at = camera_position_for_tree( + context, plantarch, plant_id, position, rig + ) + visualizer.setCameraPosition(position=camera_pos, lookAt=look_at) + visualizer.plotUpdate() + + filename = os.path.join( + OUTPUT_DIR, f"tree{plant_id}_{rig['name']}.png" + ) + visualizer.printWindow(filename) + print(f" Saved {filename}") + + print(f"\nDone. Images written to {OUTPUT_DIR}/") From 4f9885b42471d0367d3931cc9cc0dd16872ba59c Mon Sep 17 00:00:00 2001 From: Yogesh Chawla Date: Tue, 28 Jul 2026 13:26:38 -0700 Subject: [PATCH 02/42] Add Gaussian Splatting pipeline: per-tree camera rig, fruit-only, dual strategy Renders each tree from its own frontal-plane camera grid instead of an orchard-wide orbit, trains apple-only (fruit-class) Gaussians via gsplat, and runs both DefaultStrategy and MCMCStrategy on the same seed points so the two can be compared side by side. Co-Authored-By: Claude Sonnet 5 --- GAUSSIAN_SPLATTING.md | 267 +++++++++++++ apple_tree_gaussian_splatting.py | 633 +++++++++++++++++++++++++++++++ 2 files changed, 900 insertions(+) create mode 100644 GAUSSIAN_SPLATTING.md create mode 100644 apple_tree_gaussian_splatting.py diff --git a/GAUSSIAN_SPLATTING.md b/GAUSSIAN_SPLATTING.md new file mode 100644 index 0000000..e107111 --- /dev/null +++ b/GAUSSIAN_SPLATTING.md @@ -0,0 +1,267 @@ +# Apple Tree → Gaussian Splatting: automation runbook + +This documents how `apple_tree_gaussian_splatting.py` was built, so the same +project can be reproduced (or extended) as a scripted process instead of an +interactive back-and-forth. It captures the environment setup, the API facts +that had to be discovered by reading source/testing, the one serious native +bug that was hit and how it was isolated, and the design decisions that +needed a human call vs. ones that could be defaulted automatically. + +Source: `apple_tree.py` (tree growth) + `apple_tree_cameras.py` (camera rig +pattern) → `apple_tree_gaussian_splatting.py` (full pipeline). + +--- + +## 1. Objective + +Grow apple trees with PyHelios, render a multi-view dataset with known camera +poses, and train a real 3D Gaussian Splat (via the `gsplat` library) on that +dataset — seeding the Gaussians from the tree's own mesh geometry instead of +COLMAP/SfM, since PyHelios already knows exact 3D positions. Training is +restricted to the fruit class only (the semantic mask acts as a prior, see +§5), and is run once per densification strategy `gsplat` ships, so the two +can be compared on identical seed points/cameras. + +## 2. Prerequisites check (do this first, automatable) + +```bash +nvidia-smi --query-gpu=name,driver_version,memory.total --format=csv +find /usr/local -maxdepth 1 -iname "cuda*" # look for a CUDA toolkit even if nvcc isn't on PATH +/usr/local/cuda-*/bin/nvcc --version # confirm nvcc directly, not just `which nvcc` +gcc --version && g++ --version # gsplat JIT-compiles CUDA kernels, needs a working host compiler +python3 -c "import pyhelios" 2>&1 # confirm pyhelios' compiled extension already works in the base env +``` + +Key facts learned this way: +- `pyhelios/_stub.cpython-313-*.so` is compiled for a **specific Python ABI** + (3.13 here). Any new env must match that Python minor version exactly, or + the import fails. Don't assume `pyhelios` is pip-installable — it's used + via `sys.path`/cwd from the repo root, not a package registered in `pip + show`. +- `nvcc` may be installed but simply not on `PATH` (common on machines set up + for multiple CUDA versions). Check `/usr/local/cuda-*/bin` before assuming + it's missing. + +## 3. Environment setup (scripted) + +Do **not** install torch/gsplat into the base env — create a dedicated env +pinned to the same Python version as PyHelios' compiled extension: + +```bash +conda create -n gsplat python=3.13 -y +conda activate gsplat +pip install numpy pyyaml pillow scipy # pyhelios runtime deps + scipy for KDTree init + +export PATH=/usr/local/cuda-12.9/bin:$PATH # match whatever toolkit `find /usr/local` found +export CUDA_HOME=/usr/local/cuda-12.9 +pip install --index-url https://download.pytorch.org/whl/cu128 torch torchvision # cuXXX must support the GPU's arch (e.g. cu128 for Blackwell/RTX 50-series) +pip install gsplat # JIT-compiles CUDA kernels on first use, needs nvcc from above +``` + +Validation gate before writing any pipeline code — run all three in the new +env and require all to pass: + +```bash +python3 -c "import torch; assert torch.cuda.is_available(); print(torch.cuda.get_device_name(0))" +python3 -c "import gsplat; from gsplat import rasterization, DefaultStrategy, export_splats" +PYTHONPATH=/home/yogesh/PyHelios python3 -c "from pyhelios import Context; Context()" +``` +Then confirm pyhelios and torch/gsplat coexist in **one process** (headless +`Visualizer` + `torch.cuda` tensor + `gsplat` import in sequence) before +building anything on top — this is cheap insurance against a GL/CUDA context +conflict that would otherwise surface much later. + +## 4. API surface that had to be discovered (not obvious from names alone) + +Recorded here so it doesn't need re-discovery. Verify against the installed +version before trusting it — these are pinned to PyHelios (this repo, current +HEAD) and `gsplat==1.5.3`. + +**PyHelios `Visualizer`** (`pyhelios/Visualizer.py`): +- `setCameraPosition(position, lookAt)` — up vector is hardcoded to world + `+Z`, not settable. Z-up world convention throughout Helios. +- `setCameraFieldOfView(angle_FOV)` — **vertical** FOV in degrees, default 45. +- No aspect-ratio, near/far, or orthographic controls. Aspect = `width/height` + automatically, which combined with a single vertical FOV means `fx == fy` + (square pixels) always — you cannot reproduce a real camera's `fx != fy`. +- `printWindow()` output resolution == constructor `width`/`height` **only** + in `headless=True` mode (windowed mode can differ due to DPI scaling). +- No getter for view/projection matrices — must be reconstructed manually + (see §5). + +**PyHelios `Context`** (`pyhelios/Context.py`): +- `getPrimitiveVertices(uuid_list)` → `(flat_float32_array, offsets_uint32)`; + vertices for primitive `i` are `flat[offsets[i]:offsets[i+1]].reshape(-1,3)`. +- `getPrimitiveColor(uuid_list)` → Nx3 array. **Textured primitives (leaves, + bark, fruit) return `(0,0,0)`** — flat color is only meaningful for + untextured primitives. Check `getPrimitiveTextureFile(uuid)` and average + the texture image (alpha-aware) as a fallback; only 2-3 unique texture + files exist per plant species, so cache by path. +- `getDomainBoundingSphere(uuids)` / `getDomainBoundingBox(uuids)` — call + per-tree (not once for the whole orchard) to size each tree's own camera + rig off its own footprint (see §7); a shared orchard-wide box badly misses + individual trees that sit off-center in a multi-tree row. + +**`gsplat` 1.5.3**: +- `rasterization(means, quats, scales, opacities, colors, viewmats, Ks, + width, height, sh_degree, backgrounds, ...)` → `(render[C,H,W,3], + alpha[C,H,W,1], info_dict)`. `viewmats` is **world-to-camera, OpenCV + convention** (X-right, Y-down, Z-forward) — must match how the training + images were actually formed, not an idealized convention. +- Two interchangeable densification strategies, same `rasterization()` calls + but **different constructor args and hook signatures**: + - `DefaultStrategy` — clone/split/prune (original 3DGS). `params` dict + needs at least `means`, `scales` (log-space), `quats`, `opacities` + (logit-space); extra keys (`sh0`, `shN`) ride along generically as long + as one Adam optimizer per param exists in a matching `optimizers` dict. + `initialize_state(scene_scale=...)`; `step_pre_backward` and + `step_post_backward` both take `(params, optimizers, state, step, info)`. + - `MCMCStrategy` — relocation + noise injection (3DGS-MCMC paper), capped + at `cap_max` Gaussians. `initialize_state()` takes **no** `scene_scale` + arg (relocation is opacity-driven, not scale-driven). `step_pre_backward` + is a no-op (don't call it). `step_post_backward` takes an **extra + required `lr` kwarg** (the current means learning rate) — omitting it is + a `TypeError`, not a silent no-op. The paper also expects explicit L1 + opacity/scale regularization added to the loss, since there's no + split/prune step to keep them in check otherwise. +- `export_splats(means, scales, quats, opacities, sh0, shN, format="ply", + save_to=path)` — writes a standard 3DGS PLY directly, no manual PLY + encoding needed. Works identically regardless of which strategy trained + the params. + +## 5. Design decisions that were resolved manually (defaults for next time) + +These required a judgment call. Recording the answer chosen so a future +automated run can default to it without asking: + +| Decision | Choice made | Why | +|---|---|---| +| Scope of "implement gaussian splatting" | Full training pipeline via `gsplat`, not a from-scratch rasterizer, not export-only | User wanted actual trained output; `gsplat` is the standard CUDA-accelerated library, reimplementing the rasterizer adds risk with no benefit | +| Point cloud init | Sample directly from tree mesh surface (bilinear on patches, barycentric on triangles, area-proportional density) instead of COLMAP/SfM or random init | PyHelios already has exact geometry — this is strictly better-seeded than SfM sparse points, and avoids needing COLMAP at all | +| Camera pose convention | Manually derive OpenCV world-to-camera matrix (X-right, Y-down, Z-forward) from `eye`/`lookAt`/world-up=`(0,0,1)`, verified numerically (`right × down == forward`) | `gsplat` needs OpenCV-convention viewmats; Helios exposes no matrix getter, so it must be reconstructed to exactly match what Helios's internal `glm::lookAt` actually rendered | +| Real camera intrinsics (`fx,fy,cx,cy` + distortion) given by user | Derive render resolution from `(cx,cy)`, vertical FOV from `fy`; drop distortion and off-center principal point; gsplat's K matrix is a **centered** pinhole matching what was actually rendered | Helios's renderer physically cannot produce `fx != fy`, off-center principal point, or lens distortion — passing the raw values into gsplat while training on Helios pixels would be geometrically inconsistent with the images themselves | +| Loss function | `0.8 * L1 + 0.2 * (1 - SSIM)`, compact from-scratch SSIM (11×11 Gaussian window, `conv2d`) | Standard 3DGS recipe; no extra dependency needed for SSIM | +| Camera rig | Treat each tree as independent: per-tree frontal-plane grid (`NUM_GRID_COLS x NUM_GRID_ROWS` positions on one plane at a fixed distance from *that tree's own* center, always looking at it) instead of one 360° orbit around the whole orchard | A single orbit around the combined bounding box over/under-frames trees that sit off-center in a multi-tree row; a fixed per-tree frontal plane keeps every tree consistently framed regardless of where it sits, without needing to circle it | +| Target class for training | Restrict to `TARGET_CLASS = "fruit"` always (not a toggle) — seed points, dataset masking, and export are all apples-only | User wants the deliverable to be an apple-only splat, not the whole tree; the semantic mask (§ render_semantic_masks) already exists as a clean prior for this | +| Splat densification strategy | Train **once per strategy** gsplat ships (`DefaultStrategy`, `MCMCStrategy`) on the same seed points/dataset, exporting a separate `.ply` per strategy rather than picking one | Neither strategy is a strict upgrade — comparing both on the same data is cheap (just a second training pass) and avoids committing to one without evidence | + +## 6. Known bug: PyHelios headless `Visualizer` heap corruption at large resolution + +**Symptom:** `malloc(): corrupted top size`, process aborts (SIGABRT/core +dump), no Python traceback — a native heap corruption, not a Python +exception. + +**Reproduction:** occurs during multi-view rendering once a render dimension +gets large (observed failures from ~1850px up; the original camera +calibration implied 1957×1286). Always survives the *first* render in a +process; crashes on the *second or later* `plotUpdate()`/`printWindow()` +call. + +**Isolation steps that ruled out other causes** (run each independently to +confirm before assuming the fix generalizes to a different PyHelios version): +1. 3 trees / 30K primitives at large res → crashes on view 2. +2. 1 tree / ~10K primitives at the *same* large res → **still crashes** on + view 2 → not about primitive/tree count. +3. Same test with a **fresh `Visualizer` object per view** (instead of + reusing one instance across the loop) → **still crashes** → not about + object reuse/state. +4. Same test with **no `torch`/`gsplat` import anywhere in the process** → + **still crashes** → not a CUDA/allocator interaction with PyTorch; it's + internal to PyHelios's own headless rendering path. +5. Bisected resolution: 800×526, 1200×789, 1400×919, 1600×1051, 1800×1182 + all survived a 4-view sequential test; 1850×1215 and 1957×1286 crashed. + The boundary was **not perfectly monotonic** across separate process runs + (consistent with genuine heap corruption whose crash point depends on + allocator/heap layout, not a clean hardcoded limit) — treat any threshold + found this way as approximate, not exact. + +**Workaround implemented:** cap the render resolution to `MAX_RENDER_DIMENSION += 1000` (comfortably below the entire observed failure zone), scaling both +width and height from the calibration-derived resolution while preserving +aspect ratio. Validated stable across 24 sequential renders in the actual +`render_dataset()` code path before trusting it in the full pipeline. FOV is +unaffected by this scaling since it's an angle, not a pixel count. + +**If this shows up again in a different context:** re-run the bisection in +§6 rather than assuming the same numeric threshold — this is a bug in the +native `Visualizer`/Helios core (likely a framebuffer resize or PNG-encode +buffer sizing bug), not something fixable from the Python wrapper beyond +avoiding the trigger condition. Filing it upstream against Helios would be +the real fix. + +## 7. Pipeline architecture + +``` +apple_tree_gaussian_splatting.py +├── build_orchard() # reuses build_apple_tree() from apple_tree.py +├── sample_point_cloud() # mesh-surface sampling + texture-average color fallback +├── plane_camera_poses() # per-tree frontal-plane grid (NUM_GRID_COLS x NUM_GRID_ROWS), +│ # called once per tree (extends the apple_tree_cameras.py +│ # 3-shot rig; no orbiting, fixed distance per tree) +├── look_at_view_matrix() # OpenCV-convention world-to-camera matrix +├── intrinsics_matrix() # K from (width, height, vertical FOV) +├── render_dataset() # headless multi-view render + transforms.json export +├── classify_primitives() # fruit/leaf/tree split via Helios "object_label" data +├── render_semantic_masks() # flat-colored mask render -> class-index PNG per view +├── init_gaussians() # KDTree-based scale init, RGB->SH0 color init +├── train_gaussians() # gsplat rasterization + L1/SSIM loss, strategy-agnostic +├── _build_strategy() # DefaultStrategy or MCMCStrategy, per SPLAT_STRATEGIES +├── evaluate() # held-out PSNR + side-by-side comparison PNGs +└── export_ply() # gsplat.export_splats -> one .ply per strategy + +main() trains TARGET_CLASS="fruit" once per entry in SPLAT_STRATEGIES, +reusing the same rendered dataset/seed points for each. +``` + +## 8. Validation checklist before trusting a full run + +Always smoke-test at reduced scale first — this is what caught the crash in +§6 before it wasted a full 5000-iteration run: + +1. 1 tree, ~16 views, small resolution (e.g. 200×200), ~300-700 iterations. + Confirms: tree build → render → point sampling → training forward/backward + → `DefaultStrategy` density-control path (needs `iters > refine_start_iter`, + default 500) → eval → PLY export, all with no exceptions. +2. Visually inspect one ground-truth render and one eval comparison PNG + (`renders/.../eval//eval_*.png`, side-by-side GT|rendered) — + confirms camera orientation is correct (not flipped/mirrored) and colors + are sane, which numeric checks alone won't catch. Do this for **both** + strategies, not just the first — a strategy-specific bug (e.g. the + `MCMCStrategy` missing `lr` kwarg) can silently corrupt only one output. +3. Only then run the full default (3 trees, `NUM_GRID_COLS x NUM_GRID_ROWS` + views per tree, fruit class only, 5000 iterations x 2 strategies). + +## 9. Reference full run (validates the fix in §6) + +The numbers below were measured before the per-tree frontal-plane rig, +fruit-only restriction, and dual-strategy loop (§5) were added — they +validate the §6 workaround, not the current camera/strategy behavior. +**Re-run and replace these once the current pipeline has been executed +end-to-end** rather than trusting these figures for the new code path: + +``` +3 trees, 72 views (3 elevations x 24 azimuths, whole-tree class), 1000x657 +renders, 5000 iterations, DefaultStrategy only, on an RTX 5090: +30,700 seed points -> 62,767 Gaussians after density control +Held-out PSNR: 27.48 dB over 9 views +Output: renders/gaussian_splatting/apple_orchard_splats.ply (3.5 MB) +``` + +## 10. To run it + +```bash +conda activate gsplat +cd /home/yogesh/PyHelios +python3 apple_tree_gaussian_splatting.py +``` + +Trains `TARGET_CLASS="fruit"` once per entry in `SPLAT_STRATEGIES` +(`default`, `mcmc`). Output: + +``` +renders/gaussian_splatting/ +├── dataset/ # shared RGB + mask renders, one set for both strategies +├── eval/{default,mcmc}/ # per-strategy held-out comparison PNGs +├── apple_splats_fruit_default.ply +└── apple_splats_fruit_mcmc.ply +``` diff --git a/apple_tree_gaussian_splatting.py b/apple_tree_gaussian_splatting.py new file mode 100644 index 0000000..4d32214 --- /dev/null +++ b/apple_tree_gaussian_splatting.py @@ -0,0 +1,633 @@ +import functools +import json +import math +import os + +import numpy as np +import torch +import torch.nn.functional as F +from gsplat import DefaultStrategy, MCMCStrategy, export_splats, rasterization +from PIL import Image +from scipy.spatial import cKDTree +from torch import nn + +from apple_tree import build_apple_tree +from notify_slack import notify_slack +from pyhelios import Context, PlantArchitecture, Visualizer +from pyhelios.types import RGBcolor, vec3 + +DEVICE = "cuda" if torch.cuda.is_available() else "cpu" +SH0_C0 = 0.28209479177387814 # degree-0 spherical harmonic basis constant + +# Reference calibration from a real camera. PyHelios's Visualizer only exposes a +# symmetric pinhole camera (single vertical FOV, centered principal point, no +# fx != fy, no distortion), so this can't be reproduced exactly. We instead +# derive the render resolution from (cx, cy) and the vertical FOV from fy, and +# drop distortion + off-center principal point -- the gsplat training K matrix +# below is a centered pinhole matching what PyHelios actually rendered, not +# these raw values. +CAMERA_FX = 737.052978515625 +CAMERA_FY = 736.5230102539062 +CAMERA_CX = 978.52197265625 +CAMERA_CY = 643.114013671875 +CAMERA_DISTORTION = (1.1389000415802002, 1.6283799409866333, + 0.0002682750055100769, 2.426269929856062e-05, 0.17436400055885315) # k1,k2,p1,p2,k3 (unused, see note above) + +BACKGROUND_RGB = (0.70, 0.85, 1.0) +FOV_DEG = math.degrees(2 * math.atan(round(2 * CAMERA_CY) / (2 * CAMERA_FY))) + +# PyHelios's headless Visualizer has a native heap-corruption bug (crashes with +# "malloc(): corrupted top size") triggered by repeated renders once a dimension +# gets much above ~1800px -- confirmed independent of primitive count, of reusing +# vs. recreating the Visualizer, and of whether torch/gsplat are even imported. +# The full resolution implied by the camera calibration (1957x1286) crashes +# reliably within the first few views, so we render at the same aspect ratio +# scaled down to a size that's been validated stable across dozens of sequential +# renders. FOV_DEG above is an angle derived from fy/cy and stays correct +# regardless of this scaling. +MAX_RENDER_DIMENSION = 1000 +_render_scale = min(1.0, MAX_RENDER_DIMENSION / max(round(2 * CAMERA_CX), round(2 * CAMERA_CY))) +IMAGE_WIDTH = round(2 * CAMERA_CX * _render_scale) +IMAGE_HEIGHT = round(2 * CAMERA_CY * _render_scale) +NUM_GRID_COLS = 8 # horizontal camera positions within the frontal plane +NUM_GRID_ROWS = 3 # vertical camera positions within the frontal plane +GRID_SPAN_FRACTION = 0.6 # how far the grid spans across the tree's own extent/height +TEST_EVERY = 8 +DISTANCE_MARGIN = 1.05 # tight headroom so each tree nearly fills the frame + +OUTPUT_DIR = "renders/gaussian_splatting" +DATASET_DIR = os.path.join(OUTPUT_DIR, "dataset") +EVAL_DIR = os.path.join(OUTPUT_DIR, "eval") + +TRAIN_ITERS = 5000 +INIT_POINT_BUDGET = 60000 + +CLASS_NAMES = ("fruit", "leaf", "tree") # "tree" = everything else: branches, trunk, petioles, flower parts +CLASS_MASK_COLORS = {"fruit": (255, 0, 0), "leaf": (0, 255, 0), "tree": (0, 0, 255)} +MASK_BACKGROUND_RGB = (0, 0, 0) +MASK_COLOR_MATCH_THRESHOLD = 60 # pixels farther than this from every class color fall back to background + +TARGET_CLASS = "fruit" # apples only, using the semantic mask as a prior + +# Design decision: rather than committing to one gsplat densification scheme, +# train once per strategy the library ships and export a separate .ply for +# each, so the two can be compared on the same seed points/dataset: +# - "default": DefaultStrategy -- clone/split/prune (original 3DGS paper). +# - "mcmc": MCMCStrategy -- relocation + noise injection (3DGS-MCMC paper). +SPLAT_STRATEGIES = ("default", "mcmc") +MCMC_CAP_MAX = 200_000 # relocation target cap, sized for a fruit-only point budget +MCMC_OPACITY_REG = 0.01 # L1 opacity regularization weight, per the 3DGS-MCMC paper +MCMC_SCALE_REG = 0.01 # L1 scale regularization weight, per the 3DGS-MCMC paper + + +# --------------------------------------------------------------------------- +# Tree construction (reuses apple_tree.py) + point-cloud seeding +# --------------------------------------------------------------------------- + +def build_orchard(plantarch, positions, age_days): + """Build one apple tree per position, with fruit-aware collision avoidance.""" + plantarch.loadPlantModelFromLibrary("apple") + plantarch.setCollisionRelevantOrgans( + include_internodes=True, include_leaves=True, include_fruit=True + ) + plantarch.enableSoftCollisionAvoidance(enable_fruit_collision=True) + return [ + build_apple_tree(plantarch, position=position, age_days=age_days) + for position in positions + ] + + +@functools.lru_cache(maxsize=None) +def _texture_average_color(path): + """Mean RGB of a texture's opaque pixels, used when a primitive has no flat color.""" + img = Image.open(path).convert("RGBA") + arr = np.asarray(img, dtype=np.float32) / 255.0 + rgb = arr[..., :3] + mask = arr[..., 3] > 0.5 + if not mask.any(): + mask = np.ones(arr.shape[:2], dtype=bool) + return tuple(rgb[mask].mean(axis=0)) + + +def sample_point_cloud(context, uuids, target_points=INIT_POINT_BUDGET, + samples_per_unit_area=500.0, max_samples=6): + """Seed Gaussian positions/colors directly from tree surface geometry.""" + flat_verts, offsets = context.getPrimitiveVertices(uuids) + colors = context.getPrimitiveColor(uuids) + areas = context.getPrimitiveArea(uuids) + + positions, point_colors = [], [] + for i, uuid in enumerate(uuids): + verts = flat_verts[offsets[i]:offsets[i + 1]].reshape(-1, 3) + if verts.shape[0] < 3: + continue + + color = colors[i] + if color.max() < 0.03: + texture_file = context.getPrimitiveTextureFile(uuid) + if texture_file: + color = np.array(_texture_average_color(texture_file), dtype=np.float32) + + n_samples = int(np.clip(round(areas[i] * samples_per_unit_area), 1, max_samples)) + for _ in range(n_samples): + if verts.shape[0] >= 4: + u, w = np.random.rand(2) + p = ((1 - u) * (1 - w) * verts[0] + u * (1 - w) * verts[1] + + u * w * verts[2] + (1 - u) * w * verts[3]) + else: + r1, r2 = np.random.rand(2) + if r1 + r2 > 1: + r1, r2 = 1 - r1, 1 - r2 + p = verts[0] + r1 * (verts[1] - verts[0]) + r2 * (verts[2] - verts[0]) + positions.append(p) + point_colors.append(color) + + positions = np.asarray(positions, dtype=np.float32) + point_colors = np.clip(np.asarray(point_colors, dtype=np.float32), 0.02, 0.98) + if len(positions) > target_points: + idx = np.random.choice(len(positions), target_points, replace=False) + positions, point_colors = positions[idx], point_colors[idx] + return positions, point_colors + + +def classify_primitives(context, uuids): + """Split primitive UUIDs into fruit/leaf/tree using Helios's per-primitive + "object_label" data (set internally by PlantArchitecture when it builds + each organ). "tree" is everything left over: branches, trunk, petioles, + peduncles, flower parts -- there's no single "tree" label in Helios, so + it's the complement of fruit and leaf. This is a pure read (filter), so + it's safe to call at any point -- unlike render_semantic_masks below, + which mutates primitive colors.""" + uuids = list(uuids) + fruit = context.filterPrimitivesByData(uuids, "object_label", "fruit") + leaf = context.filterPrimitivesByData(uuids, "object_label", "leaf") + tree = list(set(uuids) - set(fruit) - set(leaf)) + return {"fruit": fruit, "leaf": leaf, "tree": tree} + + +def render_semantic_masks(context, class_uuids, camera_poses, out_dir, width=IMAGE_WIDTH, height=IMAGE_HEIGHT, + fov_deg=FOV_DEG): + """Render one flat-colored, unlit pass per camera pose (the same poses as + the RGB dataset, for pixel-perfect alignment) and decode it into a + single-channel class-index PNG: 0=background, 1=fruit, 2=leaf, 3=tree. + + Uses overridePrimitiveTextureColor + lighting_model="none" so every + primitive renders as its exact flat class color with no texture/shading, + then nearest-color matching recovers the class per pixel. This mutates + each primitive's stored flat color for the rest of the Context's + lifetime (only the texture-color override is restored afterward) -- call + this AFTER sample_point_cloud(), which reads primitive colors and would + otherwise pick up these synthetic mask colors instead of the tree's + actual colors. + """ + os.makedirs(out_dir, exist_ok=True) + all_uuids = [u for uuids in class_uuids.values() for u in uuids] + + for class_name in CLASS_NAMES: + color = tuple(c / 255 for c in CLASS_MASK_COLORS[class_name]) + context.setPrimitiveColor(class_uuids[class_name], RGBcolor(*color)) + context.overridePrimitiveTextureColor(all_uuids) + + reference_colors = np.array( + [MASK_BACKGROUND_RGB] + [CLASS_MASK_COLORS[c] for c in CLASS_NAMES], dtype=np.float32 + ) + + mask_paths = [] + with Visualizer(width=width, height=height, headless=True) as visualizer: + visualizer.buildContextGeometry(context, uuids=all_uuids) + visualizer.setBackgroundColor(RGBcolor(*(c / 255 for c in MASK_BACKGROUND_RGB))) + visualizer.setLightingModel("none") + visualizer.setCameraFieldOfView(fov_deg) + + for i, (eye, lookAt) in enumerate(camera_poses): + visualizer.setCameraPosition(position=eye, lookAt=lookAt) + visualizer.plotUpdate() + viz_path = os.path.join(out_dir, f"maskviz_{i:04d}.png") + visualizer.printWindow(viz_path) + + raw = np.asarray(Image.open(viz_path).convert("RGB"), dtype=np.float32) + dists = np.linalg.norm(raw[:, :, None, :] - reference_colors[None, None, :, :], axis=-1) + class_idx = np.argmin(dists, axis=-1) + class_idx[np.min(dists, axis=-1) > MASK_COLOR_MATCH_THRESHOLD] = 0 + + mask_path = os.path.join(out_dir, f"mask_{i:04d}.png") + Image.fromarray(class_idx.astype(np.uint8), mode="L").save(mask_path) + mask_paths.append(mask_path) + print(f" mask {i + 1}/{len(camera_poses)}: {mask_path}") + + context.usePrimitiveTextureColor(all_uuids) + return mask_paths + + +# --------------------------------------------------------------------------- +# Camera plane rig (extends the fixed 3-shot rig in apple_tree_cameras.py into +# dozens of views with known poses, one grid per tree, since Gaussian +# Splatting needs many views rather than a few fixed angles) +# --------------------------------------------------------------------------- + +def plane_camera_poses(center, tree_height, tree_x_extent, tree_y_extent, num_cols=NUM_GRID_COLS, + num_rows=NUM_GRID_ROWS, fov_deg=FOV_DEG, aspect_ratio=IMAGE_WIDTH / IMAGE_HEIGHT, + margin=DISTANCE_MARGIN, span_fraction=GRID_SPAN_FRACTION): + """Camera poses for ONE tree: rather than orbiting around it, the camera + stays on a single frontal plane at a fixed distance from that tree's own + center (facing along -Y) and translates across a grid of horizontal/ + vertical offsets within that plane. lookAt is pinned to the tree center + throughout, so each shot is a slightly different parallax view of the + same face rather than a new angle around the tree. + + Distance is sized off THIS tree's own height/extent (not the whole + orchard), so each tree fills its frame the same way regardless of where + it sits in the row. It must also clear the tree's extent along the + approach direction (tree_y_extent), or a wide canopy would poke through + the near clipping plane.""" + half_fov_v = math.radians(fov_deg) / 2 + half_fov_h = math.atan(aspect_ratio * math.tan(half_fov_v)) + distance_for_height = (tree_height / 2) / math.tan(half_fov_v) * margin + distance_for_width = (tree_x_extent / 2) / math.tan(half_fov_h) * margin + distance_for_clearance = (tree_y_extent / 2) * margin + distance = max(distance_for_height, distance_for_width, distance_for_clearance) + + half_span_x = (tree_x_extent / 2) * span_fraction + half_span_z = (tree_height / 2) * span_fraction + + poses = [] + for z_off in np.linspace(-half_span_z, half_span_z, num_rows): + for x_off in np.linspace(-half_span_x, half_span_x, num_cols): + eye = vec3(center.x + x_off, center.y - distance, center.z + z_off) + poses.append((eye, center)) + return poses + + +def look_at_view_matrix(eye, lookAt, world_up=(0.0, 0.0, 1.0)): + """World-to-camera matrix in OpenCV convention (X-right, Y-down, Z-forward), + matching the eye/lookAt/world-up=+Z convention the Visualizer renders with.""" + eye_np = np.array([eye.x, eye.y, eye.z], dtype=np.float32) + center_np = np.array([lookAt.x, lookAt.y, lookAt.z], dtype=np.float32) + up_np = np.array(world_up, dtype=np.float32) + + forward = center_np - eye_np + forward /= np.linalg.norm(forward) + right = np.cross(forward, up_np) + if np.linalg.norm(right) < 1e-6: + up_np = np.array([1.0, 0.0, 0.0], dtype=np.float32) + right = np.cross(forward, up_np) + right /= np.linalg.norm(right) + down = np.cross(forward, right) + + rotation_c2w = np.stack([right, down, forward], axis=1) + viewmat = np.eye(4, dtype=np.float32) + viewmat[:3, :3] = rotation_c2w.T + viewmat[:3, 3] = -rotation_c2w.T @ eye_np + return viewmat + + +def intrinsics_matrix(width, height, fov_deg): + focal = height / (2.0 * math.tan(math.radians(fov_deg) / 2.0)) + return np.array([[focal, 0, width / 2.0], [0, focal, height / 2.0], [0, 0, 1]], dtype=np.float32) + + +# --------------------------------------------------------------------------- +# Multi-view dataset rendering + export +# --------------------------------------------------------------------------- + +def render_dataset(context, all_uuids, camera_poses, out_dir, width=IMAGE_WIDTH, height=IMAGE_HEIGHT, + fov_deg=FOV_DEG, background_rgb=BACKGROUND_RGB, test_every=TEST_EVERY): + os.makedirs(out_dir, exist_ok=True) + K = intrinsics_matrix(width, height, fov_deg) + + frames = [] + with Visualizer(width=width, height=height, headless=True) as visualizer: + visualizer.buildContextGeometry(context, uuids=all_uuids) + visualizer.setBackgroundColor(RGBcolor(*background_rgb)) + visualizer.setLightingModel("phong_shadowed") + visualizer.setCameraFieldOfView(fov_deg) + + for i, (eye, lookAt) in enumerate(camera_poses): + visualizer.setCameraPosition(position=eye, lookAt=lookAt) + visualizer.plotUpdate() + filename = os.path.join(out_dir, f"view_{i:04d}.png") + visualizer.printWindow(filename) + frames.append({ + "file_path": filename, + "viewmat": look_at_view_matrix(eye, lookAt), + "split": "test" if (i % test_every == 0) else "train", + }) + print(f" rendered {i + 1}/{len(camera_poses)}: {filename}") + + transforms = { + "camera_model": "PINHOLE", + "width": width, + "height": height, + "fl_x": float(K[0, 0]), + "fl_y": float(K[1, 1]), + "cx": float(K[0, 2]), + "cy": float(K[1, 2]), + "frames": [ + { + "file_path": os.path.relpath(f["file_path"], out_dir), + "transform_matrix": np.linalg.inv(f["viewmat"]).tolist(), + "split": f["split"], + } + for f in frames + ], + } + with open(os.path.join(out_dir, "transforms.json"), "w") as fh: + json.dump(transforms, fh, indent=2) + + return frames, K + + +def load_dataset_tensors(frames, K, device=DEVICE, mask_paths=None, target_class_index=None, + background_rgb=BACKGROUND_RGB): + """Load rendered frames as training tensors. If mask_paths and + target_class_index are given (see CLASS_NAMES / TARGET_CLASS), pixels + outside that class are replaced with the background color before + training -- the semantic mask acts as a prior so the model only ever + gets photometric supervision for the target class (e.g. "fruit"), and + everywhere else is treated as empty background.""" + K_t = torch.from_numpy(K).to(device) + dataset = {"train": [], "test": []} + for i, frame in enumerate(frames): + image = np.asarray(Image.open(frame["file_path"]).convert("RGB"), dtype=np.float32) / 255.0 + if mask_paths is not None and target_class_index is not None: + class_mask = np.asarray(Image.open(mask_paths[i]).convert("L")) + keep = class_mask == target_class_index + image = np.where(keep[..., None], image, np.array(background_rgb, dtype=np.float32)) + dataset[frame["split"]].append({ + "image": torch.from_numpy(image).to(device), + "viewmat": torch.from_numpy(frame["viewmat"]).to(device), + }) + return dataset, K_t + + +# --------------------------------------------------------------------------- +# Gaussian Splatting: init, training, eval, export +# --------------------------------------------------------------------------- + +def inverse_sigmoid(x): + return math.log(x / (1 - x)) + + +def init_gaussians(positions, colors, scene_radius, device=DEVICE): + tree = cKDTree(positions) + dists, _ = tree.query(positions, k=4) + nn_dist = np.clip(dists[:, 1:].mean(axis=1), scene_radius * 1e-4, scene_radius * 0.1) + + n = positions.shape[0] + scales_init = np.repeat(nn_dist[:, None], 3, axis=1) + quats_init = np.zeros((n, 4), dtype=np.float32) + quats_init[:, 0] = 1.0 + sh0_init = (colors - 0.5) / SH0_C0 + + return { + "means": nn.Parameter(torch.tensor(positions, dtype=torch.float32, device=device)), + "scales": nn.Parameter(torch.log(torch.tensor(scales_init, dtype=torch.float32, device=device))), + "quats": nn.Parameter(torch.tensor(quats_init, dtype=torch.float32, device=device)), + "opacities": nn.Parameter(torch.full((n,), inverse_sigmoid(0.1), dtype=torch.float32, device=device)), + "sh0": nn.Parameter(torch.tensor(sh0_init, dtype=torch.float32, device=device).unsqueeze(1)), + "shN": nn.Parameter(torch.zeros((n, 0, 3), dtype=torch.float32, device=device)), + } + + +def _gaussian_window(window_size, sigma, device): + coords = torch.arange(window_size, dtype=torch.float32, device=device) - window_size // 2 + g = torch.exp(-(coords ** 2) / (2 * sigma ** 2)) + g /= g.sum() + return g.outer(g) + + +def ssim(img1, img2, window_size=11): + """Single-scale SSIM between two (H,W,3) images in [0,1].""" + channels = img1.shape[-1] + window = _gaussian_window(window_size, 1.5, img1.device).expand(channels, 1, window_size, window_size) + x = img1.permute(2, 0, 1).unsqueeze(0) + y = img2.permute(2, 0, 1).unsqueeze(0) + pad = window_size // 2 + + mu_x = F.conv2d(x, window, padding=pad, groups=channels) + mu_y = F.conv2d(y, window, padding=pad, groups=channels) + mu_x2, mu_y2, mu_xy = mu_x * mu_x, mu_y * mu_y, mu_x * mu_y + + sigma_x2 = F.conv2d(x * x, window, padding=pad, groups=channels) - mu_x2 + sigma_y2 = F.conv2d(y * y, window, padding=pad, groups=channels) - mu_y2 + sigma_xy = F.conv2d(x * y, window, padding=pad, groups=channels) - mu_xy + + c1, c2 = 0.01 ** 2, 0.03 ** 2 + ssim_map = ((2 * mu_xy + c1) * (2 * sigma_xy + c2)) / ((mu_x2 + mu_y2 + c1) * (sigma_x2 + sigma_y2 + c2)) + return ssim_map.mean() + + +def render_view(params, viewmat, K, width, height, background): + colors = torch.cat([params["sh0"], params["shN"]], dim=1) + render, alpha, info = rasterization( + means=params["means"], + quats=params["quats"], + scales=torch.exp(params["scales"]), + opacities=torch.sigmoid(params["opacities"]), + colors=colors, + viewmats=viewmat.unsqueeze(0), + Ks=K.unsqueeze(0), + width=width, + height=height, + sh_degree=0, + backgrounds=background, + packed=False, + ) + return render[0], alpha[0], info + + +def _build_strategy(strategy_name, iters): + """The two densification schemes gsplat ships take different constructor + args and different initialize_state()/step_*_backward() signatures (see + module docstring on SPLAT_STRATEGIES), so construction is centralized here + rather than branching inline in the training loop.""" + if strategy_name == "default": + return DefaultStrategy(refine_stop_iter=int(iters * 0.75), refine_every=100) + if strategy_name == "mcmc": + return MCMCStrategy(cap_max=MCMC_CAP_MAX, refine_stop_iter=int(iters * 0.75), refine_every=100) + raise ValueError(f"Unknown strategy_name: {strategy_name!r} (expected one of {SPLAT_STRATEGIES})") + + +def train_gaussians(params, dataset, K, width, height, background_rgb, strategy_name="default", + iters=TRAIN_ITERS, scene_radius=1.0, device=DEVICE): + lrs = { + "means": 1.6e-4 * scene_radius, + "scales": 5e-3, + "quats": 1e-3, + "opacities": 5e-2, + "sh0": 2.5e-3, + "shN": 2.5e-3 / 20, + } + optimizers = {name: torch.optim.Adam([params[name]], lr=lrs[name], eps=1e-15) for name in params} + means_scheduler = torch.optim.lr_scheduler.ExponentialLR( + optimizers["means"], gamma=(1.6e-6 / 1.6e-4) ** (1.0 / iters) + ) + + strategy = _build_strategy(strategy_name, iters) + strategy.check_sanity(params, optimizers) + # DefaultStrategy scales its split heuristics off scene_radius; MCMCStrategy + # has no such notion (relocation is opacity-driven, not scale-driven). + state = strategy.initialize_state(scene_scale=scene_radius) if strategy_name == "default" \ + else strategy.initialize_state() + + background = torch.tensor(background_rgb, dtype=torch.float32, device=device).unsqueeze(0) + train_views = dataset["train"] + + for step in range(1, iters + 1): + view = train_views[np.random.randint(len(train_views))] + rendered, _, info = render_view(params, view["viewmat"], K, width, height, background) + target = view["image"] + + l1 = (rendered - target).abs().mean() + loss = 0.8 * l1 + 0.2 * (1.0 - ssim(rendered, target)) + if strategy_name == "mcmc": + # 3DGS-MCMC paper regularizes opacity/scale directly since there's + # no explicit split/prune step to keep them in check otherwise. + loss = (loss + MCMC_OPACITY_REG * torch.sigmoid(params["opacities"]).abs().mean() + + MCMC_SCALE_REG * torch.exp(params["scales"]).abs().mean()) + + if strategy_name == "default": + strategy.step_pre_backward(params, optimizers, state, step, info) + loss.backward() + if strategy_name == "mcmc": + strategy.step_post_backward(params, optimizers, state, step, info, lr=means_scheduler.get_last_lr()[0]) + else: + strategy.step_post_backward(params, optimizers, state, step, info) + + for opt in optimizers.values(): + opt.step() + opt.zero_grad(set_to_none=True) + means_scheduler.step() + + if step % 200 == 0 or step == 1: + print(f" [{strategy_name}] iter {step}/{iters} loss={loss.item():.4f} l1={l1.item():.4f} " + f"gaussians={params['means'].shape[0]}") + + return params + + +def evaluate(params, dataset, K, width, height, background_rgb, out_dir, device=DEVICE): + os.makedirs(out_dir, exist_ok=True) + background = torch.tensor(background_rgb, dtype=torch.float32, device=device).unsqueeze(0) + psnrs = [] + with torch.no_grad(): + for i, view in enumerate(dataset["test"]): + rendered, _, _ = render_view(params, view["viewmat"], K, width, height, background) + rendered = rendered.clamp(0, 1) + target = view["image"] + + mse = ((rendered - target) ** 2).mean().item() + psnr = 10 * math.log10(1.0 / max(mse, 1e-10)) + psnrs.append(psnr) + + comparison = (torch.cat([target, rendered], dim=1).clamp(0, 1).cpu().numpy() * 255).astype(np.uint8) + Image.fromarray(comparison).save(os.path.join(out_dir, f"eval_{i:03d}_psnr{psnr:.1f}.png")) + + mean_psnr = float(np.mean(psnrs)) if psnrs else float("nan") + print(f"Held-out PSNR: mean={mean_psnr:.2f} dB over {len(psnrs)} views") + return mean_psnr + + +def export_ply(params, out_path): + os.makedirs(os.path.dirname(out_path), exist_ok=True) + export_splats( + means=params["means"].detach(), + scales=params["scales"].detach(), + quats=F.normalize(params["quats"], dim=-1).detach(), + opacities=params["opacities"].detach(), + sh0=params["sh0"].detach(), + shN=params["shN"].detach(), + format="ply", + save_to=out_path, + ) + + +def main(): + AGE_DAYS = 720.0 + POSITIONS = [vec3(0, 0, 0), vec3(1.5, 0, 0), vec3(3, 0, 0)] + + os.makedirs(OUTPUT_DIR, exist_ok=True) + + with Context() as context: + with PlantArchitecture(context) as plantarch: + plant_ids = build_orchard(plantarch, POSITIONS, AGE_DAYS) + tree_uuids_list = [plantarch.getAllPlantUUIDs(plant_id) for plant_id in plant_ids] + all_uuids = [u for uuids in tree_uuids_list for u in uuids] + print(f"Built {len(plant_ids)} apple trees, {len(all_uuids)} primitives total") + + camera_poses = [] + for tree_uuids in tree_uuids_list: + tree_center, _ = context.getDomainBoundingSphere(tree_uuids) + x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(tree_uuids) + tree_height = z_bounds.y - z_bounds.x + tree_x_extent = x_bounds.y - x_bounds.x + tree_y_extent = y_bounds.y - y_bounds.x + print(f" tree center=({tree_center.x:.2f}, {tree_center.y:.2f}, {tree_center.z:.2f}) " + f"height={tree_height:.2f}m x_extent={tree_x_extent:.2f}m y_extent={tree_y_extent:.2f}m") + camera_poses.extend(plane_camera_poses(tree_center, tree_height, tree_x_extent, tree_y_extent)) + + print(f"Rendering {len(camera_poses)} multi-view images...") + frames, K = render_dataset(context, all_uuids, camera_poses, DATASET_DIR) + + print("Classifying primitives by organ type...") + class_uuids = classify_primitives(context, all_uuids) + for name in CLASS_NAMES: + print(f" {name}: {len(class_uuids[name])} primitives") + + print(f"Sampling seed point cloud from tree geometry (class={TARGET_CLASS})...") + positions, colors = sample_point_cloud(context, class_uuids[TARGET_CLASS]) + print(f" {len(positions)} seed points") + # scale-init radius from the actual points being fit, not the whole + # scene -- a class filter (e.g. fruit only) covers a much smaller + # volume than the full tree, so the full bounding sphere would be + # a poor scale for the Gaussian init/LR heuristics below. + point_cloud_radius = float(np.linalg.norm(positions - positions.mean(axis=0), axis=1).max()) + + # Mutates primitive colors, so it must run after sample_point_cloud + # (see render_semantic_masks docstring) -- reuses the same camera + # poses as the RGB dataset for pixel-perfect alignment. + print(f"Rendering {len(camera_poses)} semantic mask images...") + mask_paths = render_semantic_masks(context, class_uuids, camera_poses, DATASET_DIR) + + target_class_index = CLASS_NAMES.index(TARGET_CLASS) + 1 + dataset, K_t = load_dataset_tensors(frames, K, mask_paths=mask_paths, target_class_index=target_class_index) + print(f"Train views: {len(dataset['train'])} Test views: {len(dataset['test'])}") + + # Same seed points/dataset, trained once per strategy the library ships + # (see SPLAT_STRATEGIES) so the two can be compared side by side. + results = [] + for strategy_name in SPLAT_STRATEGIES: + print(f"--- Strategy: {strategy_name} ---") + params = init_gaussians(positions, colors, point_cloud_radius) + print(f"Training {params['means'].shape[0]} Gaussians for {TRAIN_ITERS} iterations on {DEVICE}...") + params = train_gaussians(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, + strategy_name=strategy_name, iters=TRAIN_ITERS, scene_radius=point_cloud_radius) + + mean_psnr = evaluate(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, + os.path.join(EVAL_DIR, strategy_name)) + + output_ply_path = os.path.join(OUTPUT_DIR, f"apple_splats_{TARGET_CLASS}_{strategy_name}.ply") + export_ply(params, output_ply_path) + print(f"Saved trained splats to {output_ply_path}") + + results.append({ + "strategy": strategy_name, + "num_gaussians": params["means"].shape[0], + "mean_psnr": mean_psnr, + "output_ply_path": output_ply_path, + }) + + return results + + +if __name__ == "__main__": + try: + results = main() + summary = "; ".join( + f"{r['strategy']}: {r['num_gaussians']} Gaussians, PSNR={r['mean_psnr']:.2f}dB -> {r['output_ply_path']}" + for r in results + ) + notify_slack(f":white_check_mark: apple_tree_gaussian_splatting.py finished (class=fruit) — {summary}") + except Exception as e: + notify_slack(f":x: apple_tree_gaussian_splatting.py failed: {type(e).__name__}: {e}") + raise From 49282e3b0cdada8905e42ac8bfaa446e28c7cfd2 Mon Sep 17 00:00:00 2001 From: Yogesh Chawla Date: Tue, 28 Jul 2026 13:39:21 -0700 Subject: [PATCH 03/42] Sweep capture density and viewing-plane count as separate design axes CAPTURE_CONFIGS now varies image count (grid density) and number of locations (viewing planes per tree) independently from a default baseline, training/exporting a separate fruit-only .ply per (capture config, splat strategy) pair. Co-Authored-By: Claude Sonnet 5 --- GAUSSIAN_SPLATTING.md | 77 +++++++----- apple_tree_gaussian_splatting.py | 193 ++++++++++++++++++++----------- 2 files changed, 173 insertions(+), 97 deletions(-) diff --git a/GAUSSIAN_SPLATTING.md b/GAUSSIAN_SPLATTING.md index e107111..40e5280 100644 --- a/GAUSSIAN_SPLATTING.md +++ b/GAUSSIAN_SPLATTING.md @@ -19,8 +19,9 @@ poses, and train a real 3D Gaussian Splat (via the `gsplat` library) on that dataset — seeding the Gaussians from the tree's own mesh geometry instead of COLMAP/SfM, since PyHelios already knows exact 3D positions. Training is restricted to the fruit class only (the semantic mask acts as a prior, see -§5), and is run once per densification strategy `gsplat` ships, so the two -can be compared on identical seed points/cameras. +§5), and is run once per (capture config, densification strategy) pair — +sweeping image count, number of viewing planes, and `gsplat` strategy — so +these can be compared on identical seed points/geometry. ## 2. Prerequisites check (do this first, automatable) @@ -141,9 +142,16 @@ automated run can default to it without asking: | Camera pose convention | Manually derive OpenCV world-to-camera matrix (X-right, Y-down, Z-forward) from `eye`/`lookAt`/world-up=`(0,0,1)`, verified numerically (`right × down == forward`) | `gsplat` needs OpenCV-convention viewmats; Helios exposes no matrix getter, so it must be reconstructed to exactly match what Helios's internal `glm::lookAt` actually rendered | | Real camera intrinsics (`fx,fy,cx,cy` + distortion) given by user | Derive render resolution from `(cx,cy)`, vertical FOV from `fy`; drop distortion and off-center principal point; gsplat's K matrix is a **centered** pinhole matching what was actually rendered | Helios's renderer physically cannot produce `fx != fy`, off-center principal point, or lens distortion — passing the raw values into gsplat while training on Helios pixels would be geometrically inconsistent with the images themselves | | Loss function | `0.8 * L1 + 0.2 * (1 - SSIM)`, compact from-scratch SSIM (11×11 Gaussian window, `conv2d`) | Standard 3DGS recipe; no extra dependency needed for SSIM | -| Camera rig | Treat each tree as independent: per-tree frontal-plane grid (`NUM_GRID_COLS x NUM_GRID_ROWS` positions on one plane at a fixed distance from *that tree's own* center, always looking at it) instead of one 360° orbit around the whole orchard | A single orbit around the combined bounding box over/under-frames trees that sit off-center in a multi-tree row; a fixed per-tree frontal plane keeps every tree consistently framed regardless of where it sits, without needing to circle it | +| Camera rig | Treat each tree as independent: per-tree flat-plane grid at a fixed distance from *that tree's own* center, always looking at it, instead of one 360° orbit around the whole orchard | A single orbit around the combined bounding box over/under-frames trees that sit off-center in a multi-tree row; a fixed per-tree plane keeps every tree consistently framed regardless of where it sits, without needing to circle it | | Target class for training | Restrict to `TARGET_CLASS = "fruit"` always (not a toggle) — seed points, dataset masking, and export are all apples-only | User wants the deliverable to be an apple-only splat, not the whole tree; the semantic mask (§ render_semantic_masks) already exists as a clean prior for this | | Splat densification strategy | Train **once per strategy** gsplat ships (`DefaultStrategy`, `MCMCStrategy`) on the same seed points/dataset, exporting a separate `.ply` per strategy rather than picking one | Neither strategy is a strict upgrade — comparing both on the same data is cheap (just a second training pass) and avoids committing to one without evidence | +| Number of images / number of viewing planes | Sweep both as separate `CAPTURE_CONFIGS` cases (`sparse`: 1 plane, 4x2 grid; `default`: 1 plane, 8x3 grid; `multi_face`: 4 planes, 8x3 grid each), one variable changed at a time from the `default` baseline, each producing its own dataset + its own `.ply` per strategy | Answers "does more images help" (`sparse` vs `default`) and "does seeing more of the tree help" (`default` vs `multi_face`) independently, rather than guessing one fixed capture density; a full cross product of every count x every plane count wasn't worth the extra training runs for this comparison | + +**Cost note:** `CAPTURE_CONFIGS x SPLAT_STRATEGIES` = 3 x 2 = 6 full training +runs by default (each `TRAIN_ITERS` iterations), plus 3 separate render +passes (48 + 144 + 576 = 768 renders total across `sparse`/`default`/ +`multi_face`). Trim `CAPTURE_CONFIGS` or lower `TRAIN_ITERS` for a quick +smoke test rather than running the full matrix every time. ## 6. Known bug: PyHelios headless `Visualizer` heap corruption at large resolution @@ -194,23 +202,24 @@ the real fix. ``` apple_tree_gaussian_splatting.py ├── build_orchard() # reuses build_apple_tree() from apple_tree.py -├── sample_point_cloud() # mesh-surface sampling + texture-average color fallback -├── plane_camera_poses() # per-tree frontal-plane grid (NUM_GRID_COLS x NUM_GRID_ROWS), -│ # called once per tree (extends the apple_tree_cameras.py -│ # 3-shot rig; no orbiting, fixed distance per tree) +├── sample_point_cloud() # mesh-surface sampling + texture-average color fallback, once for all cases +├── plane_camera_poses() # per-tree grid across num_planes flat faces (num_cols x num_rows each), +│ # called once per tree per CAPTURE_CONFIGS case (extends the +│ # apple_tree_cameras.py 3-shot rig; no orbiting, fixed distance per tree) ├── look_at_view_matrix() # OpenCV-convention world-to-camera matrix ├── intrinsics_matrix() # K from (width, height, vertical FOV) -├── render_dataset() # headless multi-view render + transforms.json export +├── render_dataset() # headless multi-view render + transforms.json export, once per capture config ├── classify_primitives() # fruit/leaf/tree split via Helios "object_label" data -├── render_semantic_masks() # flat-colored mask render -> class-index PNG per view +├── render_semantic_masks() # flat-colored mask render -> class-index PNG per view, once per capture config ├── init_gaussians() # KDTree-based scale init, RGB->SH0 color init ├── train_gaussians() # gsplat rasterization + L1/SSIM loss, strategy-agnostic ├── _build_strategy() # DefaultStrategy or MCMCStrategy, per SPLAT_STRATEGIES ├── evaluate() # held-out PSNR + side-by-side comparison PNGs -└── export_ply() # gsplat.export_splats -> one .ply per strategy +└── export_ply() # gsplat.export_splats -> one .ply per (capture config, strategy) pair -main() trains TARGET_CLASS="fruit" once per entry in SPLAT_STRATEGIES, -reusing the same rendered dataset/seed points for each. +main() renders CAPTURE_CONFIGS once each (shared seed points across all of +them), then trains TARGET_CLASS="fruit" once per (capture config, strategy) +pair -- see §5 for what each axis controls. ``` ## 8. Validation checklist before trusting a full run @@ -218,18 +227,22 @@ reusing the same rendered dataset/seed points for each. Always smoke-test at reduced scale first — this is what caught the crash in §6 before it wasted a full 5000-iteration run: -1. 1 tree, ~16 views, small resolution (e.g. 200×200), ~300-700 iterations. - Confirms: tree build → render → point sampling → training forward/backward - → `DefaultStrategy` density-control path (needs `iters > refine_start_iter`, - default 500) → eval → PLY export, all with no exceptions. +1. 1 tree, ~16 views, small resolution (e.g. 200×200), ~300-700 iterations, + a single `CAPTURE_CONFIGS` entry. Confirms: tree build → render → point + sampling → training forward/backward → `DefaultStrategy` density-control + path (needs `iters > refine_start_iter`, default 500) → eval → PLY + export, all with no exceptions. 2. Visually inspect one ground-truth render and one eval comparison PNG - (`renders/.../eval//eval_*.png`, side-by-side GT|rendered) — - confirms camera orientation is correct (not flipped/mirrored) and colors - are sane, which numeric checks alone won't catch. Do this for **both** - strategies, not just the first — a strategy-specific bug (e.g. the - `MCMCStrategy` missing `lr` kwarg) can silently corrupt only one output. -3. Only then run the full default (3 trees, `NUM_GRID_COLS x NUM_GRID_ROWS` - views per tree, fruit class only, 5000 iterations x 2 strategies). + (`renders/.../eval///eval_*.png`, side-by-side + GT|rendered) — confirms camera orientation is correct (not + flipped/mirrored) and colors are sane, which numeric checks alone won't + catch. Do this for **both** strategies, not just the first — a + strategy-specific bug (e.g. the `MCMCStrategy` missing `lr` kwarg) can + silently corrupt only one output. For `multi_face`, check at least one + view from each of the 4 planes, not just the frontal one. +3. Only then run the full default (3 trees, all of `CAPTURE_CONFIGS`, fruit + class only, 5000 iterations x 2 strategies each — see the cost note in + §5). ## 9. Reference full run (validates the fix in §6) @@ -255,13 +268,19 @@ cd /home/yogesh/PyHelios python3 apple_tree_gaussian_splatting.py ``` -Trains `TARGET_CLASS="fruit"` once per entry in `SPLAT_STRATEGIES` -(`default`, `mcmc`). Output: +Trains `TARGET_CLASS="fruit"` once per (capture config, strategy) pair -- +`CAPTURE_CONFIGS` (`sparse`, `default`, `multi_face` by default) x +`SPLAT_STRATEGIES` (`default`, `mcmc`). See the cost note in §5 before +running the full matrix. Output: ``` renders/gaussian_splatting/ -├── dataset/ # shared RGB + mask renders, one set for both strategies -├── eval/{default,mcmc}/ # per-strategy held-out comparison PNGs -├── apple_splats_fruit_default.ply -└── apple_splats_fruit_mcmc.ply +├── dataset/{sparse,default,multi_face}/ # per-capture-config RGB + mask renders +├── eval/{sparse,default,multi_face}/{default,mcmc}/ # per-capture-config, per-strategy comparison PNGs +├── apple_splats_fruit_sparse_default.ply +├── apple_splats_fruit_sparse_mcmc.ply +├── apple_splats_fruit_default_default.ply +├── apple_splats_fruit_default_mcmc.ply +├── apple_splats_fruit_multi_face_default.ply +└── apple_splats_fruit_multi_face_mcmc.ply ``` diff --git a/apple_tree_gaussian_splatting.py b/apple_tree_gaussian_splatting.py index 4d32214..06075b1 100644 --- a/apple_tree_gaussian_splatting.py +++ b/apple_tree_gaussian_splatting.py @@ -49,12 +49,30 @@ _render_scale = min(1.0, MAX_RENDER_DIMENSION / max(round(2 * CAMERA_CX), round(2 * CAMERA_CY))) IMAGE_WIDTH = round(2 * CAMERA_CX * _render_scale) IMAGE_HEIGHT = round(2 * CAMERA_CY * _render_scale) -NUM_GRID_COLS = 8 # horizontal camera positions within the frontal plane -NUM_GRID_ROWS = 3 # vertical camera positions within the frontal plane GRID_SPAN_FRACTION = 0.6 # how far the grid spans across the tree's own extent/height TEST_EVERY = 8 DISTANCE_MARGIN = 1.05 # tight headroom so each tree nearly fills the frame +# Design decision: rather than committing to one capture density, render and +# train a separate case per entry here so "how many images" and "how many +# viewing planes" can be compared against the same baseline instead of +# guessing. Each case is fully independent (its own render pass + its own +# per-strategy training runs below), so the total cost scales as +# len(CAPTURE_CONFIGS) x len(SPLAT_STRATEGIES) full training runs -- trim +# this list (or TRAIN_ITERS) for a quick smoke test. +# - num_planes: how many faces around each tree get their own flat grid +# (1 = frontal only, matching the original single-plane rig; N = N +# evenly-spaced azimuths, e.g. front/right/back/left for N=4). This is +# the "number of locations" axis -- it does NOT orbit smoothly, it's N +# discrete flat planes. +# - num_cols x num_rows: grid density *within* each plane. This is the +# "number of images" axis at fixed camera placement. +CAPTURE_CONFIGS = ( + {"name": "sparse", "num_planes": 1, "num_cols": 4, "num_rows": 2}, # fewer images, same 1 face as baseline + {"name": "default", "num_planes": 1, "num_cols": 8, "num_rows": 3}, # baseline + {"name": "multi_face", "num_planes": 4, "num_cols": 8, "num_rows": 3}, # baseline density, 4 faces per tree +) + OUTPUT_DIR = "renders/gaussian_splatting" DATASET_DIR = os.path.join(OUTPUT_DIR, "dataset") EVAL_DIR = os.path.join(OUTPUT_DIR, "eval") @@ -225,36 +243,53 @@ def render_semantic_masks(context, class_uuids, camera_poses, out_dir, width=IMA # Splatting needs many views rather than a few fixed angles) # --------------------------------------------------------------------------- -def plane_camera_poses(center, tree_height, tree_x_extent, tree_y_extent, num_cols=NUM_GRID_COLS, - num_rows=NUM_GRID_ROWS, fov_deg=FOV_DEG, aspect_ratio=IMAGE_WIDTH / IMAGE_HEIGHT, - margin=DISTANCE_MARGIN, span_fraction=GRID_SPAN_FRACTION): - """Camera poses for ONE tree: rather than orbiting around it, the camera - stays on a single frontal plane at a fixed distance from that tree's own - center (facing along -Y) and translates across a grid of horizontal/ - vertical offsets within that plane. lookAt is pinned to the tree center - throughout, so each shot is a slightly different parallax view of the - same face rather than a new angle around the tree. - - Distance is sized off THIS tree's own height/extent (not the whole +def plane_camera_poses(center, tree_height, tree_x_extent, tree_y_extent, num_planes, num_cols, num_rows, + fov_deg=FOV_DEG, aspect_ratio=IMAGE_WIDTH / IMAGE_HEIGHT, margin=DISTANCE_MARGIN, + span_fraction=GRID_SPAN_FRACTION): + """Camera poses for ONE tree: rather than orbiting smoothly around it, the + camera sits on `num_planes` discrete flat planes facing the tree (evenly + spaced azimuths, starting from the original single-plane direction -Y) + and, within each plane, translates across a `num_cols` x `num_rows` grid + of horizontal/vertical offsets at a fixed distance. lookAt is pinned to + the tree center throughout, so each shot is a slightly different + parallax view of the same face rather than a new angle swept smoothly + around the tree. num_planes=1 reproduces the original frontal-only rig. + + Distance is sized off THIS tree's own height/footprint (not the whole orchard), so each tree fills its frame the same way regardless of where - it sits in the row. It must also clear the tree's extent along the - approach direction (tree_y_extent), or a wide canopy would poke through - the near clipping plane.""" + it sits in the row, computed per-plane the same way a single fixed + azimuth would: it must also clear the tree's footprint along the + approach direction, or a wide + canopy would poke through the near clipping plane.""" half_fov_v = math.radians(fov_deg) / 2 half_fov_h = math.atan(aspect_ratio * math.tan(half_fov_v)) + half_x, half_y = tree_x_extent / 2, tree_y_extent / 2 distance_for_height = (tree_height / 2) / math.tan(half_fov_v) * margin - distance_for_width = (tree_x_extent / 2) / math.tan(half_fov_h) * margin - distance_for_clearance = (tree_y_extent / 2) * margin - distance = max(distance_for_height, distance_for_width, distance_for_clearance) - half_span_x = (tree_x_extent / 2) * span_fraction half_span_z = (tree_height / 2) * span_fraction + base_az = -math.pi / 2 # -Y, matching the original single-plane rig this generalizes poses = [] - for z_off in np.linspace(-half_span_z, half_span_z, num_rows): - for x_off in np.linspace(-half_span_x, half_span_x, num_cols): - eye = vec3(center.x + x_off, center.y - distance, center.z + z_off) - poses.append((eye, center)) + for p in range(num_planes): + az = base_az + 2 * math.pi * p / num_planes + forward = (math.cos(az), math.sin(az)) # direction from tree center to camera + tangent = (-math.sin(az), math.cos(az)) # in-plane horizontal axis + + half_width_at_az = abs(math.sin(az)) * half_x + abs(math.cos(az)) * half_y + depth_at_az = abs(math.cos(az)) * half_x + abs(math.sin(az)) * half_y + distance_for_width = half_width_at_az / math.tan(half_fov_h) * margin + distance_for_clearance = depth_at_az * margin + distance = max(distance_for_height, distance_for_width, distance_for_clearance) + half_span_t = half_width_at_az * span_fraction + + for z_off in np.linspace(-half_span_z, half_span_z, num_rows): + for t_off in np.linspace(-half_span_t, half_span_t, num_cols): + eye = vec3( + center.x + distance * forward[0] + t_off * tangent[0], + center.y + distance * forward[1] + t_off * tangent[1], + center.z + z_off, + ) + poses.append((eye, center)) return poses @@ -547,6 +582,7 @@ def main(): POSITIONS = [vec3(0, 0, 0), vec3(1.5, 0, 0), vec3(3, 0, 0)] os.makedirs(OUTPUT_DIR, exist_ok=True) + capture_datasets = [] # one entry per CAPTURE_CONFIGS case: {"name", "frames", "K", "mask_paths"} with Context() as context: with PlantArchitecture(context) as plantarch: @@ -555,25 +591,13 @@ def main(): all_uuids = [u for uuids in tree_uuids_list for u in uuids] print(f"Built {len(plant_ids)} apple trees, {len(all_uuids)} primitives total") - camera_poses = [] - for tree_uuids in tree_uuids_list: - tree_center, _ = context.getDomainBoundingSphere(tree_uuids) - x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(tree_uuids) - tree_height = z_bounds.y - z_bounds.x - tree_x_extent = x_bounds.y - x_bounds.x - tree_y_extent = y_bounds.y - y_bounds.x - print(f" tree center=({tree_center.x:.2f}, {tree_center.y:.2f}, {tree_center.z:.2f}) " - f"height={tree_height:.2f}m x_extent={tree_x_extent:.2f}m y_extent={tree_y_extent:.2f}m") - camera_poses.extend(plane_camera_poses(tree_center, tree_height, tree_x_extent, tree_y_extent)) - - print(f"Rendering {len(camera_poses)} multi-view images...") - frames, K = render_dataset(context, all_uuids, camera_poses, DATASET_DIR) - print("Classifying primitives by organ type...") class_uuids = classify_primitives(context, all_uuids) for name in CLASS_NAMES: print(f" {name}: {len(class_uuids[name])} primitives") + # Independent of camera/capture config -- seeded once from mesh + # geometry and reused for every capture config x strategy below. print(f"Sampling seed point cloud from tree geometry (class={TARGET_CLASS})...") positions, colors = sample_point_cloud(context, class_uuids[TARGET_CLASS]) print(f" {len(positions)} seed points") @@ -583,39 +607,71 @@ def main(): # a poor scale for the Gaussian init/LR heuristics below. point_cloud_radius = float(np.linalg.norm(positions - positions.mean(axis=0), axis=1).max()) - # Mutates primitive colors, so it must run after sample_point_cloud - # (see render_semantic_masks docstring) -- reuses the same camera - # poses as the RGB dataset for pixel-perfect alignment. - print(f"Rendering {len(camera_poses)} semantic mask images...") - mask_paths = render_semantic_masks(context, class_uuids, camera_poses, DATASET_DIR) + for capture_config in CAPTURE_CONFIGS: + capture_name = capture_config["name"] + print(f"=== Capture config: {capture_name} (num_planes={capture_config['num_planes']}, " + f"grid={capture_config['num_cols']}x{capture_config['num_rows']}) ===") + + camera_poses = [] + for tree_uuids in tree_uuids_list: + tree_center, _ = context.getDomainBoundingSphere(tree_uuids) + x_bounds, y_bounds, z_bounds = context.getDomainBoundingBox(tree_uuids) + tree_height = z_bounds.y - z_bounds.x + tree_x_extent = x_bounds.y - x_bounds.x + tree_y_extent = y_bounds.y - y_bounds.x + camera_poses.extend(plane_camera_poses( + tree_center, tree_height, tree_x_extent, tree_y_extent, + num_planes=capture_config["num_planes"], + num_cols=capture_config["num_cols"], num_rows=capture_config["num_rows"], + )) + + dataset_dir = os.path.join(DATASET_DIR, capture_name) + print(f"Rendering {len(camera_poses)} multi-view images...") + frames, K = render_dataset(context, all_uuids, camera_poses, dataset_dir) + + # Mutates primitive colors, so it must run after sample_point_cloud + # (see render_semantic_masks docstring) -- reuses this capture + # config's own camera poses for pixel-perfect alignment. + print(f"Rendering {len(camera_poses)} semantic mask images...") + mask_paths = render_semantic_masks(context, class_uuids, camera_poses, dataset_dir) + + capture_datasets.append({"name": capture_name, "frames": frames, "K": K, "mask_paths": mask_paths}) target_class_index = CLASS_NAMES.index(TARGET_CLASS) + 1 - dataset, K_t = load_dataset_tensors(frames, K, mask_paths=mask_paths, target_class_index=target_class_index) - print(f"Train views: {len(dataset['train'])} Test views: {len(dataset['test'])}") - # Same seed points/dataset, trained once per strategy the library ships - # (see SPLAT_STRATEGIES) so the two can be compared side by side. + # Cross product of capture config x splat strategy, each producing its own .ply. results = [] - for strategy_name in SPLAT_STRATEGIES: - print(f"--- Strategy: {strategy_name} ---") - params = init_gaussians(positions, colors, point_cloud_radius) - print(f"Training {params['means'].shape[0]} Gaussians for {TRAIN_ITERS} iterations on {DEVICE}...") - params = train_gaussians(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, - strategy_name=strategy_name, iters=TRAIN_ITERS, scene_radius=point_cloud_radius) - - mean_psnr = evaluate(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, - os.path.join(EVAL_DIR, strategy_name)) - - output_ply_path = os.path.join(OUTPUT_DIR, f"apple_splats_{TARGET_CLASS}_{strategy_name}.ply") - export_ply(params, output_ply_path) - print(f"Saved trained splats to {output_ply_path}") - - results.append({ - "strategy": strategy_name, - "num_gaussians": params["means"].shape[0], - "mean_psnr": mean_psnr, - "output_ply_path": output_ply_path, - }) + for capture_entry in capture_datasets: + capture_name = capture_entry["name"] + dataset, K_t = load_dataset_tensors( + capture_entry["frames"], capture_entry["K"], + mask_paths=capture_entry["mask_paths"], target_class_index=target_class_index, + ) + print(f"[{capture_name}] Train views: {len(dataset['train'])} Test views: {len(dataset['test'])}") + + for strategy_name in SPLAT_STRATEGIES: + print(f"--- Capture={capture_name} Strategy={strategy_name} ---") + params = init_gaussians(positions, colors, point_cloud_radius) + print(f"Training {params['means'].shape[0]} Gaussians for {TRAIN_ITERS} iterations on {DEVICE}...") + params = train_gaussians(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, + strategy_name=strategy_name, iters=TRAIN_ITERS, scene_radius=point_cloud_radius) + + mean_psnr = evaluate(params, dataset, K_t, IMAGE_WIDTH, IMAGE_HEIGHT, BACKGROUND_RGB, + os.path.join(EVAL_DIR, capture_name, strategy_name)) + + output_ply_path = os.path.join( + OUTPUT_DIR, f"apple_splats_{TARGET_CLASS}_{capture_name}_{strategy_name}.ply" + ) + export_ply(params, output_ply_path) + print(f"Saved trained splats to {output_ply_path}") + + results.append({ + "capture": capture_name, + "strategy": strategy_name, + "num_gaussians": params["means"].shape[0], + "mean_psnr": mean_psnr, + "output_ply_path": output_ply_path, + }) return results @@ -624,7 +680,8 @@ def main(): try: results = main() summary = "; ".join( - f"{r['strategy']}: {r['num_gaussians']} Gaussians, PSNR={r['mean_psnr']:.2f}dB -> {r['output_ply_path']}" + f"{r['capture']}/{r['strategy']}: {r['num_gaussians']} Gaussians, " + f"PSNR={r['mean_psnr']:.2f}dB -> {r['output_ply_path']}" for r in results ) notify_slack(f":white_check_mark: apple_tree_gaussian_splatting.py finished (class=fruit) — {summary}") From 4f4ccc060460a9f05cb3c5dcc075493d724a89e0 Mon Sep 17 00:00:00 2001 From: Yogesh Chawla Date: Tue, 28 Jul 2026 15:08:46 -0700 Subject: [PATCH 04/42] overall md file --- active_vision_design.md | 830 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 830 insertions(+) create mode 100644 active_vision_design.md diff --git a/active_vision_design.md b/active_vision_design.md new file mode 100644 index 0000000..f1735dc --- /dev/null +++ b/active_vision_design.md @@ -0,0 +1,830 @@ +# Multi-Camera Active Vision for Apple Harvesting — Research Design + +**Stage 1: Helios simulation.** Companion to WTFRC proposal *Robotic Harvesting: Multi-Camera Fruit Detection Under Heavy Occlusion* (Vougioukas, Bailey, Kong), Objective 2 and 3. + +*Prepared July 2026. Literature current to mid-2026. Verification caveats are collected in the final section — several 2026 arXiv entries and paywalled numbers should be confirmed before they enter a thesis or paper.* + +--- + +## 0. Executive summary — six claims + +1. **Vanilla 3DGS is offline, and the reason it is the wrong primary map is not speed — it is that Gaussians represent only *occupied* space.** For a next-best-view planner, *unknown* space is the entire signal. Every serious active-3DGS system in the literature (ActiveGS, GS-Planner, HGS-Planner, ActiveGAMER, NARUTO) bolts a voxel map on the side to fix exactly this. Build the voxel map as the primary substrate and run splatting asynchronously for appearance. + +2. **In Helios you have exact camera poses, so you should not run SLAM at all.** Tracking is where most of the compute and nearly all of the fragility of GS-SLAM systems lives. You need a *mapper*, not a SLAM system. This deletes roughly half the engineering in the proposal's Objective 2 and buys you the resolution you actually need. + +3. **Exploration and exploitation are not two settings of one objective — they are structurally different optimization problems**, and conflating them is the single most common mistake in this literature. Explore is *quality-constrained* (see everything → minimize time), exploit is *budget-constrained* (fixed time → maximize value). Maximizing information gain during exploration is provably counterproductive. This split is your architectural backbone and it is a defensible novel framing. + +4. **Your five DOF are not five equal DOF.** The 2 gimbal axes are ~free; the 3 linear axes are 0.5–3 s each. Nest a cheap gradient-based gimbal refinement inside an expensive sampling-based linear-axis path planner. This decomposition matches your hardware exactly and I could not find it published. + +5. **Do not start from VGGT.** Start from a *pose-conditioned* geometry model — MapAnything (Apache-2.0 code and weights) or Depth Anything 3. Your defining asset is that you know the poses; three models released since Sept 2025 are built to consume exactly that, and there is a hard theoretical reason (15-DoF projective ambiguity under small baselines) to expect the uncalibrated path to fail on a canopy scanned through a narrow arc. + +6. **The evaluation should be built around three quantities that no field paper can compute**: the Achievable Visibility Upper Bound (what fraction of each fruit is visible from *any* reachable pose), the oracle-normalized planning score Π, and an information-gain calibration score adapted from the optical-flow sparsification literature. These convert "we detected 70% of fruits" — a number that mixes canopy difficulty, hardware limits and algorithm quality — into three separable numbers. + +--- + +## 1. Reframing the problem statement + +The proposal says, in effect: + +> multi-camera SLAM → 3DGS map → semantic fusion → NBV planner reading the 3DGS map + +I would restructure this. The pipeline as written puts a representation that cannot express ignorance at the center of a system whose whole job is to reason about ignorance. Concretely: + +| Proposal component | Recommendation | Why | +|---|---|---| +| Multi-camera SLAM for pose | **Drop in Stage 1.** Poses come free from Helios; on the robot they come from arm encoders + a one-time hand-eye calibration, with Kong's SLAM device as a correction/drift term, not the primary source | You have a *gantry*, not a handheld sensor. Encoder-derived pose is better than visual pose in a self-similar canopy, where vision-only 3DGS-SLAM catastrophically fails (see §2.3) | +| 3DGS as the map the planner reads | **Demote to a secondary, asynchronous appearance map.** Primary map = GPU occupancy octree with explicit unknown state | 3DGS cannot represent unknown space; it also systematically renders thin geometry as low-opacity haze, which is the wrong failure mode for a canopy | +| Semantic splats | **Keep the idea, move it to the voxel layer.** Per-voxel semantic posterior p(apple), p(leaf), p(branch), p(wire) | The planner needs semantics on the thing it ray-casts against | +| Single NBV planner with 3 weighted criteria | **Two planners with different objectives**, plus a nested fast local refinement | See §4 | +| "Number of fruits discovered, location/size accuracy, exploration time per fruit" | **Five-tier evaluation suite**, headlined by AVUB/NVE and Π | See §7 | + +Everything else in the proposal — the occlusion-regulation module, the Helios orchard, the physical evaluation protocol — I would keep as written. Those are strengths. + +### 1.1 The physical scales that govern every decision + +Write these on the wall, because most published performance numbers were measured 100× coarser: + +| Structure | Scale | +|---|---| +| Apple leaf lamina thickness | 0.15–0.3 mm | +| Petiole / young twig diameter | 2–10 mm | +| Trellis wire | 2–3 mm | +| Fruit diameter | 60–90 mm | + +Nearly every mapping benchmark you will read (nvblox 0.4 ms/frame, voxblox 70 ms, wavemap memory tables) was run at **5 cm voxels on indoor rooms** — 150–300× a leaf thickness. Going to 1 mm voxels is a 50× linear increase, implying 10³–10⁴× more work. **Re-derive your compute budget at your own resolution before committing to any architecture.** This is the most likely way the project's real-time claim quietly fails. + +--- + +## 2. Real-time reconstruction + +### 2.1 Yes, vanilla 3DGS is a post-process — for four separable reasons + +Kerbl et al.'s 3DGS ([arXiv:2308.04079](https://arxiv.org/abs/2308.04079)) is offline because: + +1. **Per-scene optimization.** The Gaussians *are* the free parameters. No learned prior; every scene is a fresh gradient-descent problem, 7k–30k iterations. +2. **COLMAP dependency.** Poses and the seed point cloud come from batch SfM — often slower than the splatting it feeds. +3. **Batch view sampling.** The optimizer samples randomly from *all* images every iteration. Not causal; adding an image invalidates convergence state. +4. **Adaptive density control is a global annealing heuristic.** Clone/split/prune from statistics accumulated over the whole image set, with periodic opacity resets. + +Recent work attacks 1–3 separately. Point 4 is the least solved and is quietly the source of most artifacts in incremental systems. + +### 2.2 The online variants exist, but the honest speed spread is 400× + +RGBD GS-ICP SLAM's benchmark is the most useful artifact in this literature because it puts everything on one RTX 4090: + +| System | FPS (RTX 4090, Replica) | +|---|---| +| Point-SLAM | 0.30 | +| SplaTAM ([2312.02126](https://arxiv.org/abs/2312.02126)) | 0.23 | +| GS-SLAM | 8.34 | +| Orbeez-SLAM | 24.15 | +| RGBD GS-ICP SLAM ([2403.12550](https://arxiv.org/abs/2403.12550)) | **98.11** | + +Anyone saying "real-time 3DGS SLAM" without naming a system is making a meaningless claim. And note *how* GS-ICP buys its speed: by leaning entirely on the depth sensor and skipping photometric optimization. On thin foliage, RGB-D depth is at its noisiest and most flying-pixel-prone. The speed/robustness trade is aimed away from your scene. + +Others worth knowing: MonoGS ([2312.06741](https://arxiv.org/abs/2312.06741), 3 FPS, monocular), Photo-SLAM ([2311.16728](https://arxiv.org/abs/2311.16728), runs on Jetson Orin), RTG-SLAM (SIGGRAPH 2024, [2404.19706](https://arxiv.org/abs/2404.19706)), LoopSplat (3DV 2025, loop closure causes latency spikes — bad for a fixed-rate control loop), WildGS-SLAM (CVPR 2025, [2504.03886](https://arxiv.org/abs/2504.03886)) — that last one is a cautionary tale: it identifies and *deletes* dynamic objects, which in a windy canopy means deleting the canopy. + +### 2.3 The orchard evidence is damning for vision-only splatting + +AgriGS-SLAM ([2510.26358](https://arxiv.org/abs/2510.26358), Usuelli/Rapado-Rincón/Kootstra/Matteucci) is the paper your proposal cites. Read its **baseline** table, not its headline: + +| Method (apple orchard, dormancy) | PSNR | ATE (m) | +|---|---|---| +| AgriGS-SLAM (visual + LiDAR) | 29.90 | **0.519** | +| DLO + 3DGS | 25.05 | 0.576 | +| Splat-SLAM | 19.69 | 5.26 | +| OpenGS-SLAM | 12.97 | **20.70** | +| Photo-SLAM | 8.16 | **18.95** | + +Photo-SLAM and OpenGS-SLAM do not degrade in an orchard — they *diverge*, by 19–21 metres. Two CVPR-caliber systems, total failure. AgriGS-SLAM only works because a 32-beam LiDAR is carrying the pose estimate, and even then its own ATE is 0.52 m, which is enormous for a system meant to produce per-fruit geometry. + +Also note what the paper does *not* contain: no GPU model, no runtime, no FPS anywhere in the experiments, and no thin-structure evaluation despite a densification rule explicitly motivated by "leaves and thin branches." That is a conspicuous hole and an opening for you. + +Other agricultural radiance-field work, for context: FruitNeRF (IROS 2024, [2408.06190](https://arxiv.org/abs/2408.06190)) — 12 min to 2.5 h training, authors state plainly it is "not yet suitable for real-time applications," degrades under wind. PAg-NeRF (RA-L 2023, [2309.05339](https://arxiv.org/abs/2309.05339)) — ~27 min training, 5.6 s/image inference, explicitly names thin structures as a challenge. GrowSplat (CASE 2025) — multi-session, offline alignment. + +### 2.4 The decisive argument is representational, not computational + +ActiveGS ([2412.17769](https://arxiv.org/abs/2412.17769)) states it directly: + +> "the Gaussian primitives represent only occupied space, making it hard to distinguish between unknown and free space, which are important for exploration and path planning." + +GS-Planner concurs. For an NBV planner, "I have never looked here" **is** the signal. A representation structurally incapable of expressing it cannot be the primary substrate for exploration no matter how fast it runs. + +The same argument disqualifies **surfel maps** (ElasticFusion, SuMa) entirely, and it partially disqualifies **TSDF**. In a TSDF, weight w=0 conflates three distinct conditions: never observed; observed but beyond truncation in front of the surface; and behind the surface, permanently occluded. That conflates "I should look here" with "I can never see here" — opposite conclusions for a planner. You can disambiguate by space carving, but carving is precisely the mechanism that destroys thin structures. + +**The thin-structure mechanism, precisely.** A TSDF encodes a surface as a zero crossing between a positive and a negative band of thickness τ. This requires the object be thicker than ~2τ. A leaf thinner than τ has back-face rays writing negative-then-positive into the same voxels front-face rays wrote positive-then-negative; the updates destructively average toward zero and **the surface dissolves**. Space carving makes it worse: a ray grazing a leaf edge, or a mixed pixel at a depth discontinuity (endemic on RealSense-class sensors over foliage), carves free space straight through neighbouring leaf voxels. **The failure mode is that the more views you take, the more foliage disappears** — catastrophic for a system whose entire premise is taking more views. + +The cleanest published evidence is wavemap (RSS 2023, [2306.01279](https://arxiv.org/abs/2306.01279)), which ran the controlled comparison and reported that TSDF (voxblox) "has better [surface] reconstruction performance while our approach is better at reconstructing thin objects" — illustrated by the chair that **is missing its legs** in the voxblox reconstruction. Their named target class is "thin objects such as branches, cables, or fences." That maps almost exactly onto twigs, trellis wires, and canopy support structure. + +By contrast, an **occupancy log-odds map with a proper beam model degrades gracefully**: a leaf thinner than a voxel produces a voxel converging to intermediate occupancy — which is *correct*, and which is *high entropy*, so an NBV planner naturally wants to re-observe it. A TSDF forced to place a zero crossing with insufficient thickness places nothing. For a leaf-dominated scene that asymmetry is the whole argument. + +3DGS has the analogous bias for a different reason: the photometric loss is nearly as well satisfied by a fat, low-opacity Gaussian that blurs a leaf into its background as by a correctly thin one, and the latter has weaker gradient support. **3DGS systematically represents thin high-frequency geometry as semi-transparent haze.** It renders beautifully and is geometrically wrong — exactly the failure you cannot tolerate when the downstream consumer is a grasp planner deciding whether an apple is occluded. + +### 2.5 What I would actually build + +**Primary substrate: a GPU-resident, multi-resolution occupancy map with explicit unknown state.** + +- **Three-state semantics from UFOMap** (RA-L 2020, [2003.04749](https://arxiv.org/abs/2003.04749)): occupied / free / **unknown** as first-class states, with inner nodes tagged by which states their subtree contains. This makes "does this frustum contain unknown volume?" an O(depth) query instead of O(leaves), and makes frontier extraction native. The reason every agricultural NBV paper reports 0.5–1.5 s planning cycles is that they all use OctoMap, where unknown is represented by *absence* and must be discovered by raycasting one step at a time. +- **Sensor model from wavemap**: uncertainty-aware beam model over both range σ_r *and angular* σ_θ. Their ablation shows this specific component is what recovers thin structures. It also makes voxel entropy H(p) a real quantity, so expected information gain E[ΔH] is well-posed rather than a heuristic unknown-voxel count. +- **Dual resolution**, following Freeman & Kantor's apple-fruitlet work ([2309.13669](https://arxiv.org/html/2309.13669)): coarse ~1–2 cm map over the whole canopy, fine ~3 mm maps instantiated on demand as "attention regions" around fruit clusters. GS-NBV runs a 3 mm semantic OctoMap for avocado at ~1.5 s/cycle on a laptop GPU, so millimetre-scale occupancy is demonstrably tractable. +- **Semantic layer**: per-voxel class posterior, fused from 2D detections/segmentations projected through known poses. +- **GPU-batched information-gain raycasting.** This is the single highest-leverage engineering decision in the project. In GS-NBV, viewpoint evaluation is 1.473 s of a 1.576 s cycle — **93% of the loop**. It is embarrassingly parallel. And **no agricultural NBV paper has done it**, because they are all on CPU OctoMap. That is your speedup and it is large. + +**Secondary, asynchronous: splatting for appearance and refinement scoring.** Run off the critical path, gated by an HGS-Planner-style coverage weight λ_o (the fraction of observed voxels in the local region): while a region is unexplored, drive purely on voxel IG; as λ_o → 1, let a Gaussian-derived quality term take over to decide which already-observed surfaces need a better look. HGS-Planner ([2409.17624](https://arxiv.org/abs/2409.17624)) writes this as G = G(C) + λ_o·G(Q), which is the field's own equation for "Fisher information only works once coverage is achieved." + +Two notes on the splatting branch: +- If your deliverable is canopy *surface/structure* rather than novel views, use **2DGS**, not 3DGS (cf. Sparfels, [2505.02178](https://arxiv.org/abs/2505.02178)). +- **Feed-forward 3DGS is genuinely interesting for you and under-explored in agriculture**, for one specific reason: three cameras with known 5-DOF poses give you a *calibrated multi-view rig with commandable baseline* — the ideal input regime for MVSplat ([2403.14627](https://arxiv.org/abs/2403.14627), ~45 ms/forward pass) or GS-LRM-class models. A 45 ms forward pass producing canopy Gaussians is genuinely real-time. The catch is that all of these are trained on RealEstate10K/ScanNet++/DTU — indoor real estate, not foliage — and cost-volume methods will suffer most because photo-consistency matching fails on repetitive leaf texture. Helios is exactly the fix: you can render unlimited canopy views with perfect ground truth and *retrain or fine-tune*, and you have the ground truth to evaluate it honestly. That is a legitimately strong paper on its own. + +**The trick worth stealing** from GS-Planner / HGS-Planner / MAGICIAN: rather than keeping two maps separate and fusing scores at the end, **make unknown space renderable** — inject unknown voxels as primitives into the same rasterizer as the Gaussians. Then coverage gain and quality gain both come from *one render pass per candidate pose*, at rasterizer speed, instead of a CPU raycast plus a GPU render. + +**Cheapest viable fallback:** ActiveGS's confidence accumulator, k_i = γ_i·exp(β_i) with γ_i = Σ_j (1 − d_ij/d_far)·n_i·v_ij. No gradients, no backward pass — a running per-primitive accumulator that captures the two things that matter for apples: were you close and face-on, and did you see it from diverse angles. + +### 2.6 What "real-time" actually needs to mean here + +Define it as **p99 decision latency ≤ the deadline it must meet**, and report mean/p95/p99/max. Mean latency conceals exactly what breaks closed-loop systems. + +Then note the structural slack in your favour: **arm motion dominates.** A 5-DOF viewpoint change takes 0.5–3 s (gimbal retargets in a few hundred ms; linear stages are slow). An NBV decision taking 200–500 ms is effectively *free* — it hides under motion. The reason to go faster is not loop rate; it is that you can evaluate **more candidate viewpoints** in the same wall-clock, which is what actually improves NBV quality. **Budget for candidate count, not for latency.** + +And: **do not integrate at frame rate.** Three cameras at 30 Hz is 90 frames/s, but with known poses and a mostly-static canopy you only need to integrate when a camera has *arrived somewhere new*. Frame-rate integration is a SLAM habit, not a requirement. Keyframe-on-arrival (~1–3 Hz per arm) cuts integration load by an order of magnitude and buys back the resolution you need for leaves. + +**Recommended budget:** + +| Stage | Target | Ceiling | +|---|---|---| +| Map integration, per keyframe per camera | 30–100 ms | — | +| IG evaluation (all candidates, all arms) | 200–500 ms | 1.5 s | +| Motion planning to selected viewpoint | 20–50 ms | — | +| **End-to-end per viewpoint decision** | **~250 ms** | **~2 s** | + +against 0.5–3 s of arm motion. Healthy headroom. For reference, FUEL does full UAV exploration replanning in 24 ms on a laptop CPU; HGS-Planner does a full cycle in 129–151 ms; the agricultural literature sits at 1.5–20 s. + +--- + +## 3. Moving the cameras: exploiting the 5-DOF cost asymmetry + +### 3.1 The structural insight + +Your five DOF per camera decompose into two classes with radically different cost: + +| DOF class | Axes | Motion cost | Workspace effect | +|---|---|---|---| +| **Linear** | left-right, up-down, in-out | 0.5–3 s, collision-relevant, shared corridor | Changes *position* → changes what is occluded | +| **Gimbal** | pan, tilt | ~0.1–0.3 s, collision-free, no shared resource | Changes *what is in frame* → cannot resolve occlusion | + +This asymmetry is not incidental — it maps onto a real algorithmic split. Position changes are what defeat occlusion; orientation changes only re-aim. So: + +**Outer loop (expensive, sampling-based): plan positions over the 3 linear axes.** +**Inner loop (cheap, gradient-based): refine orientation over the 2 gimbal axes.** + +For the inner loop you have two ready-made options, and they are the same idea arrived at from opposite directions: + +- **Burusa et al., Gradient-based Local NBV** (ICRA 2024, [2311.16759](https://arxiv.org/html/2311.16759)). Make the utility differentiable through the ray sampling and do gradient ascent on camera pose — visual-servoing the camera up the semantic-IG gradient. Reported vs their own sampling-based semantic NBV: **10× less computation, 28% more efficient trajectories, equivalent accuracy.** +- **Lehnert et al., 3D Move-to-See** (IROS 2019, [1809.07896](https://arxiv.org/abs/1809.07896)). Estimate the *spatial gradient* of a semantic objective by finite-differencing across a physical camera array, then servo along it. **This is structurally exactly your three-camera rig.** You can estimate the gradient without moving anything, then move all three arms along it. No map, no ray casting, one control cycle. + +Gimbal motion is nearly free relative to linear-axis motion, so this nested refinement is essentially a free 10–30% improvement. **I could not find this decomposition published, and it matches your hardware precisely. That is a paper.** + +### 3.2 Precompute a reachability/visibility roadmap + +Do not run IK and collision checking inside the planning loop. Once per cart position: + +- Sample the 5-DOF joint space densely per arm → forward kinematics → camera pose. +- Reject self-collision and canopy collision (Helios's `CollisionDetection` plugin does this; see §6). +- Nodes = surviving poses with IK cached. Edges = k-NN in joint space, **edge weight = actual execution time** under your trapezoidal velocity profiles. + +Because your linear axes are decoupled and gimbal motion is cheap, edge cost is nearly analytic — you don't need MoveIt in the loop. This is Zaenker et al.'s graph-based VMP structure (IROS 2023, [2303.03048](https://arxiv.org/abs/2303.03048)) adapted to a much simpler kinematic chain, and it turns planning into graph search. + +**Two practical findings from the literature you should honour:** +- **Do not cast rays outward from targets to generate candidate views.** Zaenker explicitly reports this fails in confined workspaces — reachability rates are too low. Sample view poses *from the reachable workspace, looking at* targets. [arXiv:2412.10515](https://arxiv.org/abs/2412.10515) independently confirms frontier-sampled candidates have poor arm reachability. +- **Normalize gain by execution time, always.** Zaenker names this as their single most important design decision. + +### 3.3 Three-arm coordination + +**Enforce non-overlapping vertical bands per arm.** With three arms in horizontally-stacked cells, this makes arm–arm collision *structurally impossible*, and coordination reduces to (i) managing the shared in-out extension corridor and (ii) deciding what each arm looks at. That trades a small loss of optimality for an enormous reduction in planning complexity, and it is how the Vougioukas hardware already works. + +**For the assignment, use sequential greedy over a submodular objective:** + +``` +A ← ∅ +for i = 1 … 3: # randomize arm order each round + a_i ← argmax_{a ∈ A_i} F(A ∪ {a}) − F(A) + A ← A ∪ {a_i} +``` + +Fifteen lines of code, and it carries a real guarantee: F(A_seq-greedy) ≥ ½·F(A*) for monotone submodular F under a matroid constraint (Fisher–Nemhauser–Wolsey). Redundancy is eliminated *structurally* — once arm 1 commits to a cluster, that cluster contributes zero marginal gain to arms 2 and 3. Randomized ordering avoids systematic bias (Corah & Michael, Auton. Robots 2019). + +**Watch out for a specific hazard here.** Three cameras choosing simultaneously is a *set* selection problem, and **modular criteria have no diminishing returns, so all three arms will converge on the same maximally-informative view.** POp-GS ([2503.07819](https://arxiv.org/abs/2503.07819)) gives direct evidence: FisherRF's criterion fails at batch selection, 18.37 dB vs D-optimality's 24.53 dB, a 6.16 dB diversity gap. Use a genuinely submodular objective — log-det is submodular, plain trace is not. I found no agricultural NBV paper handling multi-camera batch selection correctly. + +**Also exploit what three simultaneous cameras give you that one moving camera cannot:** +1. Instantaneous gradient estimation (§3.1). +2. **Instantaneous occlusion resolution** — a voxel occluded from camera 1 but visible from camera 2 is resolved *now*, not after a motion. The joint visibility V(ξ₁,ξ₂,ξ₃) is what matters. +3. Immediate stereo baseline widening for depth on thin structures where RGB-D fails. +4. **Cross-view data association at zero latency** — three simultaneous views of one apple make instance re-ID far easier than three sequential views with motion between them. + +### 3.4 The cart is an outer loop, and it makes the problem time-dependent + +A viewpoint available now may be unreachable after the cart advances. Solve over an **overlapping window** of the row rather than the whole row, re-solving as the cart moves — the MPC-style approach your own group already uses for harvest scheduling. + +Cart speed then becomes a decision variable coupling everything: slower = more views per tree = better recall = lower throughput. **Plot the recall-vs-throughput Pareto front as a function of cart speed.** That is the plot a harvesting audience actually wants, and it is not in any active-vision paper I found. + +--- + +## 4. Exploration then exploitation + +### 4.1 Why they need different objectives, not different weights + +This is the most important algorithmic idea in this document. + +Ericson, Molina & Jensfelt, *Information Gain Is Not All You Need* ([2504.01980](https://arxiv.org/abs/2504.01980)), make a sharp argument. If you are **quality-constrained** — you *must* cover everything — then the total information to be gathered is **a constant**, fixed by the requirement. Maximizing per-step gain therefore just *reorders the same work*, and does so badly: it chases distant high-entropy regions and accumulates frontier "debt" requiring backtracking. Their numbers: a simple distance-advantage heuristic gave 16% shorter paths than nearest-frontier, while IG maximization was **23% worse** than nearest-frontier, with 2× the frontier debt. + +If you are **budget-constrained** — fixed time, take the most valuable subset — then IG maximization is exactly right. + +Your Phase 1 is quality-constrained (see the whole tree). Your Phase 2 is budget-constrained (spend the remaining time where the apples are). **So they should not share an objective function.** This is a clean, defensible architectural claim and, as far as I can tell, nobody has framed agricultural active vision this way. + +### 4.2 Phase 1 — EXPLORE + +**Objective:** min T(P) subject to ∪_{ξ∈P} V(ξ) ⊇ S_reachable. + +**Define "seen" honestly.** A surface element s is *covered* iff there exists a taken view ξ with s ∈ FOV(ξ), unoccluded, at range d ∈ [d_min, d_max], and incidence angle θ_s ≤ θ_max. The range gate and incidence-angle gate are essential for RGB-D on thin structures and are routinely omitted in papers — which is one reason simulated results don't transfer. + +**You cannot guarantee 100% coverage** and should not claim it. Some canopy surface is geometrically unreachable — interior branches fully enclosed by foliage. What you *can* guarantee is coverage of the **reachable-visible subset**, and the frontier criterion determines this automatically: when no surface frontier is reachably visible, you are done. Define the guarantee that way and it is defensible. + +**Frontier definition must be re-derived** — you are covering a *surface*, not exploring free space: +- **Surface frontier:** occupied voxel adjacent to unknown → the tree surface is not yet delineated here. +- **ROI frontier:** apple-labelled voxel adjacent to unknown → there may be more apple here. +- **Free-space frontier:** the classical one, for gross workspace coverage. + +**Planner:** receding-horizon path search over the roadmap, utility + + U(ψ) = |∪_{ξ∈ψ} V_new(ξ)| / T(ψ) + +Execute the first 1–2 edges, replan. Horizon 3–5 views. This is Bircher et al.'s RH-NBVP structure with a coverage numerator and a time-normalized denominator. + +The coverage function F(A) = |∪ V(ξ)| is **monotone submodular** by construction, so greedy gives (1−1/e) ≈ 0.632 of optimal for the cardinality-constrained problem, and cost-benefit greedy with the CELF lazy-evaluation trick gives ½(1−1/e) under a time budget. Use CELF — because marginal gains are non-increasing you keep a max-heap of stale gains and only recompute the top element, typically 10–100× cheaper. + +**Be honest about the guarantee.** (1−1/e) holds for the *set* problem. Once you add ordering and travel you are in orienteering territory and the best known bound for submodular orienteering is O(log OPT). Don't overclaim. + +**Baselines you must beat**, and one you might not: +- Boustrophedon raster over the 3 linear axes at fixed gimbal angles. **This is the real commercial competitor.** +- Ericson's "distance advantage" nearest-frontier variant — five lines of code, and it beat IG maximization by ~39 percentage points of relative path length. Implement it; it may win. + +> ⚠️ **A risk you should test in week 2, not month 6.** For a fruiting-wall / V-trellis canopy with a cart moving along the row, the geometry is quasi-planar, and raster scanning is genuinely competitive. The tomato-greenhouse papers report large NBV gains partly because tomato plants are fully 3D. **If the explore-phase gap over raster turns out small, say so early and shift the emphasis to the exploit phase**, where the gains from occlusion reasoning are large regardless of canopy planarity. Reviewers will ask about this. + +### 4.3 The switch — three criteria, take the disjunction + +**(a) Frontier exhaustion.** No reachable, visible surface frontiers remain. Hard guarantee, no tuning. Downside: the tail is expensive, since the last few frontiers are deep in the canopy. Mitigate with a reachability filter and a per-frontier attempt limit. + +**(b) Marginal value theorem.** Because F is submodular, marginal gain is non-increasing in expectation. Stop when + + ΔF(ξ_t) / c(ξ_t) < η + +i.e. when *gain per unit time* falls below a threshold. Set **η equal to the predicted gain rate of the exploit phase**, and the switch happens automatically at the moment exploiting becomes more valuable per second than exploring. This is the economically correct rule and it unifies the two phases with one number. **Recommend as primary.** + +**(c) Good–Turing / Chao1 coverage estimate on apple discoveries.** Track new *apple instances* found. This is formally a species-discovery problem. With f₁ apples seen in exactly one view and f₂ seen in exactly two, the estimated number of *unseen* apples is + + f̂₀ ≈ f₁² / (2f₂), Ĉ = 1 − f₁/n + +Stop exploring when Ĉ > 0.95. **I found no use of this in robotic active vision.** It gives you a statistically principled answer to "have I found all the apples?" that depends on no voxel map, no reconstruction, and no tuning. It also directly answers the proposal's stated interest in "the rate of new fruit discovery." **This is a genuinely novel and cheap contribution — do it.** + +Also enforce a **hard cap** T_explore ≤ ρ·T_total. Start at ρ = 0.4 and **sweep it** — the explore/exploit budget split is a first-class experimental variable and one of your best plots. + +*A fourth option worth knowing:* Zaenker's VMP uses no hard switch at all — instead a **probabilistic mixture** over the three frontier types with user-set weights. That is a softer blend and arguably more elegant. Compare against it. + +*A fifth, cheap and Helios-native:* train a **learned map-completeness estimator** (Luperto et al., [2406.13482](https://arxiv.org/html/2406.13482) — 92.9% accuracy, 31–37% time saved). You have unlimited ground-truth partial-vs-complete tree pairs, so training one is nearly free. + +### 4.4 Phase 2 — EXPLOIT + +**Objective:** max F(P₁ ∪ P₂ ∪ P₃) s.t. T(P_j) ≤ B_j. This is a **submodular team orienteering problem**. Your own group already speaks this language — Vougioukas's lab formulated multi-arm harvest scheduling as a time-dependent team orienteering problem and later as a minimum-makespan VRP (Zhu & Vougioukas, [2505.10028](https://arxiv.org/html/2505.10028)). **Reuse that formulation for *views* rather than *picks*.** The mapping is exact and it gives you a decoupling template (assignment → timing → trajectory) plus their yielding rules for free. + +**Value function** — make it submodular by construction: + + F(A) = Σ_i w_i · φ( Σ_{ξ∈A} q_i(ξ) ), φ(z) = 1 − e^(−z) + +where +- **q_i(ξ)** = view quality of apple i from ξ = (visible surface fraction) × (incidence-angle term) × (range validity gate) × (predicted detector confidence) +- **w_i** = apple priority: low observation count, high current occlusion, high size-estimate uncertainty, or high prior p̂_apple +- **φ concave ⇒ F submodular ⇒ greedy guarantees apply.** The fourth view of an apple is worth less than the second, which is exactly right. + +**Include predicted-but-never-seen apples** with a discount factor. Which brings us to the strongest original idea available here. + +### 4.5 The fruit-occupancy prior — where to look for fruit you have never seen + +All existing shape-completion work (NBV-SC's superellipsoids, Pred-NBV's PoinTr-C, DM-OSVP++'s diffusion) completes *partially visible* fruit. **Nobody predicts the location of fruit that is 100% occluded.** Yet that is the dominant failure mode in apple harvesting — your own lab addresses it with physical foliage agitation precisely because vision alone can't see fully-hidden fruit. + +**Helios makes this uniquely tractable.** You have complete ground-truth canopies with known fruit positions, so you can learn + + p̂(apple at x | observed canopy geometry near x) + +conditioned on local branch structure, foliage density, height, and distance from the trunk — all strongly predictive, because apples grow on spurs at characteristic positions relative to branch architecture. Then the exploit-phase gain becomes + + I_exploit(ξ) = Σ_x P_vis(x) · p̂_apple(x) · H(x) + +i.e. you preferentially look where apples are *likely*, even where you have never seen one. No published work does this; it is Helios-native; and it directly serves the harvesting objective rather than a generic reconstruction objective. + +A useful simplification: for apples, a diffusion model is overkill. Apple *shape* is a very strong prior — near-spheres of known size distribution — so a parametric sphere/superellipsoid prior with a learned size distribution gets you 90% of the benefit at 1% of the cost. + +### 4.6 Data association: what "seen at least once" is defined over + +**This is a gap you must close explicitly.** With three cameras and a moving cart, "have I already seen this apple?" is a *data-association* problem, not a mapping problem. Count apples in a voxel map and you will double-count under registration error; track instances and you will not. + +**Define your coverage guarantee over tracked apple instances, not voxels.** Maintain a 3D multi-object track database with re-ID descriptors, per-instance observation count, confidence, estimated pose and size. Rapado-Rincón's Wageningen thesis is the reference; Wang et al. 2025 (the paper your proposal already cites) makes the point that multi-view active vision is only useful *if* you can associate detections across views — otherwise extra views produce duplicate counts rather than better recall. + +Budget real effort for this. It is unglamorous and it will otherwise silently corrupt every headline number. + +### 4.7 Summary of the planner architecture + +``` +CART (outer, MPC over overlapping row window; speed = decision variable) + │ + ├─ SHARED REPRESENTATION (one process, all 3 arms) + │ • coarse occupancy octree, ~2 cm, 3-state, + semantic posterior + │ • fine octree, ~3 mm, instantiated per fruit cluster + │ • apple instance track database ← "seen once" defined here + │ • fruit-occupancy prior p̂(apple | local canopy geometry) + │ • reachability roadmap per arm (IK + collision precomputed, + │ edge weight = real execution time) + │ • [async] 2DGS/3DGS appearance map, gated by λ_o + │ + ├─ PHASE 1: EXPLORE — quality-constrained + │ objective: min time s.t. coverage ⊇ reachable surface + │ method: submodular max-coverage + CELF + RH path search + │ coordination: sequential greedy, randomized arm order, ½ guarantee + │ switch: frontier exhaustion ∨ marginal-value-rate ∨ Good–Turing Ĉ>0.95 + │ (hard cap ρ·T_total, sweep ρ) + │ + └─ PHASE 2: EXPLOIT — budget-constrained + objective: max Σ w_i φ(Σ q_i(ξ)) s.t. T(P_j) ≤ B_j + = submodular team orienteering + method: cost-benefit greedy + CELF over roadmap paths + coordination: sequential greedy over paths + inner loop: gimbal-only gradient refinement (free) + termination: budget spent ∨ all apples above per-instance threshold +``` + +--- + +## 5. Foundation models — where to start, and what to do when it fails + +### 5.1 The headline: don't start from VGGT + +VGGT won CVPR 2025 Best Paper and it is the reference everyone knows, but it is **batch, bidirectional, and pose-solving** — three properties you specifically do not want. It has no pose-input path, memory grows quadratically (5.6 GB at 20 views → 40.6 GB at 200), and it is now beaten on accuracy by three successors. + +**Your defining asset is that you know the poses.** In Helios exactly; on the robot from arm encoders plus hand-eye calibration. Three models released since Sept 2025 accept poses, intrinsics, and/or depth as *optional conditioning inputs*: + +| Model | Pose conditioning | Streaming | License | Why it matters | +|---|---|---|---|---| +| **MapAnything** ([2509.13414](https://arxiv.org/abs/2509.13414)) | intrinsics as ray directions, poses as quat+t, depth as ray depth — 12+ input combos, per-view optional | No (batch) | **Apache-2.0 code AND `map-anything-apache` weights** | **Primary recommendation.** Factored output (depth + local raymaps + poses + metric scale) means capacity goes to depth, not pose. 2000 views on 140 GB. Training framework can fine-tune VGGT/π³/MoGe-2 as guests — one Helios loader serves every later experiment | +| **Depth Anything 3** ([2511.10647](https://arxiv.org/abs/2511.10647)) | explicit "Pose-Conditioned Depth Estimation" mode, `use_ray_pose` | **Yes** — DA3-Streaming, ultra-long video <12 GB | Apache-2.0 for BASE/SMALL/METRIC-LARGE/MONO-LARGE | Only model answering both pose-conditioning *and* streaming. +23–25% geometric accuracy over VGGT. Also finds a plain DINO transformer suffices — you can fine-tune without inventing architecture | +| **Pi3X** (Dec 2025, in [yyfz/Pi3](https://github.com/yyfz/Pi3)) | optional poses, intrinsics, depth | No | code BSD-3, weights CC-BY-NC | Base π³ is permutation-equivariant with no reference-view dependence and 56% better ATE than VGGT on Sintel, 10× more stable across runs. Use base π³ as the *unconditioned control* that pose-conditioning must beat. Note: no standalone Pi3X paper | + +Include VGGT anyway as the citation anchor (`VGGT-1B-Commercial` is license-clean), but budget about a day for it. + +### 5.2 Two theoretical reasons to expect the uncalibrated path to fail on your scene + +**(a) Projective ambiguity under small baselines.** VGGT-SLAM ([2505.12549](https://arxiv.org/abs/2505.12549)) reports: + +> "the feed-forward nature of VGGT with uncalibrated cameras introduces a **projective ambiguity**, which in addition to the Sim(3) DOF includes **shear, stretch, and perspective DOF, especially when the disparity between frames becomes small**." + +Submaps cannot be aligned by a 7-DoF similarity; they need a **15-DoF SL(4) projective homography**. For a camera on a 5-DOF arm orbiting a canopy through a narrow arc, the output is not merely noisy — it is *structurally ambiguous up to projective warp*. Your canopy comes out geometrically self-consistent and **silently metrically wrong**, corrupting canopy volume, fruit diameter, branch angle, internode length. **Known extrinsics eliminate this by construction.** This is the single best argument for the pose-conditioned path and it belongs in your introduction. + +**(b) Attention collapse under self-similar tokens.** [arXiv:2512.21691](https://arxiv.org/abs/2512.21691) formalizes VGGT's attention degeneration: global attention matrices become near rank-one, token geometry degenerates to an almost one-dimensional subspace, entropy and effective rank decay as O(1/L) in depth, and error accumulates super-linearly. The diffusion coefficient scales inversely with token count, so collapse *accelerates* with more tokens. **A dense orbit of a leaf-covered canopy is the worst case: maximal token count, minimal token diversity.** + +### 5.3 The domain evidence + +**Nobody has benchmarked this model family on apple canopies, orchards, or fruit trees.** I searched hard. The only plant-domain evaluation is [arXiv:2607.01753](https://arxiv.org/abs/2607.01753) (July 2026), testing VGGT and π³ on 26 sequences of maize, tobacco, wheat, soybean, bamboo, rapeseed, pea, broccoli. Their results are genuinely encouraging *in their regime* — initialization 6.52 min (COLMAP) → 1.58 s (π³), leaf area R² 0.936–0.944, mean leaf-angle error ≈2.04°, and COLMAP fails outright below ~30–45 views while the π³ pipeline stays usable far below that. + +But read their protocol: **closed-loop spiral orbit, 80–100 frames, single isolated potted plant.** And their own limitations section: *"repetitive textures, thin organs, occlusions, and uneven views reduce matching stability"; "Dense canopies, dynamic disturbances, cluttered backgrounds, and persistent occlusion remain challenging"; most reliable for "single-plant, close-range, closed-loop acquisition"; **"not suitable for dense field canopies."*** They fence off exactly your regime. + +Supporting out-of-domain evidence: +- **Aerial photogrammetry evaluation** ([2507.14798](https://arxiv.org/abs/2507.14798)): in the *sparse* regime these models win decisively (1 image: 0.36–0.70 m accuracy where COLMAP fails entirely; 38 images at 10% overlap: VGGT holds 35–59% completeness while COLMAP crashes from 24% to 8%). In the *standard* regime they lose badly (~70% overlap: learned 0.44–1.12 m vs COLMAP-HR 0.06–0.16 m; MASt3R camera position error 8.22–62.14 m; orientation error up to 122.6°). Effective resolution ceiling 518 px. +- **E3D-Bench** ([2506.01933](https://arxiv.org/abs/2506.01933)), 16 geometric foundation models: "current GFMs excel on simpler sub-tasks but struggle as complexity grows"; pair-view beats multi-view; **metric-scale prediction is where they fail worst**; none are real-time. Contains **no vegetation domain at all**. +- **Fruit-tree canopy reconstruction review** (Agronomy 2026, [10.3390/agronomy16131274](https://www.mdpi.com/2073-4395/16/13/1274)): SfM+MVS completeness in heavily occluded inner canopy "typically falls below 60 percent"; point-cloud networks reach >85% accuracy on primary branches but **recall below 55% for branches thinner than 10 mm** under leaf occlusion; "weak and repetitive texture of branches leads to poor feature matching stability"; wind causes "fractures in fine branches thinner than 8 mm." **DUSt3R/MASt3R/VGGT: not mentioned at all.** + +### 5.4 The reframe that opens the best toolbox + +**With exact poses and RGB-D, your problem is no longer "3D reconstruction" — it is multi-view stereo / depth fusion / depth completion.** That reframing is worth taking seriously, because it opens a much better-conditioned toolbox that cannot be corrupted by pose error at all: classical MVS, TSDF/occupancy fusion, and modern depth-completion models. + +The one to know: **PromptDA** ([2412.14015](https://arxiv.org/abs/2412.14015), Apache-2.0, ViT-S is 25 M params) takes RGB + sparse/low-res metric depth and produces 4K metric depth. That is *precisely* the primitive for cleaning up RealSense/ZED depth on a canopy. + +**Treat "does a pose-conditioned geometry model actually beat a well-tuned depth-completion + fusion pipeline?" as a live open question that your Helios stage can answer definitively.** If the classical pipeline wins, that is a more useful finding than "the foundation model worked." + +### 5.5 Complementary models, with license traps flagged + +| Role | Pick | License | Note | +|---|---|---|---| +| Metric monocular depth | **MoGe-2** ([2507.02546](https://arxiv.org/abs/2507.02546)) | **MIT** | 60 ms on A100/3090, ViT-S is 35 M. Best license/capability trade in the category | +| | Metric3D v2 | BSD-2 | ONNX export | +| Depth completion | **PromptDA** | Apache-2.0 | RGB + sparse depth → 4K metric | +| Segmentation / tracking | **SAM 2.1** | Apache-2.0 | tiny 91 FPS / large 40 FPS on A100 | +| | **EfficientSAM3** ([2511.15833](https://arxiv.org/abs/2511.15833)) | Apache-2.0 | Distills SAM 3's text-prompted concept segmentation to ~90 M params. Best route to SAM 3 capability without SAM 3's gated custom license | +| Dense features | DINOv2 (Apache) / DINOv3 (custom Meta, mandatory "Built with DINOv3" attribution that propagates) | — | See warning below | +| Matching | **XFeat** ([2404.19174](https://arxiv.org/abs/2404.19174)) Apache-2.0, or LightGlue + **DISK/ALIKED** | — | See warning below | +| Point tracking (wind) | **TAPNext++** ([2604.10582](https://arxiv.org/abs/2604.10582)) | repo Apache, weights CC-BY-SA | Causal, constant memory, 191 FPS @1024 points, explicit re-detection after occlusion | + +⚠️ **License traps that will bite you:** +- **DUSt3R, MASt3R, and anything inheriting their checkpoints (including MASt3R-SLAM) are CC-BY-NC-SA with additional MapFree restrictions.** Fast3R is FAIR non-commercial. π³ weights are CC-BY-NC. +- **SuperPoint and SuperGlue weights are academic/non-commercial only.** The ubiquitous LightGlue+SuperPoint combo is therefore non-commercial, even though LightGlue itself is Apache-2.0. Swap in DISK (Apache) or ALIKED (BSD-3). +- **YOLO-World is GPL-3.0; YOLOE is AGPL-3.0.** Copyleft traps for anything you plan to release under the project's GPLv2 or hand to a startup. +- Default **Depth Anything V2** Base/Large/Giant checkpoints are CC-BY-NC; only Small is Apache. + +⚠️ **One directly relevant empirical warning:** "DINOv3 Visual Representations for Blueberry Perception" ([2603.02419](https://arxiv.org/abs/2603.02419)) finds segmentation quality scales monotonically with backbone size but **detection does not** — ~16% mAP50 with ViT-L, and **cluster detection collapsed to below 2% mAP50.** Apple clusters are the same structure. **Use DINOv3 for dense semantic heads; use SAM-family for instance detection of clustered fruit.** + +### 5.6 Diagnostic experiments — the six that would actually tell you something + +Ordered by information per day. Each targets a specific documented failure mechanism. **This is where I would stop and wait for your Helios setup.** + +**D1 — The pose-conditioning ablation.** Same canopy, same views, four conditions: (a) images only, (b) + known intrinsics, (c) + known intrinsics and extrinsics, (d) + intrinsics, extrinsics, and simulated RGB-D depth. Report Chamfer, completeness, per-organ error. **This measurement does not exist in the literature — MapAnything's own paper never published it.** If (c)/(d) don't substantially beat (a), your central efficiency thesis is wrong and you learn it in week one. It is also a publishable figure. + +**D2 — The baseline-angle sweep.** Fix the canopy; vary arc width: 5°, 10°, 20°, 45°, 90°, 180°, 360° closed loop. Plot accuracy vs arc width, with and without pose conditioning. VGGT-SLAM's projective-ambiguity result predicts a *sharp collapse* at small disparity for uncalibrated inference and *no such collapse* with poses given. **If that prediction holds, it is the core figure of your first paper.** Nobody has published this curve. + +**D3 — Leaf-density / self-similarity sweep.** Helios lets you vary leaf area index continuously on a *fixed branch skeleton* (this is exactly what the proposal's occlusion-regulation module does). Sweep LAI from dormant to dense and measure (i) reconstruction error on the **branch geometry specifically** and (ii) attention effective rank in the aggregator. This tests the attention-collapse mechanism directly and would be the first parametric self-similarity study for this model family. + +**D4 — Thin-structure recall by diameter class.** Bin ground-truth branches by diameter (<5 mm, 5–10, 10–20, >20 mm) and report recall per bin. Benchmark to beat, from the domain literature: **<55% recall below 10 mm under leaf occlusion.** This is the metric an orchard scientist actually cares about. + +**D5 — The honest classical baseline.** With exact poses, run classical MVS / TSDF fusion / PromptDA depth-completion on the same views. **If a well-tuned classical pipeline with known poses matches the foundation model, that is the finding.** + +**D6 — Metric-scale integrity.** Check that predicted canopy volume, leaf area, fruit diameter and internode length are *metrically* correct, not just visually plausible. E3D-Bench found metric scale is where these models fail worst. A projectively-warped reconstruction can look perfect and be metrically useless. + +*Protocol note:* use a **shared global alignment** for joint camera+geometry evaluation (the UAVFF3D convention). Aligning cameras and geometry separately biases results. + +### 5.7 Decision tree for when it fails + +``` +START: D1 + D2 in Helios on MapAnything, DA3, π³ (VGGT as anchor) + +├─ Works (metric error acceptable, thin-branch recall ≥ classical MVS) +│ → move to streaming: DA3-Streaming / LONG3R / HorizonStream, +│ or the ASYNC pattern: heavy model once per settled viewpoint, +│ tiny model during motion (AsyncMDE-style). Close the loop. +│ +├─ Works with poses, fails without ← MOST LIKELY OUTCOME +│ → this IS the contribution. Commit to the pose-conditioned path. +│ Cite the SL(4) projective-ambiguity result as the theoretical reason. +│ Thesis framing: "known-pose geometry models for agricultural +│ active vision." +│ +├─ Fails on foliage (branches fine, leaves wrong) +│ TIER 1: LoRA3D-style self-supervised LoRA on Helios data +│ ([2412.07746]: ~5 min, 18 MB adapter, no labels, self-supervised). +│ Cheapest possible fix — always try first. Per-cultivar, +│ per-lighting adapters then become practical on-robot. +│ TIER 2: full fine-tune via MapAnything's training framework on +│ Helios apple / apple_fruitingwall. +│ • copy the AerialMegaDepth recipe ([2504.13157]): mixing +│ synthetic renders with real images took DUSt3R from <5% +│ to ~56% on extreme-viewpoint localization — an ~11× gain +│ • sample MANY orchard configurations, not many frames of few +│ (TartanGround: environment count, not sample count, is the +│ generalization bottleneck) +│ • include SAGE-style anti-forgetting regularization — +│ narrow fine-tuning WILL destroy general geometry priors +│ • read the Infinigen-Stereo ablation tables ([2504.16930]) — +│ the only systematic study of which randomization axes +│ matter for geometry on *vegetation* scenes +│ TIER 3: hybrid. Model as PRIOR only — initialize bundle adjustment +│ and 2DGS (not 3DGS) with joint pose refinement. +│ Add dynamic-area suppression for wind. +│ +├─ Fails sim→real (works in Helios, fails on the robot) +│ → do NOT immediately add noise to Helios. Try the INVERTED approach +│ first: Camera Depth Models ([2509.02530]) that denoise REAL depth +│ toward sim. That paper shows policies trained on raw simulated +│ depth transfer with no noise augmentation and no real fine-tuning, +│ including on "articulated, reflective, and slender objects." +│ Slender = branches. +│ Then Wat3R-style teacher–student adaptation on unlabelled real +│ orchard video (zero annotations needed). +│ +└─ Works but too slow + (1) systems-level speedups first — Speedy-MASt3R-style, free + (2) QuantVGGT 4-bit: >98% FP accuracy retained, 2.5× hardware + (3) VERIFY token-merging speedups hold at YOUR 5–20 view scale — + they are measured at ~1000 views and will badly underdeliver + (4) only then distill (eVGGT ~9×, Distill3R 5×) or go async +``` + +**One caveat on the fine-tuning literature:** nobody publishes a layer-wise LoRA ablation for VGGT-family models specifically. VGGT's own default config freezes the aggregator and trains heads only; other work suggests that for *modality* shift the problem lives inside the aggregator. This is a genuine open question and a legitimate thesis ablation. + +### 5.8 The dataset gap that justifies the whole Stage 1 + +**There is no public plant or orchard dataset providing multi-view RGB *with camera poses* plus ground-truth dense geometry** — i.e. exactly the tuple a DUSt3R/VGGT-family fine-tune consumes. Every candidate is either LiDAR/TLS point clouds with no images or poses (AgriField3D, PLANesT-3D, TreeScope), or images with 2D labels and at best stereo (AppleGrowthVision, MinneApple, PhenoBench). + +Helios fills exactly this gap, and it emits the complete tuple already: physically-based RGB + lossless float depth + exact intrinsics/extrinsics + per-pixel semantic and instance masks. **This is your Stage 1 justification, stated in one sentence.** + +--- + +## 6. Helios: what you get, what you must build + +I had the repo audited at `v1.3.78` (HEAD, 2026-07-26). Summary: **Helios covers roughly 70% of what this project needs, and covers it unusually well. Every gap is on the robotics side, plus one performance issue.** + +### 6.1 The parts that are better than expected + +**Programmatic multi-camera control is first-class.** `RadiationModel::addRadiationCamera(label, bands, position, lookat, props, AA)`, `setCameraPosition()`, `setCameraLookat()`, `setCameraOrientation()`, `getAllCameraLabels()`, and `runRadiationImaging(vector, ...)`. Cameras are stored as a map and `runBand()` iterates all of them in one dispatch — **your three-arm rig maps directly onto this.** Camera model is thin-lens with depth of field (pinhole if `lens_diameter = 0`), with real intrinsics (resolution, focal length, HFOV, sensor width) and radial/tangential distortion. + +**Real per-pixel float depth**, pushed into Context global data after every render — no file I/O needed: +```cpp +context.getGlobalData("camera_