本发明涉及磁共振成像领域，尤其涉及一种磁共振成像全自动预加重校正涡流的方法。
背景技术：
磁共振成像具有无损伤、软组织对比度高、任意方向断层等特点，目前已广泛用于医学临床诊断，而在磁共振成像过程中，系统会产生涡流，破坏主磁场和梯度场的均匀性，从而会在磁共振图像中产生伪影，尤其是快速自旋回波(fse)，弥散成像(dwi)等成像序列会对磁共振成像系统的涡流极为敏感。因此准确有效地定量并且消除磁共振成像系统中的涡流是影响磁共振图像质量最为关键的因素之一。
磁共振成像系统的涡流主要产生于施加成像梯度后，梯度电流会在梯度线圈和磁体中感应产生随时间呈自然指数形式衰减的感应电流，从而产生随时间衰减性的b0涡流场和梯度涡流场。涡流的物理形式为由单个或多个e指数组合而成的衰减函数，具体由几个e指数成分组成受系统特性、磁体和梯度线圈的设计生产、硬件设备安装等因素决定，不同的磁共振成像系统会有不同的涡流函数形式，即使在同一个磁共振成像系统中，在三个不同的磁体方向(x，y和z)的涡流衰减函数形式也会有不同。
目前校正涡流最常用且有效的两种方法：分别是利用自屏蔽线圈和预加重补偿，前者属于被动的方法，是通过在梯度线圈设计和制造上实现对涡流场的自屏蔽作用来抑制梯度线圈中的涡流，但是此方法依赖于线圈制造商的线圈设计和制造工艺，而且并不能完全消除涡流的影响，并且往往会带来负极性的涡流；而预加重补偿则是属于主动的方法，是通过对所施加的梯度波形进行预加重修正，施加与涡流幅值和衰减时间常数大小相同，但极性相反的梯度场，与系统涡流相互抵消，从而消除涡流影响。虽然该方法可以有效校正系统涡流，但在实际的磁共振成像系统调试中，涡流的极性、时间参数以及衰减指数的个数受系统的设备，硬件条件状况变化差异很大，传统的方法利用非线性拟合数学工具对涡流衰减信号无法保证能准确定义涡流衰减模型。因此实际的定量涡流调试依赖于工程师的主观经验判断，往往并不准确，无法给出合适的涡流模型，或者通常的做法是将e指数的成分个数设成为一个固定值(3个或4个)，这样做往往会造成过拟合，因此传统的涡流预加重方法无法真正做到自动化地准确调试校正系统涡流。
技术实现要素：
有鉴于此，本发明的目的是提供一种磁共振成像全自动预加重校正涡流的方法，以解决现有技术中的不足。
为了达到上述目的，本发明的目的是通过下述技术方案实现的：
提供一种磁共振成像全自动预加重校正涡流的方法，包括如下步骤：
a)分别在x、y或z方向独立施加一固定幅值位g和持续时间为d的梯度，经过延迟时间τi后，发射射频，并采集梯度回波信号；
b)利用对所有信号相位做定量计算得到各个延迟时间点τi的涡流值大小ci，得到涡流随时间衰减的曲线；
c)利用levenberg-marquardt算法，对定量得到的涡流衰减曲线ci做多e指数衰减模型的非线性拟合，并且对所有e指数的幅度值做权重正则化控制，防止出现指数项过多的过拟合发生,自动剔除掉多余的指数项后，将所有e指数衰减常数固定kj；
d)对涡流e指数幅度做svd线性拟合，根据拟合定量得到的对应于各个时间常数kj的e指数涡流幅值aj，并更新谱仪的涡流预加重参数，重新采集涡流信号；
e)反复迭代进行以上d过程，当涡流值小于某一设定的阈值t之后，程序判定为最终迭代收敛，涡流校正完毕。
上述磁共振成像全自动预加重校正涡流的方法，其中，步骤b)中利用对所有信号相位做定量计算得到各个延迟时间点τi的涡流值大小ci的，表示为公式[1]:
ci＝mi/2πγte[1]
其中，mi表示第i次采集的傅里叶变换后回波信号的斜率(单位：弧度/毫米)，γ表示物理常量旋磁比(＝0.267513弧度/毫秒/微特斯拉)，te表示回波时间(单位：微特斯拉)；
其中，步骤c)中所用的多e指数衰减的函数表示为公式[2]:
其中，n表示e指数成分的个数，τi表示梯度延迟时间，c(τi)表示第i个梯度延迟点采集计算得到的涡流大小，aj和kj分别表示第j个e指数衰减成分的幅值和时间常数；
其中，对公式[2]进行拟合时所做的正则化控制可以表示为公式[3]:
其中min{}表示最小化，λ表示正则化控制项的常数。
与已有技术相比，本发明的有益效果在于：
本发明在对涡流信号定量过程中使用了levenberg-marquardt非线性正则化最小化拟合，准确而且全自动，可以对磁共振成像系统中涡流衰减函数进行全自动地准确地判定，当涡流衰减形式得到正确的定义后，在迭代预加重校正过程中又将e指数的时间常数固定，有效地将非线性问题转化为线性模型定量，然后通过svd线性最小化方法对涡流幅值进行迭代校正，实现全过程完全自动化，不需要人工介入。
附图说明
构成本发明的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
附图1是本发明方法中定量得到的x方向涡流随时间衰减的数值定量曲线；
附图2是本发明方法中预加重校正前后系统x方向涡流曲线的效果对比图，实线所示为无预加重校正的系统x方向涡流定量曲线，虚线所示为有预加重校正的系统x方向涡流定量曲线；
附图3是本发明方法中涡流预加重校正前后水模fse图像验证的效果对比图：a为无涡流预加重校正的水模fse图像，b为有涡流预加重校正的水模fse图像。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
本发明在磁共振成像系统中的某一方向(x、y或z)施加一定量梯度，之后通过采集梯度回波信号来测量涡流值，通过测量不同延迟时间点的涡流值，得到涡流衰减曲线，然后利用levenberg-marquardt(列文伯格-马夸尔特法)的最小化算法对涡流衰减信号值做多e指数衰减的非线性拟合。并且在最小化过程中，对所有e指数的幅度值做权重正则化控制，防止出现指数项过多的过拟合发生。然后自动剔除掉多余的指数项后，将所有e指数衰减常数固定，将非线性模型转化为线性模型定量，对涡流e指数幅度做svd线性拟合，根据拟合定量得到的涡流幅值和时间常数更新谱仪的涡流预加重参数，重新采集涡流信号，反复迭代进行以上过程，迭代过程中所有e指数衰减常数固定，只做涡流e指数幅度做svd线性拟合，当涡流值小于某一设定的阈值之后，程序判定为最终迭代收敛，涡流校正完毕。全过程完全自动化，不需要人工介入。a)分别在x、y或z方向独立施加一固定幅值位g和持续时间为d的梯度，经过延迟时间τi后，发射射频，并采集梯度回波信号；
b)利用对所有信号相位做定量计算得到各个延迟时间点τi的涡流值大小ci，得到涡流随时间衰减的曲线；
c)利用levenberg-marquardt算法，对定量得到的涡流衰减曲线ci做多e指数衰减模型的非线性拟合，并且对所有e指数的幅度值做权重正则化控制，防止出现指数项过多的过拟合发生,自动剔除掉多余的指数项后，将所有e指数衰减常数固定kj；
d)对涡流e指数幅度做svd线性拟合，根据拟合定量得到的对应于各个时间常数kj的e指数涡流幅值aj，并更新谱仪的涡流预加重参数，重新采集涡流信号；
e)反复迭代进行以上d过程，当涡流值小于某一设定的阈值t之后，程序判定为最终迭代收敛，涡流校正完毕。
上述磁共振成像全自动预加重校正涡流的方法，其中，步骤b)中利用对所有信号相位做定量计算得到各个延迟时间点τi的涡流值大小ci的，表示为公式[1]:
ci＝mi/2πγte[1]
其中，mi表示第i次采集的傅里叶变换后回波信号的斜率(单位：弧度/毫米)，γ表示物理常量旋磁比(＝0.267513弧度/毫秒/微特斯拉)，te表示回波时间(单位：微特斯拉)；
其中，步骤c)中所用的多e指数衰减的函数表示为公式[2]:
其中，n表示e指数成分的个数，τi表示梯度延迟时间，c(τi)表示第i个梯度延迟点采集计算得到的涡流大小，aj和kj分别表示第j个e指数衰减成分的幅值和时间常数；
其中，对公式[2]进行拟合时所做的正则化控制可以表示为公式[3]:
其中min{}表示最小化，λ表示正则化控制项的常数。
以下实例分步介绍本发明方法的具体实现过程：本实施例所采用的全自动涡流调试技术实现于1.5t超导磁共振成像仪器。本实施例最终验证涡流校正效果的磁共振图像为快速自旋回波fse的水模图像。
a、选取一磁体方向做涡流预加重校正，本实例选取x方向，施加一固定幅值位g和持续时间为d的梯度，经过延迟时间τi后，发射射频，并采集回波时间为te的梯度回波信号。本实施例中，梯度幅值g＝10000.0微特斯拉/米，持续时间d＝20毫秒，te＝6毫秒，延迟时间τi设置成40个，分别是：1、2、3、4、5、6、8、10、15、20、25、30、40、50、60、70、80、90、100、110、120、130、140、150、160、170、180、190、200、250、300、400、500、600、700、800、1000、1200、1500和2000毫秒。
b、利用公式[1]对所有信号相位做定量计算得到各个延迟时间点τi的涡流值大小ci，得到涡流随时间衰减的曲线，如附图1所示。本实施例中，得到的40个延迟时间点的涡流值大小ci分别为：37.4620、35.5759、34.9731、34.7409、34.6033、34.7290、34.9395、34.9363、34.0809、33.7751、32.880632.0224、30.3341、28.5487、26.9940、25.4155、24.1821、22.789721.3279、19.7861、18.5141、17.5014、16.4106、15.2428、14.438613.5914、12.7890、11.8998、11.1530、8.1253、5.8067、3.3796、1.9679、1.2117、0.7434、0.4086、0.2080、0.1391、0.0813和0.0127微特斯拉。
c、利用公式[2]对涡流定量数值做levenberg-marquardt(列文伯格-马夸尔特法)多e指数衰减模型的非线性拟合，并利用公式[3]对所有e指数的幅度值做权重正则化控制，防止出现指数项过多的过拟合发生,自动剔除掉多余的指数项后，将所有e指数衰减常数固定。本实施例中初始e指数模型选用了4个指数模型，如果只用公式[2]做非线性拟合，则得到拟合结果为：a1＝-25.9644，a2＝24.0478，a3＝-5.68，a4＝71.0602(单位：微特斯拉)；k1＝12.3483，k2＝14.0353，k3＝57.8750，k4＝109.9339(单位：毫秒)。不做正则化拟合控制的话a值出现了极性交替，发生了过拟合。而采用公式[3]的正则化拟合则得到结果为：a1＝0.0281，a2＝0.0047，a3＝-0.0156，a4＝35.9238(单位：微特斯拉)；k1＝5.6944，k2＝18.9751，k3＝53.6533，k4＝176.3351(单位：毫秒)。可以发现a1，a2和a3的值均小于0.1，正则化在最小化过程中有效地防止了过拟合现象的发生，因此本实施例中涡流衰减模型确认为单指数模型，其它小幅值的e指数成分被剔除，系统时间常数被固定为k1＝176.3351毫秒。
d、对涡流e指数幅度做svd线性拟合，根据拟合定量得到的对应于各个时间常数kj的e指数涡流幅值aj，并更新谱仪的涡流预加重参数，重新采集涡流信号。本实施中磁共振系统的x方向涡流成分只有一个e指数成分k1＝176.3351毫秒，svd拟合后重新得到a1＝37.0259微特斯拉。
e、反复迭代进行以上d过程，当涡流值小于某一设定的阈值t之后，程序判定为最终迭代收敛，涡流校正完毕。本实施例中阈值设为2.0微特斯拉。收敛过程共进行了4次迭代。最终收敛时40个延迟时间点的涡流值大小分别为：0.2293、-0.6179、-0.9193、-0.7190、-0.5303、-0.3415、0.3137、0.5999、0.5898、1.0830、1.1886、1.2271、1.1131、0.8150、0.6467、0.5886、0.4130、0.3089、0.1161、-0.3292、-0.5494、-0.4621、-0.4378、-0.6785、-0.6720、-0.7268、-0.8124、-0.8818、-1.0654、-1.1561、-1.0061、-0.4635、-0.0459、0.1888、0.3444、0.3705、0.2728、0.1274、0.0856和0.1479微特斯拉。如附图2所示为预加重校正之前后系统x方向涡流曲线的效果对比图。
f、重复以上步骤(a)-(e)，用同样的方式对磁共振成像系统的y方向和z方向也进行涡流预加重校正。
g、做快速自旋回波fse序列的水模图像验证涡流校正效果。如附图3所示为涡流预加重校正之前后fse水模的效果对比图。
以上对本发明的具体实施例进行了详细描述，但本发明并不限制于以上描述的具体实施例，其只是作为范例。对于本领域技术人员而言，任何等同修改和替代也都在本发明的范畴之中。因此，在不脱离本发明的精神和范围下所作出的均等变换和修改，都应涵盖在本发明的范围内。