用于三维动态磁共振成像的球形k空间采集方法和装置与流程

本申请涉及三维动态磁共振成像技术领域，尤其涉及一种用于三维动态磁共振成像的球形k空间采集方法和装置。

背景技术：

磁共振成像(Magnetic Resonance Imaging，MRI)是利用核磁共振现象，激发人体中的自旋核，然后接收自旋核子释放的电磁信号，重建出人体组织图像信息，具有无辐射，多对比度成像和软组织对比度高等优点，目前已经成为临床医学检查的一种重要工具。磁共振成像过程中采集得到的数据叫k空间数据，所有数据整体组成k空间，所用的采集轨迹叫k空间轨迹，对k空间数据进行重建可以得到磁共振图像。动态磁共振成像是利用磁共振成像方法对人体组织器官的动态生理过程(如心脏跳动，药物代谢等)进行跟踪成像的一种技术，该技术通过对一特定的成像空间进行连续重复的扫描，得到一系列与时间相关的k空间数据，通过对这些数据筛选重建可以得到一组随时间变化的动态图像，通过对动态图像做数据分析可以获得一系列反映病变发生、发展过程中生物学和病理生理学信息的定量或半定量参数。动态磁共振成像技术已经广泛应用于如心脏电影磁共振成像(Cardiac Cine MRI)，动态对比度增强磁共振成像(Dynamic Contrast Enhanced MRI，DCE-MRI)等领域。

在这类动态磁共振成像技术中，由于受到三维k空间数据采集时间的限制，时间分辨率较低。在近似周期性运动的器官成像，如心脏电影三维成像中，虽然可以采用多运动周期重复采集和回顾性重建方法来提高时间分辨率，但受到采集方案的限制，很难保证任意时间窗内采集到的k空间数据都近似均匀分布。

目前的三维动态磁共振成像采集方法都是重复采集部分或者全部k空间数据得到的，比如三维笛卡尔采集和三维径向采集。三维笛卡尔采集方法是通过选层梯度编码和相位编码实现两个维度上的编码，然后通过频率编码实现第三维度上的数据采集，以此实现对三维k空间的填充。三维径向采集方法通过对选层，相位和频率三个方向同时编码采集实现对三维k空间的填充。对于这两类方法，重建一组三维动态图像需要重复采集全部或大量的三维k空间数据，因此时间分辨率较低，且图像重建时需要选择特定时间窗内的k空间数据，采集完一个时间窗再进行下一个时间窗采集，以保证k空间的磁共振数据均匀分布，因此图像重建的自由性受到限制。此外还有笛卡尔和径向混合采集方法，该方法在平面内采用径向采集，在选层梯度方向采用笛卡尔采集，虽然通过欠采样采集可以一定程度上降低采集时间，但是k空间数据均匀性仍然受到某一维度编码方案的限制，动态成像的时间分辨率难以提高。此外，现有动态磁共振成像技术大多是通过对后期筛选重组的数据重建来实现动态成像的，因此一般很难获得连续的动态图像。

技术实现要素：

为解决现有技术中的上述问题，本申请的一个目的在于提出一种适用于三维动态磁共振成像的球形k空间采集方法和装置，可以实现三维k空间数据的连续采集，使采集到的k空间数据在任意时间窗内都近似均匀分布，实现高时间分辨率的动态连续三维磁共振成像。

为达到上述目的，本申请实施例提出的适用于三维动态磁共振成像的球形k空间采集方法，包括：建立球形的k空间模型，确定模型中回波信号的采集轨迹，其中，所述回波信号的采集轨迹为：在所述球形内，由球心到球面的半径状的半回波信号，或者，过球心的直径状的全回波信号，回波信号采集轨迹的仰角与方位角由二维黄金分割比例计算得到；根据所述采集轨迹确定磁共振扫描的时间序列，并计算磁共振成像系统所需施加磁场的编码梯度；根据所述时间序列和编码梯度设置磁共振成像系统，并采集符合所述采集轨迹的k空间数据填充所述球形k空间。

为达到上述目的，本申请实施例提出的适用于三维动态磁共振成像的球形k空间采集装置，包括：建模模块，用于建立球形的k空间模型，确定模型中回波信号的采集轨迹，其中，所述回波信号的采集轨迹为：在所述球形内，采集由球心到球面的半径状的半回波信号，或者，过球心的直径状的全回波信号，回波信号采集轨迹的仰角与方位角由二维黄金分割比例计算得到；计算模块，用于根据所述采集轨迹确定磁共振扫描的时间序列，并计算磁共振成像系统所需施加磁场的编码梯度；采集模块，用于根据所述时间序列和编码梯度设置磁共振成像系统，并采集符合所述采集轨迹的k空间数据填充所述球形k空间。

