Simulation of the influence of magnetic field inhomogeneity and distortion correction in MR imaging

Mogneric Resonance Imaging, Vol. 8, pp. 483-489, Printed in the USA. All rights reserved.

l

0730-725X/W f3.M) + .OO Copyright 0 1990 Pergamon Press plc

1990

Original Contribution SIMULATION

OF THE INFLUENCE OF MAGNETIC FIELD INHOMOGENEITY AND DISTORTION CORRECTION IN MR IMAGING JAN WEIS AND EUBOS BUDINSK?

Institute of Measurement and Measuring Engineering, Elektrophysical Research Centre, Slovak Academy of Sciences, Dtibravska cesta 9, 842 19 Bratislava, Czechoslovakia We describe a technique for simulation and correction of the effects of an arbitrary distribution of undesired components of the static and gradient magnetic fields. This technique is applicable to direct Fourier NMR imaging. The mathematical basis and details of this technique are fully described. Computer simulation demonstrates the effectiveness of this method. Keywords: Magnetic resonance imaging; Magnetic field inhomogeneity; correction.

INTRODUCTION

computer simulation of the quadratic inhomogeneity for conventional projection methods. The magnetic field homogeneity may be disturbed by objects having magnetic susceptibility different from that of tissue. Also, differences in susceptibility between tissues may give similar effects. The principles governing the resulting changes in the image have been analysed by Ludeke et al. l2 for those simulations where the disturbance is relatively weak. For very strong disturbances produced by a steel ball, the validity of the formalism of Ludeke has been confirmed by Ericsson et a1.4 On the other hand, the literature is more rich in the methods for image restoration from nonuniform magnetic field influences. Hutchison et al.s have proposed a method to correct the phase and amplitude distortions for the spin-warp method. Lai3 has proposed a curvilinear back-projection method. Sekihara et al. ,1,6and Kawanaka and Takagi’ have proposed a restoration technique for direct Fourier imaging methods. Feig et a1.8 have proposed a combined remapping and magnetic field estimation method and Chan et a1.9 have proposed a modified transform method. In this paper we describe a simple and universal technique of simulation and correction of spin-density images for an arbitrary distribution of error magnetic fields in a direct Fourier transform imaging method. This method can be used for accurate planning of

In magnetic resonance imaging the object is placed in a highly homogeneous static magnetic field and high linear orthogonal gradient magnetic fields. Generally, a static magnetic field with a uniformity of 10-6-10-5 and linearity of gradient magnetic fields less than 3% are required. If these magnetic fields are not homogeneous, geometrical and amplitude distortions appear in the images. In our opinion, simulation of the influence of the undesired components of the magnetic field on the spin density is a very useful instrument for designing electromagnets, gradient magnetic coils, compensation coils, for theoretical study of the pulse sequences, and also for pedagogical purposes. It is surprising, therefore, that relatively few authors have taken an interest in this topic, up to the present time. For example, Sekihara et al.’ proposed an analytical solution of the imaging equation (see Eq. (2)) for direct Fourier methods. But this solution exists only if the error magnetic field distribution is in the form of a plane, rotational or elliptic functions. For image reconstruction a discrete Fourier transform must be used. Also, Wong and Rosenfeld,’ by numerically solving the Bloch equations, demonstrated the advantages of the spin-inversion imaging pulse sequence. Lai3 has proposed

RECEIVED 7/3 l/89; Acknowledgment-We

ACCEPTED

Computer simulation; Distortion

assistance in the MR measurements and for helpful discus-

2/8/90.

thank our colleagues for their

sions during the course of this work. 483

484

Magnetic Resonance Imaging 0 Volume 8, Number 4, 1990

stereotactic neurosurgery, interstitial radiotherapy and other procedures in which accurate positioning are critical, and such distortions are not negligible. MATHEMATICAL

BASIS

We shall assume that two-dimensional Fourier transform method is used to acquire data using nonideal stationary and gradient magnetic fields. The pulse sequences are shown in Fig. 1. As shown in Fig. 1, the direction of the static magnetic field is taken as the zaxis and the read-out gradient G, is placed along the x-axis and phase encoding gradient G,, is placed along the y-axis. Because the slice-select gradient is taken as G,, a 90” rf pulse inclines the magnetization which is in a planar slice perpendicular to the z-axis. When undesired components exist, the selected slices do not become a plane perpendicular to the z-axis but this effect is still negligible for the measurement of practical NMR imaging equipment. Therefore, in this paper this effect will be neglected. The distribution of magnetic fields on the slice can be expressed by means of two-dimensional coordinates (x,y). Static magnetic field B, and gradient fields G,, GY are defined as follows:

where B. is the uniform component of the static field, g, and g,, are linear components of the gradient fields and AB(x,y), Ag,(x,y), and Ag,(x,y) denote those undesired components. For the observed spinecho (SE) signal after quadrature detection at the angular frequency w = -yBo (y-gyromagnetic ratio) we can write I3: s(7,1)

