WO2019041123A1 - 一种优化平衡稳态自由进动序列的方法与装置 - Google Patents

一种优化平衡稳态自由进动序列的方法与装置 Download PDF

Info

Publication number
WO2019041123A1
WO2019041123A1 PCT/CN2017/099430 CN2017099430W WO2019041123A1 WO 2019041123 A1 WO2019041123 A1 WO 2019041123A1 CN 2017099430 W CN2017099430 W CN 2017099430W WO 2019041123 A1 WO2019041123 A1 WO 2019041123A1
Authority
WO
WIPO (PCT)
Prior art keywords
gradient
time
maximum
readout
phase
Prior art date
Application number
PCT/CN2017/099430
Other languages
English (en)
French (fr)
Inventor
朱燕杰
梁栋
刘新
郑海荣
Original Assignee
深圳先进技术研究院
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 深圳先进技术研究院 filed Critical 深圳先进技术研究院
Priority to PCT/CN2017/099430 priority Critical patent/WO2019041123A1/zh
Publication of WO2019041123A1 publication Critical patent/WO2019041123A1/zh

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves 
    • A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves  involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/565—Correction of image distortions, e.g. due to magnetic field inhomogeneities

Definitions

  • the present invention relates to the field of magnetic resonance imaging, and more particularly to a method and apparatus for optimizing a balanced steady state free precession sequence.
  • bSSFP Balanced steady state free precession
  • T2/T1 contrast is commonly used in cardiovascular magnetic resonance imaging.
  • bSSFP sequence is sensitive to the inhomogeneity of the main magnetic field, and it is easy to form "black band artifact" in the field of magnetic field inhomogeneity, and the longer the repeat time (TR) of the bSSFP sequence (for example)
  • TR repeat time
  • Shortening TR can effectively alleviate black-band artifacts while speeding up imaging and improving the signal-to-noise ratio (SNR) of the bSSFP sequence.
  • the current method of shortening TR mainly uses gradient overlap.
  • the gradients overlap the vector sum of the maximum gradient amplitude and the climb rate of the three axes of the readout phase gradient, the phase encode gradient and the layered gradient cannot exceed the maximum gradient amplitude and the climb rate allowed by the system, that is, the uniaxial climb rate.
  • the allowable gradient magnitude is constrained by the other two axes, so TR is primarily limited by the system hardware performance (maximum amplitude and rate of climb).
  • Another commonly used method for shortening TR is to use asymmetric echo: reduce the phase-to-phase moment in the readout direction, and make the echo appear early. After the signal acquisition is completed, a suitable back-concentration gradient is used in the readout direction. Maintain the sum of the gradient areas to zero and maintain the bSSFP characteristics of the sequence.
  • the ramp sampling technique can reduce the echo spacing of the sequence and has been widely used in Echo Planar Imaging (EPI) sequences.
  • EPI Echo Planar Imaging
  • the climb period sampling needs to keep the number of sampling points and the sampling bandwidth unchanged during sampling.
  • increasing the amplitude of the frequency encoding gradient can satisfy the spatial resolution requirement (ie, keeping the readout moment constant); when the sampling bandwidth is high, the gradient amplitude approaches or reaches the maximum amplitude of the readout gradient.
  • the sampling can only cover a small amount of climb time, resulting in a waste of part of the climb time.
  • Climb sampling is also used for ultra-short echo (UTR) radial sampling gradient echo (GRE) and bSSFP sequences for shortening TE time.
  • UTR ultra-short echo
  • GRE radial sampling gradient echo
  • bSSFP bSSFP sequences for shortening TE time.
  • the sampling is performed at the ramp-up of the gradient and does not need to be done at the ramp-down.
  • the phase-encoded gradient is not needed, and the readout gradient does not overlap with the selected gradient. Therefore, the climbing period sampling does not affect the gradient design of the other two axes.
  • Magnetic Resonance in Medicine 2013, 70: p657-663 also proposes a method of shortening bSSFP TR,
  • the non-layered pulse is used, so there is no need to consider the layer-selective convergence/scattering gradient;
  • the readout direction uses the asymmetric echo technique, the rising edge of the readout gradient and the phase encoding gradient overlap, and the readout decreases.
  • the climb-time sampling technique is used, because the readout gradient has gradient overlap, so the maximum gradient climb rate cannot be used; this scheme also does not perform the redistribution of the maximum gradient amplitude and the climb rate when the gradient overlaps.
  • the bSSFP sequence is highly susceptible to "black band” artifacts due to the non-uniformity of the main magnetic field, and shortening TR can reduce black band artifacts.
  • data sampling is performed only on the flattop time of the frequency encoding gradient, and the ramp time of the gradient can be overlapped with the phase encoding and the layered back/scatter gradient to shorten the bSSFP sequence. TR.
  • RF represents frequency excitation
  • the distance between two peaks on the horizontal axis represents TR
  • Slice represents the layer selection gradient
  • PE represents the phase encoding gradient
  • RO represents the frequency encoding gradient (or readout gradient)
  • ADC represents Digital analog conversion.
  • the readout phase gradient and the backscatter gradient are inconsistent, which breaks the (1st moment) balance of the 1st moment of the readout axis, so the sequence is more prone to blood flow artifacts.
  • the asymmetric echo causes some K-space data to be missing, and the missing part needs to be filled by the partial Fourier algorithm.
  • the quality of the image depends on the reconstruction algorithm, which may affect the sharpness of the image.
  • the present invention proposes a method for optimizing an equilibrium steady-state free precession sequence, which adopts a method of data acquisition in the entire frequency encoding gradient range, shortens the duration of the frequency encoding gradient, and redistributes each time.
  • the maximum amplitude and rate of climb of the gradient axis achieve the shortest effect of TR.
  • the method of the invention comprises the following steps:
  • Step a performing a sequence gradient design on the equilibrium steady-state free precession sequence so that the entire frequency encoding gradient time does not overlap with the other two-axis gradients, ie, the layer selection gradient and the phase encoding gradient;
  • Step b using the maximum gradient amplitude and rate of climb allowed by the system in the frequency encoding gradient, ie, the readout gradient;
  • Step c performing data sampling in the entire readout gradient time of the frequency encoding gradient while maintaining the image resolution unchanged, that is, the data acquisition covers the rising edge, the plateau period and the falling edge of the readout gradient;
  • Step d reallocating the maximum amplitude and the maximum climb rate of the three gradient axes of the frequency encoding gradient, the layer selection gradient and the phase encoding gradient, so that the gradient axis with the largest moment has the highest gradient maximum amplitude and the maximum climbing rate; The occupation time required for the three gradients is minimized, resulting in the shortest repetition time TR.
  • step c the sampling point time interval is kept unchanged when data is collected.
  • step c further comprises the following steps:
  • Step c1 calculating the readout moment according to the field of view of the readout direction set in the equilibrium steady-state free precession sequence and the size of the readout matrix
  • the FOV is the field of view in the readout direction
  • the unit is millimeter
  • the Readout data points are the readout matrix size
  • Step c3 According to the calculation result of step c2, the readout time is obtained: 2*Ramp Time+Flattop time;
  • Step c4 Adjust the read bandwidth so that the number of collected points is not less than the original matrix size.
  • step d an optimal allocation scheme of three axial gradient performances is obtained by an exhaustive method, including the following steps:
  • Step d2 setting intermediate parameters, including: Gro, Gpe, and Gss are respectively the maximum gradient amplitudes allowed by the readout phase gradient, the phase encoding gradient, and the layer selection gradient; Rro, Rpe, and Rss are readout phase gradients, respectively.
  • Step d3 Set the search range of Gro and Rro, and calculate the required dephasing time (Ddephase) for each point in the search range, as follows:
  • Dpe Mpe/Gpe+Gpe ⁇ Rpe
  • Step d4 obtaining a series of phase-spreading times by step c3, selecting a minimum of the phase-spreading time, wherein the minimum phase-to-phase corresponding to Gro, Gpe, Gss, Rro, Rpe, and Rss is the maximum gradient of the three gradient axes Amplitude and maximum rate of climb.
  • the method is also combined with an asymmetric echo.
  • the asymmetry factor is defined as the actual number of points/full number of points.
  • the equilibrium steady state free precession sequence can be replaced with a multi-echo gradient echo sequence.
  • step e is further included, and according to the symmetry of the bSSFP sequence, the same gradient performance allocation and processing are also performed on the frequency encoding back-concentration gradient, the layer-selected phase-diffusion gradient and the phase-encoding recovery gradient after the readout gradient.
  • the present invention also provides a computer readable medium having a program stored therein that is computer executable to cause a computer to perform processing including the steps described above.
  • the present invention also provides an apparatus for optimizing the equilibrium steady state free precession sequence bSSFP for the above method.
  • the device includes different modules for implementing the various steps mentioned in the above method.
  • the apparatus may comprise: a sequence gradient design module that performs a sequence gradient design on the equilibrium steady-state free precession sequence such that the entire frequency encoding gradient time does not overlap with the other two-axis gradients, ie, the layer selection gradient and the phase encoding gradient Maximum gradient amplitude and climb rate determination mode Block, using the maximum gradient amplitude and rate of climb allowed by the system in the frequency encoding gradient, ie, the readout gradient; the data sampling module performs the entire readout gradient time of the frequency encoding gradient while maintaining the image resolution unchanged Data sampling, that is, data acquisition covers the rising edge, plateau period and falling edge of the readout gradient; and the maximum amplitude and maximum rate of climb reallocation module, redistributing the three gradient axes of frequency encoding gradient, layer selection gradient and phase encoding gradient The maximum amplitude and maximum rate of climb make the gradient axis with the largest moment have the highest gradient maximum amplitude and maximum climb rate; thus the time required for the three gradients is minimized, thus obtaining the
  • the data sampling module keeps the sampling point time interval unchanged when collecting data.
  • the data sampling module further comprises the following modules:
  • the FOV is the field of view in the readout direction, the unit is millimeter, and the Readout data points are the readout matrix size;
  • Reading time acquisition module according to the above calculation result, obtaining readout time: 2*Ramp Time+Flattop time;
  • the bandwidth adjustment module adjusts the read bandwidth so that the number of collection points is not less than the original matrix size.
  • an optimal allocation scheme for the three axis gradient performance is obtained by an exhaustive method, and the following modules are also included:
  • the intermediate parameter setting module sets intermediate parameters, including: Gro, Gpe, and Gss are respectively the maximum gradient amplitude allowed by the readout phase gradient, the phase encoding gradient, and the layer selection gradient; Rro, Rpe, and Rss are readout densities, respectively.
  • the dephasing time calculation module sets the search range of Gro and Rro. For each point in the search range, the required dephasing time (Ddephase) is calculated separately.
  • the calculation method is as follows:
  • Dpe Mpe/Gpe+Gpe ⁇ Rpe
  • the maximum gradient amplitude and the maximum rate of climb obtaining module obtains a series of phase-spreading times by the phase-dissipation time calculation module, and selects the smallest of the phase-spreading times, and the minimum phase-spreading time corresponds to Gro, Gpe, Gss, Rro, Rpe, Rss is the maximum gradient amplitude and maximum climb rate of the three gradient axes.
  • the device is also combined with an asymmetric echo.
  • the asymmetry factor is defined as the actual number of points/full number of points.
  • the equilibrium steady state free precession sequence can be replaced with a multi-echo gradient echo sequence.
  • the apparatus further includes a processing module that performs the same gradient performance allocation and processing on the frequency encoding back-concentration gradient, the layer-selective phase-diffusion gradient, and the phase-encoding recovery gradient after the readout gradient according to the symmetry of the bSSFP sequence.
  • a processing module that performs the same gradient performance allocation and processing on the frequency encoding back-concentration gradient, the layer-selective phase-diffusion gradient, and the phase-encoding recovery gradient after the readout gradient according to the symmetry of the bSSFP sequence.
  • the method of the invention performs data acquisition in the time range of the entire frequency encoding gradient, and can effectively shorten the duration of the frequency encoding gradient.
  • the gradient of this axis can use the system's fastest climb rate and the highest gradient amplitude, which speeds up the acquisition and shortens the gradient time.
  • the method of the invention can also be combined with asymmetric echo techniques to further shorten the TR.
  • FIG. 1 is a schematic diagram of a prior gradient pre-optimization gradient scheme.
  • FIG. 2 is a schematic illustration of an optimized gradient overlap scheme of the present invention.
  • Figure 3 is a TR comparison diagram of three methods of pre-optimization (std-bSSFP), optimization (rs-bSSFP), and optimization of the pre-sequence using asymmetric echo (asym-bSSFP).
  • FIG. 6 is a flow chart of a method of optimizing an equilibrium steady state free precession sequence of the present invention.
  • step a sequence gradient design of the equilibrium steady-state free precession sequence (bSSFP), so that the gradient of the entire frequency encoding gradient time (ie, the readout gradient time RO) does not overlap with the other two axes, ie, the layer selection gradient (Slice) And the phase encoding gradient (PE) overlaps.
  • the distance between the two dashed lines is the readout gradient time.
  • FIG. 2 is a schematic diagram of an optimized gradient overlapping scheme of the present invention. It can be seen from Fig. 2 that during the readout gradient time, the ramp up time, flattop time and ramp down time of the readout gradient do not occur with the gradients of the other two axes. overlapping.
  • Step b Use the maximum gradient amplitude and rate of climb allowed by the system in the frequency encoding gradient, ie the readout gradient. These two value systems are automatically given.
  • Step c In the case of keeping the image resolution unchanged, the whole of the frequency encoding gradient
  • the data is sampled within the readout gradient time, that is, the data acquisition covers the rising edge, the plateau period, and the falling edge of the readout gradient.
  • the sampling point time interval is maintained at the time of acquisition.
  • performing data collection in the above step c further includes the following steps:
  • Step c1 calculating the readout moment according to the field of view of the readout direction set in the equilibrium steady-state free precession sequence and the size of the readout matrix
  • the FOV is the field of view in the readout direction
  • the unit is millimeter (mm)
  • the Readout data points is the size of the readout matrix
  • Step c3 According to the calculation result of step c2, the readout time is obtained: 2*Ramp Time+Flattop time; wherein the time of the rising edge and the falling edge are equal.
  • Step c4 Adjust the read bandwidth so that the number of collected points is not less than the original matrix size.
  • the readout gradient time can be obtained.
  • step d redistribute the maximum amplitude and maximum climb rate of the three gradient axes of the frequency encoding gradient, the layer selection gradient and the phase encoding gradient, so that the gradient axis with the largest moment has the highest gradient, the maximum amplitude and the maximum climb. rate.
  • the optimal allocation scheme of the three axial gradient performances can be obtained by the exhaustive method, that is, the maximum gradient amplitude and the maximum climbing rate of the three gradient axes are obtained, which specifically includes the following steps.
  • Step d2 setting intermediate parameters, including: Gro, Gpe, and Gss are respectively the maximum gradient amplitudes allowed by the readout phase gradient, the phase encoding gradient, and the layer selection gradient; Rro, Rpe, and Rss are readout phase gradients, respectively.
  • Dpe Mpe/Gpe+Gpe ⁇ Rpe
  • Step d4 obtaining a series of phase-spreading times by step c3, selecting a minimum of the phase-spreading time, wherein the minimum phase-to-phase corresponding to Gro, Gpe, Gss, Rro, Rpe, and Rss is the maximum gradient of the three gradient axes Amplitude and maximum rate of climb.
  • step d the maximum amplitude and the rate of climb of the three gradient axes are re-divided, so that the gradient axis with the largest moment has a higher gradient maximum amplitude and the rate of climb, so that the occupation time required for the three gradients is minimized, thereby Let TR be minimized. Because the three gradients occupy a period of time together, the longest gradient required determines the length of time.
  • Step d may further include step e, according to the symmetry of the bSSFP sequence, the frequency encoding back-concentration gradient, the layer-selected phase-diffusion gradient, and the phase-encoding recovery gradient after the readout gradient are also performed.
  • step e the frequency encoding back-concentration gradient, the layer-selected phase-diffusion gradient, and the phase-encoding recovery gradient after the readout gradient are also performed.
  • This method effectively reduces the TR of the entire bSSFP sequence with fixed system hardware conditions (maximum gradient amplitude and rate of climb allowed by the system).
  • the present invention draws the data in the table as a graph, as shown in FIG. 3, and after optimization, the TR is significantly shortened.
  • the TR is significantly shortened.
  • the readout matrix size/resolution (mm) is 256/1.17
  • TR 3.32 before optimization
  • TR 2.92 after optimization
  • the TR is obviously shortened.
  • the optimized TR of the present invention is slightly lower than the TR of the standard sequence using asymmetric echo techniques.
  • the asymmetric echo technique can also be added to the optimized sequence of the present invention to further shorten the TR, that is, the asymmetric echo technique is added to the foregoing embodiment 1, as Embodiment 2, which will be described in detail below.
  • Embodiment 2 which will be described in detail below.
  • the same points as in Embodiment 1 are not described herein again.
  • the asymmetric factor needs to be defined as the actual number of points taken/the number of points taken.
  • RO moment asymmetric RO moment ⁇ asymmetric factor.
  • the K-space data is multiplied by a linear phase with a phase slope of G ramp is the readout gradient amplitude after the ramp sampling method, and ⁇ t is the time interval of the sampling point, which is equivalent to shifting the image in the readout direction - ⁇ s;
  • Multiply the meshed data by a linear phase with a phase slope of G std is the amplitude of the readout gradient required for the unsampled ramp sampling method, which is equivalent to moving the image back to the initial position in the readout direction.
  • the TR can be further shortened by incorporating the asymmetric echo technique in the optimized sequence of the present invention as described above.
  • the method of the invention can also be applied to multi-echo gradient echo sequences, reducing echo spacing, scanning trajectories in Cartesian coordinates, and sampling methods for radial trajectories. .
  • This method is most effective in two-dimensional imaging or non-isotropic imaging.
  • the present invention provides a computer readable medium having a program stored therein, the program being computer executable to cause a computer to execute including The processing of each step is described.
  • the present invention also provides an apparatus for optimizing the equilibrium steady state free precession sequence bSSFP for the above method.
  • the device includes different modules for implementing the various steps mentioned in the above method.
  • the apparatus may include: a sequence gradient design module that performs a sequence gradient design on the equilibrium steady-state free precession sequence such that the entire frequency encoding gradient time does not overlap with the other two-axis gradients, ie, the layer selection gradient and the phase encoding gradient; Maximum gradient amplitude and climb rate determination module, using the maximum gradient amplitude and climb rate allowed by the system in the frequency coding gradient, ie, the readout gradient; the data sampling module, in the case of maintaining the image resolution unchanged, in the frequency encoding gradient The data is sampled during the entire readout gradient time, that is, the data acquisition covers the rising edge, the plateau period and the falling edge of the readout gradient; and the maximum amplitude and maximum rate of climb reallocation module, reallocating the frequency encoding gradient, the layer selection gradient and The maximum amplitude and maximum rate of climb of the three gradient axes of the phase encode gradient, so that the gradient axis with the largest moment has the highest gradient maximum amplitude and maximum climb rate; thus the time required for the three gradient
  • the data sampling module keeps the sampling point time interval unchanged when collecting data.
  • the data sampling module further comprises the following modules:
  • the FOV is the field of view in the readout direction, the unit is millimeter, and the Readout data points are the readout matrix size;
  • the reading time obtaining module obtains the reading time according to the calculation result of step c2: 2*Ramp Time+Flattop time;
  • the bandwidth adjustment module adjusts the read bandwidth so that the number of collection points is not less than the original matrix size.
  • the optimal allocation scheme for the three axis gradient performance can be obtained by the exhaustive method, and the following modules are also included:
  • the intermediate parameter setting module sets intermediate parameters, including: Gro, Gpe, and Gss are respectively the maximum gradient amplitude allowed by the readout phase gradient, the phase encoding gradient, and the layer selection gradient; Rro, Rpe, and Rss are readout densities, respectively.
  • the dephasing time calculation module sets the search range of Gro and Rro. For each point in the search range, the required dephasing time (Ddephase) is calculated separately.
  • the calculation method is as follows:
  • Dpe Mpe/Gpe+Gpe ⁇ Rpe
  • the maximum gradient amplitude and the maximum rate of climb obtaining module obtains a series of phase-spreading times by the phase-dissipation time calculation module, and selects the smallest of the phase-spreading times, and the minimum phase-spreading time corresponds to Gro, Gpe, Gss, Rro, Rpe, Rss is the maximum gradient amplitude and maximum climb rate of the three gradient axes.
  • the device is also combined with an asymmetric echo.
  • the asymmetry factor is defined as the actual number of points/full number of points.
  • the equilibrium steady state free precession sequence may be replaced with a multi-echo gradient echo sequence.
  • the apparatus can also include a processing module that performs the same gradient performance allocation and processing on the frequency encoded back-concentration gradient, the layer-selected phase-diffusion gradient, and the phase-coded recovery gradient after the readout gradient based on the symmetry of the bSSFP sequence.

