In magnetic resonance imaging, the ideal setup features imaging gradients that create a linear spatial field variation while the main B₀ field remains static. However, in practice, the gradients deviate from linearity, and B₀ drifts and other sources of field variation affect the data acquisition. The scanner’s gradient system often fails to produce the gradient waveforms exactly as programmed in the pulse sequence because of hardware imperfections and inherent physical effects, such as limited coil and amplifier bandwidth, eddy currents and concomitant fields induced by the gradients, mutual coupling of the gradient channels, mechanical vibrations, and thermal variations.

Deviations from the expected gradient fields significantly alter the nominal k-space sampling trajectory, leading to violations of the encoding model used for reconstruction. If left uncorrected, these errors result in various image artifacts such as ghosting, blurring, shifts, geometric distortions, and signal drop-offs. B₀ drifts and other time-varying field perturbations complicate MRI techniques that depend on temporal stability for phase-sensitive measurements (e.g., Dixon method, phase-contrast MRI, QSM), consistency across multiple time series (e.g., functional MRI), or for preventing errant phase accrual over multiple TRs (e.g., balanced SSFP sequence).

Fourier Transform MRI, which constitutes 99% of current MRI technology, assumes data is sampled on a perfectly defined grid. However, k-space samples often deviate from the expected trajectory, leading to image artifacts if not addressed in the image reconstruction model.

Examples of nominal EPI and spiral k-space trajectories, along with actual trajectories measured concurrently or predicted from the gradient impulse response function captured with a Field Camera.

Several methods exist for measuring and correcting for field deviations. Most of these require separate calibration measurements, modifications to the pulse sequence, or specific modeling assumptions. Consequently, the corrections often fail to accurately represent the actual field evolutions during the imaging experiment.

A hardware-based solution that requires fewer modeling assumptions is to directly measure the magnetic field with a dynamic field camera using NMR probes. Assuming the field components can be accurately expanded in terms of a series of spatial harmonic functions, the temporal and spatial field variations can be calculated from the signals measured by the NMR probe array. These measurements are used for determining the k-space phase evolutions of linear, second, and up to third order in space. These k-space phase coefficients can then be integrated within an expanded signal model to enable more accurate image reconstruction.

Field Probes

The image below shows the basic design of an NMR field probe, and an illustration of NMR probes distributed over the 20mm-diameter spherical head of a Dynamic Field Camera (DFC).

In the DFC, 16 ¹H probes are used, while for a Clip-on Camera, up to 16 ¹⁹F probes can be mounted around the coil frame. Each probe encapsulates an NMR active liquid droplet approximately 0.8 mm in diameter, enabling k-space excursions up to the equivalent of a 0.8 mm image resolution. The lifetime of the probe signal is approximately 70-100 ms and depends on the composition of the NMR active material and the field strength.

¹⁹F-based probes of a Clip-on Camera mounted on a head coil.

NMR probes serve as effective sensors, directly measuring phase accrual induced by the magnetic field with high temporal resolution, limited only by the spectrometer bandwidth. However, their primary limitation is the short signal lifetime, which restricts monitoring for TRs shorter than 110 ms. For TRs less than 110 ms, common in fast GRE sequences, interpolation is necessary for the TRs that occur between successive probe measurements.

During the imaging experiment, each NMR probe in the array measures the phase accrual due to the magnetic field at its specific location. The precise spatial coordinates of the probes are determined from a separate calibration sequence, as these positions are necessary to solve for the k-space phase coefficients corresponding to each spatial basis function.

Dynamic Field Model

The magnetic field strength, |B(r,t)| is modeled as the sum of a dynamic and a static component. The dynamic field component is represented as a linear combination of N_L spatial basis functions, h_l(r), each of which is scaled with a corresponding set of temporal coefficients c_l(t). Since the NMR probes can only measure field changes, the static field component is added as a separate term. The static field, B_ref(r), can be measured at a reference calibration time using a separate field mapping sequence, such as a multi-echo GRE. This field includes the main B₀ field and its spatial inhomogeneities, including those caused by variations in magnetic susceptibility within the imaged object.

Thus, the magnetic field strength as a function of space and time can be written as

The phase accrual, which is directly measured by the NMR probes, is proportional to the time integral of the magnetic field strength and is given by

Given an appropriate set of spatial basis functions h_l(r) and the probe signals, the dynamic phase coefficients, k_l(t), can be determined.

Expansion of the Field in Terms of Spherical Harmonic Functions

Any set of N_L spatial basis functions can be chosen to model the spatial distribution of the magnetic field. However, one of the most effective is the set of real spherical harmonics, which is the standard choice of basis set for shimming applications on modern scanners. For example, with N_L = 16, the field can be represented by spherical harmonics up to 3^rd order. The first 4 terms correspond to the 0^th and linear order (x, y, and z) components of the dynamic field. The 0^th order term reflects changes in the global background field, while the h₁(r) = x, h₂(r) = y, and h₃(r) = z terms correspond to the usual k-space definitions for k_x(t), k_y(t), and k_z(t), respectively.

Spherical harmonics up to 3^rd order. This example shows surface plots of 16 spherical harmonic functions.

Finding the Positions and Resonance Frequencies of the Field Probes

