一种磁感应磁声内窥图像的建模与仿真方法
【专利摘要】一种磁感应磁声内窥图像的建模与仿真方法，所述方法首先建立生物腔体组织的横截面模型和生物腔体组织横截面的电磁特性参数模型；然后应用有限元分析的方法，对腔体组织在磁激励场中产生感应电流的过程进行仿真，再由感应电流仿真腔体组织在静磁场中产生声源的过程；之后根据声源的分布仿真出腔体组织产生的磁声信号；最后根据磁声信号重建出生物腔体的横截面图像。本发明可灵活地调整生物腔体组织横截面模型，改变模型中腔体组织的变异类型及程度，设置多层腔体壁组织的电磁特性参数和声学参数，精确地仿真声速不均匀的多层生物腔体组织的声场，因此可以为成像算法和图像后处理算法的研究和性能测试等提供足够的样本图像数据。
【专利说明】
-种磁感应磁声内窥图像的建模与仿真方法
技术领域
[0001] 本发明设及一种医学成像技术，特别是对生物腔体的磁感应磁声内窥横截面图像 进行建模与仿真。
【背景技术】
[0002] 生物组织磁感应磁声内窥化ndoscopic Magnetoacoustic Tomography with Magnetic Induct ion, EMAT-MI)成像属于功能成像，其成像的物理基础是生物组织的磁声 效应，即生物组织变异部位与非变异部位的电导率不同，导致产生的感应满流不同，从而在 静磁场中产生不同强度的超声信号（即磁声信号KEMAT-MI结合了生物电阻抗成像的高对 比度和超声扫描成像的高空间分辨率，可得到生物腔体组织（含变异组织）内部的形态结 构、组织类型和成分等信息。然而EMAT-MI成像技术尚处于实验室研究阶段，未大规模应用 于实际，而改善成像仪器的结构和校准参数、优化图像分析和处理算法等都是建立在分析 大量样本图像数据的基础之上的。因此，为了获取足够多的样本图像数据，有必要寻找一种 对生物腔体组织磁感应磁声内窥图像进行建模与仿真的方法。
【发明内容】
[0003] 本发明的目的在于针对现有技术之弊端，提供一种磁感应磁声内窥图像的建模与 仿真方法，为成像和图像后处理的研究和性能测试等提供数据源。
[0004] 本发明所述问题是W下述技术方案解决的：
[0005] -种磁感应磁声内窥图像的建模与仿真方法，所述方法首先建立生物腔体组织的 横截面模型和生物腔体组织横截面的电磁特性参数模型;然后应用有限元分析的方法，对 腔体组织在磁激励场中产生感应电流的过程进行仿真，再由感应电流仿真腔体组织在静磁 场中产生声源的过程;之后根据声源的分布仿真出腔体组织产生的磁声信号；最后根据磁 声信号重建出生物腔体的横截面图像。
[0006] 上述磁感应磁声内窥图像的建模与仿真方法，具体处理按W下步骤进行：
[0007] a.建立生物腔体组织的横截面模型：
[000引成像导管位于生物腔体组织横截面模型的中屯、，接收磁声信号的超声探测器位于 成像导管顶端，将超声探测器看作理想的点换能器，其扫描轨迹为平行于成像平面、半径趋 近于0的圆形轨迹；
[0009] b.建立生物腔体组织横截面的电磁特性参数模型：
[0010] W生物腔体组织横截面模型的中屯、为起始点，将模型等角度划分为m份，每一份近 似为多层组织，对模型施加脉冲磁激励，并对每一份组织接收到的磁声信号进行仿真，成像 导管上的超声探测器所处的角度为
[0011] 目 i = 360(i-l)/m
[001^ 其中，i = l，2, . . .，m，对应的成像区域的角度范围为[0ia，0ib]，其中0ia=0广180/ m，目化=目 i+180/m。
[0013] 确定每个成像区域的每层组织的电导率和厚度参数，形成生物腔体组织横截面的 电磁特性参数模型。W血管的横截面模型为例，多层血管壁组织的参数设置示例如表1所 列。
[0014] C.仿真EMAT-MI成像中声源产生的过程：
[0015] ①应用仿真软件构造亥姆霍兹化elmholtz)线圈，并对亥姆霍兹线圈施加高斯脉 冲电流作为激励源，产生交变磁场Bi,沿亥姆霍兹线圈轴向施加一个稳恒匀强磁场Bo,多层 生物腔体组织的横截面模型同轴放置于亥姆霍兹线圈中间;稳恒磁场和亥姆霍兹线圈的参 数设置示例如表2所列。
[0016] ②W空气为背景域，W自由四面体为单元对亥姆霍兹线圈W及生物腔体组织的横 截面模型进行划分，进而采用有限元分析的方法仿真得到生物腔体组织中感应满流J的分 布；
[0017] ③根据感应满流J仿真出生物腔体组织的声源▽?(/><巧,)的分布，其中，"X"表示 向量积；
[0018] d.根据声源的分布仿真生物腔体组织产生的磁声信号，得到超声探测器在时刻t、 角度9i、位置r处接收到的多层生物腔体组织产生的磁声信号的声压pi(r，t)a = l，2，...， m)；
[0019] e.重建极坐标系下的EMAT-MI图像：
[0020] ①根据步骤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] 本发明可灵活地调整生物腔体组织横截面模型，改变模型中腔体组织的变异类型 及程度，设置多层腔体壁组织的电磁特性参数和声学参数，精确地仿真声速不均匀的多层 生物腔体组织的声场，因此可W为成像算法和图像后处理算法的研究和性能测试等提供足 够的样本图像数据。
【附图说明】
[0036] 图1是含有脂质斑块的血管横截面模型示例；
[0037] 图2是EMAT-MI成像导管上的超声探测器在角度01处接收磁声信号的示意图；
[0038] 图3是将图1中的血管横截面模型等角度划分后，将其中的一份近似为多层组织的 示意图；
[0039] 图4是0-1坐标系中的多层生物腔体组织的示意图；
[0040] 图5是亥姆霍兹化eImholtz)线圈激励示意图，图中R是线圈半径。
[0041 ]表1是多层血管壁组织参数设置举例。
[0042] 表2是稳恒磁场和亥姆霍兹线圈的参数设置示例。
[0043] 文中各符号为:X、Y、腔体横截面模型所在的X-Y平面直角坐标系的横轴和纵轴，其 中，成像导管的中屯、为坐标原点，水平向右的方向为X轴正方向，垂直于X轴向上的方向为Y 轴正方向;m、横截面模型被等角度分割的总份数;01、成像导管的第i个成像角度，其中，i = 1，2,...，m;目13、目lb、成像导管在角度目1处进行成像时对应的角度范围的上、下限；目、1、目-1坐 标系的横轴和纵轴;Bi、激励磁场的磁通密度;Bo、稳恒磁场的磁通密度;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处的激励磁场和稳 恒磁场的叠加磁通密度，其中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] 下面结合附图对本发明作进一步详述。
[0045] 本发明方法的步骤包括：
[0046] (1)建立多层生物腔体组织横截面的形态模型：
[0047] W血管横截面形态模型为例，如附图1所示，模型包括成像导管(接收磁声信号的 超声探测器位于成像导管顶端）、管腔、血管壁内膜/中膜、外膜和斑块(脂质斑块、纤维斑 块、巧化斑块或混合斑块)五部分。其中成像导管位于模型中屯、，周围由内向外依次为管腔、 斑块、血管壁内膜/中膜和外膜。可根据组织变异的类型（如巧化、纤维化或者脂质斑块)和 大小，W及血管内腔、血管壁内膜/中膜、外膜的厚度建立不同的血管横截面形态模型。本发 明方法忽略超声探测器的孔径效应，将其看作理想的点换能器，其扫描轨迹为平行于成像 平面、半径趋近于0的圆形轨迹。模型所在的坐标系为X-Y平面直角坐标系，其中坐标原点是 成像导管中屯、，水平向右的方向为X轴正方向，垂直于X轴向上的方向为Y轴正方向。
[004引（2)建立多层生物腔体组织横截面的电磁特性参数模型：
[0049] 首先，如附图2所示，W生物腔体组织横截面形态模型的中屯、为起始点，将模型等 角度划分为m份，对模型施加脉冲磁激励，并对每一份组织接收到的磁声信号进行仿真。成 像导管上的超声探测器所处的成像角度为
[0050] 目 i = 360(i-l)/m (1)
[0051] 其中，i = l，2，. . .，m。对应的成像区域的角度范围为[0ia，0ib]，其中0ia=0广180/ m，目ib =目i+180/m。例如，当m=360时，将附图1中含有脂质斑块的血管横截面形态模型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轴向上的方向，表示多层组织 的厚度。确定每个成像角度对应的组织中每层的电导率和厚度参数，形成多层生物腔体组 织的电磁特性参数模型。W血管的横截面模型为例，多层血管壁组织的参数设置示例如表1 所列。不同成像角度对应的多层生物腔体组织仅厚度参数不同，图4中的生物腔体组织分为 a、b、c、d、e五层，若形态模型中不包含斑块，则c、d层的厚度参数为0。
[0化3] (3)仿真EMAT-MI成像中声源产生的过程：
[0054] 本发明采用有限元分析的方法进行声源的仿真，具体步骤如下：
[0055] 首先，应用COMSOL软件构造亥姆霍兹化elmholtz)线圈，如附图5所示，并施加高斯 脉冲电流作为激励源产生交变磁场Bi,沿线圈轴向施加一个稳恒匀强磁场Bo,多层生物腔体 组织模型置于线圈中间，其轴向与亥姆霍兹线圈的轴向同向；然后，W空气为背景域，W自 由四面体为单元对亥姆霍兹线圈W及生物腔体组织模型进行划分，进而仿真得到生物腔体 组织中感应满流J的分布;最后，根据感应满流J仿真出生物腔体组织的声源▽?(/>< A)的分 布，其中，"X"表示向量积。稳恒磁场和亥姆霍兹线圈的参数设置示例如表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处稳恒磁场的磁通密度。
[0064] 生物组织的声学特性与流体相似，在理想流体中，声波的传播满足牛顿第二定律、 质量守恒定律和物态方程。设随时间变化的激励电流信号为n(t)，则由声波的连续性方程、 运动方程和物态方程(张海欄.理论声学.北京:高等教育出版社.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] 根据式(4)即可仿真得到超声探测器在时刻t、角度01、位置r处接收到的多层生物 腔体组织产生的磁声信号的声压91^，*)(1 = 1，2，...，111)。
[OOW (5)重建极坐标下的EMAT-MI图像：
[0072] 图像重建的实质是由超声探测器接收到的磁声信号计算出生物组织的电导率分 布。具体方法如下：
[0073] 首先，根据步骤(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. 一种磁感应磁声内窥图像的建模与仿真方法，其特征是，所述方法首先建立生物腔 体组织的横截面模型和生物腔体组织横截面的电磁特性参数模型;然后应用有限元分析的 方法，对腔体组织在磁激励场中产生感应电流的过程进行仿真，再由感应电流仿真腔体组 织在静磁场中产生声源的过程;之后根据声源的分布仿真出腔体组织产生的磁声信号；最 后根据磁声信号重建出生物腔体的横截面图像。2. 根据权利要求1所述的磁感应磁声内窥图像的建模与仿真方法，其特征是，所述方法 包括以下步骤： a. 建立生物腔体组织的横截面模型： 成像导管位于生物腔体组织的横截面模型的中心，接收磁声信号的超声探测器位于成 像导管顶端，将超声探测器看作理想的点换能器，其扫描轨迹为平行于成像平面、半径趋近 于O的圆形轨迹； b. 建立生物腔体组织横截面的电磁特性参数模型： 以生物腔体组织的横截面模型的中心为起始点，将模型等角度划分为m份，每一份近似 为多层腔体壁组织，对模型施加脉冲磁激励，并对每一份组织接收到的磁声信号进行仿真， 成像导管上的超声探测器所处的成像角度为： 9i = 360(i-l)/m 其中，1 = 1，2，...，111，对应的成像区域的角度范围为[0^，0^]，其中0^ = 01-18〇/111，0^ = 0i+18〇/m, 确定每个成像区域的每层组织的电导率和厚度参数，形成生物腔体组织横截面的电磁 特性参数模型； c. 仿真EMAT-MI成像中声源产生的过程： ① 应用仿真软件构造亥姆霍兹线圈，并对亥姆霍兹线圈施加高斯脉冲电流作为激励 源，产生交变磁场B1，沿亥姆霍兹线圈轴向施加一个稳恒勾强磁场B ο，多层生物腔体组织的 横截面模型同轴放置于亥姆霍兹线圈中间； ② 以空气为背景域，以自由四面体为单元对亥姆霍兹线圈以及生物腔体组织的横截面 模型进行划分，进而采用有限元分析的方法仿真得到生物腔体组织中感应涡流J的分布； ③ 根据感应涡流J仿真出生物腔体组织的声源▽·(·/><與)的分布，其中，"X"表示向量 积； d. 根据声源的分布仿真生物腔体组织产生的磁声信号，得到超声探测器在时刻t、角度 Q1、位置r处接收到的多层生物腔体组织产生的磁声信号的声压pi(r，t)(i = l，2,. . .，m); e. 重建极坐标系下的EMT-MI图像： ① 根据步骤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日
【发明人】孙正, 马真, 毛娟
【申请人】华北电力大学（保定）