Iterative methods for image reconstruction problems in MRI are essential for advanced techniques including: compressed sensing, low rank reconstruction, SENSE reconstruction¹, non-Cartesian trajectories^2,3, field correction for off-resonance distortion², and joint-reconstruction based parameter estimation⁴. Iterative reconstructions can enable high resolution imaging with high SNR, such as submillimeter isotropic diffusion weighted acquisitions on a 3T clinical scanner⁵, by taking advantage of parallel imaging and physics-based field correction. However, all of these advanced reconstruction techniques typically have high computational demands and may not be available to clinicians and researchers.

Fostering robust research and exchange in MR reconstruction techniques requires sharing open and extensible code for others to refine, extend, and improve. Previous work with graphics processing units (GPUs) as hardware accelerators for image reconstruction, such as the IMPATIENT toolkit^6-8, produced very high speedup factors over CPU reconstructions for non-Cartesian field corrected SENSE reconstructions. However, its implementation is unfriendly for MR physicists to extend and adapt, particularly to accommodate additional problem formulations such as motion correction models or alternative regularization schemes. In this work, we utilize compiler directives (available through OpenACC 2.0) to target GPU accelerators with intuitive, high level C/C++ syntax to build extensible reconstruction software and tools. The result, PowerGrid, is an open source toolkit that makes accelerators and iterative reconstruction as accessible to MR physicists as the commonly-adapted Image Reconstruction Toolbox (IRT)⁹.

Building upon a set of hardware-accelerated operators such as the discrete Fourier transform (DFT) and the nonuniform fast Fourier transform (NUFFT) via gridding, objects implementing reconstruction techniques were coded in C++ to take advantage of the accelerated transforms. Figure 1 shows an example of how these objects maintain similar syntax and structure to those based in MATLAB from the IRT⁹ including a demonstration of a complete gradient descent solver in 6 lines of PowerGrid/C++ code with intuitive syntax using freely available libraries¹⁰. Table 2 shows a list of features implemented as part of PowerGrid.

The NUFFT and the DFT were written in C++ for CPUs and annotated for use with GPU accelerators (PGI/NVIDIA OpenACC Toolkit 15.7) as shown in Figure 2. The OpenACC compiler directives enable the use of CPU based C code on accelerators such as the GPU with suggestions to the compiler about parallelism and data movement.

Test data sets were acquired on a Siemens 3 T Trio scanner using spiral readout acquisitions and either a 32-channel or 12-channel head coil. The 12-channel data was compressed to 4 virtual coils using Siemens hardware coil compression. Reconstructions were performed on a quad-core Intel Xeon E5620 Linux workstation, 12 GB RAM, with a Tesla K40 GPU and Ubuntu Linux 14.04.

Table 2 compares the IRT/MATLAB reconstructions on CPU only versus GPU enabled reconstructions utilizing PowerGrid. For the smaller 2D dataset (192x192x1x32 coils, R=2), the speed up factor was much more modest (1.69x) due to the dominance of CPU-GPU transfers and CPU-based linear algebra vs. the small Fourier transforms transforms on GPU. For larger datasets, such as 3D datasets, the NUFFT on GPU via gridding performed better than the IRT/MATLAB and the DFT on GPU.

For the smaller 3D test dataset (64x64x16x4 coils, R=1), PowerGrid had a speed up factor of 8.96x via GPU acceleration versus the fastest MATLAB implementation in the IRT. The larger 3D test dataset (120x120x16x32 coils, R=1), Figure 3, had a slower speedup of (2.23x) over the IRT implementation due to the increased number of coils as the SENSE operator is not currently GPU accelerated and additional coils results in additional CPU-GPU transfers and overhead.

Finally, the modular nature of PowerGrid was demonstrated with a nonlinear motion induced phase correction algorithm¹¹ derived from the SENSE operator. The phase corrected example sped up by 4.0x over the IRT.

PowerGrid provides access to GPU acceleration for iterative reconstructions that is both adaptable and extensible for the implementation of additional reconstruction algorithms. By maintaining an object and operator based syntax, implementing novel algorithms in C++ are simplified in terms of programming time and required expertise. In our initial observations, reasonably well-experienced IRT users can perform an initial port of their MATLAB-based code to C++ using PowerGrid in less than a day. The combination of an open-source, freely available (http://mrfil.github.io/PowerGrid) toolkit with field inhomogeneity correction and accelerated non-uniform Fourier transforms makes the practical use of longer readout, non-Cartesian trajectories more accessible and more computationally tractable.