阅读说明：本技术 一种磁共振系统b0场图的校正方法 (Correction method of B0 field map of magnetic resonance system ) 是由 杨刚 张龙江 张志强 卢光明 许强 戚荣丰 张运明 张其锐 程晓青 邱连丽 于 2021-08-18 设计创作，主要内容包括：本发明提出了一种磁共振系统B0场图的校正方法,首先获取磁共振系统的三维B0场图,并沿每个维度划分为若干层；然后沿三维B0场图至少其中两个维度分别逐层进行2D相位去卷绕,获取中间B0场图；最后选定所述至少两个维度的其中一个作为基准方向,获取沿所述基准方向上B0场图中相邻两层的相位差,且该相位差满足全局优化规则。该方法基于对3D相位每个层面逐一做2D去卷的结果再做3D相位去卷绕计算,并引入了全局优化规则,因此匀场结果具有优良的稳定性和准确性；而且,本发明的方法计算效率高,更加符合实际应用中的快速响应需求,具有更好的临床应用前景。(The invention provides a correction method of a magnetic resonance system B0 field pattern, which comprises the steps of firstly obtaining a three-dimensional B0 field pattern of the magnetic resonance system, and dividing the field pattern into a plurality of layers along each dimension; then respectively carrying out 2D phase unwrapping layer by layer along at least two dimensions of the three-dimensional B0 field map to obtain a middle B0 field map; and finally, selecting one of the at least two dimensions as a reference direction, and acquiring the phase difference of two adjacent layers in the B0 field pattern along the reference direction, wherein the phase difference meets the global optimization rule. The method is based on the results of performing 2D deconvolution on each layer of the 3D phase one by one, then performing 3D phase deconvolution calculation, and introducing a global optimization rule, so that shimming results have excellent stability and accuracy; moreover, the method of the invention has high calculation efficiency, better meets the rapid response requirement in practical application and has better clinical application prospect.)

1. A correction method of a magnetic resonance system B0 field image is characterized by comprising the following steps:

s1, acquiring a three-dimensional B0 field map of the magnetic resonance system, and dividing the map into a plurality of layers along each dimension;

s2, respectively carrying out 2D phase unwrapping layer by layer along at least two dimensions of the three-dimensional B0 field map to obtain a middle B0 field map;

s3, selecting one of the at least two dimensions as a reference direction, and acquiring phase differences of two adjacent layers in the B0 field pattern along the reference direction as shimming parameters to realize unwrapping, wherein the phase differences of the two adjacent layers meet the global optimization rule.

2. The method for correcting the B0 field pattern of the magnetic resonance system as set forth in claim 1, wherein the global optimization rule is: using other directions different from the reference direction for the phase differenceVerifying, and taking the phase difference if the verification result is idealShimming is carried out; otherwise, the objective function Z is optimized, and the optimized phase difference is taken for shimming.

3. The method of claim 2, wherein the phase difference is verified in Y and Z directions with the X direction as a reference direction.

4. The method for correcting the B0 field pattern of the magnetic resonance system according to claim 3, wherein the verification method is: if:

（1）

and is

（2）

The verification result is ideal;

wherein FX is the intermediate phase diagram along the X direction; FY is the intermediate phase plot along the Y direction; FZ is the intermediate phase map along the Z direction; i is the coordinate in the X direction; j is the coordinate in the Y direction and k is the coordinate in the Z direction.

5. The method of claim 4, wherein if the verification result is not ideal, i.e. the equations (1) and (2) cannot be satisfied exactly at the same time, the following objective function Z is optimized:（3）

solving for minimising objective function ZAs the shimming parameters after correction.

6. A method as claimed in claim 5, wherein the modified shim parameters are used to correct a B0 field map of the magnetic resonance systemComprises the following steps:

（4）。

Technical Field

The invention relates to a correction method of a B0 field diagram, wherein a 3D phase unwrapping method is applied; belongs to the technical field of magnetic resonance imaging.

Background

In the magnetic resonance imaging system, the magnet itself can provide a relatively uniform static magnetic field B0 in a certain space, and the uniformity of the magnetic field is a key factor of the quality of the magnetic resonance image. However, due to the difference in magnetic susceptibility of the tissues, local magnetic field variation (i.e., inhomogeneity) may occur at the boundary, resulting in problems such as magnetic sensitivity artifacts, poor fat-pressing effect, and image deformation. Therefore, in order to improve the image quality, it is necessary to reduce the inhomogeneity of the magnetic field (i.e., shimming) so that the magnetic field becomes a spatially continuous distribution and uniform image without color jump. In the prior art, two-dimensional or three-dimensional phase unwrapping technology is generally adopted for B0 shimming, namely, B0 field images are corrected. Generally, compared with a two-dimensional phase unwrapping algorithm, three-dimensional phase unwrapping has the advantages of being more difficult, more difficult to guarantee stability, longer in calculation time and difficult to apply in practice.

