一种用于分布式麦克风阵列的声源定位系统及定位方法

摘要

本发明公开了一种用于分布式麦克风阵列的声源定位系统及其定位方法，所述的系统包括信号能量比计算模块A、数据库建立模块、GMM训练模块、信号能量比计算模块B和声源坐标计算模块；信号能量比计算模块A的输出端与数据库建立模块输入端连接；数据库建立模块输出端与GMM训练模块输入端连接，GMM训练模块输出端与声源坐标计算模块输入端连接；信号能量比计算模块B的输出端也与声源坐标计算模块输入端连接；声源坐标计算模块输出端输出声源坐标。由于本发明采用麦克风接收信号能量作为主要参数，仅与该麦克风和声源间的距离有关，故无需考虑阵列几何结构与尺寸信息，也无需考虑声音信号传输延迟时间、平稳性、宽带或窄带等因素。

说明书

技术领域

本发明涉及一种麦克风阵列的声源定位技术，特别是一种用于分布式麦克风阵列的声源定位系统及定位方法。

背景技术

麦克风阵列由多个麦克风排列组成，各麦克风同时采集并联合处理声音信号，在得到信号时域与频域信息的同时，还能够获得声音的空域(空间方位)信息。用麦克风阵列和目标的声音信息来获取目标位置是一种有效的定位方法。由于麦克风拾取声音属于被动探测方式，基于麦克风阵列的声源定位技术在军事探测、安全监控、音/视频会议系统、人机交互等领域应用广泛。

近年来，随着网络通信与移动计算技术的发展，麦克风阵列由规则几何拓扑结构发展为无特定几何形状和尺寸信息的分布式结构，这使得阵列的搭建、更改与应用更加灵活，但也带来了许多新的问题，导致已有的规则麦克风阵列定位理论和方法无法直接应用于分布式麦克风阵列中。

苏淳于2012年在公开号为CN102707262A的专利《一种基于麦克风阵列的声源定位系统》中，采用基于到达时间差的声源定位算法搭建了真实声场环境下基于麦克风阵列的声源定位系统。该系统先用广义互相关函数方法估计麦克风间接收信号的延迟时间，之后根据阵列的几何结构，建立几何关系的方程组求得声源位置。基于到达时间差的声源定位方法应用于传统的、结构规则且尺寸较小的麦克风阵列中，声源定位效果较好。但对于麦克风间距离较大的分布式阵列，受声音信号波长的影响，尤其对于（准）周期声音信号，直接估计麦克风间信号延迟时间会出现周期模糊现象，导致估计结果出现偏差甚至完全错误。此外，用计算两个信号互相关系数的方法来获得麦克风信号间的时延，需要大量数据参与运算，其计算复杂度较高，进行实时处理较难。

Carlos T I,Olivier C.发表的《Evaluation of a MUSIC-based real-time soundlocalization of multiple sound sources in real noisy environments 》 (InternationalConference on Intelligent Robots and Systems,Saint Louis,USA,2009:2027-2032.)基于空间谱估计技术来实现麦克风阵列声源定位。该方法建立各麦克风接收信号的相关矩阵，对其进行空间谱分解，划分出互相正交的信号子空间和噪声子空间，再估计空间谱峰来确定波达方向，进而确定声源位置。该技术虽然是阵列信号处理领域研究最广的定位方法，但其要求阵列中麦克风等间隔线性排列，而对于几何结构不规则的分布式阵列，该方法不再适用，并且该方法多适用于远场、窄带、平稳的声音信号，应用局限性较大。

发明内容

为解决现有技术存在的上述问题，本发明要实现一种算法简单、定位过程计算量小、易于实时实现，适用于几何结构不规则的分布式麦克风阵列的声源定位系统及定位方法。