本申请实施例提供的技术方案，能够在任一时刻采集得到近似均匀的k空间数据，进而可以实现连续动态的三维成像，实现高时间分辨率的动态成像，具体优点如下：

1、均匀性方面：可以实现任意长度时间窗内、任意位置时间窗内、任意组合时间窗内采集的数据，在三维球形k空间内都近似均匀分布。因此，在进行图像重建时数据的选择具有较高的自由性，通过适当的图像重建方法可以获得具有较高时间分辨率的动态影像。

2、适用范围方面：分别提出了适用于全回波采集和半回波采集的k空间采集方法，都可以实现近似均匀的k空间分布，因此既可应用在采集全回波信号的序列(如稳态自由进动序列等)，也可应用在采集半回波信号的序列(如超短回波时间序列等)，应用范围广。

3、应用方面：有利于对动态生理过程进行三维磁共振成像(如动态对比度增强磁共振成像，DCE-MRI)，还有利于对近似周期性运动的器官进行三维磁共振成像(如心脏，胃，肺等)。

本申请附加的方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本申请的实践了解到。

附图说明

为了更清楚地说明本申请实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本申请一实施例的适用于三维动态磁共振成像的球形k空间采集方法的流程示意图；

图2是本申请一实施例的采集半回波信号的回波线示意图；

图3是本申请一实施例的采集全回波信号的回波线示意图；

图4是本申请一实施例的球形k空间数据采集的流程示意图；

图5是本申请一实施例采集到的球形k空间的数据分布情况示意图；

图6是本申请一实施例的动态磁共振图像重建中k空间时间窗的三种选择模式的示意图；

图7是本申请一具体实施例的相应不同的时间窗选择模式下的k空间数据均匀性对比示意图；

图8是本申请一实施例的适用于三维动态磁共振成像的球形k空间采集装置的结构示意图；

图9是本申请另一实施例的适用于三维动态磁共振成像的球形k空间采集装置的结构示意图。

具体实施方式

本申请实施例提供一种用于三维动态磁共振成像的球形k空间采集方法和装置。

为了使本技术领域的人员更好地理解本申请中的技术方案，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都应当属于本申请保护的范围。

图1是本申请一实施例提出的适用于三维动态磁共振成像的球形k空间采集方法的流程示意图。如图1所示，该方法包括：

步骤101，建立球形的k空间模型，确定模型中回波信号的采集轨迹，其中，所述回波信号的采集轨迹为：在所述球形内，由球心到球面的半径状的半回波信号，或者，过球心的直径状的全回波信号，回波信号采集轨迹的仰角θ与方位角由二维黄金分割比例计算得到。

步骤102，根据所述采集轨迹确定磁共振扫描的时间序列，并计算磁共振成像系统所需施加磁场的编码梯度。

步骤103，根据所述时间序列和编码梯度设置磁共振成像系统，并采集符合所述采集轨迹的k空间数据填充所述球形k空间。

具体地，本实施例采集的三维k空间整体是一个球体。在三维球形k空间中，每次采集回波的方向即采集轨迹的仰角θ和方位角可以分别由二维黄金分割比例γ₁，γ₂计算得到，然后根据所述采集轨迹确定磁共振扫描的时间序列，并计算磁共振成像系统所需施加磁场的编码梯度，然后根据扫描时间序列及编码梯度编写软件程序输入到磁共振成像系统，即可采集符合上述采集轨迹的回波信号实现对三维球形k空间的填充。这样可以使采集到的k空间数据更均匀，提高动态成像的时间分辨率。

进一步地，二维黄金分割比例系数为γ₁和γ₂，是根据广义斐波那契数列的特征向量得到的，取小数点后四位后，γ₁≈0.6823和γ₂≈0.4656。所述回波信号的仰角θ依据γ₁计算，方向角依据γ₂计算，或者，仰角θ依据γ₂计算，方向角依据γ₁计算。

基于以上，本发明实施例可以实现三维k空间的连续采集，并且在任意采集时间窗内都可以得到近似均匀的k空间分布，大大提高了重建数据时的选择自由性，同时结合欠采样重建技术可以实现高时间分辨率的动态磁共振图像重建，适用于临床的应用需求。