Landscapes

  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • General Health & Medical Sciences (AREA)
  • Radiology & Medical Imaging (AREA)
  • Engineering & Computer Science (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Signal Processing (AREA)
  • Medical Informatics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Biophysics (AREA)
  • Pathology (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • General Physics & Mathematics (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

一种优化平衡稳态自由进动序列的方法，包括以下步骤：步骤a：在频率编码梯度的整个读出梯度时间内进行数据采样；步骤b：对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度发生重叠；步骤c：在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；步骤d：重新分配频率编码散相梯度、层选回聚梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；使三个梯度共同占用的时间降到最短，相应的，频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理，从而获得最短的重复时间。

Description

一种优化平衡稳态自由进动序列的方法与装置 技术领域
本发明涉及磁共振成像领域，尤其涉及一种优化平衡稳态自由进动序列的方法与装置。
背景技术
平衡稳态自由旋进(balanced steady state free precession,bSSFP)是心血管成像中最常用的序列之一，具有较高的信噪比和良好的T2/T1对比度，常用于心血管磁共振成像上，如心脏电影(cardiac cine)，单次激发成像(single-shot imaging)等。然而，bSSFP序列对主磁场的不均匀性较为敏感，在磁场不均匀处易形成“黑带”伪影(dark band artifact)，而且，bSSFP序列的重复时间(repeat time，TR)越长(例如，在高分辨率成像的情况下)，对主磁场的不均匀性越敏感。这个问题在高场磁共振成像中尤为突出。缩短TR可以有效地减轻黑带伪影，同时加快成像速度，提高bSSFP序列的信噪比(SNR)。
目前缩短TR的方法主要采用梯度重叠。梯度重叠时，读出散相梯度、相位编码梯度和层选梯度三个轴的最大梯度幅值和爬升率的矢量和不能超过系统允许的最大梯度幅值和爬升率，即单轴的爬升率和允许的梯度幅值会受到另外两轴的制约，因此，TR主要受系统硬件性能(梯度的最大幅值和爬升率)的限制。
另外一个常用的缩短TR的方法是采用不对称回波(asymmetric echo)：将读出方向的散相矩减少，令回波提早出现，信号采集完成后，在读出方向采用合适的回聚梯度，维持梯度面积总和为零，保持序列的bSSFP特性。
另外，爬升期采样(ramp sampling)技术可以减小序列的回波间距(echo spacing)，已经广泛应用于平面回波成像(Echo Planar Imaging,EPI)序列中。在EPI中，爬升期采样需要保持采样点数和采样带宽不变，在采样 时间固定的情况下，提高频率编码梯度的幅值可以满足空间分辨率的要求(即保持读出矩不变)；当采样带宽较高时，梯度幅值接近或达到读出梯度的最大幅值，无法进一步提高，为满足读出矩的要求，采样仅能覆盖少量爬升时间，造成部分爬升时间的浪费。爬升期采样也用于超短回波(ultra-short TE，或UTE)的径向采样的梯度回波(GRE)和bSSFP序列，用于缩短TE时间。在这两个应用中，由于重点在于减短TE的时间，采样在梯度上升沿(ramp-up)时进行，不需要在梯度下降沿(ramp-down)进行。在UTE采用的超短径向采样轨迹中，不需要相位编码的梯度，读出梯度也不会与选层梯度重叠，因此爬升期采样没有影响另外两轴的梯度设计。
Oliver Bieri在2013年第70期《医学磁共振》第657～663页发表的《超快速稳态自由进动及其在1.5特斯拉的体内1H形态和功能性肺成像中的应用》(Oliver Bieri,Ultra-Fast Steady State Free Precession and Its Application to In Vivo 1H Morphological and Functional Lung Imaging at 1.5Tesla.Magnetic Resonance in Medicine2013,70:p657–663)中也提出了一种缩短bSSFP TR的方法，此应用中，采用非层选脉冲，因此无需考虑层选回聚/散相梯度；读出方向采用了不对称回波技术，读出梯度的上升沿和相位编码梯度重叠，而读出的下降沿采用了爬升期采样技术，因为读出梯度存在梯度重叠，因此无法采用最大梯度爬升率；该方案也未进行梯度重叠时最大梯度幅值和爬升速率的重新分配。
如上所述，在高场磁共振成像中，由于主磁场的不均匀性，bSSFP序列极易形成“黑带”伪影，而缩短TR可以减少黑带伪影。在通常的bSSFP序列中，只在频率编码梯度的平台期(flattop time)进行数据采样，梯度的爬升期(ramp time)可以和相位编码及层选回聚/散相梯度重叠，来缩短bSSFP序列的TR。实际中，由于频率编码方向的图像分辨率通常较高，读出方向的散相梯度所需时间要大于相位编码梯度和层选回聚梯度的时间，因此频率编码梯度的爬升时间是浪费掉的，如图1所示，RF表示频率激发，横轴上两个峰之间的距离代表TR，Slice表示层选梯度，PE表示相位编码梯度，RO表示频率编码梯度(或读出梯度)，ADC表示数字模拟转换。
在不对称回波bSSFP中，读出散相梯度和回聚梯度不一致，这打破了读出轴1阶矩的(1st moment)平衡，因此序列较容易产生血流伪影。 同时，不对称回波造成部分K空间数据缺失，缺少的部分需要使用部分傅里叶算法进行填充，图像的质量依赖于重建算法，可能会影响图像的清晰度。
发明内容
针对现有技术的上述问题，本发明提出了一种优化平衡稳态自由进动序列的方法，采用在整个频率编码梯度范围内进行数据采集的方式，缩短频率编码梯度的时长，并重新分配每个梯度轴的最大幅值和爬升率，达到TR最短的效果。
本发明方法包括以下步骤：
步骤a：对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；
步骤b：在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；
步骤c：在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；
步骤d：重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的占用时间降到最短，从而获得最短的重复时间TR。
优选地，步骤c中，采集数据时保持采样点时间间隔不变。
优选地，步骤c还包括以下步骤：
步骤c1：根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
Figure PCTCN2017099430-appb-000001
其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
步骤c2：为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
步骤c3：根据步骤c2的计算结果，获得读出时间：2*Ramp Time+Flattop time；
步骤c4：调整读出的带宽，使得采集点数不小于原始矩阵大小。
优选地，步骤d中，采用穷举法获得三个轴梯度性能的最优分配方案，包括如下步骤：
步骤d1：设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
步骤d2：设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
步骤d3：设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
Figure PCTCN2017099430-appb-000002
Dro＝Mro/Gro+Gro×Rro；
Dpe＝Mpe/Gpe+Gpe×Rpe；
Mss＝Grf×(Drf+Grf×Rss)/2；
Dss＝Mss/Gss+(Gss+Grf)×Rss；
Ddephase＝max(Dro,Dpe,Dss)；
步骤d4：通过步骤c3得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
优选地，步骤d3中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
优选地，所述方法还与不对称回波结合。
优选地，与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。
优选地，在计算读出矩时，计算公式为RO moment asymmetric＝RO moment×asymmetric factor。
优选地，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
优选地，步骤d之后还包括步骤e，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
本发明还提供了一种计算机可读介质，该计算机可读介质具有存储在其中的程序，该程序是计算机可执行的，以使计算机执行包括上述各步骤的处理。
本发明还提供了一种用于上述方法的优化平衡稳态自由进动序列bSSFP的装置。该装置包括了不同模块用于实现上述方法中提及的各个步骤。
优选地，该装置可以包括：序列梯度设计模块，对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；最大梯度幅值和爬升率确定模 块，在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；数据采样模块，在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；以及最大幅值和最大爬升率重新分配模块，重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的时间降到最短，从而获得最短的重复时间TR。
优选地，所述数据采样模块在采集数据时保持采样点时间间隔不变。
优选地，所述数据采样模块还包括以下模块：
读出矩计算模块，根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
Figure PCTCN2017099430-appb-000003
其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
读出梯度爬升时间获得模块，为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
读出时间获得模块，根据上述计算结果，获得读出时间：2*Ramp Time+Flattop time；和
带宽调整模块，调整读出的带宽，使得采集点数不小于原始矩阵大小。
优选地，最大幅值和最大爬升率重新分配模块中，采用穷举法获得三个轴梯度性能的最优分配方案，还包括如下模块：
输入参数设置模块，设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO  moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
中间参数设置模块，设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
散相时间计算模块，设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
Figure PCTCN2017099430-appb-000004
Dro＝Mro/Gro+Gro×Rro；
Dpe＝Mpe/Gpe+Gpe×Rpe；
Mss＝Grf×(Drf+Grf×Rss)/2；
Dss＝Mss/Gss+(Gss+Grf)×Rss；
Ddephase＝max(Dro,Dpe,Dss)；
最大梯度幅值和最大爬升率获得模块，通过散相时间计算模块得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
优选地，散相时间计算模块中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
优选地，所述装置还与不对称回波结合。
优选地，与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。
优选地，在计算读出矩时，计算公式为 RO moment asymmetric＝RO moment×asymmetric factor。
优选地，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
优选地，该装置还包括处理模块，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
有益效果
本发明方法在整个频率编码梯度的时间范围内进行数据采集，可有效的缩短频率编码梯度的时长。
由于频率编码梯度不再与其他两轴的梯度重叠，此轴的梯度可以采用系统的最快爬升速率和最高的梯度幅值，加快了采集速度，缩短了梯度时间。
在读出梯度前面的三个梯度(即读出散相梯度，相位编码梯度，层选回聚梯度)，也可以通过重新分配每个梯度轴的最大梯度幅值和爬升率，达到TR最短的效果；由于bSSFP序列的对称性，读出梯度后面的三个梯度(即读出回聚梯度，相位编码回聚梯度，层选散相梯度)可采用和之前完全相同的方法分配梯度轴的最大梯度幅值和爬升率。该方法可在系统硬件条件(系统允许的最大梯度幅值和爬升率)固定的情况下，有效地减短bSSFP的TR。
本发明方法还可以和不对称回波技术结合起来，进一步缩短TR。
附图说明
图1是现有的优化前的梯度重叠方案的示意图。
图2是本发明的优化后的梯度重叠方案的示意图。
图3是优化前(std-bSSFP)，优化后(rs-bSSFP)，以及优化前序列采用不对称回波(asym-bSSFP)的三种方法的TR比较图。
图4(a)是优化后的bssfp采集的图像，TR＝3.2ms，图4(b)是未优化的bssfp采集的图像，TR＝3.6ms。
图5(a)是优化后的bssfp采集的图像，TR＝2.7ms，图5(b)是未优化的bssfp采集的图像，TR＝3.16ms。
图6是本发明的优化平衡稳态自由进动序列的方法的流程图。
具体实施方式
在下列说明中，为了提供对本发明的彻底了解而提出许多具体细节。本发明可在不具有部分或所有这些具体细节的情况下实施。在其他情况下，为了不对本发明造成不必要的混淆，不详述众所周知的过程操作。虽然本发明将结合具体实施例来进行说明，但应当理解的是，这并非旨在将本发明限制于这些实施例。
实施例1
下面结合图6来说明本发明实施例1的优化平衡稳态自由进动序列的方法的具体步骤。
首先，步骤a：对平衡稳态自由进动序列(bSSFP)进行序列梯度设计，使得整个频率编码梯度时间(即读出梯度时间RO)内不和其他两轴的梯度即层选梯度(Slice)及相位编码梯度(PE)发生重叠。RO所在的轴上，两条虚线之间的距离即为读出梯度时间。如图2所示为本发明的优化后的梯度重叠方案的示意图。从图2中可以看出，在读出梯度时间内，读出梯度的上升沿(ramp up time)、平台期(flattop time)和下降沿(ramp down time)均未与其他两轴的梯度发生重叠。
步骤b：在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率。这两个值系统会自动给出。
步骤c：在保持图像分辨率不变的情况下，在频率编码梯度的整 个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿。优选地，在采集时保持采样点时间间隔不变。
具体地，上述步骤c中进行数据采集还包括以下步骤：
步骤c1：根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
Figure PCTCN2017099430-appb-000005
其中，FOV为读出方向的视野，单位为毫米(mm)，Readout data points为读出矩阵大小；
步骤c2：为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
步骤c3：根据步骤c2的计算结果，获得读出时间：2*Ramp Time+Flattop time；其中，上升沿和下降沿的时间是相等的。
步骤c4：调整读出的带宽，使得采集点数不小于原始矩阵大小。
通过上述步骤，可以获得读出梯度时间。
接下来，步骤d：重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度、最大幅值和最大爬升率。
步骤d中，可以采用穷举法获得三个轴梯度性能的最优分配方案，即获得三个梯度轴最大梯度幅值和最大爬升率，具体包括如下步骤。
步骤d1：设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
步骤d2：设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率散相编码时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选回聚梯度所需时间的加和；
步骤d3：设置Gro，Rro的搜索范围，例如，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
Figure PCTCN2017099430-appb-000006
Dro＝Mro/Gro+Gro×Rro；
Dpe＝Mpe/Gpe+Gpe×Rpe；
Mss＝Grf×(Drf+Grf×Rss)/2；
Dss＝Mss/Gss+(Gss+Grf)×Rss；
Ddephase＝max(Dro,Dpe,Dss)；
步骤d4：通过步骤c3得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
通过步骤d，重新分了三个梯度轴的最大幅值和爬升率，使得矩最大的梯度轴具有较高的梯度最大幅值和爬升率，使三个梯度需要的占用时间降到最短，从而令TR减至最短。因为这三个梯度共同占用一段时间，所需时间最长的梯度决定了这段时间的长短。
步骤d之后还可以包括步骤e，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行 相同的梯度性能分配和处理。该方法可在系统硬件条件(系统允许的最大梯度幅值和爬升率)固定的情况下，有效地减短整个bSSFP序列的TR。
如下表所示，可以比较优化前(std-bSSFP)，优化后(rs-bSSFP)，以及优化前序列用不对称回波(asym-bSSFP)三种方法的TR。
Figure PCTCN2017099430-appb-000007
为直观显示，本发明将表中的数据画成了曲线图，如图3所示,经过优化，TR明显缩短。例如在读出矩阵大小/分辨率(mm)为256/1.17的情况下，优化前TR＝3.32，优化后TR＝2.92，不对称回波的情况下TR＝2.94。可见，经过优化，TR明显缩短了。另外，本发明的优化的TR略低于采用不对称回波技术的标准序列的TR。
另外，本发明方法经过试验验证，具有缩短TR的有益效果。如图4所示为水膜实验，(a)为优化后的bssfp采集的图像，TR＝3.2ms；(b)为未优化的bssfp采集的图像，TR＝3.6ms。如图5所示为人体实验，(a)为优化后的bssfp采集的图像，TR＝2.7ms；(b)为未优化的bssfp采集的图像TR＝3.16ms。
实施例2
本发明优化后的序列中也可以加入不对称回波技术，进一步缩短TR，即在前述实施例1的基础上加入不对称回波技术，作为实施例2，下面将进行详细说明。与实施例1相同之处，在此不再赘述。
本发明优化后的序列和不对称回波(asymmetric echo)结合 时，计算方法修改如下：
需要定义不对称因子(asymmetric factor)为实际采点数/全采点数。在步骤c中计算读出矩时，RO moment asymmetric＝RO moment×asymmetric factor。在步骤d中，频率编码方向散相矩仍为全采RO moment的一半，RO dephase moment＝RO moment/2。其余步骤不变。如果用了不对称回波，RO moment(读出矩)会变小，但是本文中仍采用全采的读出矩来计算。
采用上述方法采样，采样点距离不再均匀，因此数据需要采用网格化(grdding)方法将其插值为均匀采样的数据[1,2]。
如果沿读出方向有平移，平移距离为△s，需进行校正，具体方法如下：
根据傅立叶移位原理，对K空间数据乘以一个线性相位，相位斜率为
Figure PCTCN2017099430-appb-000008
Gramp为采用ramp sampling方法后的读出梯度幅值，Δt为采样点的时间间隔，相当于将图像在读出方向平移-△s；
对相位校正后的数据做regridding；
将网格化后的数据乘以一个线性相位，相位斜率为
Figure PCTCN2017099430-appb-000009
Gstd为未采样ramp sampling方法时所需的读出梯度幅值，相当于将图像在读出方向移回初始位置。
后续的图像重建方案等同于未做ramp sampling的重建方案。
通过上述在本发明优化后的序列中可以加入不对称回波技术，可以进一步缩短TR。
另外，本发明方法也可用于多回波(multi-echo)的梯度回波序列，减少回波间距(echo spacing)，可以用在笛卡尔坐标的扫描轨迹，也适用于径向轨迹的采样方法。此方法在二维成像或非各向同性成像中效果最明显。
此外，本发明还提供了一种计算机可读介质，该计算机可读介质具有存储在其中的程序，该程序是计算机可执行的，以使计算机执行包括上 述各步骤的处理。
本发明还提供了一种用于上述方法的优化平衡稳态自由进动序列bSSFP的装置。该装置包括了不同模块用于实现上述方法中提及的各个步骤。
例如，该装置可以包括：序列梯度设计模块，对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；最大梯度幅值和爬升率确定模块，在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；数据采样模块，在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；以及最大幅值和最大爬升率重新分配模块，重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的时间降到最短，从而获得最短的重复时间TR。
优选地，所述数据采样模块在采集数据时保持采样点时间间隔不变。
优选地，所述数据采样模块还包括以下模块：
读出矩计算模块，根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
Figure PCTCN2017099430-appb-000010
其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
读出梯度爬升时间获得模块，为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
读出时间获得模块，根据步骤c2的计算结果，获得读出时间：2*Ramp Time+Flattop time；
带宽调整模块，调整读出的带宽，使得采集点数不小于原始矩阵大小。
最大幅值和最大爬升率重新分配模块中，可以采用穷举法获得三个轴梯度性能的最优分配方案，还包括如下模块：
输入参数设置模块，设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
中间参数设置模块，设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
散相时间计算模块，设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
Figure PCTCN2017099430-appb-000011
Dro＝Mro/Gro+Gro×Rro；
Dpe＝Mpe/Gpe+Gpe×Rpe；
Mss＝Grf×(Drf+Grf×Rss)/2；
Dss＝Mss/Gss+(Gss+Grf)×Rss；
Ddephase＝max(Dro,Dpe,Dss)；
最大梯度幅值和最大爬升率获得模块，通过散相时间计算模块得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
散相时间计算模块中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
在另外的实施例中，所述装置还与不对称回波结合。与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。在计算读出矩时，计算公式为RO moment asymmetric＝RO moment×asymmetric factor。
在另外的实施例中，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
该装置还可以包括处理模块，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
显然，本技术领域中的普通技术人员应当认识到，以上的实施例仅是用来说明本发明，而并非用作为对本发明的限定，只要在本发明的实质精神范围内，对以上所述实施例的变化、变型都将落在本发明的权利要求书范围内。

Claims (30)

  1. 一种优化平衡稳态自由进动序列bSSFP的方法，其特征在于，包括以下步骤：
    步骤a：对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；
    步骤b：在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；
    步骤c：在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；以及
    步骤d：重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的时间降到最短，从而获得最短的重复时间TR。
  2. 根据权利要求1所述的方法，其特征在于，步骤c中，采集数据时保持采样点时间间隔不变。
  3. 根据权利要求2所述的方法，其特征在于，步骤c还包括以下步骤：
    步骤c1：根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
    Figure PCTCN2017099430-appb-100001
    其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
    步骤c2：为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
    读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
    读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
    步骤c3：根据步骤c2的计算结果，获得读出时间：2*Ramp Time+ Flattop time；
    步骤c4：调整读出的带宽，使得采集点数不小于原始矩阵大小。
  4. 根据权利要求1所述的方法，其特征在于，步骤d中，采用穷举法获得三个轴梯度性能的最优分配方案，包括如下步骤：
    步骤d1：设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
    步骤d2：设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
    步骤d3：设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
    Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
    Figure PCTCN2017099430-appb-100002
    Dro＝Mro/Gro+Gro×Rro；
    Dpe＝Mpe/Gpe+Gpe×Rpe；
    Mss＝Grf×(Drf+Grf×Rss)/2；
    Dss＝Mss/Gss+(Gss+Grf)×Rss；
    Ddephase＝max(Dro,Dpe,Dss)；
    步骤d4：通过步骤c3得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
  5. 根据权利要求4所述的方法，其特征在于，步骤d3中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
  6. 根据权利要求1所述的方法，其特征在于，所述方法还与不对称回波结合。
  7. 根据权利要求6所述的方法，其特征在于，与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。
  8. 根据权利要求7所述的方法，其特征在于，在计算读出矩时，计算公式为RO moment asymmetric＝RO moment×asymmetric factor。
  9. 根据权利要求1～8中任一项所述的方法，其特征在于，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
  10. 根据权利要求1～8中任一项所述的方法，其特征在于，步骤d之后还包括步骤e，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
  11. 一种优化平衡稳态自由进动序列bSSFP的装置，其特征在于，包括：
    序列梯度设计模块，对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；
    最大梯度幅值和爬升率确定模块，在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；
    数据采样模块，在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；以及
    最大幅值和最大爬升率重新分配模块，重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的时间降到最短，从而获得最短的重复时间TR。
  12. 根据权利要求11所述的装置，其特征在于，所述数据采样模块在采集数据时保持采样点时间间隔不变。
  13. 根据权利要求12所述的装置，其特征在于，所述数据采样模块还包括以下模块：
    读出矩计算模块，根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
    Figure PCTCN2017099430-appb-100003
    其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
    读出梯度爬升时间获得模块，为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
    读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
    读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
    读出时间获得模块，根据步骤c2的计算结果，获得读出时间：2*Ramp Time+Flattop time；
    带宽调整模块，调整读出的带宽，使得采集点数不小于原始矩阵大小。
  14. 根据权利要求11所述的装置，其特征在于，最大幅值和最大爬升率重新分配模块中，采用穷举法获得三个轴梯度性能的最优分配方案，还包括如下模块：
    输入参数设置模块，设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
    中间参数设置模块，设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss 分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
    散相时间计算模块，设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
    Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
    Figure PCTCN2017099430-appb-100004
    Dro＝Mro/Gro+Gro×Rro；
    Dpe＝Mpe/Gpe+Gpe×Rpe；
    Mss＝Grf×(Drf+Grf×Rss)/2；
    Dss＝Mss/Gss+(Gss+Grf)×Rss；
    Ddephase＝max(Dro,Dpe,Dss)；
    最大梯度幅值和最大爬升率获得模块，通过散相时间计算模块得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
  15. 根据权利要求14所述的装置，其特征在于，散相时间计算模块中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
  16. 根据权利要求11所述的装置，其特征在于，所述装置还与不对称回波结合。
  17. 根据权利要求16所述的装置，其特征在于，与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。
  18. 根据权利要求17所述的装置，其特征在于，在计算读出矩时，计算公式为RO moment asymmetric＝RO moment×asymmetric factor。
  19. 根据权利要求11～18中任一项所述的装置，其特征在于，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
  20. 根据权利要求11～18中任一项所述的装置，其特征在于，还包括处理模块，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
  21. 一种计算机可读介质，该计算机可读介质具有存储在其中的程序，该程序是计算机可执行的以使计算机执行包括以下步骤的处理：
    步骤a：对平衡稳态自由进动序列进行序列梯度设计，使得整个频率编码梯度时间内不和其他两轴的梯度即层选梯度及相位编码梯度发生重叠；
    步骤b：在频率编码梯度即读出梯度中使用系统允许的最大梯度幅值和爬升率；
    步骤c：在保持图像分辨率不变的情况下，在频率编码梯度的整个读出梯度时间内进行数据采样，即数据采集覆盖读出梯度的上升沿、平台期和下降沿；以及
    步骤d：重新分配频率编码梯度、层选梯度及相位编码梯度这三个梯度轴的最大幅值和最大爬升率，使得矩最大的梯度轴具有最高的梯度最大幅值和最大爬升率；这样使三个梯度需要的时间降到最短，从而获得最短的重复时间TR。
  22. 根据权利要求21所述的计算机可读介质，其特征在于，步骤c中，采集数据时保持采样点时间间隔不变。
  23. 根据权利要求22所述的计算机可读介质，其特征在于，步骤c还包括以下步骤：
    步骤c1：根据平衡稳态自由进动序列中设定的读出方向的视野和读出矩阵大小计算读出矩
    Figure PCTCN2017099430-appb-100005
    其中，FOV为读出方向的视野，单位为毫米，Readout data points为读出矩阵大小；
    步骤c2：为保证读出时间为最短，将读出梯度幅值和爬升率设置为系统允许的最大梯度幅值Gmax和最大爬升率，令Rmax＝1/最大爬升率，则获得
    读出梯度爬升时间Ramp Time(us)＝Gmax(mT/m)*Rmax((m*us)/mT)，
    读出梯度的平台期时间Flattop Time＝RO moment/Gmax-Ramp Time；
    步骤c3：根据步骤c2的计算结果，获得读出时间：2*Ramp Time+Flattop time；
    步骤c4：调整读出的带宽，使得采集点数不小于原始矩阵大小。
  24. 根据权利要求21所述的计算机可读介质，其特征在于，步骤d中，采用穷举法获得三个轴梯度性能的最优分配方案，包括如下步骤：
    步骤d1：设置输入参数Mro，Mpe，Gmax，Rmax，Grf和Drf，其中，Mro表示频率编码方向散相矩(RO dephase moment＝RO moment/2)，Mpe表示相位编码方向最大的相位矩(Maximum PE moment)，Gmax表示系统允许的最大梯度幅值，Rmax表示1/系统允许的最大爬升率，Grf表示层选梯度幅值，Drf表示层选梯度的平台期时间；
    步骤d2：设置中间参数，包括：Gro，Gpe，Gss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度幅值；Rro，Rpe，Rss分别为读出散相梯度、相位编码梯度和层选梯度所允许的最大梯度爬升率的倒数；Dro为频率编码散相时间；Dpe为相位编码所需时间；Dss为层选梯度下降沿和层选梯度回聚所需时间的和；
    步骤d3：设置Gro，Rro的搜索范围，对于处于搜索范围内的每个点，分别计算所需的散相时间(Ddephase)，计算方式如下：
    Gpe＝Gro×Mpe/Mro；Rpe＝Rro×Mro/Mpe；
    Figure PCTCN2017099430-appb-100006
    Dro＝Mro/Gro+Gro×Rro；
    Dpe＝Mpe/Gpe+Gpe×Rpe；
    Mss＝Grf×(Drf+Grf×Rss)/2；
    Dss＝Mss/Gss+(Gss+Grf)×Rss；
    Ddephase＝max(Dro,Dpe,Dss)；
    步骤d4：通过步骤c3得到一系列散相时间，选出其中最小的散相时间，所述最小的散相时间对应的Gro，Gpe，Gss，Rro，Rpe，Rss即为三个梯度轴最大梯度幅值和最大爬升率。
  25. 根据权利要求24所述的计算机可读介质，其特征在于，步骤d3中，设置Gro Range＝[Gmax/2,Gmax],Rro Range＝[Rmax,Rmax*3]，搜索间隔为0.2mT/m。
  26. 根据权利要求21所述的计算机可读介质，其特征在于，所述方法还与不对称回波结合。
  27. 根据权利要求26所述的计算机可读介质，其特征在于，与不对称回波结合的情况下，需要定义不对称因子asymmetric factor为实际采点数/全采点数。
  28. 根据权利要求27所述的计算机可读介质，其特征在于，在计算读出矩时，计算公式为RO moment asymmetric＝RO moment×asymmetric factor。
  29. 根据权利要求21～28中任一项所述的计算机可读介质，其特征在于，所述平衡稳态自由进动序列可以替换成多回波的梯度回波序列。
  30. 根据权利要求21～28中任一项所述的计算机可读介质，其特征在于，步骤d之后还包括步骤e，根据bSSFP序列的对称性，对读出梯度后面的频率编码回聚梯度、层选散相梯度和相位编码恢复梯度也进行相同的梯度性能分配和处理。
PCT/CN2017/099430 2017-08-29 2017-08-29 一种优化平衡稳态自由进动序列的方法与装置 WO2019041123A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
PCT/CN2017/099430 WO2019041123A1 (zh) 2017-08-29 2017-08-29 一种优化平衡稳态自由进动序列的方法与装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/CN2017/099430 WO2019041123A1 (zh) 2017-08-29 2017-08-29 一种优化平衡稳态自由进动序列的方法与装置

Publications (1)

Publication Number Publication Date
WO2019041123A1 true WO2019041123A1 (zh) 2019-03-07

Family

ID=65524555

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2017/099430 WO2019041123A1 (zh) 2017-08-29 2017-08-29 一种优化平衡稳态自由进动序列的方法与装置

Country Status (1)

Country Link
WO (1) WO2019041123A1 (zh)

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1711965A (zh) * 2004-06-16 2005-12-28 西门子公司 在磁共振断层造影中显示具很短t2弛豫时间组织的方法
CN101299061A (zh) * 2007-05-03 2008-11-05 巴塞尔大学 用于与磁化率相关的磁场畸变的检测和成像的磁共振方法
US8283924B2 (en) * 2008-10-20 2012-10-09 University Of Southern California Eddy-current artifact reduction in balanced steady-state free precession magnetic resonance imaging
CN104297708A (zh) * 2014-10-22 2015-01-21 浙江大学 磁共振成像系统中的旋转倾斜梯度线圈组件
CN104523274A (zh) * 2014-12-25 2015-04-22 中国科学院深圳先进技术研究院 一种利用稳态自由进动序列的磁共振成像方法
US20160116559A1 (en) * 2014-10-24 2016-04-28 The General Hospital Systems and methods for steady-state magnetic resonance fingerprinting

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1711965A (zh) * 2004-06-16 2005-12-28 西门子公司 在磁共振断层造影中显示具很短t2弛豫时间组织的方法
CN101299061A (zh) * 2007-05-03 2008-11-05 巴塞尔大学 用于与磁化率相关的磁场畸变的检测和成像的磁共振方法
US8283924B2 (en) * 2008-10-20 2012-10-09 University Of Southern California Eddy-current artifact reduction in balanced steady-state free precession magnetic resonance imaging
CN104297708A (zh) * 2014-10-22 2015-01-21 浙江大学 磁共振成像系统中的旋转倾斜梯度线圈组件
US20160116559A1 (en) * 2014-10-24 2016-04-28 The General Hospital Systems and methods for steady-state magnetic resonance fingerprinting
CN104523274A (zh) * 2014-12-25 2015-04-22 中国科学院深圳先进技术研究院 一种利用稳态自由进动序列的磁共振成像方法

Similar Documents

Publication Publication Date Title
Polak et al. Highly‐accelerated volumetric brain examination using optimized wave‐CAIPI encoding
Robison et al. Correction of B0 eddy current effects in spiral MRI
Polak et al. Nonlinear dipole inversion (NDI) enables robust quantitative susceptibility mapping (QSM)
Tsai et al. Reduced aliasing artifacts using variable‐density k‐space sampling trajectories
US8120358B2 (en) Magnetic resonance imaging with high spatial and temporal resolution
Brodsky et al. Generalized k‐space decomposition with chemical shift correction for non‐cartesian water‐fat imaging
JP7075420B2 (ja) 可変コントラストのスタック・オブ・スター収集を使用したｍｒイメージング
US9588208B2 (en) Methods, systems and apparatuses for rapid segmented, accelerated, and simultaneous multi-slice echo planar imaging
Stehning et al. Fast isotropic volumetric coronary MR angiography using free‐breathing 3D radial balanced FFE acquisition
CN105143904A (zh) 用于多信道发射机的多频带射频/磁共振成像脉冲设计
KR20020086289A (ko) 자기장의 불균일성을 고려하여 위상 부호화된 개별영상으로부터 화학적 이동이 상이한 스핀 스펙트럼을추출하기 위한 방법 및 장치
Benkert et al. Dynamically phase‐cycled radial balanced SSFP imaging for efficient banding removal
Schirda et al. Rosette spectroscopic imaging: optimal parameters for alias‐free, high sensitivity spectroscopic imaging
CN107110940A (zh) 用于mri径向或螺旋成像的方法
Robison et al. Three‐dimensional ultrashort echo‐time imaging using a FLORET trajectory
JPWO2011007691A1 (ja) 磁気共鳴イメージング装置および磁気共鳴イメージング方法
CN107728090B (zh) 一种优化平衡稳态自由进动序列的方法与装置
Moran et al. High‐resolution 3D radial bSSFP with IDEAL
Sprenger et al. NeuroMix—A single‐scan brain exam
US20150137811A1 (en) System and method for magnetic resonance imaging using highly accelerated projection imaging
JP5846450B2 (ja) 磁気共鳴イメージング装置、計測空間座標補正方法、及び、画像再構成方法
van der Plas et al. Time‐encoded golden angle radial arterial spin labeling: simultaneous acquisition of angiography and perfusion data
CN103513204B (zh) 一种磁共振成像中k空间数据的轨迹校正方法和装置
US20130076355A1 (en) Fast, Low Energy Deposition and Homogeneous T2 Weighted Variable Amplitude PSIF (T2 VAPSIF) Imaging in the Presence of B0inhomogeneities
WO2021123384A1 (en) A method and controller for magnetic resonance imaging

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 17923934

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205A DATED 23/09/2020)

122 Ep: pct application non-entry in european phase

Ref document number: 17923934

Country of ref document: EP

Kind code of ref document: A1