基于双高斯衰减模型的超声回波处理方法*

王大为 王召巴 陈友兴 李海洋 王浩坤

1) (中北大学信息与通信工程学院,太原 030051)

2) (山西师范大学物理与信息工程学院,临汾 041000)

信号降噪与特征提取是超声检测数据处理的关键技术.基于超声信号有特定结构而噪声和超声信号的结构无关,本文提出一种旨在解决强噪声背景下超声回波的参数估计和降噪问题的方法.该方法将超声回波的参数估计和降噪问题转换为函数优化问题,首先根据工程经验建立超声信号的双高斯衰减数学模型,然后根据观测回波和建立的超声信号模型确定目标函数,接着选择人工蜂群算法对目标函数进行优化从而得到参数的最优估计值,最后由估计出的参数根据建立的超声信号数学模型重构出无噪的超声估计信号.通过仿真和实验表明本文方法可以准确估计出信噪比大于—10 dB的含噪超声回波中的无噪信号,且效果优于基于自适应阈值的小波降噪方法和经验模态分解方法; 此外相比常用的指数模型和高斯模型,本文提出的双高斯衰减超声信号模型与实测超声信号更接近,其均方误差为9.4 × 10—5,波形相似系数为0.98.

1 引言

超声检测是重要的无损检测技术之一,也是目前国内外使用频率最高且发展较快的一种无损检测技术,在工业生产等实践中得到了广泛应用[1,2].信号降噪与特征提取是超声检测数据处理的关键技术,同时也是超声无损检测的核心环节,其性能优劣直接影响着无损检测结果的表征,因此对其进行研究有重要的学术意义和实用价值[3].Wang等[4]为评价点焊接头强度,首先利用快速傅里叶变换和小波包变换对检测到的超声信号进行处理,然后利用数理统计方法提取超声检测数据的时域、频域和小波域特征,再用支持向量机的分类器对提取的特征进行分类.Meng等[5]为分类碳纤维增强复合材料中的孔洞和分层,首先对采集到的超声信号进行小波分解,然后利用深度卷积神经网络从小波系数中学习每个信号紧凑而有效的表示,最后利用支持向量机的方法对学习结果进行分类.孙灵芳等[6]针对换热管薄层污垢超声检测信号降噪问题提出了一种基于改进完备总体经验模态分解的小波收缩阈值信号处理方法.该方法首先引入夹角余弦计算原始信号和固有模态函数的相似程度,判断信号和噪声主导模态的分界点,然后利用小波收缩阈值方法提取噪声主导模态中的细节信息,最后重构得到降噪后的超声信号.这些方法的共同特点是先对检测到的超声回波降噪,然后再进行特征提取和表征.对于信噪比相对较高的超声回波,这些方法可以取得期望的处理效果,但对于信噪比较低的超声回波,例如回波中信号被噪声淹没时,这些方法难以对超声回波中的特征进行提取和表征.为解决该问题,本文作者在前期工作中提出了一种强噪声背景下微弱超声信号提取方法[7],然而该方法对实测超声信号处理结果的优劣依赖于所建立的超完备字典中是否包含能够和实测超声信号相匹配的原子.为增强对实测超声信号的处理能力,本文提出双高斯衰减超声信号模型,并基于该模型提出一种超声回波处理方法.该方法将超声信号处理问题转换为对含参函数的优化问题,通过利用人工智能优化算法对目标函数求解而实现从观测到的含噪回波中估计无噪信号参数.

2 问题描述

在超声脉冲回波检测中,观测到的超声回波y(t)是测试样本中各界面反射回波sk(t)的线性组合[8],即

其中sk(t)为发射脉冲遇到第k个界面后的反射回波信号,ak为对应的反射系数,ξ(t)表示来自于测试系统和被测材料中的噪声.由于实际测量到的观测回波y(t)中通常包含噪声,本文主要解决如何从包含噪声的观测回波中提取超声信号特征参数及重构无噪超声信号.根据工程经验每一个超声回波sk(t)都是由超声传感器调制发射的,故其有特定的结构并可以表示为一组含参数的非线性函数;而噪声ξ(t)通常是一个和超声信号结构无关的高斯随机过程[9-11]因此,在超声无损检测中由观测回波估计出无噪原始信号的问题可以转化成函数在参数集上的优化问题.本文通过四个环节来解决该问题,即建立超声信号模型、确定目标函数、优化目标函数并估计参数和重构信号.流程图如图1所示.

