JPH0647021A - Motion correction method in mr imaging - Google Patents

Abstract

PURPOSE:To provide a motion correcting method by detecting a motion in parallel movability in a slice surface generated suddenly by only one of the upper half or the lower half of measuring data. CONSTITUTION:Measuring data are divided into a center low zone and a high zone and their motions are detected and corrected by different methods, respectively. After the motion correction processing of the center low zone is executed, as for the motion correction to the high zone, when the measuring data are divided into an upper half and a lower half, the data of a half of the side in which the motion is generated are estimated by utilizing a complex conjugate point symmetric property from the data of a half of the side in which the motion does not occur, a phase difference between the estimated data and the measuring data of a half of the side in which the motion is generated is taken, and by estimating a phase deviation caused by the motion from the phase difference data by the minimum square method to detect the motion, the motion correction processing is executed. In such a manner, the motion correction can be executed without deteriorating S/N of an image.

Description

Detailed Description of the Invention

[0001]

BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention relates to in-vivo tomography utilizing the nuclear magnetic resonance phenomenon, and in particular, the movement of an object to be inspected during radiography, in particular, the artefact caused by parallel movement in a slice plane. Regarding the method of suppressing.

[0002]

2. Description of the Related Art Since MRI has a long imaging time, artifacts due to the movement of an object to be inspected are likely to occur. Therefore, many methods of correcting the movement have been proposed. The following are typical examples of parallel movement.

(1) Using the Gerchberg-Saxton algorithm SRM, Book of Abstracts (SMRM) (1990) p561 SRM, Book of Abstracts (1991) p743 (2) Using projection data Using navigation echo Radiology 173, p225-263
(1989) without navigation echo Radiology 179, p139 (199)
1) Here, the influence of the parallel movement on the measurement signal will be summarized.
First, all of the methods proposed below are directed to what motion has occurred during each phase encoded signal. This is because the measurement time of each phase-encoded signal is several milliseconds, but the measurement waiting time is as long as several hundreds of milliseconds, so it can be assumed that movement occurs during this measurement waiting time. Based on the −to−View effect. In this case the original signal S (k _x, k _y) and the signal S including the effect of translation in the slice plane '(k _x,
The k _y), the following relation.

[0004]

[Number 1] _{S '(k x, k y} ) = exp _{[φ (k x, k y)} i ] _{S (k x, k y)} ... (1)

[0005]

[Number 2] _{φ (k x, k y)} = - 2π (p (k y) k x / N x + q (k y) k y / N y) ... (2) where, k _x = -N _{x / 2~N x / 2-1, k x} = -N y /
2 to N _y / 2-1, N _x , N _y is the number of sample points, p
_{(k y), q (k} y) , the line k _y -1 and k happened readout direction between the _y, a positional shift amount of the phase encode direction, i is the imaginary unit.

That is, the parallel movement appears as a linear phase error in the measurement signal.

(1) utilizes the fact that the inspection object is in a bounded area. First, the boundary line of the inspection object is determined in the image and divided into the inside and the outside. Then perform the following procedure.

In the image, a portion protruding from the boundary to the outside is cut off as a false image.

The image data obtained in (1) is subjected to inverse Fourier transform to be converted into signal data.

The phase difference between the signal data obtained in step 1 and the original measurement data is taken, and a linear phase error is estimated by applying the least squares method to the phase difference data.

A phase correction process is performed on the original measurement data based on the phase error estimated at.

The signal data obtained in (4) is Fourier transformed into an image.

The above operation is repeated until the phase error of the signal data is sufficiently corrected.

This method is effective for small movements of 1 pixel or less.

(2) The hybrid space data obtained by one-dimensional Fourier transforming the data of each line in the reading direction shows the position information of the inspection object along the reading direction at the time of measuring each line. It uses what you have, and there are the following two methods.

Using a navigation echo During the measurement of each phase encode signal, a signal called a navigation echo with a zero phase encode (zero encode signal) is measured, and from the absolute value of the hybrid space data obtained from the navigation echo, The amount of positional deviation in the reading direction is detected by detecting the edge.

Without navigation echo In the −32 to +32 lines at the center of the measurement space where the S / N is high, even if the navigation echo is not used, from the absolute value of the hybrid space data for each line in the reading direction, By detecting the edge, it is possible to accurately detect the positional deviation amount in the reading direction between the respective measurements. However, since the S / N of the measurement signal becomes lower as it goes to the line in the high frequency range, the edge detection error becomes large, and it becomes difficult to accurately detect the positional deviation amount in the high frequency range. Therefore, when measuring the data required to reproduce the image, a small marker that clearly appears in the hybrid space data in the high frequency part is set outside the inspection object, and this marker is measured. The amount of positional deviation in the reading direction is detected.

The item (2) is effective for a large motion of 1 pixel or more.

The problems of the above methods (1) and (2) are as follows.

(1) Utilization of Gerchberg-Saxton algorithm With respect to a large parallel movement of 1 pixel or more, the effect does not appear unless the algorithm is repeated a considerable number of times. Therefore, since the Fourier transform is included twice in the algorithm, a lot of calculation time is required until the effect appears.

The boundary line of the inspection object on the image must be determined manually or semi-automatically by a technique such as region growing, which is troublesome.

(2) Using projection data (i) Using navigation echo In addition to the echo for reproducing a normal image, the navigation echo is measured, which increases the measurement time.

Only parallel movement in the reading direction can be corrected.

(Ii) Those that do not use navigation echo Even if a marker is attached, a detection error of several tens of pixels often occurs in the high frequency region.

Only parallel movement in the reading direction can be corrected.

As described above, each of the parallel movement correction methods usually has advantages and disadvantages.

[0027]

DISCLOSURE OF THE INVENTION The present invention is to perform a motion correction process by detecting a motion from the measurement data itself, without using a special sequence, for a sudden parallel movement in the slice plane. It is also effective for movements of 1 pixel or more. Sudden movements of the patient usually occur within a limited time during the measurement, such as the first half of the measurement until the patient calms down or the second half of the measurement when the patient becomes tired and can't stay still. I think that it occurs in only half. For example, when parallel movement occurs during measurement of the upper half of the measurement space, if the lowest line in the body motion generation line is j, the data (ky ₎
It is considered that <j) does not include the influence of the motion, and since it includes half of all the measured data, the complex conjugate point symmetry is used from this data to include the influence of the remaining motion. The data (k _y ≧ j) can be estimated and a motion compensated image is obtained. However, since the part including the influence of motion is not used from the original measurement data, S /
N deteriorates. Therefore, in order to obtain the motion-corrected image while suppressing the deterioration of the S / N, the data (k _y
It is necessary to correct the motion for ≧ j).

