CN103487794A - 一种基于小波包变换的水底混响抑制方法 - Google Patents

Abstract

本发明涉及一种基于小波包变换的水下混响抑制方法，该方法首先对接收到的阵元域信号进行滤波预处理与常规波束形成后，然后，采用小波包变换方法提取波束域数据的时频特征，获得混响与目标回波的差异特征，并据此进行分类识别，再根据特征的分布情况计算接收信号与发射脉冲的整体相似度，最后，对匹配滤波过程进行修正，达到混响抑制的目的；采用该方法可以避免对混响空间分布特性的建模与分析，具有更好的普适性与鲁棒性，应用前景可观。

Description

一种基于小波包变换的水底混响抑制方法

技术领域

本发明属于主动声纳数字信号处理领域，特别涉及一种基于小波包变换的水底混响抑制方法。

背景技术

混响是影响主动声纳检测性能的重要因素之一，它是发射信号通过水中颗粒、水体以及水面的散射叠加而成。其中，来自于水面和水底的界面混响最为严重。混响的产生机制使得它们在时频域上的特性与目标回波非常相似，难以通过常规的能量检测方法进行区分。由于在常规的主动声纳系统中，主要采用匹配滤波进行检测，因此，混响的存在会带来大量的虚警，使得目标难以被发现。

虽然，可以通过合理设计发射波束的指向性以避免其与水面或者水底接触，或是利用自适应波束形成提高接收波束的主瓣分辨率，或是增大基阵孔径以提高角度-距离分辨率从，或是设计具有针对特性的主动发射信号以抑制混响，但是这些前端处理方法的局限性太大，不具备足够的普适性。对于像前视声纳、扫雷声纳这种有针对性地向水面、水底发射信号的主动声纳而言，常规的前端处理方法是不够的。基于这些原因，一些后端处理算法成为近些年来的研究热点。

常用的混响抑制算法主要包含如下几种：空间处理方法，利用自适应波束形成的方法提高接收波束主瓣分辨率，以达到混响抑制的目的，如文献1“J.N.Maksym,M.S.Wunsch.Adaptive beamforming against reverberation for a three-sensor array.J.Acoust.Soc.Am vol 102,pp:3433-3438,1997”；图形学方法，通过提高主动声纳图像中的混响背景与目标回波的对比度来抑制混响，如文献2“G.Ginolhac,J.Chanussot,andC.Hory.Morphological and statistical approaches to improve detection in the presence ofreverberation.IEEE.J.Oceanic engineering,vol.30,pp:881-899,2005”；PCI方法，利用子空间的方法将回波信号分解为混响子空间与回波加白噪声子空间，通过设置一些门限准则去除混响相关的分量，从而抑制混响，如文献3“Guillaume Ginolhac,Genevieve Jourdain.Principal Component Inverse Algorithm for Detection in thePresence of Reverberation.IEEE jouranl of Oceanic Engineering,vol27,No.2,pp:310-321,April 2002”；预白化方法，将原本非平稳、有色的混响噪声进行白化处理，利用修正协方差方法估计白化滤波器的参数，从而抑制混响，如文献4B.W.Choi,E.H.Bae,J.S.Kim and K.K.Lee.Improved prewhitening method for linearfrequency modulationreverberation using dechirping transformation.J.Acoust.Soc.Am,vol 123,pp:EL21-EL25,2008”。

文献5“S.Kay，J.Salisbury.Improved active sonar detection using autoregressiveprewhiteners.J.Acouts.Soc.Amer.,vol.87,pp.1603-1611,Arr.1990.”针对混响的非平稳、有色特征提出了次优的分帧自适应预白化检测器，作者用修正协方差方法估计白化滤波器的参数，用广义似然比检测器（GLRT）检测，提高了主动声纳在混响中的检测能力。

文献6“V. Carmillet,G.Jourdain.Detection of phase-or frequency-modulated signalsin reverberatioin noise.J.Acoust.Soc.Amer.,vol.105,no.6,pp.3375-3389,Jun.1999.”将混响的预白化方法推广到调频调相信号，验证了该算法对于这些主动声纳发射信号的适用性。

