一种基于特征结构的广义旁瓣相消超声成像波束合成方法

摘要

本发明涉及一种基于特征结构的广义旁瓣相消超声成像波束合成方法，该方法首先对接收阵元的采样信号进行延时和前向平滑处理，得到前向空间协方差矩阵估计，并对其进行对角加载后与方向向量结合，同时基于最小方差准则，计算出最优的自适应加权矢量，并构造与最小方差等效解的广义旁瓣相消器，获得自适应与非自适应两部分加权矢量；然后对接收信号的前向空间协方差矩阵估计进行特征分解，构建特征阈值信号子空间，并将整体的加权矢量投影到特征阈值信号子空间的左奇异矢量空间中，得到新的自适应波束加权矢量；最后利用投影所获得的加权矢量与超声阵元接收数据进行点乘运算得到超声波束合成数据；本方法能够去除干扰较大的特征矢量，使得超声图像的抗噪能力更强、分辨率更好。

说明书

技术领域

本发明属于超声成像技术领域，涉及一种基于特征结构的广义旁瓣相消超声成像波束合 成方法。

背景技术

超声具有对人体无害，实时成像，成本较低等优点，目前广泛应用于血管、腹腔、心脏 等器官的医学检测领域。超声成像算法是超声成像的关键技术。延时叠加是目前最简单，最 传统，也是应用最为广泛的一种成像技术。这种方法的成像速度快，但分辨率较低，且栅瓣 严重。

目前广泛应用的是采用最小方差自适应波束形成方法来克服延时叠加算法的这一缺点， 它是根据保持期望方向增益不变，并且使阵列输出能量最小的原则，计算出聚焦延时后信号 的加权矢量。由于该方法是实时根据回波数据计算加权值，所以该方法可以有效地降低旁瓣 等级,从而提高横向分辨率，同时可以显著的提高图像的对比度，但是该方法的缺点是稳健性 远不如传统的延时叠加算法，当存在噪声时，抗噪声能力很差。为了提高最小方差波束合成 算法的稳定性，很多学者提出了对角加载和特征空间法。但是该方法在噪声较强时便失去了 其优势。

在超声成像中，根据最小方差准则计算出最优加权矢量然后进行成像时，虽然图像分辨 率有所提升，但其对比度和鲁棒性均有不足，同时也会使得部分有用信号相消，而且分辨率 也有很大的提升空间。广义旁瓣相消实质是最小方差的一种等效结构，利用普通的广义旁瓣 相消成像时图像质量提升不多。特别是当超声回波数据中存在噪声时，最小方差和广义旁瓣 相消算法的鲁棒性下降较大。

综上所述，有必要提出一种自适应算法，在能提高超声图像横向分辨率的同时，还能保 持算法的稳定性和较好的抗噪能力。

发明内容

有鉴于此，本发明的目的在于提供一种基于特征结构的广义旁瓣相消超声成像波束合成 方法，该方法基于最小方差准则，构造广义旁瓣相消器，将自适应加权矢量分解成自适应与 非自适应两部分加权矢量，然后根据接收信号的前向空间协方差矩阵估计构建特征阈值信号 子空间，并将加权矢量投影到此空间的左奇异矢量空间中，获得新的加权矢量。

为达到上述目的，本发明提供如下技术方案：

一种基于特征结构的广义旁瓣相消超声成像波束合成方法，包括以下步骤：

步骤一：对超声阵元接收的回波信号进行放大处理与AD转换，以获得超声成像所需的 回波数据；

步骤二：将接收阵列依次划分为具有一个重叠阵元的子阵，然后对相应划分子阵接收的 回波信号进行前向平滑处理，以获得样本协方差矩阵；

步骤三：对样本协方差矩阵进行对角加载，以获得前向空间协方差矩阵估计；

步骤四：将空间协方差矩阵估计与方向向量结合，计算得到自适应波束加权矢量；

步骤五：根据最小方差准则，构造广义旁瓣相消器，将加权矢量进行分解以获得自适应 与非自适应两部分加权矢量；

