Abstract:
A method for correcting MR data elements acquired in the presence of an imperfectly linear gradient field, each element comprising nominal first and second quadrature components. The quadrature components are multiplied by an artifact correction factor to provide respective first and second quantities. A Gradwarp geometric correction operation is applied to respective quantities to provide a corrected first and second quadrature component, corresponding to each data element. For a given data element, the arctangent function is applied to the result obtained by dividing the first corrected component by the second corrected component to provide a corrected phase component for use in forming an image.

Description:
BACKGROUND OF THE INVENTION 
     The invention disclosed and claimed herein generally pertains to phase and phase difference magnetic resonance (MR) imaging. More particularly, the invention pertains to a method for correcting the effect of distortion in a set of MR data to be used in phase or phase difference imaging, wherein the distortion results from non-linearity in one or more gradient magnetic fields used in acquiring the data. 
     As is well known, acquired MR imaging data can be represented in quadrature. That is, the MR signal at each pixel location can be represented by two complex components which are 90 degrees out of phase, commonly referred to as the I and Q components, respectively. Thus, the MR signal at each pixel location has a magnitude equal to √(I 2  +Q 2 ), and a phase component equal to the arctangent of (Q/I), or tan -1  (Q/I). 
     In the past, it has been most common to employ signal magnitude in forming MR images. However, it is now recognized that for certain applications, such as imaging of flowing fluids (such as blood or cerebral spinal fluid) and shimming, there are distinct benefits in using the phase component at respective pixel locations to form MR phase and phase difference images. As is known by those of skill in the art, a phase difference image can be constructed by first obtaining two phase images of a subject. An imaging parameter is varied between the two phase images, so that the respective phases thereof are different. For example, in phase contrast angiography, also known as phase mapping, wherein flow velocity along an axis is to be measured, the value of the first gradient moment along the axis is varied to produce two different data sets. One of the phase images is then subtracted from the other, on a pixel-by-pixel basis, to provide the phase difference image. Thus, given two phase images φ A  =tan -1  (Q A  /I A ) and φ B  =tan -1  (Q B  /I B ), the phase difference Δφ may be obtained from subtraction, i.e., Δφ=φ B  -φ A . As known to those skilled in the art, the phase difference is more commonly calculated with a single arctangent via the expression: 
     
         Δφ=tan.sup.-1 {(I.sub.A Q.sub.B -I.sub.B Q.sub.A)/(I.sub.A I.sub.B +Q.sub.B Q.sub.A)}                                Equation 1 
    
     The single arctangent function ensures Δφ in Equation 1 covers a range of 2π radians or 360 degrees. The invention herein and the teachings set forth with respect thereto are considered to apply to both phase images and phase difference images. Hereinafter, the term &#34;phase information image&#34; is used to refer to both phase and phase difference images. 
     As is further well known, a set of MR image data may be acquired under conditions wherein the X-, Y- and/or Z- gradient magnetic fields are non-linear relative to their respective axes, across at least a portion of the field of view. Imperfect linearity, or nonuniformity, of the gradient, which results in image distortion, is of particular concern if the field of view is comparatively large, which often is the case, for example, for imaging the spine in the sagittal plane. Accordingly, techniques have been developed to geometrically correct an image for the distorting effects of gradient nonuniformity, as exemplified by commonly assigned U.S. Pat. No. 4,591,789, issued May 27, 1986 to Glover and Pelc. Such patent is directed to a technique hereinafter referred to as &#34;Gradwarp&#34;. 
     The Gradwarp technique applies a two step process to a set of MR data elements, wherein the elements comprise nominal values of an MR signal parameter at respective locations in a set of pixel locations. In accordance with the first step, a magnification/minification, pixel displacement or warp factor, which is determined by the gradient nonuniformity, is applied to displace or remap respective data elements to their actual locations. In the second step, an interpolation is performed, using the displaced data elements, to determine the actual or corrected values of the MR signal parameter at respective pixel locations. The interpolation step commonly uses polynomial interpolation, such as a cubic or quartic function, such as fitting a cubic spline curve through four displaced data elements close to a pixel location, and then sampling the curve at the pixel location. 
     Techniques such as Gradwarp have generally been quite effective in correcting distortion of the above type in images constructed from MR signal magnitudes. However, if the Gradwarp technique is applied to the respective phase values of a phase image, the corrected phase image, or phase map, will not be unwrap-compatible. Unwrap-compatibility is discussed in more detail hereinafter, in connection with FIG. 2 of the drawings. In addition, overshoot and ringing, also described hereinafter, can occur in applying the interpolation step to compute corrected phase values at pixel locations which are close to a boundary or interface between materials of very different MR properties, such as body tissue and air. Overshoot causes artifacts in the form of bright and dark spots which lie along the boundary. 
     SUMMARY OF THE INVENTION 
     A method is provided for correcting MR data, acquired for use in forming a phase information image, for the effects of an associated imperfectly linear gradient magnetic field. The method includes the step of reconstructing MR image elements from data acquired in the presence of the imperfectly linear gradient field, each of the image elements comprising first and second quadrature components corresponding to a pixel location in a set of pixel locations. The quadrature components corresponding to a pixel location are multiplied by an artifact correction factor, likewise corresponding to the pixel location, to provide respective first and second quantities. A geometric correction (Gradwarp) having a specified relationship to gradient non-linearity is applied to each of the first and second quantities to provide respective corrected first and second components, one of the corrected first components and one of the corrected second components corresponding to each of the pixel locations. For a given pixel location, the arctangent function is applied to the result obtained by dividing the first corrected component corresponding thereto by the second corrected component corresponding thereto, to provide a corrected phase component for use in forming a phase information image. 
     In a preferred embodiment, the artifact correction factor has the form M n , where, for a particular pixel location, M is the square root of the sum of the squares of the nominal first and second quadrature components respectively corresponding to the particular location. Preferably, n is on the order of -0.5. 
     Preferably also, the geometric correction step comprises displacing respective first and second quadrature components to their respective actual locations, and then employing displaced first quadrature components in performing interpolations to provide the corrected first components, and likewise employing displaced second quadrature components in performing interpolations to provide the corrected second components. 
     OBJECTS OF THE INVENTION 
     An object of the invention is to provide an improved method for correcting geometric distortion, resulting from gradient field non-linearity, in MR phase and phase difference imaging. 
     Another object is to provide a method of the above type which does not introduce other distortions or artifacts, such as those resulting from flawed unwrapping, or overshoot associated with Gradwarp interpolation. 
     Another object is to provide an improved method of the above type for angiography and shimming applications. 
     Another object is to provide a method of the above type for use in phase mapping or phase contrast imaging to image flow of blood or cerebral spinal fluid. 
     These and other objects of the invention will become more readily apparent from the ensuing specification, taken together with the accompanying drawings. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     FIG. 1 is a simplified block diagram showing structure for implementing an embodiment of the invention. 
     FIGS. 2a and 2b are a set of curves respectively illustrating seamless and flawed unwrapping. 
     FIGS. 3a and 3b are a set of graphs respectively illustrating overshoot, and the substantial reduction in overshoot achieved by an embodiment of the invention. 
     FIGS. 4a and 4b show a schematic set of images respectively depicting the overshoot artifact, and the correction thereof achieved by an embodiment of the invention. 
    
    
     DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT 
     Referring to FIG. 1, there is shown an MR imaging system 10 which, for purposes of emphasizing the invention, is shown to comprise a data acquisition system 12, reconstruction processor 14 and a data storage or display device 16. As is very well known, data acquisition system 12 generally comprises a configuration of gradient coils, main magnet and RF transmit and receive coils (not shown), together with components for driving the coils in accordance with specified imaging waveforms. Magnetic fields respectively generated by the driven coils relative to a subject of imaging (not shown) enable the receive coil to acquire a set of MR data for constructing an image of a section taken through the subject. The data set comprises the MR signals S. As stated above, the acquired MR data is commonly represented in the form of two quadrature components. This may be accomplished with a device 18 known as a quadrature splitter, coupled between acquisition system 12 and processor 14. The quadrature data are then reconstructed into a complex image with a complex Fourier transform. The Fourier transformation on the data provides the resultant complex image with I and Q channels, there being a 90 degree phase shift therebetween. The Fourier transform operation is carried out by processor 14, and is represented therein as a block 20. Thus, data for each pixel in the form of quadrature components I m  (x, y, z) and Q m  (x, y, z) are provided for further use by processor 14. 
     As also stated above, the phase value φ m  for the pixel (x 0 , Y 0 , z 0 ) can be found from the relation φ m  =tan -1  [Q m  (x 0 , y 0 , z 0 )/I m  (x 0 , y 0 , z 0 )]. However, in the event of, for example, X-gradient non-linearity, the values of Q m  (x 0 , y 0 , z 0 ) and I m  (x 0 , y 0 , z 0 ) would only be nominal with regard to pixel location (x 0 , y 0 , z 0 ), in that the MR signal represented thereby would have actually been acquired, for example, at a position (x 1 , y 0 , z 0 ). Accordingly, an image formed from respective phase values φ m  would be geometrically distorted. 
     The reconstruction processor 14 could generate phase values φ m  directly from I m  and Q m . The Gradwarp technique described above could then be applied directly to respective phase values φ m . However, it has been found that an image formed of phase data constructed thereby would not be unwrap compatible. 
     To illustrate unwrap compatibility, reference is made to FIG. 2(a), which shows a plot of signal intensity versus phase. The discussion which follows in connection with FIG. 2 assumes that the expression in Equation 1, above, or one which is mathematically equivalent thereto, was used. It will be readily apparent that all phase values derived with the arctangent function are limited to a range of 2π radians, such as -π to π (or -180 degrees to 180 degrees). However, it is frequently desirable to go beyond such limits, in order to extend the range over which signal intensity is linearly proportional to phase. This is accomplished by a technique known in the art as unwrapping or phase unwrapping. In accordance with such technique, the signal-phase curve of FIG. 2(a) includes a first branch B1 extending between -π and π, and a second branch B2 extending between π and 3π, both branches depicting linear relationships between signal intensity or contrast and φ. A phase value which exceeds π or 180 degrees by the amount δ, such as a δ of 10 degrees so that phase is 190 degrees, is represented as a phase of -180 degrees +δ, i.e., -170 degrees. Accordingly, as φ increases to 190 degrees, the signal intensity jumps to branch B2. The signal intensity value on branch B2 corresponding to 190 degrees is the same as the signal intensity value on branch B1 corresponding to a phase of -170 degrees. Thus, branch B2 extends the phase range by 2π. For FIG. 2(a), branch B2 could be displaced upwardly to align with linear branch B1, as shown by dotted line B2&#39;. FIG. 2(a) thereby depicts seamless unwrapping or unwrap compatibility, i.e., a continually linear relationship between signal and φ even as φ varies across π. 
     Referring to FIG. 2(b), there are shown branches B1 and B2 for respective corrected phase values provided by applying the Gradwarp process directly to respective phase values of a geometrically distorted phase image. Corrected phase values which are close to π will tend to be non-linear with respect to signal. Thus, branches B1 and B2 of FIG. 2(b) do not match proximate to π, since the upward displacement B2&#39; shown therein does not align with branch B1. Accordingly, the image provided by such corrected values is not unwrap-compatible. If a phase image is not unwrap compatible, significant artifacts can occur therein. Such artifacts are caused by non-linearities in the unwrapped signal intensity as a function of phase. 
     To provide a phase image which is corrected for gradient non-linearity, and at the same time is unwrap-compatible, the Gradwarp geometric correction process described above is independently applied to the I m  and Q m  quadrature components corresponding to each pixel location, by selected operation of processor 14. Thus, for a given location (x 0 , y 0 , z 0 ), the uncorrected or nominal value I m  (x 0 , y 0 , z 0 ) is multiplied by a warp factor to displace it to the actual location of the signal to which it corresponds. An interpolation is then performed, using the respective displaced I m  values, to determine a corrected quadrature value of I c  at the location (x 0 , y 0 , z 0 ). Usefully, such operation is expressed as I c  (x 0 , y 0 , z 0 )=GW{I m  (x 0 , y 0 , z 0 )}. Similarly, the Gradwarp technique is applied to Q m  (x 0 , y 0 , z 0 ) to compute a corrected value Q c  at pixel location (x 0 , y 0 , z 0 ). Such operation is usefully expressed as Q c  (x 0 , y 0 , z 0 )=GW{Q m  (x 0 , y 0 , z 0 )}. 
     After respective values of I c  and Q c  have been determined, the distortion corrected phase value φ c  (x 0 , y 0 , z 0 ) can be computed by operating processor 14 according to the following expression: 
     
         φ.sub.c =tan.sup.-1 [GW{Q.sub.m }/GW{I.sub.m }]        Equation 2 
    
     After computing a set of corrected phase values φ c , such values may be coupled to device 16 for storage, or to construct and display a phase image, according to requirements. 
     The technique of Equation 2 thus provides a phase information image which is both unwrap-compatible and corrected for non-linear gradient distortion. However, such image may still be affected by significant artifacts. As stated above, the interpolation step of the Gradwarp process commonly uses a cubic or quartic function. As a result, conditions referred to as over-shoot and ringing can occur at the boundary between tissue and air. This condition is depicted in FIG. 3(a), which shows signal strength versus spatial position along the X-axis. For purposes of illustration, the X-axis passes through the air-tissue boundary. The tissue produces signal, while the air does not. The cubic or quartic function causes the signal to &#34;over-shoot&#34; or &#34;ring&#34; proximate to the boundary, generating signal intensity peaks P 1  and P 2 , respectively, on opposite sides thereof. The signal peaks cause respective dark and light spots along the boundary, as shown by the arrows A l  and A d , respectively, in FIG. 4(a). 
     It has been found that the artifacts resulting from signal intensity overshoot can be substantially eliminated by multiplying the I m  and Q m  values corresponding to a given pixel location (x 0 , y 0 , z 0 ) by an artifact correction factor M n , prior to applying the Gradwarp correction process thereto. Thus, a modified set of corrected phase values φ c  &#39; may be generated for a phase information image, wherein: 
     
         φ.sub.c &#39;=GW [{Q.sub.m M.sup.n }/GW {I.sub.m M.sup.n }]Equation 3 
    
     For a pixel location having respective quadrature components I m  and Q m , the value of M to be used for Equation 3 is M =√(I m   2  +Q m   2 ). It has been found that a very useful value for the exponent n in Equation 3 is n=-0.5 at all pixel locations, although other values of n may be used as well. A value of n substantially equal to -0.5 has been found to be optimal in conjunction with cubic interpolation. Other forms of interpolation may benefit from different values of n. As shown by FIG. 3(b), the effect of the artifact correction factor M n  is to significantly reduce the MR signal levels between tissue and air during the intermediate step of Gradwarp application. Since the correction factor is applied to both the numerator and denominator of Equation 3, the over-all contrast in the phase image is not disturbed, as shown by FIG. 4(b) of the drawings. It is also shown by FIG. 4(b), that over-shoot and ringing resulting from Gradwarp interpolation in computing respective values φ c  &#39; produces substantially reduced peaks P 1  &#39; and P 2  &#39;, which are less apparent in a phase information image formed from respective φ c  &#39; values. 
     Referring once again to FIG. 1, computation of respective values φ c  &#39; by processor 14, in accordance with Equation 3, is represented by a block 22. 
     Obviously, many modifications and variations of the present invention are possible in light of the above teachings. It is therefore to be understood that within the scope of the disclosed concept, the invention may be practiced otherwise than as specifically described.