为了实现上述目的，本发明的技术方案如下：一种用于分布式麦克风阵列的声源定位系统，包括信号能量比计算模块A、数据库建立模块、GMM训练模块、信号能量比计算模块B和声源坐标计算模块；所述的信号能量比计算模块A的输入端输入指纹采集点处各麦克风采集到的信号，输出端与数据库建立模块输入端连接；数据库建立模块输出端与GMM训练模块输入端连接，GMM训练模块输出端与声源坐标计算模块输入端连接；所述的信号能量比计算模块B的输入端输入未知位置处各麦克风采集到的信号，输出端也与声源坐标计算模块输入端连接；声源坐标计算模块输出端输出声源坐标；所述的GMM为高斯混合模型Gaussian Mixture Model的简称。

信号能量比计算模块A根据在各个指纹采集点处采集到的信号生成训练矢量，数据库建立模块将该训练矢量存储起来构成数据库，GMM训练模块根据数据库中的数据进行训练，获得GMM参数，该参数供声源坐标计算模块使用；信号能量比计算模块B采集未知位置的信号生成观测矢量，声源坐标计算模块根据该矢量和GMM参数计算出该未知位置处的声源坐标。

一种用于分布式麦克风阵列的声源定位方法，包括以下步骤：

A、在已部署好的麦克风阵列工作场地，将声源移动至事先确定好的各个指纹采集点处，将采集到的信号输入到信号能量比计算模块A，进行信号能量比计算，输出计算结果至数据库建立模块；所述的信号能量比计算对每个指纹采集点处采集的每帧信号进行一次计算；所述的信号能量比计算模块A的工作方法包括以下步骤：

首先检测输入信号的有效性，定义分布式麦克风阵列中麦克风数目为M，第1个麦克风接收到的输入为s₁(n)，当前帧在原始信号序列中的起始位置为start，信号帧长为L=qf_s，qf_s为信号采样频率，q为帧长对应的等效时间，其取值范围是0.05<q<0.5；则第1个麦克风的当前帧短时对数能量E₁为

$>E_1=\underset{n=\mathrm{start}}{\overset{\mathrm{start}+L-1}{\Sigma }}{s}_1^2\left(n\right)$>

将E₁与信号/噪声判定阈值E_thr进行比较，若E₁<E_thr，判定当前帧为噪声帧，否则，判定当前帧为信号帧；若当前帧为噪声帧，则估计背景噪声能量并更新阈值E_thr;更新阈值E_thr的方法为

if  E_min＜E₁

  E_min＝E₁

  E_thr＝aE_min

end if

其中，E_min表示记录的第1个麦克风最小短时能量，a为一个门限因子，其取值范围1<a<100；

定义第m个麦克风接收到输入为s_m(n)，若当前帧为噪声帧，则第m个麦克风的背景噪声E_m估计方法为

$>E_m=\underset{n=\mathrm{start}}{\overset{\mathrm{start}+N-1}{\Sigma }}{s}_m^2\left(n\right),m=\mathrm{1,2},···,M$>

ε_m,curr＝(1-b)·ε_m,old+b·E_m

其中，ε_m,curr为第m个麦克风当前帧背景噪声能量估计值，ε_m,old为第m个麦克风前一次背景噪声能量估计值，b为平滑因子，表示平滑过程中当前帧能量E_m所占比重，其取值范围0<b<1；

若当前帧为信号帧，去除第1个和第m个麦克风背景噪声能量的影响，之后以第1个麦克风信号能量作参考，计算第m个麦克风与第1个麦克风的修正后的信号能量比er_m，

$>{\mathrm{er}}_m=\frac{E_m-{ϵ}_{m,\mathrm{curr}}}{E_1-{ϵ}_1}$>

然后将er_m，m=2,3,…,M，共M-1个信号能量比组成能量比矢量ER，即

ER＝[er₂,er₃,…,er_M]^T

B、将信号能量比计算模块A计算得到的信号能量比以及对应的指纹采集点坐标输入到数据库建立模块，并将其存储起来用于建立数据库，数据库建立模块只需离线执行一次，GMM训练模块输出的GMM参数输出至声源坐标计算模块；所述的数据库建立模块的工作方法包括以下步骤：