步骤六：对空间协方差矩阵估计进行特征分解，构建特征阈值信号子空间；

步骤七：将整体的加权矢量投影到特征阈值信号子空间的左奇异矢量空间中，得到新的 自适应波束加权矢量；

步骤八：最后利用投影所获得的加权矢量与阵元回波数据进行点乘运算得到超声成像的 波束合成数据。

进一步，在步骤二中，具体包括以下步骤：

1)令步骤一中得到的聚焦延时后的超声回波信号为x(k),表示为x(k)＝[x₁(k),x₂(k),…, x_N(k)]，其中N表示超声阵列的阵元个数，k为对应于深度的采样时刻；

2)把N阵元阵列依次划分为具有一个重叠阵元的L个子阵，分别计算各个子阵的样本 协方差矩阵R_l(k)，然后计算所有子阵的样本协方差矩阵的平均公式如下：

$\tilde{R}\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{R}_l\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{x}_d^l\left(k\right)x_d^l{\left(k\right)}^H$

公式中$x_d^l\left(k\right)={[x_l\left(k\right),x_{l+1}\left(k\right),...,x_{l+L-1}\left(k\right)]}^T$为第l个子阵的输出向量，表示的共 轭转置。

进一步，在步骤三中，通过下列公式对平均的样本协方差矩阵进行对角加载得到 对角加载后的前向空间协方差矩阵估计

${\tilde{R}}_{DL}=\tilde{R}\left(k\right)+\gamma I$

其中，△为加入的空间噪声与信号功率比，为信号等效功 率，I表示单位矩阵。

进一步，在步骤五中，根据如下步骤构造广义旁瓣相消器，获得分解的自适应与非自适 应加权矢量：

1)先通过最小方差原理计算自适应加权矢量w为：

$w=\frac{{\tilde{R}}_{DL}^{-1}a}{a^H{\tilde{R}}_{DL}^{-1}a}$

其中，a为方向向量,表示前向空间协方差矩阵估计的逆运算；

2)构造与此解等效的广义旁瓣相消器，将加权矢量w分解成自适应与非自适应加权矢 量部分：

w＝w_q-Bw_a s.tB^Hw_q＝0

其中，w_q为非自适应权矢量，位于约束子空间中，对应的加权系数全部为1； 为自适应权矢量，B矩阵为L×(L-1)阻塞矩阵，根据约束条件 B^Hw_q＝0计算获得，B^H是B矩阵的共轭转置。

进一步，在步骤六中，根据如下步骤得到特征阈值信号子空间：

对前向空间协方差矩阵估计进行特征分解，并且将得到的特征值λ_i从大到小排列， 计算相邻两个特征值的比值λ_i/λ_i+1，当比值大于设定阈值T_h时，则构成特征阈值信号子空间 R_s＝[a₁,a₂,…a_i,a_i+1]，否则，R_s＝[a₁,a₂,…a_i]，其中a_i为相应特征值λ_i对应的特征向量，其 余特征向量构成干扰信号子空间。

进一步，在步骤七中，根据如下步骤将整体的加权矢量投影到特征阈值信号子空间的左 奇异矢量空间中，以获得新的加权矢量：

1)根据如下公式对特征阈值信号子空间进行奇异值分解：

R_s＝UDV^H

其中，U是左奇异矩阵，D为奇异值在对角线上的对角矩阵，V为右奇异矩阵，V^H是V 的共轭转置；

2)将步骤五中获得的权矢量w按如下公式投影到R_s较大特征值所对应的左奇异矢量空 间中：

$w_{ES}=U_s{U}_s^{H}w=U_s{U}_s^H{w}_q-U_s{U}_s^H{\mathrm{Bw}}_a$

其中，w_ES为投影后新的自适应加权矢量，U_s为较大特征值对应的左奇异矩阵，为 较大特征值对应的左奇异矩阵的共轭转置。

进一步，在步骤八中，通过如下公式利用投影所获得的加权矢量与阵元回波数据进行点 乘运算得到超声成像的波束合成数据：

