Discussion

In this paper we present a perturbation method for calculating the $\mathrm{B}_0$ field for an object with varying spatial susceptibility. A fast, first-order calculation is presented for voxel-based objects, using the analytical voxel solution. The accuracy of this method is tested using the analytical solution for a sphere as well as phantom and human in-vivo data. These results indicate that highly localised errors of less than 1 ppm are achieved generally, which is very similar to other calculation methods, and sufficient for most MR imaging purposes.

There are several main contributions from this work. The first is the use of a principled, perturbation method for arriving at the field approximation. This is useful in that it allows the magnitude of the error terms (second-order and higher) to be estimated which then permits the relevant applicability of the method to be assessed. For instance, the method cannot be used for metallic objects where the susceptibility difference, $\delta$, is large, but can be used for some slightly higher susceptibility substances like graphite [11]. Without knowledge of how these errors scale it is not possible to know when to have confidence in applying an approximation without substantial experimental testing and validation. It is also possible, although potentially analytically intractable, to extend the approximation to higher orders to increase the accuracy. In addition, the formulation of the perturbation equations is separate from the object model specification and could be used with other object models, such as boundary element methods.

The other significant contribution of this work is the ability to calculate more than just the $z$ component of the field. In particular, the $x$ and $y$ components can be calculated just as easily (although separately) as well as the gradients of the fields (evaluated at the voxel centres), and formulations are provided for all these cases. More interestingly, it is possible to calculate the field at different object orientations by linearly combining `basis' images. This allows the field to be determined, without further approximation, at any orientation in a very efficient manner, if the basis images have been precalculated and stored. Such calculations will allow the interaction between susceptibility fields and motion artefacts to be explored more easily, a current research interest of the authors.

In the field calculations used here there are two main sources of approximation: (1) neglecting all terms beyond the first-order term and (2) representing the object by a voxel-based model. The first approximation limits the range of objects for which this method could be applied. For instance, it is not useful for metallic objects which have very large $\delta$ and potentially non-zero currents, but is applicable for the typical range of biological tissues encountered.

The second approximation is potentially more limiting, as the use of a voxel-based model for the object will cause errors that are not as easily estimated as the perturbation approximation errors. In particular, voxel-based models are likely to cause greater errors in the calculation for large voxel sizes, especially at the boundaries, as indicated in the sphere results (see Figure 2). By reducing the size of the voxels the spatial extent of this error can be reduced. Alternative models such as boundary element methods [4,1,3,2] are likely to be physically accurate in capturing the object shape, but have two main disadvantages. One is that boundary meshes are more difficult to instantiate from images and the second is that they require more computation for the field calculation as each element (triangle of the mesh) is potentially unique and requires separate calculations. In contrast, voxel-based models [12,5] are easy to instantiate and very efficient to calculate (using Fast Fourier Transforms). Furthermore, the numerical results on the spherical object indicate that similar errors are obtained, regardless of the method chosen. Finally, both of these object models have an advantage over finite Fourier representations [6] since they can ensure that the object has finite spatial extent, which is not possible with the Fourier method.