Fast Unwrapping using Discrete Gradient Evaluation (FUDGE): an analytical correction to the Laplacian-based phase unwrapping technique for discrete data.
Amanda Ching Lih Ng1, Meei Pyng Ng2, Sonal Josan3, Shawna Farquharson4, Claire Mulcahy4, and Roger J Ordidge1

1Dept of Anatomy & Neuroscience, The University of Melbourne, Parkville, Australia, 2Dept of Mathematics & Statistics, The University of Melbourne, Parkville, Australia, 3Siemens Healthcare, Melbourne, Australia, 4Imaging, The Florey Institute of Neuroscience and Mental Health, Melbourne, Australia


Laplacian-based phase unwrapping is commonly used to pre-process phase for methods such as Quantitative Susceptibility Mapping (QSM). However, the formulation was derived with the assumption of a continuous signal and a continuous Fourier transform. When applied to discrete MRI phase data, serious errors in phase can occur, resulting in substantial errors in QSM estimates. We present a mathematically correct Laplacian-based phase unwrapping formula, based on the assumption of the discrete nature of MRI phase data and processing. Our results reflect the mathematical predictions of the old and new formulations.


The phase of the complex-valued MRI gradient echo data contains information on local magnetic field changes. The phase is implicitly wrapped into the range $$$\small{\left[-\pi,\pi\right)}$$$. In order to recover the field information, the phase needs to be unwrapped. Laplacian-based phase unwrapping [1-3] is a commonly used technique, particularly as a pre-processing step in Quantitative Susceptibility Mapping [4]. The method unwraps phase by calculating the Laplacian of the unwrapped phase, $$$\small{\phi_u}$$$, using the wrapped phase, $$$\small{\phi_w}$$$:

$$\phi_u = \nabla^{-2} \left( \cos \phi_w \nabla^2 \left( \sin \phi_w \right) - \sin \phi_w \nabla^2 \left( \cos \phi_w \right) \right), \tag{1}$$

where $$$\small{\nabla^2}$$$ and $$$\small{\nabla^{-2}}$$$ are the Laplacian and inverse Laplacian operations, computed using the time derivative property of the Fourier transform. Whilst this is analytically correct, it was derived with the assumption of a continuous signal and employing the continuous Fourier transform. These assumptions do not translate well to discrete signals and the use of the discrete Fast Fourier Transform (FFT), which are present in the context of MRI phase data and computation. In fact, the FFT of the discrete Laplacian operator is commonly used to calculate the Laplacian operations. Yet, even when employing this modification, the resulting unwrapped phase (for 1D) becomes

$$\phi\left[x\right] = \nabla^{-2} \left[ \sin\left(\phi_{w}\left[x+1\right]-\phi_{w}\left[x\right]\right) - \sin\left(\phi_{w}\left[x\right]-\phi_{w}\left[x-1\right]\right) \right] ,\tag{2}$$

where $$$\small{x}$$$ is the discrete image space coordinate. From this, it can be seen that a change of 0 is indistinguishable from a change of $$$\small{\pi}$$$,resulting in incorrect unwrapped phase values. We present an analytically correct alternative to this method, by correctly assuming the discrete nature of the data and computations.


We derived a discrete version of Laplacian-based unwrapping, where the second order partial derivative is

$$\frac{\delta^2\phi_{u}}{\delta t^2} = \tan^{-1}\frac{\sin\left(\phi_{w}\left[t+1\right]-\phi_{w}\left[t\right]\right)}{\cos\left(\phi_{w}\left[t+1\right]-\phi_{w}\left[t\right]\right)}-\tan^{-1}\frac{\sin\left(\phi_{w}\left[t\right]-\phi_{w}\left[t-1\right]\right)}{\cos\left(\phi_{w}\left[t\right]-\phi_{w}\left[t-1\right]\right)},\tag{3}$$

the 3D Laplacian is defined as$$ \nabla^2 \phi_{u} = \frac{\delta^2\phi_{u}}{\delta x^2} + \frac{\delta^2\phi_{u}}{\delta y^2} + \frac{\delta^2\phi_{u}}{\delta z^2}\tag{4}$$and the inverse Laplacian is calculated using the Discrete Cosine Transform:

$$\phi_{u}\left[x,y,z\right] = \text{DCT}^{-1}\left(\left[2\cos\left(\frac{\pi}{N}\left(k-1\right)\right)-2\right]^{-1}\text{DCT}\left(\nabla^2 \phi_{u}\right)\right).\tag{5}$$

