CN109884592A - 一种面向低频高斯噪声源的声源定位仿真方法 - Google Patents

Abstract

本发明公开了一种面向低频高斯噪声源的声源定位仿真方法，包括：将声源所在平面作为声源扫描平面，将声源扫描平面划分成若干个网格；设置麦克风阵列采集声信号，输出麦克风阵列信号作为阵列输出信号，对阵列输出信号进行频谱分析以及时频谱分析，得到声信号的特征频率；对阵列输出信号表达式的等号两边的式子分别求协方差矩阵得到等号左侧协方差矩阵和等号右侧协方差矩阵；将等号左侧协方差矩阵按行展开，并转置成列向量形式，将得到的列向量作为压缩感知的测量值；将等号右侧协方差矩阵按行展开，得到压缩感知的矩阵表达式；通过测量值和矩阵表达式计算声源的能量向量，定位声源位置。上述方法具有高分辨率和极好的抗噪性。

Description

一种面向低频高斯噪声源的声源定位仿真方法

技术领域

本发明涉及信号处理领域，具体涉及一种面向低频高斯噪声源的声源定位仿真方法。

背景技术

声源定位，就是用非接触、遥感方式表征声音的位置和能量信息，找到异常噪声源的位置，以可视化图像的形式给出声音的空间分布特征。

目前高分辨率声源定位方法主要是谱估计类算法，例如MUSIC算法和ESPRIT算法。

但是，现有技术中至少存在以下缺点和不足：

MUSIC算法和ESPRIT算法都是基于特征值分解的算法，所以必须预知信源数目才能准确估计声源的位置。但是在实际应用中，根本无法准确预估声源的数目，再考虑到声音传播过程中的折射，反射现象造成的混响效应，对于该算法的性能也会造成极大的影响。除此之外，MUSIC算法还无法对声源的绝对能量信息进行准确的估计，只能预估能量的相对值大小。

压缩感知技术是指一种寻找欠定线性系统的稀疏解的技术。它可以突破奈奎斯特采样定理的极限，用较少的测量值复原出原始信号，从而应用于获取和重构稀疏信号。将压缩感知技术用在阵列信号处理领域，通过稀疏波束图整形的方法限制波束图中阵列增益较大的元素个数，同时鼓励较大的阵列增益集中在波束主瓣中，从而达到降低旁瓣水平同时，提高主瓣中阵列增益水平，最终有效提高声源定位算法的分辨率。

现有的声源定位算法大多存在分辨率低，抗噪性差的问题，尤其在低信噪比条件下无法达到较好的分辨率。因此，亟需一种能够有效提高定位的分辨率，同时极大的提高算法的抗噪性，在低信噪比、低频条件下也能够得到高分辨率的声源定位方法。

发明内容

针对本领域存在的不足之处，本发明提供了一种面向低频高斯噪声源的声源定位仿真方法，引入了压缩传感技术，能够有效提高定位的分辨率，同时极大的提高了算法的抗噪性，在低信噪比、低频条件下也能够得到高分辨率。

为了提高声源定位算法的分辨率和稳定性，本发明提出将压缩感知技术应用在阵列信号处理中，极大的提高了波束图的主瓣增益，从而获得高精度和高分辨率。

一种面向低频高斯噪声源的声源定位仿真方法，包括：

(1)将声源所在平面作为声源扫描平面，将声源扫描平面划分成若干个网格；

(2)设置麦克风阵列采集声信号，输出麦克风阵列信号作为阵列输出信号，对阵列输出信号进行频谱分析以及时频谱分析，得到声信号的特征频率；

(3)对阵列输出信号表达式的等号两边的式子分别求协方差矩阵得到等号左侧协方差矩阵和等号右侧协方差矩阵，所述的阵列输出信号表达式为：

X＝AS+N

其中，X为阵列输出信号矩阵，A为麦克风阵列的方向矩阵，S为声信号矩阵，N为噪声矩阵；

(4)将等号左侧协方差矩阵按行展开，并转置成列向量形式，将得到的列向量作为压缩感知的测量值；

(5)将等号右侧协方差矩阵按行展开，得到压缩感知的矩阵表达式；

(6)通过测量值和矩阵表达式计算声源的能量向量，定位声源位置。

步骤(1)中，所述的声源扫描平面指的是一个假想平面，即在空间中约定的一个平面区域，将声源扫描平面划分为若干个网格，每个网格都是一个声源可能存在的位置点，通过逐网格扫描的方式遍历每一个网格进行计算，找出声源所在位置。

步骤(2)中，优选地，所述的麦克风阵列所在平面和声源扫描平面同轴设置，方便标定声源的位置坐标。

优选地，所述的麦克风阵列中相邻麦克风阵元之间的距离不大于待测声源的声信号的半波长。

所述的频谱分析可采用快速傅里叶变换(FFT)。

