一種磁感應磁聲內窺圖像的建模與仿真方法
【專利摘要】一種磁感應磁聲內窺圖像的建模與仿真方法，所述方法首先建立生物腔體組織的橫截面模型和生物腔體組織橫截面的電磁特性參數(shù)模型；然后應用有限元分析的方法，對腔體組織在磁激勵場中產生感應電流的過程進行仿真，再由感應電流仿真腔體組織在靜磁場中產生聲源的過程；之后根據(jù)聲源的分布仿真出腔體組織產生的磁聲信號；最后根據(jù)磁聲信號重建出生物腔體的橫截面圖像。本發(fā)明可靈活地調整生物腔體組織橫截面模型，改變模型中腔體組織的變異類型及程度，設置多層腔體壁組織的電磁特性參數(shù)和聲學參數(shù)，精確地仿真聲速不均勻的多層生物腔體組織的聲場，因此可以為成像算法和圖像后處理算法的研究和性能測試等提供足夠的樣本圖像數(shù)據(jù)。
【專利說明】
-種磁感應磁聲內窺圖像的建模與仿真方法
技術領域
[0001] 本發(fā)明設及一種醫(yī)學成像技術，特別是對生物腔體的磁感應磁聲內窺橫截面圖像 進行建模與仿真。
【背景技術】
[0002] 生物組織磁感應磁聲內窺化ndoscopic Magnetoacoustic Tomography with Magnetic Induct ion, EMAT-MI)成像屬于功能成像，其成像的物理基礎是生物組織的磁聲 效應，即生物組織變異部位與非變異部位的電導率不同，導致產生的感應滿流不同，從而在 靜磁場中產生不同強度的超聲信號（即磁聲信號KEMAT-MI結合了生物電阻抗成像的高對 比度和超聲掃描成像的高空間分辨率，可得到生物腔體組織（含變異組織）內部的形態(tài)結 構、組織類型和成分等信息。然而EMAT-MI成像技術尚處于實驗室研究階段，未大規(guī)模應用 于實際，而改善成像儀器的結構和校準參數(shù)、優(yōu)化圖像分析和處理算法等都是建立在分析 大量樣本圖像數(shù)據(jù)的基礎之上的。因此，為了獲取足夠多的樣本圖像數(shù)據(jù)，有必要尋找一種 對生物腔體組織磁感應磁聲內窺圖像進行建模與仿真的方法。

【發(fā)明內容】