From equations (1) and (2), the positional shift amount p in the reading direction and the phase encoding direction when k _y ≧ j
(k _y), if the estimated q (k _y) is, by performing the phase correction processing, it is possible to correct the translation. Therefore,
An object of the present invention, such _{p (k y), q (} k y) (k y ≧
j) is to be detected.

[0029]

[Means for Solving the Problems] In order to achieve the above object, the following procedure is performed.

[I] As shown in FIG. 3, it is divided into a low-frequency center 302 and high-frequency regions 301 and 303. And the low frequency center 302
The steps are divided into the step of correcting the motion of and the step of correcting the motion of the high range 301.

[1] Correction of motion in low-frequency center The parallel movement generated in the low-frequency center 302 is performed using a known technique. For example, the motion correction of the low frequency center 302 can be performed by the following method. First of all, when measuring
Only when measuring a low center 302, to detect motion in the phase encoding direction, k _y axis 304 toward zero encoded signals (navigation echo) {NAV
_ky (y)} is measured.

The hybrid space data {h _ky (x)} is obtained by one-dimensional Fourier transforming the measurement data of each line k _{y in} the reading direction of the low frequency band 302. The edge E _ky is detected from | h _ky (x) | of each line k _y . Based on the zero encoded data of the edge, positional deviation amount ΔX in the direction reading for each line k _y (k _y) =
_Eky −E ₀ is calculated. _{p (k y) = ΔX (} k y), as q (k _y) = 0, the number
Phase correction processing is performed based on (1) and (2). Similarly, the navigator echo signals {NAV _ky (y)}, performs a procedure of the above-obtained positional shift amount ΔY in the phase encoding direction (k _y). As _{p (k y) = 0,} q (k y) = ΔY (k y), the number
Phase correction processing is performed on the measurement data based on (1) and (2).

[2] High-frequency motion correction The high-frequency 301 is estimated by the half-scan imaging method from the low-frequency center 302 corrected in [1] and the high-frequency 303 of the original measurement data. The phase difference between the high-frequency part 301 estimated in step 1 and the data of the high-frequency part 301 of the original measurement data is calculated. Smoothing processing is performed in order to remove the noise of the phase difference data obtained in. With respect to the smoothed phase difference data in the high frequency band 301 obtained in step 1, p of the equation (2) is applied to each line in the reading direction.
_{(k y), r (k} y) = q (k y) k y a least-squares estimation. In estimated p (k _y), the q (k _y), the number (1), (2)
Based on this, the phase correction processing is performed on the data corrected in [1].

[II] In [2] of [I], the least squares method is applied to the phase difference data, but since the phase cyclically takes the same value every 2π, it is 1 in the reading direction.
For movements of more than one pixel, the number (2) at each k _y is folded back and does not become a straight line, as shown in FIG. Therefore, in order to apply the least-squares method, it is necessary to perform the unwrap processing that eliminates aliasing and restores a straight line. However, it is difficult to unwrap the actual measurement data due to the influence of noise. Therefore, in each line in the reading direction, a section that does not include a turnaround point is determined from the low frequency central portion [-c, c] of each line as follows, and the least squares method is applied in this section.

Moving dispersion is applied to the phase difference data of each line in the reading direction. Here, the moving dispersion is to obtain the dispersion at each point k _x of the line based on the phase difference data between that point and points b before and after. Threshold value 6 in the section [-c, c] of each line
Points a ₁ <a ₂ <... <a _k having the above dispersion value are _obtained . _{-C, a 1, a 2,} ..., a k, seeking the length of the adjacent two points to each other in the c, the largest interval between application section of the minimum square method.

[0036]

[Action]

[1] Motion correction in the center of the low range The motion correction in the center of the low range is a known technique, and as described above, in the −32 to +32 line in the center of the measurement space where the S / N is high, no navigation echo is used. By examining the hybrid space data for each line in the reading direction, it is possible to accurately detect the positional deviation amount in the reading direction between each measurement. Further, movement correction in the phase encoding direction is performed by using a navigation echo, but since the navigation echo is measured only when measuring the low-frequency center 302, the increase in measurement time can be small.

[2] High-frequency motion correction Since parallel movement occurs only in the upper half of the measurement space, it is considered that the measurement data in the lower half of the measurement space does not include the influence of motion. In the ideal case where the image takes real values, the measurement data are in a complex conjugate relationship with respect to the point-symmetrical positions with respect to the origin, so by using the measurement data of the lower half that does not affect the movement, We can estimate the upper half of the data, which does not include the impact. However, in reality, due to blood flow, device distortion, etc., phase distortion occurs in the image and a complex component appears, so the symmetry of the measurement data does not hold,
Unable to estimate the upper half data with good accuracy. The known half-scan imaging method approximates the phase distortion of an image with a phase map obtained from the low-frequency center of the measurement data, and corrects this to restore the symmetry of the areas 301 and 302 of the measurement data. Is.

Therefore, by using the half scan imaging method with the data of the low-frequency center 102 and the data of the high-frequency region 303 that have no motion influence in [1], it is possible to obtain the accuracy without the influence of motion. Good high range 301
The data can be estimated. Here, taking the phase difference between the data of the high band 301 estimated by the high band 301 and half-scan imaging method of the original measurement data, by applying the least squares method in the phase difference data in the out direction each line k _y read , The phase error of the equation (2) can be estimated, and the parallel movement can be corrected by performing the phase correction processing. As described above, since the obtained phase-corrected image uses all the measurement data, the S / N is superior to the motion-corrected image obtained by the half scan imaging method.

As described above, when applying the least squares method,
Consider limiting the application of the least-squares method in order to avoid the influence of phase folding. If there is no movement in the phase encoding direction of the inspection object and there is no offset in the phase, the applicable section of the least squares method is read out constantly for each line in the direction, and moreover, it is taken in the line center low range. Good. However, when there is movement in the phase encoding direction, as can be seen from equation (2), the phase offset becomes larger or smaller as the line goes to higher frequencies, so that the folding point as shown in FIG. 501 comes in the center. For this reason, there are 5 turnaround points for each line.
It is necessary to detect 01 and obtain a section that does not include the turning point, and apply the least squares method in this section.