所述的快速傅里叶变换是离散傅里叶变换的快速运算方法，能够极大程度的提高运算速度。本发明可调用MATLAB中的FFT算法对原始声信号进行处理。

所述的时频谱分析可采用短时傅里叶变换(STFT)。

所述的短时傅里叶变换解决了快速傅里叶变换没有时间分辨率的问题，通常用来分析非平稳信号。该算法的原理是对一段信号进行加窗，将其分割成多段，每一段内分别进行傅里叶变换，得到该段内的局部频谱，最后按照时间顺序把各段拼接起来，得到信号的时频谱。

本发明可调用MATLAB中的STFT算法对原始声信号进行处理。

步骤(6)的具体步骤为：

(6-1)分别用L1范数和L2范数对所述测量值和矩阵表达式进行最小化约束；

(6-2)调用MATLAB的凸优化(CVX)工具包计算声源的能量向量，定位声源位置。

优选地，在步骤(6)之后，可将定位结果以图像化的形式展示在声源扫描平面上，非常直观。

本发明的面向低频高斯噪声源的声源定位仿真方法可用于大多数机械设备的声源定位中，例如汽车，石化管道等，给出机械设备的噪声分布空间特征，最终服务于机械设备的状态检测与故障诊断领域。

本发明与现有技术相比，主要优点包括：本发明提出的面向低频高斯噪声源的声源定位仿真方法能够实现对声源进行有效的定位，尤其是具有极好的抗噪性，在低信噪比条件下依然保持较高的分辨率，同时该算法能够较准确的复原声源的能量信息。除此之外，在中高频率段该算法展现出了优良的性能，在低频段该算法也依旧保持着较好的精度。

附图说明

图1为实施例的面向低频高斯噪声源的声源定位仿真方法的流程示意图；

图2为实施例的MATLAB仿真声源示意云图；

图3为实施例的导入的麦克风阵列示意图；

图4为实施例的声源信号中心频率为300Hz情况下的定位结果云图；

图5为实施例的声源信号中心频率为500Hz情况下的定位结果云图；

图6为实施例的声源信号中心频率为1000Hz情况下的定位结果云图。

具体实施方式

下面结合附图及具体实施例，进一步阐述本发明。应理解，这些实施例仅用于说明本发明而不用于限制本发明的范围。下列实施例中未注明具体条件的实验方法，通常按照常规条件，或按照制造厂商所建议的条件。

如图1所示，本实施例的面向低频高斯噪声源的声源定位仿真方法，包括以下步骤：

S01，用麦克风阵列采集声信号。

S02，对阵列输出信号表达式的等号两边分别求协方差矩阵，具体步骤如下所示：

阵列信号的传播模型如图2所示，假设有M个声源，N个麦克风，则第i个麦克风的输出信号可以表示为：

其中：

S_k(t)表示声源信号；

n_i(t)表示干扰噪声；

a_ik表示相位延迟，具体表示如下：

其中：

表示第k个声源到第i个阵元传播过程中的相位信息；

Δr_i(k)表示第k个声源到第i个阵元的间距，代表幅值衰减；

f为对声信号进行频谱分析和时频谱分析得到的特征频率，也可称为声源的信号中心频率；

c表示声音传播速度；

j为虚数单位，j²＝-1；

e为自然对数的底数。

用矩阵形式可以表示为：

X＝AS+N (3)

其中：

X＝[x₁(t),x₂(t),...,x_N(t)]^T；将X称为麦克风阵列输出信号；

将A称为麦克风阵列的方向矩阵；

S＝[S₁(t),S₂(t),…,S_M(t)]^T；将S称为声源信号矩阵或声信号矩阵；

N＝[n₁(t),n₂(t),...,n_N(t)]^T；将N称为噪声向量或噪声矩阵。

对公式(3)两边同时求协方差矩阵，等式左侧协方差矩阵可以表示为：

其中，X^H表示X的共轭转置。

将矩阵R逐行抽取转化为列向量：

由公式(3)可以得到：

E[XX^H]＝E[ASS^HA^H] (6)。

将各项代入公式(6)右侧表达式，按照矩阵形式运算得到：

同样的，将该矩阵按行展开，写成矩阵相乘的形式，可以得到下式：

S03，将公式(5)中的R_V代入公式(8)中，写成如下形式：

R_V＝GP+Q (9)

其中：

调用凸优化(CVX)工具包求解该问题：

S04，求解出声源的能量和位置，通过在扫描平面上逐点扫描的方式以图像化的形式展示出来。

为突显本发明方法的优越性，将借助MATLAB软件对本发明提出的方法进行仿真，具体步骤包括：

(1)如图2所示，设定1.4×1.4m²的正方形网格扫描平面为假想的声源所在平面。其中，每个小网格为0.05×0.05m²(未在图中显示)。