$y\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{w}_{ES}^H\left(k\right)x_d^l\left(k\right)$

其中，y(k)表示计算采样点k时刻的自适应波束信号，表示w_ES的共轭转置，表 示第l个子阵的输出向量。

本发明的有益效果在于：本发明所述方法引入了期望信号的方向矢量，并且能去除干扰 较大的特征矢量，从而在期望信号功率较小或大致等于噪声功率时，可以保持本算法的性能， 对噪声具有良好的鲁棒性，使得超声图像的抗噪能力更强、分辨率和对比度更好，从而显著 提高超声图像的质量。

附图说明

为了使本发明的目的、技术方案和有益效果更加清楚，本发明提供如下附图进行说明：

图1为广义旁瓣相消结构示意图；

图2为本方法的点目标成像与传统的DAS，以及MV-CF算法的成像效果图；

图3为点目标在纵向深度为60mm处的横向截面图；

图4为本发明所做的斑散射目标成像与传统的DAS以及MV-CF算法的成像效果图；

图5为本发明所述方法的流程示意图。

具体实施方式

下面将结合附图，对本发明的优选实施例进行详细的描述。

图5为本发明所述方法的流程示意图，如图所示，本方法包括以下步骤：步骤一：对超 声阵元接收的回波信号进行放大处理与AD转换，以获得超声成像所需的回波数据；步骤二： 将接收阵列依次划分为具有一个重叠阵元的子阵，然后对相应划分子阵接收的回波信号进行 前向平滑处理，以获得样本协方差矩阵；步骤三：对样本协方差矩阵进行对角加载，以获得 前向空间协方差矩阵估计；步骤四：将空间协方差矩阵估计与方向向量结合，计算得到自适 应波束加权矢量；步骤五：根据最小方差准则，构造广义旁瓣相消器，将加权矢量进行分解 以获得自适应与非自适应两部分加权矢量；步骤六：对空间协方差矩阵估计进行特征分解， 构建特征阈值信号子空间；步骤七：将整体的加权矢量投影到特征阈值信号子空间的左奇异 矢量空间中，得到新的自适应波束加权矢量；步骤八：最后利用投影所获得的加权矢量与阵 元回波数据进行点乘运算得到超声成像的波束合成数据。

具体步骤如下：

步骤1、令经过延时后的采样信号为x(k)＝[x₁(k),x₂(k),…,x_N(k)]，其中，N表示传感器阵 列的阵元个数。

步骤2、把N阵元的阵列依次划分为具有一个重叠阵元的L个子阵，分别计算各个子阵 的样本协方差矩阵R_l(k)，然后计算所有子阵的样本协方差矩阵的平均公式如下：

$\tilde{R}\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{R}_l\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{x}_d^l\left(k\right)x_d^l{\left(k\right)}^H$

其中，$x_d^l\left(k\right)={[x_l\left(k\right),x_{l+1}\left(k\right),...,x_{l+L-1}\left(k\right)]}^T$为第l个子阵的输出向量，表示的共轭转置；

步骤3、对样本协方差矩阵进行对角加载，以获得前向空间协方差矩阵估计

具体的对角加载公式为：

${\tilde{R}}_{DL}=\tilde{R}\left(k\right)+\gamma I$

其中，△为加入的空间噪声与信号功率比，为信号等效功率，I 表示单位矩阵。

步骤4、根据最小方差的约束方程$\begin{array}{cc}w={\mathrm{argminw}}^H{\tilde{R}}_{DL}w& s.{\mathrm{tw}}^{H}a=1\end{array}>$

解得自适应加权矢量为：

其中为接收信号的前向空间协方差矩阵估计，a为方向向量，为前向空间协方差矩阵 的逆运算。

步骤5、如图1所示，构造广义旁瓣相消器，根据等效解原则将权矢量分解为自适应权 矢量和非自适应权矢量两部分，图中上半支路w_q为非自适应权矢量部分，下半支路中，w_a为自适应权矢量部分，矩阵B为阻塞矩阵，只允许非期望信号进入。此时加权向量可以写成：

