qMRI mOOC · Chapter 5

B₀ Mapping

Dual-echo mapping · Phase unwrapping · Multi-echo methods

Slides from the qMRLab mOOC chapter by Alexandre d'Astous
NeuroPoly Lab, Polytechnique Montréal

📖 B0 Mapping

Why map B₀?

Image reconstruction assumes a perfectly uniform B₀. It never is.

A non-uniform B₀ causes

Signal loss (T₂^*)
Geometric distortion (EPI)
Failed fat saturation
Broad MRS lines

A B₀ map lets you

Shim the field
Unwarp EPI (FUGUE)
Mitigate T₂^* signal loss
Compute susceptibility, χ

Measured at 3 T, after shimming: one axial slice, dual-echo GRE.

B0 sets the precession frequency and the bulk magnetization; the RF coils are tuned to it. Reconstruction assumes it is perfectly homogeneous. Inhomogeneity leads to signal loss, distortions (Jezzard 1995), poor fat saturation (Anzai 1992) and broader MRS linewidths; in extreme cases no image at all. B0 maps are used for active shimming, retrospective distortion correction (FSL FUGUE, Jenkinson 2012), recovery of T2* signal decay (An 2002) and QSM. The map: the 2013 Stikov 3 T session (the B1 deck's dataset), the qMT protocol's B0 map (mtspgr_b0, Boudreau, Stikov & Pike, MRM 2018): a two-point phase-difference GRE, TR 25 ms, 7°, one 5 mm slice, 2 mm pixels. Complex phase difference divided by 2π·ΔTE with ΔTE = 2.24 ms (the data and the original pipeline agree on 2.24 ms; the paper reports 4/8.48 ms). Even after shimming the field varies across the slice by roughly −30 to +10 Hz, lowest frontally toward the sinuses. The scalp is masked: fat's chemical shift is not field.

What a B₀ map shows

Per voxel: actual − expected field
Units: Hz, µT or ppm

$$\delta_{ppm} = \frac{\Delta f}{f_0}\,10^6 = \frac{\Delta B}{B_0}\,10^6$$

B₀	1 ppm =
64 mT	2.7 Hz	0.064 µT
1.5 T	64 Hz	1.5 µT
3 T	128 Hz	3 µT
7 T	298 Hz	7 µT

A 3 T brain, field map in ppm. Menu: magnitude, phase, or the map in Hz or µT.

Where ΔB₀ comes from

The scanner

Magnet imperfections: < 1 ppm, fixed by passive shimming
Slow drift from heating

The subject

Tissue magnetic susceptibility χ

$$\Delta B_z = B_0\,(\Delta\chi \otimes d)$$

The dipole kernel d (image space) and D (k-space, as used in QSM).

From χ to ΔB₀, and back

If χ is known, ΔB₀ is its convolution with the dipole kernel d. B₀ maps simulated this way are useful for designing shim coils and evaluating shimming algorithms.

$$\textcolor{#CC79A7}{\Delta B_0(\mathbf{r})}\,[\text{ppm}] = \textcolor{#009E73}{\int} \textcolor{#0072B2}{\Delta\chi(\mathbf{r}')}\, \textcolor{#D55E00}{d(\mathbf{r}-\mathbf{r}')}\, \textcolor{#009E73}{\mathrm{d}^3\mathbf{r}'}$$

Speed 1×

$\textcolor{#0072B2}{\Delta\chi(\mathbf{r}')}$

$\textcolor{#D55E00}{d(\mathbf{r})}$

B₀ ↑ · ◯ the r′ summed

$\textcolor{#D55E00}{d(\mathbf{r}-\mathbf{r}')}$ on $\textcolor{#0072B2}{\Delta\chi(\mathbf{r}')}$

→

$\textcolor{#0072B2}{\Delta\chi(\mathbf{r}')}\,\textcolor{#D55E00}{d(\mathbf{r}-\mathbf{r}')}$

$\textcolor{#009E73}{\int\!\cdot\;\mathrm{d}^3\mathbf{r}'}$

$\textcolor{#009E73}{\displaystyle\sum_{\mathbf{r}'}}$→

$\textcolor{#CC79A7}{\Delta B_0(\mathbf{r})}$

r: the voxel where ΔB₀ is computed, which moves across the map. r′: every voxel of Δχ, which the sum runs over.

In practice the convolution is computed in k-space, as a multiplication: $\Delta B_0 = \mathcal{F}^{-1}\{D \cdot \mathcal{F}\{\Delta\chi\}\}$ (Marques & Bowtell, Concepts Magn Reson B 2005;25B:65–78).

Quantitative Susceptibility Mapping (QSM) asks the inverse question: from ΔB₀ (the phase), compute χ, which tells us about the tissue.

We come back to it on slides 19–20 (unwrapping and multi-echo mapping for QSM), and in the T2 lecture, slide 26.

Forward problem: knowing the susceptibility distribution, the field is its convolution with the dipole kernel d = (3cos²θ − 1)/(4π r³), θ from B0. The demo does it voxel by voxel: the kernel is centred on one voxel, multiplied by Δχ everywhere (middle), and summed, in 3D, to give the field at that voxel (right); then the next voxel. Here the kernel is cut off at 10 voxels, so the field is an approximation, but what is shown is exactly what is summed. In practice the convolution is done in k-space as a multiplication by D = 1/3 − kz²/k² (the convolution theorem; Marques & Bowtell 2005), as in Shimming Toolbox's susceptibility-to-fieldmap-fft. The yellow voxel is r; the kernel is cut off at 10 voxels, and the sum is over the 3D volume. One voxel: the field is the kernel. A uniform sphere: uniform field inside, a dipole outside. An interface: the field changes most where B0 crosses it (⊥), less when it runs along it (‖). QSM is the inverse problem: from the measured field (phase), recover χ; ill-posed where D = 0 (the magic angle cone), hence regularization. χ tells us about iron, myelin, calcium.

Air–tissue interfaces dominate

Air (+χ) next to tissue (−χ): largest offsets
ΔB₀ scales with B₀ (here 7 T)
Every subject needs its own shim

Menu: cylinders or brain. χ (left), field (middle), field with the slow background removed (right).

* This is the ΔB₀ map used in QSM to compute χ.

ΔB₀ effect 1:
dephasing

Different frequencies in a voxel → phases spread → faster decay
Demodulating at the wrong frequency → oscillations

$$\frac{1}{T_2^{*}} = \frac{1}{T_2} + \frac{1}{T_2'}$$

Menu: one, two or many frequencies (isochromats).

t = 0.000 s Speed 1×

ΔB₀ effect 2: wrong frequency, wrong place

Frequency encodes position: off-resonance shifts voxels
Errors build up along the readout: EPI is worst
Also: slice shifts, failed fat sat

Red: a 0.1 mT/m parasitic gradient (0.25 % of the 40 mT/m encoding gradients).

How to measure ΔB₀: phase changes linearly with time

$$\phi(\mathbf r,t) = \phi_0(\mathbf r) + \gamma\,\Delta B_0(\mathbf r)\,t$$

Slope = γΔB₀: the field map
φ₀ (coil, RF) is unknown → at least two echoes

Play; the menu switches between laboratory and rotating frame. In the rotating frame the on-resonance spin (blue) stands still.

B₀ mapping technique: dual echo method

The phase of an MRI pixel at TE is:

$$\phi(TE) = \phi_0 + \gamma\,\Delta B_0\,TE$$

where φ₀ is the initial phase after the excitation pulse, and is independent of TE.

$$\Delta B_0 = \frac{\phi(TE_2) - \phi(TE_1)}{\gamma\,\Delta TE}$$

Complex difference keeps Δφ in [−π, π): $\angle\big(e^{i\phi_2}\,e^{-i\phi_1}\big)$

A cylinder at one constant frequency offset.

From the scanner: DICOM typically stores phase as unsigned integers (e.g. 0 to 4095) with a rescale slope and intercept (e.g. 2 and −4096). Depending on how you read it, you get one of two ranges, both spanning [−π, π):

0 to 4095: $\phi = \dfrac{2\pi}{4096}\,\text{value} - \pi$ −4096 to 4094: $\phi = \dfrac{\pi}{4096}\,\text{value}$

Data ranges differ by vendor, converter and reader, so check yours. It caught us too: qMRLab/mooc-data#1 ↗

⏸ Exercise 5-1

Phase at longer echo times

Measured phase lives in [−π, π)
Sharp jumps are phase wraps
Longer TE → more wraps

Phase images come from real and imaginary data: ∠ = atan2(imag, real)

B₀ mapping technique: dual echo in the brain

Each echo has wraps
The difference often removes them…
…not always: large offsets (sinuses) leave residual wraps

qMRLab mOOC dataset: qMRLab/mooc-data ↗

QSM consensus dataset (Chan et al., CC BY 4.0): Zenodo ↗

Choosing the echo times

Phase noise: $\sigma_{\phi} = 1/\mathrm{SNR}_{mag}$ → short TE₁
Longer ΔTE: more precise field
Shorter ΔTE: fewer wraps
Fat and water in phase at 3 T: ΔTE ≈ 2.34 ms
GRE, monopolar readouts, high bandwidth (not EPI)

No wraps in Δφ if

$$\Delta TE < \frac{1}{2\,\Delta f_{max}}$$

B₀

Max ΔB₀1.0 ppm

Δf_max128 Hz

ΔTE_max3.91 ms

Benefits and pitfalls. At least two phase images, since φ0 is unknown. GRE, or spin echo with an asymmetric second echo shifted 1-2 ms. High bandwidth, thin slices and multi-shot reduce B0 artifacts in the field map itself, so EPI is generally avoided (Akcakaya 2022). Monopolar readouts avoid misregistration between odd and even echoes. σphase = 1/SNRmag (Eq. 5.7, Brown 2014), so TE1 is short. ΔTE: fat has a ~3.35 ppm shift, which is mistaken for field; choose ΔTE so fat and water are in phase (~2.34 ms at 3 T). Longer ΔTE: better estimate if SNR allows; shorter: fewer wraps. Table 5.1 (calculator): with offsets up to ±Δf_max, the phase difference 2πΔf·ΔTE stays within ±π if ΔTE < 1/(2Δf_max). At 1 ppm this reproduces the book's table (3.91 ms at 3 T, 7.82 ms at 1.5 T). Eddy currents from fast switching between echoes: one echo per RF pulse avoids them, at twice the scan time. Note: this calculator also answers the chapter's first exercise.

Phase wrapping

We measure phase modulo 2π
True phase = measured + 2πn: find n

Two families

Temporal: across echoes
Spatial: across neighbouring voxels

Wrapped phase (blue) and true phase (red) along one line.

Temporal unwrapping

Assume the phase changes by < π between echoes
Add or remove 2π where it jumps
Fit a line: slope = γΔB₀

One voxel at four echo times: measured (blue), unwrapped (red), linear fit (green).

Spatial unwrapping

Neighbours vary smoothly: region growing from a seed
The seed sets the answer: solutions differ by 2nπ

Drag the divider: wrapped phase (left), the three possible unwrapped solutions (right).

Unwrapping in 2D

Needs one connected region, and noise masked out
Pick n so the mean phase in the ROI is ≈ 0 (after shimming)
Δφ maps unwrap more reliably than single echoes

Simulated wrapped (left) and unwrapped (right) phase. More wraps = more inhomogeneity.

Phase singularities

Two paths from A to B cross a different number of wraps
No unwrapper can give one right answer
Cause: poor coil combination

Fix the coil combination, not the unwrapper

⏸ Exercise 5-2

Phase unwrapping software

Tool	Approach	Notes
PRELUDE	Spatial, region growing	Robust to noise; slow
SEGUE	Spatial, like PRELUDE	Much faster
BEST PATH	Spatial 3D, quality-guided	Unwraps high-quality regions first
ROMEO	Spatial + temporal	Fast; beats PRELUDE and BEST PATH with noise
UMPIRE	Temporal, ≥ 3 uneven echoes	Handles ΔTE beyond π; noise-sensitive
Laplacian	Spatial, global	Robust but inaccurate: fine for QSM, not for shimming or qMRI

Multi-echo data: unwrap by fitting the slope

Unwrap each echo spatially, then across echoes
Fit a line through all echoes
A 2π offset moves the line, not its slope

More echoes → more accurate maps (QSM, shimming)

Two unwrapped solutions of three echoes, exactly 2π apart.

Multi-echo data: choosing echo times

TE₁, TE₂ close: no wraps, a first field estimate
It predicts the wraps of a late TE₃
Refit with all three: more precise
Late echoes have less SNR: weight by magnitude

⏸ Exercise 5-3

Shimming and real-time mapping

Shimming

A B₀ map tells us the field to cancel
Currents in coils make known field patterns: 0th order (f₀), 1st order (X, Y, Z), 2nd order (Z², ZX, ZY, X²−Y², XY)
Choose the currents that make the field most uniform, over a volume or a slice (i.e. minimizing the spatial SD of the field)

Real-time (breathing)

A map in much less than a breath
2D, slice by slice: each slice has its own time
Dual-echo GRE, minimum TR and TEs (EPI is fast, but distorted)

shimming-toolbox.org ↗

Shimming: the B0 map is the field to cancel. The scanner's shim coils (and gradients, for 1st order) each make a known field pattern, the spherical harmonics: 0th order is the centre frequency f0, 1st order the linear X, Y, Z gradients, 2nd order Z², ZX, ZY, X²−Y² and XY. The currents are chosen so that the sum of these patterns cancels the measured field as well as it can over the region of interest, by least squares: over a whole volume (volume-wise), or slice by slice with the currents changed for each slice (slice-wise, or dynamic, shimming). Exercise 5-3 does this by hand, then compares with the optimal currents. Real-time mapping of respiration: acquire a map in much less than a respiratory cycle; 2D so slice timing can be associated with each slice's field; EPI is faster but distorted, so a 2D dual-echo GRE with minimum TR, short TEs and short RF pulses.

Key takeaways

ΔB₀ grows with B₀ and is dominated by air–tissue interfaces
Phase grows linearly with TE: φ = φ₀ + γΔB₀t, so two echoes at least
ΔTE trades precision against wraps: no wraps if ΔTE < 1/(2Δf_max)
Unwrapping (temporal or spatial) leaves a 2nπ ambiguity; singularities cannot be unwrapped
More echoes: fit the slope; a close echo pair helps unwrap a late echo