CN103837858A - 一种用于平面阵列的远场波达角估计方法及*** - Google Patents

Abstract

本发明涉及一种用于平面阵列的远场波达角估计方法及***，所述的方法包含：步骤101）将平面阵列作为麦克风对组成的线性子阵列的几何组合，且每一个线性子阵列决定一个子波达角；步骤102）在假定波达方向x已知的情况下：采用全局波达角来计算各子阵列的期望子波达角θ_i；通过各子阵列的时间差计算子阵列的估计子波达角步骤103）基于估计子波达角和期望子波达角构造代价函数为：步骤104）将代价函数收敛时的波达角作为最终确定的波达角的值，完成声源定位。所述τ_i采用如下策略获得：提取麦克风阵列的输出数字信号，对每一帧的数字化声音信号做加窗预处理，傅立叶变换并且在频域白化信号；计算预处理后信号的交叉相关，求取所有麦克风对之间的时间延迟τ_i。

Description

一种用于平面阵列的远场波达角估计方法及***

技术领域

本发明涉及宽带信号、远场的麦克风阵列声源定位技术领域，具体的说，本发明涉及一种用于平面阵列的远场波达角估计方法。

背景技术

远场波达角估计在麦克风阵列的应用中占居重要位置，它可用于远程会议，为麦克风阵列指示波束聚焦的方向，为会议摄像头提供指向信息。它能指示声源目标所在的空间方位，为后续的信息采集与处理提供重要的空间信息。

传统的宽带波达角估计算法分为两类，第一类是以基于相位变换的联合可控响应功率(SRP-PHAT)算法为代表的空间谱估计算法，另外一类是以基于相位变换的广义交叉相关（GCC-PHAT）为代表的时间差算法。前者以鲁棒性著称，具有较好的抗噪声、抗混响能力，但计算复杂度高；后者以计算效率著称，但鲁棒性差，容易遭受噪声和混响的干扰。这两类算法的代价函数存在两个共同的问题。其中，波达角是指声源入射的方向，相对于平面阵列的方位角和仰角。

其一，他们的代价函数通常拥有多个极值，而且根据最优准则，无法导出波达角的最优解，因此，不得不采用计算复杂度较高的格点搜索方法，寻求一个最优估计。格点搜索是导致这些算法高复杂度的主要原因。尽管有些算法提出了一些格点搜索的优化方法，但优化后的格点搜索的计算复杂度仍然不可忽略。

其二，基于交叉相关的时间差估计方法容易受到噪声和混响的影响，由于各麦克风所处于的位置不同，他们受到噪声和混响干扰的程度不同，计算精度差异较大。此外，时间差计算的精度，也和麦克风与声源的相对位置有关，在垂直于麦克风对的方向，时间延迟估计较准，然而在接近于平行方向，时间延迟估计精度较差。因此，我们必须在代价函数中考虑这种差异。不幸的是，大多数算法往往忽略了这一点，它们平等对待每一时间延迟信息。对于某一时间差的大幅度扰动，往往导致较大的波达角估计误差。

以上两个问题根源于代价函数的定义，因此，解决上述问题必须从代价函数的定义入手。

发明内容

本发明目的在于，为克服现有技术的上述缺陷，本发明提供了一种用于平面阵列的远场波达角估计方法及***。

为实现上述目的，本发明提供了一种用于平面阵列的远场波达角估计方法，所述的方法包含：

步骤101）将平面阵列作为麦克风对组成的线性子阵列的几何组合，且每一个线性子阵列决定一个子波达角；

步骤102）在假定波达方向x已知的情况下：

首先，采用全局波达角来计算各子阵列的期望子波达角θ_i；然后采用通过各子阵列的时间差计算子阵列的估计子波达角

步骤103）基于针对各子阵列估计波达角和期望波达角构造代价函数为：

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

