CN116559779A - 一种可迭代式频率聚焦变换的宽带声源定位方法 - Google Patents

Abstract

本发明公开了一种可迭代式频率聚焦变换的宽带声源定位方法，过程如下：通过麦克风阵列采集环境中的声音信号，将声音信号经过滤波等处理后通过波束形成方法得到预定位位置，然后通过可迭代式频率聚焦变换算法将宽带信号聚焦为窄带信号，并利用MUSIC算法对窄带信号进行定位估计，再将定位结果作为预定位位置迭代一次得到更准确的定位结果。本发明提供的可迭代式频率聚焦变换的宽带声源定位方法，解决了高分辨率谱估计算法MUSIC算法只能对窄带信号进行定位，而原有频率聚焦变换算法(如RSS算法)定位准确度对预定位准确度依赖程度过高的问题，本发明进行频率聚焦时对预定位位置准确度完全没有要求，且保证了频率聚焦的有效性。

Description

一种可迭代式频率聚焦变换的宽带声源定位方法

技术领域

本发明涉及阵列信号处理技术领域，尤其是涉及宽带声源定位方式，具体涉及一种可迭代式频率聚焦变换的宽带声源定位方法。

背景技术

在电力设备故障检测中，声学检测技术作为无损检测的典型代表，具有采集便捷、定位方便、检测精度高等优势，其中，使用声学成像技术可以实现声音可视化，可以让电力设备检测人员更加直观地看到电力设备内部的声源分布情况。

声学成像技术，又叫做声学相机，指的是使用麦克风阵列和声源定位技术计算声源分布情况与相机采集的图像进行融合并显示的一种实现声音可视化的方法。声源定位技术有传统的基于麦克风阵列的声源定位方法，也有微型仿生机械声源定位、利用光学传感方式定位等方法，这里研究基于麦克风阵列的声源定位方法。

基于麦克风阵列的声源定位技术有三大类算法，分别是基于波束形成方法、基于高分辨率谱估计方法、基于到达时间差方法，其中基于高分辨率谱估计方法有多重信号分类方法MUSIC和借助旋转不变性的信号参数估计(ESPRIT)两大类方法，多重信号分类方法MUSIC定位准确率高、分辨率高、可以实现多声源定位，是实现电力设备故障检测定位的一种很好的选择，但是MUSIC只针对于窄带信号进行计算，而电力设备故障时产生的异响大多数都是宽带信号，将宽带信号当作窄带信号进行定位时会造成有用信号的一个浪费，并且当有多个不同频率范围的声源需要定位时需要多次计算才能定位到多个声源，故研究人员又进而研究宽带声源定位算法。

后来Wang和Kaveh首先提出了一种相干信号子空间处理方法(Coherent Singal-Subspace Method，CSM)，这是一种有效的宽带信号估计方法，主要思想就是用传统波束形成算法对声源数和声源位置进行预估计，将宽带信号在频域上划分成多个窄带频率信号，再使用聚焦矩阵把各个窄带频率聚焦到选定的中心频率上，再对相关矩阵进行特征值分解，应用窄带声源定位算法进行声源定位来实现相干声源的定位。此外，还有双边相关变换(TCT)、旋转信号子空间(RSS)、信号子空间变换(SST)等聚焦矩阵算法也能实现宽带声源定位，但以上大多数算法最后的准确度都十分依赖预定位位置的准确度，而预定位位置本就是比较粗略的估计，所以这类算法的准确度无法保证。

发明内容

本发明的目的是为了解决现有技术中的上述缺陷，提供一种可迭代式频率聚焦变换的宽带声源定位方法，通过高分辨率谱估计实现更加准确的宽带声源定位，本发明中频率聚焦变换算法对预定位位置的准确度完全没有要求，通过迭代可自己提高定位准确度，实现更有效的频率聚焦。

本发明的目的可以通过采取如下技术方案达到：

一种可迭代式频率聚焦变换的宽带声源定位方法，所述宽带声源定位方法包括以下步骤：

S1、建立空间直角坐标系，确定传声器阵列中各传声器的位置、摄像头的位置、声源聚焦平面的位置，并使用传声器阵列采集声音信号X；

S2、选定目标频率范围的上界f_max和下界f_min，将采集的声音信号进行带通滤波得到目标频率范围的滤波声音信号X_f；

S3、将声源聚焦平面S划分为若干个子区域集合D₀，设D₀中的某个子区域为d_ij，其中i和j分别为子区域的行下标和列下标，对D₀中的所有d_ij采用波束形成方法计算声功率分布矩阵；