w＝w_q-Bw_a S.t B^Hw_q＝0

步骤6、对获得的前向空间协方差矩阵估计进行特征值分解，获取相应的特征值和特征 向量，并且将其特征值从大到小排列。计算相邻的两个特征值之比λ_i/λ_i+1，假如此比值大于 某个阈值，那么就构成期望阈值信号子空间R_s：

R_s＝[a₁,a₂,…a_i,a_i+1]，否则，R_s＝[a₁,a₂,…a_i]。

步骤7、对获取的期望阈值信号子空间R_s进行奇异值分解，得到：

R_s＝UDV^H

其中，U是左奇异矩阵，D为奇异值在对角线上的对角矩阵，V为右奇异矩阵，V^H是V的共 轭转置。

步骤8、将整体的权矢量w投影到R_s的较大特征值对应的左奇异矢量空间中，得 到本方法的自适应权矢量：

$w_{ES}=U_s{U}_s^{H}w=U_s{U}_s^H{w}_q-U_s{U}_s^H{\mathrm{Bw}}_a$

其中，w_ES为投影后新的自适应加权矢量，U_s为较大特征值对应的左奇异矩阵，为较大特 征值对应的左奇异矩阵的共轭转置。

步骤9、最后利用投影所获得的加权矢量与阵元数据进行点乘运算得到成像回波数据， 最后其波束合成输出为：

$y\left(k\right)=\frac{1}{N-L+1}\underset{l=1}{\overset{N-L+1}{\Sigma }}{w}_{ES}^H\left(k\right)x_d^l\left(k\right)$

其中，y(k)表示计算采样点k时刻的自适应波束信号，表示w_ES的共轭转置，表 示第l个子阵的输出向量。

为了验证该发明的有效性，利用Field II对医学成像中常用的点散射目标和斑散射目标进 行成像。Field II是基于线性系统空间响应的原理，它的仿真效果与实际成像很接近，已被国 际上广泛认同为仿真超声系统的标准。设置14个目标点，间隔5mm，均匀分布在40mm-70mm 的深度范围内，同时在接收数据中加入10dB的噪声，超声图像的动态范围均为40dB。同时 设定半径为3mm的圆形吸声斑，深度在20-35mm范围之间，设定图像对比度范围为60dB， 并在接收的回波数据中加入20dB的噪声。使用3种算法分别进行成像，并比较3种成像算 法的抗噪性能，对比度和分辨率。

图2给出了本发明点目标成像与原始的算法DAS，MV-CF，ES-GSC的对比图像。图3给 出了本发明算法生成的点目标成像图中在深度60mm处的横向截面图与现有的算法相比.由图 分析可知：从成像结果图2、图3中可以看出传统DAS算法成像质量最差，旁瓣等级最高， 抗噪性能最差，MV-CF旁瓣等级降低，横向分辨率有所提高，但是抗造性能基本没有提高， 本发明提出ES-GSC算法在分辨率、对比度和鲁棒性上均有所提升，ES-GSC算法成像质量最 好，而且近场和远场区域的图像分辨率和对比度都得到了较大的改善。

图4给出了本发明算法生成的斑散射目标成像图与DAS，MV-CF，ES-GSC的比较图。 表1列出了各种算法的对比度(CR)。由表1分析可知，ES-GSC算法中心平均功率相对于 DAS降低了7dB，与MV-CF算法相近；ES-GSC算法对比度有略微提高，与MV-CF相当。 相对于DAS算法提高了4.4dB。背景区域方差表征算法的鲁棒性，可以看出DAS算法鲁棒 性最好，本专利提出的算法相对于MV-CF背景区域方差下降，因此，本专利提出的ES-GSC 算法具有更好的鲁棒性。

最后说明的是，以上优选实施例仅用以说明本发明的技术方案而非限制，尽管通过上述 优选实施例已经对本发明进行了详细的描述，但本领域技术人员应当理解，可以在形式上和 细节上对其作出各种各样的改变，而不偏离本发明权利要求书所限定的范围。