其中，M为总的麦克风对的数目，g_i为第i对麦克风连线的方向，τ_i为第i对麦克风估计的时间延迟，d_i为第i对麦克风的间距，x为波达方向，w_i为给予第i对麦克风的权重，c为声波在空气中的传递速度；

步骤104）将上述代价函数收敛时的波达角的取值作为最终确定的波达角的值，完成声源定位。

上述期望子波达角θ_i的计算公式为：

θ_i=arccos(cτ_i/d_i)

所述估计子波达角的计算公式为：

${\widehat{\mathrm{\θ}}}_i=\arccos \left(x^T{g}_i\right).$

上述τ_i采用如下策略获得：

步骤1：提取麦克风阵列的输出数字信号，对每一帧的数字化声音信号做加窗预处理，进行傅立叶变换，并且在频域白化信号；

步骤2：计算预处理后信号的交叉相关，对于来自每对麦克风的傅立叶谱进行点乘，求取所有麦克风对之间的时间延迟τ_i。

上述步骤104）进一步包含：

步骤104-1）假定权重系数相等，即w_i=1/M，初始化指向x=[1,0,0]；对于所有麦克风对，计算它们的连线表示的单位方向矢量，g_i＝[g_i，1,g_i，2，0]^T，其中，g_i′=[g_i，1，g_i,2]^T；0表示所有的麦克风位于同一个平面内；

步骤104-2）将所有的g_i′和τ_i代入公式 $\left(\begin{array}{c}{\hat{x}}_1\\ {\hat{x}}_2\end{array}\right)={\left[{\mathrm{\Σ}}_{i=1}^M{w}_i{g}_i^{\′}{g}_i^{\′T}\right]}^{-1}{\mathrm{\Σ}}_{i=1}^M{\mathrm{cw}}_i{\mathrm{\τ}}_i{g}_i^{\′}/d_i$ 和

中，得到波达方向的单位矢量

步骤104-3）如果

无限接近于1，则结束迭代，

即为最终的波达方向，并跳转到步骤104-5）；否则并跳转到步骤104-4）继续执行迭代；

步骤104-4）根据波达方向

求取新的权重系数， $w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2),$ 其中 ${\mathrm{\δ}}_i=\arccos \left({\hat{x}}^T{g}_i\right)-\arccos ({\mathrm{c\τ}}_i/d_i),$ ${\mathrm{\σ}}^2={\mathrm{\Σ}}_i{\mathrm{\δ}}_i^2/M;$ 跳转到步骤104-2）；

步骤104-5）将波达方向矢量

转化为方位角和仰角，完成声源定位。

上述步骤104-5）采用反正切三角函数将波达方向矢量

转化为方位角和仰角，其中，方位角为

当

时，仰角表示为 $\arccos (\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2}/\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2+{{\hat{x}}_3}^2})$ 和，当

时，仰角表示为

基于上述方法本发明提供了一种用于平面阵列的远场波达角估计***，所述的***包含：

设置模块，用于将平面阵列作为麦克风对组成的线性子阵列的几何组合，且每一个线性子阵列决定一个子波达角；

波达角获取模块，用于在给定一个未知的DOA情况下采用如下策略获取波达角：

首先，采用全局波达角来计算各子阵列的期望子波达角

然后采用通过测量子阵列的时间差计算各子阵列的估计子波达角θ_i；

代价函数生成模块，用于基于针对各子阵列的上述估计子波达角和期望子波达角构造代价函数为：

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

M为总的麦克风对的数目，g_i为第i对麦克风连线的方向，τ_i为第i对麦克风估计的时间延迟，d_i为第i对麦克风的间距，x为声源入射方位。w_i为给予第i对麦克风的权重；和

波达角确定输出模块，用于将上述代价函数收敛时的波达角的取值作为最终确定的波达角的值，完成声源定位。

上述期望子波达角θ_i的计算公式为：

θ_i=arccos(cτ_i/d_i)

所述估计子波达角

的计算公式为：