3 本文方法

3.1 超声回波模型

图1 超声信号处理流程图Fig.1.Flow chart of ultrasonic signal processing.

在超声检测中,超声回波通常是一个被超声换能器中心频率调制的宽带信号.当前,超声回波数学模型有指数模型[10]和高斯模型[11],其包络分别按照指数规律和高斯规律变化.指数模型和高斯模型在超声换能器特性模拟和超声信号数值分析方面有重要作用,并得到了广泛应用.但实际中超声换能器发射的超声脉冲包络并非严格的指数包络或高斯包络,而是产生了非线性失真.为提高本文方法对实测超声信号参数估计的准确性,希望建立的超声信号数学模型包络尽可能地接近超声换能器发射脉冲的包络.基于通过两个不同中心位置的高斯窗函数叠加来模拟高斯模型和发射脉冲包络之间存在的非线性失真,本文提出超声信号的双高斯衰减模型,其数学表达式如下：

式中参数为 θ=[A,a,b,T,τ,fc,φ].其中 A 是超声回波幅度; a,b和T与超声换能器特性有关,决定超声回波的包络,因此定义为包络参数; 包络参数有直观的物理意义,a,b分别是第一个和第二个高斯包络的中心位置,T是高斯包络的尺度参数,其决定着高斯包络的时宽; τ为渡越时间; fc是超声发射脉冲的中心频率,φ是初相位.双高斯衰减模型是高斯模型的拓展,当a=b时,该模型退化为高斯模型.不同包络参数[a,b,T ]对应的双高斯包络如图2所示.

从图2中可以看出,本文提出的双高斯衰减模型的包络结构具有可塑性,本文正是利用该特性模拟超声换能器发射脉冲相对于高斯模型的非线性失真.

3.2 目标函数

目标函数的选择直接影响模型参数估计的准确性.文献[8]中采用内积作为目标函数评价超声回波和字典中原子的相似度; 文献[10]中以均方误差函数作为目标函数采用指数模型估计超声回波渡越时间.在本文的应用背景中以内积或均方误差函数作为目标函数存在目标函数最优收敛值随被处理回波信噪比变化的问题,不利于判断优化算法是否正确收敛.余弦相似度是一种有效衡量数据序列相似性的方法[12],本质上是对内积的归一化; 此外,以余弦相似度作为目标函数,则目标函数值是否为1可成为判断优化算法是否正确收敛的标准.设f,g分别为两个长度相同的离散时间序列,则余弦相似度定义为

图2 超声回波模型包络 (a) [a,b,T ]=[1.5,2.5,1]; (b) [a,b,T ]=[3,3,1]; (c) [a,b,T ]=[1.7,3.2,1]; (d) [a,b,T ]=[1,3,1]Fig.2.Envelope of ultrasonic echo model： (a) [a,b,T ]=[1.5,2.5,1]; (b) [a,b,T ]=[3,3,1]; (c) [a,b,T ]=[1.7,3.2,1]; (d) [a,b,T ]=[1,3,1].

当两个序列有相同的形状时余弦相似度为1,当两个序列无关时余弦相似度为0.在[0,1]区间上,余弦相似度值越大表明两个序列形状越接近.设由测试系统观测到的回波为y(t),待估计的超声信号为s(θ; t),则观测回波和待估计超声信号的余弦相似度可表示为

为分析问题简单起见,假设y(t),s(θ; t)是经过范数归一化处理的,同时利用内积运算的线性性质得

由于超声信号有特定的结构,而噪声和超声信号的结构无关,所以理论上噪声和超声信号不相关,故等式右边第二项为零,即

其中bk是范数归一化系数.因此由观测回波y(t)求超声信号s(θ; t)可以转化为在参数空间θ=[A,a,b,T,τ,fc,φ]上对(8)式所示目标函数的最大值优化,

