Abstract:
In a magnetic resonance imaging system, object movement along an axis associated with readout magnetic gradient pulses is detected, quantized and compensated for in reconstructing an image slice. Phase-encoding magnetic gradient pulses in the sequence are arranged so that low-order phase-encoding pulses are distributed substantially uniformly during scanning time. In an N-slice acquisition, object movement is determined y arranging each of the N sets of phase-encoding gradient pulses comprising a multi-slice sequence such that at each scan time, at least one of the N sets contains a low-order phase-encoding gradient pulse. During data acquisition, echoes returned from low-order phase-encoding gradient pulses have sufficient signal-to-noise ratios to identify object position along the readout axis. Identification of object movement from one slice may be valid for other slices relatively adjacent physically in the object, so that all N slices are movement corrected during image reconstruction. Two-dimensional object position may also be identified and corrected by suitably exchanging the axes of N sets of pulses within a sequence. For a single slice acquisition, object position may be identified when scans associated with low-order phase-encoding gradient pulses are distributed evenly during imaging. After one-dimensional Fourier transformation, object position relative to a reference box function is determined using cross-correlation.

Description:
FIELD OF THE INVENTION 
     This invention relates to magnetic resonance imaging, and more specifically to methods for compensating for patient motion in magnetic resonance imaging. 
     BACKGROUND OF THE INVENTION 
     Magnetic resonance imaging (&#34;MRI&#34;) is a useful diagnostic procedure that analyzes changes in magnetic orientations to produce an image of an object under examination. Unfortunately, movement of the object during MRI relative to a reference position adversely affects the quality of the image produced. 
     FIG. 1 shows a typical MRI system 2, a system with which the present invention may be practiced. MRI system 2 includes a powerful typically superconducting magnet 4 that creates a longitudinal magnetic field h o  along the z-axis when power supply 6 is energized. To promote homogeneity of field h o , magnet 4 includes shim coil mechanisms (not depicted) that are coupled via bus 8 to a computer system 10 that generates appropriate corrective shimming gradient signals. 
     The object 12 to be imaged is placed within a radio frequency (&#34;RF&#34;) coil 14 within magnet 4, and is thus exposed to the magnetic field h o . Object 12 includes atoms whose nuclei have typically randomly oriented magnetic dipoles. During MRI, the h o  field polarizes the population of such nuclei (or &#34;spins&#34;), causing the spins to align parallel to the magnetic field direction, here along the z-axis. A net magnetization M z  results. 
     Computer 10 via bus 8 and RF coil 14 subjects the spins to a brief pulse of radio frequency (&#34;RF&#34;) energy at a resonant frequency called the Larmor frequency (RF source not depicted), thus exciting the oriented spins. This excitation temporarily reorients the spin axis of the nuclei, tipping the nuclei magnetic moment towards a plane transverse to the z-axis, i.e., the x-y plane. After termination of the RF excitation pulse, the net transverse magnetization is detected (or acquired) in real time as RF signals through RF coil 14 and a receiver (not depicted). 
     These signals are digitized (digitizing mechanism not depicted) and coupled to computer 10 for processing. The acquired raw RF signals include magnetic field gradient information, and represent the spatial frequency information domain, or so-called raw k-space data matrix (or array). By using Fourier transformation upon the raw k-space data matrix, a reconstructed image 16 of object 12 may be obtained and displayed on monitor. However, if object 12 moved during scanning, the quality of image 16 can be degraded. 
     The amplitude of the detected nuclear magnetic resonance (&#34;NMR&#34;) signal depends upon the degree of RF excitation and the relaxation properties of the spins. During imaging, computer system 10 controls magnetic field gradients G x , G y , G z  along the x, y, z axes that are needed to spatially correlate locations of the spins in object 12. 
     Typically 256 y-axis gradient levels are used during an MRI scan, which levels (or phase encoding steps) range from -128 through 0 to +128, where the 0th phase encoding for a slice refers to the line (or &#34;echo&#34;) having no phase encoding gradient. Each line or echo will comprise 256 points that correspond to one step along the y-axis in the so-called k-space (or momentum space) domain. Essentially these gradient levels allow system 2&#39;s sensitivity to be limited to a single x-y plane slice along the z-axis. 
     By suitably controlling the gradient strengths G x , G y , G z  during each MRI image cycle, the spatial distribution of nuclei spin excitation may be controlled, and the location of the resulting NMR signals can be identified. The NMR data for reconstructing images may be collected using several known techniques, including multiple angle projection reconstruction and Fourier transform. Such techniques typically comprise a pulse sequence made up of a plurality of sequentially implemented views. Each view may include one or more NMR experiments, each using an RF excitation pulse and a magnetic field gradient pulse to encode spatial information into the resulting NMR signal. As is understood by those skilled in the relevant art, the NMR signal may be a free induction decay or, preferably, a spin-echo signal. 
     As further understood by those skilled in the relevant art, low frequency information including motion information is present only in the center of the raw k-space data, with higher frequency information appearing at the beginning and end of this data matrix. 
     Computer system 10 uses the gradient magnetic field information to identify where spatially the intensity values represented by the acquired k-space radio data belong. Such techniques are commonly found in clinical MRI systems using two dimensional multi-slice imaging sequences, wherein the acquisition of the data for the slices is interleaved in time. 
     Computer system 10 transforms the k-space data to a real space image using conventional spatial Fourier transform methods. The real space image 16 may then be displayed on a monitor 18, to show contrast changes as a function of x-axis and y-axis positions. Typically image 16 depicts a &#34;slice&#34; of object 12 in the x-y plane, at chosen different locations along the z-axis. 
     If object 14 moves during imaging (depicted by dashed lines 20), so-called motion artifacts result that can degrade the quality of the reconstructed image 16. However if object 12 motion could be measured, it would be possible to compensate for that movement and generate an improved reconstruction image. 
     Although the k-space data includes information that can identify object movement, it has been difficult in the prior art to measure such movement directly from the k-space data at all times during the scan. This difficulty arises because, as depicted in FIG. 2, the signal to noise (&#34;S/N) ratio is poor during most of the image scan time. 
     FIG. 2 depicts a normal prior art collection of k-space data for a single slice along the z-axis. Assume that a total of 256 lines are to be scanned from time 0 to time 256, wherein one scan acquires data for one line. As noted, each line will have 256 data points that correspond to one k-space step along the y-axis, and collectively the 256 lines represent 256 magnetic field gradient levels along the y-axis. According to the prior art, the same acquisition order is used for each slice, e.g., here -128, through 0 to +127. 
     According to FIG. 2, although the scan starts at time t=0 (or line -128) when possibly object 10 movement occurs, the S/N is relatively low (e.g., poor) and artifact movement signals cannot be readily detected. It is not until about time t=96 to time t=159 that the S/N ratio improves adequately, with time t=128 (corresponding to line 0) representing peak S/N ratio. Thus t=128, which corresponds to a low order phase gradient (e.g., 0) representing low frequency data, when the center portion of the k-space data is acquired represents the best time to detect artifact signals representing motion, because at all other acquisition times the S/N is relatively inadequate. 
     Within this region of high S/N, object 12 movement could be detected, the movement displacement could be calculated, and a suitable correction could then be made computer system 10, using well known techniques, to improve the quality of the displayed reconstructed image 16. Unfortunately, however, object movement that occurred during time t=0 to time t=95, as well as movement occurring between time t=159 to t=256 will go unrecognized because the artifact signal is too submerged in the noise, due to poor S/N. Thus no displacement can be ascertained, no correction can be generated for such movement, and the resultant display image quality remains degraded. Stated differently, prior art approaches such as depicted in FIG. 2 are &#34;blind&#34; to object movement occurring other than at the period of optimum S/N. 
     One prior art attempt to reduce motion artifacts in a magnetic resonance imaging system is set forth in U.S. Pat. No. 4,937,526 issued to Ehman, et al. (1990). Applicants incorporate Ehman, et al. herein by reference, as this patent provides a very detailed description of MRI theory, accompanied by numerous explanatory figures. 
     Ehman, et al. discloses a method wherein additional RF and magnet gradient pulses are generated to create a second data matrix (in k-space) having a maximum S/N at all lines. This second matrix is then used to provide positional information identifying artifact movement. The positional change information is then used to improve the resultant image quality. 
     While Ehman, et al.&#39;s system improves image quality, the improvement requires considerable overhead. For example, Ehman, et al.&#39;s technique requires additional imaging time due to the necessary generation of additional RF and gradient pulses. Further, Ehman, et al. generates essentially twice the amount of initial k-space data to be stored in memory 14, and presents additional data to be processed by computer system 10, which processing requires additional time. 
     What is needed for MRI is a method that can recognize and compensate for object 12 motion, utilizing the k-space data itself generated during normal MRI scanning. Preferably such method should not require that memory 14 be increased beyond what is needed to store the image k-space data itself, should not require the generation of additional RF and gradient pulses (thereby not extending imaging time), and should not impose substantial additional processing time upon computer system 10. Such a method preferably should measure object translational movement in the readout direction by self-calibration, using only high S/N echoes from different slices. The present invention discloses such a method. 
     SUMMARY OF THE PRESENT INVENTION 
     The present invention provides methods for compensating for artifact motion in magnetic resonance imaging, without requiring additional computer memory beyond that needed to store the conventional k-space data, without requiring additional RF and gradient pulses, and without imposing substantial additional computational requirements. As such, the present invention may be used with existing MRI systems. As noted, although the raw k-space data matrix includes information needed to identify artifact motion, the signal/noise (&#34;S/N&#34;) ratio throughout the matrix is not sufficiently to allow accurate object position determination. 
     The present invention provides a measure of translational object motion in the readout direction by self-calibration, wherein only the high S/N ratio echoes are used. Essentially the present method offsets the starting phase encoding number of each slice in a two-dimensional multi-slice imaging system. This permits acquisition of high S/N echoes (from different slices) to be spread out reasonably evenly during the entire imaging time period. Assuming that each slice is displaced by the same amount in the readout direction, object displacement data from one slice may be used to correct displacement in other slices. The method works especially well when a small number of relatively thin slices are used to cover a small physical volume. 
     A first embodiment is directed to single axis motion correction. In this embodiment, the order in which the raw k-space data (e.g., phase encodings) for each slice is acquired is rearranged to distribute low-order phase encodings uniformly throughout the time of the scan. This rearrangement allows estimates of object displacement in one dimension, using low-order (e.g., high S/N ratio) phase encoding echoes. 
     Other preferred embodiments are directed to two-axis motion correction for multi-slice imaging. For example, a second embodiment uses a different phase encoding order for each slice. Data acquisition for each slice starts with a different phase encoding number, using interleaving, and within each slice, the phase encoding number is sequentially incremented by one (modulo the total number of echoes) for each echo. Again, object position is estimated from each slice&#39;s low-order (e.g., high S/N) phase encodings. A final displacement curve (displacement vs. time) is obtained by combining the position information obtained from all slices. This final motion data will be used to correct motion for all slices. Time interpolation of these position measurements, or displacement data, is used to obtain the shift values for high-order phase encoding echoes that were acquired between two neighboring points. 
     Furthermore, different readout directions for different slices may be used, to allow object positional measurement in two dimensions. For example, the first slice may have an x-axis readout direction and a y-axis phase encode direction. The second slice would have y-axis readout direction and x-axis phase encode direction. 
     Object position is then estimated from each radio frequency (&#34;RF&#34;) echo by using edge information obtained from the modulus of the one-dimensional Fourier transform of that echo. As used herein, &#34;modulus&#34; means the vector sum of a complex value, for example, the modulus of A+iB being (A 2  +B 2 ) 1/2 . The positional information is then used to correct the raw data using Fourier shift methods, so that the reconstructed image corresponds to one position of the object, with minimal motion artifacts being displayed. In contrast to the prior art, according to the present invention, object movement can be detected and measured at several sample time intervals, and essentially no &#34;blind&#34; regions are encountered. 
     Other features and advantages of the invention will appear from the following description in which the preferred embodiments have been set forth in detail, in conjunction with the accompanying drawings. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     FIG. 1 depicts a conventional magnetic resonance imaging system, a system with which the present invention may be practice; 
     FIG. 2 depicts a raw data (or k-space) matrix, showing acquisition order for two-dimensional imaging in a conventional magnetic resonance imaging system wherein the same acquisition order is used for each slice for multi-slice imaging using interleaving; 
     FIG. 3A depicts a k-space data matrix for slices 1, 3, 5 and 7, showing acquisition order for recording motion in the x-axis direction for an MRI signal processing system, according to the present invention; 
     FIG. 3B depicts k-space data matrix for slices 2, 4, 6 and 8, showing acquisition order for recording motion in the y-axis direction for an MRI signal processing system, according to the present invention; 
     FIG. 4 depicts a spin-echo pulse sequence for multi-slice two-dimensional imaging, wherein object motion is detected and measured along the x-axis and y-axis according to a second embodiment of the present invention; 
     FIG. 5 depicts a preferred phase encoding order used with the image pulse sequence of FIG. 4, wherein object motion is detected and measured along two axes, according to a second preferred embodiment of the present invention; 
     FIG. 6 depicts the relationship between image slices and readout gradients for the embodiment of FIG. 5, according to the present invention; 
     FIG. 7 depicts readout direction so orientated that the region of the object being imaged is totally enclosed within the field of view, according to the present invention; 
     FIG. 8 depicts a box-car function superimposed on the one-dimensional Fourier transform of the object echo, where the box-car function is used as a reference for motion detection, according to the present invention; 
     FIG. 9 is a flow chart depicting motion detection and measurement of object motion along two axes, according to the present invention; 
     FIGS. 10A and 10B depict measured displacement data along the x-axis and y-axis, using the present invention; 
     FIGS. 11A-B, and 11C-D depict, respectively, a one-dimension Fourier transform of 0th and 32nd phase coding echoes, and the corresponding cross-correlation function with a reference box-car (displacement expressed in video pixel units); 
    
    
     DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT 
     Applicants have discovered that to meaningfully use the information contained in a conventional multiple slice, k-space data array to identify artifact motion, it is necessary to rearrange the order in which data for each slice are acquired, or collected. 
     According to the present invention, rearrangement is such that the low-order phase encodings (e.g., high S/N, corresponding to center of k-space) are distributed reasonably uniformly throughout the time of the scan. As a result, high S/N information is available at any time during imaging. This advantageously permits motion recording using high S/N ratios echoes to be made for the entire imaging period. Applicants, invention with its re-ordered phase encodings is compatible with any of the standard magnetic resonance imaging pulse sequences, including multi-slice spin echo, fast spin echo, and gradient echo. 
     Object position is then estimated from each radio frequency (&#34;RF&#34;) echo by using the modulus of the one-dimensional Fourier transform of that echo, using edge detection information, in one or in two dimensions. The positional information is then used to correct the MRI image data so that the reconstructed image corresponds to one position of the object, with minimal motion artifacts being displayed. 
     FIG. 4 schematically depicts the spin-echo pulse sequence for four image slices, for a second preferred embodiment wherein motion detection and measurement is provided a both x and y axes. While only four image slices are depicted, thereafter additional slices may be added by repeating these four slice as long as the cycle time (TR) permits. As will be understood by those skilled in the art of MRI, the RF pulse train shows a relatively small 90° pulse, followed by a larger 180° pulse. The various gradient signals are denoted G x , G y , G z , and A/D denotes the time when data are collected. In FIG. 4, alternate slices have their phase encoding directions in the y-axis and x-axis directions. While FIG. 4 depicts four slices, it is understood that a different number of slices could be used if desired. 
     For the embodiment depicted in FIG. 4, the collection order of phase encodings in a single slice is rearranged to provide a more reasonably uniform time distribution of low order phase encodings. For example, one collection order for phase encodings of a 256 phase encoding image, wherein one in every four phase encodings lies in the optimal range -32 through +32: -128, 0, -64, +64, -127, +1, -63, +65, -126, +2, -62, +66, and so forth. 
     In the above sequence, [-128, 0, -64, +64] are data for the first four echoes, [-127, +1, -63, +65] represent data for the second set of echoes, and so forth. Although other orderings are possible, it is desirable to maintain the phase encoding order monotonic (modulo the total number of phase encodings), to avoid introducing strong artifacts from uncorrected object 12 motion. 
     FIG. 5 depicts the preferred phase encoding order used for the image pulse sequence shown in FIG. 4, while FIG. 6 depicts the relationship between image slices and readout gradients for the same embodiment. Preferably a monotonic phase encoding order is maintained with multi-slice images by using a different ordering for the phase encodings in each slice. In FIG. 5, the collection order of phase encodings in a single slice is rearranged to provide a more reasonably uniform time distribution of low order phase encodings. In FIG. 5, the cross-hatched &#34;rectangles&#34; schematically denote regions of high S/N, for example during acquisition of k-space regions -32 through 0 through +32. FIG. 5 depicts one collection order for phase encodings of a 256 phase encoding image, wherein one in every four phase encodings lies in the optimal range -32 through +32: -96, -32, +32, +96, -95, -31, +33, +97, -94, -30, +34, +98, -93. -29, +35, and so on. 
     It should be appreciated that eight images are simultaneously acquired, wherein four images are acquired from the x-axis, and four images are acquired from the y-axis. With respect to the above phase encodings, [-96, -32, +32, +96, -95, -31, +33, +97] represent, respective, four x-dimension and four y-dimension slice of the same respective four images, and so forth. Note that the preferred encoding is such that the first four slices are acquired 64 &#34;units&#34; apart, e.g., during high S/N, and that +1 is added to the phase encoding in going from an image&#39;s first slice to that image&#39;s second slice, e.g., in going from -96 to -95, then to -94, and so on. 
     Note from FIG. 5 that at any given time, at least one slice is being acquired that has optimum S/N. For example, object 12 motion occurring during time 0≦t≦64 may be recognized and quantized because slice 3 and slice 4 each experience high S/N during that time interval. However if movement occurred, for example, during 64≦t≦96, then detection and measurement of such movement is possible because during this time interval, slices 1 and 2 each experience high S/N. Thus, unlike the prior art approach of FIG. 2, wherein there were &#34;blind&#34; regions at 0≦t≦96 and at 156≦t≦256 that precluded movement detection and measurement, the present invention has no blind spots. 
     Although many other orderings are possible, it is nonetheless desirable to maintain the phase encoding order monotonic (modulo the total number of phase encodings) to avoid introducing strong artifacts from uncorrected object 10 motion. Note from FIG. 5 that preferably there are no discontinuities that could introduce artifacts in acquiring data from one slice to another slice. 
     Note from FIG. 6 that slices 1, 3, 5, and 7 are readout using the G x  gradient, and slices 2, 4, 6, and 8 are readout using the G y  gradient. 
     It is understood that the readout direction is to be oriented such that the region of the object being imaged is totally enclosed within the field of view, this requirement being depicted in FIG. 7. In FIG. 7, the cross-section slice of object 10 is denoted as 30, whereas the field of view of the image acquired is denoted as the bounding box 32. 
     As noted, other readout arrangements are possible as it is not necessary that the slice data be collected in sequential order, or that the starting phase encodings be arranged in order. For example, for an eight-slice image with 256 phase encodings (i.e., phase encoding numbers ranging from -128 to +127), the following acquisition sequence could be used wherein a 32 echo sequence is employed: 
     
         ______________________________________Slice 1, phase encoding              0Slice 2, phase encoding              -32Slice 3, phase encoding              -64Slice 4, phase encoding              -96Slice 5, phase encoding              -128Slice 6, phase encoding              96Slice 7, phase encoding              64Slice 8, phase encoding              32Slice 1, phase encoding              1Slice 2, phase encoding              -31Slice 3, phase encoding              -63etc.______________________________________ 
    
     For example, for the same number of slices and phase encoding as discussed above, another useful acquisition sequence could be: 
     
         ______________________________________Slice 1, phase encoding              -112Slice 2, phase encoding              -80Slice 3, phase encoding              -48Slice 4, phase encoding              -16Slice 5, phase encoding              16Slice 6, phase encoding              48Slice 7, phase encoding              80Slice 8, phase encoding              112Slice 1, phase encoding              -111Slice 2, phase encoding              -79etc.______________________________________ 
    
     In this manner, the present invention avoids artifacts that arise when the first phase encoding collected on a slice is phase encoding number zero. 
     However the above scheme could result in one slice causing weak artifacts on the next slice, due to imperfect slice profiles, especially where a zero order phase encoding is followed by an echo. However those skilled in the relevant art will appreciate that standard MRI techniques exist to reduce this effect. For example, one such known technique is simply not to collect slices in the same order as their physical positions in space. Another technique is to use tailored RF pulses and slice selection gradient waveforms that produce better slice profiles. Yet another known technique is to tailor the gradient waveform shapes to minimize the generation of stimulation echoes. Those skilled in the art further recognize that combinations of these techniques may be used as well. 
     The modulus of the one-dimensional Fourier transform of each echo must be computed to measure the position of the object under examination, here object 10. A reference echo is needed for each slice, which reference echo may normally be selected as the 0th phase encoding for that slice, i.e., the echo with no phase encoding gradient. The position of object 10 has to be measured from each reference echo, which measurement can be done with a simple threshold, such as an edge detection algorithm. 
     To either side of the object 12, the one-dimensional transform will be just noise, but within object 12 there will be finite signal with noise superposed. This S/N relationship is depicted in FIG. 8. In FIG. 8, the modulus of the one-dimensional Fourier transform of the reference echo (e.g., 0th phase encoding) 34 is plotted against the box-car function 36 corresponding to the width (W) of object 12. The object width (W) is determined by the left-most and right-most positions where the transform crosses a threshold value, the box-car height indicating the threshold value chosen. 
     The position of the object 12 as a function of time can be determined by comparing each of the moduli of the transforms of the echoes with the modulus of the transform of the reference echo. Prior art comparison techniques, such as that proposed by Felmlee, et al. in Radiology 179, pp. 139-142 (1991) wherein threshold detection determines the boundaries of the non-zero signal region of each transformed echo, suffer from several drawbacks. Felmlee, et al.&#39;s technique does not work well for in vivo measurements, as signals from all phase encodings are weak except for the 0th phase encoding. Also, simple threshold techniques tend to be sensitive to noise as the threshold level must be lowered to match the reduced signal levels. 
     By contrast, applicants&#39; position measuring algorithm uses a box-car function derived from the 0th phase encoding, as depicted in FIG. 8. The box-car function is zero outside the object, and has a constant non-zero value (conveniently taken to equal one) inside the object. This box-car function is cross-correlated with the moduli of the transforms of all the echoes, except the reference. The positions of the maxima of the cross-correlations indicate the positions of the object 10 at the times of those echoes. This position information can be used to reduce the effects of motion by correcting the phase of the data. Various publications by Felmlee, et al. disclose such phase correction, including U.S. Pat. No. 4,937,526 (1990), Radiology 179, pp. 139-142 (1991), and Magnetic Resonance Imaging, 8S1, p. 25 (1990). 
     As noted, high-order phase encodings do not have sufficient S/N ratio to allow accurate estimation of displacement at the time of those echoes, and thus the present invention uses information from one slice to correct other slices. Because of the order of acquisition of the phase encodings, at any time at least one of the slices has echoes with high S/N ratios. 
     According to the present invention, using position information from one slice to correct another requires knowledge of the relative positions of the edges of the two slices. This information can be obtained by using position measurements on the two slices made at similar times. If extra 0 phase encoding echoes are included at the start of the imaging sequence, these can be used to make the relative position measurements. Alternatively, regular phase encodings may be used for the measurements. 
     Optionally, if desired one could smooth out position measurements as a function of time since in practice object 10 movements are relatively infrequent. 
     The present invention may be extended to measure position in two dimensions by interchanging the directions of the readout and phase encoding directions on alternate slices. For example, on slices 1, 3, 5 and 7 the readout gradient can be in the x direction and phase encoding gradient can be in the y direction, whereas on slices 2, 4, 6 and 8 the phase encoding gradient can be in the x direction and the readout gradient in the y direction. The above-described method is then applied separately for the odd and even-numbered slices. The y positions of the odd-numbered slices are interpolated from the measurements made on the even-numbered slices, and the x positions of the even-numbered slices are interpolated from the odd-numbered slices. 
     Motion measurement and correction according to the present invention may advantageously be used in conjunction with other numerical techniques, for example to obtain initial estimates of the position of object 12. Further numerical techniques such as the algorithm described by R. W. Gerchberg and W. O. Saxton in Optik, 35, pp. 237-248 (1972) may be employed to improve those position estimates. In many cases such numerical algorithms do not work unless the initial position estimates are reasonably accurate. 
     FIG. 9 is a flow chart depiction of the present invention showing determination of displacement along a given axis, here the x-axis. In FIG. 9, it is understood that the nomenclature &#34;1D FT&#34; denotes the one-dimensional Fourier transform, and that by practicing the procedure shown, displacement of object 12 along the x-y plane may be ascertained. In FIG. 9, using the hybrid data sets and reference box car functions, cross-correlation occurs to yield the displacement data sets. The displacement data sets are aligned and the four data sets are combined to yield the displacement data for motion in the x-axis direction. Of course, it is understood that the same procedure is followed using echoes from slices 2, 3, 4 and 8 to obtain displacement data for y-axis direction motion. 
     FIGS. 10A and 10B depict actual test data showing determination of object 12 displacement along the x-axis and y-axis, according to the present invention. In obtaining the data depicted, a test subject&#39;s knee was examined, and the subject was instructed to move his knee alternatively along one axis, and then along the other axis. However simultaneous movement of object 12 along the x-y plane could be similarly detected and measured, according to the present invention. 
     FIGS. 11A-11D depict various measurements according to the present invention. In FIG. 11A, the one-dimensional Fourier transform is shown for the 0th gradient, in a format similar to what is depicted in FIG. 8. The cross-correlation of these h o  (x) with reference box-car functions is depicted in FIG. 11B, wherein zero x-axis displacement is determined from maxima. FIG. 11C is similar to FIG. 11A, except that a lower S/N phase encoding echo is shown, namely phase encoding number -32. The cross-correlation of what is shown in FIG. 11C is depicted in FIG. 11D. FIG. 11D is similar to what is shown in FIG. 11B except that there is now depicted an x-axis object displacement or about 20 pixels, in a rightward direction. 
     Modifications and variations may be made to the disclosed embodiments without departing from the subject and spirit of the invention as defined by the following claims.