将声源移动到事先确定好的各个指纹采集点处，为保证有足够多的训练用矢量，降低指纹参数的随机误差，需在每个采集点多次计算输入信号能量比，并将能量比矢量与对应采集点坐标依次记录下来；设置一个平面区域，将每个能量比矢量ER_i与对应采集点坐标SL_i组成联合矢量，存储至数据库，即

$>x_i=\left(\begin{array}{c}{\mathrm{ER}}_i\\ {\mathrm{SL}}_i\end{array}\right)---\left(6\right)$>

其中，i=1，2…，N，N表示数据库中共有N组数据；

C、将数据库中的数据输入到GMM训练模块，将其作为训练矢量，训练GMM，输出GMM参数至声源坐标计算模块，在数据库不变的情况下，GMM训练模块只需离线执行一次；所述的GMM训练模块的工作过程如下：

GMM用多个高斯分布函数的组合逼近任意形状的概率分布；假设GMM中具有K个高斯混合分量，那么对于D维的观测矢量x_i，其联合概率分布表示为

$>p(x_i;\Theta )=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(x_i;{\mu }_k,{\Sigma }_k)---\left(7\right)$>

其中，w_k是第k个高斯分量的权值，满足f_k(x_i;μ_k,Σ_k)为D维高斯概率分布函数，表达式为

$>f_k(x_i;{\mu }_k,{\Sigma }_k)=\frac{1}{{\left(2\pi \right)}^{D/2}{\left|{\Sigma }_k\right|}^{1/2}}\exp [-\frac{1}{2}{(x_i-{\mu }_k)}^T{\Sigma }_k^{-1}(x_i-{\mu }_k)]---\left(8\right)$>

其中，μ_k是第k个高斯分量的均值矢量，Σ_k是第k个高斯分量的协方差矩阵，而表示全部高斯混合分量的参数，也就是整个GMM的参数；

对于含有N组观测数据的数据集X＝[x₁,x₂,…,x_N]，其对数似然函数表示为

$>{\log }_{e}p(X;\Theta )=\underset{i=i}{\overset{N}{\Sigma }}{\log }_{e}p(x_i;\Theta )---\left(9\right)$>

用数据集X来训练GMM，就是找到一组参数Θ，使得对数似然函数log_ep(X;Θ)的期望最大，求解这一组参数采用的是期望最大化算法，即EM算法；EM算法是一种迭代算法，主要包括两个步骤：求期望的EXPECT步骤和最大化的MAXIMIZE步骤；对于GMM，其具体过程如下：

C1、采用K均值聚类方法确定参数初始值

C2、求期望步骤，即EXPECT步骤：求训练数据集中第i个数据在第k个高斯分量状态下的概率

$>p(k_i=k;x_i,{\Theta }^{\left(0\right)})=\frac{w_k^{\left(0\right)}{f}_k(x_i;{\mu }_k^{\left(0\right)},{\Sigma }_k^{\left(0\right)})}{\underset{j=1}{\overset{K}{\Sigma }}{w}_j^{\left(0\right)}{f}_j(x_i;{\mu }_j^{\left(0\right)},{\Sigma }_j^{\left(0\right)})}---\left(10\right)$>

C3、最大化步骤，即MAXIMAZE步骤：求使log_ep的期望取最大值的Θ：各混合分量的均值矢量、方差以及权值，即

$>{\mu }_k^{\left(1\right)}=\frac{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})x_i}{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})}---\left(11\right)$>

$>{\Sigma }_k^{\left(1\right)}=\frac{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})(x_i-{\mu }_k^{\left(1\right)}){(x_i-{\mu }_k^{\left(1\right)})}^T}{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})}---\left(12\right)$>

$>w_k^{\left(0\right)}=\frac{1}{N}\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})---\left(13\right)$>

其中，上角标数字表示当前迭代次数；

C4、重复步骤C2和C3，直到达到收敛条件，迭代停止；

