基于梯度电感约束和图形处理单元加速算法的磁共振双平面自屏蔽梯度线圈设计方法

2020-04-02黄清明

中国医学影像技术 2020年3期

黄清明,郑刚

(1.上海健康医学院医学影像学院，上海 201318；2.上海理工大学光电信息与计算机工程学院，4.医学影像工程研究所，上海 200093；3.上海市分子影像重点实验室，上海 201318)

MR成像过程中，梯度线圈在灵敏区内产生期望的梯度磁场强度，但涡流使灵敏区实际梯度磁场为激发梯度磁场和涡流产生磁场的合成，导致脉冲上升和下降沿拖长，抑制梯度磁场快速变化，影响梯度切换率[1]；空间定位编码会使样品空间定位线性变差，导致图像畸变和伪影等失真现象，延长系统回波时间，减缓MR成像时间。为加快MR成像速度、获得高质量图像，必须减弱甚至消除涡流的影响[2-4]。

梯度线圈设计的实质是根据已知磁场分布计算梯度线圈的绕线模式，使其通电后产生期望磁场强度。1981年，BANGERT等[5]提出梯度线圈最简单模型和制作梯度线圈的基本方法，使MRI得以实现；1993年，TURNER等[6]总结梯度线圈设计的矩阵求逆法[7]、流函数法[8]和目标场点法[9]等，利用目标场点法预设圆柱表面磁场大小，通过计算机辅助软件计算圆柱形电流密度分布函数，经离散处理得到梯度线圈的绕线模型。之后，LIU等[10]提出双平面梯度线圈设计方法，计算出梯度电感能量分布和占电感总值的比例；GREEN等[11]提出单平面梯度线圈设计理论，以傅里叶变换电流密度函数得到梯度线圈功耗系数最小值。本研究采用目标场点法，增加梯度电感最小化和梯度磁场自屏蔽等约束条件建立数学模型，利用数值计算软件计算梯度线圈的理论参数性能和闭环绕线模式，以图形处理单元(graphic processing unit， GPU)加速算法有效解决了上述复杂数学模型中电流密度系数矩阵运算时间长的问题，并以电磁场仿真软件对设计方案进行正演验证，可在一定程度上减小梯度线圈电感值和线圈能耗，使梯度线圈具有良好梯度切换率[12]。

1 梯度线圈性能指标

MR通过梯度线圈产生的梯度磁场实现信号空间定位，主要性能指标包括梯度强度、梯度线性度、梯度切换率和梯度电感等[13]。

1.1 梯度线圈性能指标

1.1.1 梯度强度指梯度磁场强度大小，即每单位长度磁场强度的变化。快速成像和获取高分辨率MRI需要较高梯度强度[14]，梯度强度越高，选择脉冲序列越灵活。

(1)

式中γ为原子核的磁旋比，N为某个方向采集的信号数，BW为射频带宽，Gmax为最大梯度磁场强度，χ为分辨像素的最小限度。选定射频带宽BW后，梯度强度Gmax越大，则切面层厚越薄，图像空间分辨率越高，扫描视野越小。

1.1.2 梯度线性度指梯度强度变化率的大小，用梯度磁场强度变化的斜率表示。

(2)

式中γmax、γmin和γave分别为灵敏区内各点梯度强度变化斜率最大值、最小值和平均值。梯度线性度决定图像空间定位的准确性，λ值越小表示梯度磁场越精确，空间定位越准确图像质量越好，反之图像几何失真度越大。

1.1.3 梯度切换率指梯度强度从零增加到某一预定梯度强度值的速度。梯度磁场变化时，扫描序列处于等待状态，梯度切换率增加，则回波时间和重复时间减小，梯度变化速度越快，成像时间越短[15-16]。

1.1.4 梯度电感指梯度线圈电感量的大小，主要取决于线圈匝数、绕制方式和材料等。

(3)

式中μ0为真空磁导率，I为电流，J表示电流密度，v和v′表示电流密度分布的线圈空间。梯度线圈通电时，线圈导线周围产生电磁场，处于电磁场范围内的线圈导线发生“自感”和“互感”作用，在电导体中感应出涡流，抑制梯度磁场快速变化。

