超声二维面阵的三维宽波束小区域快速空化成像方法
技术领域
本发明属于超声检测与超声成像技术领域，具体涉及超声二维平面换能器 阵列，以及基于二维面阵的三维宽波束和复合三维宽波束的小区域快速空化成像 方法。
背景技术
三维超声技术的发展和进步，为非侵入性的诊断技术又开辟了一个新的领 域，它能够克服二维超声空间显像的不足，成为二维超声技术的重要辅助手段。 三维超声的进步体现在能够迅速地对容积图像数据进行储存、处理和显示其三维 立体图像，并且能够得到多平面的图像。
目前三维超声成像系统实现的三维成像通常是由一维线阵或相控阵旋转而 实现，即：先获得传统的二维图像，并记录这些图像；随后重建三维图像；产生 矢状面、管状面和横截面图像，最终显示多个三维成像图。换言之，三维超声是 将连续不同平面的二维图像进行计算机处理，得到一个重建的有立体感的图形。 早期的三维重建一次必须采集大量的二维图像(10～50幅)，并将其存在计算机 内，进行脱机重建和联机显示，单次三维检查的图像数据所需的存储空间达数十 兆字节，成像需要数小时甚至数天时间。而采用二维阵列换能器可对声束实现三 维控制，可以使成像探头不需要机械运动而能够实现实时的三维成像。
近年来三维超声与高速的计算机技术的联合使其具备了更多的临床实用 性：表面成像在80年代首次应用于胎儿；90年代初期开始了切面重建和三个互 交平面成像；容积成像则开始于1991年；1994发展了散焦成像；到1996年开 始了实时超声束跟踪技术；而基于二维阵列的三维超声成像系统最早由美国 Duke大学的Von Ramm和Smith教授研制成功，主要用于心脏的成像，并于1998 年推出了世界上第一台商用实时三维超声成像系统；飞利浦公司于2002年推出 xMatrix矩阵容积探头技术，该技术结合“子阵”和“微通道”可用于大规模阵元超 声平面阵成像。但其每一个阵元后面都附加有一“微通道”模拟电路进行延时控 制，技术复杂，成本较高。
二维阵列与一维阵列相比较，在功能上具有以下优点：能够获得任意平面 的图像，并标明其在空间的方向和位置，有利于对图像进行仔细分析，减少主观 因素干扰；具有精确的体积计算功能，基于一维阵列的超声只能获取一个组织结 构的三个切面，通过三个切面的径线粗略地估测体积，当目标形态不规则时则无 法估计，基于二维阵列的超声可处理多平面资料，模拟出组织的形状，利用特定 的容积计算公式得出体积大小，使体积的测量更为精确，尤其对不规则形器官或 病灶体积的测量更具优越性；能够对感兴趣结构重建三维立体图像，使结果直观。
值得注意的是，二维阵列成像探头受面阵规模限制，其视角范围通常较小， 在大区域成像时会受到一定的限制。增加阵元数目可以增大成像范围，但通常以 大规模增加物理数据通道为代价，其成像技术方案及系统复杂、成本较大。
但是，高强度聚焦超声(HIFU)作用时，对其感兴趣的监控成像目标通常 主要集中在焦点区域，即在一个较小的范围内进行成像，因而可以发挥二维阵列 对于小区域成像的优势。
另一方面，针对HIFU作用时焦点区域空化微泡等具有瞬态物理特性信息 的活动，基于一维线阵的传统平移扫描方式，或切面重建的三维成像方法使得不 同平面的图像存在时间差；此外，传统B模式成像不同扫描线之间也存在时间 差，其帧率一般较低(<100Hz)，所以成像的异步性问题在瞬态物理信息成像 中是不能接受的。法国的Fink小组作者成功地利用平面波发射的概念，即利用 线阵所有阵元同时发射三维宽波束，使得超快速成像的帧频速度可以达到 5000Hz。虽然单次平面波的声透射提供了快速的成像方案，但是，随机多样性 噪声、回波伪像、阴影等的存在，却容易造成图像的质量下降。随后，该小组又 提出基于一维线阵的平面波复合快速超声成像法用于空化成核的检测成像，其以 不同角度相干波复合来提高成像信噪比，但以牺牲了部分成像帧率作为此成像方 案的代价。
因此，如何能在现有超声成像设备基础上，并结合射频数据采集装置，针 对具有瞬态物理信息的空化微泡等活动，进行小区域快速、高信噪比的时空分布 三维成像，仍是本领域的难题之一。
发明内容
本发明的目的在于提供一种超声二维面阵的三维宽波束小区域快速空化成 像方法。
为了实现上述目的，本发明采用了以下技术方案：
1)二维面阵向空化区域发射与二维面阵垂直的单一方向三维宽波束或依次 发射多个互相存在方向差异的三维宽波束；
2)根据目标点位置以及三维宽波束的方向，计算对应回波信号的时间延迟， 然后对回波信号分别进行时间延迟校正，得到校正后的回波信号；
3)采用基于特征空间的最小方差自适应波束合成算法对校正后的回波信号 进行二维波束合成；对于二维面阵发射所述单一方向三维宽波束的情况，经过二 维波束合成得到所述单一方向下的二维波束合成结果；对于二维面阵依次发射所 述多个互相存在方向差异的三维宽波束的情况，经过二维波束合成得到对应方向 下的二维波束合成结果；
4)对所述单一方向下的二维波束合成结果根据空间坐标信息进行插值得到 体数据，根据体数据利用视觉化工具函式库进行面绘制三维成像；或者，将所述 多个互相存在方向差异的三维宽波束对应的二维波束合成结果连续叠加并平均 化，得最终的改进的二维波束合成结果，然后进行插值得到体数据，并利用视觉 化工具函式库进行面绘制三维成像。
对于二维面阵发射所述单一方向三维宽波束的情况，假设二维面阵阵元个 数为M×N，及空间目标点坐标(x,y,z)，其中x为方位角方向坐标，y为俯仰 角方向坐标，z为轴向深度坐标，则接收阵元的回波信号相对于参考阵元的时间 延迟为：
τ mn ( r p ) = ( z + z 2 + ( x - x m ) 2 + ( y - y n ) 2 ) / c ]]>
其中，m＝1,...,M，n＝1,...,N，c为声速；
则校正后的回波信号的矩阵形式为：
![]()
其中，chmn(t)是(xm,yn,0)位置处阵元的回波信号的时域射频数据。
对于二维面阵依次发射所述多个互相存在方向差异的三维宽波束的情况， 假设二维面阵阵元个数为M×N，及空间目标点坐标(x,y,z)，其中x为方位角 方向坐标，y为俯仰角方向坐标，z为轴向深度坐标，接收阵元的回波信号相对 于参考阵元的时间延迟为：
τ x ( m , n , i ) = ( z cos α i + x sin α i + y ) / c + z 2 + ( x - x m ) 2 + ( y - y n ) 2 / c ]]>
τ y ( m , n , j ) = ( z cos β j + x + y sin β j ) / c + z 2 + ( x - x m ) 2 + ( y - y n ) 2 / c ]]>
其中，αi表示相对于x轴第i个方向上三维宽波束的倾斜角度，βj表示相对 于y轴第j个方向上三维宽波束的倾斜角度，m＝1,...,M，n＝1,...,N，c为声速；
则校正后的回波信号的矩阵形式为：
![]()
或
![]()
其中，chmn(t)是(xm,yn,0)位置处阵元的回波信号的时域射频数据。
所述基于特征空间的最小方差自适应波束合成算法包括以下步骤：
1)设目标点沿x轴方向对应的有效阵元数目为K，则校正后的回波信号为：
X x - axis ( t ) = ch 1 ( t ) ch 2 ( t ) . . . ch K ( t ) ]]>
2)将K个有效阵元分为K-L+1个子阵，每个子阵包含L个阵元，并计算 Xx-axis(t)的空间平滑自相关矩阵Rx(t)：
R x ( t ) = 1 K - L + 1 Σ l = 1 K - L + 1 X x - axis l ( t ) X x - axis l ( t ) H ]]>
其中，![]()
表示沿x轴的第l个子阵；
3)用![]()
代替所述自相关矩阵Rx(t)，其中I为单位矩阵， γ＝Δ·trace{Rx(t)}，trace{Rx(t)}为所述自相关矩阵的迹，Δ为空间噪声与信号 功率比；
4)计算最优加权系数wx(t)：
w x ( t ) = R ^ x ( t ) - 1 a a H R ^ x ( t ) - 1 a ]]>
其中，a为单位向量，a＝[1,1,...,1]T；
5)对![]()
进行奇异值分解：![]()
其中S为对角阵，U、V分别 是分解后的左、右矩阵，选取对角阵中特征值大于阈值所对应的向量重构矩阵 Us，计算基于特征空间的最优加权系数![]()
![]()
6)计算x轴目标点的最优输出yx(t)：
![]()
7)重复步骤1至步骤6，沿x轴遍历计算所有目标点的最优输出；
8)根据步骤1至步骤7，沿y轴遍历计算上述所有目标点的最优输出。
优选的，所述面绘制三维成像选用matching cube算法获取等值面而完成。
优选的，二维面阵阵元均匀分布在11行12列的栅格上，阵元总数为128。
优选的，二维面阵尺寸为9.17mm×8.4mm。
对于二维面阵依次发射所述多个互相存在方向差异的三维宽波束的情况， 优选采用五个互相存在方向差异的三维宽波束，分别为相对于二维面阵的x轴方 向以及y轴方向发射倾斜角度为±5°以及0°的三维宽波束。即采用空间的5个角 度三维宽波束进行复合空化成像。
本发明的有益效果体现在：
本发明摒弃了常规超声成像技术所采用的线性扫描方式，而采用了三维宽 波束成像方法，可以克服在具有瞬态物理信息的成像中所存在的时间异步性和运 动模糊问题，其成像帧率可以达到5kHz。本发明为了进一步克服单次发射三维 宽波束快速成像信噪比不足的问题，利用多角度连续三维宽波束进行复合成像。 并举例采用空间5个连续不同角度三维宽波束发射进行空化复合成像，其成像帧 率仍可达到1kHz。
进一步的，本发明针对小区域空化成像，设计一种二维面阵，其通道数为 128，二维面阵阵元数与现有超声成像设备物理通道数一致，可方便与现有常用 超声成像设备直接匹配而不需要额外数据通道，因而不会大幅增加成像方案的复 杂度，以及不会产生较大的成本。且其尺寸大小方便与HIFU探头进行同轴配置， 可以与HIFU探头同时移动，有利于HIFU作用时空化活动的监控成像。
附图说明
图1为举例设计的超声二维面阵结构示意图；(a)二维平面图；(b)三维 立体图；
图2为举例设计的二维平面超声换能器阵列聚焦发射场俯仰角和方位角波 束仿真示意图；(a)为发射声场方位角的声场波束；(b)为发射场俯仰角的声场 波束；
图3为二维平面超声换能器阵列接收回波，其目标点处回波信号的时间延 迟计算示意图；
图4为超声二维面阵三维宽波束小区域快速空化成像方法示意图；
图5为感兴趣的3组随机分布点散射子目标超声二维面阵三维宽波束小区 域快速空化成像仿真结果；(a)视角为正前方斜视-45°，(b)视角为正前方斜视 45°；
图6为超声二维面阵复合三维宽波束改变发射角度进行小区域空化成像示 意图；(a)为与二维面阵垂直的单一方向三维宽波束发射，(b)为偏转角度三维 宽波束发射；
图7为超声二维面阵复合三维宽波束改变角度发射后，其回波信号时间延 迟计算示意图；(a)为沿x轴(方位角方向)偏转角度回波信号延时计算示意图， (b)为沿y轴(俯仰角方向)偏转角度回波信号延时计算示意图；
图8为感兴趣的空化区域，沿x轴(方位角方向)3次不同角度单次发射三 维宽波束成像及其复合成像仿真结果；
图9为感兴趣的空化区域，沿y轴(俯仰角方向)3次不同角度单次发射三 维宽波束成像及其复合成像仿真结果；
图10为感兴趣的空化区域，沿x轴、y轴方向发射5个不同角度三维宽波 束共同复合成像仿真结果。
具体实施方式
下面结合附图和实施例对本发明做详细说明。
本发明提供一种用于小区域空化三维成像的二维平面超声换能器阵列设 计，以及在此基础上结合三维宽波束发射方式，针对具有瞬态物理信息成像的一 种超声二维面阵的三维宽波束小区域快速空化成像方法和复合三维宽波束小区 域快速空化成像方法。
1.超声三维宽波束小区域快速空化成像的二维平面超声换能器阵列(即二维 面阵)设计
二维平面超声换能器阵列由M×N个阵元组成，阵元中心距离为d，焦点与 面阵几何中心垂直距离为F。以图1所示二维面阵为例，其单阵元尺寸为 700μm×700μm，阵元间距为70μm，阵元均匀分布在11行12列栅格上。阵元 总数和数据通道数都为128，以方便和现有超声成像设备的数据通道直接匹配。 二维面阵发射中心频率为2MHz，阵元中心间隔与波长比为1，有效带宽为60％， 焦点距离面阵几何中心30mm。该面阵尺寸为9.17mm×8.4mm，可与HIFU方便 完成同轴位置配置，有利于对HIFU激励空化过程等焦点或其他目标区域进行小 区域的三维空化监控成像。
参阅图2所示，对所设计二维平面超声换能器阵列发射场声场进行仿真。 仿真表明在焦点处(30mm深度)二维面阵可形成聚焦，并且，发射采用二维汉 明窗进行变迹，以降低旁瓣影响。
2.超声二维面阵的三维宽波束小区域快速空化成像方法
二维平面超声换能器阵列的所有阵元发射一次超声波，即二维面阵发射三 维宽波束，可以覆盖空化目标所在的大部分区域，回波信号经处理从而获取整个 空化区域的体数据。在垂直于声束方向上，每一部分的图像数据都是同时采集到 的，因此，极大避免了常规线扫描方式中扫描线之间的采样时间差，进而提高了 成像速度，保证了成像的瞬时性，其成像帧率可以达到5000Hz。
首先二维面阵向成像目标区域，即空化区域，发射周期为125纳秒的单脉 冲超声三维宽波束，其发射方向与面阵垂直，三维宽波束在遇到成像目标后发生 散射，产生与发射波方向相反的回波。这些回波被二维平面超声换能器阵列接收， 利用多通道射频数据采集设备将回波信号采样并转换为数字信号并存储在计算 机硬盘中。随后根据目标点的位置计算有效孔径数目进行接收聚焦，并计算相对 于参考阵元的时间延迟，根据时间延迟对接收聚焦后的射频数据校正。
参阅图3所示，假设二维平面超声换能器阵列阵元个数为M×N，及空间 目标点坐标(x,y,z)，其中x为方位角方向坐标，y为俯仰角方向坐标，z为轴向 深度坐标，则接收阵元的回波信号相对于参考阵元的时间延迟为：
τ mn ( r p ) = ( z + z 2 + ( x - x m ) 2 + ( y - y n ) 2 ) / c - - - ( 1 ) ]]>
其中，m＝1,...,M，n＝1,...,N，c为声速。对于校正后的回波信号，其矩 阵形式为：
![]()
其中chmn(t)是阵元(xm,yn,0)接收信号的时域射频数据。
对经过时间延迟校正后的射频数据进行基于二维特征空间的最小方差自适 应波束合成，实现目标点的二维最优幅度变迹。遍历成像目标区域所有目标点， 得到成像目标区域的二维波束合成结果。然后，对二维波束合成结果根据空间坐 标信息进行插值得到规则体数据，最后利用视觉化工具函式库(Visualization Toolkit,VTK)进行超声面绘制三维成像显示，参阅图4所示。
所述基于二维特征空间的最小方差自适应波束合成和超声面绘制三维成 像，具体步骤为：
1)设目标点沿x轴方向对应的有效孔径(有效阵元)数目为K，经时间延 迟校正后的接收信号为：
X x - axis ( t ) = ch 1 ( t ) ch 2 ( t ) . . . ch K ( t ) - - - ( 3 ) ]]>
2)将K个阵元分为K-L+1个子阵，每个子阵包含L个阵元，并计算 Xx-axis(t)的空间平滑自相关矩阵Rx(t)：
R x ( t ) = 1 K - L + 1 Σ l = 1 K - L + 1 X x - axis l ( t ) X x - axis l ( t ) H - - - ( 4 ) ]]>
其中，![]()
表示沿x轴的第l个子阵，[ ]H表示共轭转置；
3)用![]()
来代替自相关矩阵Rx(t)，从而实现对角线加载。其 中I为单位矩阵，γ＝Δ·trace{Rx(t)}，trace{Rx(t)}为所述自相关矩阵的迹，Δ为 算法加入的空间噪声与信号功率比；
4)计算最优加权系数wx(t)：
w x ( t ) = R ^ x ( t ) - 1 a a H R ^ x ( t ) - 1 a - - - ( 5 ) ]]>
其中，a为单位向量，即a＝[1,1,...,1]T；[ ]T表示转置；
5)对![]()
进行奇异值分解![]()
其中S为对角阵，U、V分别是 分解后的左、右矩阵。选取对角阵中特征值大于阈值(如最大特征值的0.5倍) 所对应的向量重构矩阵Us，计算基于特征空间的最优加权系数![]()
![]()
6)计算x轴目标点的最优输出yx(t)：
![]()
7)重复步骤1-6，沿x轴遍历计算所有目标点的最优输出；
8)沿y轴方向，重复步骤1-7，从而得到基于二维特征空间的最小方差自 适应波束合成后的二维波束合成结果；
9)对二维波束合成结果，根据空间坐标信息进行插值得到规则体数据，然 后设定要抽提的等值面的值v1；
10)在体素(由空间相邻的8个数据点组成的小立方体)中抽提等值面， 比较v1与这8个顶点的数据值的大小，例如对于体素小立方体的一边顶点上的两 个数据v2、v3，若v2<v1并且v3>v1，则可以确定该边上即含有等值点，假设 体数据符合线性，则可以通过线性插值得到等值点的坐标，依次比较体数据的 12条边的8个顶点数据就可以得到该体素内的等值面的形状；
11)基于Marching Cube算法，遍历三维空间体数据内的所有体素就可以得 到要抽提的等值面；
12)在VTK编程环境中，建立可视化管道，将上一步骤得到的等值面多边 形数据输入到渲染环境进行三维显示。
参阅图5所示，仿真为深度25～40mm处随机分布的3组点散射子目标， 利用前述二维面阵对目标进行三维宽波束二维面阵小区域高分辨率高速三维成 像，发射频率2MHz。从图5可以看出，与二维成像相比，三维宽波束二维面阵 小区域高分辨率三维成像能够直观地描述空化区域的空间分布情况。
3.二维面阵的复合三维宽波束小区域快速空化成像方法
参阅图6(a)所示，能够达到二维面阵的快速空化三维成像的最好方法是 向媒介发射与二维面阵垂直的单一方向三维宽波束，然而，该单一方向三维宽波 束的回波成像信噪比较差，会影响成像模式的鲁棒性。因此，提出一种改善三维 宽波束快速空化成像的方法，参阅图6(b)所示，即利用不同倾斜角度、混合 三维宽波束的连续再组合提高成像质量。
参阅图7(a)所示，假设用相对于x轴(方位角方向)的角度α发射三维 宽波束(规定向正半轴方向偏转角度为正)，则发射波到目标点(x,y,z)的时间为：
τex(α,x,y,z)＝(zcosα+xsinα+y)/c (8)
返回位置(xm,yn,0)的阵元时间为：
τ rec ( x m , y n , 0 , x , y , z ) = z 2 + ( x - x m ) 2 + ( y - y n ) 2 / c - - - ( 9 ) ]]>
则总的时间延迟为：
τ(xm,yn,0,x,y,z)＝τex+τrec (10)
其中，c为介质中的声速。
参阅图7(b)所示，相对于y轴(俯仰角方向)的角度β发射三维宽波束， 其计算过程与前述相似：
τex(β,x,y,z)＝(zcosβ+x+ysinβ)/c
τ rec ( x m , y n , 0 , x , y , z ) = z 2 + ( x - x m ) 2 + ( y - y n ) 2 / c - - - ( 11 ) ]]>
τ(xm,yn,0,x,y,z)＝τex+τrec
选择n个三维宽波束，倾斜角度γi,i＝1,2,...,n。则二维面阵的复合三维宽波 束快速空化成像的具体步骤为：
1)二维面阵向空化区域发射倾斜角度为γi的三维宽波束，并利用多通道射 频数据采集设备采集并储存回波信号；
2)根据当前倾斜角度γi计算时间延迟，对回波数据进行校正处理；
3)对校正后的回波数据，利用前述基于二维特征空间的最小方差自适应波 束合成方法得到倾斜角度为γi时的二维波束合成结果；
4)重复步骤1-4，将所有角度的二维波束合成结果连续叠加并平均化来获 得最终的改进的二维波束合成结果；
5)对改进的二维波束合成结果进行超声面绘制三维成像。
请参阅图8所示，假设深度为30mm处存在感兴趣的空化区域(随机点散射 子)，回波信号加入信噪比为60dB的随机噪声，利用前述二维面阵沿x轴(方位 角方向)发射3次不同角度三维宽波束，其角度分别为-5°,0°,5°，中心频率皆为 2MHz，得到沿x轴每个单一角度的三维宽波束快速空化成像结果和复合三维宽 波束空化成像结果(视角为正前方斜视45°)。
请参阅图9所示，在相同仿真条件下，利用二维面阵沿y轴(俯仰角方向) 发射三次不同角度三维宽波束，其角度分别为-5°,0°,5°，得到沿y轴每个单一角 度的三维宽波束快速空化成像结果和复合三维宽波束空化成像结果(视角为正前 方斜视45°)。图8和图9表明，沿某个轴方向的多角度复合三维宽波束空化成像 较单一三维宽波束空化成像能够适当提高成像的信噪比，改进成像质量。
参阅图10所示，针对x轴、y轴5种不同角度三维宽波束共同复合成像的 结果(视角为正前方斜视45°)。在图8、图9所示的关于x轴或y轴不同角度发 射的三维宽波束成像的基础上，图10将关于x轴和y轴不同角度发射的三维宽 波束所得结果进行复合叠加，其增加了信号的信噪比，从而提高成像的质量，但 帧率有所下降。
单一三维宽波束易受随机噪声影响而降低成像质量，由不同角度的几个三 维宽波束得到的复合体数据可以显著提高信噪比，从而降低随机噪声的影响，改 善三维成像的质量。发射n个三维宽波束需要更多的时间，因而其帧率会降低n 倍。如前述，采用单一三维宽波束的帧率可以达到5kHz，那么5个不同角度的 复合三维宽波束发射的帧频仍可保证为1kHz，依然高于B模式的成像帧率。因 此，质量和帧频的制衡可以针对具体成像需求进行优化选取。