由此估计得到GMM参数从而建立GMM；

D、在线定位阶段，声源处于某个未知位置，将采集到的信号输入到信号能量比计算模块B，进行信号能量比计算，输出计算结果至声源坐标计算模块，信号能量比计算每帧执行一次；所述的信号能量比计算模块B的工作过程与信号能量比计算模块A的工作过程相同；

E、将GMM参数与信号能量比计算模块B输出的信号能量比输入到声源坐标计算模块，进行声源坐标的计算，输出声源坐标；所述的声源坐标计算模块的工作过程如下：

高斯混合回归是一种多变量非线性回归模型，其在对数据的联合概率密度构造GMM之后，导出条件概率密度和回归函数；对于联合矢量$>\left(\begin{array}{c}{\mathrm{ER}}_i\\ {\mathrm{SL}}_i\end{array}\right),$>由训练GMM模块训练所得GMM参数可构成其联合概率密度：

$>p(\mathrm{ER},\mathrm{SL})=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER},\mathrm{SL};{\mu }_k,{\Sigma }_k)$>

$>=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k\left(\mathrm{SL}\right|\mathrm{ER};m_k\left(\mathrm{ER}\right),{\sigma }_k^2)f_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})$>    (14)

式中

$>\underset{k=1}{\overset{K}{\Sigma }}{w}_k=1,{\mu }_k=\left(\begin{array}{c}{\mu }_{\mathrm{kE}}\\ {\mu }_{\mathrm{kS}}\end{array}\right),{\Sigma }_k=\left(\begin{array}{cc}{\Sigma }_{\mathrm{kEE}}& {\Sigma }_{\mathrm{kES}}\\ {\Sigma }_{\mathrm{kSE}}& {\Sigma }_{\mathrm{kSS}}\end{array}\right)---\left(15\right)$>

$>m_k\left(\mathrm{ER}\right)={\mu }_{\mathrm{kS}}+{\Sigma }_{\mathrm{kSE}}{\Sigma }_{\mathrm{kEE}}^{-1}(\mathrm{ER}-{\mu }_{\mathrm{kE}})---\left(16\right)$>

$>{\sigma }_k^2={\Sigma }_{\mathrm{kSS}}-{\Sigma }_{\mathrm{kSE}}{\Sigma }_{\mathrm{kEE}}^{-1}{\Sigma }_{\mathrm{kES}}---\left(17\right)$>

其中，ER表示能量比矢量，SL表示与之对应的采集点坐标，ER|SL表示在SL采集点处采集得能量比为ER的情况；μ_kE与μ_kS分别表示ER和SL在第k个高斯分量上的均值矢量，Σ_kEE和Σ_kSS分别表示ER和SL在第k个高斯分量上的自相关矩阵；Σ_kES和Σ_kSE分别表示ER和SL、SL和ER在第k个高斯分量上的协方差矩阵；m_k(ER)表示ER的期望，表示方差；

计算边缘密度

$>p\left(\mathrm{ER}\right)=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})---\left(18\right)$>

得到条件概率密度