= K

=

p(x,y)exp[-iy(AB(x,y) (2)

+ G,(x,Y))~

-~r~‘Gy(-wVld~dy.

In Eq. (2) we neglected the term describing the relaxation effect, K is a constant expressing the sensitivity of the receiver, etc., I= -L, -L + 1, . . . - l,O, 1,. . . L - 1, L, p(x,y) is spin density, and 7 = f - 2F (Fig. 1). Defining the new coordinate system (x’,y’) as transformation @-‘(x’,y’) : x=x’

- AB(x,y)&

- &x(x,y)4tx

(3) Y = Y’ - Ag,

ky)&

for the observed signal we can write: s(~,f) = K

B(x,y)

ss

Bo + AB(x,y)

ss

P’(X’,Y’)

(4)

x exp( -iyg,x’T

G,(x,Y) = gxx + &x(x,y)

(1)

where p’(x’,y’)

G,(x,Y) = guy+ &y(w)

- iyg,T’,y’l)dx’dy

= p(~-‘(x’,y’))/v(~-‘(x:yl))

and V(x,y) = V(+-‘(x’,y’) transformation3:

n’““’

B, n90.

v(a-‘(x’,y’))

=

(5)

is the Jacobian of the

6x’/6x

6x’/6y

6y’/6x

6y’/6y

.

(6)

I

---__ ---_

Equations (3) and (5) indicate that the undesired components of the static magnetic field AB(x,y) and gradient magnetic field Ag,(x,y) , Agv (x, y) caused geometrical and amplitude distortions of the spin-density images.

; 1

I ,I,

Gx ’

I

Tl

/

I

s

I

COMPUTER

Ti

Ifi I

0”

!’ 0 Fig. 1. Pulse sequences

and received

AT

signal.

SIMULATION

In practice, spin-density images can only be obtained at discrete points on an image matrix consisting of N * N pixels. Therefore, as a first step of simulation we change the continuous variables (x,y) and (x’,y’) into pixel units (1, J), (I’, J’), respectively:

Simulation of influence of magnetic field inhomogeneity and distortion 0 J. WEISSAND

Z = IFIX(x/Ax

+ N/2)

J = IFIX(N/2

- y/Ay)

(7)

I’ = IFIX(x’/Ax’) J’ = IFIX(y’/Ay’) where function IFIX indicates the maximum integer not greater than A; Ax,Ay are the pixel widths in the x and y directions, and Ax’ = Ax, Ay’ = Ay are the pixel widths in the x’ and y’ directions and as follows from Eq. (3): + Ag,(Z,J))/(g,Ax) (9)

y’/Ay’ = J + Ag,(Z,J)/(g,Ay).

We assume that the distributions of the original spin density p (x, y), undesired components AB (x, y), Ag, (x, y), and Ag,,(x,y) are known at image matrix (Z,J). The amplitude distortion5 we compute at the second step of our simulation: P,(Z,J)

= P(Z,J)/V(Z,J)

where V(Z,J) is a correction matrix calculated by numerical derivations of the AB (Z, J) , Ag, (Z, J) , Ag,(Z,J) at Eq. (6). In our case we performed these derivations using a cubic spline algorithm. The following interpolation was used for geometric distortion: P’(Z;J’)

= (1 - Q)(l

- DJ)P,(Z,J)

+ D,(l - DJ)P,(Z_ + D,D,p,(Z-

P(Z,J) = P,U,J)W,J) To our knowledge, Sekihara did not mention this possibility of simulation. He used for simulation the method briefly described in the introduction of this paper. The main disadvantage of the above-described methods is the use of numerical derivation for computing the matrix V(Z, J). Another source of error is the geometrical interpolation procedure. This error results from the implicit assumption of the square form of the pixel (Z,J) at raster (I’, J’). Therefore we propose the second procedure for simulation and distortion correction, which follows from strictly geometrical interpolation of the raster (Z, J) at raster (I’, J’). For simplicity we assume that undesired components of gradient magnetic fields Ag, (x,y) = Ag,,(x, y) = 0, e.g., J = J’. Again, at the first step of simulation we change the continuous variable (x, y), (x’, y’ ) into rasters (Z,J), (I’, J’) by Eq. (7), (8), (9). Then by a suitable interpolation method we compute for each pixel Z at row J = J’ the magnetic inhomogeneity B,, B, on the left and right edge of the pixel, respectively. According to Eq. (8) and (9) we compute the position of the left and right edge of the pixel (Z,J) at raster (I’, J’): L = Z - 0.5 + B,/(g,Ax) R = Z + 0.5 + B,(g,YAx)

Z,’ = IFIX

1,J)

1,J - 1)