文献7“A.Trucco.Experimental results on the detection of embeded objects by aprewhitening filter.IEEE,Joural of Oceanic Engineering,vol 26,No 4,pp.783-794,OC T,2001.”认为混响噪声并不一定服从高斯分布，利用高阶统计的方法估计白化滤波器的参数，取得了比传统的修正协方差方法略好的效果，但是，作者也认为这一结果是由于混响的非高斯统计特性造成的，同时作者也指出高阶统计的方法需要较多的样本才能获得较好的效果。

总的来说，这些方法需要对混响数据帧之间的平稳特性进行分析，并且需要充分考虑水声信道的非线性特性，其建模过程较复杂，而且需求大量的执行资源。实际应用中需要一种普适性更强的方法抑制混响，对随时变化的混响环境具有较强的鲁棒性。

发明内容

本发明的目的在于，为了解决以往混响抑制算法中的建模过程复杂、计算开销过大、受混响分布特性影响较大等局限性，提出了一种使用小波包变换的特征提取算法，获得混响与目标回波的差异特征，并进行分类识别。根据识别结果，对匹配滤波过程进行修正，从而达到混响抑制的目的。

为了实现上述发明目的，本发明提供的一种基于小波包变换的水底混响抑制方法，该方法首先对接收到的阵元域信号进行滤波预处理与常规波束形成后，然后，采用小波包变换方法提取波束域数据的时频特征，获得混响与目标回波的差异特征，并据此进行分类识别，再根据特征的分布情况计算接收信号与发射脉冲的整体相似度，最后，对匹配滤波过程进行修正，达到混响抑制的目的；

所述的水底混响抑制方法的步骤包括：

1）对接收阵所有接收到的阵元域信号进行预处理，得到波束域数据信号；

2）采用小波包变换方法提取波束域数据的时频特征；

3）根据发射信号的中心频率与带宽，选取小波包变换的分解层数与小波包系数组，对所获取的离散小波分解系数进行多维样本矩阵重构；

4）依据重构的多维样本矩阵计算接收信号的小波包系数矩阵和发射脉冲的小波包系数矩阵的相关特征，获得彼此之间的多维相似度特征，并据此计算最终的整体相似度，该相似度的计算公式为：

$\mathrm{Sim}={\left(\underset{i=1}{\overset{\dim }{\mathrm{\Π}}}{(c-{\mathrm{\λ}}_i^2)}^{-2}\right)}^d$

其中，c与d是规整化参数，作用是控制相似度值的变化范围，c的取值需要满足

d为经验参数；dim为重构样本矩阵的维度，与选取的小波包系数组数量保持一致；d的一般取值为5，为经验参数，具体情况下可以适当改变。

5）根据波束域接收信号的相似度特征与整体相似度，对匹配滤波过程进行修正，达到混响抑制的目的。

作为上述技术方案的一种改进，所述的步骤1）的预处理是根据发射脉冲信号的中心频率与带宽，对第K帧与K+1帧的阵元域接收信号进行滤波处理，获得对应带宽内的时域信号，然后，对时域信号进行快速傅里叶变换FFT，选择带宽内对应的谱线数据，最后获得频域波束形成的结果。

作为上述技术方案的一种改进，所述的步骤2）的小波包变换方法是对各个接收波束数据与发射信号进行小波包变换，获得小波包系数组；包括：

21）根据尺度函数φ与小波函数ψ，产生一组能够用于分解与重构信号的函数族；

选取合适的参数j，对于离散的数字信号而言，使得2^j大于等于信号的奈奎斯特采样率；

设接收信号为

那么，将接收信号f_j分解为f_j=ω_j-1+f_j-1，其中：

${\mathrm{\ω}}_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{b}_k^{j-1}\mathrm{\ψ}(2^{j-1}x-k)$

$f_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{a}_k^{j-1}\mathrm{\φ}(2^{j-1}x-k)$

$b_k^{j-1}=\frac{a_{2k}^j-a_{2k+1}^j}{2},$ $a_k^{j-1}=\frac{a_{2k}^j+a_{2k+1}^j}{2},$

通过有限次分解，得到：f_j＝ω_j-1+ω_j-2+...+ω₀+f₀；

式中，ω表示小波函数分量，f表示小波尺度分量；k表示区间[1,2^j]内的整数，j为频率抽选参数，x表示时间轴，φ表示尺度函数，ψ表示小波函数，a表示尺度函数分量参数，b表示小波函数分量系数，Z表示整数集合；

22）对每一组ω分量进行进一步细分，所有的ψ均通过两尺度方程进行进一步分解，其表达式为：

${\mathrm{\ψ}}_{2n}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{h}_k\mathrm{\ψ}(2t-k)$

${\mathrm{\ψ}}_{2n+1}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{g}_k\mathrm{\ψ}(2t-k)$

其中，

{h_k}_k∈Z与{g_k}_k∈Z分别为一组低通滤波器和一组高通滤波器，获得的a_k,b_k,h_k,g_k系数组就是小波包分解对应的小波包系数组。

作为上述技术方案的一种改进，所述的步骤3）根据发射信号的中心频率与带宽选取小波包变换的分解层数，是对于第i层的小波包分解而言，一共可以获得2ⁱ组系数，按照从低频到高频的顺序进行系数组排列，这1~2ⁱ组系数实现了对于整个信号频带的均匀划分。

作为上述技术方案的一种改进，所述的步骤3）根据发射信号的中心频率与带宽选取小波包系数组，如果是分解n层，那么只选取第n层中对应带宽的系数组。

作为上述技术方案的一种改进，所述的步骤4）依据重构的多维样本矩阵计算接收信号的小波包系数矩阵和发射脉冲的小波包系数矩阵的相关特征的矩阵方程为：

C_xx ^-1C_xyC_yy ^-1C_yxw_x=λ²w_x

C_yy ^-1C_yxC_xx ^-1C_xyw_y=λ²w_y

其中，λ_i为样本矩阵第i维数据对应的相关系数，C_xx,C_xy,C_yx,C_yy分别表示多维样本矩阵的自协方差与互协方差矩阵。

作为上述技术方案的一种改进，所述的步骤5）的对匹配滤波过程进行修正，是匹配滤波以帧长作为步长，将常规匹配滤波的结果与相似度系数值相乘，获得最终的输出结果。

作为上述技术方案的一种改进，所述的接收阵为半圆弧接收阵。

本发明根据发射阵与接收阵的参数选择合适的发射脉冲信号，并持续接收回声信号。对所有接收到的阵元域信号进行滤波预处理与常规波束形成，采用小波包变换方法提取波束域数据的时频特征，并据此识别混响与回波信号。根据特征的分布情况计算接收信号与发射脉冲的相似度，并对匹配滤波过程进行修正，最终达到混响抑制的目的。

本发明的优点是：

1)采用离散小波变换提取目标回波与混响的差异特征，利用了两者在时频域中能量分布的差异特性，不用对混响进行建模，也不用充分考虑混响的时空分布特性，具有良好的鲁棒性与普适性。

2)可以根据发射信号的带宽特性灵活选择离散小波变换参数，获得不同时频分辨率下的特征，适应不同类型目标的探测与识别任务。

3)计算过程较为简洁，容易实现实时处理。

本发明利用小波包变换灵活的时频分析手段，获得混响与目标回波的能量在时频空间中的分布特性差异，并据此计算它们与发射脉冲之间的结构相似度，根据获得的结构相似度，对匹配滤波过程进行修正，实现混响与目标回波的区分，达到抑制混响的目的。

采用该方法可以避免对混响空间分布特性的建模与分析，具有更好的普适性与鲁棒性，应用前景可观。

附图说明

图1是本发明中采用的小波包变换对于数字信号时频结构的划分示意图；

图2是本发明中具体实施的湖试条件说明；

图3是本发明基于小波包变换的水底混响抑制方法的流程示意图；

图4a是当目标位于弱混响区域时，采用本发明算法处理前的结果示意图；

图4b是当目标位于弱混响区域时，采用本发明算法处理后的结果示意图；

图5a是当目标位于强混响区域时，采用本发明算法处理前的结果示意图；

图5b是当目标位于强混响区域时，采用本发明算法处理后的结果示意图。

具体实施方式

下面结合附图对本发明进行进一步说明。

如图1-3所示，本发明基于小波包变换的水底混响抑制方法的具体步骤如下：

1)用线阵发射调频脉冲信号，并用半圆弧阵接收回声信号。根据发射脉冲信号的中心频率与带宽，对第K帧与K+1帧的阵元域接收信号进行滤波处理，获得对应带宽内的时域信号。然后对时域信号进行快速傅里叶变换(FFT)，选择带宽内对应的谱线数据，获得频域波束形成的结果。

2)对各个接收波束数据与发射信号进行小波包变换，获得小波包系数组。选择合适的尺度函数φ与小波函数ψ，产生一组能够用于分解与重构信号的函数族。选取合适的参数j，使得2^j大于信号的奈奎斯特采样率。对于离散的数字信号而言，只要使得2^j大于等于信号的采样率即可。假设接收信号表达式为

那么f_j可以分解为f_j=ω_j-1+f_j-1，其中：

${\mathrm{\ω}}_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{b}_k^{j-1}\mathrm{\ψ}(2^{j-1}x-k)$ (1)

$f_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{a}_k^{j-1}\mathrm{\φ}(2^{j-1}x-k)$

$b_k^{j-1}=\frac{a_{2k}^j-a_{2k+1}^j}{2},$ $a_k^{j-1}=\frac{a_{2k}^j-a_{2k+1}^j}{2}---\left(2\right)$

通过有限次分解，可以得到：f_j=ω_j-1+ω_j-2+...+ω₀+f₀。ω表示小波函数分量，f表示小波尺度分量。

对每一组ω分量进行进一步细分。所有的ψ均可以通过两尺度方程进行进一步分解，其表达式如下：

${\mathrm{\ψ}}_{2n}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{h}_k\mathrm{\ψ}(2t-k)$ (3)

${\mathrm{\ψ}}_{2n+1}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{g}_k\mathrm{\ψ}(2t-k)$

其中

{h_k}_k∈Z与{g_k}_k∈Z可以分别看成一组低高通滤波器。获得的a_k,b_k,h_k,g_k系数组就是小波包分解对应的系数组。

3)根据发射信号的中心频率与带宽，决定小波包变换的分解层数与系数组选取方式。小波包分解对于信号频带的划分，如图1所示：接收信号的采样率为fs，那么对于第i层的小波包分解而言，一共可以获得2ⁱ组系数。按照从低频到高频的顺序进行系数组排列，那么这1~2ⁱ组系数实现了对于整个频带的均匀划分。例如对于第1层小波包系数而言，它们对应的频带分别为0~fs/4与fs/4~fs/2；而对第2层小波包系数而言，它们对应的频带分别为0~fs/8、fs/8~fs/4、fs/4~3fs/8与3fs/8~fs/2。小波包分解层数需要选择合适的值，要求发射信号的带宽尽可能严格对应到其中几组小波包系数。不采用多层的小波包系数组，如果是分解n层，那么只选取第n层中对应带宽的系数组。

4)对所获取的离散小波分解系数进行多维样本矩阵重构，并计算它们与发射脉冲对应的离散小波分解系数矩阵的相似度。要求相似度计算法则既能保证混响与目标回波之间具有显著差异，同时要求不同类型目标回波之间的相似度比值不能过大。相似度计算公式如下：

$\mathrm{Sim}={\left(\underset{i=1}{\overset{\dim }{\mathrm{\Π}}}{(c-{\mathrm{\λ}}_i^2)}^{-2}\right)}^d---\left(4\right)$

其中c与d是规整化参数，作用是控制相似度值的变化范围；dim为重构样本矩阵的维度，与选取的小波包系数组数量保持一致；

5)根据波束域接收信号的特征提取与相似度计算结果，对匹配滤波过程进行修正，从而使混响数据的相关峰不再明显。λ_i为样本矩阵第i维数据对应的相关系数，其值可以通过如下矩阵方程求解：

C_xx ^-1C_xyC_yy ^-1C_yxw_x=λ²w_x    (5)

C_yy ^-1C_yxC_xx ^-1C_xyw_y=λ²w_y

λ_i即为式(5)中的λ第i个根。其中C_xx,C_xy,C_yx,C_yy分别表示多维样本矩阵的自协方差与互协方差矩阵。

根据发射信号的脉宽与接收信号的单帧长度，确定联合处理的帧数。例如采用8ms的发射信号，接收帧长为1ms，那么需要以1ms为步长，每次处理8帧数据。根据计算获得的相似度系数值，对匹配滤波过程进行修正。匹配滤波也以帧长作为步长，将常规匹配滤波的结果与相似度系数值相乘，获得最终的输出结果。

图2为本发明一实施例湖试中主动声纳基阵与目标的吊放情况，以及水下环境的示意图。其中基阵与目标的吊放深度均为15m，混响源主要是水面与水下山体，目标采用一个小型球状金属材质目标。目标的方位可以任意移动，但是吊放深度保持15m不变。采用HFM信号作为发射脉冲，其中心频率为30kHz，带宽为8kHz。

图3为接收信号的整体处理流程示意图。其中小波包变换中的分解层数以及系数组的选择由发射信号的带宽与中心频率决定，且所选的系数组数决定了多维样本矩阵的维度。接收信号的小波包系数矩阵需要与发射脉冲的小波包系数矩阵进行相关特征计算，获得彼此之间的多维相似度特征，并据此计算最终的整体相似度。整体相似度将被用于匹配滤波过程的修正。

图4a和图4b为球目标位于弱混响区时的处理前后结果对比示意图。可以看到，直接采用匹配滤波后，混响的强度要高于目标回波，因此将带来大量的虚警，严重影响接收机的检测性能。采用本发明中的算法处理后，混响得到了极大的抑制，通过设置合适的判决门限可以消除其带来的虚警，信混比提升了约13.0dB。

图5a和图5b为球目标位于强混响区时的处理前后结果对比示意图。可以看到，直接采用匹配滤波后，混响的强度要远高于目标回波，且无法通过设置门限获得目标的准确方位与距离。采用本发明中的算法处理后，不仅有效降低了混响，而且球目标的检测结果也相当理想，信混比提升了约12.3dB。

总而言之，本发明能够极大地提升信混比，在不影响目标回波检测性能的前提下，达到混响抑制的目的。即使目标淹没在强混响区域内，也能够获得理想的检测结果。

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (8)

1.一种基于小波包变换的水底混响抑制方法，该方法首先对接收到的阵元域信号进行滤波预处理与常规波束形成后，然后，采用小波包变换方法提取波束域数据的时频特征，获得混响与目标回波的差异特征，并据此进行分类识别，再根据特征的分布情况计算接收信号与发射脉冲的整体相似度，最后，对匹配滤波过程进行修正，达到混响抑制的目的；所述的水底混响抑制方法的步骤包括：

1）对接收阵所有接收到的阵元域信号进行预处理，得到波束域数据信号；

2）采用小波包变换方法提取波束域数据的时频特征；

3）根据发射信号的中心频率与带宽，选取小波包变换的分解层数与小波包系数组，对所获取的离散小波分解系数进行多维样本矩阵重构；

4）依据重构的多维样本矩阵计算接收信号的小波包系数矩阵和发射脉冲的小波包系数矩阵的相关特征，获得彼此之间的多维相似度特征，并据此计算最终的整体相似度，该相似度的计算公式为：

$\mathrm{Sim}={\left(\underset{i=1}{\overset{\dim }{\mathrm{\Π}}}{(c-{\mathrm{\λ}}_i^2)}^{-2}\right)}^d$

其中，c与d是规整化参数，作用是控制相似度值的变化范围，c的取值需要满足

d为经验参数；dim为重构样本矩阵的维度，与选取的小波包系数组数量保持一致；

5）根据波束域接收信号的相似度特征与整体相似度，对匹配滤波过程进行修正，达到混响抑制的目的。

2.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤1）的预处理是根据发射脉冲信号的中心频率与带宽，对第K帧与K+1帧的阵元域接收信号进行滤波处理，获得对应带宽内的时域信号，然后，对时域信号进行快速傅里叶变换FFT，选择带宽内对应的谱线数据，最后获得频域波束形成的结果。

3.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤2）的小波包变换方法是对各个接收波束数据与发射信号进行小波包变换，获得小波包系数组；包括：

21）根据尺度函数φ与小波函数ψ，产生一组能够用于分解与重构信号的函数族；

选取合适的参数j，对于离散的数字信号而言，使得2^j大于等于信号的奈奎斯特采样率；

设接收信号为

那么，将接收信号f_j分解为f_j=ω_j-1+f_j-1，其中：

${\mathrm{\ω}}_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{b}_k^{j-1}\mathrm{\ψ}(2^{j-1}x-k)$

$f_{j-1}=\underset{k\∈Z}{\mathrm{\Σ}}{a}_k^{j-1}\mathrm{\φ}(2^{j-1}x-k)$

$b_k^{j-1}=\frac{a_{2k}^j-a_{2k+1}^j}{2},$ $a_k^{j-1}=\frac{a_{2k}^j+a_{2k+1}^j}{2}$

通过有限次分解，得到：f_j＝ω_j-1+ω_j-2+...+ω₀+f₀；

式中，ω表示小波函数分量，f表示小波尺度分量；k表示区间[1,2^j]内的整数，j为频率抽选参数，x表示时间轴，φ表示尺度函数，ψ表示小波函数，a表示尺度函数分量参数，b表示小波函数分量系数，Z表示整数集合；

22）对每一组ω分量进行进一步细分，所有的ψ均通过两尺度方程进行进一步分解，其表达式为：

${\mathrm{\ψ}}_{2n}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{h}_k\mathrm{\ψ}(2t-k)$

${\mathrm{\ψ}}_{2n+1}\left(t\right)=\sqrt{2}\underset{k}{\mathrm{\Σ}}{g}_k\mathrm{\ψ}(2t-k)$

其中，{h_k}_k∈Z与{g_k}_k∈Z分别为一组低通滤波器和一组高通滤波器，获得的a_k,b_k,h_k,g_k系数组就是小波包分解对应的小波包系数组。

4.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤3）根据发射信号的中心频率与带宽选取小波包变换的分解层数，是对于第i层的小波包分解而言，一共可以获得2ⁱ组系数，按照从低频到高频的顺序进行系数组排列，这1~2ⁱ组系数实现了对于整个信号频带的均匀划分。

5.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤3）根据发射信号的中心频率与带宽选取小波包系数组，如果是分解n层，那么只选取第n层中对应带宽的系数组。

6.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤4）依据重构的多维样本矩阵计算接收信号的小波包系数矩阵和发射脉冲的小波包系数矩阵的相关特征的矩阵方程为：

C_xx ^-1C_xyC_yy ^-1C_yxw_x=λ²w_x

C_yy ^-1C_yxC_xx ^-1C_xyw_y=λ²w_y

其中，λ_i为样本矩阵第i维数据对应的相关系数，C_xx,C_xy,C_yx,C_yy分别表示多维样本矩阵的自协方差与互协方差矩阵。

7.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的步骤5）的对匹配滤波过程进行修正，是匹配滤波以帧长作为步长，将常规匹配滤波的结果与相似度系数值相乘，获得最终的输出结果。

8.根据权利要求1所述的基于小波包变换的水底混响抑制方法，其特征在于，所述的接收阵为半圆弧接收阵。

Images