$>p\left(\mathrm{SL}\right|\mathrm{ER})=\frac{p(\mathrm{ER},\mathrm{SL})}{p\left(\mathrm{ER}\right)}=\underset{k=1}{\overset{K}{\Sigma }}\frac{w_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{f}_k\left(\mathrm{SL}\right|\mathrm{ER};m_k\left(\mathrm{ER}\right),{\sigma }_k^2)---\left(19\right)$>

根据式(19)，在已知观测矢量ER′的情况下，计算其条件概率期望m(ER′)，得到声源坐标估计值即

$>\hat{S}L=m\left({\mathrm{ER}}^{\prime }\right)=E\left[\mathrm{SL}\right|\mathrm{ER}={\mathrm{ER}}^{\prime }]=\underset{k=1}{\overset{K}{\Sigma }}\frac{w_k{f}_k(E{R}^{\prime };{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(E{R}^{\prime };{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{m}_k\left({\mathrm{ER}}^{\prime }\right)---\left(20\right).$>

本发明所述的帧长对应的等效时间q的最佳值为0.2；所述的门限因子a的最佳值为a=4；所述的平滑因子b的最佳值为b=0.05。

本发明带来的有益效果如下：

1、本发明采用麦克风接收信号能量作为主要参数，它仅与该麦克风和声源间的距离有关，故无需考虑阵列几何结构与尺寸信息，也无需考虑声音信号传输延迟时间、平稳性、宽带或窄带等因素。

2、本发明采用声音信号能量比作为位置指纹参数，并进行有效信号输入检测以及背景噪声能量估计，从而消除背景噪声能量对真实信号能量比值的影响。

3、本发明用高斯混合回归(简称：GMR)实现信号能量比矢量到声源坐标的非线性映射，具有较高的定位精度与较快的映射速度。仿真测试结果显示，在本仿真条件下，本发明定位误差在0.5m之内的概率可达60%，在1m之内的概率达96%；定位速度与传统的基于位置指纹与概率的定位方法相当。

附图说明

本发明共有附图4张，其中：

图1是基于位置指纹与GMR的声源定位系统功能框图。

图2是信号能量比计算流程图。

图3是分布式麦克风阵列、指纹采集点以及声源分布示意图。

图4是GMR定位误差的累计分布示意图。

图中：1、信号能量比计算模块A，2、数据库建立模块，3、GMM训练模块，4、声源坐标计算模块，5、信号能量比计算模块B。

具体实施方式

下面结合附图对本发明进行进一步地描述。图1是基于位置指纹与GMR的声源定位系统功能框图。具体步骤如下：

A、在已部署好的麦克风阵列工作场地，将声源移动至事先确定好的各个指纹采集点处，将采集到的信号输入到信号能量比计算模块A，进行信号能量比计算，输出计算结果至数据库建立模块；所述的信号能量比计算对每个指纹采集点处采集的每帧信号进行一次；信号能量比计算流程如图2所示。所述的信号能量比计算模块A的工作方法包括以下步骤：

首先检测输入信号的有效性，定义分布式麦克风阵列中麦克风数目为M，第1个麦克风接收到输入为s₁(n)，当前帧在原始信号序列中的起始位置为start，信号帧长为L，则第1个麦克风的当前帧短时对数能量E₁为

$>E_1=\underset{n=\mathrm{start}}{\overset{\mathrm{start}+L-1}{\Sigma }}{s}_1^2\left(n\right)$>

将E₁与信号/噪声判定阈值E_thr进行比较，若E₁<E_thr，判定当前帧为噪声帧，否则，判定当前帧为信号帧；若当前帧为噪声帧，则估计背景噪声能量并更新阈值E_thr;更新阈值E_thr的方法为

if  E_min＜E₁

  E_min＝E₁

  E_thr＝aE_min

end if

其中，E_min表示记录的第1个麦克风最小短时能量，a为一个门限因子，取值为a=4；

定义第m个麦克风接收到输入为s_m(n)，若当前帧为噪声帧，则第m个麦克风的背景噪声E_m估计方法为

$>E_m=\underset{n=\mathrm{start}}{\overset{\mathrm{start}+N-1}{\Sigma }}{s}_m^2\left(n\right),m=\mathrm{1,2},···,M$>

ε_m,curr＝(1-b)·ε_m,old+b·E_m

其中，ε_m,curr为第m个麦克风当前帧背景噪声能量估计值，ε_m,old为第m个麦克风前一次背景噪声能量估计值，b为平滑因子，表示平滑过程中当前帧能量E_m所占比重，取值为b=0.05；

若当前帧为信号帧，去除第1个和第m个麦克风背景噪声能量的影响，之后以第1个麦克风信号能量作参考，计算第m个麦克风与第1个麦克风的修正后的信号能量比er_m，

$>{\mathrm{er}}_m=\frac{E_m-{ϵ}_{m,\mathrm{curr}}}{E_1-{ϵ}_1}$>

然后将er_m，m=2,3,…,M，共M-1个信号能量比组成能量比矢量ER，即

ER＝[er₂,er₃,…,er_M]^T

B、将信号能量比计算模块A计算得到的信号能量比以及对应指纹采集点坐标输入到数据库建立模块，并将其存储起来用于建立数据库，数据库建立模块只需离线执行一次，GMM训练模块输出的GMM参数输出至声源坐标计算模块；所述的数据库建立模块的工作方法包括以下步骤：

将声源移动到事先确定好的各个指纹采集点处，为保证有足够多的训练用矢量，降低指纹参数的随机误差，需在在每个采集点多次计算输入信号能量比，并将能量比矢量与对应采集点坐标依次记录下来；设置一个平面区域，将每个能量比矢量ER_i与对应采集点坐标SL_i组成联合矢量，存储至数据库，即

$>x_i=\left(\begin{array}{c}{\mathrm{ER}}_i\\ {\mathrm{SL}}_i\end{array}\right)---\left(6\right)$>

其中，i=1，2…，N，N表示数据库中共有N组数据。

C、将数据库中的数据输入到GMM训练模块，将其作为训练矢量，训练GMM，输出GMM参数至声源坐标计算模块，在数据库不变的情况下，GMM训练模块只需离线执行一次；所述的GMM训练模块的工作过程如下：

GMM用多个高斯分布函数的组合逼近任意形状的概率分布；假设GMM中具有K个高斯混合分量，那么对于D维的观测矢量x_i，其联合概率分布表示为

$>p(x_i;\Theta )=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(x_i;{\mu }_k,{\Sigma }_k)---\left(7\right)$>

其中，w_k是第k个高斯分量的权值，满足f_k(x_i;μ_k,Σ_k)为D维高斯概率分布函数，表达式为

$>f_k(x_i;{\mu }_k,{\Sigma }_k)=\frac{1}{{\left(2\pi \right)}^{D/2}{\left|{\Sigma }_k\right|}^{1/2}}\exp [-\frac{1}{2}{(x_i-{\mu }_k)}^T{\Sigma }_k^{-1}(x_i-{\mu }_k)]---\left(8\right)$>

其中，μ_k是第k个高斯分量的均值矢量，Σ_k是第k个高斯分量的协方差矩阵，而表示了全部高斯混合分量的参数，也就是整个GMM的参数；

对于含有N组观测数据的数据集X＝[x₁,x₂,…,x_N]，其对数似然函数表示为

$>{\log }_{e}p(X;\Theta )=\underset{i=i}{\overset{N}{\Sigma }}{\log }_{e}p(x_i;\Theta )---\left(9\right)$>

用数据集X来训练GMM，就是找到一组参数Θ，使得对数似然函数log_ep(X;Θ)的期望最大，求解这一组参数采用的是期望最大化算法，即EM算法；EM算法是一种迭代算法，主要包括两个步骤：求期望的EXPECT步骤和最大化的MAXIMIZE步骤；对于GMM，其具体过程如下：

C1、采用k均值聚类方法确定参数初始值

C2、EXPECT步骤：求训练数据集中第i个数据在第k个高斯分量状态下的概率

$>p(k_i=k;x_i,{\Theta }^{\left(0\right)})=\frac{w_k^{\left(0\right)}{f}_k(x_i;{\mu }_k^{\left(0\right)},{\Sigma }_k^{\left(0\right)})}{\underset{j=1}{\overset{K}{\Sigma }}{w}_j^{\left(0\right)}{f}_j(x_i;{\mu }_j^{\left(0\right)},{\Sigma }_j^{\left(0\right)})}---\left(10\right)$>

C3、MAXIMIZE步骤：求使log_ep的期望取最大值的Θ：各混合分量的均值矢量、方差以及权值，即

$>{\mu }_k^{\left(1\right)}=\frac{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})x_i}{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})}---\left(11\right)$>

$>{\Sigma }_k^{\left(1\right)}=\frac{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})(x_i-{\mu }_k^{\left(1\right)}){(x_i-{\mu }_k^{\left(1\right)})}^T}{\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})}---\left(12\right)$>