In order to detect the turnaround point, first, the movement dispersion defined above is applied to the phase difference data for each line. The phase value near the turning point is because the amplitude from −π to π due to the turning is added to the vibration due to the influence of noise.
The dispersion value near the turnaround point (value close to 2π) is larger than the dispersion value at other points. Therefore, an appropriate threshold value is determined, and a point having a variance value equal to or larger than this threshold value is a turning point candidate. The high-frequency part of each line is strongly influenced by noise, so it is excluded from the application section of the least-squares method. Therefore, the section [-c, c] at the center of the line is set and the section to which the least square method is applied is determined from within this section. In this section, the points at which the variance is greater than or equal to the threshold are a ₁ , a ₂ , ...,
If a _l is set, such a point is localized near the turning point, so that the distance between two adjacent points in this vicinity is short.
Therefore, if you choose _{-c <a 1 <a 2 <} ... <a l <c adjacent two points to each other largest section seeking distance of the turning point in this segment are not included. Therefore, by setting this section as the section to which the least squares method is applied, it is possible to avoid the influence of the turning point when the least squares method is applied.

[0041]

【Example】

<First Embodiment> FIG. 2 shows an example of an MRI system to which the present invention can be applied. The static magnetic field generation system 201 generates a uniform static magnetic field, and the transmission system 202 generates a high frequency pulse magnetic field for exciting spins. The gradient magnetic field generation system 203
It is possible to generate each gradient magnetic field in the x, y, and z directions and change the magnetic field in each direction. The reception system 204 is the subject 200.
The electromagnetic wave radiated from is received, detected, converted into a digital signal, and then supplied to the processing device 205. The processing device 205 performs various calculations on the data from the reception system 204 to generate image data, and displays the image on the CRT display device 206. Pulse sequence file 207
Holds control information that regulates the operation sequence of the systems 201 to 204, and the sequence control unit 208 reads the control information from the pulse sequence file 207.
In accordance therewith, the operations of total 201 to 204 are controlled.

FIG. 13 shows a pulse sequence in the system of FIG. 2, which is in accordance with the control information from the pulse sequence file 207.
It is generated under the control of 8. RF 1301 indicates the timing of the high frequency pulsed magnetic field generated by the transmission system 202, and Gz 1302, Gy 1303, Gx 1304.
Indicates the timing at which the gradient magnetic fields in the z, y, and x directions are generated by the gradient magnetic field generation system 203, respectively.
Reference numeral 1305 indicates the timing at which the reception system 204 measures the measurement data signal 1307. RF1301 frequency and Gz1
Select the slice plane in the z direction at the level of 302, and
Position separation in the y direction is performed at the level of 303, and Gx130
Position separation in the x direction is performed at the level of 4. Since position separation in the y direction cannot be performed at one time, the level of Gy1303 is changed and the measurement data signal 130 is changed for each level.
Measure 7. Gy1303 typically varies over 256 levels, so 25
Six measurements are taken. In this way, the gradient magnetic field whose level is changed during scanning (in this example, the y-direction gradient magnetic field)
Is the phase encode pulse 1306. In the pulse sequence of FIG. 13, data measured by sequentially increasing or decreasing the phase encode level is S (k _x , k _y ), and k _x = −N _x / 2 to N _x / 2-1,
Let k _y = −N _y / 2 to N _y / 2-1.

FIG. 1 shows the flow of the motion correction process when the inspection object is translated in the slice plane in either the upper half or the lower half of the measurement space.
In step 101, the measurement space is set to the low center 302 and the high range 3
It is divided into 01 and 302. Here, the low center 302 is set to -N
~ + N lines. Normally, N is usually set to 8, 6, 32. Next, in step 102, the motion of the center of the low frequency band 302 is corrected, and then in the step 103, the motion of the high frequency band is corrected. From the measurement data itself, whether the generated motion occurred in the upper half of the measurement data, Since it is not possible to determine whether the movement occurred in the half, the movement correction processing as if the movement occurred in the upper half of the measurement space and the movement correction processing as that in the lower half are performed in parallel. Therefore, in step 104, the two data that have been subjected to the motion correction processing are Fourier transformed to obtain
Two reproduced images are obtained. First, the low-frequency motion correction step 102 will be described with respect to three examples of different motion detection methods.

[Low-Band Center Motion Correction Processing Example 1] A low-band center motion correction processing example 1 will be described with reference to FIG. In this example, first the edge of the hybrid space data obtained by one-dimensional Fourier transforming the measurement data in the center of the low frequency band in the reading direction is detected, and the position shift amount in the reading direction is detected to detect the movement in the reading direction. Make a correction. By this processing, the slope of the phase component of the measurement data due to the movement in the reading direction becomes flat, and the change in the phase offset due to the movement in the phase encoding direction for each line is calculated from the data that has been subjected to the movement correction processing in the reading direction. It can be detected by examining the average of the phase for each line. That is, the positional shift amount in the phase encoding direction can be detected by detecting the average change amount of the phase of each line with reference to the average phase of the zero encode line. In order for this method to be effective, it is necessary to remove in advance the linear phase distortion that is uniformly on the phase component of the measurement data due to the apparatus distortion.

[0045] for the data to remove the linear phase distortion due to [Step 601] device distortion, each line k _y belonging to low center 302 (k _y = -N~N) by Fourier transform one-dimensional, hybrid Obtain space data {h _ky (x)}.

[0046] [step 602] _{{| h ky (x) |} } by detecting edges from obtains the position deviation amount [Delta] X (k _y) in each _{k y (k y = -N~N)} .

[Step 603] p (k _y ) = ΔX (k _y ), q (k _y ) = 0 k _y ε [−N, N] where Interval [1], (2) -N, measured in N] data S (k _x, k _y) data subjected to phase correction processing in _{S '(k x, k y} ) is determined.

[Step 604] For each line k _y belonging to the low-pass center 302, the phase component arg (S ′ (k _x , k _y )) of the data S ′ (k _x , k _y ) obtained in step 603 is calculated. to calculate the average w (k _y).

[0049] For each line k _y belonging to [Step 605] low center 302, the positional shift amount of the phase encode direction _{ΔY (k y) = (w} (k y) -w (0)) / k y Calculate as.

[Step 605] p (k _y ) = 0, q (k _y ) = ΔY (k _y ) k _y ε [−N, N] is set based on the equations (1) and (2). -N, N]
Data S ′ obtained in step 603 for phase correction processing
Apply to (k _x , k _y ).

