Imaging voxel distortions resulting from gradient non-linearities are minimal near the magnet’s isocenter, but are generally always present and increase with distance from isocenter. Troublesome in abdominal imaging, these distortions are generally subtle in neuro-applications, but still pose problems as they may introduce head-size bias into volumetric studies e.g. larger crania producing larger distortions. Hence, the correction of the gradient distortion, or gradient unwarping [1-2], is an important step for certain quantitative image analyses methods; e.g. volumetric measures of neuroanatomical structures for atrophy assessment. Gradient unwarping is typically applied as a black-box procedure, slice-wise (2D) post-processing step to magnitude images near the end of the manufacturer's reconstruction pipeline. Research applications, however, often source data at a more primitive step (raw k-space) and, in practice, usually neglect the gradient warp due to the lack of effective unwarping methods. These distorted images represent the basis of many sophisticated image reconstructions, such as phase MRI, non-Cartesian reconstructions, or compressed sensing. [3-6] In particular, when resulting images are to be compared with images reconstructed by the scanner software (commonly T1w images), it is essential to unwarp the research images to avoid bias. [1-2,7-8] While some manufacturers may provide information regarding their unwarping process, others do not, requiring a laborious measurement of the field non-linearities using special hardware.

In this work we present a simple method to reverse-engineer the manufacturer’s correction of the gradient distortions, and use this information to unwarp images.

The foundations of gradient unwarping and the use of spherical harmonics to describe gradient field non-linearities are detailed in references. [2,7,9]

Phantom: Our methodology requires a phantom that: (a) is large enough to fill the entire possible imaging volume of the scanner, and (b) contains sufficiently detailed MR-visible structure. We constructed such a phantom by stacking the magnet bore with 48 consumer water bottles (500ml) (Figure 1a).

Data Acquisition: A GE-system was used because the spherical harmonic coefficients used by the GE platform are contained on the console and thus can be used as “ground truth” to determine the accuracy of our method. We used a T₂* gradient echo (GRE) sequence, TR:2.7 ms, TE:0.74 ms, flip angle:3^o, field of view: 440x440x367 mm³, slice thickness:1.7 mm, acquisition matrix 256x256x216, voxel size: 1.72x1.72x1.7 mm³. GRE images were acquired twice, once each with the manufacturer's gradient unwarping processing turned on or off. As the scanner's unwarping only operates in 2D, multiple slice orientations (transversal and sagittal) were required to measure the 3D warp information. We operated our GE 3T Signa Excite HD 12.0 whole-body scanner (General Electric, Milwaukee, WI) in body-coil transmit and receive mode. Figure 1b-c shows representative slices of warped (1b) and GE-unwarped (1c) images.

Warpfield determination: Non-linear 2D registration was performed with ANTS [10] on a slice-wise basis to match the warped and unwarped images and the resulting slice-wise warpfields were merged into a single displacement field.

Spherical harmonics fitting: Spherical harmonics were fitted to the ANTS displacement field using in-house developed code (Python using Scipy and Matplotlib libraries [11-12]) solving for each spatial direction the following system: $$$W=SC$$$, where $$$W$$$ is a vector concatenating the displacement in a certain direction, $$$C$$$ concatenates the spherical harmonic coefficients and $$$S$$$ is the spherical harmonic expansion system matrix. The system was solved in a least-squares sense using the pseudo-inverse of $$$S$$$.

Demonstration: The gradient displacement fields from the fitted and system coefficients were generated, and applied to warped images of a structure phantom. Details are in reference [8].

Figure 2 contains the generated displacement fields from our fit, the manufacturer and the absolute difference.

Figure 3 shows the results of gradient unwarping in a phantom. Note that in this case the phantom was intentionally positioned far from isocenter in order to exacerbate the effects of the gradient distortions.

The fitted harmonics were in very good agreement with GE's supplied coefficients (see Figure 2), especially near the isocenter where the majority of imaging applications apply. Figure 3 indicates only minor deviations even at the extreme situation presented here.

A further strength of this methodology is its repeatability and ease of setup – different phantoms can be created, imaged and the resulting coefficients can be averaged to obtain even better estimates of the spherical harmonic coefficients.

The displacement field generated from our method produced acceptable unwarping results, demonstrating the algorithm's efficacy in reverse engineering a gradient displacement field. This method allows researchers to reverse-engineer the gradient non-linearities inherent in their research images, without laborious measurements or expensive hardware.