如图2所示，随机设置了6个仿真声源，初始坐标分别为(-0.3,0)，(0.3,0.2)，(0,0.4)，(0.6,-0.3)，(-0.5,-0.3)，(0.2,-0.4)；右侧的色标显示了声源的能量大小，左侧和下侧的柱状图表示的是在该坐标轴方向上的所有声源的能量的叠加。

(2)导入麦克风阵列，采集声信号，并进行频谱分析和时频谱分析，得到声信号的特征频率。

仿真的麦克风阵列的示意图如图3所示，为7×5的阵列，阵列面积2.4×1.6m²，阵列平面和声源所在平面平行，距离为1m。设置声信号的信噪比为0dB，信号和噪声均采用高斯分布。

(3)分别测试信号中心频率为300Hz，500Hz，1000Hz情况下的定位结果，结果如图4，5，6所示。

为了进一步定量表示定位结果的位置精度和能量精度，引出相对误差和相关系数两个参量。

相对误差主要用来表征定位结果的能量精确度，公式如下：

其中：

x_t表示真值，即仿真声源各个网格点的幅值；

x_r表示重构值，即仿真结果上各个网格点对应的幅值；

‖.‖₂的下标2表示取L2范数。

相关系数主要用来表征定位结果的位置精确度，公式如下：

其中：

x_t表示真值，即仿真声源各个网格点的幅值；

D(x_t)表示真值的方差；

x_r表示重构值，即仿真结果上各个网格点对应的幅值；

D(x_r)表示重构值的方差；

Cov(x_t,x_r)表示真值和重构值的协方差，具体表示如下：

Cov(x_t,x_r)＝E(x_t*x_r)-E(x_t)*E(x_r) (13)

其中：

E(x_t*x_r)表示真值和重构值之积的数学期望；

E(x_t)表示真值的数学期望；

E(x_r)表示重构值的数学期望。

根据公式(11)以及公式(12)，计算得出声源中心频率为300HZ，500HZ，1000HZ情况下定位结果的相对误差分别为38.2％，13.9％和4.35％；相关系数分别为0.8997，0.9894和0.9981。

由仿真结果可见，本实施例所述的面向低频高斯噪声源的定位方法能够有效对低频声源进行定位，能同时保证定位的位置精度和能量精度的准确性，而且还能够保持很高的分辨率。对比图4，5，6可以看出，随着信号中心频率的增高，定位的精度和分辨率都逐渐增高，中心频率增大到1000Hz时，几乎可以完全复原声源的位置和能量信息。

此外应理解，在阅读了本发明的上述描述内容之后，本领域技术人员可以对本发明作各种改动或修改，这些等价形式同样落于本申请所附权利要求书所限定的范围。

Claims (6)

1.一种面向低频高斯噪声源的声源定位仿真方法，包括：

(1)将声源所在平面作为声源扫描平面，将声源扫描平面划分成若干个网格；

(2)设置麦克风阵列采集声信号，输出麦克风阵列信号作为阵列输出信号，对阵列输出信号进行频谱分析以及时频谱分析，得到声信号的特征频率；

(3)对阵列输出信号表达式的等号两边的式子分别求协方差矩阵得到等号左侧协方差矩阵和等号右侧协方差矩阵，所述的阵列输出信号表达式为：

X＝AS+N

其中，X为阵列输出信号矩阵，A为麦克风阵列的方向矩阵，S为声信号矩阵，N为噪声矩阵；

(4)将等号左侧协方差矩阵按行展开，并转置成列向量形式，将得到的列向量作为压缩感知的测量值；

(5)将等号右侧协方差矩阵按行展开，得到压缩感知的矩阵表达式；

(6)通过测量值和矩阵表达式计算声源的能量向量，定位声源位置。

2.根据权利要求1所述的面向低频高斯噪声源的声源定位仿真方法，其特征在于，所述的麦克风阵列所在平面和声源扫描平面同轴设置。

3.根据权利要求1或2所述的面向低频高斯噪声源的声源定位仿真方法，其特征在于，所述的麦克风阵列中相邻麦克风阵元之间的距离不大于声信号的半波长。

4.根据权利要求1所述的面向低频高斯噪声源的声源定位仿真方法，其特征在于，所述的频谱分析采用快速傅里叶变换，时频谱分析采用短时傅里叶变换。

5.根据权利要求1所述的面向低频高斯噪声源的声源定位仿真方法，其特征在于，步骤(6)的具体步骤为：

(6-1)分别用L1范数和L2范数对所述测量值和矩阵表达式进行最小化约束；

(6-2)调用MATLAB的凸优化工具包计算声源的能量向量，定位声源位置。

6.根据权利要求1或5所述的面向低频高斯噪声源的声源定位仿真方法，其特征在于，在步骤(6)之后，将定位结果以图像化的形式展示在声源扫描平面上。