[Low-Band Center Motion Correction Processing Example 2] A low-frequency center motion correction processing example 2 will be described with reference to FIG. In this example, the navigation echo is used to detect the movement in the phase encoding direction in the low center 302. Use pulse sequence is phase encoding level indicated by viewing Figure 14 when -N～N, in addition to the low-center data necessary for image reproduction, to measure the k _y axis 304 toward zero encoded signals (navigation echo). Here, a navigation echo signal _{{NAV k (k y)}} (k = -N~
N). Regarding the movement in the reading direction, by detecting the edge from the hybrid space data obtained by one-dimensional Fourier transforming the measurement data of each line in the reading direction at the center of the low range, the movement in the phase encoding direction Detects an edge from hybrid space data obtained by one-dimensional Fourier transform of a navigation echo, thereby detecting the amount of positional deviation.

[Step 701] The measurement data of each line k _y belonging to the center 302 of the low frequency band is subjected to one-dimensional Fourier transform to obtain {h _ky (x)}.

[0053] [step 702] _{{| h ky (x) |} } by detecting edges from obtains the position deviation amount [Delta] X (k _y) in each _{k y (k y = -N~N)} .

[Step 703] Each k (k = -N to N)
On the other hand, NAV _k (ky) is one-dimensional Fourier transformed to obtain hybrid space data {HNAV _k (y)}.

[Step 704] {| HNAV _k (y)
|} Detects an edge from, obtain each k _y (k _y = -N~N) at position shift amount ΔY (k _y).

[0056] As [step 705] _{p (k y) = ΔX (} k y), q (k y) = ΔY (k y) k y ∈ [-N, N], the number (1), (2) based interval [-N, N] in the measurement data S (k _x, k _y) performs a phase correction process.

[Low-Band Center Motion Correction Processing Example 3] In this example, the case of multi-slice is considered. FIG. 15 is an example of a multi-slice spin echo sequence. In this example, the phase encoding direction and the reading direction during the second slice imaging are replaced with the phase encoding direction and the reading direction during the other slice imaging. As a result, with respect to the parallel movement that occurs during the echo measurement of the low-frequency center of each slice, the movement in the k _x axis 305 direction in the measurement space of FIG. 3 is obtained from the hybrid space data obtained from the echo of the first slice, and k _y With respect to the movement in the direction of the axis 304, the amount of positional deviation can be detected from the hybrid space data obtained from the echo of the second slice. The spin echo signals obtained from the slice _{S1 (k x, k y)} , S2 (k x, k y), ..., Sn (k x, k y) and. Hereinafter, a low-frequency motion correction processing example 3 will be described with reference to FIG.

[Step 801] Each k _{y of} the first slice
= -N～N echo by the measurement data S1 (k _x, k _y) is Fourier transformed one dimension to obtain {h1 _ky (x)}.

[Step 802] An edge is detected from {| h1 _ky (x) |}, and k in each k _y (k _y = −N to N) is detected.
_x-axis 305 determine the direction of the positional deviation amount ΔX (k _y).

[Step 803] Each k _{x of} the second slice
= -N～N echo by the measurement data S2 (k _x, k _y) is Fourier transformed one dimension to obtain {h2 _kx (y)}.

[Step 804] An edge is detected from {| h2 _kx (y) |}, and k in each k _y (k _y = −N to N) is detected.
_y-axis 304 determine the direction of the positional deviation amount [Delta] Y (k _y).

[0062] As [step 805] _{p (k y) = ΔX (} k y), q (k y) = ΔY (k y) k y ∈ [N, N], the number (1), based on (2) phase correction processing section [-N, N] Te and _{S1 (k x, k y)} , S2 (k x, k y), ...,
Apply to Sn (k _x , k _y ).

[0063] Hereinafter, data subjected to motion compensation processing of the _{low-S1 (k x, k y)} , S2 (k x, k y), ..., Sn (k x, k y) , respectively, the movement of the high frequency Perform correction processing.

Next, the high frequency motion correction step 103 will be described. Detection of high-frequency movement is performed in high frequencies 301 and 30.
It is utilized that the measurement data of 3 are in a complex conjugate point symmetry relationship with each other. First, the data in the high range 301 is converted to the high range 30.
Estimated from the measurement data of No. 3, and similarly, the data of the high frequency region 303 is estimated from the measurement data of the high frequency region 301. Then, the phase difference between the estimated high frequency data and the high frequency measurement data is calculated, and the phase error due to the parallel movement is obtained by applying the least square method to the phase difference data. Here, when the least squares method is applied to the phase difference data, noise removal processing is performed on the phase difference data in order to improve estimation accuracy. However, different methods are used in the following two high-frequency motion correction processing examples. Apply.

[High Frequency Motion Correction Processing Example 1] In high frequency motion correction processing example 1, noise removal processing of phase difference data
The moving average processing is performed on the phase difference data of each line.
Hereinafter, a high-frequency motion correction processing example 1 will be described with reference to FIG.

[Step 901] Estimated data IU (k _x of the high frequency band 301 based on the data of the low frequency band 302 and the high frequency band 303, which are motion-corrected in the low frequency band correction step 102, by a known half-can imaging method. , K _y ) (k _y
≧ N + 1) is calculated. Similarly, the estimated data I of the high band 303 is calculated from the data of the low band center 302 and the data of the high band 301.
D (k _x, k _y) determining the (k _y ≦ -N-1).

[Step 902] High band 3 of S (k _x , k _y ).
01 and the estimated data IU (k _x , k _y ) are calculated.

[0068] _{PU (k x, k y)} = arctan (S (k x, k y) * conj (IU (k x, k y))) (k y ≧ N + 1) in the same manner, S (k _x, high region 303 of k _y ) and estimated data I
The phase difference with D (k _x , k _y ) is taken.