综上,采用余弦相似度作目标函数时,从包含多个回波的观测回波中提取超声信号特征参数及重构无噪超声信号的问题可以转化为对多个单个回波参数估计问题的线性组合.而对从单个回波中提取超声信号特征参数的问题即是在θ张成参数空间中寻找一组参数使得目标函数(9)式取最大值,

显然,该目标函数的理论最大值为1.

3.3 目标函数优化及信号重构

很多优化算法都可以用来求解本文的目标函数,例如模拟退火算法[13-15]、粒子群优化算法[16-18]、人工蜂群算法[19-21]等.不失一般性,选择人工蜂群(artificial bee colony,ABC)算法优化本文目标函数并以目标函数作为评价每个解的适应度函数,其基本过程如下.

首先在由θ张成的参数空间中随机生成初始解xi(i=1,2,···,SN ),SN是解的个数,每个解xi是一个与θ结构相同的d维向量.根据(10)式计算每个解的适应度值.初始化之后整个种群将进行引领蜂、跟随蜂和侦查蜂搜寻过程的重复循环,直到达到设定的最大迭代次数,

1) 引领蜂阶段

引领蜂对当前解进行邻域搜索,产生新解并贪婪选择较优解.每个引领蜂由(11)式产生一个新解,

其中,i∈{ 1,2,···,SN },j∈{ 1,2,···,d },k 是随机产生的整数,k∈{ 1,2,···,SN }且 k≠i; φij为[—1,1]之间的随机数.检验新解的适应度值并根据贪心原则保留优质解,即若新解的适应度值比旧解大,则引领蜂保留新解舍弃旧解; 反之,保留旧解.在所有引领蜂完成搜索过程之后会把解的信息及适应度值信息分享给跟随蜂.

2) 跟随蜂阶段

跟随蜂根据(12)式计算每个解的选择概率,

然后在[0,1]区间内产生一个随机数,如果解的概率大于该随机数,则跟随蜂由(11)式产生一个新解,并根据贪心原则保留适应度值大的解.

3) 侦查蜂阶段

在所有跟随蜂完成搜索过程之后,如果一个解通过有限次循环不能被进一步改良,则该解被舍弃.设解xi被舍弃,则此解对应的引领蜂将变成一个侦查蜂.此侦查蜂由(13)式产生一个新解代替被舍弃解,

其中lj和uj分别是变量xij的下界和上界.然后返回引领蜂搜索过程,开始重复循环.

利用ABC算法对本文目标函数进行优化的基本步骤如下：

Step 1随机初始化种群xi,i=1,2,···,SN ;

Step 2根据(10)式计算种群中各初始解xi的适应度值fit(xi);

Step 3引领蜂根据(11)式产生新解vi,并计算适应度值fit(vi),根据贪心原则保留适应度值大的解;

Step 4跟随蜂根据(12)式计算解xi的选择概率p(xi),依轮盘赌方式选择解xi;

Step 5若解xi被选择,根据(11)式产生新解vi,并计算新解的适应度; 跟随蜂根据贪心原则保留优质解;

Step 6记录最优解,即适应度值最大的解θbest;

Step 7侦查蜂决定是否存在需要放弃的解,如果存在,根据(13)式随机产生一个新解替代旧解;

Step 8重复Step 3至Step 7,直到达到最大循环次数.

经人工蜂群算法优化后得到最优解θbest=[A,a,b,T,τ,fc,φ],不仅能直接获取超声回波的特征参数,还可以利用特征参数根据(3)式重构出超声信号s(θbest; t).

4 数值仿真与分析

目前,广泛用于评价仿真信号的指标有均方误差(mean square error,MSE)、估计信噪比(estimated signal-to-noise ratio,ESNR)和波形相似系数(normalized correlation coefficient,NCC)[22,23].其定义如下：

4.1 超声信号仿真

为定量分析本文算法对超声回波的处理能力,以(3)式中的超声信号与3个不同强度的脉冲序列进行卷积模拟发射的超声脉冲通过三层不同介质界面后产生的超声回波.无噪超声信号及加入噪声后的超声回波如图3所示.

