任自艷 謝德馨 李會香
(沈陽工業(yè)大學(xué)電氣工程學(xué)院 沈陽 110870)
醫(yī)用磁共振成像技術(shù)(MRI)自20世紀(jì)80年代發(fā)展以來,由于其獲得的圖像清晰、分辨率高,并具有對軟組織層次顯示好的優(yōu)點,因此成為現(xiàn)今臨床上不可缺少的五大醫(yī)療影像設(shè)備(超聲、放射學(xué)、MRI、X-CT和核醫(yī)學(xué))之一。臨床上獲得高質(zhì)量圖像的前提是,MRI系統(tǒng)的主磁體在成像區(qū)能夠產(chǎn)生足夠強度和高均勻度的靜磁場。然而,在磁體裝配完成之后,由于加工誤差以及溫度、環(huán)境等因素的影響,主磁場的均勻度往往與設(shè)計指標(biāo)相差較大,不能滿足成像的要求。這就需要憑借特定的工藝過程來校正成像區(qū)域的磁場,以提高主磁場的均勻度,這一過程稱為墊補,也叫勻場。
教育部高等學(xué)校博士學(xué)科點專項科研基金資助項目(20050142003)。
目前的勻場方法主要分為有源勻場和無源勻場兩種[1]。有源勻場主要采用通電線圈產(chǎn)生的場來補償成像區(qū)域磁場不均勻的諧波分量[2-3]。無源勻場則多采用在主磁體極面上粘貼用永磁材料(釤鈷等)制成的小磁片,利用小磁片陣列產(chǎn)生的磁場抵消成像區(qū)域的不均勻磁場分量[4-5]。本文研究的是開放式永磁型 MRI主磁體的勻場技術(shù),采用無源勻場方法。
將整數(shù)規(guī)劃方法引入MRI磁體的勻場,這一思想最初由B.Dorri等人[5]提出,它為后續(xù)研究奠定了基礎(chǔ)。文獻(xiàn)[6]建立了一個基于整數(shù)規(guī)劃的無源勻場數(shù)學(xué)模型,并將相應(yīng)計算軟件的求解結(jié)果以圖形的方式顯示出來,直接指導(dǎo)勻場過程。此方法所用墊片數(shù)目少,勻場時間短,減小了勻場對經(jīng)驗的依賴性。由于未考慮負(fù)性墊片(即以主磁體產(chǎn)生磁場的方向為參考,若小磁片充磁方向與主磁場的方向相反即為負(fù)片,充磁方向與主磁場的方向相同即為正片)的影響,致使主磁場的均勻度在提高了20%以后很難再進(jìn)一步改善。
本文結(jié)合文獻(xiàn)[6]方法的優(yōu)點,并考慮其不足,將整數(shù)規(guī)劃算法與實際人工勻場的模式相結(jié)合,以開放式永磁型MRI磁體系統(tǒng)的研究為基礎(chǔ),提出一種基于整數(shù)規(guī)劃思想的中心勻場方法。在這一方法的優(yōu)化模型中,以成像區(qū)域縱向中軸線的磁感應(yīng)強度值為基礎(chǔ),優(yōu)化勻場所用小磁片放置的位置和數(shù)目,使得小磁片在空間產(chǎn)生的磁場能抵消成像區(qū)域不均勻的分量,縮小成像區(qū)的磁場值與中心值的差距。本文詳細(xì)介紹了這一方法的數(shù)學(xué)模型與實施步驟,最后通過勻場實驗驗證所提出的勻場方法的有效性。
本文研究的成像區(qū)域為主磁體上、下極面之間直徑為 40cm的球域,在該球域內(nèi),應(yīng)能產(chǎn)生豎直方向的均勻磁場。主磁體產(chǎn)生磁場的均勻度可用成像空間內(nèi)磁感應(yīng)強度的正反向峰峰值偏差與磁感應(yīng)強度的平均值之比[7-8]來衡量,即
式中 η ——磁場的不均勻度,通常用10-6表示;
Bmax——成像區(qū)域內(nèi)磁感應(yīng)強度的最大值;
Bmin——成像區(qū)域內(nèi)磁感應(yīng)強度最小值。
目前,勻場工程師憑借以往的經(jīng)驗試探性地確定用于補償不均勻場的小磁片的片數(shù)及位置,一般需要兩個工程師調(diào)試3~4天完成,最后所得墊片的數(shù)目較多,其位置也沒有規(guī)律(實際墊片分布如圖1所示,圖中小圓圈代表小磁片)。這種人工勻場工作量大、效率低,不能擺脫經(jīng)驗的束縛,因此,改進(jìn)勻場技術(shù),提高勻場效率是本文研究的重點。
圖1 無源勻場示意圖Fig.1 Sketch map of passive shimming
成像區(qū)球域是勻場的目標(biāo)區(qū)域,在這一區(qū)域中,由于不存在電流源,因而標(biāo)量磁位φ 滿足拉普拉斯方程
亦即φ 為調(diào)和函數(shù)。由調(diào)和函數(shù)的性質(zhì)可知,其最大值和最小值均位于所研究區(qū)域的邊界上,因此實際的勻場目標(biāo)區(qū)域可縮減到成像區(qū)域的邊界球面。
磁體裝配后,首先需要測定成像區(qū)球面上的磁場,這些磁場的采樣值稱為基礎(chǔ)場?;A(chǔ)場值的采集可利用核磁共振磁場強度檢測儀PT2025來完成?;A(chǔ)場采集示意圖如圖2所示,成像區(qū)球面自上而下分為13層,每層球面上平均設(shè)定12個采樣點,共得到156個數(shù)據(jù),以此為基礎(chǔ)對無源勻場工作進(jìn)行規(guī)劃。在圖2中,Bm為成像球域縱向中軸線上的磁感應(yīng)強度,在球域內(nèi)測定,下文稱為中心場值;Bi為成像球面上i點的基礎(chǔ)場值。
圖2 基礎(chǔ)場樣本采集示意圖Fig.2 Sample selection picture of basic magnetic field
在實際的人工勻場過程中,首先將成像區(qū)域的中心場值 Bm與成像區(qū)球面上的磁感應(yīng)強度 Bi進(jìn)行比較。在同一層面上,如果 Bi低于 Bm,則在與極面相應(yīng)的位置上貼正性磁片,從而提高周圍場強;反之,在與極面相應(yīng)的位置貼上負(fù)性磁片,從而降低周圍場強。以球的赤道面為界,上半球的勻場在上極面貼片,下半球的勻場在下極面貼片,球赤道面的勻場則在上下極面共同完成。
對應(yīng)于圖2所示的樣本采集圖,在實際成像空間內(nèi)自上而下測量13中心場值,通常這些測量值的均勻度較高,例如由式(1)算得其不均勻度可以僅為7.699×10-6。如果球面上各點的磁場值通過勻場均能逼近中心場值,則成像區(qū)的均勻度將符合MRI的成像要求。中心場值的均勻性為模型的建立提供了事實依據(jù)。本文將人工以中心場值為準(zhǔn)的勻場過程稱為中心勻場策略。
本文將上述人工勻場思路與優(yōu)化理論相結(jié)合,提出了基于整數(shù)規(guī)劃的中心勻場方法,實現(xiàn)了理論研究與實際勻場過程的結(jié)合。建模時選擇最小化磁場偏差的峰峰值為優(yōu)化目標(biāo),在一系列不等式約束下實現(xiàn)小磁片放置的位置和數(shù)目的優(yōu)化。小磁片可能放置的位置作預(yù)先設(shè)定,其最終放置的位置由該位置上的磁片數(shù)決定,當(dāng)磁片數(shù)為零時,該位置即不放置磁片。這一方法的數(shù)學(xué)模型為
式中 T——磁場偏差的峰峰值;
λ ——與主磁體產(chǎn)生的磁場大小相關(guān)的系數(shù),為常值;
Xi——位置i處墊補片的數(shù)目,取整數(shù);
N——墊補位置的總數(shù);
M——磁場測量點的總數(shù);
Bj——當(dāng)前勻場周期開始之前空間點j磁場感應(yīng)強度的測量值;
?Bij——在位置 i上放置單個墊片從而導(dǎo)致點 j磁感應(yīng)強度幅值的變化;
Bm——球面上每個測量點對應(yīng)的中心場值。
表達(dá)式-ni≤Xi≤ni確保每一個墊補位置放置小磁片的數(shù)目在合適的范圍內(nèi),ni表示在位置 i上最多可以放置墊片的數(shù)目。式(4)中的前兩個約束條件表示:在空間點 j處的磁感應(yīng)強度與極面上所有墊片對點j的貢獻(xiàn)之和,減去對應(yīng)的中心場值Bm,得到點j處磁場偏差值的大小應(yīng)在-T/2~+T/2的范圍內(nèi)。
置于極面上確定位置的一定大小的小磁片對球面上采樣點處磁場的貢獻(xiàn)稱為“影響值”,即將放置單個小磁片后的球面上的磁場值與放置前的磁場值相減得到的差值。為了避免解析求解的困難,所有的影響值均在恒定條件(溫度、環(huán)境等)下通過實測得到的。
第3節(jié)所述數(shù)學(xué)模型的求解借助Lingo 8.0軟件提供的優(yōu)化函數(shù)采用分支定界法[9-10]來實現(xiàn)。獲得勻場方案后,本文采用逐步軟件仿真驗證的方法,在上述數(shù)學(xué)模型的軟件實現(xiàn)中,分若干步驟進(jìn)行,最后獲得校正后的勻場方案。理論研究過程中,采用磁場的數(shù)值分析方法對小磁片墊補方案的總體影響進(jìn)行校驗。
(1)極面上等分12份間隔30°的勻場模型
選擇首次采集到的156個數(shù)據(jù)為基礎(chǔ)場;按照式(3)和式(4)建立整數(shù)規(guī)劃數(shù)學(xué)模型;利用Lingo軟件對其求解。為了便于建立極面上小磁片位置分布的模型,用 Matlab軟件編制程序,實現(xiàn)數(shù)字信息的圖形化。下極面上間隔 30°徑向射線處的小磁片分布如圖3所示,圖上的小圓圈代表勻場用的小磁片,旁邊的數(shù)字代表該位置放置磁片的數(shù)目。
圖3 極面墊片分布示意圖Fig.3 Distribution of shims on the polar
(2)極面上從15°開始等分12份的勻場模型
如圖4所示,實線處表示步驟1勻場校正后的磁場值,虛線為步驟(2)勻場采用的值。以步驟(1)的勻場方案校正后的156個值為基礎(chǔ),通過相鄰實線刻度的平均插值,獲得成像球面上15°處的磁感應(yīng)強度值。然后求解整數(shù)規(guī)劃模型,獲得15°徑向射線處的墊片分布。
(3)極面上從7.5°開始等分24份的勻場模型與步驟(2)相似,以步驟(2)勻場校正后的結(jié)果為基礎(chǔ),進(jìn)行插值、模型求解、勻場,最終獲得墊補信息。
圖4 插值獲得15°處的磁場強度示意圖Fig.4 Interpolation chart of obtaining magnetic intensity on the 15-degree positions
實際勻場按照最后獲得的墊片陣列圖在上下極面貼片,然后用磁場強度檢測儀測量成像球域面上的磁感應(yīng)強度值就可得到勻場后的結(jié)果。在本文中,研究的是勻場軟件算法,為了節(jié)省時間,驗證上述勻場數(shù)學(xué)模型的有效性,故采用數(shù)值分析方法計算小磁片陣列的影響值。數(shù)值分析方法就是基于數(shù)值仿真技術(shù),通過對實際磁體工作環(huán)境的模擬,將不同位置需要的小磁片放置在磁體極面上,采用ANSOFT靜磁場求解軟件包得到一次勻場方案中的磁片陣列對空間成像區(qū)域的影響值。
一個小磁片在空間產(chǎn)生的磁場值約為 10-7T左右,相對于MRI主磁體0.3T的磁場來說,兩者數(shù)量級相差很大。再考慮到有限元分析方法的舍入誤差和截斷誤差,若選擇包括永磁主磁體和鐵軛在內(nèi)的系統(tǒng)進(jìn)行磁場分析,顯然不能保證對磁場影響值計算的精度。綜合以上因素,并考慮到永磁主磁體產(chǎn)生的磁場在 0.3~0.5T的低場因素,本文只對小磁片陣列進(jìn)行建模,并將磁體系統(tǒng)的成像球域設(shè)為求解域。仿真后采集求解域的磁場值,進(jìn)行勻場,獲得每一步驟的勻場結(jié)果。
表1和表2分別列出了兩臺MRI主磁體勻場過程的計算結(jié)果。表中的不均勻度采用式(1)計算獲得,提高幅度是勻場前后不均勻度的差值與勻場前不均勻度的百分比值。
從表1和表2可以看出,兩臺MRI主磁體勻場后的總體均勻度分別提高了52.56%和41.08%。
表1 勻場前后不均勻度的比較Tab.1 Inhomogeneity comparison of fore and after shimming
表2 勻場前后不均勻度的比較Tab.2 Inhomogeneity comparison of fore and after shimming
本文研究的勻場方法最后在 MRI磁體生產(chǎn)企業(yè)進(jìn)行了實踐。選用的實驗磁體基礎(chǔ)場不均勻度為2002.75×10-6,勻場后磁體成像空間不均勻度降低到了902.57×10-6,均勻度提高54.93%。圖5和圖6分別為勻場前后成像樣品區(qū)采樣點的磁通密度分布情況。
圖5 基礎(chǔ)場的磁通密度分布Fig.5 Distribution of basic magnetic flux intensity
圖6 勻場后磁通密度分布Fig.6 Distribution of magnetic flux intensity after shimming
以上理論分析與勻場實驗均表明了以中心場值為基礎(chǔ)的勻場思路是正確的,它可以方便快捷地獲得勻場墊片的分布信息,并大幅度提高成像區(qū)域的均勻度。
本文針對開放式永磁型 MRI主磁體的勻場問題,提出了一種以成像區(qū)域內(nèi)中軸線上磁感應(yīng)強度值為基礎(chǔ)的無源勻場方法,即中心勻場方法。所獲得的勻場方案經(jīng)實驗驗證可將成像區(qū)的均勻度提高50%以上,更符合實際情況,可以更直接地指導(dǎo)人工勻場。
[1]謝德馨, 王恩龍, 張艷麗.基于動態(tài)規(guī)劃的永磁MRI磁體的無源勻場[J]. 沈陽工業(yè)大學(xué)學(xué)報, 2007,29(4): 396-399.Xie Dexin, Wang Enlong, Zhang Yanli. Passive shimming of permanent magnet for MRI based on dynamic programming[J]. Journal of Shenyang University of Technology, 2007, 29(4): 396-399.
[2]吳海澄, 周荷琴, 劉正敏. 有源勻場電流優(yōu)化算法研究[J]. 北京生物醫(yī)學(xué)工程, 2007, 26(5): 466-469.Wu Haicheng, Zhou Heqin, Liu Zhengmin. Research on the optimization algorithm of active shimming current[J]. Beijing Biomedical Engineering, 2007,26(5): 466-469.
[3]Yamamoto S, Yamada T, Morita M, et al. Field correction of a high-homogeneous field superconducting magnet using a least squares method[J]. IEEE Transactions on Magnetics, 1985, 21(2): 689-701.
[4]Hector Sanchez Lopez, Liu Feng, Ewald Weber, et al.Passive shim design and a shimming approach for biplanar permanent magnetic resonance imaging magnets[J]. IEEE Transactions on Applied Superconductivity, 2008, 44(3): 394-402.
[5]Dorri B, Vermilyea M E. Passive shimming of MR magnets: algorithm, hardware, and results[J]. IEEE Transactions on Applied Superconductivity, 1993,3(1): 254-257.
[6]王恩龍. 一種特殊永磁型磁體的勻場技術(shù)研究[D].沈陽: 沈陽工業(yè)大學(xué), 2007.
[7]張艷麗, 謝德馨, 王恩龍, 等. 開放式永磁機構(gòu)無源勻場方法[J]. 電工技術(shù)學(xué)報, 2007, 22(6): 1-5.Zhang Yanli, Xie Dexin, Wang Enlong, et al. Passive shimming method for open permanent magnet structure[J]. Transactions of China Electrotechnical Society, 2007, 22(6): 1-5.
[8]王恩龍, 謝德馨, 張艷麗, 等. 基于整數(shù)規(guī)劃的永磁型磁共振成像裝置主磁體的無源勻場技術(shù)研究[J].中國生物醫(yī)學(xué)工程學(xué)報, 2008, 27(4): 533-536.Wang Enlong, Xie Dexin, Zhang Yanli, et al. A passive shimming method based on integer programming for permanent magnet of mri devices[J].Chinese Journal of Biomedical Engineering, 2008,27(4): 533-536.
[9]袁新生, 邵大宏, 郁時煉. LINGO和Excel在數(shù)學(xué)建模中的應(yīng)用[M]. 北京: 科學(xué)出版社, 2007.
[10]甘應(yīng)愛, 田豐, 李維錚, 等. 運籌學(xué)[M]. 北京: 清華大學(xué)出版社, 1995.