Optimize TRK reader

.github/workflows/codespell.yml

@@ -28,4 +28,4 @@ jobs:
         uses: codespell-project/actions-codespell@v2
         with:
           path: packages/niivue
-          ignore_words_list: aNormal,inout,aLo,hAx,LOD,mater
+          ignore_words_list: aNormal,inout,aLo,hAx,LOD,mater,lod,nd

AGENTS.md

@@ -175,7 +175,7 @@ the same options CI uses (skip generated dirs that aren't in CI's checkout):
 
 ```bash
 codespell --skip './.git,*/dist,*/node_modules,*/coverage,bun.lock' \
-  --ignore-words-list aNormal,inout,aLo,hAx,LOD,mater \
+  --ignore-words-list aNormal,inout,aLo,hAx,LOD,mater,lod,nd \
   packages/niivue/src packages/niivue/examples packages/niivue/*.md
 ```
 

packages/niivue/examples/index.html

@@ -128,6 +128,7 @@ <h1>Voxel Scenes</h1>
           <li><a href="./vox.basic.html" target="_blank" rel="noopener noreferrer">Basic</a></li>
           <li><a href="./vox.clip.html" target="_blank" rel="noopener noreferrer">Clip planes</a></li>
           <li><a href="./vox.affine.html" target="_blank" rel="noopener noreferrer">Affine Manipulation</a></li>
+          <li><a href="./transform.html" target="_blank" rel="noopener noreferrer">Transform Gizmo</a></li>
           <li><a href="./vox.stats.html" target="_blank" rel="noopener noreferrer">Statistical Overlay</a></li>
           <li><a href="./vox.torso.html" target="_blank" rel="noopener noreferrer">Labeled Torso</a></li>
           <li><a href="./vox.atlas.html" target="_blank" rel="noopener noreferrer">Atlas</a></li>

packages/niivue/examples/transform.html

@@ -0,0 +1,48 @@
+<!DOCTYPE html>
+<html lang="en">
+  <head>
+    <meta charset="utf-8" />
+    <meta name="viewport" content="width=device-width,initial-scale=1.0" />
+    <title>NiiVue - Transform Gizmo</title>
+    <link rel="icon" href="data:,">
+    <link rel="stylesheet" href="/niivue.css" />
+    <style>
+      #canvas-container {
+        display: block;
+        position: relative;
+        flex: 1;
+      }
+      #gizmo {
+        position: static;
+        flex: none;
+        width: 200px;
+        height: 600px;
+        cursor: grab;
+        touch-action: none;
+      }
+    </style>
+  </head>
+  <body>
+    <header>
+      <label for="target">Volume</label>
+      <select id="target">
+        <option value="0">background (mni152)</option>
+        <option value="1" selected>overlay (mpld_asl)</option>
+      </select>
+      &nbsp; Drag a knob around its ring. One turn = 90 mm or 3.6x scale; rotation is limited to &plusmn;30 deg.
+      &nbsp;
+      <button id="resetAffine" type="button">Reset</button>
+      &nbsp;
+      <label for="webgpuCheck">WebGPU</label>
+      <input type="checkbox" id="webgpuCheck" checked />
+      <br />
+      <span id="dials"></span>
+    </header>
+    <main>
+      <div id="canvas-container"><canvas id="gl1"></canvas></div>
+      <canvas id="gizmo" width="200" height="600"></canvas>
+    </main>
+    <footer id="location">NiiVue</footer>
+    <script type="module" src="./transform.js"></script>
+  </body>
+</html>

packages/niivue/examples/transform.js