${\widehat{\mathrm{\θ}}}_i=\arccos \left(x^T{g}_i\right).$

当获取τ_i时采用如下模块：

预处理模块，提取麦克风阵列的输出数字信号，对每一帧的数字化声音信号做加窗预处理，进行傅立叶变换，并且在频域白化信号；

计算预处理后信号的交叉相关，对于来自每对麦克风的傅立叶谱进行点乘，求取所有麦克风对之间的时间延迟τ_i。

上述波达角确定输出模块进一步包含：

初始化子模块，用于假定权重系数相等，即w_i=1/M，初始化指向x=[1,0,0]；对于所有麦克风对，计算它们的连线表示的单位方向矢量，g_i′=[g_i，1,g_i,2]^T；

第一处理子模块，用于将所有的g_i′和τ_i代入公式 $\left(\begin{array}{c}{\hat{x}}_1\\ {\hat{x}}_2\end{array}\right)={\left[{\mathrm{\Σ}}_{i=1}^M{w}_i{g}_i^{\′}{g}_i^{\′T}\right]}^{-1}{\mathrm{\Σ}}_{i=1}^M{\mathrm{cw}}_i{\mathrm{\τ}}_i{g}_i^{\′}/d_i$ 和中，得到波达方向的单位矢量

判决决策子模块，用于判断如果其中η表示一个小于且接近于1的常数，则结束迭代，

即为最终的波达方向，并将结果输出给波达角输出子模块；否则

并将该设置结果输入迭代处理子模块继续执行迭代；

迭代处理子模块，用于根据波达方向求取新的权重系数， $w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2),$ 其中 ${\mathrm{\δ}}_i=\arccos \left({\hat{x}}^T{g}_i\right)-\arccos ({\mathrm{c\τ}}_i/d_i),$ ${\mathrm{\σ}}^2={\mathrm{\Σ}}_i{\mathrm{\δ}}_i^2/M;$ 跳转到第一处理子模块；和

波达角输出子模块，用于将波达方向矢量

转化为方位角和仰角，完成声源定位。

上述波达角输出子模块采用反正切三角函数将波达方向矢量

转化为方位角和仰角，其中，方位角为当时，仰角表示为 $\arccos (\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2}/\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2+{{\hat{x}}_3}^2})$ 和，当

时，仰角表示为

与现有技术相比本发明的技术优势在于：

1．计算速度快。传统的宽带信号波达角估计方法中，由于代价函数有多个极值，不得不采取全局搜索的方法，导致很高的计算复杂度。本算法通过定义一个凹函数，直接获取波达角的解析解，避免了复杂的搜索过程。因而具有很好的计算效率。它的计算速度是经典的SRP-PHAT的30倍，是基于时间延迟估计方法的5倍。

2．鲁棒性好。传统的基于时间延迟求解波达角的类似方法中，每一对相关函数或者时间延迟被赋予同等的权重，那些不可靠的延迟或相关函数导致较大的估计误差。本方法中，通过赋予不同的权重，避免不可靠的时间延迟信息对最终估计结果的干扰。因此它的鲁棒性好，它的性能接近于以鲁棒型著称的SRP-PHAT，它的鲁棒性大幅度超出以计算速度著称的GCC-PHAT。

附图说明

图1估计波达角的流程图

图2用于算法评价的声学仿真环境

图3-a为设定房间的混响时间为200毫秒，得到的“错误率——误差门限”曲线；

图3-b为设定混响时间为300毫秒得到的“错误率——误差门限”曲线；

图3-c为设定混响时间为200毫秒，外加10dB的白噪声的“错误率——误差门限”曲线。

具体实施方式

下面结合附图对本发明的内容进行详细阐述。

本发明提出了一种快速的波达角估计方法，下面以宽带语音声源为例，说明如何使用该算法进行定位。