S4、设定阈值γ，使用谱峰值搜索方法筛选出聚焦平面的峰值点，其中超过阈值γ的峰值点认为可能为声源，筛选出所有超过阈值γ的峰值点的区域集合D₁作为预定位位置，预定位位置总数记为N；

S5、将滤波声音信号X_f和预定位位置D₁通过频率聚焦变换算法进行聚焦，从而将频率为(f_min，f_max)范围内的滤波声音信号X_f聚焦成中心频率为f₀的窄带信号上；

S6、利用多重信号分类MUSIC算法对窄带信号进行定位；

S7、将步骤S6的定位结果重复步骤S4、S5、S6再进行一次定位得到二次定位结果。

进一步地，所述步骤S1中，由于采用的是单目摄像头，单目摄像头采集的图像数据为平面图，所以设计传声器阵列为二维平面阵列，其中，二维均匀圆形阵列最为常见且适用于大多数声源定位算法，本发明传声器阵列为L通道的圆形阵列，圆形阵列所在平面与地面垂直，摄像头置于圆形阵列的几何中心位置，以摄像头光心坐标为原点，以垂直于阵列平面往声源方向为z轴、以垂直地面向上方向作为y轴建立空间直角坐标系，圆形阵列的半径为r，每两个传声器之间与摄像头所成的角度为则第l，l＝1，2，...，L个传声器在空间直角坐标系中的坐标(x_l，y_l，z_l)为：

由于阵列尺寸相对于电力设备与阵列之间的距离不能忽略，考虑到电力设备故障检测的应用场景，当在室内情况下进行故障检测时检测距离为1～5m，在室外情况下进行故障检测时检测距离应为3～50m，因此声源模型设置为远场平面波模型，设定声源聚焦平面为S，声源聚焦平面划分子区域的行数和列数分别为N_x和N_y，声源聚焦平面S距离阵列平面的距离为d。

传声器阵列采集声音信号的长度可根据传声器的采样率f_s计算，摄像头采集图像实现流畅播放的采集速度为30帧/秒，声源定位结果要与摄像头图像同步融合则与摄像头采集速度一致，则一次声音数据采集的长度为则L通道的传声器阵列采集的声音信号X为/>的矩阵。

进一步地，所述步骤S3中，波束形成方法是一种复杂度低、定位带宽宽的高速低分辨率的声源定位方法，所以用波束形成方法计算预定位位置不会给后续的定位增加复杂度又能得到宽带声源的粗略的位置，波束形成方法计算声功率矩阵的步骤为：

S3.1、计算波束形成器在期望方向上的加权向量的构成为 其中ω＝2π*f₀，τ_(i，j，l)表示声源聚焦平面S中第i行第j列的子区域到第l个传声器的声音时延，声音时延用空间距离s除以声速v计算，

S3.2、计算声音信号X的协方差矩阵R：

上式中，X^H为声音信号X的共轭转置矩阵；

S3.3、计算波束形成器的输出功率，公式如下：

P_CBF(i，j)＝ω^H(i，j)*R*ω(i，j) (3)

