Skip to content

Build the cortical ribbon with a winding number instead of mris_volmask - #880

Open
m-reuter wants to merge 4 commits into
Deep-MI:devfrom
m-reuter:volmask
Open

m-reuter wants to merge 4 commits into
Deep-MI:devfrom
m-reuter:volmask

Conversation

@m-reuter

@m-reuter m-reuter commented Sep 13, 2026

Copy link
Copy Markdown
Member

recon-all -cortribbon was 8.07 of the 81 serial minutes in the reference 1.0mm run, more than
topology correction. Its mris_volmask decides for every voxel whether it is inside the white
and pial surfaces; recon_surf/volmask.py answers the same question with a generalized winding
number in about ten seconds.

Not a speed trick: white and pial coincide exactly at 8,653 of 141,883 lh vertices along the
medial wall, so this is not two clean closed surfaces. A winding number is a solid angle
integral per point rather than a fill, so coincident sheets and self-intersections change its
value smoothly instead of letting a region escape.

Results

Against the ribbon.mgz, lh.ribbon.mgz and rh.ribbon.mgz from each subject's CI run:

subject ribbon lh.ribbon rh.ribbon time
1.0mm 119 0 0 9.7 s
0.8mm 86 1 2 14.8 s

Every differing voxel is accounted for, in two groups.

Where inside is undefined (0 and 3 voxels). Inside is a discrete predicate, so two correct
implementations must agree wherever it is defined. The three at 0.8mm lie at distance 0.0000 mm
from a triangle. Nowhere else is it close: 0 query points fall within 0.1 of the 0.5 cut.

The hemisphere overlap (119 and 83), a bias this PR removes. A few hundred voxels near the
medial wall are inside both hemispheres. mris_volmask gives all of them to the left, inflating
left cortex in the same direction in every subject. Each now goes to the hemisphere whose pial
it lies deeper inside:

1.0mm overlap 231 112 lh / 119 rh (was 231 / 0) lh cortex -117, rh cortex +117
0.8mm overlap 168 85 lh / 83 rh (was 168 / 0) lh cortex -80, rh cortex +81

Deltas mirror, so nothing is gained or lost, only reassigned. Four rules were compared; they
agree with each other only 45 to 73% against a 50% baseline, so these voxels are genuinely
ambiguous, and depth was picked for being the only one near even by construction. It costs under
half a second, running on a few hundred points. The per-hemisphere files are unchanged, since
they record each hemisphere's claim before arbitration, as mris_volmask does.

Changes

  • recon_surf/volmask.py, new. Checks the surfaces match the reference grid, takes face
    orientation from the mesh, writes nothing until everything is checked, honours the thread
    budget, and evaluates in 4M point slabs so memory does not grow with resolution.
  • recon_surf/recon-surf.sh, replaces the recon-all -cortribbon call.
  • recon_surf/long_compat_segmentHA.py, was still building the ribbon the slow way, about 5 of
    its 8.5 minutes; now uses the same step.
  • tools/build/, mris_volmask moves from copy_files to link_files, saving 2.7 MB. It
    cannot simply be dropped: every recon-all run probes about 70 binaries with -all-info and
    exits 1 if one is missing, 17 times per run for this one. not-here.sh answers the probe and
    fails loudly on a real call.
  • pyproject.toml, requirements.txt, add libigl.

High resolution

Measured, since 0.7mm is already in use: cost is linear in query count, 21.8 s at 0.7mm and
102.9 s at 0.4mm, still four times faster than mris_volmask is at 1.0mm. Memory was the real
constraint, and slabbing caps a term that would otherwise reach 1.2 GB per surface at 0.4mm.
Slabbing never triggers at 1.0 and 0.8mm, so it was exercised by forcing small slabs: identical
output at 4M, 1M, 250k and 50k on both subjects.

Dependency

libigl, abi3 wheels for every platform we
ship, 13.9 MB on linux. The wheels use std::thread rather than OpenMP on both platforms
(0 of 956 dynamic symbols match GOMP_, omp_, __kmpc) and read IGL_NUM_THREADS;
OMP_NUM_THREADS is set alongside in case of a source build. Method from Jacobson et al.,
SIGGRAPH 2013, with the fast evaluation from Barill et al., SIGGRAPH 2018.

Caveats

  • The embedded colour table is lost, a 12,947 byte footer becomes 20 bytes. Nothing downstream
    reads it; it affects label names in freeview only.
  • Depth is the most symmetric tie-break, not a demonstrably correct one. It removes a fixed
    direction, it does not establish per-voxel truth.
  • Two subjects, no unit tests, and the orientation branch is exercised by nothing in CI.
    Comparing this branch against the dev run at b688279c is the meaningful check, since the
    stored reference tarballs are stale for unrelated reasons.
  • The segmentHA path is checked by reading and compilation, not by a run.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant