基于梯度電感約束和圖形處理單元加速算法的磁共振雙平面自屏蔽梯度線圈設計方法

黃清明,鄭 剛

(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并行運算加速算法，縮短了系數矩陣的計算時間。梯度線圈性能分析結果表明，采用最小化梯度電感能有效降低渦流影響，實現磁場自屏蔽，具有良好的梯度切換率，改善了梯度線圈整體性能。

圖3 梯度線圈參數化設計的GUI A.縱向梯度線圈設計的GUI； B.橫向梯度線圈設計的GUI

圖5 梯度線圈品質因素Q與線圈性能指標的關系 A.線圈品質因素Q與梯度非線性度的關系； B.線圈品質因素Q與梯度切換率的關系； C.線圈品質因素Q與梯度電感的關系； D.線圈品質因素Q與回路電阻的關系