The field pattern unwrapping technique is the core of the correction, and existing 3D unwrapping algorithms include path-guided methods, minimum p-norm methods, and the like. In the minimum p-norm method, a target phase diagram is preset firstly, three-direction difference is obtained point by point for an original phase diagram (a wrapped phase diagram), then the difference is subtracted by the difference corresponding to the preset target phase diagram to make difference, so that three difference values are obtained for each pixel point, then p-th power of the difference values in the three directions of all the points is added to construct an objective function, and the problem is converted into a problem of optimizing the objective function. The method has better accuracy and stability, a square matrix of N x N is constructed in the calculation process, N is the number of pixel points participating in calculation, and the calculation amount of three-dimensional de-winding is far larger than that of two-dimensional de-winding. In the path-oriented method, an evaluation function of the phase diagram quality is selected firstly, full wrapping starts from one point of three-dimensional data, the evaluation function is calculated along each direction, the direction with the highest evaluation value is selected, and then one-dimensional phase unwrapping is carried out along the direction. However, the path-oriented method lacks global calculation of the phase diagram, and when a certain path is wrong, calculation errors of subsequent pixel points on the path may be caused, and the method also lacks practical application value.

In view of the above, it is necessary to provide a better correction method for the B0 field pattern.

Disclosure of Invention

In order to overcome the defects of the prior art, the invention aims to provide a stable, accurate and high-calculation-efficiency B0 field pattern correction method.

In order to achieve the above object, the present invention adopts the following technical solutions:

the invention discloses a correction method of a B0 field pattern of a magnetic resonance system, which comprises the following steps:

s1, acquiring a three-dimensional B0 field map of the magnetic resonance system, and dividing the map into a plurality of layers along each dimension;

s2, respectively carrying out 2D phase unwrapping layer by layer along at least two dimensions of the three-dimensional B0 field map to obtain an intermediate B0 field map, namely an intermediate phase map;

s3, selecting one of the at least two dimensions as a reference direction, and acquiring phase differences of two adjacent layers in the B0 field pattern along the reference direction as shimming parameters to realize unwrapping, wherein the phase differences of the two adjacent layers meet the global optimization rule.

Preferably, the global optimization rule is as follows: using other directions different from the reference direction for the phase differenceVerifying, and taking the phase difference if the verification result is idealShimming is carried out; otherwise, the objective function Z is optimized to obtain the optimized phase difference for shimming.

Preferably, the X direction is used as a reference direction, and the Y and Z directions are used to verify whether the phase difference is ideal. Of course, the Y direction may be used as the reference direction, and the X and Z directions may be used to verify whether the phase difference is ideal.

More preferably, the verification method is as follows: if:

and is

（2）

It indicates that the verification result is ideal, and the unwrapping of the phase is ended and the phase difference is outputtedAs shimming parameters, then performing unwrapping operation on the next phase data until all calculations are completed;

wherein FX is the intermediate phase diagram along the X direction; FY is the intermediate phase plot along the Y direction; FZ is the intermediate phase map along the Z direction; i is a coordinate in the X direction; j is the coordinate in the Y direction and k is the coordinate in the Z direction. Thus, the coordinates of FX (i, j, K + 1) and FY (i, K +1, j) are the same, only the phase values of the intermediate phase map in different directions.

Further preferably, if the verification result is not ideal, i.e. the equations (1) and (2) cannot be strictly established at the same time, the objective function Z is optimized:（3）

solving for minimising objective function ZAnd the obtained value is the shimming parameter after correction.

More preferably, the specific optimization process is as follows:

finally, obtaining the corrected shimming parametersComprises the following steps:

（4）。

the invention has the advantages that:

(1) the invention provides a correction method of a B0 field map, wherein a 3D phase unwrapping method is applied, shimming can be carried out on a B0 field map of a magnetic resonance system, 3D phase unwrapping calculation is carried out on the basis of the result of carrying out 2D unwrapping on each layer of a 3D phase one by one, the information of the 2D unwrapping result is fully utilized, the calculation result of 2D unwrapping on each layer is not changed, 1D phase unwrapping is carried out on a third dimension by utilizing the existing 2D unwrapping result, and meanwhile, a global optimization rule is introduced, so that the correction result has excellent stability and accuracy, and the imaging quality is effectively improved;

