concept-collection / mri-spins
MRI spin physics demo
Jeremy Magland <jmagland@flatironinstitute.org> committed commit 158a71609817 Browse files
2 changed files+999−0
README.mdadded+26−0View file
@@ -0,0 +1,26 @@
1+# mri-spins
2+
3+Interactive 3D demonstration of MRI spin physics.
4+
5+**Live: https://concept-collection.github.io/mri-spins/**
6+
7+A sample of isochromats (sphere or cube) evolves under the Bloch equations in
8+the rotating frame. Apply RF pulses, gradients, and a Larmor frequency offset;
9+adjust T1/T2 relaxation; and view the magnetization as vector arrows or
10+colored balls. A live trace shows the signal (real and imaginary parts of the
11+net transverse magnetization).
12+
13+A built-in spoiled gradient echo sequence plays out 16 phase-encoded TRs,
14+driving the controls automatically and acquiring signal during the readout
15+windows. Optional gradient sounds click like a real scanner.
16+
17+## Things to try
18+
19+- Apply a 90° pulse: the spins precess at Δf₀ while T2 attenuates the signal
20+ and Mz recovers with T1.
21+- Turn up Gx to dephase the spins, then apply a 180° pulse to refocus them
22+ into a spin echo.
23+- Turn on sound and run the spoiled gradient echo.
24+
25+Single static page, no build step — three.js is loaded from a CDN. Serve
26+locally with any static file server, e.g. `npx serve .`
index.htmladded+973−0View file
@@ -0,0 +1,973 @@
1+<!DOCTYPE html>
2+<html lang="en">
3+<head>
4+<meta charset="utf-8">
5+<meta name="viewport" content="width=device-width, initial-scale=1">
6+<title>MRI spin physics</title>
7+<style>
8+ html, body {
9+ margin: 0;
10+ height: 100%;
11+ overflow: hidden;
12+ background: #0d1117;
13+ color: #e6edf3;
14+ font-family: -apple-system, BlinkMacSystemFont, "Segoe UI", Helvetica, Arial, sans-serif;
15+ font-size: 13px;
16+ }
17+ #scene { position: fixed; inset: 0; }
18+ #panel {
19+ position: fixed;
20+ top: 12px;
21+ left: 12px;
22+ bottom: 12px;
23+ width: clamp(270px, 26vw, 430px);
24+ overflow-y: auto;
25+ background: rgba(22, 27, 34, 0.92);
26+ border: 1px solid #30363d;
27+ border-radius: 10px;
28+ padding: 14px 16px;
29+ box-sizing: border-box;
30+ backdrop-filter: blur(4px);
31+ }
32+ #panel h1 { font-size: 16px; margin: 0 0 4px; }
33+ #panel .sub { color: #8b949e; margin: 0 0 10px; line-height: 1.35; }
34+ section { border-top: 1px solid #30363d; padding: 10px 0; }
35+ section h2 {
36+ font-size: 11px;
37+ text-transform: uppercase;
38+ letter-spacing: 0.06em;
39+ color: #8b949e;
40+ margin: 0 0 8px;
41+ }
42+ .row {
43+ display: flex;
44+ align-items: center;
45+ gap: 8px;
46+ margin: 6px 0;
47+ }
48+ .row label { flex: 0 0 92px; color: #c9d1d9; }
49+ .row input[type=range] { flex: 1; min-width: 0; accent-color: #58a6ff; }
50+ .row .val {
51+ flex: 0 0 78px;
52+ text-align: right;
53+ font-variant-numeric: tabular-nums;
54+ color: #8b949e;
55+ }
56+ button {
57+ background: #21262d;
58+ color: #e6edf3;
59+ border: 1px solid #30363d;
60+ border-radius: 6px;
61+ padding: 5px 10px;
62+ font-size: 13px;
63+ cursor: pointer;
64+ }
65+ button:hover { background: #30363d; }
66+ button.primary { background: #1f6feb; border-color: #1f6feb; }
67+ button.primary:hover { background: #388bfd; }
68+ button:disabled { opacity: 0.45; pointer-events: none; }
69+ input:disabled { opacity: 0.5; }
70+ button.flash, button.primary.flash { background: #e05555; border-color: #e05555; }
71+ .modes { display: flex; gap: 0; }
72+ .modes button {
73+ flex: 1;
74+ border-radius: 0;
75+ padding: 5px 4px;
76+ font-size: 12px;
77+ }
78+ .modes button:first-child { border-radius: 6px 0 0 6px; }
79+ .modes button:last-child { border-radius: 0 6px 6px 0; }
80+ .modes button.active { background: #1f6feb; border-color: #1f6feb; }
81+ .row .modes { flex: 1; }
82+ .hint { color: #8b949e; line-height: 1.4; margin: 10px 0 0; }
83+ #timeReadout { font-variant-numeric: tabular-nums; color: #8b949e; margin-left: auto; }
84+ .legendRow { display: flex; align-items: center; gap: 10px; margin-top: 8px; }
85+ .legendRow span { color: #8b949e; line-height: 1.3; }
86+ #signalPanel {
87+ position: fixed;
88+ left: calc(clamp(270px, 26vw, 430px) + 24px);
89+ right: 12px;
90+ bottom: 12px;
91+ height: 130px;
92+ background: rgba(22, 27, 34, 0.85);
93+ border: 1px solid #30363d;
94+ border-radius: 10px;
95+ overflow: hidden;
96+ }
97+ #signalPanel canvas { position: absolute; inset: 0; width: 100%; height: 100%; }
98+ #signalLegend {
99+ position: absolute;
100+ top: 6px;
101+ left: 12px;
102+ color: #8b949e;
103+ pointer-events: none;
104+ }
105+ #signalLegend b { font-weight: 600; }
106+ #axesLegend {
107+ position: fixed;
108+ right: 14px;
109+ top: 12px;
110+ color: #8b949e;
111+ background: rgba(22, 27, 34, 0.7);
112+ border-radius: 6px;
113+ padding: 4px 8px;
114+ }
115+ #axesLegend b { font-weight: 600; }
116+</style>
117+<script type="importmap">
118+{
119+ "imports": {
120+ "three": "https://cdn.jsdelivr.net/npm/three@0.170.0/build/three.module.js",
121+ "three/addons/": "https://cdn.jsdelivr.net/npm/three@0.170.0/examples/jsm/"
122+ }
123+}
124+</script>
125+</head>
126+<body>
127+<div id="scene"></div>
128+
129+<div id="panel">
130+ <h1>MRI spin physics</h1>
131+ <p class="sub">Bloch simulation of a sample of isochromats, viewed in the rotating frame. B<sub>0</sub> is along z.</p>
132+
133+ <section>
134+ <h2>Simulation</h2>
135+ <div class="row">
136+ <button id="playBtn" class="primary">Pause</button>
137+ <button id="resetBtn">Reset M</button>
138+ <span id="timeReadout">t = 0.0 s</span>
139+ </div>
140+ <div class="row">
141+ <label for="speed">Speed</label>
142+ <input type="range" id="speed" min="-1" max="0.6" step="0.01" value="-0.6">
143+ <span class="val" id="speedVal">0.25×</span>
144+ </div>
145+ </section>
146+
147+ <section>
148+ <h2>Display</h2>
149+ <div class="row">
150+ <label>Sample</label>
151+ <div class="modes">
152+ <button id="sampleSphere" class="active">Sphere</button>
153+ <button id="sampleCube">Cube</button>
154+ </div>
155+ </div>
156+ <div class="row">
157+ <label>Show</label>
158+ <div class="modes">
159+ <button id="shapeBalls">Balls</button>
160+ <button id="shapeArrows" class="active">Arrows</button>
161+ </div>
162+ </div>
163+ <div class="row">
164+ <label>Color</label>
165+ <div class="modes">
166+ <button id="colorTransverse" class="active" title="Transverse magnetization: hue = phase, saturation = magnitude">M<sub>xy</sub></button>
167+ <button id="colorLongitudinal" title="Longitudinal magnetization">M<sub>z</sub></button>
168+ <button id="colorNone" title="Uniform color">None</button>
169+ </div>
170+ </div>
171+ <div class="row">
172+ <label for="resolution">Grid N³</label>
173+ <input type="range" id="resolution" min="2" max="32" step="1" value="16">
174+ <span class="val" id="resolutionVal">16³</span>
175+ </div>
176+ <div class="legendRow" id="hueLegend">
177+ <canvas id="hueRing" width="64" height="64"></canvas>
178+ <span>Hue = phase of M<sub>xy</sub><br>Saturation = |M<sub>xy</sub>|</span>
179+ </div>
180+ <div class="legendRow" id="mzLegend" style="display:none">
181+ <canvas id="mzBar" width="120" height="12"></canvas>
182+ <span>M<sub>z</sub>: −1 … +1</span>
183+ </div>
184+ </section>
185+
186+ <section>
187+ <h2>RF pulse</h2>
188+ <div class="row">
189+ <label for="flip">Flip angle</label>
190+ <input type="range" id="flip" min="0" max="180" step="1" value="90">
191+ <span class="val" id="flipVal">90°</span>
192+ </div>
193+ <div class="row">
194+ <button id="rfBtn" class="primary">Apply RF pulse (about x′)</button>
195+ </div>
196+ </section>
197+
198+ <section>
199+ <h2>Pulse sequence</h2>
200+ <div class="row">
201+ <button id="seqBtn" class="primary" style="flex:1">Run spoiled gradient echo</button>
202+ </div>
203+ <p class="hint" id="seqStatus" style="display:none"></p>
204+ <p class="hint">16 TRs of 30° RF → stepped Gy phase-encode → Gx readout (gradient echo) → Gy rewinder → Gz spoiler. The signal is acquired only during the readout windows, shaded in the trace. Raise the speed to play it out faster.</p>
205+ </section>
206+
207+ <section>
208+ <h2>Field</h2>
209+ <div class="row">
210+ <label for="df0">Δf₀ offset</label>
211+ <input type="range" id="df0" min="-10" max="10" step="0.1" value="3">
212+ <span class="val" id="df0Val">3.0 Hz</span>
213+ </div>
214+ <div class="row">
215+ <label for="gx">Gx</label>
216+ <input type="range" id="gx" min="-10" max="10" step="0.1" value="0">
217+ <span class="val" id="gxVal">0.0 Hz/cm</span>
218+ </div>
219+ <div class="row">
220+ <label for="gy">Gy</label>
221+ <input type="range" id="gy" min="-10" max="10" step="0.1" value="0">
222+ <span class="val" id="gyVal">0.0 Hz/cm</span>
223+ </div>
224+ <div class="row">
225+ <label for="gz">Gz</label>
226+ <input type="range" id="gz" min="-10" max="10" step="0.1" value="0">
227+ <span class="val" id="gzVal">0.0 Hz/cm</span>
228+ </div>
229+ <div class="row">
230+ <button id="zeroGradBtn">Zero gradients &amp; offset</button>
231+ <button id="soundBtn">🔇 Sound: off</button>
232+ </div>
233+ <p class="hint">The sample spans ±1 cm on each axis; a gradient of G Hz/cm shifts a spin at position x by G·x Hz.</p>
234+ </section>
235+
236+ <section>
237+ <h2>Relaxation</h2>
238+ <div class="row">
239+ <label for="t1">T1</label>
240+ <input type="range" id="t1" min="1" max="4" step="0.01" value="3">
241+ <span class="val" id="t1Val">1000 ms</span>
242+ </div>
243+ <div class="row">
244+ <label for="t2">T2</label>
245+ <input type="range" id="t2" min="1" max="4" step="0.01" value="2.4771">
246+ <span class="val" id="t2Val">300 ms</span>
247+ </div>
248+ <p class="hint">Defaults are soft-tissue-like; try T2 ≈ 2000 ms (water/CSF) to keep M<sub>xy</sub> around longer.</p>
249+ </section>
250+
251+ <p class="hint">Try: apply a 90° pulse and watch the spins precess at Δf₀ while T2 slowly attenuates them and M<sub>z</sub> recovers (T1). Turn up Gx and the spins fan out (dephase); a 180° pulse refocuses them into a spin echo.</p>
252+</div>
253+
254+<div id="axesLegend"><b style="color:#e05555">x′</b> <b style="color:#4caf50">y′</b> <b style="color:#4d9fff">z (B₀)</b></div>
255+
256+<div id="signalPanel">
257+ <canvas id="signalCanvas"></canvas>
258+ <div id="signalLegend">Signal ⟨M<sub>xy</sub>⟩ &nbsp; <b style="color:#e05555">Re (x′)</b> &nbsp; <b style="color:#4caf50">Im (y′)</b></div>
259+</div>
260+
261+<script type="module">
262+import * as THREE from 'three';
263+import { OrbitControls } from 'three/addons/controls/OrbitControls.js';
264+import { mergeGeometries } from 'three/addons/utils/BufferGeometryUtils.js';
265+
266+// ---------- simulation state ----------
267+
268+let N = 16;
269+let count = 0;
270+let pos = null; // Float32Array, xyz per isochromat, coordinates in [-1, 1]
271+let M = null; // Float32Array, (Mx, My, Mz) per isochromat
272+
273+let playing = true;
274+let simTime = 0;
275+let speed = 0.25;
276+let df0 = 3, gx = 0, gy = 0, gz = 0; // df0 in Hz; gradients in Hz/cm (cube coords are cm, spanning ±1)
277+let T1 = 1, T2 = 0.3; // seconds
278+let flipDeg = 90;
279+let sample = 'sphere';
280+let shape = 'arrows';
281+let colorMode = 'transverse';
282+
283+// pulse-sequence playout state
284+let seq = null; // active sequence {name, duration, events, status?}
285+let seqTime = 0;
286+let seqIndex = 0;
287+let acquiring = false;
288+let savedControls = null;
289+
290+function resetM() {
291+ for (let i = 0; i < count; i++) {
292+ M[3 * i] = 0;
293+ M[3 * i + 1] = 0;
294+ M[3 * i + 2] = 1;
295+ }
296+ simTime = 0;
297+ sigSamples.length = 0;
298+ acqWindows.length = 0;
299+}
300+
301+// One Bloch step in the rotating frame: exact z-rotation for the local
302+// off-resonance, plus exponential T1/T2 relaxation. Unconditionally stable.
303+function stepPhysics(dt) {
304+ const e1 = Math.exp(-dt / T1);
305+ const e2 = Math.exp(-dt / T2);
306+ const twoPiDt = 2 * Math.PI * dt;
307+ for (let i = 0; i < count; i++) {
308+ const k = 3 * i;
309+ const f = df0 + gx * pos[k] + gy * pos[k + 1] + gz * pos[k + 2];
310+ const phi = -twoPiDt * f;
311+ const c = Math.cos(phi), s = Math.sin(phi);
312+ const mx = M[k], my = M[k + 1];
313+ M[k] = (mx * c - my * s) * e2;
314+ M[k + 1] = (mx * s + my * c) * e2;
315+ M[k + 2] = 1 + (M[k + 2] - 1) * e1;
316+ }
317+}
318+
319+// Instantaneous RF rotation about the x' axis (z tips toward -y', NMR convention).
320+function applyRF() {
321+ const a = flipDeg * Math.PI / 180;
322+ const c = Math.cos(a), s = Math.sin(a);
323+ for (let i = 0; i < count; i++) {
324+ const k = 3 * i;
325+ const my = M[k + 1], mz = M[k + 2];
326+ M[k + 1] = my * c - mz * s;
327+ M[k + 2] = my * s + mz * c;
328+ }
329+}
330+
331+// ---------- three.js scene ----------
332+
333+const container = document.getElementById('scene');
334+const renderer = new THREE.WebGLRenderer({ antialias: true });
335+renderer.setPixelRatio(Math.min(window.devicePixelRatio, 2));
336+renderer.setSize(window.innerWidth, window.innerHeight);
337+container.appendChild(renderer.domElement);
338+
339+const scene = new THREE.Scene();
340+scene.background = new THREE.Color(0x0d1117);
341+
342+// default view: x right, y up, z (B0) toward the observer
343+const camera = new THREE.PerspectiveCamera(45, window.innerWidth / window.innerHeight, 0.1, 100);
344+camera.up.set(0, 1, 0);
345+camera.position.set(0, 0, 5.2);
346+
347+const controls = new OrbitControls(camera, renderer.domElement);
348+controls.enableDamping = true;
349+controls.dampingFactor = 0.1;
350+
351+scene.add(new THREE.AmbientLight(0xffffff, 0.85));
352+const dirLight = new THREE.DirectionalLight(0xffffff, 0.6);
353+dirLight.position.set(2, -3, 5);
354+scene.add(dirLight);
355+
356+// sample boundary: cube edges, or three great circles for the sphere
357+const boundaryMat = new THREE.LineBasicMaterial({ color: 0x384049 });
358+const boxEdges = new THREE.LineSegments(
359+ new THREE.EdgesGeometry(new THREE.BoxGeometry(2.2, 2.2, 2.2)),
360+ boundaryMat
361+);
362+scene.add(boxEdges);
363+
364+const sphereEdges = new THREE.Group();
365+for (const plane of ['xy', 'xz', 'yz']) {
366+ const pts = [];
367+ for (let a = 0; a <= 64; a++) {
368+ const t = a / 64 * 2 * Math.PI;
369+ const u = 1.1 * Math.cos(t), v = 1.1 * Math.sin(t);
370+ if (plane === 'xy') pts.push(new THREE.Vector3(u, v, 0));
371+ else if (plane === 'xz') pts.push(new THREE.Vector3(u, 0, v));
372+ else pts.push(new THREE.Vector3(0, u, v));
373+ }
374+ sphereEdges.add(new THREE.Line(new THREE.BufferGeometry().setFromPoints(pts), boundaryMat));
375+}
376+scene.add(sphereEdges);
377+boxEdges.visible = sample === 'cube';
378+sphereEdges.visible = sample === 'sphere';
379+
380+// axes with labels
381+function addAxis(dir, color, label) {
382+ const g = new THREE.BufferGeometry().setFromPoints([
383+ new THREE.Vector3(0, 0, 0),
384+ new THREE.Vector3(dir[0], dir[1], dir[2]).multiplyScalar(1.7),
385+ ]);
386+ scene.add(new THREE.Line(g, new THREE.LineBasicMaterial({ color })));
387+ const canvas = document.createElement('canvas');
388+ canvas.width = canvas.height = 64;
389+ const ctx = canvas.getContext('2d');
390+ ctx.font = '42px sans-serif';
391+ ctx.textAlign = 'center';
392+ ctx.textBaseline = 'middle';
393+ ctx.fillStyle = '#' + new THREE.Color(color).getHexString();
394+ ctx.fillText(label, 32, 34);
395+ const sprite = new THREE.Sprite(new THREE.SpriteMaterial({
396+ map: new THREE.CanvasTexture(canvas),
397+ depthTest: false,
398+ }));
399+ sprite.position.set(dir[0], dir[1], dir[2]).multiplyScalar(1.88);
400+ sprite.scale.setScalar(0.24);
401+ scene.add(sprite);
402+}
403+addAxis([1, 0, 0], 0xe05555, "x′");
404+addAxis([0, 1, 0], 0x4caf50, "y′");
405+addAxis([0, 0, 1], 0x4d9fff, 'z');
406+
407+// isochromat meshes (rebuilt when resolution changes)
408+let spheres = null;
409+let arrows = null;
410+
411+const sphereGeoBase = new THREE.SphereGeometry(1, 10, 8);
412+const arrowGeoBase = mergeGeometries([
413+ new THREE.CylinderGeometry(0.025, 0.025, 0.72, 6).translate(0, 0.36, 0),
414+ new THREE.ConeGeometry(0.065, 0.28, 8).translate(0, 0.86, 0),
415+]); // unit-length arrow along +Y, base at origin
416+
417+function buildGrid() {
418+ const coords = [];
419+ for (let i = 0; i < N; i++) {
420+ for (let j = 0; j < N; j++) {
421+ for (let l = 0; l < N; l++) {
422+ const x = -1 + 2 * i / (N - 1);
423+ const y = -1 + 2 * j / (N - 1);
424+ const z = -1 + 2 * l / (N - 1);
425+ if (sample === 'sphere' && x * x + y * y + z * z > 1.000001) continue;
426+ coords.push(x, y, z);
427+ }
428+ }
429+ }
430+ if (coords.length === 0) coords.push(0, 0, 0);
431+ count = coords.length / 3;
432+ pos = Float32Array.from(coords);
433+ M = new Float32Array(count * 3);
434+ resetM();
435+
436+ for (const old of [spheres, arrows]) {
437+ if (old) {
438+ scene.remove(old);
439+ old.material.dispose();
440+ old.dispose();
441+ }
442+ }
443+
444+ const spacing = 2 / (N - 1);
445+ const radius = Math.min(0.08, spacing * 0.11);
446+
447+ spheres = new THREE.InstancedMesh(
448+ sphereGeoBase,
449+ new THREE.MeshBasicMaterial(),
450+ count
451+ );
452+ const mtx = new THREE.Matrix4();
453+ const white = new THREE.Color(0xffffff);
454+ for (let i = 0; i < count; i++) {
455+ mtx.makeScale(radius, radius, radius);
456+ mtx.setPosition(pos[3 * i], pos[3 * i + 1], pos[3 * i + 2]);
457+ spheres.setMatrixAt(i, mtx);
458+ spheres.setColorAt(i, white);
459+ }
460+ scene.add(spheres);
461+
462+ arrows = new THREE.InstancedMesh(
463+ arrowGeoBase,
464+ new THREE.MeshLambertMaterial(),
465+ count
466+ );
467+ for (let i = 0; i < count; i++) arrows.setColorAt(i, white);
468+ scene.add(arrows);
469+
470+ spheres.visible = shape === 'balls';
471+ arrows.visible = shape === 'arrows';
472+}
473+
474+// ---------- per-frame coloring / orientation ----------
475+
476+const tmpColor = new THREE.Color();
477+const tmpQuat = new THREE.Quaternion();
478+const tmpMtx = new THREE.Matrix4();
479+const tmpDir = new THREE.Vector3();
480+const tmpPos = new THREE.Vector3();
481+const tmpScale = new THREE.Vector3();
482+const yAxis = new THREE.Vector3(0, 1, 0);
483+
484+function transverseColor(mx, my, out, lightBase, lightGain) {
485+ const mag = Math.min(1, Math.hypot(mx, my));
486+ let h = Math.atan2(my, mx) / (2 * Math.PI);
487+ if (h < 0) h += 1;
488+ out.setHSL(h, mag, lightBase + lightGain * mag);
489+}
490+
491+// Mz colormap: blue (-1) -> gray (0) -> orange (+1); arrows get a lighter
492+// mid-gray so near-zero-Mz arrows stay visible against the dark background.
493+const negC = new THREE.Color(0x2f6fed);
494+const midBall = new THREE.Color(0x3a3f46);
495+const midArrow = new THREE.Color(0x848d97);
496+const posC = new THREE.Color(0xff9d3c);
497+const noneBall = new THREE.Color(0x7d8590);
498+const noneArrow = new THREE.Color(0xadb6c0);
499+
500+function setDataColor(i, out, forArrows) {
501+ const k = 3 * i;
502+ if (colorMode === 'none') {
503+ out.copy(forArrows ? noneArrow : noneBall);
504+ } else if (colorMode === 'transverse') {
505+ if (forArrows) transverseColor(M[k], M[k + 1], out, 0.55, 0);
506+ else transverseColor(M[k], M[k + 1], out, 0.18, 0.37);
507+ } else {
508+ const mz = Math.max(-1, Math.min(1, M[k + 2]));
509+ const mid = forArrows ? midArrow : midBall;
510+ if (mz >= 0) out.lerpColors(mid, posC, mz);
511+ else out.lerpColors(mid, negC, -mz);
512+ }
513+}
514+
515+function updateVisuals() {
516+ if (shape === 'arrows') {
517+ const spacing = 2 / (N - 1);
518+ const maxLen = Math.min(0.55, spacing * 0.95);
519+ for (let i = 0; i < count; i++) {
520+ const k = 3 * i;
521+ const mag = Math.min(1.2, Math.hypot(M[k], M[k + 1], M[k + 2]));
522+ if (mag < 1e-4) {
523+ tmpScale.setScalar(0);
524+ tmpQuat.identity();
525+ } else {
526+ tmpDir.set(M[k], M[k + 1], M[k + 2]).normalize();
527+ tmpQuat.setFromUnitVectors(yAxis, tmpDir);
528+ tmpScale.setScalar(mag * maxLen);
529+ }
530+ tmpPos.set(pos[k], pos[k + 1], pos[k + 2]);
531+ tmpMtx.compose(tmpPos, tmpQuat, tmpScale);
532+ arrows.setMatrixAt(i, tmpMtx);
533+ setDataColor(i, tmpColor, true);
534+ arrows.setColorAt(i, tmpColor);
535+ }
536+ arrows.instanceMatrix.needsUpdate = true;
537+ arrows.instanceColor.needsUpdate = true;
538+ } else {
539+ for (let i = 0; i < count; i++) {
540+ setDataColor(i, tmpColor, false);
541+ spheres.setColorAt(i, tmpColor);
542+ }
543+ spheres.instanceColor.needsUpdate = true;
544+ }
545+}
546+
547+// ---------- gradient sounds ----------
548+// Every step change of a gradient produces a click (like scanner gradient
549+// noise): different pitch per axis, louder for bigger amplitude jumps.
550+
551+let audioCtx = null;
552+let soundOn = false;
553+let noiseBuf = null;
554+const GRAD_FREQ = { x: 140, y: 240, z: 400 }; // body resonance per axis
555+const lastClickAt = { x: 0, y: 0, z: 0 };
556+
557+function makeNoiseBuffer() {
558+ const len = Math.floor(audioCtx.sampleRate * 0.05);
559+ const buf = audioCtx.createBuffer(1, len, audioCtx.sampleRate);
560+ const d = buf.getChannelData(0);
561+ for (let i = 0; i < len; i++) d[i] = Math.random() * 2 - 1;
562+ return buf;
563+}
564+
565+// A gradient switch sounds like a mechanical knock: a broadband impact
566+// transient (band-passed noise) plus a damped body resonance (sine with a
567+// fast pitch drop). Louder for bigger amplitude jumps.
568+function gradClick(axis, delta) {
569+ if (!soundOn || !audioCtx) return;
570+ const now = audioCtx.currentTime;
571+ if (now - lastClickAt[axis] < 0.05) return; // throttle manual slider drags
572+ lastClickAt[axis] = now;
573+ const f = GRAD_FREQ[axis];
574+ const amp = 0.12 + 0.45 * Math.min(1, Math.abs(delta) / 10);
575+
576+ const osc = audioCtx.createOscillator();
577+ const og = audioCtx.createGain();
578+ osc.type = 'sine';
579+ osc.frequency.setValueAtTime(f * 1.8, now);
580+ osc.frequency.exponentialRampToValueAtTime(f, now + 0.012);
581+ og.gain.setValueAtTime(amp, now);
582+ og.gain.exponentialRampToValueAtTime(0.001, now + 0.09);
583+ osc.connect(og).connect(audioCtx.destination);
584+ osc.start(now);
585+ osc.stop(now + 0.1);
586+
587+ const noise = audioCtx.createBufferSource();
588+ noise.buffer = noiseBuf;
589+ const bp = audioCtx.createBiquadFilter();
590+ bp.type = 'bandpass';
591+ bp.frequency.value = f * 5;
592+ bp.Q.value = 0.9;
593+ const ng = audioCtx.createGain();
594+ ng.gain.setValueAtTime(amp * 0.9, now);
595+ ng.gain.exponentialRampToValueAtTime(0.001, now + 0.02);
596+ noise.connect(bp).connect(ng).connect(audioCtx.destination);
597+ noise.start(now);
598+ noise.stop(now + 0.03);
599+}
600+
601+// ---------- UI wiring ----------
602+
603+function bindSlider(id, fmt, cb) {
604+ const input = document.getElementById(id);
605+ const val = document.getElementById(id + 'Val');
606+ const update = () => {
607+ const v = parseFloat(input.value);
608+ val.textContent = fmt(v);
609+ cb(v);
610+ };
611+ input.addEventListener('input', update);
612+ update();
613+}
614+
615+const fmtHz = v => v.toFixed(1) + ' Hz';
616+const fmtHzCm = v => v.toFixed(1) + ' Hz/cm';
617+const fmtT = v => Number(Math.pow(10, v).toPrecision(3)) + ' ms'; // slider value is log10(ms)
618+
619+// speed slider is log10(multiplier), 0.1x .. 4x
620+bindSlider('speed', v => Number(Math.pow(10, v).toPrecision(2)) + '×', v => { speed = Math.pow(10, v); });
621+bindSlider('resolution', v => v + '³', v => {
622+ if (v !== N) { N = v; buildGrid(); }
623+});
624+bindSlider('flip', v => v + '°', v => { flipDeg = v; });
625+bindSlider('df0', fmtHz, v => { df0 = v; });
626+bindSlider('gx', fmtHzCm, v => { if (v !== gx) gradClick('x', v - gx); gx = v; });
627+bindSlider('gy', fmtHzCm, v => { if (v !== gy) gradClick('y', v - gy); gy = v; });
628+bindSlider('gz', fmtHzCm, v => { if (v !== gz) gradClick('z', v - gz); gz = v; });
629+bindSlider('t1', fmtT, v => { T1 = Math.pow(10, v) / 1000; });
630+bindSlider('t2', fmtT, v => { T2 = Math.pow(10, v) / 1000; });
631+
632+// set a slider programmatically, updating its readout and bound variable
633+function setControl(id, v) {
634+ const input = document.getElementById(id);
635+ input.value = v;
636+ input.dispatchEvent(new Event('input'));
637+}
638+
639+const playBtn = document.getElementById('playBtn');
640+playBtn.addEventListener('click', () => {
641+ playing = !playing;
642+ playBtn.textContent = playing ? 'Pause' : 'Play';
643+});
644+document.getElementById('resetBtn').addEventListener('click', () => {
645+ if (seq) stopSequence();
646+ resetM();
647+});
648+const rfBtn = document.getElementById('rfBtn');
649+rfBtn.addEventListener('click', applyRF);
650+document.getElementById('zeroGradBtn').addEventListener('click', () => {
651+ for (const id of ['df0', 'gx', 'gy', 'gz']) setControl(id, 0);
652+});
653+const soundBtn = document.getElementById('soundBtn');
654+soundBtn.addEventListener('click', () => {
655+ soundOn = !soundOn;
656+ if (soundOn && !audioCtx) audioCtx = new (window.AudioContext || window.webkitAudioContext)();
657+ if (soundOn && !noiseBuf) noiseBuf = makeNoiseBuffer();
658+ if (soundOn && audioCtx.state === 'suspended') audioCtx.resume();
659+ soundBtn.textContent = soundOn ? '🔊 Sound: on' : '🔇 Sound: off';
660+});
661+
662+function setSample(s) {
663+ if (s === sample) return;
664+ sample = s;
665+ document.getElementById('sampleCube').classList.toggle('active', s === 'cube');
666+ document.getElementById('sampleSphere').classList.toggle('active', s === 'sphere');
667+ boxEdges.visible = s === 'cube';
668+ sphereEdges.visible = s === 'sphere';
669+ buildGrid();
670+}
671+document.getElementById('sampleCube').addEventListener('click', () => setSample('cube'));
672+document.getElementById('sampleSphere').addEventListener('click', () => setSample('sphere'));
673+
674+function setShape(s) {
675+ shape = s;
676+ document.getElementById('shapeBalls').classList.toggle('active', s === 'balls');
677+ document.getElementById('shapeArrows').classList.toggle('active', s === 'arrows');
678+ spheres.visible = s === 'balls';
679+ arrows.visible = s === 'arrows';
680+}
681+function setColorMode(c) {
682+ colorMode = c;
683+ document.getElementById('colorTransverse').classList.toggle('active', c === 'transverse');
684+ document.getElementById('colorLongitudinal').classList.toggle('active', c === 'longitudinal');
685+ document.getElementById('colorNone').classList.toggle('active', c === 'none');
686+ document.getElementById('hueLegend').style.display = c === 'transverse' ? 'flex' : 'none';
687+ document.getElementById('mzLegend').style.display = c === 'longitudinal' ? 'flex' : 'none';
688+}
689+document.getElementById('shapeBalls').addEventListener('click', () => setShape('balls'));
690+document.getElementById('shapeArrows').addEventListener('click', () => setShape('arrows'));
691+document.getElementById('colorTransverse').addEventListener('click', () => setColorMode('transverse'));
692+document.getElementById('colorLongitudinal').addEventListener('click', () => setColorMode('longitudinal'));
693+document.getElementById('colorNone').addEventListener('click', () => setColorMode('none'));
694+
695+// legends
696+{
697+ const ctx = document.getElementById('hueRing').getContext('2d');
698+ for (let deg = 0; deg < 360; deg += 3) {
699+ ctx.beginPath();
700+ ctx.strokeStyle = `hsl(${deg}, 100%, 50%)`;
701+ ctx.lineWidth = 9;
702+ ctx.arc(32, 32, 24, -deg * Math.PI / 180, -(deg + 4) * Math.PI / 180, true);
703+ ctx.stroke();
704+ }
705+ ctx.fillStyle = '#8b949e';
706+ ctx.font = '13px sans-serif';
707+ ctx.textAlign = 'center';
708+ ctx.textBaseline = 'middle';
709+ ctx.fillText('φ', 32, 33);
710+}
711+{
712+ const ctx = document.getElementById('mzBar').getContext('2d');
713+ const grad = ctx.createLinearGradient(0, 0, 120, 0);
714+ grad.addColorStop(0, '#2f6fed');
715+ grad.addColorStop(0.5, '#3a3f46');
716+ grad.addColorStop(1, '#ff9d3c');
717+ ctx.fillStyle = grad;
718+ ctx.fillRect(0, 0, 120, 12);
719+}
720+
721+// ---------- signal trace ----------
722+
723+const sigCanvas = document.getElementById('signalCanvas');
724+const sigCtx = sigCanvas.getContext('2d');
725+const SIG_WINDOW = 4; // seconds of sim time shown
726+const sigSamples = []; // {t, re, im, gap}, re/im are means over all isochromats
727+const acqWindows = []; // {t0, t1} acquisition spans, shaded in the trace
728+let sigPendingGap = false; // next sample starts a new stroke (acq gap)
729+let sigDpr = 1;
730+
731+function resizeSignalCanvas() {
732+ sigDpr = Math.min(window.devicePixelRatio, 2);
733+ sigCanvas.width = sigCanvas.clientWidth * sigDpr;
734+ sigCanvas.height = sigCanvas.clientHeight * sigDpr;
735+}
736+resizeSignalCanvas();
737+
738+function sampleSignal() {
739+ let re = 0, im = 0;
740+ for (let i = 0; i < count; i++) {
741+ re += M[3 * i];
742+ im += M[3 * i + 1];
743+ }
744+ sigSamples.push({ t: simTime, re: re / count, im: im / count, gap: sigPendingGap });
745+ sigPendingGap = false;
746+ if (acquiring && acqWindows.length) acqWindows[acqWindows.length - 1].t1 = simTime;
747+ const cutoff = simTime - SIG_WINDOW;
748+ while (sigSamples.length && sigSamples[0].t < cutoff) sigSamples.shift();
749+}
750+
751+function drawSignal() {
752+ const w = sigCanvas.width, h = sigCanvas.height;
753+ if (!w) return;
754+ sigCtx.clearRect(0, 0, w, h);
755+ const xOf = t => w - (simTime - t) / SIG_WINDOW * w;
756+ const yOf = v => h / 2 - v * h * 0.42;
757+
758+ // shade acquisition windows
759+ while (acqWindows.length && acqWindows[0].t1 < simTime - SIG_WINDOW) acqWindows.shift();
760+ sigCtx.fillStyle = 'rgba(88, 166, 255, 0.10)';
761+ for (const wdw of acqWindows) {
762+ const x0 = Math.max(0, xOf(wdw.t0));
763+ const x1 = Math.min(w, xOf(wdw.t1));
764+ if (x1 > x0) sigCtx.fillRect(x0, 0, x1 - x0, h);
765+ }
766+
767+ // gridlines: one per second of sim time, plus 0 and ±1 levels
768+ sigCtx.lineWidth = 1;
769+ sigCtx.strokeStyle = 'rgba(48, 54, 61, 0.7)';
770+ for (let ts = Math.ceil(simTime - SIG_WINDOW); ts <= simTime; ts++) {
771+ sigCtx.beginPath();
772+ sigCtx.moveTo(xOf(ts), 0);
773+ sigCtx.lineTo(xOf(ts), h);
774+ sigCtx.stroke();
775+ }
776+ sigCtx.setLineDash([4 * sigDpr, 4 * sigDpr]);
777+ for (const v of [-1, 1]) {
778+ sigCtx.beginPath();
779+ sigCtx.moveTo(0, yOf(v));
780+ sigCtx.lineTo(w, yOf(v));
781+ sigCtx.stroke();
782+ }
783+ sigCtx.setLineDash([]);
784+ sigCtx.strokeStyle = '#414954';
785+ sigCtx.beginPath();
786+ sigCtx.moveTo(0, yOf(0));
787+ sigCtx.lineTo(w, yOf(0));
788+ sigCtx.stroke();
789+
790+ sigCtx.lineWidth = 1.5 * sigDpr;
791+ for (const [key, color] of [['re', '#e05555'], ['im', '#4caf50']]) {
792+ sigCtx.strokeStyle = color;
793+ sigCtx.beginPath();
794+ for (let i = 0; i < sigSamples.length; i++) {
795+ const s = sigSamples[i];
796+ if (i === 0 || s.gap) sigCtx.moveTo(xOf(s.t), yOf(s[key]));
797+ else sigCtx.lineTo(xOf(s.t), yOf(s[key]));
798+ }
799+ sigCtx.stroke();
800+ }
801+}
802+
803+// ---------- pulse sequences ----------
804+
805+const seqBtn = document.getElementById('seqBtn');
806+const seqStatusEl = document.getElementById('seqStatus');
807+const SEQ_RUN_LABEL = 'Run spoiled gradient echo';
808+
809+// Spoiled gradient echo: each TR is RF -> Gx dephase lobe + stepped Gy
810+// phase-encode -> Gx readout (echo at TR center of readout) -> Gy rewinder
811+// -> Gz spoiler.
812+function makeSpoiledGradientEcho() {
813+ const TR = 0.9, nTR = 16, alpha = 30;
814+ const events = [];
815+ for (let r = 0; r < nTR; r++) {
816+ const t0 = r * TR;
817+ const pe = -7.5 + 15 * r / (nTR - 1); // Gy phase-encode step, -7.5..7.5 Hz/cm
818+ events.push(
819+ { t: t0, type: 'rf', flip: alpha },
820+ { t: t0 + 0.02, type: 'grad', gx: -5, gy: pe }, // dephase + PE lobes, 0.15 s
821+ { t: t0 + 0.17, type: 'grad', gx: 5, gy: 0 }, // readout, 0.3 s; echo at center
822+ { t: t0 + 0.17, type: 'acqOn' },
823+ { t: t0 + 0.47, type: 'acqOff' },
824+ { t: t0 + 0.47, type: 'grad', gx: 0 },
825+ { t: t0 + 0.49, type: 'grad', gy: -pe }, // PE rewinder, 0.15 s
826+ { t: t0 + 0.64, type: 'grad', gy: 0 },
827+ { t: t0 + 0.66, type: 'grad', gz: 10 }, // spoiler, 0.2 s
828+ { t: t0 + 0.86, type: 'grad', gz: 0 },
829+ );
830+ }
831+ events.sort((a, b) => a.t - b.t);
832+ return {
833+ name: 'Spoiled gradient echo',
834+ duration: nTR * TR,
835+ events,
836+ status: t => `TR ${Math.min(nTR, Math.floor(t / TR) + 1)} / ${nTR}`,
837+ };
838+}
839+
840+// everything the sequence drives or that would rebuild the sample mid-run
841+const SEQ_LOCKED = ['flip', 'gx', 'gy', 'gz', 'df0', 't1', 't2', 'resolution',
842+ 'rfBtn', 'zeroGradBtn', 'sampleCube', 'sampleSphere'];
843+
844+function setSeqUIDisabled(dis) {
845+ for (const id of SEQ_LOCKED) document.getElementById(id).disabled = dis;
846+}
847+
848+function flashRF() {
849+ rfBtn.classList.add('flash');
850+ setTimeout(() => rfBtn.classList.remove('flash'), 250);
851+}
852+
853+function applyEvent(ev) {
854+ if (ev.type === 'rf') {
855+ setControl('flip', ev.flip);
856+ applyRF();
857+ flashRF();
858+ } else if (ev.type === 'grad') {
859+ for (const k of ['gx', 'gy', 'gz']) if (k in ev) setControl(k, ev[k]);
860+ } else if (ev.type === 'acqOn') {
861+ acquiring = true;
862+ acqWindows.push({ t0: simTime, t1: simTime });
863+ } else if (ev.type === 'acqOff') {
864+ acquiring = false;
865+ sigPendingGap = true;
866+ }
867+}
868+
869+function startSequence(makeFn) {
870+ seq = makeFn();
871+ seqTime = 0;
872+ seqIndex = 0;
873+ acquiring = false;
874+ savedControls = { flip: flipDeg, gx, gy, gz };
875+ resetM();
876+ for (const id of ['gx', 'gy', 'gz']) setControl(id, 0);
877+ sigPendingGap = true;
878+ setSeqUIDisabled(true);
879+ seqBtn.textContent = 'Stop sequence';
880+ seqStatusEl.style.display = 'block';
881+ playing = true;
882+ playBtn.textContent = 'Pause';
883+}
884+
885+function stopSequence() {
886+ seq = null;
887+ acquiring = false;
888+ sigPendingGap = true;
889+ setControl('flip', savedControls.flip);
890+ setControl('gx', savedControls.gx);
891+ setControl('gy', savedControls.gy);
892+ setControl('gz', savedControls.gz);
893+ setSeqUIDisabled(false);
894+ seqBtn.textContent = SEQ_RUN_LABEL;
895+ seqStatusEl.style.display = 'none';
896+}
897+
898+seqBtn.addEventListener('click', () => {
899+ if (seq) stopSequence();
900+ else startSequence(makeSpoiledGradientEcho);
901+});
902+
903+// Advance sim time by dt. During a sequence, step exactly to each event
904+// boundary so event timing does not depend on frame rate.
905+function advanceSim(dt) {
906+ if (!seq) {
907+ stepPhysics(dt);
908+ simTime += dt;
909+ return;
910+ }
911+ let remaining = dt;
912+ while (remaining > 1e-9 && seq) {
913+ const ev = seq.events[seqIndex];
914+ const tStop = ev ? ev.t : seq.duration;
915+ const step = Math.min(remaining, Math.max(0, tStop - seqTime));
916+ if (step > 0) {
917+ stepPhysics(step);
918+ seqTime += step;
919+ simTime += step;
920+ remaining -= step;
921+ }
922+ if (ev) {
923+ if (seqTime >= ev.t - 1e-9) {
924+ applyEvent(ev);
925+ seqIndex++;
926+ }
927+ } else if (seqTime >= seq.duration - 1e-9) {
928+ stopSequence(); // sequence finished; free evolution resumes
929+ if (remaining > 1e-9) {
930+ stepPhysics(remaining);
931+ simTime += remaining;
932+ }
933+ return;
934+ }
935+ }
936+}
937+
938+// ---------- main loop ----------
939+
940+const timeReadout = document.getElementById('timeReadout');
941+const clock = new THREE.Clock();
942+
943+buildGrid();
944+
945+function animate() {
946+ requestAnimationFrame(animate);
947+ const rawDt = Math.min(0.05, clock.getDelta());
948+ if (playing) {
949+ const dt = rawDt * speed;
950+ advanceSim(dt);
951+ if (!seq || acquiring) sampleSignal();
952+ }
953+ if (seq) {
954+ seqStatusEl.textContent =
955+ (seq.status ? seq.status(seqTime) + ' · ' : '') + 'sequence t = ' + seqTime.toFixed(2) + ' s';
956+ }
957+ timeReadout.textContent = 't = ' + simTime.toFixed(1) + ' s';
958+ updateVisuals();
959+ drawSignal();
960+ controls.update();
961+ renderer.render(scene, camera);
962+}
963+animate();
964+
965+window.addEventListener('resize', () => {
966+ camera.aspect = window.innerWidth / window.innerHeight;
967+ camera.updateProjectionMatrix();
968+ renderer.setSize(window.innerWidth, window.innerHeight);
969+ resizeSignalCanvas();
970+});
971+</script>
972+</body>
973+</html>