其中，P_CBF(i，j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度。

进一步地，波束形成方法计算所得的声源聚焦平面S的声功率强度的分布中，声功率强度值越大的子区域越有可能存在声源，假设超过阈值的峰值点认为可能为声源，低于阈值的峰值点认为是环境中的噪声或者混响，所以需要设定阈值筛选出可能为声源的区域集合。所述步骤S4中阈值γ由下式确定：

γ＝λ(P_max-P_min)+P_min (4)

上式中，P_max为D₀中子区域的声功率最大值，P_min为D₀中子区域的声功率最小值，λ为缩放系数，λ的参考取值范围是0.5＜λ＜1；

其中，谱峰值搜索方法的搜索依据为

上式中，P(i，j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度，同理，P(i-1，j)、P(i+1，j)、P(i，j-1)、P(i，j+1)分别为声源聚焦平面S中与第i行第j列子区域相邻的左、右、下、上子区域的声功率强度，对波束形成方法的声功率强度结果进行谱峰值搜索时P(i，j)＝P_CBF(i，j)，对其他方法的声功率强度或谱函数估计结果进行谱峰值搜索也同理，将筛选出来的所有子区域d_ij的区域集合记为D₁。

进一步地，由于MUSIC算法是窄带声源定位算法，即MUSIC算法进行定位时选取的定位频率为f₀时，所能定位到的带宽B＜＜f₀，所以使用MUSIC算法实现宽带声源定位最重要的是将滤波声音信号的X_f的其他频率值上的分量聚焦到频率f₀上来，由于已有的RSS算法等使用酉矩阵进行频率聚焦，定位效果直接受到预定位位置的极大影响，故在此发明中，提出了一种非酉矩阵进行频率聚焦的方法。

步骤S5中，使用非酉矩阵进行频率聚焦的频率聚焦变换算法过程如下：

S5.1、首先将频率范围为(f_min，f_max)范围内的滤波声音信号X_f划分为K个子带，各个频率分量子带为X(f_k)，k＝1，……，K；

S5.2、已知在步骤S4中得到预定位位置区域集合D₁，则采用频率聚焦矩阵为

上式中，为频率为f₀的L*N_pre维的初始导向矩阵，E表示L*(L-N_pre)的单位矩阵，/>为频率为f_k的L*N_pre维的初始导向矩阵，计算方法为

上式中，τ_(l，n)表示第l个传声器到第n个预定位位置的时延，n＝1，……，N_pre，已知D₁为步骤S4中筛选出的可能为声源的子区域集合，经过式(6)可以对集合D₁中各个子区域位置的声音功率在频率为f₀的初始导向矩阵进行加成，在频率为f_k的初始导向矩阵进行消去，即在声源发出时，将频率为f_k的声音信号分量按照频率为f₀的声音信号分量的传播时延来计算，以此来达到频率聚焦到f₀的目的；

S5.3、使用频率聚焦矩阵对各个频率分量子带进行频率聚焦得到中心频率为f₀的窄带信号

S5.4、聚焦变换后，各频率点下的方向矩阵所包含的频率信息相等，对聚焦后阵列各频率点下的协方差矩阵求和、平均，得到协方差矩阵R_y：

上式中，R_y表示L*L维的对称矩阵。

进一步地，多重信号分类MUSIC算法分辨率高，针对窄带信号定位效果好，所述步骤S6中，使用多重信号分类MUSIC对聚焦后的窄带信号进行谱函数估计步骤为：

S6.1、根据窄带信号计算协方差矩阵R_y

S6.2、对上式得到的协方差矩阵进行特征值分解

R_y＝U∑U^H (12)

上式中，U表示矩阵R_y的特征向量组成的L*L维的矩阵，∑表示矩阵R_y的特征值组成的L*L维的对角矩阵；

S6.3、按特征值的大小顺序，把与预定位位置数量N_pre相等的最大特征值对应的特征向量看成信号子空间，把剩下的L-N_pre个特征值对应的特征向量看成噪声子空间，即

上式中，∑_S表示前N_pre个最大特征值组成的L*N_pre维的对角矩阵，U_S表示∑_S中的特征值对应的特征向量组成的L*N_pre维的矩阵，∑_N表示剩下L-N_pre个特征值组成的L*(L-N_pre)维的对角矩阵，U_N表示∑_N中的特征值对应的特征向量组成的L*(L-N_pre)维的矩阵；

S6.4、对声源聚焦平面S重新扫描计算谱函数估计P_MUSIC(i，j)，估计公式为

上式中，P_MUSIC(i，j)表示声源聚焦平面S中第i行第j列的子区域的谱函数估计值，a(i，j)表示声源聚焦平面S中第i行第j列的子区域与传声器阵列的导向向量，公式为a(i，j)＝w(i，j)。

本发明相对于现有技术具有如下的优点及效果：

(1)本发明公开的可迭代式频率聚焦变换的宽带声源定位方法，提出了一种可迭代式频率聚焦变换算法构造聚焦矩阵，进行频率聚焦时不再受预定位位置准确度的影响，并且可以通过迭代提高准确度。

(2)本发明使用频率聚焦方法对宽带声源进行聚焦，最大限度地利用声音信号中有用信息。当有多个不同频率范围的声源需要进行定位时，使用常规的窄带多重信号算法需要对不同频率范围的声源选择其中心频率进行多次定位，而进行频率聚焦后只需进行一次定位就能定位到多个不同频率范围的声源。

附图说明

此处所说明的附图用来提供对本发明的进一步理解，构成本申请的一部分，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：

图1是本发明公开的可迭代式频率聚焦变换的宽带声源定位方法总体流程图；

图2是本发明实施例公开的传声器阵列-声源聚焦平面模型图；

图3是本发明实施例公开的本发明、RSS、SST聚焦算法的信噪比-相对误差比对图；

图4是本发明实施例公开的宽带声源定位方法、RSS-MUSIC、传统MUSIC、常规波束形成法的窄带声源中心频率-相对误差对比图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例1

本实施例公开了一种可迭代式频率聚焦变换的宽带声源定位方法，提出了一种可迭代式频率聚焦变换算法构造聚焦矩阵，进行频率聚焦时不再受预定位位置准确度的影响，并且可以通过迭代提高准确度，且使用频率聚焦方法对宽带声源进行聚焦之后可以最大限度地利用到声音信号里的有用信息，当有多个不同频率范围的声源需要进行定位时，使用常规的窄带多重信号算法需要对不同的频率范围的声源选择其中心频率进行多次定位，而进行频率聚焦后只需进行一次定位就能定位到多个不同频率范围的声源。

图1是本发明实施例公开的可迭代式频率聚焦变换的宽带声源定位方法总体流程图，具体包含如下步骤：

S1、建立空间直角坐标系，确定传声器阵列中L个传声器的位置、摄像头的位置、声源聚焦平面的位置，具体情况如附图2所示，并使用传声器阵列采集声音信号X；

本实施例中，传声器阵列为L通道的圆形阵列，圆形阵列所在平面与地面垂直，摄像头置于圆形阵列的几何中心位置，以摄像头光心坐标为原点，以垂直于阵列平面往声源方向为z轴、以垂直地面向上方向作为y轴建立空间直角坐标系，圆形阵列的半径为r，每两个传声器之间与摄像头所成的角度为则第l(l＝1，2，...，L)个传声器在空间直角坐标系中的坐标(x_l，y_l，z_l)为：

声源模型设置为远场平面波模型，设定声源聚焦平面为S，聚焦平面网格点的行数和列数分别为N_x和N_y，声源聚焦平面距离阵列平面的距离为d；

根据传声器的采样率f_s计算声音信号X的长度，设置一次声音数据采集的长度为则L通道的传声器阵列采集的声音信号X为/>的矩阵。

S2、选定目标频率范围的上界f_max和下界f_min，将采集的声音信号进行带通滤波得到目标频率范围的滤波声音信号X_f；

S3、将声源聚焦平面S划分为若干个子区域集合D₀，设D₀中的某个子区域为d_ij，其中i和j分别为子区域的行下标和列下标，对D₀中的所有d_ij采用波束形成方法计算声功率分布矩阵；

在本实施例中，波束形成方法计算声功率矩阵的步骤为：

S3.1、计算波束形成器在期望方向上的加权向量的构成为 其中ω＝2π*f₀，τ_(i，j，l)表示声源聚焦平面S中第i行第j列的子区域到第l个传声器的声音时延，i＝1，…，N_x，j＝1，…，N_y，声音时延用空间距离s除以声速c计算，/>

S3.2、计算声音信号X的协方差矩阵R：

上式中，X^H为声音信号X的共轭转置矩阵；

S3.3、计算波束形成器的输出功率，公式如下：

P_CBF(i，j)＝ω^H(i，j)*R*ω(i，j) (3)

其中，P_CBF(i，j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度。

S4、设定阈值γ，使用谱峰值搜索方法筛选出声源聚焦平面的峰值点，其中超过阈值γ的峰值点认为可能为声源，筛选出所有超过阈值γ的峰值点的区域集合D₁作为预定位位置，预定位位置总数记为N_pre；

在本实施例中，阈值γ由下式确定：

γ＝λ(P_max-P_min)+P_min (4)

上式中，P_max为D₀中子区域的声功率最大值，P_min为D₀中子区域的声功率最小值，λ为缩放系数，λ的参考取值范围是0.5＜λ＜1；

其中，谱峰值搜索方法的搜索依据为

上式中，P(i，j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度，同理，P(i-1，j)、P(i+1，j)、P(i，j-1)、P(i，j+1)分别为声源聚焦平面S中与第i行第j列子区域相邻的左、右、下、上子区域的声功率强度，对波束形成方法的声功率强度结果进行谱峰值搜索时P(i，j)＝P_CBF(i，j)，对其他方法的声功率强度或谱函数估计结果进行谱峰值搜索也同理，将筛选出来的所有子区域d_ij的区域集合记为D₁。

S5、将滤波声音信号X_f和预定位位置D₁通过频率聚焦变换算法进行聚焦，从而将频率为(f_min，f_max)范围内的滤波声音信号X_f聚焦成中心频率为f₀的窄带信号上；

在本实施例中，频率聚焦变换算法过程如下：

S5.1、首先将频率范围为(f_min，f_max)范围内的滤波声音信号X_f划分为K个子带，各个频率分量子带为X(f_k)，k＝1，……，K；

S5.2、已知在步骤S4中得到预定位位置区域集合D₁，则采用频率聚焦矩阵为

上式中，为频率为f₀的L*N_pre维的初始导向矩阵，E表示L*(L-N_pre)的单位矩阵，/>为频率为f_k的L*N_pre维的初始导向矩阵，计算方法为

上式中，τ_(l，n)表示第l个传声器到第n个预定位位置的时延，n＝1，……，N_pre；

S5.3、使用频率聚焦矩阵对各个频率分量子带进行频率聚焦得到中心频率为f₀的窄带信号

S5.4、聚焦变换后，各频率点下的方向矩阵所包含的频率信息相等，对聚焦后阵列各频率点下的协方差矩阵求和、平均，得到协方差矩阵R_y：

上式中，R_y表示L*L维的对称矩阵。

S6、利用多重信号分类MUSIC算法对窄带信号进行定位；

在本实施例中，使用多重信号分类MUSIC对聚焦后的窄带信号进行谱函数估计步骤为：

S6.1、根据窄带信号计算协方差矩阵R_y

S6.2、对上式得到的协方差矩阵进行特征值分解

R_y＝U∑U^H (12)

上式中，U表示矩阵R_y的特征向量组成的L*L维的矩阵，∑表示矩阵R_y的特征值组成的L*L维的对角矩阵；

S6.3、按特征值的大小顺序，把与预定位位置数量N_pre相等的最大特征值对应的特征向量看成信号子空间，把剩下的L-N_pre个特征值对应的特征向量看成噪声子空间，即

上式中，∑_S表示前N_pre个最大特征值组成的L*N_pre维的对角矩阵，U_S表示∑_S中的特征值对应的特征向量组成的L*N_pre维的矩阵，∑_N表示剩下L-N_pre个特征值组成的L*(L-N_pre)维的对角矩阵，U_N表示∑_N中的特征值对应的特征向量组成的L*(L-N_pre)维的矩阵；

S6.4、对声源聚焦平面S重新扫描计算谱函数估计P_MUSIC(i，j)，估计公式为

上式中，P_MUSIC(i，j)表示声源聚焦平面S中第i行第j列的子区域的谱函数估计值，a(i，j)表示声源聚焦平面S中第i行第j列的子区域与传声器阵列的导向向量，公式为a(i，j)＝w(i，j)。

S7、将步骤S6的定位结果重复步骤S4、S5、S6再进行一次定位得到二次定位结果。

实施例2

为了验证本发明提出的方法在定位宽带声源的有效性和准确性，现通过将算法运行在AMD Ryzen 7 5800H with Radeon Graphics 3.20GHz计算机和matlab2021a仿真平台上。

分别比较几种频率聚焦变换的宽带声源定位算法的准确度，下面结合实施例1的步骤，将已有的频率聚焦算法RSS聚焦算法和SST聚焦算法与本发明提出的聚焦算法进行比较，本发明最大的优点在于不依赖预定位准确度，本实施例采用随机选取的位置作为预定位位置，即在预定位位置极为不准确的情况下如果本发明的频率聚焦算法仍能定位到准确位置则可验证本发明定位方法的准确性。

可迭代式频率聚焦变换的宽带声源定位方法与已有的宽带声源RSS-MUSIC方法和SST-MUSIC方法进行比较仿真包括以下步骤：

S1、建立空间直角坐标系，确定传声器阵列中各传声器的位置、声源聚焦平面的位置，在声源聚焦平面内的设置随机声源位置P_t，并根据传声器阵列位置和声源位置计算出虚拟的接收声音信号X；

本实施例中，传声器阵列为L＝16通道的圆形阵列，圆形阵列所在平面与地面垂直，以圆形阵列的几何中心位置为原点，以垂直于阵列平面往声源方向为z轴、以垂直地面向上方向作为y轴建立空间直角坐标系，圆形阵列的半径为r＝0.07m，每两个传声器之间与原点所成的角度为则第l个传声器在空间直角坐标系中的坐标(x_l，y_l，z_l)为式(1)：

声源模型设置为远场平面波模型，设定声源聚焦平面为S，声源聚焦平面为S距离阵列平面的距离d根据声源的远近而定，可设d为1～50m的随机整数，则声源聚焦平面为S的大小为2d×2d的矩形平面，具体如附图2所示，聚焦平面网格点的行数和列数分别为N_x和N_y，均设为100。

设传声器阵列的采样率为f_s＝48kHz，设声源数为M＝1，单通道的数据采集长度为length＝1600，通过远场宽带声源信号模型生成频率范围为(6000，12000)Hz的宽带声源信号，设宽带声源信号划分次数W＝11，则划分的子频率带f_w∈[6500，7000，7500，…，11500]，宽带声源信号的位置P_t可表示为

则生成宽带声源声音信号的公式如下：

X(f_w)＝A(f_w)S(f_w)+N，w＝1，……，W (16)

上式中，f_w表示范围为(6000，12000)Hz中划分的子带频率，A(f_w)表示频率为f_w的L*M维的初始导向矩阵，τ_l表示的是第l个传声器到声源位置P_test的时延，S(f_w)表示中心频率为f_w的M*length维度的声音信号矩阵，N表示L*length维度的噪声信号矩阵，定义信噪比为0dB。

S2、选定目标频率范围的上界f_max＝6000Hz和下界f_min＝12000Hz，将采集的声音信号进行带通滤波得到目标频率范围的滤波声音信号X_f；

S3、在聚焦平面内选择M个随机位置作为预定位位置，预定位位置记为P_pre，可能的声源数为M；

S4、同实施例1中S5计算出聚焦后窄带信号的协方差矩阵R_y-本发明，并分别计算RSS聚焦算法和SST聚焦算法的协方差矩阵R_y-RSS、R_y-SST；

S5、同实施例1中S6步骤计算本发明的MUSIC算法的谱函数估计P_{MUSIC-本发明}以及RSS聚焦算法和SST聚焦算法的MUSIC算法的谱函数估计P_MUSIC-RSS、P_MUSIC-SST；

S6、使用实施例1中步骤S4的谱峰值搜索方法分别计算P_{MUSIC-本发明}、P_MUSIC-RSS、P_MUSIC-SST的声源定位位置P_{estimate-本发明}、P_estimate-RSS、P_estimate-SST；

S7、将步骤S6的定位结果重复步骤S4、S5、S6再进行一次定位得到二次定位的结果，并与P_t相比计算出相对误差；

S8、重复50～100组以上步骤，统计误差并求取平均值；

S9、设置信噪比snr值为snr＝-15，-10，-5，5，10，15,20,25,30再重复以上步骤，对比不同信噪比下三种算法的准确度。

得到的比较结果如附图3所示，由对比图可知，在预定位位置极为不准确的情况下，本发明的可迭代式频率聚焦算法的定位结果的相对误差比较RSS聚焦算法和SST聚焦算法都要低，且在任何信噪比下都有较好的准确性，而RSS聚焦算法和SST聚焦算法由于非常依赖预定位位置的准确度，在进行一次迭代后仍然没有提高准确度，而本发明的可迭代式频率聚焦算法对预定位位置的准确度没有要求，所以在预定位位置极为不准确的情况下仍然能表现出较好的定位准确性。

实施例3

为了验证本发明提出的方法在定位宽带声源的有效性和准确性，现通过将算法运行在AMD Ryzen 7 5800H with Radeon Graphics 3.20GHz计算机和matlab2021a仿真平台上。

当有多个不同频率范围的声源需要进行定位时，使用常规的窄带多重信号分类算法需要对不同频率范围的声源选择其中心频率进行多次定位，而进行频率聚焦后只需进行一次定位就能定位到多个不同频率范围的声源，为了验证本发明的算法对宽频带内多声源的宽频性，将传统MUSIC算法、常规波束形成算法、RSS聚焦算法的MUSIC算法和基于本发明的可迭代式频率聚焦变换的MUSIC算法进行比较，来验证频率聚焦的有效性。

可迭代式频率聚焦变换的宽带声源定位方法仿真包括以下步骤：

S1、建立空间直角坐标系，确定传声器阵列中各传声器的位置、声源聚焦平面的位置，在声源聚焦平面内的设置随机声源位置P_test，并根据传声器阵列位置和声源位置计算出虚拟的接收声音信号X；

本实施例中，传声器阵列为L＝16通道的圆形阵列，圆形阵列所在平面与地面垂直，以圆形阵列的几何中心位置为原点，以垂直于阵列平面往声源方向为z轴、以垂直地面向上方向作为y轴建立空间直角坐标系，圆形阵列的半径为r＝0.07m，每两个传声器之间与摄像头所成的角度为则第l个传声器在空间直角坐标系中的坐标(x_l，y_l，z_l)为式(1)：

由于阵列尺寸相对于电力设备与阵列之间的距离不能忽略，考虑到电力设备故障检测的应用场景，当在室内情况下进行故障检测时检测距离为1～5m，在室外情况下进行故障检测时检测距离应为3～50m，因此声源模型设置为远场平面波模型，设定声源聚焦平面为S，声源聚焦平面为S距离阵列平面的距离d根据声源的远近而定，可设d为1～50m的随机整数，则声源聚焦平面为S的大小为2d×2d的矩形平面，具体如附图所示，聚焦平面网格点的行数和列数分别为N_x和N_y，均设为100。

设传声器阵列的采样率为f_s＝48kHz，设声源数为2或者3个声源，即M∈[2，3]，单通道的数据采集长度为length＝1600，本实例使用中心频率f_m∈[3200Hz，3400Hz，3600Hz，…，5800Hz]的M个窄带声源信号进行测试，声源信号的位置P_test表示随机生成的M个声源位置，第m个声源位置(x_m，y_m，z_m)可表示为

窄带信号的定义为带宽远小于中心频率，即B＜＜f，第m个窄带声源信号模型公式如下：

X_m＝A(f_m)S(f_m)+N (19)

上式中，A(f_m)表示频率为f_m的L*1维的初始导向矩阵，τ_(l，m)表示第l个传声器到第m个声源位置的时延，S(f_m)表示中心频率为f_m的长度为length的声音信号向量，N表示L*length维度的噪声信号矩阵，定义信噪比为0dB，则多个声源的声音信号X的公式为

S2、选定目标频率范围的上界f_max＝3000Hz和下界f_min＝6000Hz，将采集的声音信号进行带通滤波得到目标频率范围的滤波声音信号X_f；

S3、同实施例1中S3；

S4、同实施例1中S4；

S5、同实施例1中S5计算出聚焦后窄带信号的协方差矩阵R_y-本发明，以及传统MUSIC算法的协方差矩阵R_y-tra和RSS聚焦算法的协方差矩阵R_y-RSS，传统MUSIC算法的协方差矩阵计算公式同式(2)；

S5、同实施例1中S6步骤计算本发明的谱函数估计P_{MUSIC-本发明}，用传统MUSIC算法的谱函数估计P_MUSIC-tra，以及RSS聚焦算法的MUSIC算法的谱函数估计P_MUSIC-RSS，三种算法的定位频率均为

S6、使用实施例1中步骤S4的谱峰值搜索方法分别计算P_{MUSIC-本发明}、P_MUSIC-tra和P_MUSIC-RSS的声源定位位置P_{estimate-本发明}、P_estimate-tra、P_estimate-RSS，并与P_test相比计算出相对误差；

S8、重复100～200组以上步骤，根据窄带声源的中心频率和不同的估计方法，统计误差并求取平均值，记录δ_本发明(f)，δ_MUSIC-tra(f)，δ_MUSIC-RSS(f)，δ_波束形成(f)，f＝3200Hz，3400Hz，3600Hz，…，5800Hz。

得到的比较结果如附图4所示，由对比图可知，在四种算法的定位频率都在中心频率f₀＝4500Hz的情况下，在窄带声源的中心频率为f₀时各算法的相对误差最小，而窄带声源的中心频率越远离定位频率时相对误差逐渐升高，四种算法的定位结果都呈现出这一趋势。其中，本发明方法定位结果的相对误差最低，波束形成法的定位结果为预定位位置，RSS聚焦算法的MUSIC算法由于太依赖于预定位位置的准确度而受到影响，而本发明的可迭代式频率聚焦变换的宽带声源定位方法不受预定位结果影响，即使有多个中心频率不同的声源同时出现，只需要一个定位频率就能定位到某个更宽的频率范围的多个声源，呈现出更好的聚焦效果和定位准确度，验证了本发明方法的有效性。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (6)

1.一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在于，所述宽带声源定位方法包括以下步骤：

S1、建立空间直角坐标系，确定传声器阵列中L个传声器的位置、摄像头的位置、声源聚焦平面的位置，并使用传声器阵列采集声音信号X；

S2、选定目标频率范围的上界f_max和下界f_min，将采集的声音信号进行带通滤波得到目标频率范围的滤波声音信号X_f；

S3、将声源聚焦平面S划分为若干个子区域集合D₀，设D₀中的某个子区域为d_ij，其中i和j分别为子区域的行下标和列下标，对D₀中的所有d_ij采用波束形成方法计算声功率分布矩阵；

S4、设定阈值γ，使用谱峰值搜索方法筛选出声源聚焦平面的峰值点，其中超过阈值γ的峰值点认为可能为声源，筛选出所有超过阈值γ的峰值点的区域集合D₁作为预定位位置，预定位位置总数记为N_pre；

S5、将滤波声音信号X_f和预定位位置D₁通过频率聚焦变换算法进行聚焦，从而将频率为(f_min,f_max)范围内的滤波声音信号X_f聚焦成中心频率为f₀的窄带信号上；

S6、利用多重信号分类MUSIC算法对窄带信号进行定位；

S7、将步骤S6的定位结果重复步骤S4、S5、S6再进行一次定位得到二次定位结果。

2.根据权利要求1所述的一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在于，所述步骤S1中，传声器阵列为L通道的圆形阵列，圆形阵列所在平面与地面垂直，摄像头置于圆形阵列的几何中心位置，以摄像头光心坐标为原点，以垂直于阵列平面往声源方向为z轴、以垂直地面向上方向作为y轴建立空间直角坐标系，圆形阵列的半径为r，每两个传声器之间与摄像头所成的角度为则第l个传声器在空间直角坐标系中的坐标(x_l,y_l,z_l)为：/>

其中，l＝1,2,…,L，将声源模型设置为远场平面波模型，设定声源聚焦平面为S，聚焦平面网格点的行数和列数分别为N_x和N_y，声源聚焦平面距离阵列平面的距离为d；

根据传声器的采样率f_s计算声音信号X的长度，设置一次声音数据采集的长度为则L通道的传声器阵列采集的声音信号X为/>的矩阵。

3.根据权利要求2所述的一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在于，所述步骤S3中，波束形成方法计算声功率矩阵的步骤为：

S3.1、计算波束形成器在期望方向上的加权向量的构成为 其中ω＝2π·f₀，τ_(i,j,l)表示声源聚焦平面S中第i行第j列的子区域到第l个传声器的声音时延，i＝1,…,N_x,j＝1,…,N_y，声音时延用空间距离s除以声速c计算，/>

S3.2、计算声音信号X的协方差矩阵R：

上式中，X^H为声音信号X的共轭转置矩阵；

S3.3、计算波束形成器的输出功率，公式如下：

P_CBF(i,j)＝ω^H(i,j)*R*ω(i,j) (3)

其中，P_CBF(i,j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度。

4.根据权利要求3所述的一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在于，所述步骤S4中阈值γ由下式确定：

γ＝λ(P_max-P_min)+P_min (4)

上式中，P_max为D₀中子区域的声功率最大值，P_min为D₀中子区域的声功率最小值，λ为缩放系数，λ的参考取值范围是0.5＜λ＜1；

其中，谱峰值搜索方法的搜索依据为

上式中，P(i,j)表示声源聚焦平面S中第i行第j列的子区域的声功率强度，P(i-1,j)、P(i+1,j)、P(i,j-1)、P(i,j+1)分别为声源聚焦平面S中与第i行第j列子区域相邻的左、右、下、上子区域的声功率强度，对波束形成方法的声功率强度结果进行谱峰值搜索时P(i,j)＝P_CBF(i,j)，将筛选出来的所有子区域d_ij的区域集合记为D₁。

5.根据权利要求1所述的一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在于，所述步骤S5中，频率聚焦变换算法过程如下：

S5.1、首先将频率范围为(f_min,f_max)范围内的滤波声音信号X_f划分为K个子带，各个频率分量子带为X(f_k),k＝1,……,K；

S5.2、已知在步骤S4中得到预定位位置区域集合D₁，则采用频率聚焦矩阵为

上式中，为频率为f₀的L*N_pre维的初始导向矩阵，E表示L*(L-N_pre)的单位矩阵，为频率为f_k的L*N_pre维的初始导向矩阵，计算方法为

上式中，τ_(l,n)表示第l个传声器到第n个预定位位置的时延，n＝1,……,N_pre；

S5.3、使用频率聚焦矩阵对各个频率分量子带进行频率聚焦得到中心频率为f₀的窄带信号

S5.4、聚焦变换后，各频率点下的方向矩阵所包含的频率信息相等，对聚焦后阵列各频率点下的协方差矩阵求和、平均，得到协方差矩阵R_y：

上式中，R_y表示L*L维的对称矩阵。

6.根据权利要求5所述的一种可迭代式频率聚焦变换的宽带声源定位方法，其特征在进行于，所述步骤S6中，使用多重信号分类MUSIC算法对聚焦后的窄带信号进行谱函数估计的过程如下：

S6.1、根据窄带信号计算协方差矩阵R_y：

S6.2、对上式得到的协方差矩阵进行特征值分解

R_y＝U∑U^H (12)

上式中，U表示矩阵R_y的特征向量组成的L*L维的矩阵，∑表示矩阵R_y的特征值组成的L*L维的对角矩阵；

S6.3、按特征值的大小顺序，把与预定位位置数量N_pre相等的最大特征值对应的特征向量看成信号子空间，把剩下的L-N_pre个特征值对应的特征向量看成噪声子空间，即

上式中，Σ_S表示前N_pre个最大特征值组成的L*N_pre维的对角矩阵，U_S表示∑_S中的特征值对应的特征向量组成的L*N_pre维的矩阵，∑_N表示剩下L-N_pre个特征值组成的L*(L-N_pre)维的对角矩阵，U_N表示∑_N中的特征值对应的特征向量组成的L*(L-N_pre)维的矩阵；

S6.4、对声源聚焦平面S重新扫描计算谱函数估计P_MUSIC(i,j)，估计公式为

上式中，P_MUSIC(i,j)表示声源聚焦平面S中第i行第j列的子区域的谱函数估计值，a(i,j)表示声源聚焦平面S中第i行第j列的子区域与传声器阵列的导向向量，并且，a(i,j)＝w(i,j)。