+ (1 - D,P,P,(Z,J

485

Step 2

and

x’/Ax’ = Z+ (AB(Z,J)

L. BUDINSKI

- 1)

where DI = x//Ax - Z’, DJ = y’/Ay - J’. We note that this method of simulation is a reverse procedure to the image restoration technique, which was originally proposed by Sekihara et al.’ This restoration algorithm uses the following relations:

Z; = IFIX D/ = L - Zi D, = R - Z; .

This transform is shown in Fig. 2 for two possible cases a, 6. Then we compute the number of the pixels p’(Z’, J’) by pixel p(Z, J):

Step 1 P,(Z, J) = (1 - &I(1 - D,)P’(~‘, J’)

+ D1(l - D,)p’(Z’ + (1 - D,)D,p’(Z’,J’ + D,D,p’(Z’

+ 1, J’) + 1)

+ l,J’ + 1)

n = z; - z; + 1.

If n = 1 (Fig. 2(a)) for distorted spin density, p’( Z’,J’) can be written: p’(Z’, J’) = p’(Z’, J’) + p(Z, J).

Magnetic Resonance Imaging 0 Volume 8, Number 4, 1990

486

?‘(I;=I;,J'l

degradation by magnetic field inhomogeneity was estimated by Sekihara et al.’ and these artifacts were experimentally demonstrat$l%j+rYoung et al.‘O RESULTS AND DISCUSSION

I-J ?(I,Jl

?‘(I;,J')

?'U',,J'l 4

I-

F

Fig. 2. Two kinds of the transform

-I

of the raster

(I,J)

into

(Z‘,J’).

If n > 1 (Fig. 2(b)) then: p’(K,J’)

= P’(K,J’)

+ W(WP(Z,J)

where K=Z;,z;+ 1,...z; W(Z/) = (1 - Z&)/X, W(Z;) = 0,./X,, W(Z/ + 1) = W(Z/ + 2) = . . . = wtz: = l/X, X,=D,-D,+n-1.

- 1)

We present here the computer simulation by the second above-mentioned method. We point out that differences between the first and second methods in the presented examples are less than the variation between neighbouring gray levels of the intensity scale (0,255), e.g. less than 100/255 = 0.4%. For simulation we used the theoretical head phantom” shown in Figs. 5(A) and 6(A), and the real image of the human hand shown in Fig. 7. This axial human hand image (10 mm thick) was measured by a home-built experimental MR tomograph with a resistive magnet, having 400 mm room-temperature bore, and 0.0855 T field strength. The images consist of a (128, 128) image matrix. The pixel widths Ax = Ay = Ax’ = Ay’ were equal to 1 mm. For simulation and experimental measuring of human hand images the read-out gradient g, = 0.9175 lOA T/m was used. Also, we assume linear gradient magnetic fields, e.g., Ag,(x,y) = Ag,(x,y) = 0. The error field distribution was assumed to take the form of a hyperbolic paraboloid AB(x,y) = ABm,,(x2 - y2) shown in Fig. 3 and the form of real magnetic field at our experimental tomograph shown in Fig. 414 measured in selected points by the NMR gaussmeter. We note that hyperbolic paraboloid was used, because for this form of inho-

As follows from the above-mentioned equations, the weighting coefficients W(K) express a percentual part of the pixel area p (Z, J) which is superimposed on the pixels p( Z’,J' ). It should be pointed out that these weighting coefficients are calculated under the implicit assumption of linear dependence of inhomogeneity between values B,, B,. By analogy, the image restoration from a nonuniform magnetic field influence can be written: p(Z, J) = (1 - D,)p’(Z;, J) + p’(Zi+

1, J)

+ . . . + p’(Z: - 1, J) + D,p’(Z:, J)

for n > 1 (Fig. 2(B)) and ~(4 J) = (D, - D/)p’(Z;

= Z:, J)

for n = 1 (Fig. 2(A)). We point out that in our simulation and distortion correction we do not consider degradation of the images by decreasing the signal-to-noise ratio. This

Fig. 3. Error paraboloid.

field distribution

of the form of hyperbolic

Simulation

of

influence of magnetic field inhomogeneity and distortion 0

J. WEISS AND L. BUDINSK*

487

(D), and Figs. 6(B), (C), (D)) are not presented, because they are almost entirely free from the influence of field nonuniformities and therefore equal to the original images (Figs. 5(A), 6(A)). The results, applying the proposed correction method to experimental data, are shown in Fig. 7. The distorted image of a human hand is shown in Fig. 7(A). This image was measured by our experimental tomograph without the use of the compensation magnetic coils system. Inhomogeneity in a place of the measured slice was approximately +50 ppm. Restored data are shown in Fig. 7(B). CONCLUSION

Methods have been developed for simulation and correction of the spin density images for an arbitrary distribution of error magnetic fields. The advantage of the proposed methods of simulation is the fact, that they do not require application of a discrete Fourier transformation and magnetic field inhomogeneities only in matrix form need to be known. Therefore, the proposed correction method can be easily implemented in existing NMR imaging systems.

Fig. 4. Real error field distribution tomograph.

of an experimental

mogeneity there is no analytical solution of the imaging Eq. (2) as performed by Sekihara’ for rotational

and elliptic paraboloid. Figures 5(B), (C), (D) show images resulting from this simulation when a hyperbolic paraboloid (Fig. 3) was used for ti(x,u) and the maximum field nonuniformity in the place of the phantom was +25, +50, +-75 ppm, respectively. Figures 6(B), (C), (D) show images distorted by the real form of the magnetic field nonuniformity .(Fig. 4), but the maximum field nonuniformity in the place of the phantom was chosen to be +25, +50, f75 ppm, respectively. Maximum geometrical distortion introduced by inhomogeneity +25 ppm is approximately 2 mm and maximum relative error of the brightness (Fig. 5(B)) is up to - 15%. Maximum geometrical distortion introduced by inhomogeneity +75 ppm is approximately 6 mm and maximum relative error of the brightness (Fig. 5(D)) is up to -37%. Images restored from a distortion (Figs. 5(B), (C),

Fig. 5. Simulation of the influence of the hyperbolic paraboloid inhomogeneity. (A) Original phantom, (B)-(D) distorted images by nonuniformity k25, *SO, f75 ppm, respectively.

488

Magnetic Resonance Imaging 0 Volume 8, Number 4, 1990

D Fig. 6. Simulation of the influence of the real form inhomogeneity. nonuniformity *25, f50, +75 ppm, respectively.

(A) Original phantom,

(B)-(D) distorted images by

Simulation of influence of magnetic field inhomogeneity and distortion 0

J. WEISS AND L. BUDINSK~

489

REFERENCES 1. Sekihara, K.; Kuroda, M.; Kohno, H. Image restoration from non-uniform magnetic field influence for direct Fourier NMR imaging. Phys. Med. Biol. 29(l): 15-24; 1984. 2. Wong, T.S.; Rosenfeld, D. Spin-inversion imaging: A technique for NMR imaging under magnetic fields with high field nonuniformities. IEEE Trans. Med. Imaging MI-6(2):148-156;

1987.

3. Lai, Ch. Reconstructing NMR images under magnetic fields with large inhomogeneities. J. Phys. E: Sci. Znstrum. 15:1093-1100; 1982. 4. Ericsson, A.; Hemmingsson, A.; Jung, B.; Sperber, G.O. Calculation of MRI artifacts caused by static field disturbances. Phys. Med. Biol. 33(10):1103-l 112; 1988. 5. Hutchison, J.M.S.; Sutherland, R.J.; Mallard, J.R. NMR imaging: Image recovery under magnetic fields with large inhomogeneities. J. Phys. E: Sci. Znstrum. 11: 217-221; 1978. 6. Sekihara, K.; Matsui, S.; Kohno, H. NMR imaging for magnets with large nonuniformities. IEEE Trans. Med. Imaging MI-4(4):193-199;

1985.

7. Kawanaka, A.; Takagi, M. Correction method of image distortion due to non-uniformity of static magnetic field in NMR imaging. Proc. IEEE Comp. Sot. Znt. Sympos. Med. Images Icons 173-177; 1985.

8. Feig, E.; Greenleaf, F.; Perlin, M. Magnetic resonance imaging with non-uniform fields. Phys. Med. Biol. 31(10):1091-1099; 1986. 9. Chan, S.K.; Wong, T.S.; Rosenfeld, D. A correction technique for NMR imaging with large magnetic field School. Elect. Engg., Univ. Sydney, non-uniformities. Tech. Rep.; 1984.

10. Young, I.R.; Cox, I.J.; Bryant, D.J.; Bydder, G.M. The benefits of increasing spatial resolution as a means of reducing artifacts due to field inhomogeneities. Magn. Reson. Imaging 6:585-590;

1988.

11. Shepp, L.A.; Logan, B.F. The Fourier reconstruction of a head section. IEEE Trans. Nucl. Sci. NS-21:21-43; 1974.

12. Liideke, K.M.; Roschmann, P.; Tischler, R. Susceptibility artifacts in NMR imaging. Magn. Reson. Imaging 3:329; 1985.

Fig. 7. Transversal human hand image: (A) Real measurement. Inhomogeneity in a place of measured slice was k50 ppm. (B) Corrected image.

13. Hinshaw, W.S.; Lent, A.H. An introduction to NMR imaging: From the Bloch equation to the imaging equation. Proc. IEEE 71:328; 1983. 14. Weis, J.; Frollo, I.; Budinsky, B. Magnetic field distribution measurement by the modified FLASH method. Z. Naturforsch. 44a (in press).