$>w_k^{\left(0\right)}=\frac{1}{N}\underset{i=1}{\overset{N}{\Sigma }}p(k_i=k;x_i,{\Theta }^{\left(0\right)})---\left(13\right)$>

其中，上角标数字表示当前迭代次数。

C4、重复步骤C2和C3，直到达到收敛条件，迭代停止；

由此就估计得到GMM参数从而建立GMM；

D、在线定位阶段，声源处于某个未知位置，将采集到的信号输入到信号能量比计算模块B，进行信号能量比计算，输出计算结果至声源坐标计算模块，信号能量比计算每帧执行一次；所述的信号能量比计算模块B的工作过程与信号能量比计算模块A的工作过程相同；

E、将GMM参数与信号能量比计算模块B输出的信号能量比输入到声源坐标计算模块，进行声源坐标的计算，输出声源坐标；所述的声源坐标计算模块的工作过程如下：

高斯混合回归GMR是一种多变量非线性回归模型，其是在对数据的联合概率密度构造GMM之后，导出条件概率密度和回归函数；对于联合矢量$>\left(\begin{array}{c}{\mathrm{ER}}_i\\ {\mathrm{SL}}_i\end{array}\right),$>由训练GMM模块训练所得GMM参数可构成其联合概率密度：

$>p(\mathrm{ER},\mathrm{SL})=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER},\mathrm{SL};{\mu }_k,{\Sigma }_k)$>