1.2 梯度线圈的品质因素梯度线圈品质因数Q与梯度线圈效率η、梯度电感L和梯度均匀性δ之间关系如下[17]：

(4)

2 梯度线圈的设计理论

梯度线圈设计方法有分离导线法和分布电流密度法。本研究采用分布电流密度法(即目标场点法)，加入梯度电感最小化和梯度磁场自屏蔽等约束条件进行优化。

2.1 梯度线圈目标场点法电流密度设双平面梯度线圈位于平面z=±a，梯度线圈半径ρ满足ρmin≤ρ≤ρmax，ρmin为线圈最小半径，ρmax为线圈最大半径。极坐标系上，通电梯度线圈的电流密度J(ρ,φ)分解成径向分量Jρ(ρ，φ)和切向分量Jφ(ρ，φ)，对电流密度J(ρ,φ)进行傅里叶变换[18]：

(5)

其中，c=π/(ρmax-ρmin)，Uq为电流密度系数，Q为展开的级数，q和k是整数。k=0时，为纵向梯度线圈z，产生纵向梯度磁场，电流密度表达式：

(6)

k=1时，为横向梯度线圈x或y，产生横向梯度磁场，电流密度表达式：

(7)

k≥2为对应阶次的匀场线圈。

2.2 自屏蔽约束双平面梯度线圈梯度线圈开合使电导体产生涡流，致梯度磁场不稳定，引起图像畸变[19]。为抑制涡流产生，在梯度线圈外增加一组屏蔽线圈，梯度线圈产生MR梯度磁场，与屏蔽线圈通反向电流产生的磁场叠加，抵消梯度系统涡流产生的磁场，这种结构称为自屏蔽线圈。设屏蔽线圈位置为z=±b(b>a)，屏蔽线圈电流密度函数表达式类似于梯度线圈的电流密度表达式，根据毕奥-萨伐尔定律可知：

(8)

Dq和Eq代表P(x,y,z)的功能函数，预设的目标场内的点Bz在目标场点的灵敏区和屏蔽区内求解(8)式，得Q和P2个未知变量。给出目标点和对应磁场强度值Bi(i=1， 2， …，Q+P)，将通过(8)式求解出的DiQ和EiQ写成矩阵形式：

(9)

求得梯度线圈和屏蔽线圈的电流密度函数，得到两组线圈绕线模式。

2.3 构建梯度电感最小约束函数设计梯度线圈时，既要找到最佳电流密度分布，满足灵敏区磁场强度要求，更应减轻梯度电感对涡流的影响[20]。梯度线圈为电感L和电阻R组成的电磁场回路，梯度电感L影响电流从零值上升到最大值63%的速度，即时间常数τ=L/R。τ越小越好，故应增大电阻R或减小电感L，而增大电阻R会导致线圈功率消耗增加，因此选择梯度电感L尽可能小。增加梯度电感L最小化约束条件的磁能表示为[21]：

(10)

(11)

其中，U为电流密度表达式的系数。用拉格朗日乘子法构造梯度电感最小约束函数：

(12)

对(12)式电流密度函数求偏导，推导出电流密度函数为：

U=W-1D(DTW-1D)-1B

(13)

根据电流密度散度为零的特点，双平面线圈的流函数I(ρ,φ)满足如下关系：

(14)

图1 电流密度系数矩阵GPU加速运算并行运算流程

对于横向梯度线圈y或x的流函数表示为：

I=-∑Usin[qc(ρ-ρmin)]cosφ

(15)

对于纵向梯度线圈z的流函数表示为：

(16)

电流密度的等量线图为：

I(ρ,φ)=Imin+(i+1/2)I0， (i=0， 1， 2， …，N-1)

(17)

式中I0=(Imax-Imin)/N，Imax为线圈平面内流函数最大值，Imin为线圈平面内流函数最小值，N为离散的线圈匝数，绘制(17)式等高线的曲线分布，保证每匝线圈中电流为I0，绕线位置与等高线重合，即得到满足场点要求的梯度线圈绕线[22]。

3 电流密度系数矩阵GPU加速运算

GPU具有高显存带宽、良好的浮点计算和并行运算能力。对双平面梯度线圈目标场点法加入梯度电感最小和自屏蔽线圈等约束条件，反演运算矩阵的求解属于大数据量并行处理计算过程，执行相同流程的并行化浮点运算。利用MATLAB平台在多核CPU基础上搭载多块GPU计算异构，即以分布式并行计算方式实现计算梯度线圈电流密度系数矩阵[23]。