进一步地，本发明分别提出了基于全回波采集和基于半回波采集两种不同的采集实施例。本发明实施例中的三维动态磁共振成像的球形k空间采集轨迹如图2和图3所示，整体采集空间为一个球形，图中θ_n和分别表示第n次采集回波信号线的仰角和方位角。为了得到任意层面的空间信息，磁共振成像系统在x，y，z三个坐标方向均使用物理梯度磁场。本发明提出的球形k空间采集轨迹在实际应用中可以与磁共振成像系统中的三个物理梯度磁场任意旋转对应，都属于本专利的保护范围。

在本申请的一个实施例中，针对半回波采集的信号(如采集FID信号：自由感应衰减(free induction decay,FID)是磁共振成像(MRI)中最简单的信号形式)，如图2中所示，每次采集回波线从三维球形k空间中心开始，沿球的径向方向向外采集。

在本申请的另一个实施例中，球形k空间采集的全回波信号如图3中球体直径所示，每次采集的回波从球面一点开始穿过球心采集，采集轨迹为球的一根直径。

无论是全回波采集还是半回波采集，在球形k空间中，每次采集回波的方向均可由θ和两个参数确定。因此，设计适用于高时间分辨率的三维动态磁共振球形k空间采集方案只需优化θ和这两个参数。因此如果要实现高时间分辨率的三维k空间采集，只需优化θ和这两个参数。

假定三维k空间的采集矩阵为：R_x×R_y×R_z，其中R_x，R_y，R_z分别表示x，y，z方向的编码数。在三维球形k空间中，R_x＝R_y＝R_z＝2R。针对半回波采集和全回波采集，可分别采用不同的采集方案。

根据本申请的一个实施例，如图4所示，回波信号的采集轨迹通过以下步骤计算：

1)设定第n次采集的回波信号在球坐标系下对应的参数为：仰角θ_n与方位角初始回波的编号为n＝n₀。

2)根据二维黄金分割比例系数，计算第n次采集的回波信号的仰角θ_n与方位角

具体地，当所述回波信号是由球心出发的半回波信号时，通过下式计算其仰角和方位角：

θ_n＝arccos(2mod(n·γ1_,1)-1)

当所述回波信号是直径状的全回波信号时，通过下式计算其仰角和方位角：

θ_n＝arccos(mod(n·γ₁,1))

或者，根据本申请另一实施例，通过下式计算其仰角和方位角：

θ_n＝arccos(2mod(n·γ₁,1)-1)

其中，arccos(i)为求i的反余弦值，mod(a,b)为求a/b的余数。

3)根据θ_n和计算第n次采集的回波信号上的各个采集点在笛卡尔坐标系下的坐标x_n,y_n,z_n：

z_n＝r·cos(θ_n)

对于半回波信号：r＝0,1,...,R，其中r从0到R表示半回波线上各个采集点；

对于全回波信号：r＝-R,-R+1,...,R，其中r从-R到R表示全回波线上各个采集点。

4)根据预设截止条件判断采集是否结束，即第n次采集的回波信号是否为采集的最后一条回波信号。

其中，预设的截止条件由人为设定，例如是采集时间为20分钟，或者，采集10万次回波信号等。

5)若是，则采集结束，否则使n＝n+1，重复1)到4)直至采集结束。

根据本申请的一个实施例，上述回波采集的所述初始回波的编号n₀可以为任意自然数，能够从任意位置一条轨迹开始采集均匀的k空间数据对球形k空间进行填充。

根据本申请的一个实施例，计算磁共振成像系统所需施加磁场的编码梯度，可以根据第n次采集轨迹上各个采集点坐标x_n，y_n，z_n计算第n次采集时所需施加的三个方向磁场的编码梯度。在不考虑实际硬件条件等限制时，第n次采集时所施加的编码梯度G_x,y,_z与x_n，y_n，z_n的关系可由下式描述：

$G_{x,y,z}=\frac{k_{x,y,z}}{\mathrm{\γ}t}$

其中，k_x,y,_z为第n次采集轨迹在k空间内的坐标，对应于x_n，y_n和z_n；γ为旋磁比，t梯度施加的时间。

在所述根据所述编码采集相应的k空间数据之后，对所述k空间数据进行处理，重建得到连续的动态图像。