$>=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k\left(\mathrm{SL}\right|\mathrm{ER};m_k\left(\mathrm{ER}\right),{\sigma }_k^2)f_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})$>  (14)

式中，

$>\underset{k=1}{\overset{K}{\Sigma }}{w}_k=1,{\mu }_k=\left(\begin{array}{c}{\mu }_{\mathrm{kE}}\\ {\mu }_{\mathrm{kS}}\end{array}\right),{\Sigma }_k=\left(\begin{array}{cc}{\Sigma }_{\mathrm{kEE}}& {\Sigma }_{\mathrm{kES}}\\ {\Sigma }_{\mathrm{kSE}}& {\Sigma }_{\mathrm{kSS}}\end{array}\right)---\left(15\right)$>

$>m_k\left(\mathrm{ER}\right)={\mu }_{\mathrm{kS}}+{\Sigma }_{\mathrm{kSE}}{\Sigma }_{\mathrm{kEE}}^{-1}(\mathrm{ER}-{\mu }_{\mathrm{kE}})---\left(16\right)$>

$>{\sigma }_k^2={\Sigma }_{\mathrm{kSS}}-{\Sigma }_{\mathrm{kSE}}{\Sigma }_{\mathrm{kEE}}^{-1}{\Sigma }_{\mathrm{kES}}---\left(17\right)$>

其中，ER表示能量比矢量，SL表示与之对应的采集点坐标，ER|SL表示在SL采集点处采集得能量比为ER的情况；μ_kE与μ_kS分别表示ER和SL在第k个高斯分量上的均值矢量，Σ_kEE和Σ_kSS分别表示ER和SL在第k个高斯分量上的自相关矩阵；Σ_kES和Σ_kSE分别表示ER和SL、SL和ER在第k个高斯分量上的协方差矩阵；m_k(ER)表示ER的期望，表示方差。

计算边缘密度

$>p\left(\mathrm{ER}\right)=\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})---\left(18\right)$>

得到条件概率密度