[0069] _{PD (k x, k y)} = arctan (S (k x, k y) * conj ((ID (k x, k y))) (k y ≦ -N -1) where, conj ( ) Represents a complex conjugate.

[Step 903] Each line k _y (k _y ≧ N
_{+ 1, k y ≦ -N-} 1) in the phase difference data P _ky (k _{x)
= PU (k x, k y} ) or PD (k _x, k _y) with respect to, M
PU (k _x, k _y) or _{MPD (k x, k y)} = (P ky (k x -b) + P ky (k x -b + 1) + ... + P ky (k x) + ... + P ky (k x as _{+ b-1) + P ky} (k x + b)) / (2b + 1) (k x = -N x + b~N x -b), performs smoothing the phase difference data.

[Step 904] The process goes to the step of determining the least squares method applicable section described below, and each k of the high frequency band 301
In _{y (k y ≧ N + 1} ), MPU (k x, k y) least square method applied section Request _{[su (k y), eu (} k y) ]. Similarly, in each k _y (k _y ≦ −N−1) of the high frequency band 303, the least square method applied section [sd of MPD (k _x , k _y ) is applied.
(k _y ), ed (k _y )] is obtained.

[Step 905] Each k _y (k
In _y ≧ N + 1), the interval _{[su (k y), eu (} k y) ]
Applying the least squares method to MPU (k _x , k _y ) in
Of p (k _y), to estimate _{r (k y) = q (} k y) k y (k y ≧ N + 1).

Similarly, each k _y (k _y ≤−) of the high frequency band 303
N-1 in), M in the interval _{[sd (k y), ed (} k y) ]
The least squares method is applied to PD (k _x , k _y ) and p of the equation (2) is applied.
(k _y), to estimate _{r (k y) = q (} k y) k y (k y ≦ -N-1).

[0074] p (k _y) estimated in [Step 906] Step 905, the _{q (k y) (k y} ≧ N + 1), number (1),
Based on (2), if k _y ≧ N + 1, then S (k _x , k _y )
Then, the data SU (k _x , k _y ) subjected to the phase correction processing is obtained.

Similarly, p estimated in step 905
(k _y), q (k _y) by (k _y ≦ -N-1), number (1), (2)
Obtained in basis _{k y ≦ -N-1, S} (k x, k y) to the phase correction processing performed data SD (k _x, k _y) to.

[Least Square Method Applied Section Determining Step] The least square method applied section determining step will be described with reference to FIG.

[Step 1001] As input, phase difference data PU (k _x , k _y ) or P on each line k _y
Receive D (k _x , k _y ).

[0078] Phase difference data Pk _y (k _x) in [Step 1002] Each line _{k y = PU (k x,} k y) or P
For D (k _x , k _y ), the moving variance σ (k _x ) (k _x = −N _x
+ B~N _x -b) take.

[0079] _{μ (k x) = (Pk} y (k x -b) + Pk y (k x -b + 1) + ... + Pk y (k x) + ... + Pk y (k x + b-1) + Pk y (k x + b)) / (2b + 1) σ (k x) = ((Pk y (k x -b) -μ (k x)) 2 + (Pk y (k x -b + 1) -μ (k x))
² _{+ ... + (Pk y (k} x) -μ (k x)) 2 + ... + (Pk y (k x + b-1) -μ (k x)) 2 + (Pk y (k x + b) -μ (k _x )) ² ) / (2b + 1) [Step 1003] The moving variance value σ (k _x ) is [−c,
c], k _x exceeding the threshold value 6 is obtained, and −c <
Let a ₁ <a ₂ <... <a _l <c.

[Step 1004] length (0) = a ₁ +
c, length (i) = _{ai + 1- ai} (i = 1 to l-1), length
(1) = c−a _l . At this time, length (i) (i = 1 to
The value of i that maximizes l) is calculated and set as imax.

[Step 1005] As an output, [a
_imax , a _{imax + 1} ] is returned. Where imax = 0, imax + 1
= L + 1, a _imax = -c and a _{imax + 1} = c.

[High-Band Motion Correction Processing Example 2] In the high-band motion correction processing example 2, a method of averaging the phase difference data in the phase encoding direction is adopted in the noise removal processing.
More specifically, in a certain line, a block consisting of this line and g-1 lines on the high frequency side is considered, and by averaging the phase difference data in the phase encoding direction from the phase difference data in this block, Smooth the phase difference data. Then, the phase difference data for each line in which the smoothing process in this way, by applying the least square method, the number p (k _y) of (2), but to estimate q (k _y),
P in the motion generation line and g lines on the low frequency side
The value of (k _y) or q (k _y) increases or decreases monotonically. This is because the contribution to the average of the front and the rear of the motion generation line in the block changes monotonically as the block described above crosses the motion generation line from the low frequency side to the high frequency side. The sudden change in the slope or offset of the phase difference data of each line appears as a monotonous change in the slope or offset of the phase difference data of each line after smoothing. This monotonous change continuously occurs g times the width of the block.

Therefore, g- on the low frequency side from the motion generation line
P according to the above method in one line (k _y), q (k
Since estimates of _y) is not correct, p (k _y in g-th line from the motion generating line low frequency side), it embeds the value of q (k _y). If the two motion generation lines are not separated by more than g lines, p between the two generation lines
(k _y), the value of q (k _y) can not be estimated by the above method, the method of motion compensation processing example 1 of the high frequency, the lines each p (k _y), the value of q (k _y) To ask. Finally,
If the motion is occurring in the -N-g to -N-1 lines and the N + 1 to N + g lines near the boundary between the center of the low frequency range and the high frequency range, the generation line cannot be specified, so the above processing cannot be performed. . Therefore, the low-frequency motion correction step 10
In 2, by performing the motion correction in a region in addition to the low center 302 by the vertical g-line _{(-N-g ≦ k y ≦} N + g), to avoid this.

Hereinafter, a high-frequency motion correction processing example 2 will be described with reference to FIG.

[Steps 1101 to 1102] Same as steps 901 to 902 in the high frequency motion correction processing example 1.

[Step 1103] In line k _y ((N _y / 2-1) -g + 1 ≧ k _y ≧ N + 1) of high band 301, EPU (k _x , k _y ) = PD (k _x , k _y ). _{+ PD (k x, k y} +1) + ... + PD (k x, k y + g-1) by a / g, smoothing the phase difference data.