本发明涉及一种快速而鲁棒的声源定位算法，该算法把平面阵列看作麦克风对组成的线性子阵列的几何组合，每一个线性子阵列决定一个子波达角。传统的波达角由方位角和仰角描述，为了方便起见，这里波达角表示为一个三维的单位向量x=[x₁,x₂,x₃]^T，其中^T表示矩阵的转置。假定第i个子阵列所在的位置表示为g=[g_i，1,g_i，2,0]^T，由于所有的麦克风都在一个平面内，所以第三个维度为零表示。那么假定波达方向x已知的情况下，期望的子阵列波达角可以表示为：

θ_i=arccos(cτ_i/d_i)

通过广义交叉相关可以估计第i个子阵列的时间差，从而可以计算该子阵列的估计子波达角：

${\widehat{\mathrm{\θ}}}_i=\arccos \left(x^T{g}_i\right)$

其中，c表示空气在声波中的传播速度，d_i表示子阵列的孔径。在理想的声学条件下，

应当无限接近θ_i。在实际环境中，由于外界的干扰，

通常偏离θ_i。本发明采用如下的代价函数描述这种偏离关系，

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

这里w_i表示第i个子阵列的权重系数，反映他的可靠性。最优波达角估计通过最小化该权重均方误差获得，即

$\underset{x}{\min }\mathrm{\ϵ}\left(x\right)\mathrm{subjectto}:x^{T}x=1$

这是一个基于Kuhn-Tucker条件的最小化约束问题，对应的拉格朗日梯度方程可以表示为

L(x,μ)=ε(x)+μ(x^Tx-1)

从

我们可以获得DOA的最优估计

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

这里g_i′=[g_i，1，g_i,2]^T，w_i通过θ_i与

之间的差值求解，

δ_i=arccos(x^Tg_i)-arccos(cτ_i/d_i)

假定δ_i遵从零均值的高斯分布，它的方差可以表示为

那么归一化的权重系数可以表示为：

$w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)$

w_i可以通过迭代的方式获取。首先假设所有子阵列的权重系数相同，从而求取一个最优波达角的估计，用此估计求取新的权重系数，用新的权重系数求取新的波达角估计。此迭代不断循环，直至波达角收敛。

如图1所示，包括下列步骤：

步骤1：提取麦克风输出的数字信号，对这一帧的数字化声音信号做预处理，对数字信号进行加窗处理，设每帧长度为F点。先补零到N点，N≥F，N=2^j，j为整数，且j≥8，进行N点离散傅里叶变换，得到离散谱

其中，

表示缓存中第个麦克风采集信号的第n个采样点，

表示缓存中第

个采集信号傅里叶变换值（k=0,1,...,N-1）；然后白化幅度谱，

所述的预处理包括加窗、或/和预加重；所述的加窗函数采用汉宁窗或哈宁窗。

步骤2：进行交叉相关处理，对于来自每对麦克风的白化傅立叶幅度谱进行点乘。对于由第p个和第q个麦克风组成的第i个对，求其反傅利叶变换其实部real(x_i,n)表示两个信号的相关函数，对real(x_i,n)循环移动N/2个采样点，得到相关函数序列z_i，k，根据麦克风的间距d_i计算最大时间延迟τ_max=fix(f_sd_i/c)，其中f_s表示信号的采样率，c表示声速，fix表示取整数运算。在从k₁=N/2-τ_max到k₂=N/2+τ_max的范围内，在

序列中搜索最大值，该值对应的时间索引第p个和第q个麦克风信号间的时间延迟。以此类推，求取所有麦克风对之间的时间延迟τ_i。

步骤3：假定权重系数相等，即w_i=1/M，初始化指向对于所有麦克风对，计算它们的连线表示的单位方向矢量，g_i′=[g_i，1,g_i，2]^T。