图3 仿真超声回波 (a)理想超声信号; (b)含噪超声回波Fig.3.Simulation of ultrasonic echo： (a) Ideal ultrasonic signal; (b)noisy ultrasonic echo.

图3(a)所示超声信号由三个回波组成,第一个回波渡越时间为 1 00 μs,第二个回波渡越时间为500 μs,第三个回波渡越时间为 8 00 μs.但在工程实际中获取的观测回波通常是如图3(b)所示的被噪声污染的信号,经计算其信噪比为5 dB,从图中已经无法准确分辨出第二个回波和第三个回波的渡越时间.

4.2 渡越时间估计

利用本文方法可以直接提取超声回波渡越时间.对图3(a)所示的信号,给其加入不同大小的噪声,利用本文方法分别对不同信噪比的超声回波渡越时间进行估计,结果列于表1.

表1 渡越时间τ估计误差Table 1. Estimation error of transit time τ.

从表1可以看出,在相同信噪比下估计第一个回波渡越时间的绝对误差小于第二个回波,第二个回波小于第三个回波.这是因为三个回波幅值依次减小,第三个回波最先被噪声淹没(见图3(b)).随信噪比的降低估计误差增大,但在信噪比低至—10 dB时,估计第一个回波渡越时间的绝对误差为1.3809 μs,相对误差为1.3809%,估计第二个回波渡越时间的绝对误差为 4 .4061 μs,相对误差为0.8812%,由于噪声太大,无法估计出第三个回波的渡越时间.

4.3 降噪性能分析

为了验证本文方法的降噪能力,将本文方法和当前超声信号处理中常用的小波阈值降噪方法、经验模态分解(empirical mode decomposition,EMD)方法进行对比.为方便研究,不失一般性地选择图3(a)中的第二个回波作为研究对象对其加入噪声,当信噪比SNR=—10 dB时,本文方法、小波阈值降噪方法以及EMD方法的降噪结果如图4所示.

如图4所示,当SNR=—10 dB时,超声信号已完全被噪声淹没,采用本文方法可以不失真地重构出超声回波,该回波的渡越时间为 5 00 μs ; 而小波方法和EMD方法已产生失真,从图中无法准确分辨出回波的渡越时间.为进一步验证本文方法对超声信号降噪处理的有效性,对上述原始信号加入不同信噪比的噪声后,分别采用本文方法、小波方法和EMD方法进行降噪处理,得到的降噪性能评价指标结果如表2所列.

图4 不同方法降噪结果对比 (a)黑色实线为含噪回波,红色虚线为原始信号; (b)本文方法处理结果; (c)小波方法处理结果;(d) EMD方法处理结果Fig.4.Comparison of de-noising by different methods： (a) Noisy echo plotted in a black solid line,original signal plotted in a red dotted line; (b) signal de-noised by our proposed method; (c) signal de-noised by wavelet method; (d) signal de-noised by EMD method.

均方误差MSE反映了原始信号和去噪后信号之间的平均误差,波形相似系数NCC衡量原始信号和去噪后信号之间的波形相似度,估计信噪比ESNR反映了降噪后信号的凸显程度.MSE越小,ESNR越大,NCC越接近与1,说明降噪效果越好[23].从表2中可以看出,信噪比较高时本文方法、小波方法和EMD方法都可以取得较好的降噪结果,在信噪比较低时本文方法降噪结果明显优于小波方法和EMD方法.

5 实测信号验证

在实验室中用RITEC公司生产的RITECRAM-5000-SNAP超声检测系统采集了用透射法测试金属构件拉伸疲劳实验的超声回波进行处理,从而验证本文方法是否可以有效处理实测超声信号.如图5所示,超声回波采集装置由RAMSNAP系统、匹配电阻、衰减器、低通滤波器、超声换能器、示波器及计算机等组成.实验中设置RAM-SNAP系统激发频率为5 MHz、周期数为8的脉冲信号,经过匹配电阻和低通滤波器传输到超声换能器上,再由换能器将电信号转换成超声波发射到被测金属试件上表面,超声波在被测试件内部传播到达试件下表面后发生透射,透射的超声波被接收换能器接收并转换为电信号,最后将接收到的电信号在示波器上显示.