We shall refer to the original algorithm as Laplacian-based Uwrapping using Continuous K-space (LUCK) and the new algorithm as Fast Unwrapping using Discrete Gradient Evaluation (FUDGE). We applied both methods to two numerical simulation datasets and an experimental in vivo MRI gradient echo data. The numerical simulations comprised: a phase image with spheres of varying phase values; a phase image containing a sphere (diameter=5mm, susceptibility=0.3ppm, B0=3T, TE=20ms). In vivo MRI brain data from a healthy volunteer (35yo, female) was acquired on a 7T research MRI scanner (Siemens, Erlangan, Germany) in accordance with local ethics (fully flow-compensated GRE, TE=7.65ms,15.3ms, TR=18ms, FA=13$$$^{\circ}$$$, iPAT=3, matrix=366$$$\scriptsize{\times}$$$316$$$\scriptsize{\times}$$$224, voxel=0.6mm isotropic, TA=7:42min). iLSQR [5] was used to calculate QSM for the second numerical simulation.


LUCK produced erroneous phase in the first numerical simulation (Fig. 1), whilst FUDGE produced phase values exact to within numerical floating point error. The second numerical simulation (Fig. 2) demonstrated similar results. The susceptibility of the sphere was estimated by the iLSQR method at 0.281ppm, 0.239ppm and 0.281ppm for the true phase, LUCK unwrapped phase and FUDGE unwrapped phase, respectively. For LUCK, this equated to 85.2% of the susceptibility calculated using the true phase. The unwrapped in vivo phase images (Fig. 3) demonstrated much less error in the FUDGE image than the LUCK image, particularly in the vicinity of the periventricular vasculature. FUDGE was 4.9 times faster to compute than LUCK on the in vivo data.


The simulation results reflected the errors predicted by the analysis of the mathematical formulation of LUCK. Similarly, the lack of error in the FUDGE results reflects the expected accuracy when using a mathematical formulation based on the correct assumptions of discrete data and computation. The inaccuracies in unwrapped phase produced with LUCK have expected effects on subsequent QSM maps, with the 15% error in susceptibility being substantial, particularly for a quantitative method. The in vivo results demonstrated greater accuracy with FUDGE compared to LUCK. The presence of open-ended fringe lines, a phase feature common at high field (e.g. 7T), resulted in some residual error in the phase unwrapped with FUDGE, although to a much lesser extent than in the LUCK image. FUDGE was also much faster, due to there being only 2 FFTs in the algorithm, compared to 6 FFTs in the LUCK algorithm.


The analytical, numerical and experimental results demonstrate that relying on LUCK will result in erroneous phase data and QSM maps, however data processed with FUDGE produces numerically accurate results.


No acknowledgement found.


[1] Schofield, Marvin A., and Yimei Zhu. 2003. “Fast Phase Unwrapping Algorithm for Interferometric Applications.” Optics Letters 28 (14): 1194–96. doi:10.1364/OL.28.001194.

[2] Li, Wei, Alexandru V. Avram, Bing Wu, Xue Xiao, and Chunlei Liu. 2014. “Integrated Laplacian-Based Phase Unwrapping and Background Phase Removal for Quantitative Susceptibility Mapping.” NMR in Biomedicine 27 (2): 219–27. doi:10.1002/nbm.3056.

[3] Li, Wei, Bing Wu, and Chunlei Liu. 2011. “Quantitative Susceptibility Mapping of Human Brain Reflects Spatial Variation in Tissue Composition.” NeuroImage 55 (4): 1645–56. doi:10.1016/j.neuroimage.2010.11.088.

[4] Wang, Yi, and Tian Liu. 2015. “Quantitative Susceptibility Mapping (QSM): Decoding MRI Data for a Tissue Magnetic Biomarker.” Magnetic Resonance in Medicine 73 (1): 82–101. doi:10.1002/mrm.25358.

[5] Li, Wei, Nian Wang, Fang Yu, Hui Han, Wei Cao, Rebecca Romero, Bundhit Tantiwongkosi, Timothy Q. Duong, and Chunlei Liu. 2015. “A Method for Estimating and Removing Streaking Artifacts in Quantitative Susceptibility Mapping.” NeuroImage 108 (March): 111–22. doi:10.1016/j.neuroimage.2014.12.043.


Figure 1. Numerical simulation images and line profiles demonstrate substantial errors in the LUCK-processed phase and exact results from the FUDGE processing

Figure 2. Errors produced in unwrapped phase by LUCK translate to inaccurate susceptibility after QSM processing.

Figure 3. In vivo brain data: FUDGE produced less error in the unwrapped phase compared to LUCK.

Proc. Intl. Soc. Mag. Reson. Med. 24 (2016)