Similarly, the line k _y (−N) of the high frequency band 303
_{y / 2 + g-1 ≦} k y ≦ -N-1 _{in), EPD (k x, k} y) = PU (k x, k y) + PU (k x, k y -1) + ... + PU (k x, with _{k y -g + 1) / g} , smoothing the phase difference data.

[Step 1104] Go to the least squares method applicable section determination step in the high frequency motion correction processing example 1, and set the high frequency range 30
The least-squares-applied section [su of EPU (k _x , k _y ) in the line k _y ((N _y / 2-1) -g + 1 ≧ k _y ≧ N + 1) of 1 [su
(k _y ), eu (k _y )] is obtained.

Similarly, the line k _y (-N of the high frequency band 303)
EPD (k _x , in _y / 2 + g−1 ≦ k _y ≦ −N−1)
k _y) least square method applied section Request _{[sd (k y), ed (} k y) ].

[Step 1105] In the line k _y ((N _y / 2-1) -g + 1 ≧ k _y ≧ N + 1) of the high frequency band 301,
Section applying the least square method in _{[su (k y), eu (} k y) ] EPU (k _x, k _y) in, p (k _y) of the number _{(2), r (k y} ) = q
(k _y) estimates the _{k y ((N y / 2-1} ) -g + 1 ≧ k y ≧ N + 1).

Similarly, each line k _y (-
N _y / 2 + g−1 ≦ k _y ≦ −N−1), the interval [s
d (k _y), applying the least square method by ed (k _y)] EPD (k _x, k _y) in, p (k _y) of the number _{(2), r (k y} ) = q (k y ) k _y (-
Estimate N _y / 2 + g−1 ≦ k _y ≦ −N−1).

[Step 1106] In the partial region of the high band 301 ((N _y / 2-1) -g + 1 ≧ k _y ≧ N + 1),
_{p (k y), q (} k y) is N + ₁ ≦ j 1 the start line to increase or decrease continuously g times in the phase encoding positive direction <
j ₂ <... <j _n <N _y , similarly, in the partial region (−N _y / 2 + g−1 ≦ k _y ≦ −N−1) of the high frequency region 303, p
(N _y <j ₁ ′ <j is a starting line where (k _y ), q (k _y ) continuously increases or decreases g times in the negative direction of phase encoding.
_{2, '<... <j n' ≦} -N-1.

[Step 1107] Each section [j _k , j _{k + 1}
−1] (k = 1 to n) (where j _{n + 1} = N _y ), if j _{k + 1} −j _k ≧ g, then [j _k , j _{k + 1} −g] in, p estimated in step _{1105 (k y), q (} k y)
Accordingly, in the _{[j k + 1 -g + 1,} j k + 1 -1 ], p
(j _{k +} 1 -g), the _{q (j k + 1 -g)} , obtaining S (k _x, k _y) data was performed phase correction process to _{SU (k x, k y)} .

If j _{k + 1} −j _k <g, then [j _k , j _{k + 1}
In -1] obtained by performing the step 905 p (k _y), q by _{(k y), S (k} x, k y) data was performed phase correction process to SU (k _x, k _y) To get

Similarly, each section [j _k ′ −1,
j _{k + 1} ′] (k = _{0 to} n−1) (where j ₀ = N _y ), if j _{k + 1} ′ −j _k ′ ≧ g, then [j _k ′ +
g, 'in], p estimated in step 1105 (k _y), the q (k _{y), [j k' j k + 1 +} g-1,
j _k ′], p (j _{k + 1} ′ + g), q (j _{k + 1} ′ +
The g), to obtain the S (k _x, k _y) data was performed phase correction process to _{SD (k x, k y)} .

If j _{k + 1} ′ −j _k ′ <g, then
_{[J k ', j k + 1} ' -1 ] in obtained by executing steps 905 p (k _y), the q (k _y), S
Data SD (k _x , k) obtained by performing phase correction processing on (k _x , k _y ).
_y ).

<Embodiment 2> In Embodiment 2, the case of high-speed spin echo imaging will be considered. FIG. 16 (a) shows that in the system of FIG. 2, when the 512 × 256 matrix measurement space is divided into five blocks as shown in FIG. It is an example of a fast spin echo sequence in which each block is measured by assigning one echo. In this case, the imaging time is about ⅕ of that of the conventional spin echo imaging method, so that the inspection object does not move much during the imaging, and even if a sudden movement of the inspection object occurs, the measurement can be performed. The situation where it happens in either the first half or the second half becomes even more realistic.

In the first embodiment, the half scan imaging method was used to detect the motion generated in the high frequency region. For this reason, it was necessary to correct the motion in the low frequency central portion of the measurement data. In the example, a multi-echo sequence having a number of echoes larger than the number of blocks is used, and two multi-echoes are assigned to the block 403 corresponding to the central part of the low frequency range, and the measurement is performed in either the first half or the second half of the measurement. A high-speed spin echo imaging method capable of obtaining data in the central portion of the low frequency band that is not affected by the motion even if the motion occurs will be described. After that, a motion correction method from data obtained by the high speed spin echo imaging method will be described.

FIG. 16B is a high-speed spin echo sequence using a multi-echo sequence of 6 echoes, in which the first echo and the second echo are respectively blocked 401.
And block 402, and the fifth and sixth echoes are assigned to blocks 404 and 405, respectively.
Also, the third echo and the fourth echo are assigned to the block 403. For the third echo, the block 403 is measured from top to bottom, and for the fourth echo, the block 403 is measured from bottom to top.
If the block 403 is measured in this way, and if the movement occurs only in either the first half or the second half of the measurement, the upper half of the block 403 measured by the third echo,
The lower half of the block 403 measured by the fourth echo,
Alternatively, the lower half of the block 403 measured by the third echo and the block 403 measured by the fourth echo
The data in block 403 from either combination in the top half does not include motion effects.

The motion detecting method in the data measured by the high speed spin echo imaging method will be described below. Motion detection is performed using the complex conjugate point symmetry of the measurement data of blocks 402 and 404. If the motion can be detected in one block, the motion correction process may be performed from the detected motion in the other blocks as well. If the measurement data measured by the i-th echo (i = 1, 2, ..., 6) is S _i (k _x , k _y ) (i = 1, 2, ..., 6), the high-speed spin as described above is performed. Since it was measured by the echo imaging method, there is the following correspondence.