表2 不同方法降噪结果对比Table 2. Comparison of denoising results using different methods.

图5 超声回波采集装置Fig.5.Ultrasonic echo acquisition device.

5.1 超声回波模型性能对比

为验证本文提出的双高斯衰减模型的有效性,将目前常用的指数模型、高斯模型及双高斯衰减模型分别和采集到的实测超声回波进行匹配,结果如图6所示.

图6 不同超声信号模型对比Fig.6.Comparison of different ultrasonic signal models.

从图6可以看出,相比指数模型和高斯模型,本文提出的双高斯衰减模型与实测超声回波匹配度更高; 此外,经计算本文模型和实测回波的均方误差 MSE=9.4 × 10—5,波形相似系数 NCC=0.98; 指数模型和实测回波的均方误差MSE=3.5 × 10—4,波形相似系数 NCC=0.93; 高斯模型和实测回波的均方误差 MSE=2.9 × 10—4、波形相似系数NCC=0.94.因此,本文提出的双高斯衰减模型和实验室中采集到的实际超声回波信号更接近.

5.2 实测回波处理

由于和工程实际相比实验室环境中的噪声很小,为验证本文方法对实测超声回波的降噪处理能力,将无激励时测试到的系统噪声放大后叠加到实测回波中以增加实测回波中的噪声.实测超声回波、加入噪声的实测超声回波及本文方法、小波方法和EMD方法处理结果对比如图7所示.

图7 实测超声回波处理结果对比 (a)实测超声回波;(b)混入系统噪声后的实测超声回波,信噪比为—10 dB;(c)本文方法结果; (d)小波方法结果; (e) EMD方法结果Fig.7.Comparison of measured ultrasonic echo processing results： (a) Measured ultrasonic echo; (b) measured ultrasonic echo contaminated by system noise; (c) echo processed by our proposed method; (d) echo processed by wavelet method; (e) echo processed by EMD method.

经计算,小波处理结果的MSE=1.6 × 10—3,NCC=0.62,ESNR=1.64 dB; EMD方法的MSE=2.6 × 10—3,NCC=0.65,ESNR=—0.49 dB;本文方法的 MSE=0.2 × 10—3,NCC=0.97,ESNR=11.58 dB.此外,从图7(c)中可以分辨出本文处理结果回波渡越时间为 7 .29 μs,和图7(a)中结果一致.因此在信噪比低至—10 dB时本文方法仍可以提取出实测含噪超声回波中的超声信号及对应的渡越时间,而此时小波方法和EMD方法已经失真.

6 结 论

基于超声信号有特定结构而噪声和超声信号的结构无关,本文提出了一种超声回波处理方法,并通过对仿真信号和实测回波的处理验证了该方法能有效解决强噪声中的超声回波参数估计和降噪问题且处理效果优于小波方法和经验模态分解方法.本文主要创新如下：

1) 将超声回波的参数估计和降噪问题转换为函数优化问题,通过建立超声信号模型、确定目标函数、优化目标函数和根据优化结果重构无噪信号四步来实现信号参数估计和降噪;

2) 提出了双高斯衰减超声信号模型,相比于目前广泛应用的指数模型和高斯模型,该模型产生的仿真超声信号与实际测量到的超声信号更接近;

3) 基于超声信号有特定结构而噪声和超声信号的结构无关这一事实,在理论上证明了以余弦相似度作为目标函数可以将多回波信号估计问题转化为多个单回波信号估计问题的线性组合,对于单回波估计问题可通过目标函数值是否为1判断优化算法是否正确收敛.

同时还存在一些问题,例如在实验中发现信噪比低于—10 dB时本文方法的稳定性和准确性都变差.这可能是由于采用的优化算法不能正确收敛导致的.开发适合解决本文目标函数的智能优化算法以增强本文方法处理效果将是我们下一步要重点研究解决的问题.