(2) the method has high calculation efficiency and better meets the requirement of quick response in practical application. The size of the 3D phase data to be unwrapped is M × M, the calculated time to unwrapp a 2D phase map of M × M is t0, and the total time to complete unwrapping the 3D phase in the prior art is about 3 × M × t 0. And this application develops a new way, does not consider all pixel points to go on 3D to go to convolute from beginning on the whole, but carries out 2D to go to convolute layer by layer, utilizes 2D information to go on 1D to go to convolute on this basis at last, and this kind of 2+1 method has very big advantage at computational efficiency.

Drawings

FIG. 1 is a schematic diagram of data after a 3D phase data matrix is unwrapped from a 2D phase layer by layer along an X direction;

FIG. 2 is a schematic diagram of data after 2D phase unwrapping of a 3D phase data matrix layer by layer in the Y direction;

FIG. 3 is a schematic diagram of data after 2D phase unwrapping of a 3D phase data matrix layer by layer in the Z direction;

FIG. 4 is a graph of unwind data for the k-th layer in FX;

FIG. 5 is a graph of unwind data for the k +1 th layer in FX;

fig. 6 is a wrapped phase plot of a 32X 32 data phase matrix along the X direction for an embodiment of the present invention;

fig. 7 is a wrapped phase plot of a 32 x 32 data phase matrix along the Y direction for an embodiment of the present invention;

fig. 8 is a wrapped phase plot of a 32 x 32 data phase matrix along the Z direction of an embodiment of the present invention;

fig. 9 is a diagram showing FX along the X direction after the 32X 32 data phase matrix is unwrapped along the X direction slice 2D in accordance with an embodiment of the present invention;

fig. 10 is a diagram showing FX along the Z direction after the 32X 32 data phase matrix is unwrapped along the X direction slice 2D in accordance with an embodiment of the present invention;

FIG. 11 is a layer diagram of the FX-based 3D deconvolution calculation results of FIG. 9 along the X direction;

FIG. 12 is a layer diagram of the FX-based 3D deconvolution calculation of FIG. 9 along the Z-direction.

Detailed Description

The invention is described in detail below with reference to the figures and the embodiments.

The B0 field map correction method is suitable for B0 shimming of a magnetic resonance system, and specifically comprises the following steps:

s1, acquiring a three-dimensional B0 field map of the magnetic resonance system, and dividing the map into a plurality of layers along each dimension.

Usually, we acquire B0 field map by directly acquiring raw data, and in this embodiment, we can acquire three-dimensional B0 field phase map M by using fast gradient echo sequence. The method comprises the steps of collecting two echoes when the frequency coding gradient is in the same direction by using the fast gradient echo, wherein the echo time is respectively TE1 and TE2, filling the two echoes in K space 1 and K space 2 respectively, and obtaining complex images of two different echo times through Fourier reconstruction for calculating a B0 field phase diagram.

S2, performing 2D phase unwrapping on the 3D phase data matrix layer by layer along three directions X, Y and Z, respectively, to obtain layer unwrapped phase data in three directions, the results are shown in fig. 1, fig. 2 and fig. 3, respectively. The square broken line frame represents the 2D unwrapped layer, the hollow circle represents the data after the 2D unwrapping, and it can be seen that the layer directions of the three unwrapping layers are perpendicular to each other. For convenience, the phase diagrams in three directions after decoiling are sequentially denoted as FX, FY and FZ, specifically, taking the X direction as an example, wherein the data of the ith row and the jth column of the k layer is denoted as FX (i, j, k), and the recording of other phase points is the same.

Specifically, the layer-by-layer 2D phase unwrapping in this step includes the following sub-steps:

(1) acquiring a phase diagram to be processed;

(2) preprocessing the phase diagram to obtain a phase mask (mask) corresponding to the phase diagram;

(3) segmenting the phase diagram at least along a first direction and a second direction respectively according to the phase mask to obtain at least two segmented phase diagrams, wherein the segmented phase diagrams correspond to the segmented directions;

(4) respectively carrying out phase unwrapping processing on the at least two segmented phase diagrams to obtain at least two intermediate phase diagrams;

(5) and determining a unwrapped phase map corresponding to the phase map according to the at least two intermediate phase maps.

For the two-dimensional phase map of the present embodiment, one dimension corresponds to the phase encoding direction, and the other dimension corresponds to the frequency encoding direction.

S3, selecting one of the at least two dimensions as a reference direction, and acquiring phase differences of two adjacent layers in the B0 field pattern along the reference direction as shimming parameters to realize unwrapping, wherein the phase differences of the two adjacent layers meet the global optimization rule.

Herein the followingTo illustrate the example of interlayer unwinding in the X direction, starting with X first layer data, the phase difference to be added to the 2D unwinding result of the second layer is calculated by calculating the difference between FX second layer and first layer and performing contrast verification in the Y and Z directionsAnd sequentially calculating until the phase difference of the last layer is obtained, and finishing the calculation.