@@ -0,0 +1,188 @@
+import { mat4, vec3 } from 'gl-matrix'
+import NiiVue from '../src/index.ts'
+
+const RADII = [82, 64, 46] // concentric rings for X, Y, Z
+const CELL = 200 // per-dial height
+const AXES = ['X', 'Y', 'Z']
+const COLORS = ['#e44', '#4c4', '#48f']
+// one knob revolution = 360 units of that dial
+const DIALS = [
+  { name: 'Translate (mm)', unit: 0.25, fmt: (v) => v.toFixed(1) },
+  // ponytail: Euler XYZ gimbal-locks near 90; clamp instead of switching to quaternions
+  { name: 'Rotate (deg)', unit: 1, fmt: (v) => v.toFixed(0), max: 30 },
+  { name: 'Scale', unit: 0.01, fmt: (v) => (1 + v).toFixed(2) },
+]
+
+const nv1 = new NiiVue({ backgroundColor: [0.2, 0.2, 0.2, 1] })
+const ctx = gizmo.getContext('2d')
+const C = gizmo.width / 2
+
+// per volume: dial values, affine at load (rows), centre in world mm at load
+const state = []
+let vol = 0
+let values = null
+let pivot = null
+let hover = null
+let drag = null
+
+const toRows = (m) =>
+  [0, 1, 2, 3].map((r) => [0, 1, 2, 3].map((c) => m[c * 4 + r]))
+const toMat = (rows) =>
+  mat4.fromValues(...[0, 1, 2, 3].flatMap((c) => rows.map((row) => row[c])))
+const rad = (d) => (d * Math.PI) / 180
+
+// world-space T = translate(t) * about pivot( Rx Ry Rz S )
+function buildTransform() {
+  const [t, r, s] = values
+  const m = mat4.fromTranslation(
+    mat4.create(),
+    vec3.add(vec3.create(), t, pivot),
+  )
+  mat4.rotateX(m, m, rad(r[0]))
+  mat4.rotateY(m, m, rad(r[1]))
+  mat4.rotateZ(m, m, rad(r[2]))
+  mat4.scale(m, m, [1 + s[0], 1 + s[1], 1 + s[2]])
+  mat4.translate(m, m, vec3.negate(vec3.create(), pivot))
+  return m
+}
+
+async function apply() {
+  const m = mat4.multiply(
+    mat4.create(),
+    buildTransform(),
+    toMat(state[vol].original),
+  )
+  await nv1.setVolumeAffine(vol, toRows(m))
+  document.getElementById('dials').textContent = DIALS.map(
+    (d, i) => `${d.name}: ${values[i].map(d.fmt).join(', ')}`,
+  ).join('   ')
+}
+
+// 12 o'clock is zero; clockwise increases
+const knobAngle = (dial, axis) =>
+  rad(-90 + values[dial][axis] / DIALS[dial].unit)
+const knobPos = (dial, axis) => {
+  const a = knobAngle(dial, axis)
+  const r = RADII[axis]
+  return [C + r * Math.cos(a), dial * CELL + CELL / 2 + r * Math.sin(a)]
+}
+
+function draw() {
+  ctx.clearRect(0, 0, gizmo.width, gizmo.height)
+  ctx.font = '13px system-ui'
+  ctx.textAlign = 'center'
+  for (const [dial, d] of DIALS.entries()) {
+    const cy = dial * CELL + CELL / 2
+    ctx.lineWidth = 2
+    for (const axis of [0, 1, 2]) {
+      ctx.strokeStyle = COLORS[axis]
+      ctx.globalAlpha = 0.4
+      ctx.beginPath()
+      ctx.arc(C, cy, RADII[axis], 0, Math.PI * 2)
+      ctx.stroke()
+    }
+    ctx.globalAlpha = 1
+    ctx.fillStyle = '#ddd'
+    ctx.fillText(d.name, C, cy - RADII[0] - 8)
+    for (const axis of [0, 1, 2]) {
+      const [x, y] = knobPos(dial, axis)
+      const hot = hover && hover.dial === dial && hover.axis === axis
+      ctx.fillStyle = COLORS[axis]
+      ctx.beginPath()
+      ctx.arc(x, y, hot ? 11 : 8, 0, Math.PI * 2)
+      ctx.fill()
+      ctx.fillStyle = '#000'
+      ctx.fillText(AXES[axis], x, y + 4)
+    }
+  }
+}
+
+function hit(x, y) {
+  for (const dial of [0, 1, 2]) {
+    for (const axis of [0, 1, 2]) {
+      const [kx, ky] = knobPos(dial, axis)
+      if (Math.hypot(x - kx, y - ky) < 12) return { dial, axis }
+    }
+  }
+  return null
+}
+
+const pos = (e) => {
+  const b = gizmo.getBoundingClientRect()
+  return [e.clientX - b.left, e.clientY - b.top]
+}
+const angleAt = (dial, x, y) => Math.atan2(y - (dial * CELL + CELL / 2), x - C)
+
+gizmo.onpointerdown = (e) => {
+  const [x, y] = pos(e)
+  const h = hit(x, y)
+  if (!h) return
+  drag = { ...h, angle: angleAt(h.dial, x, y) }
+  gizmo.setPointerCapture(e.pointerId)
+}
+gizmo.onpointermove = (e) => {
+  const [x, y] = pos(e)
+  if (!drag) {
+    hover = hit(x, y)
+    draw()
+    return
+  }
+  const angle = angleAt(drag.dial, x, y)
+  let delta = angle - drag.angle
+  if (delta > Math.PI) delta -= 2 * Math.PI
+  if (delta < -Math.PI) delta += 2 * Math.PI
+  drag.angle = angle
+  const d = DIALS[drag.dial]
+  const v = values[drag.dial][drag.axis] + ((delta * 180) / Math.PI) * d.unit
+  values[drag.dial][drag.axis] = d.max
+    ? Math.max(-d.max, Math.min(d.max, v))
+    : v
+  draw()
+  apply()
+}
+gizmo.onpointerup = () => {
+  drag = null
+}
+
+function selectVolume(i) {
+  vol = i
+  ;({ values, pivot } = state[i])
+  draw()
+  apply()
+}
+target.onchange = function () {
+  selectVolume(Number(this.value))
+}
+resetAffine.onclick = async () => {
+  for (const v of values) v.fill(0)
+  draw()
+  await apply()
+}
+webgpuCheck.onclick = function () {
+  nv1.reinitializeView({ backend: this.checked ? 'webgpu' : 'webgl2' })
+}
+
+nv1.addEventListener('locationChange', (e) => {
+  document.getElementById('location').textContent = e.detail.string
+})
+await nv1.attachToCanvas(gl1)
+nv1.sliceType = 3
+await nv1.loadVolumes([
+  { url: '/volumes/mni152.nii.gz', calMin: 30, calMax: 80, colormap: 'gray' },
+  { url: '/volumes/mpld_asl.nii.gz', colormap: 'red', opacity: 0.5 },
+])
+for (const [i, v] of nv1.volumes.entries()) {
+  const original = nv1.getVolumeAffine(i)
+  const d = v.hdr.dims
+  state.push({
+    original,
+    values: DIALS.map(() => [0, 0, 0]),
+    pivot: vec3.transformMat4(
+      vec3.create(),
+      [d[1] / 2, d[2] / 2, d[3] / 2],
+      toMat(original),
+    ),
+  })
+}
+nv1.setCrosshairPos([41, -16, 58])
+selectVolume(Number(target.value))