[0102] Measurement data S ₁ block _{401 (k x, k y)} 77 ≦ k y measurement data S ₂ of ≦ 127 block _{402 (k x, k y)} 26 ≦ k measurement data S ₃ of _y ≦ 76 block 403 (k _x, k _y), the measurement of _{S 4 (k x, k y} ) -25 ≦ k y ≦ 25 block 404 measurement data _{S 5 (k x, k y} ) -76 ≦ k y ≦ -26 block 405 data _{S 6 (k x, k y} ) -127 ≦ k y ≦ -78 Further, each block other than the block 403, is measured from top to bottom by the assigned echo. With the above preparations, the motion detection method will be described with reference to FIG.

[Step 1201] The upper half data S ₃ (k _x , k _y ) (0 ≦ k _y ≦ 25) and the lower half data S ₄ (k _x , k _y ) (−1 ≦ k) of the block 403. a central low-frequency measurement data by _y ≦ -25), half of the data S ₅ (k _x on the block 404, k _y) (- from 26 ≦ k _y ≦ -51), the half scan imaging method, block 402 estimation data ID ₂ under half (k _x, k _y) to calculate the _{(26 ≦ k y ≦ 51)} .

Similarly, the upper half data S ₄ (k _x , k _y ) (0 ≦ k _y ≦ 25) of the block 403 and the lower half data S
₃ (k _x , k _y ) (− 1 ≦ k _y ≦ −25) in the central low frequency band, and the lower half of the block 404, S ₅ (k _x ,
From k _y ) (− 76 ≦ k _y ≦ −51), the estimated data IU ₂ (k _x , k _y ) (51 ≦ k _y ≦ 76) of the upper half of the block 402 is calculated by the half scan imaging method.

[Step 1202] The phase difference between the upper half data S ₂ (k _x , k _y ) (51 ≦ k _y ≦ 76) of the block 402 and the estimated data IU ₂ (k _x , k _y ) is calculated.

PU ₂ (k _x , k _y ) = arctan (S ₂ (k _x , k _y ) * conj (IU ₂ (k _x , k _y ))) (51 ≦ k _y ≦ 76) Similarly, Data S in the lower half of block 402
₂ (k _x , k _y ) (26 ≦ k _y ≦ 51) and estimated data ID ₂ (k
_x, take the phase difference between the k _y).

PD ₂ (k _x , k _y ) = arctan (S ₂ (k _x , k _y ) * conj ((ID ₂ (k _x , k _y ))) (26 ≦ k _y ≦ 51) where: conj () represents a complex conjugate.

[Step 1203] Each line k _y (26 ≦
k _y ≦ 76 phase difference) data PU ₂ (k _x, smooth, k _y) and _{PD 2 (k x, k y} ) The method of step 903 in the high range of the motion compensation processing example 1 of Example 1 with respect to Conversion to MPU
₂ (k _x , k _y ) (51 ≦ k _y ≦ 76) and MPD ₂ (k _x , k _y ).
Calculate (26 ≦ k _y ≦ 51).

[Step 1204] Go to the least-squares method applicable section determination step of the high-frequency motion correction processing example 1 of the first embodiment, and in each k _y (51 ≦ k _y ≦ 76), MPU ₂ (k _x ,
k _y) least square method applied section Request _{[su (k y), eu (} k y) ].

Similarly, in each k _y (26 ≦ k _y ≦ 51), the least-squares-applied interval [s of MPD ₂ (k _x , k _y ) [s]
d (k _y), determine the ed (k _y)].

[Step 1205] Each k _y (51 ≦ k _y ≦
In 76), MP in the interval _{[su (k y) 0eu (k} y) ]
The least squares method is applied to U ₂ (k _x , k _y ) to obtain p ₂ (k
_y), estimating the _{r 2 (k y) = q} 2 (k y) k y (51 ≦ k y ≦ 76).

Similarly, for each k _y (26 ≦ k _y ≦ 51), MPD ₂ (k _x ,), in the interval [sd (k _y ), ed (k _y )].
k _y) applying the least squares method to, p ₂ (k _y of number (2)), r
_{2 (k y) = q 2 (} k y) k y (26 ≦ k y ≦ 51) is estimated.

[Step 1206] Based on p ₂ (k _y ), q ₂ (k _y ) (51 ≦ k _y ≦ 76) estimated in step 1205, the upper half of the measured data for each block is added to the upper half of the block 401. S ₁ (k _x , k _y ) 102 ≦ k _y ≦ 127 p (k _y ) = p
_{2 (k y -51), q} (k y) = upper half of q ₂ (k _y -51) block _{402 S 2 (k x, k} y) 51 ≦ k y ≦ 76 p (k y) = p
_{2 (k y), q (} k y) = upper half of q ₂ (k _y) blocks _{403 S 3 (k x, k} y) 0 ≦ k y ≦ 25 p (k y) = p
_{2 (k y +51), q} (k y) = upper half of q ₂ (k _y +51) block _{404 S 5 (k x, k} y) -51 ≦ k y ≦ -26 p (k y) = p
_{2 (k y + 2 × 51} ), q (k y) = q 2 (k y + 2 × 51) upper half S ₆ block _{405 (k x, k y)} -102 ≦ k y ≦ -77 p (k y ) = P
_{2 (k y + 3 × 51} ), by the _{q (k y) = q 2} (k y + 3 × 51), after the phase correction processing, the entire measurement space obtained by collecting these data The data SU (k _x , k _y ) is obtained.

Similarly, based on p ₂ (k _y ), q ₂ (k _y ) (26 ≦ k _y ≦ 51) estimated in step 1205, the lower half of the measured data for each block is added to the lower half of the block 401. _{S 1 (k x, k y} ) 77 ≦ k y ≦ 102 p (k y) = p
_{2 (k y -51), q} (k y) = q 2 (k y -51) the lower half of the block _{402 S 2 (k x, k} y) 26 ≦ k y ≦ 51 p (k y) = p
_{2 (k y), q (} k y) = q 2 (k y) the lower half of the block _{403 S 3 (k x, k} y) -25 ≦ k y ≦ 0 p (k y) = p
_{2 (k y +51), q} (k y) = q 2 (k y +51) the lower half of the block _{404 S 5 (k x, k} y) -76 ≦ k y ≦ -51 p (k y) = p
_{2 (k y + 2 × 51} ), q (k y) = q 2 (k y + 2 × 51) the lower half of the block _{405 S 6 (k x, k} y) -127 ≦ k y ≦ -102 p (k y ) = P
_{2 (k y + 3 × 51} ), by the _{q (k y) = q 2} (k y + 3 × 51), after the phase correction processing, the entire measurement space obtained by collecting these data Data SD (k _x , k _y ) is obtained.