本申请实施例可以实现在任意时间窗内采集的数据在三维球形k空间内都近似均匀分布，因此有利于提高动态成像的时间分辨率。下面通过具体模拟实验评估采用本申请实施例的方法采集到的k空间数据的均匀性。

以k空间分辨率为100×100×100，k空间采集的是半回波信号为例进行阐述，全回波信号的情况与此类似，不再赘述。在三维球形k空间内连续采集500次半回波，采集的球形k空间的数据在球面上的分布情况如图5所示。

为了评估k空间数据均匀性，可计算球面上每个点与之邻域内所有点的平均距离，然后统计所有这些距离值的标准差，标准差越趋近于0，表示每个点与之邻域内点之间的距离越相近，即表示分布越均匀。每个点的邻域范围为在球面上以该点为中心，立体角为Ω＝4π²/T的球面圆域，其中T为采集总次数。图6所示为动态磁共振图像重建中k空间时间窗的三种选择模式，其中(a)为不同长度的时间窗，(b)为不同位置的时间窗，(c)为不同组合的时间窗。图7给出了相应不同的时间窗选择模式下的k空间数据均匀性对比，包括不同的采集次数(a)、不同的时间窗位置(b)以及不同的时间窗组合次数(c)三种条件下三维球形k空间数据分布均匀性对比，其中(b)和(c)的采集次数为5000。

由以上评估结果可以看出，本发明实施例的采集方法采集的k空间数据在三维球形k空间内均匀性很好，在三种不同时间窗选择模式下的统计标准差均接近于0，证明采用本发明的方法得到的k空间数据具有较好的空间均匀性。

综上所述，本发明实施例的球形k空间采集方法具备如下优点：

1、均匀性方面：可以实现任意长度时间窗内、任意位置时间窗内、任意组合时间窗内采集的数据，在三维球形k空间内都近似均匀分布。因此，在进行图像重建时数据的选择具有较高的自由性，通过适当的图像重建方法可以获得具有较高时间分辨率的动态影像。

2、适用范围方面：分别提出了适用于全回波采集和半回波采集的k空间采集方法，都可以实现近似均匀的k空间分布，因此既可应用在采集全回波信号的序列(如稳态自由进动序列等)，也可应用在采集半回波信号的序列(如超短回波时间序列等)，应用范围广。

3、应用方面：有利于对动态生理过程进行三维磁共振成像(如动态对比度增强磁共振成像，DCE-MRI)，还有利于对近似周期性运动的器官进行三维磁共振成像(如心脏，胃，肺等)。

基于同一发明构思，本申请实施例还提供了一种适用于三维动态磁共振成像的球形k空间采集装置，可以用于实现上述实施例所描述的方法，如下面的实施例所述。由于适用于三维动态磁共振成像的球形k空间采集装置解决问题的原理与适用于三维动态磁共振成像的球形k空间采集方法相似，因此适用于三维动态磁共振成像的球形k空间采集装置的实施可以参见适用于三维动态磁共振成像的球形k空间采集装置的实施，重复之处不再赘述。以下所使用的，术语“单元”或者“模块”可以实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现，但是硬件，或者软件和硬件的组合的实现也是可能并被构想的。

图8是本申请一实施例的适用于三维动态磁共振成像的球形k空间采集装置的结构示意图。本实施例的装置可以为实现相应功能的逻辑部件构成，也可以为运行有相应功能软件的电子设备。如图8所示，该适用于三维动态磁共振成像的球形k空间采集装置包括：建模模块10、计算模块20和采集模块30。

具体地，建模模块10用于建立球形的k空间模型，确定模型中回波信号的采集轨迹，其中，所述回波信号的采集轨迹为：在所述球形内，采集由球心到球面的半径状的半回波信号，或者，过球心的直径状的全回波信号，回波信号采集轨迹的仰角与方位角由二维黄金分割比例计算得到。

计算模块20用于根据所述采集轨迹确定磁共振扫描的时间序列，并计算磁共振成像系统所需施加磁场的编码梯度；

采集模块30用于根据所述时间序列和编码梯度设置磁共振成像系统，并采集符合所述采集轨迹的k空间数据填充所述球形k空间。

图9所示是本申请另一实施例的适用于三维动态磁共振成像的球形k空间采集装置的结构示意图。如图9所示，在图8的基础上，建模模块10还包括第一计算单元11、第二计算单元12和第三计算单元13。

具体地，所述建模模块10还用于通过以下步骤确定回波信号的采集轨迹：