步骤4：将所有的g_i′和τ_i代入公式 $\left(\begin{array}{c}{\hat{x}}_1\\ {\hat{x}}_2\end{array}\right)={\left[{\mathrm{\Σ}}_{i=1}^M{w}_i{g}_i^{\′}{g}_i^{\′T}\right]}^{-1}{\mathrm{\Σ}}_{i=1}^M{\mathrm{cw}}_i{\mathrm{\τ}}_i{g}_i^{\′}/d_i$ 和

中，得到波达方向的单位矢量

步骤5：如果

无限接近于1，则结束迭代，

即为最终的波达方向，并跳转到步骤7；否则

并跳转到步骤6继续执行迭代。

步骤6：根据波达方向

求取新的权重系数， $w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2),$ 其中 ${\mathrm{\δ}}_i=\arccos \left({\hat{x}}^T{g}_i\right)-\arccos ({\mathrm{c\τ}}_i/d_i),$ ${\mathrm{\σ}}^2={\mathrm{\Σ}}_i{\mathrm{\δ}}_i^2/M,$ 跳转到步骤4。

步骤7：其中采用反正切三角函数可以将波达方向矢量

转化为方位角和仰角。其中，方位角为当

时，仰角表示为 $\arccos (\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2}/\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2+{{\hat{x}}_3}^2})$ 和，当

时，仰角表示为

基于上述实施例的算法，对声源定位算法的性能进行评估。我们采用image声学仿真工具，模拟一个普通房间的声学环境。房间中放置一个9元的麦克风平面阵列，8个麦克风均匀放置在圆周上，一个放置在圆心。话者距离麦克风中心位置大约1.8米处。房间的几何平面图如图2所示。

采用TIMIT数据库中的60秒的语音信号。算法性能的评价指标为给定门限下的错误率，具体定义为：在给定一个方位角误差门限的情况下，误差大于门限的样本所占总体样本的百分比定义为错误率。图2为实例比较了本发明的权重最小均方误差算法WMMSE、不考虑权重的最小均方误差算法MMSE、经典的基于相位变换的广义交叉相关算法（GCC-PHAT）、基于相位变换的联合可控响应功率（SRP-PHAT）效果比对图。图3-a、3-b和3-c表示算法评价的数据指标对比。首先通过仿真工具，设定房间的混响时间为200毫秒，得到的“错误率——误差门限”曲线如图3-a所示；然后设定混响时间为300毫秒，如图3-b所示；最后设定混响时间为200毫秒，外加10dB的白噪声，如图3-c所示。该实验表明，算法的性能接近于SRP-PHAT，但远远超出GCC-PHAT。同时我们在普通台式电脑上对比了算法的计算效率。本算法分别为GCC-PHAT的5倍，SRP-PHAT的30倍。

总之，本发明提供了一种快速而鲁棒的声源定位算法，该算法把平面阵列看作麦克风对组成的线性子阵列的几何组合，每一个线性子阵列决定一个子波达角。该子波达角既可以通过交叉相关来估计，也可以通过给定的全局波达角来计算。在理想的声学条件下，估计的子波达角和计算得到的子波达角相等。我们采用这两个子波达角余弦差的平方作为代价函数，波达角的最优估计应该使代价函数最小。根据该准则，可以推导出最优波达角为所有麦克风对时间差的线性函数。同时，我们在该函数中引入权重系数，衡量每一子阵列的可靠度。我们在仿真环境和真实环境下测试了该算法，它的计算速度比SRP-PHAT高出30倍，与以计算效率著称的GCC-PHAT相比，甚至高出5倍。它的鲁棒性大幅超出GCC-PHAT，接近于以鲁棒性闻名的SRP-PHAT。该算法适用于低功耗硬件平台上的声源定位。

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

Claims (10)

1.一种用于平面阵列的远场波达角估计方法，所述的方法包含：

步骤101）将平面阵列作为麦克风对组成的线性子阵列的几何组合，且每一个线性子阵列决定一个子波达角；

步骤102）在假定波达方向x已知的情况下：