packages/niivue/package.json

@@ -42,7 +42,7 @@
     "codegen:assets": "node scripts/generate-assets.js",
     "gen-fixtures": "bun scripts/generate-range-poc-fixture.ts && bun scripts/generate-tile-range-poc-fixture.ts",
     "ensure-samples": "bun scripts/ensure-samples.ts",
-    "dev": "bun run codegen:assets && bun run gen-fixtures && bun run ensure-samples && bunx --bun vite --port 5273 --strictPort --open /examples/index.html",
+    "dev": "bun run codegen:assets && bun run gen-fixtures && bunx --bun vite --port 5273 --strictPort --open /examples/index.html",
     "dev:perf": "bun run codegen:assets && bunx --bun vite --port 5273 --strictPort --open /examples/benchmark.html",
     "build": "bun run codegen:assets && bunx --bun vite build --config vite.config.lib.ts",
     "build:examples": "bun run codegen:assets && bun run gen-fixtures && bunx --bun vite build --config vite.config.examples.ts",

packages/niivue/src/mesh/tracts/readers/trk.ts

@@ -99,31 +99,41 @@ export async function read(buffer: ArrayBufferLike): Promise<NVTractData> {
   m[14] = vox2ras[14]
   m[15] = vox2ras[15]
 
-  // Parse streamline data
-  const i32 = new Int32Array(buffer.slice(hdr_sz))
-  const f32 = new Float32Array(i32.buffer)
-  const ntracks = i32.length
-  if (ntracks < 1) throw new Error('Empty TRK file')
+  // Parse streamline data: views on the original buffer, no copy
+  const nWords = (buffer.byteLength - hdr_sz) >> 2
+  if (nWords < 1) throw new Error('Empty TRK file')
+  const i32 = new Int32Array(buffer, hdr_sz, nWords)
+  const f32 = new Float32Array(buffer, hdr_sz, nWords)
+
+  // Pass 1: count streamlines and points for exact allocation
+  let nStreamlines = 0
+  let nPoints = 0
+  let idx = 0
+  while (idx < nWords) {
+    const n_pts = i32[idx]
+    if (n_pts < 1) break
+    const next = idx + 1 + n_pts * (3 + n_scalars) + n_properties
+    if (next > nWords) break
+    nStreamlines++
+    nPoints += n_pts
+    idx = next
+  }
+  const nWordsValid = idx
 
-  // Over-provision arrays
-  let vertices = new Float32Array(ntracks)
-  let offsets = new Uint32Array(ntracks / 4)
+  const vertices = new Float32Array(nPoints * 3)
+  const offsets = new Uint32Array(nStreamlines + 1)
+  const dpvVals = dpvNames.map(() => new Float32Array(nPoints))
+  const dpsVals = dpsNames.map(() => new Float32Array(nStreamlines))
+
+  // Pass 2: read
   let npt = 0
   let npt3 = 0
   let noffset = 0
-
-  // Temporary arrays for scalars (converted to Float32Array after)
-  const dpvVals: number[][] = dpvNames.map(() => [])
-  const dpsVals: number[][] = dpsNames.map(() => [])
-
-  let idx = 0
-  while (idx < ntracks) {
+  idx = 0
+  while (idx < nWordsValid) {
     const n_pts = i32[idx]
     idx++
-    if (!Number.isFinite(n_pts) || n_pts < 1) break
-
-    offsets[noffset++] = npt
-
+    offsets[noffset] = npt
     for (let j = 0; j < n_pts; j++) {
       const px = f32[idx],
         py = f32[idx + 1],
@@ -135,33 +145,27 @@ export async function read(buffer: ArrayBufferLike): Promise<NVTractData> {
       vertices[npt3++] = px * m[8] + py * m[9] + pz * m[10] + m[11]
       // Per-vertex scalars
       for (let s = 0; s < n_scalars; s++) {
-        dpvVals[s].push(f32[idx])
-        idx++
+        dpvVals[s][npt] = f32[idx++]
       }
       npt++
     }
     // Per-streamline properties
     for (let s = 0; s < n_properties; s++) {
-      dpsVals[s].push(f32[idx])
-      idx++
+      dpsVals[s][noffset] = f32[idx++]
     }
+    noffset++
   }
 
   // Fence-post: final offset
-  offsets[noffset++] = npt
-
-  // Trim arrays
-  vertices = vertices.slice(0, npt3)
-  offsets = offsets.slice(0, noffset)
+  offsets[noffset] = npt
 
-  // Convert scalar arrays to Record<string, Float32Array>
   const dpv: Record<string, Float32Array> = {}
   for (let i = 0; i < dpvNames.length; i++) {
-    dpv[dpvNames[i]] = Float32Array.from(dpvVals[i])
+    dpv[dpvNames[i]] = dpvVals[i]
   }
   const dps: Record<string, Float32Array> = {}
   for (let i = 0; i < dpsNames.length; i++) {
-    dps[dpsNames[i]] = Float32Array.from(dpsVals[i])
+    dps[dpsNames[i]] = dpsVals[i]
   }
 
   return { vertices, offsets, dpv, dps, groups: {}, dpvMeta: {}, dpsMeta: {} }