The values of the spatial basis functions h_l(r_j) at each probe location r_j are necessary. Additionally, we must determine the resonance frequencies ω_ref(r_j) of the probes to obtain an accurate phase from each one. Thus, we first establish the positions and resonance frequencies of the probes through a calibration scan (shown below).

Left: Pulse sequence schematic for the calibration scan. The first free induction decay (FID) acquisition, conducted without active gradients, measures the probes’ frequencies at their respective locations. The following three scans, each activating one gradient axis while turning off the others, determine the probes’ spatial positions. The transistor-transistor logic (TTL) trigger activates the field camera’s data sampling. Right: Example of probe positions mapped onto a 3D grid of the scanner frame.

To accurately measure the positions, the field probe array must be placed at the isocenter of the magnet, where the gradients are linear.

The probe signals and the MR signal from imaging experiments must share a common demodulation frequency. This ensures that the phase coefficients from the probes accurately reflect the signal evolution during the imaging scan.

Solving for the Dynamic Field Component – the Phase Coefficients k_l(t)

To ensure a well-conditioned inverse problem for determining the phase coefficients, k_l(t), the number of probes N_P must be at least equal to N_L, the number of spatial basis functions (i.e., N_P ≥ N_L).

As mentioned above, the position r_j of each NMR probe (j = 1, 2, … N_P) is found in a separate calibration scan. The phase measured by the j^th probe (located at r_j) is

where γ_P is the gyromagnetic ratio of the NMR probe sample. The probes’ reference frequencies, ω_ref(r_j), are determined from a calibration scan conducted at a time that also establishes the reference state for the static field B_ref(r). Consequently, the field mapping sequence for estimating B_ref(r) should occur immediately after the calibration scan, as ideally ω_ref(r_j) = γ_PB_ref(r_j).

Given the probe phase signals ϕ_j​(t), the reference frequencies ω_ref(r_j), the probe positions r_j, and the values of the spatial basis functions at those positions h_l(r_j), the problem is to solve for the phase coefficients k_l(t). The above model can be represented in matrix form as

where the probing matrix P and the probes’ phase and reference frequency vectors are

This problem is posed mathematically as a least-squares problem to solve for k(t). Using the pseudoinverse of the probing matrix P, defined by P^†=(P^TP)⁻¹P^T , the solution is

Evolution of the Phase Coefficients k_l(t)

Using Spherical Harmonics up to 1st Order

An instructive example is to consider modeling the spatial component of the field with spherical harmonics up to 1st order, which corresponds to a set of N_L = 4 basis functions. In this case, the probing matrix simplifies to

The corresponding phase expansion is then

Neglecting the global phase term k₀(t) reveals the standard k-space representation. Here, only the linear spatial field variations induced by the applied imaging gradients contribute to the signal’s overall phase.

k-space Trajectory Measured in a Balance SSFP Scan

Balance steady-state free precession (bSSFP) sequences depend on the transverse magnetization remaining consistent from TR to TR. Field monitoring can detect violations of this steady-state condition. The figure below shows errant phase oscillations, particularly in the G_x and G_y gradients, that accumulate over the TR periods. This accumulation of phase errors leads to artifacts, such as banding, commonly seen in bSSFP images.

Top: k-space measurements from a balanced SSFP acquisition for the B0 and the X, Y and Z terms. Bottom: Zoomed version shows potential violation of balanced behavior due to oscillations and eddy currents.

Field Monitoring in Diffusion-weighted EPI

The figure below shows an example of the evolution of the phase coefficients k_l(t) during a diffusion-weighted Cartesian EPI scan. There are N_L = 16 phase coefficients corresponding to 16 spherical harmonic functions h_l(t) needed to model spatial phase variations up to 3rd order.

The first row in the figure displays the phase coefficients of the EPI scan without diffusion-sensitizing gradients. The second, third, and fourth rows illustrate the differences in phase coefficients compared to the top row, resulting from diffusion-sensitizing gradients applied along the read, phase-encoding, and slice-selected axes, respectively.

(a) Evolution of phase coefficients during a 45.8 ms single-shot EPI readout without prior diffusion gradients. All plots are scaled to reflect the maximum phase contribution that occurs within the imaging volume. Change in these coefficients induced by diffusion weighting in (b) the frequency-encoding direction, (c) the phase-encoding direction, and (d) the slice direction. Image taken from: Wilm et al., “Higher order reconstruction for MRI in the presence of spatiotemporal field perturbations”, Magn. Reson. Med. 2011; 65:1690-1701, doi: 10.1002/mrm.22767, link: Wilm et al. (2011).

Field Monitoring and the Expanded Signal Model

In common Fourier imaging, it is assumed that the first-order phase coefficients align perfectly on the nominal k-space grid defined by the imaging gradients. However, this assumption does not hold in practice, and severe image artifacts occur if there are large deviations from the expected trajectory. Alternatively, field monitoring can enhance image reconstruction by integrating the actual measurements of the phase coefficients into the signal model used for reconstruction.

The expanded signal model becomes

In this case, the phase coefficients measured by the field camera are directly integrated into the encoding matrix. Solving the inverse imaging problem using this more general, algebraic framework results in more accurate reconstruction.