首先，采用全局波达角来计算各子阵列的期望子波达角θ_i；然后采用通过各子阵列的时间差计算子阵列的估计子波达角

步骤103）基于针对各子阵列估计波达角和期望波达角构造代价函数为：

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

其中，M为总的麦克风对的数目，g_i为第i对麦克风连线的方向，τ_i为第i对麦克风估计的时间延迟，d_i为第i对麦克风的间距，x为波达方向，w_i为给予第i对麦克风的权重，c为声波在空气中传递的速度；

步骤104）将上述代价函数收敛时的波达角的取值作为最终确定的波达角的值，完成声源定位。

2.根据权利要求1所述的用于平面阵列的远场波达角估计方法，其特征在于，所述期望子波达角θ_i的计算公式为：

θ_i=arccos(cτ_i/d_i)

所述估计子波达角

的计算公式为：

${\widehat{\mathrm{\θ}}}_i=\arccos \left(x^T{g}_i\right).$

3.根据权利要求1或2所述的用于平面阵列的远场波达角估计方法，其特征在于，所述τ_i采用如下策略获得：

步骤1：提取麦克风阵列的输出数字信号，对每一帧的数字化声音信号做加窗预处理，进行傅立叶变换，并且在频域白化信号；

步骤2：计算预处理后信号的交叉相关，对于来自每对麦克风的傅立叶谱进行点乘，求取所有麦克风对之间的时间延迟τ_i。

4.根据权利要求1所述的用于平面阵列的远场波达角估计方法，其特征在于，所述步骤104）进一步包含：

步骤104-1）假定权重系数相等，即w_i=1/M，初始化指向x=[1,0,0]；对于所有麦克风对，计算它们的连线表示的单位方向矢量，g_i′＝[g_i，1,g_i,2]^T，其中g_i，1和g_i,2分别表示表示g_i的第一维和第二维；

步骤104-2）将所有的g_i′和τ_i代入代价函数公式 $\left(\begin{array}{c}{\hat{x}}_1\\ {\hat{x}}_2\end{array}\right)={\left[{\mathrm{\Σ}}_{i=1}^M{w}_i{g}_i^{\′}{g}_i^{\′T}\right]}^{-1}{\mathrm{\Σ}}_{i=1}^M{\mathrm{cw}}_i{\mathrm{\τ}}_i{g}_i^{\′}/d_i$ 和

中，得到波达方向的单位矢量

其中，

和

表示

的第一维和第二维分量；

步骤104-3）如果

无限接近于1，则结束迭代，

即为最终的波达方向，并跳转到步骤104-5）；否则

并跳转到步骤104-4）继续执行迭代；

步骤104-4）根据波达方向

求取新的权重系数， $w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2),$ 其中 ${\mathrm{\δ}}_i=\arccos \left({\hat{x}}^T{g}_i\right)-\arccos ({\mathrm{c\τ}}_i/d_i),$ ${\mathrm{\σ}}^2={\mathrm{\Σ}}_i{\mathrm{\δ}}_i^2/M;$ 跳转到步骤104-2）；

步骤104-5）将波达方向矢量

转化为方位角和仰角，完成声源定位。

5.根据权利要求4所述的用于平面阵列的远场波达角估计方法，其特征在于，所述步骤104-5）采用反正切三角函数将波达方向矢量

转化为方位角和仰角，其中，方位角为

当

时，仰角表示为 $\arccos (\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2}/\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2+{{\hat{x}}_3}^2})$ 和，当时，仰角表示为

6.一种用于平面阵列的远场波达角估计***，所述的***包含：

设置模块，用于将平面阵列作为麦克风对组成的线性子阵列的几何组合，且每一个线性子阵列决定一个子波达角；

波达角获取模块，用于在给定一个未知的DOA情况下采用如下策略获取波达角：

首先，采用全局波达角来计算各子阵列的期望子波达角

然后采用通过测量子阵列的时间差计算各子阵列的估计子波达角θ_i；