[Step 1207] Phase correction data SU
(k _x, k _y), the image playback SD (k _x, k _y) the Fourier transform.

[0115]

According to the present invention, the parallel movement in the slice plane, which suddenly occurs only in the first half or the second half until all the data necessary for image reproduction is measured, does not require a special sequence. By detecting the motion from the measurement data itself and correcting the motion from the detected motion to the measurement data, it is possible to suppress the parallel motion motion artifact without degrading the S / N.

[Brief description of drawings]

FIG. 1 is a processing flowchart of motion correction processing according to a first embodiment of the present invention.

FIG. 2 is a diagram showing an example of a system configuration for implementing the present invention.

FIG. 3 is a diagram showing each region when the measurement space is divided into a low frequency center and a high frequency region.

FIG. 4 is a diagram showing allocation of a measurement space of each echo in high-speed spin echo imaging of 6 echoes.

FIG. 5 is a diagram illustrating folding back of phase difference data in each line of a measurement space.

FIG. 6 is a processing flowchart of an example of motion correction processing of the first center of the low frequency range according to the first embodiment of the present invention.

FIG. 7 is a processing flowchart of a second low-frequency center motion correction processing example according to the first embodiment of the present invention.

FIG. 8 is a processing flowchart of a third low-frequency center motion correction processing example according to the first embodiment of the present invention.

FIG. 9 is a processing flow chart of a first high-frequency motion correction processing example according to the first embodiment of the present invention.

FIG. 10 is a processing flow chart for determining an applicable section when the least squares method is applied to the phase difference data of each line in the measurement space in the first embodiment of the present invention.

FIG. 11 is a processing flow chart of a second high-frequency motion correction processing example according to the first embodiment of the present invention.

FIG. 12 is a processing flowchart of an example of motion correction processing according to the second embodiment of the present invention.

FIG. 13 is a diagram showing an example of a sequence for capturing measurement data for which motion correction processing is performed according to the first embodiment of the present invention.

FIG. 14 is a diagram showing an example of a sequence for photographing the low frequency center of measurement data for which motion correction processing is performed in the first embodiment of the present invention.

FIG. 15 is a diagram showing an example of a sequence for capturing measurement data for which motion correction processing is performed in the first embodiment of the present invention.

FIG. 16 is a diagram showing a conventional sequence and a sequence of the present invention in a high speed spin echo imaging method.

[Explanation of symbols]

Reference numeral 200 ... Subject, 201 ... Static magnetic field generation system, 202 ... Transmission system, 203 ... Gradient magnetic field generation system, 204 ... Reception system, 205
... processing device, 206 ... CRT display device, 207 ... pulse sequence file, 208 ... sequence control unit, 30
4 ... k _y axis, 305 ... k _x axis, 501 ... turn-around point, 1
301 ... RF pulse, 1302 ... z-direction gradient magnetic field, 13
03 ... y-direction gradient magnetic field, 1304 ... x-direction gradient magnetic field, 1
305 ... Signal, 1306 ... Phase encode pulse, 1
307 ... Measurement signal, 1401 ... Navigation echo.

Front page continuation (72) Inventor Hideaki Koizumi 882, Ige, Katsuta-shi, Ibaraki Hitachi Ltd. Measuring Instruments Division

Claims (9)

[Claims] 1. A magnetic resonance imaging apparatus comprising a static magnetic field, a gradient magnetic field, a high-frequency magnetic field generator, a detector for extracting a nuclear magnetic resonance signal from an inspection object, and a processor for performing various calculations including an image reproduction configuration. , Two or more different types of motion detection means are prepared, the measurement data is divided into two or more blocks, and each block detects the movement by any one type of motion detection means, and the movement of the measurement data is corrected. A motion correction method characterized by performing. 2. The motion correction method according to claim 1, wherein the measurement data is divided into low-frequency data in the central portion and high-frequency data other than that, and the high-frequency data is subjected to motion detection based on phase information of the measurement data. A characteristic motion compensation method. 3. The low-frequency center correction processing according to claim 2 includes processing for correcting movement in the line direction based on position information of data obtained by one-dimensional Fourier transform of measurement data of each line. A motion compensation method characterized by the above. 4. The low-frequency center correction processing according to claim 2, wherein the phase components of each line are averaged, and the phase offset is detected by the deviation amount of the average value of the phase of each line,
A motion correction method comprising a process of correcting motion in a direction perpendicular to a line (phase encoding direction). 5. A magnetic resonance imaging apparatus having a static magnetic field, a gradient magnetic field, a high frequency magnetic field generator, a detector for extracting a magnetic resonance signal from an object to be inspected, and a processor for performing various calculations including an image reproducing configuration, A motion correction method characterized by dividing measurement data into an upper half and a lower half, estimating a phase error due to movement of the upper half measurement data and the lower half measurement data, and performing phase correction processing based on the phase error. 6. The motion correction method according to claim 5, wherein the phase error is estimated by estimating half data on one side from half measurement data on the other side by utilizing complex conjugate point symmetry. A motion correction method characterized in that a phase difference between data and the estimated measurement data is calculated, and a phase error due to motion is estimated by applying a least square method to the phase difference data. 7. The method of applying the least squares method to phase difference data according to claim 6, wherein a block in which a few lines adjacent to the line are added to the line of interest is considered as a line within the block. A motion correction method comprising a process of smoothing phase difference data by averaging phase data in a vertical direction (phase encoding direction). 8. The process of applying the method of least squares to phase difference data according to claim 6, obtains a section which does not include a point where the phase difference values are −π and π in the central low frequency section of each line. And a method of applying a least squares method in this section. 9. A high-speed method for obtaining measurement data necessary for image reproduction by dividing a measurement space into a plurality of blocks, allocating spin echo signals measured at different echo times (TE) using a multi-echo sequence for each block. In the spin echo imaging method, the number of echoes is increased by one more than the number of blocks, and the block in the center of the measurement space including the line with zero phase encoding is measured with two multi-echoes, and each echo measures the block in the opposite direction. A high-speed spin echo imaging method characterized by: