基于循环维纳滤波器的声学成像优化方法

何 通，徐赋民，任豪杰，张二亮

（郑州大学机械与动力工程学院，河南 郑州 450001）

1 引言

循环平稳源是旋转机械工作时的主要噪声源，来自齿轮、轴承等关键部件，其辐射的振动与噪声信号表现为循环平稳性［1］。循环平稳性为信号的处理与分析提供一个强大的工具，已有学者将其用于轴承故障诊断［2］、疲劳损伤检测［3］等领域，引起了广泛关注。旋转机械的噪声信号类型通常很复杂，包含各种平稳源和循环平稳源所辐射信号的混合。考虑到循环平稳源实际上蕴含了大多数重要信息的情况下，对其辐射的循环平稳信号进行盲提取并应用声学成像技术，可以准确地识别、定位和可视化重建主要噪声源，进而有针对性的指导旋转机械的设计及改进，从而减振降噪。

盲分离混合源是通信、生物医学等领域的一类无监督学习问题。盲源分离通常依赖于信号相互统计独立性的假设。相互统计独立性假设对许多复杂机械系统可能不再满足。考虑到混合源中的重要信息往往由一种信号携带，其余的则是干扰或噪声。因此，只需恢复一个感兴趣的信号，盲源分离便简化为盲信号提取。基于循环平稳性的盲信号提取更适合被应用于机械系统，尤其在旋转或往复机械研究方面得到广泛关注。循环平稳信号处理方法最早在文献［4］中被提及，被成功应用于通信领域。后来，它被很多学者引入到机械领域，相继出现了多种盲信号提取方法，例如多元循环回归［5］，子空间盲提取［6］，降秩循环回归［7］和角时循环维纳滤波器［8］等。

盲信号提取方法目前大多被用于机械的振动信号分析，而非声学信号处理。但振动信号需要采用加速度计等进行传统的表面接触式测量，这样会造成两个弊端，一是由于加速度计本身的重量可能使其附着的机械结构固有特性发生改变，二是在高温或腐蚀环境等情况下不能适用。随着基于麦克风阵列测量的声学成像技术［9］的快速发展，声学信号处理取得很大进展。近场声全息作为其中一种重要的噪声源识别、定位及可视化重建技术，不仅具有非接触式测量、空间分辨率高的优点，还包含丰富的声场信息，可以为噪声控制提供指导依据。贝叶斯近场声全息是基于贝叶斯理论对近场声全息方法的优化，与传统的近场声全息相比具有更好的重建精度，最早在文献［10］中被提出。然而，这些声学成像技术仅能用于声源重建，而不能用于声源分离。目前而言，尽管其可能的应用价值巨大，但有关方面的研究却极少，仅在文献［11-12］中尝试去分离不同噪声源。

综上所述，针对噪声环境中循环平稳源的重建问题，提出了一种基于循环维纳滤波器的声学成像优化方法。首先对麦克风阵列测量的信号进行滤波，盲提取出由感兴趣的循环平稳源辐射的循环平稳信号，然后引入贝叶斯近场声全息重建循环平稳源。实验结果表明，所提方法能够很好的消除其它噪声源的干扰，仅重建感兴趣的循环平稳源。

2 问题描述

2.1 盲信号提取问题

盲信号提取问题是从观测信号xi(n)中提取循环平稳信号xci(n)，可以写成：

式中：i=1，…，M，M—传感器个数；pi(n)—一阶循环平稳信号；xci(n)—二阶循环平稳信号；ni(n)—“噪声”部分，包含所有与pi(n)和xci(n)均不相关的干扰或噪声。

2.2 循环平稳信号及其二阶描述

循环平稳信号指的是信号本身不具有周期性，但信号统计量（如均值、自相关函数等）呈周期性变化的信号。均值是周期性的信号为一阶循环平稳信号（CS1），更一般的循环平稳信号是二阶循环平稳信号（CS2），即其自相关函数是周期性的。在此，主要研究CS2的盲提取算法，因为CS1通过同步平均法［1］不难得到。

CS2信号xci(n)的时变自相关函数为

式中：n—时间；τ—时滞；αk≠0，k=1，…，K—循环频率。循环谱密度由(τ)的傅里叶变换得到，即：

式中：f—谱频率。

根据CS2信号的二阶循环平稳性，其循环谱密度(f) ≠0，即以循环频率αk∈A隔开的两个频谱分量之间存在非零相关性。

3 理论研究

3.1 循环维纳滤波器

在信号平稳的情况下，盲提取问题的解决方案是利用线性时不变滤波器，众所周知，它是维纳滤波器。有趣的是，在信号循环平稳的情况下，盲提取问题的最优解具有线性周期时变结构，因此，可以将其命名为循环维纳滤波器。本节主要基于线性周期时变结构介绍了盲提取CS2信号的工作原理，以及循环维纳滤波器的计算问题。

3.1.1 盲提取CS2信号的工作原理

对于循环频率αk∈A，线性周期时变滤波器可写成以下形式

式中：gi，k(n)—线性时不变滤波器。输入-输出关系可写为：

式中：*—卷积算符；zi，k(n)=xi(n)⋅e-2πjαkn—输入信号xi(n)的频移，αk可通过计算频谱相关性估算得到［13］。

可以看出，对信号的线性时变滤波等效于对信号的频移进行线性时不变滤波的总和。基于线性周期时变原理盲提取CS2信号的结构框图，如图1所示。

图1 盲提取CS2信号的结构框图Fig.1 Block Diagram of Blind Extraction of CS2 Signal

3.1.2 循环维纳滤波器的计算

如上所述，已经显示了如何基于线性周期时变原理提取CS2信号，现在剩下的事情就是找到滤波器gi，k(n)，或等效地求出其频域形式Gi，k(f)。

由节2.2可知，盲提取CS2信号正是基于两个频谱分量之间的非零相关性，因此，一种更相关的方法是从频域上进行计算。那么，式（5）的频域形式为：

循环维纳滤波器是在最小均方误差意义上最佳地提取CS2分量，即：

由此，求得其矩阵形式的最优滤波器为：

最终，将最优滤波器代入式（6）可求得CS2分量(f)，再进行傅里叶逆变换得到(n)。遍历麦克风阵列上所有声压传感器，便计算出

3.2 贝叶斯近场声全息

如上所述，经过循环维纳滤波器进行滤波后，问题就变成：已知，重建感兴趣的循环平稳源s1。为此，采用贝叶斯近场声全息方法。

设声源为s1(r)，r∈Γ是感兴趣的源面区域，声源到麦克风阵列的正向传播模型为：

式中：G(r)—传递矩阵；v—测量噪声。声源分布可以写成以下形式：

式中：M—空间基函数ϕm(r)的数量（等于麦克风阵列上传感器的个数）；cm—系数。

因此，估计声源的问题就归结为寻找最优空间基函数Φ及其系数c的问题。可以通过后验概率分布最大化，获得最优参数和，这就是所谓的最大后验估计（［］表示随机变量的概率密度函数）。

基于贝叶斯公式，后验概率分布可表示为：

式中：[1(c，Φ)]—似然函数，表示给定源分布s1(c，Φ)下观察到测量数据xc的概率；s1(c，Φ)—先验概率分布，表示源分布的空间先验信息，可通过引入孔径函数来表示，满足

最终，重建感兴趣的循环平稳源为：

3.3 声学成像优化方法

将提出的基于循环维纳滤波器的声学成像优化方法进行梳理总结，具体流程，如图2所示。

图2 声学成像优化方法流程图Fig.2 Flow Chart of Acoustic Imaging Optimization Method

4 盲信号提取算法评估

在声学成像优化方法中，循环维纳滤波器的性能对声源重建结果影响较大，而许多因素都会影响其降噪质量，因此需要对盲信号提取算法进行评估，探究各参数对盲提取精度的影响。

4.1 仿真信号模型

根据式（1）对含噪声ni(n)的CS2信号xci(n)进行建模，一阶成分pi(n)已由同步平均法滤除。循环周期由T表示，且长度为L的合成信号模型为：

式中：δ(n)—离散狄拉克脉冲；Al—遵循高斯分布N(μA，σ)的随机振幅；τl—遵循高斯分布N(0，σ)的一系列独立且均匀分布的随机变量；ni(n)—高斯白噪声。信号模型式（14）可以很好地对许多旋转机械信号进行建模，因此将该信号用于盲信号提取算法的仿真验证。

4.2 算法评估

为了评估循环维纳滤波器的降噪效果，可通过计算提取的CS2 信号(n)和真实的CS2 信号xci(n)之间的均方误差来实现，即：