For a clearer understanding and for a more clear implementation of the invention, the interlayer unwinding of the k-th layer and the k +1 layer in FX is described below as an example. As shown in fig. 4 and 5, the points corresponding to the asterisks are the ith row and the jth column, and are respectively located in the k layer and the k +1 layer, the phase values thereof are respectively marked as FX (i, j, k) and FX (i, j, k + 1), the difference between the phase values of the k +1 layer and the k layer is calculated, and the offset value is obtained。

If the phase difference is not constantIf the verification result is ideal, the shim parameters can be directly used as the shim parameters; if the verification result is not ideal, global optimization is required. Therefore, the phase difference is corrected by the following methodAnd carrying out comparison and verification.

During verification, other directions different from the reference direction are adopted for comparison and verification. In this embodiment, if the interlayer unwinding is performed with the X direction as the reference direction, the unwinding result is verified in both the Y and Z directions, and if:

（1）

and is

（2）

It indicates that the verification result is ideal, at which point the unwrapping of the phase is ended and the calculation result is output, followed by the unwrapping operation of the next phase data.

If the unwrapped phase data is not ideal, i.e. equations (1) and (2) above cannot be strictly satisfied at the same time, then the objective function is optimized, and the objective function Z to be optimized is:

（3）

in the above equation (3), the optimized shimming parameters are obtained by solving the delta (k) that minimizes the objective function Z.

To solve for delta (k) that minimizes Z, the derivative of the objective function Z to delta (k) is required and is equal to 0:

can find out：

（5）

According to the method, delta (k) between adjacent layers after global optimization is obtained, and finally, delta (k) is added to the k +1 th layer data in each FX to obtain final 3D unwrapped phase data, so that B0 shimming is completed.

For better understanding and implementing the present invention, this embodiment constructs a 3D phase data 32 x 32 (i.e., M is 32), and fig. 6 to 8 are winding phase diagrams of the matrix along X, Y and Z directions, respectively, and it can be seen that the color jump is obvious and the winding degree is severe.

First, the slice unwrapping phase data after unwrapping along the X, Y and Z directions are obtained by using the 2D unwrapping algorithm, and fig. 9 is a diagram showing the data matrix of the present embodiment after 2D unwrapping along the X direction along the slice. As can be seen from fig. 9, there is a jump in brightness (flare) between different 2D unwrapped slices along the slice 2D unwrapped data in the X-direction because the slice 2D unwrapping calculations are independent of each other. This jump can be seen more clearly if FX is displayed in the Y-direction or Z-direction, and fig. 10 is a display of FX in the Z-direction after unwinding in the X-direction layer 2D, in which bright and dark jumps can be seen clearly. Therefore, the B0 shimming requirements are not yet met after 2D decoiling.

Then, 3D unwrapping calculation based on X direction is performed on the matrix data, where the global optimization rule mentioned above should be followed. FIG. 11 is a diagram showing FX-based 3D unwrapped slice representations along the X-direction, from which it is clear that the phase along the slice direction between slices has a gradual trend, but there are no distinct slice intensity jumps in FIG. 9.

Finally, in order to verify the shimming effect, the FX-based 3D unwrapping effect is shown along the Z direction, as shown in fig. 12, there is also a clear gradual trend between layers, the space is continuous, and there is no color jump, so the method of the present invention has an excellent B0 shimming effect.

In summary, the correction method of the B0 field pattern provided by the invention fully utilizes the information of the 2D unwrapping result, carries out 1D unwrapping on the third dimension, and utilizes the global optimization rule, so that the calculation result has excellent stability and accuracy. Moreover, the efficiency of the calculation method of the present invention is greatly improved, for example, the calculation time for unwrapping a 2D phase map of M × M is t0, the total time for completing unwrapping a 3D phase is about 3 × M × t0, if the calculation time/pixel point number is defined, the calculation time/pixel point number for unwrapping a 3D phase is only 3 times (independent of the pixel point) that for unwrapping a 2D phase, and for many 2D/3D general algorithms, the ratio of the 3D and 2D calculation time/pixel point number tends to increase linearly with the pixel point number. However, the algorithm develops a new method, all pixel points are not considered globally at the beginning to perform 3D unwinding, only 2D unwinding is performed layer by layer, and finally 1D unwinding is performed by using 2D information on the basis, so that the 2+1 method has great advantages in calculation efficiency.

The foregoing illustrates and describes the principles, general features, and advantages of the present invention. It should be understood by those skilled in the art that the above embodiments do not limit the present invention in any way, and all technical solutions obtained by using equivalent alternatives or equivalent variations fall within the scope of the present invention.