电流密度系数矩阵GPU加速运算并行运算流程见图1。CPU按照预设梯度线圈参数选择目标场点，对GPU进行显存容量分配并复制原始数据，调用Matrix函数执行流函数和绕线模式运算处理，循环进行上述运算流程直至结束。将GPU运算所得结果数据复制到计算机内存并保存，释放显存中的数据[24]。

4 实现梯度线圈设计算法

基于梯度电感约束和GPU加速算法的梯度线圈的算法实现见图2。算法流程：①设定真空磁导率、梯度线圈平面的最小半径等参数值；②设定梯度线圈平面的最大直径、线圈两平面间距、灵敏区直径、目标场点个数、流函数展开级数、线圈匝数和预期梯度强度等参数；③选取目标场点；④给定极坐标下电流密度函数和磁场强度的函数表达式；⑤GPU调用矩阵函数加速运算后，将计算结果传回CPU，得到系数矩阵的求解结果；⑥在梯度电感最小化约束条件下计算电流密度函数的系数；⑦比较实际磁场强度与预期梯度强度，计算目标场点的梯度均匀性，若在允许范围5%以内跳转到流程⑧，否则跳转至流程③，重新选取目标场点；⑧离散化处理流函数，得到线圈绕线模型；⑨结束计算[25]。

由于梯度线圈的空间对称性，流程③选取目标场点时，应将目标场点设置在第1象限，且所取任意3个点不可在同一平面上，否则运算电流密度系数矩阵时将出现无解。

在梯度线圈目标场点法电流密度数学模型基础上，加入梯度电感最小化和自屏蔽线圈约束条件，采用GPU加速算法实现梯度线圈的流函数和计算绕线分布。利用数值计算工具MATLAB软件，实现编译算法梯度线圈参数化设计的图形用户界面(graphical user interface， GUI)如图3。

运用上述基于GUI的梯度线圈参数化设计软件得到梯度线圈绕线模型，通过电磁仿真软件ANSOFT正演计算梯度线圈各项性能指标，并对线圈边界进行修正，得到实际应用意义上的线圈绕线分布(图4)。

图2 基于梯度电感约束和GPU加速算法的梯度线圈算法实现流程图

5 实验与讨论

梯度电感是涡流产生的施予方。梯度线圈工作电流为脉冲电流，为使其上升前沿陡峭、加快MR成像速度，在满足其他性能指标的前提下，梯度电感值要尽可能小。根据梯度线圈的品质因素公式(4)，梯度线圈品质因数与梯度线圈效率、梯度电感和梯度均匀性之间存在一定函数关系，故以梯度线圈品质因素作为评价其性能的依据。理想梯度线圈要求梯度线圈效率尽可能高，梯度电感尽可能小，梯度磁场均匀性较好。事实上，梯度线圈绕线模式的关键在于选择梯度线圈品质因素Q值，Q值不同，得到的梯度绕线模式不同。以不同梯度线圈品质因素Q模拟仿真，对测得的梯度非线性度、梯度切换率、梯度电感和梯度电路电阻相关数据进行分析，发现梯度线圈品质因素Q与线圈性能指标的关系如图5。

随着梯度线圈品质因素Q值增大，梯度线性度指标改善，梯度切换率变小，但线圈能量损耗增大、系统效率降低，使梯度线圈整体性能变差。因此，适当选择Q值对梯度线圈性能指标而言非常重要。

6 结论

梯度电感引起的涡流导致梯度磁场性能指标变差，进而影响MR成像速度和图像质量。为缓解涡流对MRI系统的不利因素，本研究基于梯度线圈目标场点法，建立梯度电感最小和梯度磁场自屏蔽等约束条件下的数学模型，将电磁场计算的物理反问题转换为约束条件下数学模型的优化数值计算问题，采用GPU并行运算加速算法，缩短了系数矩阵的计算时间。梯度线圈性能分析结果表明，采用最小化梯度电感能有效降低涡流影响，实现磁场自屏蔽，具有良好的梯度切换率，改善了梯度线圈整体性能。