选择采样频率fs=8192Hz，数据长度L=219。将信号按循环频率αk=k/T，k=1，…，K成比例的频移，由于二阶循环平稳信号xci(n)与其频移形式相关，而噪声ni(n)与其频移形式不相关。因此，在对频移的信号进行滤波后，噪声减小。

噪声信号比和用于频移的循环频率数量的影响，如图3 所示。可以看出，增加循环频率的数量K会增强降噪效果。但随着循环频率数量的增加，计算量将变得非常大，这将降低算法的计算效率。因此，算法中使用的循环频率数量应在计算速度和降噪效果之间进行折衷。另外，噪声信号比越低，盲提取精度越高。

图3 噪声信号比和循环频率数量的影响Fig.3 The Influence of the Noise-to-Signal Ratio and the Number of Cycle Frequencies

采样频率和数据长度与上述保持一致，选择循环频率个数K=4，考虑噪声信号比为-5dB的情形。滤波前后的结果和残差，如图4（a）、图4（b）所示。其中，图4（a）中，部分颜色是使用循环维纳滤波器之前的测量信号，另一部分颜色是使用循环维纳滤波器之后盲提取的CS2 信号。图4（b）中，残差图符合平稳白噪声的形式。因此，这验证了盲提取CS2信号的有效性。

图4 CS2信号的盲提取Fig.4 Blind Extraction of CS2 Signal

5 声源重建实验

为了进一步验证所提的声学成像优化方法的有效性，进行了双音响声源重建实验，该实验在混响室中进行。实验装置如图5所示，主要由三部分组成：（1）双音响声源（按一定距离放置的两个白色音响，音响1为循环平稳信号源，由CS2信号驱动，音响2为干扰源，由白噪声信号驱动）；（2）麦克风阵列（由半径为0.5m的圆形阵列支架及64个G.R.A.S.麦克风传感器组成），用于采集由双音响声源所辐射的声压信号；（3）PXIe-1085信号采集器，用于麦克风阵列上的声压传感器所采集信号的接收和存储。

图5 声源重建实验Fig.5 Sound Source Reconstruction Experiment

声源重建实验的具体流程如下：

（1）实验场景布置：将麦克风传感器与信号采集器对应相连，且阵列面与声源面（xy平面）的距离为z=0.5m，音响1在xy平面上距离麦克风阵列中心的位置为（-0.1，0.05）m，音响2在xy平面上距离麦克风阵列中心的位置为（0.25，0.05）m。

（2）信号的产生：利用Matlab分别生成循环平稳信号xci(n)=和白噪声信号ni(n)，并将其转换成可播放格式，从而作为双音响的信号源。

（3）采集声压数据：同时播放音响1和音响2，利用麦克风传感器对其辐射信号的声压进行采集，采样频率fs=8192Hz，采样时间为40s，然后将数据存储到信号采集器中。

（4）数据处理：分别采用贝叶斯近场声全息方法和这里所提的声学成像优化方法处理数据，进行声源重建。

其中，在2422Hz处声源的重建结果，如图6所示。图6（a）是直接采用贝叶斯近场声全息方法得到的重建结果，重建声源的平面位置中心分别为（-0.11，0.06）m和（0.27，0.04）m，与两个音响的实际位置非常接近。由结果可知，贝叶斯近场声全息方法虽然可以较好的实现声源的定位，但并不能抑制其它干扰源。图6（b）则是用基于循环维纳滤波器的声学成像优化方法得到的结果，重建声源的平面位置中心为（-0.1，0.06）m，与音响1的实际位置非常接近。由结果可知，声学成像优化方法既能较好的实现声源的定位，又可以从含有其它干扰源的噪声环境中仅重建感兴趣的循环平稳源。

图6 2422Hz处声源重建结果Fig.6 Sound Source Reconstruction Result at 2422Hz

6 结论

基于循环维纳滤波器提出了一种声学成像优化方法，通过实验验证了所提方法的有效性和可行性。与传统的振动信号分析方法相比，该方法具有非接触式测量和可视化重建的优点。与传统的声学成像技术相比，又解决了在含有其它干扰源的噪声环境中准确地定位和重建感兴趣的循环平稳源的问题。因此，该方法具有一定的应用价值，将为机械系统，尤其是旋转机械系统中关键部件的噪声控制和结构优化等方面提供指导依据。