Abstract:
A magnetic resonance signal is detected at a plurality of different periods of time in accordance with the phase shift of a plurality of chemical shift components while varying the intensity of a gradient magnetic field with time, a plurality of magnetic resonance signal thus obtained are used for estimating a magnetic resonance signal due to only a predetermined chemical shift component, the phase rotation of the predetermined chemical shift component is corrected in the estimated resonance signal, and further a blur due to chemical shift components other than the predetermined one is corrected to obtain a tomogram formed of only the predetermined chemical shift component. Alternatively, one measured resonance signal is modified by an image reconstruction algorithm peculiar to each of a plurality of chemical shift components, to obtain tomograms corresponding to the chemical shift components.

Description:
CROSS-REFERENCE TO RELATED APPLICATIONS 
     This application is a continuation-in-part of copending application Ser. No. 262,282, filed Oct. 25, 1988, now U.S. Pat. No. 4,901,021. 
    
    
     BACKGROUND OF THE INVENTION 
     The present invention relates to an apparatus for obtaining a tomogram of a living body by utilizing the nuclear magnetic resonance phenomenon, and more particularly to an image reconstructing device suited to separate chemical shift components from measured resonance signals and to form a tomogram of each chemical shift component, by a high-speed imaging method in which a resonance signal is measured in a state that the strength of a gradient magnetic field is varied with time. Incidentally, the image reconstruction means that a tomogram is reconstructed from a measured resonance signal. 
     A first method of forming a tomogram of each chemical shift component in the magnetic resonance imaging (hereinafter referred to as &#34;MRI&#34;) is described in, for example, an article entitled &#34;Chemical-Shift Imaging with Large Magnetic Field Inhomogeneity&#34; (Magnetic Resonance in Medicine, Vol. 4, No. 5, 1987, pages 452 to 460). In the first method, resonance signals are measured at some different periods of time in a state that the strength of a gradient magnetic field is kept constant, tomograms of an object are reconstructed from measured resonance signals by the two-dimensional Fourier transformation method, and tomograms corresponding to chemical shift components are separated from the reconstructed tomograms. In more detail, each of the chemical shift components is slightly different in resonance frequency from each other. Thus, when a resonance signal is measured at different periods of time, tomograms corresponding to the chemical shift components are different from each other in phase, that is, in the state of nuclear spin. Further, image processing is carried out so that tomograms corresponding to undesired chemical shift components cancel each other on the basis of the above phase difference, to obtain tomograms each formed of a single chemical shift component. 
     While, a novel imaging method is described in, for example, an article entitled &#34;A Novel Fast Scanning System&#34; (Proceedings of 5th Society of Magnetic Resonance in Medicine, 1986, pages 156 and 157). In this method, a resonance signal is measured in a state that the strength of a gradient magnetic field is varied with time. As a result, the resonance signal can be continuously measured, and thus a time necessary for imaging can be shortened in a great degree. 
     A general method of reconstructing a tomogram for a resonance signal which is obtained in the abovementioned manner, that is, the so-called correlation method is discussed on page 781 of the Proceedings of the 6th Society of Magnetic Resonance in Medicine held in 1987. 
     A second method of forming tomograms of a plurality of chemical shift components at the same time from a resonance signal which is measured only once, is discussed on page 230 of the Proceedings of the 6th Society of Magnetic resonance in Medicine held in 1987. The second method utilizes a fact that when the strength of two gradient magnetic fields in different directions is varied with time to obtain two-dimensional tomograms, the difference in resonance frequency between chemical shift components generates blurred tomograms. Let us pay attention to the atomic nucleus of hydrogen, by way of example. When a tomogram is reconstructed in accordance with the resonance frequency of the atomic nucleus of hydrogen contained in water, a resonance signal due to the atomic nucleus of hydrogen contained in fat is subjected to additional phase rotation, and thus a blurred tomogram is superposed on the tomogram due to water. When a tomogram is reconstructed in accordance with the resonance frequency of the atomic nucleus of hydrogen contained in fat, a tomogram due to water is blurred. That is, when a tomogram is reconstructed in accordance with the resonance frequency of the atomic nucleus of hydrogen contained in one of water and fat, a blurred tomogram due to the other substance is superposed on the tomogram due to one substance. When a blurred tomogram can be neglected, the tomogram due to water and the tomogram due to fat can be obtained from the single resonance signal. 
     The first method has given no consideration to the matter that tomograms due to chemical shift components may be reconstructed at high speed, by using a time-varying gradient magnetic field. The reason for this is as follows. 
     Now, let us express the tomogram formed of the i-th chemical shift component by ρ i  (x, y), where x and y designate two-dimensional coordinate values on the tomogram, with the origin at the center of a field of view. Further, let us express the deviation of a resonance frequency of each chemical shift component from a reference frequency by Δω i . It is to be noted that the reference frequency means the resonance frequency corresponding to a specific chemical shift component under a static magnetic field. Then, the resultant resonance signal S(t) of respective resonance signals of chemical shift components is given by the following equation: ##EQU1## where ##EQU2## and j=√1. Further, G x  and G y  indicate the strength of x- and y-gradient magnetic fields, respectively, γ: a nuclear gyromagnetic ratio, and C: a constant. 
     The inhomogeneity of static magnetic field and the T 2  effect (namely, transverse relaxation time effect) have no connection with the following explanation, and hence explanation thereof will be omitted. Now, let us consider a case where G x  and G Y   are constant as in the two-dimensional Fourier transformation method, for example, a case where G x  (t)=G x   o  (namely, a constant) and G y  (t)=0. The resonance signal S(t) is rewritten as follows: ##EQU3## Thus, tomograms formed of different chemical shift components are shifted from each other by an amount equal to Δω i  /(YG x   o ). As a result, an overlapping tomogram is formed. In a case where G x  and G y  vary with time, however, signal components due to different chemical shift components do not appear as signal components which are merely shifted from each other, but interact on each other in a complicated manner. Thus, it is impossible to reconstruct tomograms corresponding to chemical shift components at the same time. That is, the first chemical shift imaging method is inapplicable to a case where both G x  and G y  vary with time. 
     While, according to the second method, the point spread function (hereinafter referred to as &#34;PSF&#34;) at a focused tomogram is equal in average power (that is, average brightness) to the PSF at a tomogram which is out of focus, and hence it is necessary to remove the influence of the tomogram which is out of focus, on the focused image. In other words, the shape of the PSF is a function for indicating the degree of blur, and a tomogram due to a desired chemical shift component is equal in power to a tomogram due to an undesired chemical shift component. Hence, it is impossible to obtain focused tomograms of a plurality of chemical shift components at the same time. 
     Various methods of reconstructing tomograms from the resonance signal of the equation (1) have been known. The most general one of the above methods is the correlation method. Accordingly, the correlation method which is described in IEEE Transactions on Medical Imaging, Vol. 7, No. 1, pp. 26-31, March 1988, for example, will be used in the following explanation. It is needless to say that various methods other than the correlation method can be used for reconstructing tomograms. 
     Now, let us express a tomogram which is reconstructed in accordance with the k-th chemical shift component, by ρ k  (x, y). Then, ρ k  (x, y) is generally given by a complex function as follows: 
     
         ρ.sub.k (x, y)=∫dts(t)exp[+jγ{kx(t)x+ky(t)y}·exp[-jΔω.sub.k t]·W(t)                                     (3) 
    
     where W(t) indicates a weight coefficient for controlling a PSF. 
     When the equation (1) is substituted in the equation (3), we can obtain the following equation: ##EQU4## where a mathematical sign * indicates two-dimensional convolution which is given by the following equation: ##EQU5## Further, in the equation (4), h ki  (x, y) indicates a PSF which is given by the following equation: 
     
         h.sub.ki (x, y)=C∫dt exp[jγ{kx(t)x+ky(t)y]·exp[+j(Δω.sub.i -Δω.sub.k)t]·W(t)                    (6) 
    
     In more detail, h kk  indicates a focused PSF, and h ki  (k≠i) a PSF which is out of focus. The PSF h ki  can be obtained from the time-variation pattern of each gradient magnetic field and the weight coefficient W(t), and hence is considered to be a known function. 
     From the equation (6), we can obtain the following equation: 
     
         ∫dxdy|h.sub.kk (x, y)|.sup.2 =∫dxdy|h.sub.ki (x, y)|.sup.2 for all i (7) 
    
     As is evident from the equation (7), it is impossible to neglect the influence of chemical shift components other than the k-th chemical shift component on the tomogram ρ k  given by the equation (4). This is true of a case where only two chemical shift components exist as in the atomic nucleus of hydrogen. 
     SUMMARY OF THE INVENTION 
     It is accordingly a first object of the present invention to provide an image reconstructing device which can solve the above problems and can form a tomogram of each chemical shift component directly by the high-speed imaging method using time-varying gradient magnetic fields, with the aid of a tomogram reconstruction algorithm peculiar to each chemical shift component. 
     It is a second object of the present invention to provide an image reconstructing device, in which a chemical shift component included in the equation (4) and put out of focus is estimated and removed, to reconstruct tomograms of a plurality of chemical shift components from a resonance signal measured only once. 
     In order to attain the first object, according to an aspect of the present invention, a resonance signal due to a desired chemical shift component is produced from a plurality of resonance signals measured at different periods of time, and an image reconstruction algorithm peculiar to the desired chemical shift component is applied to the produced resonance signal, to form a tomogram of the desired chemical shift component. 
     A resonance signal S.sub.ΔT (t) obtained at a time which is deviated from a predetermined time by ΔT, is given by the following equation: ##EQU6## The final factor on the right-hand side of the equation (8) indicates that each chemical shift component contributing to the resonance signal S.sub.ΔT (t) has phase shift of exp[jΔω i  ΔT]. The resonance frequency ω i  of each chemical shift component is known, and hence the above phase shift is determined depending upon the time deviation ΔT. By measuring resonance signals for a plurality of values of ΔT and performing an appropriate arithmetic operation for the resonance signals, a resonance signal S i  (t) due to only the i-th chemical shift component can be obtained as follows: 
     
         S.sub.i (t)=C∫dxdyρ.sub.i (x, y)exp[-jγ{k.sub.x (t)x+k.sub.y (t)y}+jΔω.sub.i t]                            (9) 
    
     In a case where a tomogram of the i-th chemical shift component is reconstructed, a resonance signal S i  (t) which is obtained by removing the phase rotation due to a change in resonance frequency (that is, the phase shift from a reference frequency at a measuring time) from the equation (9), is expressed as follows: 
     
         S.sub.i (t)=S.sub.i (t)exp[-jΔω.sub.i t]       (10) 
    
     A tomogram due to the i-th chemical shift component is reconstructed by using the resonance signal S i  (t). 
     In order to attain the second object, according to another aspect of the present invention, there is provided an image reconstructing device, in which a plurality of sampling operations are performed for a resonance signal measured only once, to obtain necessary data, a tomogram ρ k  is reconstructed from the above data while paying attention to the k-th chemical shift component, and the influence of chemical shift components other than the k-th one on the tomogram ρ k  is estimated and removed. 
     In more detail, the initial estimated value ρ k .sup.(o) of the tomogram due to the k-th chemical shift component is set as follows: 
     
         ρ.sub.k.sup.(o) =ρ.sub.k 
    
     The (n+1)th estimated value ρ k .sup.(n+1) of the tomogram is calculated from the n-th estimated value ρ k  .sup.(n) (where n&gt;0) as follows: ##EQU7## where Λn is given by the following equation: ##EQU8## Further, R e  (C) in the equation (12) indicates the real part of a complex number C, and θ n  (≧0) indicates a threshold parameter dependent upon a number n. 
     In actual processing, the estimation based upon the equation (11) is repeated till an appropriate end condition is satisfied, to obtain an excellent tomogram of the k-th chemical shift component. 
     Now, let us consider the n-th estimated value ρ k .sup.(n) as an approximation to the true tomogram ρ k  of the k-th chemical shift component. Then, the influence of chemical shift components other than the k-th component on the estimated value ρ k .sup.(n) can be removed by using the equation (11). Thus, a better approximation to the true tomogram ρ k  is obtained. By repeating the above processing, the estimated value ρ k .sup.(n) can approach to the true tomogram ρ k . 
     The tomogram ρ k  is proportional to the spin density of the k-th chemical shift component. Accordingly, ρ k  is a real number, and has a value which is not negative. Further, the function Λn included in the equation (11) corresponds to an operation for taking out a component which is a real number and has a positive or zero value, from the n-th estimated value ρ k .sup.(n) at a time the threshold parameter θ n  is equal to zero. Accordingly, a value Λn [ρ k .sup.(n) ] comes closer to the true tomogram ρ k , as compared with the n-th estimated value ρ k .sup.(n). Thus, the accuracy of the estimated value, due to the equation (11) is improved. 
     Further, the speed of convergence can be increased by decreasing the threshold parameter θ n  as the number n is larger. This is because the second term on the right-hand side of the equation (4) contributes greatly to ρ k  for the number n having a small value, and a threshold parameter θ n  having a large value contributes only a little to the estimated value due to the equation (11). 
     As has been explained in the above, according to the present invention, chemical shift components are, from the first, separated from each other in a high-speed imaging method using time-varying gradient magnetic fields, to form a tomogram of each chemical shift component. Thus, a tomogram due to each chemical shift component can be constructed at high speed. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     FIG. 1 is a flow chart showing an example of the processing in an embodiment of the present invention. 
     FIG. 2 is a block diagram showing the circuit configuration of an MRI apparatus. 
     FIG. 3 is a schematic diagram showing an example of the pulse sequence used in the present invention. 
     FIG. 4 is a graph showing the distribution of signal measuring point in a spatial frequency region. 
     FIG. 5 is a schematic diagram showing another example of the pulse sequence used in the present invention. 
     FIG. 6 is a flow chart showing the processing in a first embodiment of the present invention. 
     FIG. 7 is a flow chart showing the processing in a second embodiment of the present invention. 
    
    
     DESCRIPTION OF THE PREFERRED EMBODIMENTS 
     Now, a first embodiment of the present invention will be explained, with reference to FIGS. 1 to 6. 
     FIG. 2 is a block diagram showing the circuit configuration of an MRI apparatus, in which the present invention is realized. Referring to FIG. 2, the MRI apparatus includes a sequence controller 201 for controlling each part of the apparatus in accordance with a predetermined procedure, a transmitter 202 for sending out a high-frequency pulse to generate magnetic resonance, a gradient magnetic field driver 204 for generating a gradient magnetic field, a magnetic field controller 203 for controlling the gradient magnetic field driver, a receiver 205 for receiving a magnetic resonance signal which is generated by an object to be inspected, a processor 206 for reconstructing a tomogram and for performing various arithmetic operations, a CRT display 207 for displaying a tomogram, and an external memory 208 for storing therein tomogram data and others. 
     FIG. 3 shows an example of the pulse sequence which is used in a first embodiment of the present invention. Referring to FIG. 3, a 90°-RF pulse 302 is generated at the same time as a pulse 301 for generating a gradient mangetic field in a Z-direction is applied to a coil. Thus, spins in a desired slice plane (namely, an X-Y plane in this example) are excited. Next, a 180°-RF pulse 303 for generating a spin echo and a pulse 304 for generating a gradient magnetic field in the Z-direction are simultaneously generated. Thereafter, a resonance signal is measured for a period T (indicated by reference numeral 307) in a state that a gradient magnetic field 305 in a Y-direction and a gradient magnetic field 306 in an X-direction are both generated. The strength of each of the gradient magnetic field 305 and the gradient magnetic field 306 varies with time. The gradient magnetic field 305 and the gradient magnetic field 306 are applied to obtain a resonance signal at various points in the spatial frequency region corresponding to the X-Y plane. Let us express the time interval between the generation of the 90°-RF pulse 302 and the generation of the 180°-RF pulse 303 by T E  /2. When a time T E  /2 has elapsed after the 180°-RF pulse 303 was generated, the phase shift of every atom is returned to zero, and a spin echo signal is generated. Further, when a time ΔT has elapsed after the spin echo signal was generated, the gradient magnetic field 305 in the Y-direction and the gradient magnetic field 306 in the X-direction are both generated, and the measurement of resonance signal is started. In the first embodiment of the present invention, a set of necessary data for reconstructing a tomogram is obtained from a plurality of times of measurements. The value of ΔT is varied in each of the measurements. The resonance data is obtained by a plurality of samplings in each measurement. 
     FIG. 4 shows an example of the distribution of k x  and k y  which indicate the integrated values of the gradient magnetic fields 306 and 305, respectively. 
     Now, let us suppose that G x  (t) and G y  (t) are given by the following equations: ##EQU9## where A and B are constants. 
     Then, k x  (t) and k y  (t) are given by the following equations: ##EQU10## 
     Thus, the combination of k x  (t) and k y  (t) will produce the locus shown in FIG. 4. In FIG. 4, reference numeral 401 indicates one of measuring points which are arranged at a time interval in the period 307. 
     FIG. 5 shows a sequence in a case where a 180°-RF pulse is not used, but an echo signal due to the reversal of a gradient magnetic field is measured. Referring to FIG. 5, a pulse 501 for generating a gradient magnetic field in the Z-direction and a 90°-RF pulse 502 are simultaneously generated to excite only spins in a desired slice plane, as in the sequence of FIG. 3. When a time ΔT has elapsed after the 90°-RF pulse 502 was generated, a resonance signal is measured for a period 505 in a state that a gradient magnetic field 503 in the Y-direction and a gradient magnetic field 504 in the X-direction are generated. 
     The present invention can be realized in the same manner for the sequence of FIGS. 3 and 5. 
     FIG. 1 is a flow chart showing the steps of procedure for reconstructing tomograms in accordance with the present invention, for a case where only two chemical shift components are detected. 
     In step 101 of FIG. 1, a resonance signal S.sub.ΔT.sbsb.1 (t) for ΔT=ΔT 1  and a resonance signal S.sub.ΔT.sbsb.2 (t) for ΔT=ΔT 2  (ΔT 2  ≠ΔT 1 ) are measured (where 0≦t≦T). It is to be noted that the above-mentioned ΔT indicates ΔT in FIGS. 3 and 5. 
     In step 102, a resonance signal S 1  (t) due to only the first chemical shift component and a resonance signal S 2  (t) due to only the second chemical shift component are determined by the following equations: 
     
         S.sub.1 (t)=(S.sub.ΔT.sbsb.1 (t)e.sup.-jΔω 2.sup.ΔT 1-S.sub.ΔT.sbsb.1 (t)e.sup.-jΔω 2.sup.ΔT 2) 
    
     
          /{e.sup.j(Δωw 1.sup.-Δω 2.sup.)ΔT 1-e.sup.j(Δω 1.sup.-Δω 2.sup.)ΔT 2} 
    
     
         S.sub.2 (t)=(S.sub.ΔT.sbsb.1 (t)e.sup.-jΔω 1.sup.ΔT 1.sup.- S.sub.ΔT.sbsb.2 (t)e.sup.-jΔω 1.sup.ΔT 2) 
    
     
          /{e.sup.j(Δω 2.sup.-Δω 1.sup.)ΔT 1-e.sup.j(Δω 2.sup.-Δω 1.sup.)ΔT 2} 
    
     In order to determine the resonance signals S 1  (t) and S 2  (t) by the above equations, it is required to satisfy the following formula: 
     
         (Δω.sub.2 -Δω.sub.1)×(ΔT.sub.2 -ΔT.sub.1)≠2nπ 
    
     where n is an integer. 
     Accordingly, the values of ΔT 1  and ΔT 2  are previously selected so as to satisfy, for example, the following equation: 
     
         (Δω.sub.2 -Δω.sub.1)×(ΔT.sub.2 -ΔT.sub.1)=π 
    
     In step 103, the phase orientation due to a change in resonance frequency is corrected as indicated by the following equations: ##EQU11## Thus, resonance signals S 1  (t) and S 2  (t), each of which is a resonance signal due to a specific chemical shift component, are obtained from the resonance signals S 1  (t) and S 2  (t). 
     In step 104, a tomogram ρ 1  (x, y) due to the first chemical shift component is determined from the resonance signal S 1  (t) by the following equation: ##EQU12## where W(t) is equal to √k x  (t) 2  +k y  (t) 2  for a case where the locus of FIG. 4 is used, and Re[] is a sign for indicating a real part. 
     In step 105, a tomogram ρ 2  (x, y) due to the second chemical shift component is determined, as in the step 104, from the resonance signal S 2  (t) by the following equation: ##EQU13## That is, the tomogram due to the first chemical shift component and the tomogram due to the second chemical shift component can be reconstructed by the above steps. 
     FIG. 6 is a flow chart showing the steps of procedure for reconstructing tomograms in accordance with the present invention, for a case where n chemical shift components are included. 
     In step 601 of FIG. 6, a resonance signal S.sub.ΔT.sbsb.1 (t) for ΔT=ΔT 1 , a resonance signal S.sub.ΔT.sbsb.2 (t) for ΔT=ΔT 2  and a resonance signal S.sub.ΔT.sbsb.n (t) for ΔT=ΔT n  are measured (where 0≦t≦T). It is to be noted that the above-mentioned T and ΔT indicate T and ΔT shown in FIGS. 3 and 5, respectively. 
     In step 602, a resonance signal S i  (t) due to only the i-th chemical shift component (where i=1, 2, . . . , and n) is determined by the following equation: ##EQU14## where a matrix A is given by the following equation: ##EQU15## Values of ΔT 1  to ΔT n  are selected so that an inverse matrix A -1  can be formed. 
     In step 603, a resonance signal S i  (t), in which the phase rotation due to a change in resonance frequency has been corrected, is obtained in the following manner: 
     
         S.sub.i (t)=S.sub.i (t)e.sup.-jΔω.sbsp.i.sup.t 
    
     where i=1, 2, 3, . . . , and n. 
     In step 604, a tomogram ρ i  (x, y) due to the i-th chemical shift component is determined from the resonance signal S i  (t) by the following equation: ##EQU16## That is, n tomograms corresponding to n chemical shift components can be separately reconstructed by the above steps. 
     As can be seen from the above explanation, the first embodiment is independent of how the strength of gradient magnetic field varies with time, and is also independent of the algorithm for reconstructing a tomogram. Further, the present embodiment is applicable to a case where a three-dimensional spin-density image is formed. In the first and second objects of the invention mentioned above, a resonance signal necessary for reconstructing a tomogram is measured by generating a single resonance. The present invention can be applied in the situation in which a resonance signal necessary for reconstructing a tomogram consists of a plurality of part signals. Each of the part signals is measured by generating a single resonance. In the situation in which the tomogram consists of a plurality of part signals, the value of ΔT shown in FIG. 3 or FIG. 5 is fixed to a certain value and measurements are carried out by generating a plurality of resonances. A set of resonance signals necessary for reconstructing a tomogram is obtained from the measurements. Then, the value of ΔT is set to another value and the same measurements are again taken, so that another set of resonance signals necessary for reconstructing a tomogram is obtained. 
     In the present invention, all the steps of procedure for reconstructing a tomogram include only linear operations. Accordingly, the order of the step for obtaining S i  (t) from S.sub.ΔT (t) and the step for forming a tomogram of a chemical shift component, may be inverted, that is, an arithmetic operation for obtaining S i  (t) from S.sub.ΔT (t) may be performed after an arithmetic operation for obtaining ρ i  (x, y) from S i  (t). 
     Next, a second embodiment of the present invention will be explained, with reference to FIG. 7. 
     FIG. 7 is a flow chart showing the steps of procedure for separating a tomogram due to water and a tomogram due to fat from each other by the second embodiment in a case where the atomic nucleus of hydrogen is detected by magnetic resonance imaging techniques. 
     Referring to FIG. 7, a resonance signal S(t) is measured on the basis of the pulse sequence of FIG. 3 (step 701). In step 702, a tomogram ρ 1 .sup.(o) due to water and a tomogram ρ 2 .sup.(o) due to fat are reconstructed by the following equation: 
     
         ρ.sub.i.sup.(o) (x, y)=ρdtS(t)exp[+jγ{kx(t)x+ky(t)y}]exp[-jΔω.sub.i t]W(t) i=1, 2 
    
     where the weight coefficient W(t) is given as follows: 
     
         W(t)={kx(t).sup.2 +ky(t).sup.2 }.spsp.1/2 
    
     In step 703, the value of a variable n is made equal to zero. 
     In step 704, ρ 2 .sup.(n+1) is calculated from ρ 1 .sup.(n) and ρ 2 .sup.(o) by the following equation: 
     
         ρ.sub.2.sup.(n+1) =ρ.sub.2.sup.(o) -h.sub.21 *Λ.sub.n.sup.[ ρ.sub.1.sup.(n)] 
    
     where h 21  indicates a PSF given by the equation (6) and Λ n  is a function defined by the equation (12). 
     In step 705, ρ 1 .sup.(n+1) is calculated from ρ 2 .sup.(n+1) and ρ 1 .sup.(o) by the following equation: ##EQU17## 
     In step 706, the value of the variable n is incremented by one. 
     In step 707, it is checked whether or not a predetermined constant N (for example, N=10) is greater than the value of the variable n. When it is judged that the constant N is greater than the value of n, the processing in the step 704 is again carried out. When it is judged that the constant N is less than or equal to the value of n, the processing in step 708 is carried out. In the step 708, ρ 1 .sup.(N) and ρ 2 .sup.(N) are delivered to terminate the separation processing. 
     The threshold parameter θ n  used in the step 704 is given by, for example, the following equation: ##EQU18## where max{ } indicates a maximum value. 
     As can be seen from the above explanation, the second embodiment is independent of how the strength of gradient magnetic field varies with time and is also imdependent of the algorithm for reconstructing a tomogram. Further, the second embodiment is applicable to a case where a three-dimensional spin-density image is formed. 
     In the steps 704 and 705, convolution is carried out on a tomogram. Alternatively, corresponding processing may be carried out for the resonance signal. That is, the processing in the step 704 may be replaced by the following processing. A corrected signal which is obtained by subtracting a component corresponding to ##EQU19## from the resonance signal S(t) is formed as follows: ##EQU20## Then, ρ 2 .sup.(n+1) is reconstructed from the corrected signal by the equation (3). The processing in the step 705 is replaced by processing similar to the above-mentioned. 
     Further, in the steps 704 and 705, only a real part on the left-hand side of each equation is used in the next arithmetic operation. Accordingly, in order to reduce the amount of calculation, the right-hand side of each equation may be calculated so that only a real part is obtained.