$>p\left(\mathrm{SL}\right|\mathrm{ER})=\frac{p(\mathrm{ER},\mathrm{SL})}{p\left(\mathrm{ER}\right)}=\underset{k=1}{\overset{K}{\Sigma }}\frac{w_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(\mathrm{ER};{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{f}_k\left(\mathrm{SL}\right|\mathrm{ER};m_k\left(\mathrm{ER}\right),{\sigma }_k^2)---\left(19\right)$>

根据式(19)，在已知观测矢量ER′的情况下，计算其条件概率期望m(ER′)，得到声源坐标估计值即

$>\hat{S}L=m\left({\mathrm{ER}}^{\prime }\right)=E\left[\mathrm{SL}\right|\mathrm{ER}={\mathrm{ER}}^{\prime }]=\underset{k=1}{\overset{K}{\Sigma }}\frac{w_k{f}_k(E{R}^{\prime };{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{\underset{k=1}{\overset{K}{\Sigma }}{w}_k{f}_k(E{R}^{\prime };{\mu }_{\mathrm{kE}},{\Sigma }_{\mathrm{kEE}})}{m}_k\left({\mathrm{ER}}^{\prime }\right)---\left(20\right).$>

下面通过测试实验来说明本发明带来的有益效果。

1、测试环境与条件

包含7个麦克风的分布式麦克风阵列工作于区域为一个10m×10m的矩形，矩形的四个顶点坐标分别为(0,0)m、(0,10)m、(10,0)m和(10,10)m。在区域中设置100个指纹采集点。选取不同的声音信号作为声源发出的信号，包括语音、乐器音、歌曲、乐曲、鸣笛声、按键音以及白噪声，在区域中设置20个声源位置。分布式麦克风阵列、指纹采集点以及声源的分布如图3所示。分别用确定性方法、概率法、神经网络(Artificial Neural Network，ANN)方法以及基于GMR的方法对声源位置依次进行估计，为满足需要，共建立四个数据库(如表1所示)。假设在训练、定位之前，各麦克风的增益按MADHU N,MARTIN R.发表的《Low-complexity,robust algorithm for sensor anomaly detection andself-calibration of microphone arrays.IET Signal Processing》[2011,5(1):97-103]的方法已校准好。采集声音信号时的采样频率为16KHz。

表1数据库说明

2、基于位置指纹与GMR的声源定位结果

分别基于数据库DB1、DB2与DB3，使用GMR方法实现声源定位。图4给出基于三种数据库、高斯分量为16个和24个情况下，GMR方法的误差累计分布。由图中可见，基于数据库DB1时，定位误差小于1m的概率低于使用数据库DB2和DB3情况，但除了一种情况外，其余5种情况的定位误差小于1m的概率均可达90%以上。

3、GMR声源定位与已有位置指纹声源定位效果比较

将GMR方法与已有的基于确定性方法、概率法、神经网络的位置指纹定位效果进行比较，结果如表2所示，其中定位所用时间是在主频为2.26GHz的HP笔记本电脑上测得的处理长度为140秒的信号数据时，程序执行所用时间，存储数据类型为双精度浮点型（即C语言中的double型）。

从表2可以看出，四类方法均能实现声源的有效定位。从定位误差角度评价，概率法定位误差最大，确定性方法居中，而BP神经网络法和GMR方法定位误差较小，定位误差小于0.5m的情况可达60%以上，明显优于其他方法。从定位速度角度评价，确定性方法计算量最小、速度最快，概率法和GMR方法定位速度相当，耗时约为确定性方法的3倍，而神经网络法定位速度非常慢，耗时为确定性方法的15倍或30倍。从存储数据量角度评价，广义回归神经网络(Generalized Regression Neural Network，GRNN)方法存储数据量较大（这是由于输入训练样本较多），确定性方法、径向基函数神经网络(Radial Basis FunctionNeural Network，RBF-NN)方法与反向传播神经网络(Back Propagation NeuralNetwork，BP-NN)方法存储量较小，概率法和GMR方法存储量居中。

综上所述，GMR方法能够在保证较快定位速度的同时，具有较高的定位精度，整体性能优于其余三类方法。

表2四类定位方法性能比较