代价函数生成模块，用于基于针对各子阵列的上述估计子波达角和期望子波达角构造代价函数为：

$\mathrm{\ϵ}\left(x\right)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[{\cos \mathrm{\θ}}_i-\cos {\widehat{\mathrm{\θ}}}_i]}^2$

$=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{w}_i{[x^T{g}_i-{\mathrm{c\τ}}_i/d_i]}^2$

M为总的麦克风对的数目，g_i为第i对麦克风连线的方向，τ_i为第i对麦克风估计的时间延迟，d_i为第i对麦克风的间距，x为声源入射方位。w_i为给予第i对麦克风的权重；c为声波在空气中的传递速度；和

波达角确定输出模块，用于将上述代价函数收敛时的波达角的取值作为最终确定的波达角的值，完成声源定位。

7.根据权利要求6所述的用于平面阵列的远场波达角估计***，其特征在于，所述期望子波达角θ_i的计算公式为：

θ_i=arccos(cτ_i/d_i)

所述估计子波达角

的计算公式为：

${\widehat{\mathrm{\θ}}}_i=\arccos \left(x^T{g}_i\right).$

8.根据权利要求6或7所述的用于平面阵列的远场波达角估计***，其特征在于，获取τ_i时采用如下模块：

预处理模块，提取麦克风阵列的输出数字信号，对每一帧的数字化声音信号做加窗预处理，进行傅立叶变换，并且在频域白化信号；和

计算预处理后信号的交叉相关，对于来自每对麦克风的傅立叶谱进行点乘，求取所有麦克风对之间的时间延迟τ_i。

9.根据权利要求6所述的用于平面阵列的远场波达角估计***，其特征在于，所述波达角确定输出模块进一步包含：

初始化子模块，用于假定权重系数相等，即w_i=1/M，初始化指向

对于所有麦克风对，计算它们的连线表示的单位方向矢量，g_i′=[g_i，1,g_i,2]^T；

第一处理子模块，用于将所有的g_i′和τ_i代入代价函数公式 $\left(\begin{array}{c}{\hat{x}}_1\\ {\hat{x}}_2\end{array}\right)={\left[{\mathrm{\Σ}}_{i=1}^M{w}_i{g}_i^{\′}{g}_i^{\′T}\right]}^{-1}{\mathrm{\Σ}}_{i=1}^M{\mathrm{cw}}_i{\mathrm{\τ}}_i{g}_i^{\′}/d_i$ 和中，得到波达方向的单位矢量

其中，

和

分别表示

的第一维和第二维分量；

判决决策子模块，用于判断如果

其中η表示一个小于且接近于1的常数，则结束迭代，

即为最终的波达方向，并将结果输出给波达角输出子模块；否则

并将该设置结果输入迭代处理子模块继续执行迭代；

迭代处理子模块，用于根据波达方向

求取新的权重系数， $w_i=\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2)/{\mathrm{\Σ}}_i\exp ({\mathrm{\δ}}_i^2/{\mathrm{\σ}}^2),$ 其中 ${\mathrm{\δ}}_i=\arccos \left({\hat{x}}^T{g}_i\right)-\arccos ({\mathrm{c\τ}}_i/d_i),$ ${\mathrm{\σ}}^2={\mathrm{\Σ}}_i{\mathrm{\δ}}_i^2/M;$ 跳转到第一处理子模块；和波达角输出子模块，用于将波达方向矢量

转化为方位角和仰角，完成声源定位。

10.根据权利要求9所述的用于平面阵列的远场波达角估计***，其特征在于，所述波达角输出子模块采用反正切三角函数将波达方向矢量

转化为方位角和仰角，其中，方位角为当

时，仰角表示为 $\arccos (\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2}/\sqrt{{{\hat{x}}_1}^2+{{\hat{x}}_2}^2+{{\hat{x}}_3}^2})$ 和，当

时，仰角表示为

Images