[0003] 本發(fā)明的目的在于針對現(xiàn)有技術之弊端，提供一種磁感應磁聲內窺圖像的建模與 仿真方法，為成像和圖像后處理的研究和性能測試等提供數(shù)據(jù)源。
[0004] 本發(fā)明所述問題是W下述技術方案解決的：
[0005] -種磁感應磁聲內窺圖像的建模與仿真方法，所述方法首先建立生物腔體組織的 橫截面模型和生物腔體組織橫截面的電磁特性參數(shù)模型;然后應用有限元分析的方法，對 腔體組織在磁激勵場中產生感應電流的過程進行仿真，再由感應電流仿真腔體組織在靜磁 場中產生聲源的過程;之后根據(jù)聲源的分布仿真出腔體組織產生的磁聲信號；最后根據(jù)磁 聲信號重建出生物腔體的橫截面圖像。
[0006] 上述磁感應磁聲內窺圖像的建模與仿真方法，具體處理按W下步驟進行：
[0007] a.建立生物腔體組織的橫截面模型：
[000引成像導管位于生物腔體組織橫截面模型的中屯、，接收磁聲信號的超聲探測器位于 成像導管頂端，將超聲探測器看作理想的點換能器，其掃描軌跡為平行于成像平面、半徑趨 近于0的圓形軌跡；
[0009] b.建立生物腔體組織橫截面的電磁特性參數(shù)模型：
[0010] W生物腔體組織橫截面模型的中屯、為起始點，將模型等角度劃分為m份，每一份近 似為多層組織，對模型施加脈沖磁激勵，并對每一份組織接收到的磁聲信號進行仿真，成像 導管上的超聲探測器所處的角度為
[0011] 目 i = 360(i-l)/m
[001^ 其中，i = l，2, . . .，m，對應的成像區(qū)域的角度范圍為[0ia，0ib]，其中0ia=0廣180/ m，目化=目 i+180/m。
[0013] 確定每個成像區(qū)域的每層組織的電導率和厚度參數(shù)，形成生物腔體組織橫截面的 電磁特性參數(shù)模型。W血管的橫截面模型為例，多層血管壁組織的參數(shù)設置示例如表1所 列。
[0014] C.仿真EMAT-MI成像中聲源產生的過程：
[0015] ①應用仿真軟件構造亥姆霍茲化elmholtz)線圈，并對亥姆霍茲線圈施加高斯脈 沖電流作為激勵源，產生交變磁場Bi,沿亥姆霍茲線圈軸向施加一個穩(wěn)恒勻強磁場Bo,多層 生物腔體組織的橫截面模型同軸放置于亥姆霍茲線圈中間;穩(wěn)恒磁場和亥姆霍茲線圈的參 數(shù)設置示例如表2所列。
[0016] ②W空氣為背景域，W自由四面體為單元對亥姆霍茲線圈W及生物腔體組織的橫 截面模型進行劃分，進而采用有限元分析的方法仿真得到生物腔體組織中感應滿流J的分 布；
[0017] ③根據(jù)感應滿流J仿真出生物腔體組織的聲源▽?(/><巧,)的分布，其中，"X"表示 向量積；
[0018] d.根據(jù)聲源的分布仿真生物腔體組織產生的磁聲信號，得到超聲探測器在時刻t、 角度9i、位置r處接收到的多層生物腔體組織產生的磁聲信號的聲壓pi(r，t)a = l，2，...， m)；
[0019] e.重建極坐標系下的EMAT-MI圖像：
[0020] ①根據(jù)步驟d得到的m個磁聲信號的聲壓pi(r，t)(i = l，2,. . .，m)計算出
[0021]
，
[0022] 式中，Cs是超聲波的波速，ri是0-1平面中超聲探測器所在位置(與成像導管在X-Y 平面中的成像角度01相對應）；
[0023] ②由下式重建出在角度01處的聲源分布：
[0024]
;
[0025] 式中，q是ri處的單位矢量；
[00%]③將上式代入下式，得到腔體組織在角度01處的電導率分布：
[0027]
[0028] 將〇i(r)作為角度01、位置r處的腔體組織橫截面極坐標視圖的灰度值，從而得到極 坐標下的EMAT-MI圖像；
[00巧]f.圖像的坐標轉換：
[0030] 將極坐標下的EMAT-MI圖像轉換為X-Y直角坐標系下的橫截面視圖，具體方法如 下：
[0031] 建立X-Y平面直角坐標系，坐標原點是成像導管中屯、，水平向右的方向為X軸正方 向，垂直于X軸向上的方向為Y軸正方向，設極坐標系中的一點（j，k)的灰度值為f (j，k)，該 點在X-Y坐標系中的對應點的坐標為（j/，k/ )，灰度值為g(j/，k/ )，其中j G [O，231 ]，k G [O， d]，j/ G[-d，d]，k/ ￡[-(1，(1]，(1為極坐標視圖中極徑的最大值，則有：
[0032] g(j'，k')=f(j，k)，
[0033] 其中
[0034]
[0035] 本發(fā)明可靈活地調整生物腔體組織橫截面模型，改變模型中腔體組織的變異類型 及程度，設置多層腔體壁組織的電磁特性參數(shù)和聲學參數(shù)，精確地仿真聲速不均勻的多層 生物腔體組織的聲場，因此可W為成像算法和圖像后處理算法的研究和性能測試等提供足 夠的樣本圖像數(shù)據(jù)。
【附圖說明】
[0036] 圖1是含有脂質斑塊的血管橫截面模型示例；
[0037] 圖2是EMAT-MI成像導管上的超聲探測器在角度01處接收磁聲信號的示意圖；
[0038] 圖3是將圖1中的血管橫截面模型等角度劃分后，將其中的一份近似為多層組織的 示意圖；
[0039] 圖4是0-1坐標系中的多層生物腔體組織的示意圖；
[0040] 圖5是亥姆霍茲化eImholtz)線圈激勵示意圖，圖中R是線圈半徑。
[0041 ]表1是多層血管壁組織參數(shù)設置舉例。
[0042] 表2是穩(wěn)恒磁場和亥姆霍茲線圈的參數(shù)設置示例。
[0043] 文中各符號為:X、Y、腔體橫截面模型所在的X-Y平面直角坐標系的橫軸和縱軸，其 中，成像導管的中屯、為坐標原點，水平向右的方向為X軸正方向，垂直于X軸向上的方向為Y 軸正方向;m、橫截面模型被等角度分割的總份數(shù);01、成像導管的第i個成像角度，其中，i = 1，2,...，m;目13、目lb、成像導管在角度目1處進行成像時對應的角度范圍的上、下限；目、1、目-1坐 標系的橫軸和縱軸;Bi、激勵磁場的磁通密度;Bo、穩(wěn)恒磁場的磁通密度;J、生物腔體組織中 的感應滿流；X、向量積;V、哈密頓算子;j、k、0-1坐標系沿0方向和1方向的單位向量;t、時 間；r、0-l坐標系中的一點;pi(r，t)、超聲探測器在時刻t角度01位置r處接收到的多層生物 腔體組織產生的磁聲信號的聲壓，其中i = l，2, . . .，m;Ji(r)、生物腔體組織在角度01位置r 處產生的感應滿流，其中i = l，2,. . .，m;B〇i(r)、生物組織在角度01位置r處的激勵磁場和穩(wěn) 恒磁場的疊加磁通密度，其中i = 1，2，. . .，m; n (t)、隨時間變化的激勵電流信號；Cs、超聲波 的波速;P、生物組織的密度;vii(r)和vei(r)分別是角度01方向位置r處的質點在1方向和0方 向的振動速度，其中i = l，2,. . .，m;(j，k)、位于0-1平面內腔體組織上一點的坐標；A 0、A 1、坐標軸0和坐標軸1方向的的離散空間間距；A t、離散時間間距;n、離散時間點;pn(j，k)、 位置(j，k)的質點在時刻n產生的磁聲信號的聲壓;v；；0'，A)、<(/，&)、時刻n位置(j，k)的質 點在e方向和1方向的振動速度;P(j，k)、生物組織在點（j，k)處的密度;cs(j，k)、磁聲信號 在位置（j，k)的傳播速度;ri、0-l平面中超聲探測器所在位置(與成像導管在X-Y平面中的 成像角度相對應）；q、超聲探測器所在位置ri處的單位矢量;〇i(r)、生物腔體組織在角度 目i位置r處的電導率，其中i = I，2，. . .，m; (j/，k/ )、0-1坐標系中的一點（j，k)在X-Y直角坐 標系中對應點的坐標;；^〇，1〇、邑(^，4/)、點〇，1〇和點〇/，4/)的灰度值;(1、極坐標視圖極 徑的最大值。
【具體實施方式】
[0044] 下面結合附圖對本發(fā)明作進一步詳述。
[0045] 本發(fā)明方法的步驟包括：
[0046] (1)建立多層生物腔體組織橫截面的形態(tài)模型：
[0047] W血管橫截面形態(tài)模型為例，如附圖1所示，模型包括成像導管(接收磁聲信號的 超聲探測器位于成像導管頂端）、管腔、血管壁內膜/中膜、外膜和斑塊(脂質斑塊、纖維斑 塊、巧化斑塊或混合斑塊)五部分。其中成像導管位于模型中屯、，周圍由內向外依次為管腔、 斑塊、血管壁內膜/中膜和外膜?？筛鶕?jù)組織變異的類型（如巧化、纖維化或者脂質斑塊)和 大小，W及血管內腔、血管壁內膜/中膜、外膜的厚度建立不同的血管橫截面形態(tài)模型。本發(fā) 明方法忽略超聲探測器的孔徑效應，將其看作理想的點換能器，其掃描軌跡為平行于成像 平面、半徑趨近于0的圓形軌跡。模型所在的坐標系為X-Y平面直角坐標系，其中坐標原點是 成像導管中屯、，水平向右的方向為X軸正方向，垂直于X軸向上的方向為Y軸正方向。
[004引（2)建立多層生物腔體組織橫截面的電磁特性參數(shù)模型：
[0049] 首先，如附圖2所示，W生物腔體組織橫截面形態(tài)模型的中屯、為起始點，將模型等 角度劃分為m份，對模型施加脈沖磁激勵，并對每一份組織接收到的磁聲信號進行仿真。成 像導管上的超聲探測器所處的成像角度為
[0050] 目 i = 360(i-l)/m (1)
[0051] 其中，i = l，2，. . .，m。對應的成像區(qū)域的角度范圍為[0ia，0ib]，其中0ia=0廣180/ m，目ib =目i+180/m。例如，當m=360時，將附圖1中含有脂質斑塊的血管橫截面形態(tài)模型Wx-Y 坐標系的原點為中屯、等分為360份，超聲探測器對每一份組織采集磁聲信號的過程中所處 的角度為目 1 = 0°，目2 = 1°，...，目i = 360(i-l)/360,...，目360 = 3 5 9°。
[0052] 然后，如附圖4所示，在0-1坐標系中，每一份多層組織的表面平行于0軸且垂直于1 軸。其中9軸正方向為水平向右的方向，1軸正方向為垂直于9軸向上的方向，表示多層組織 的厚度。確定每個成像角度對應的組織中每層的電導率和厚度參數(shù)，形成多層生物腔體組 織的電磁特性參數(shù)模型。W血管的橫截面模型為例，多層血管壁組織的參數(shù)設置示例如表1 所列。不同成像角度對應的多層生物腔體組織僅厚度參數(shù)不同，圖4中的生物腔體組織分為 a、b、c、d、e五層，若形態(tài)模型中不包含斑塊，則c、d層的厚度參數(shù)為0。
[0化3] (3)仿真EMAT-MI成像中聲源產生的過程：
[0054] 本發(fā)明采用有限元分析的方法進行聲源的仿真，具體步驟如下：
[0055] 首先，應用COMSOL軟件構造亥姆霍茲化elmholtz)線圈，如附圖5所示，并施加高斯 脈沖電流作為激勵源產生交變磁場Bi,沿線圈軸向施加一個穩(wěn)恒勻強磁場Bo,多層生物腔體 組織模型置于線圈中間，其軸向與亥姆霍茲線圈的軸向同向；然后，W空氣為背景域，W自 由四面體為單元對亥姆霍茲線圈W及生物腔體組織模型進行劃分，進而仿真得到生物腔體 組織中感應滿流J的分布;最后，根據(jù)感應滿流J仿真出生物腔體組織的聲源▽?(/>< A)的分 布，其中，"X"表示向量積。穩(wěn)恒磁場和亥姆霍茲線圈的參數(shù)設置示例如表2所列。
[0化6] 表1 [0化7]
[0化引
[0化9]
[0060] (4)仿真生物腔體組織產生的磁聲信號：
[0061] 生物組織產生的磁聲信號的本質是超聲波，描述磁聲信號在聲學均勻介質中傳播 的物理模型為：
[006^
P)
[0063] 其中，i = l，2，...，m; V為哈密頓算子；"?"為向量的點積;j和k分別為0-1坐標系 的0方向和1方向的單位向量;t為時間；r為0-1坐標系中的一點，如附圖4所示；Cs為超聲波 的波速;Pi(r，t)為超聲探測器在時刻t角度01位置r處接收到的多層生物腔體組織產生的磁 聲信號的聲壓;Ji(r)為生物組織在角度01位置r處產生的感應滿流;Boi(r)為生物組織在角 度0i位置r處穩(wěn)恒磁場的磁通密度。
[0064] 生物組織的聲學特性與流體相似，在理想流體中，聲波的傳播滿足牛頓第二定律、 質量守恒定律和物態(tài)方程。設隨時間變化的激勵電流信號為n(t)，則由聲波的連續(xù)性方程、 運動方程和物態(tài)方程(張海欄.理論聲學.北京:高等教育出版社.2007.)可將式(2)改寫成：
[0065]
(3)
[OOW 其中，i = l ,2, . . . ,m;pi(;r,t)是聲波的聲壓;P是生物組織的密度;Vii(T)和V0i(;r) 分別是角度0i方向位置r處的質點在1方向和0方向的振動速度。
[0067]義用時域有限差分法（Yee K S.Numerical solution of initial boundary value problems involving Maxwell's equations in isotropic media.IEEE Trans .Antennas Propag，1966，14(3): 302-307 ?)和Ca:rtesian網格劃分法（Sakaguchi T， Hirano T,Watanabe Y，et al. Inner head acoustic field for bone-conducted sound calculated by finite-difference time-domain method. Japanese journal of applied physics，2002，41 (5S): 3604-3608.)對方程(3)中的各個物理量離散化得到：
[0068]
(4)
[0069] 式中，（j，k)表示位于0-1平面內腔體組織上一點的坐標；A 0和A 1分別表示坐標 軸0和坐標軸1方向的離散空間間距；A t表示離散的時間間距;n表示離散時間點；pn( j，k) 是位置（j，k)的質點在時刻n產生的磁聲信號的聲壓；v；：aA)和畔化均是時刻n位置（j，k) 的質點在e方向和1方向的振動速度;P(j，k)是生物組織在點（j，k)處的密度；cs(j，k)是磁 聲信號在位置(j，k)的傳播速度。
[0070] 根據(jù)式(4)即可仿真得到超聲探測器在時刻t、角度01、位置r處接收到的多層生物 腔體組織產生的磁聲信號的聲壓91^，*)(1 = 1，2，...，111)。
[OOW (5)重建極坐標下的EMAT-MI圖像：
[0072] 圖像重建的實質是由超聲探測器接收到的磁聲信號計算出生物組織的電導率分 布。具體方法如下：
[0073] 首先，根據(jù)步驟(4)得到的m個磁聲信號pi(r，t)(i = l，2, . . .，m)計算出
[0074]
(5)
[0075] 然后，由下式重建出在角度01處的聲源分布：
[0076]
(6)
[0077] 式中，i = l，2, ...，m;ri是0-1平面中超聲探測器所在位置(與成像導管在X-Y平面 中的成像角度01相對應）；q是ri處的單位矢量。
[0078] 最后，將式(6)代入下式重建出腔體組織在角度01處的電導率分布：
[0079]
(7)
[0080] 其中，i = l，2，. . .，m。將〇i(r)作為角度01、位置r處的腔體組織橫截面極坐標視圖 的灰度值。
[0081 ] (6)圖像的坐標轉換：
[0082] 將步驟(5)中獲得的0-1坐標系中的極坐標視圖轉換為X-Y直角坐標系下的橫截面 視圖。具體方法如下：
[0083] 設0-1坐標系中的一點（j，k)的灰度值為f( j，k)，該點在X-Y坐標系中的對應點的 坐標為0'，k')，灰度值為g(j'，k')，其中jG[0,23T]，kG[0，d]，j'G[-d，d]，k'G[-d，d]，d 為極坐標視圖縱坐標的最大值。那么有：
[0084] g(j',k')=f(j,k) (8)
[0085] 其中
[0086] 餅
【主權項】
1. 一種磁感應磁聲內窺圖像的建模與仿真方法，其特征是，所述方法首先建立生物腔 體組織的橫截面模型和生物腔體組織橫截面的電磁特性參數(shù)模型;然后應用有限元分析的 方法，對腔體組織在磁激勵場中產生感應電流的過程進行仿真，再由感應電流仿真腔體組 織在靜磁場中產生聲源的過程;之后根據(jù)聲源的分布仿真出腔體組織產生的磁聲信號；最 后根據(jù)磁聲信號重建出生物腔體的橫截面圖像。2. 根據(jù)權利要求1所述的磁感應磁聲內窺圖像的建模與仿真方法，其特征是，所述方法 包括以下步驟： a. 建立生物腔體組織的橫截面模型： 成像導管位于生物腔體組織的橫截面模型的中心，接收磁聲信號的超聲探測器位于成 像導管頂端，將超聲探測器看作理想的點換能器，其掃描軌跡為平行于成像平面、半徑趨近 于O的圓形軌跡； b. 建立生物腔體組織橫截面的電磁特性參數(shù)模型： 以生物腔體組織的橫截面模型的中心為起始點，將模型等角度劃分為m份，每一份近似 為多層腔體壁組織，對模型施加脈沖磁激勵，并對每一份組織接收到的磁聲信號進行仿真， 成像導管上的超聲探測器所處的成像角度為： 9i = 360(i-l)/m 其中，1 = 1，2，...，111，對應的成像區(qū)域的角度范圍為[0^，0^]，其中0^ = 01-18〇/111，0^ = 0i+18〇/m, 確定每個成像區(qū)域的每層組織的電導率和厚度參數(shù)，形成生物腔體組織橫截面的電磁 特性參數(shù)模型； c. 仿真EMAT-MI成像中聲源產生的過程： ① 應用仿真軟件構造亥姆霍茲線圈，并對亥姆霍茲線圈施加高斯脈沖電流作為激勵 源，產生交變磁場B1，沿亥姆霍茲線圈軸向施加一個穩(wěn)恒勾強磁場B ο，多層生物腔體組織的 橫截面模型同軸放置于亥姆霍茲線圈中間； ② 以空氣為背景域，以自由四面體為單元對亥姆霍茲線圈以及生物腔體組織的橫截面 模型進行劃分，進而采用有限元分析的方法仿真得到生物腔體組織中感應渦流J的分布； ③ 根據(jù)感應渦流J仿真出生物腔體組織的聲源▽·(·/><與)的分布，其中，"X"表示向量 積； d. 根據(jù)聲源的分布仿真生物腔體組織產生的磁聲信號，得到超聲探測器在時刻t、角度 Q1、位置r處接收到的多層生物腔體組織產生的磁聲信號的聲壓pi(r，t)(i = l，2,. . .，m); e. 重建極坐標系下的EMT-MI圖像： ① 根據(jù)步驟d得到的m個磁聲信號的聲壓pi(r，t)(i = l，2,. . .，m)計算出式中，cs是超聲波的波速，^是0-1平面中超聲探測器所在位置(與成像導管在X-Y平面 中的成像角度Gi相對應）； ② 由下式重建出在角度I處的聲源分布： 式中，q是^處的單位矢量；③將上式代入下式，得到腔體組織在角度91處的電導率分布：將〇1&)作為角度01、位置r處的腔體組織橫截面極坐標視圖的灰度值，從而得到極坐標 下的EMAT-MI圖像； f.圖像的坐標轉換： 將極坐標下的EMT-MI圖像轉換為X-Y直角坐標系下的橫截面視圖，具體方法如下： 建立X-Y平面直角坐標系，坐標原點是成像導管中心，水平向右的方向為X軸正方向，垂 直于X軸向上的方向為Y軸正方向，設極坐標系中的一點（j，k)的灰度值為f (j，k)，該點在X-Y坐標系中的對應點的坐標為(jW)，灰度值為g( jW)，其中j e [〇，231]，k e [〇，d]，Jv e [-d，d] ,k7 e [_d，d]，d為極坐標視圖中極徑的最大值，則有：
【文檔編號】G06T11/00GK106023277SQ201610332698
【公開日】2016年10月12日
【申請日】2016年5月18日
【發(fā)明人】孫正, 馬真, 毛娟
【申請人】華北電力大學（保定）