1)设第n次采集的回波信号在球坐标系下对应的参数为：仰角θ_n与方位角初始回波的编号为n＝n₀；

2)根据二维黄金分割比例系数，计算第n次采集的回波信号的仰角θ_n与方位角

3)根据θ_n和计算第n次采集的回波信号上的各个采集点在笛卡尔坐标系下的坐标x_n,y_n,z_n；

4)根据预设截止条件判断采集是否结束，即第n次采集的回波信号是否为采集的最后一条回波信号；

5)若是，则采集结束，否则使n＝n+1，重复1)到4)直至采集结束。

其中，所述二维黄金分割比例系数为γ₁和γ₂，γ₁≈0.6823，γ₂≈0.4656；所述回波信号的仰角θ依据γ₁计算，方向角依据γ₂计算，或者，仰角θ依据γ₂计算，方向角依据γ₁计算。

其中，所述初始回波的编号n₀可以为任意自然数。

建模模块10还包括第一计算单元11，用于当所述回波信号是由球心出发的半回波信号时，通过下式计算其仰角和方位角：

θ_n＝arccos(2mod(n·γ₁,1)-1)

其中，arccos(i)为求i的反余弦值，mod(a,b)为求a/b的余数。

建模模块10还包括第二计算单元12，用于当所述回波信号是直径状的全回波信号时，通过下式计算其仰角和方位角：

θ_n＝arccos(mod(n·γ₁,1))

其中，arccos(i)为求i的反余弦值，mod(a,b)为求a/b的余数。

建模模块10还包括第三计算单元13，用于当所述回波信号是直径状的全回波信号时，通过下式计算其仰角和方位角：

θ_n＝arccos(2mod(n·γ₁,1)-1)

其中，arccos(i)为求i的反余弦值，mod(a,b)为求a/b的余数。

本实施例的装置可以实现三维k空间数据的连续采集，提高三维动态磁共振成像的时间分辨率，具体体现为如下优点：

1、均匀性方面：可以实现任意长度时间窗内、任意位置时间窗内、任意组合时间窗内采集的数据，在三维球形k空间内都近似均匀分布。因此，在进行图像重建时数据的选择具有较高的自由性，通过方法可以获得具有较高时间分辨率的动态影像。

2、适用范围方面：分别提出了适用于全回波采集和半回波采集的k空间采集方法，都可以实现近似均匀的k空间分布，因此既可应用在采集全回波信号的序列(如稳态自由进动序列等)，也可应用在采集半回波信号的序列(如超短回波时间序列等)，应用范围广。

3、应用方面：有利于对动态生理过程进行三维磁共振成像(如动态对比度增强磁共振成像，DCE-MRI)，还有利于对近似周期性运动的器官进行三维磁共振成像(如心脏，胃，肺等)。

需要说明的是，在本申请的描述中，术语“第一”、“第二”等仅用于描述目的，而不能理解为指示或暗示相对重要性。此外，在本申请的描述中，除非另有说明，“多个”的含义是两个或两个以上。

流程图中或在此以其他方式描述的任何过程或方法描述可以被理解为，表示包括一个或更多个用于实现特定逻辑功能或过程的步骤的可执行指令的代码的模块、片段或部分，并且本申请的优选实施方式的范围包括另外的实现，其中可以不按所示出或讨论的顺序，包括根据所涉及的功能按基本同时的方式或按相反的顺序，来执行功能，这应被本申请的实施例所属技术领域的技术人员所理解。

应当理解，本申请的各部分可以用硬件、软件、固件或它们的组合来实现。在上述实施方式中，多个步骤或方法可以用存储在存储器中且由合适的指令执行系统执行的软件或固件来实现。例如，如果用硬件来实现，和在另一实施方式中一样，可用本领域公知的下列技术中的任一项或他们的组合来实现：具有用于对数据信号实现逻辑功能的逻辑门电路的离散逻辑电路，具有合适的组合逻辑门电路的专用集成电路，可编程门阵列(PGA)，现场可编程门阵列(FPGA)等。

本技术领域的普通技术人员可以理解实现上述实施例方法携带的全部或部分步骤是可以通过程序来指令相关的硬件完成，所述的程序可以存储于一种计算机可读存储介质中，该程序在执行时，包括方法实施例的步骤之一或其组合。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本申请的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。

尽管上面已经示出和描述了本申请的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本申请的限制，本领域的普通技术人员在本申请的范围内可以对上述实施例进行变化、修改、替换和变型。