文章目录

一. 关于麦克风阵列
二. 关于声源定位
三. 基于广义互相关(GCC)计算时延
四. 基于时延差的声源定位法
- 1. 近场模型
- 2. 远场模型
五. 三维空间阵列的声源定位系统实现
- 1. 推导过程
六. 六元圆形麦克风阵列声源定位
七. 相关链接

一. 关于麦克风阵列

麦克风阵列： 麦克风阵列是由一定数目的声学传感器（麦克风）按照一定规则排列的多麦克风系统，而基于麦克风阵列的声源定位是指用麦克风拾取声音信号，通过对麦克风阵列的各路输出信号进行分析和处理，得到一个或者多个声源的位置信息。

麦克风阵列系统的声源定位技术研究意义在于： 输入的信息只有两个方向难以确定声源的位置，人类的听觉系统主要取决于头和外耳气压差声波实现声源定位。假使没有这个压力差，只能定位在平面上声源的位置，但就无法知道声音是从前面，或从后面传来的。因此，由人的听觉系统，科技研发人员得到了灵感，使用多个麦克风系统可以实现在三维空间中的声源位置的定位，麦克风的数量越多，所接收到的信息量也越多。声源的声源定位和声源增强是实现智能处理的两个关键问题，而声源定位是实现语音增强的前提和基础。一个麦克风的信息量较少，使得声源定位所需的信息缺乏，而麦克风阵列克服了上述缺点，充分利用每个麦克风信号之间的数据相关性，并加以融合，可以实现声源定位。

麦克风阵列声源定位技术的应用： 广泛应用于国防、智能机器人、视频会议及语音增强等众多领域,尤其在当下以智能办公和智能家居为主要室内场景的远场语音交互系统中。

二. 关于声源定位

目前基于麦克风阵列的声源定位方法主要有三种：基于最大输出功率的可控波束成形的定位方法、基于高分辨谱估计的定位方法、基于到达时延差估计的定位方法(Time Difference of Arrival，TDOA)。

基于最大输出功率的可控波束成形定位法。 波束形成法的原理是将麦克风接收到的信号进行滤波加权求和来形成波束，按照一定的规律对声源位置进行搜索，当麦克风达到最大输出功率时，为时搜索到的声源位置即为真实的声源方位。波束形成可分为常规的波束形成CBF（Conventional Beam Forming）、CBF+Adaptive Filter和自适应波束形成ABF（Adaptive Beam Forming）。
基于高分辨谱估计的定位法。 基于高分辨率谱估计的定位方法通过分解协方差矩阵估计声源方位。适合多个声源的情况，且声源的分辨率与阵列尺寸无关，突破了物理限制。该方法的优点是不受采样频率限制，且在一定程度下可以实现任意程度的定位，但是该方法计算复杂度较高，抗噪和抗混响性能较差，因此该方法适合在一些特定的环境下使用。这类方法可以拓展到宽带处理，但是对误差十分敏感（如麦克风单体误差，通道误差），适合远场模型，且矩阵运算量巨大。
基于到达时延差估计的定位法。 TDOA（time difference of arrival）是先后估计声源到达不同麦克风的时延差，通过时延来计算距离差，再利用距离差和麦克风阵列的空间几何位置来确定声源的位置。可分为TDOA估计（估计信号到达各麦克风的时间差）和TDOA定位（运用几何关系确定声源位置）两步。

三. 基于广义互相关(GCC)计算时延

时延估计有很多种，比较经典就是广义互相关函数 (Generalized Cross Correlation, GCC) 估计时延，这里简单介绍基于广义互相关函数估计时延的方法。

在噪声存在情况下，一个由远处声源发出的，并且被两个不同空间中的麦克风监听的信号可以数学建模为：

x1=s1(t)+n1(t)x_1=s_1(t)+n_1(t)x1=s1(t)+n1(t) (1) x2=αs1(t−D)+n2(t)x_2=αs_1(t-D)+n_2(t)x2=αs1(t−D)+n2(t) (2)

其中，s(t)s(t)s(t) 是声音信号，n1(t)、n2(t)n_1(t)、n_2(t)n1(t)、n2(t)是两个声音传感器检测噪声。三者是稳定的随机过程，且互不相关。

计算 x1x_1x1 与 x2x_2x2 的互相关函数：

Rx1x2(τ)=E[x1(t)∗x2(t−τ)]R_{x_1x_2}(τ)=E [ x_1(t)*x_2(t-τ) ]Rx1x2(τ)=E[x1(t)∗x2(t−τ)] (3) R^x1x2(τ)=1T−τ∫τTx1(t)x2(t−τ)dt\hat R_{x_1x_2}(τ)=\frac{1}{T-τ}\int_τ^T {x_1(t)x_2(t-τ)} \,{\rm d}tR^x1x2(τ)=T−τ1∫τTx1(t)x2(t−τ)dt (4)

其中估计的时延 DDD 为互相关函数值达到最大值时取得的 τττ 值，即：

D^=argmaxτRx1x2(τ)\hat D=\underset {\ τ}{argmax} R_{x_1x_2}(τ)D^= τargmaxRx1x2(τ) (5)

MATLAB 例程：

% 导入两个麦克风的音频数据
[y_0,Fs] = audioread('音轨-0.wav');
[y_1] = audioread('音轨-1.wav');fprintf('采样频率：%d\n ······\n', Fs);% 取出两段音频前2048个采样点，并作互相关处理，绘制曲线。
A = y_0(1:2048);
B = y_1(1:2048);[value,delay] = xcorr(A,B);subplot(1,2,1);
plot(delay, value);D = zeros(1,926); % 以2048个点为一帧，计算互相关得到的时延，绘制出两段音频的时延变化。
for a = 1:2048:size(y_0,1)-2048A = y_0(a:a+2048);B = y_1(a:a+2048);[value,delay]=xcorr(A,B);value_max_idx = find(value==max(value));D1 = delay(value_max_idx);D((a-1)/2048+1) = D1;
endsubplot(1,2,2);
plot(D);

使用ReSpeaker的树莓派六麦克风套件，拿手机放歌作为声源围绕麦克风阵列移动，录制音频，并保存为6个wav音频文件，采用率为44100Hz。对其中两个麦克风音频文件做处理，如果取初始的2048个采样点做互相关处理就能得到下面图一，图一的峰值的横坐标即为估计时延。如果以2048个采样点作为一帧，对两段音频的每一帧进行互相关，并找到每帧的估计时延，绘制出来就能得到时延的变化曲线，如下面图二。

注意： 这里的时延并不是指具体的时间，而是采样点数。想要得到具体时间还需要除以采样频率。以2048采样点为一帧，就是以0.464s为一帧（204844100≈0.464\frac{2048}{44100}≈0.464441002048≈0.464s）。

四. 基于时延差的声源定位法

在对麦克风阵列进行建模之前，我们需要分清楚什么近场与远场。顾名思义，离麦克风近则符合近场模型 ，离得远则符合远场模型 。

假设LLL 为阵列间距，λλλ 为声波波长，MMM为声源与麦克风的距离，我们定义 2L2λ\frac{2L^2}{λ}λ2L2 为远近场临界值。

当 M<2L2λM<\frac{2L^2}{λ}M<λ2L2 时，符合近场模型，此时声源到达麦克风阵列的波形视为球面波。

当 M<2L2λM<\frac{2L^2}{λ}M<λ2L2 时，符合远场模型，此时声源到达麦克风阵列的波形视为平面波。

可听声的机械波频带为20Hz ~20000Hz，机械波波长大约在1.7cm ~ 17m（声速取340m/s）。然而在现实生活中，过高频率或过低频率的声波都是非常少量的，以人的声音为例，人语音频带的范围大概为300Hz至3400Hz，波长范围为0.1m~ 1.12m，当阵列间隔取4cm时，远近场临界值范围为 0m~3.2m。若取中间值，则可以认为1.6m内符合近场模型，1.6m外符合远场模型。

1. 近场模型

当 M<2L2λM<\frac{2L^2}{λ}M<λ2L2 时，符合 近场模型，此时声源到达麦克风阵列的波形视为 球面波。

近场模型至少需要3个麦克风，以最简单的3麦克风模型为例（如图：y1、y2、y3y_1、y_2、y_3y1、y2、y3）。假设 τ12、τ13τ_{12} 、τ_{13}τ12、τ13分别表示第一个麦克风与第二和第三个麦克风之间的时延，那么有：

τ12=r2−r1cτ_{12}=\frac{r_2-r_1}{c}τ12=cr2−r1 (6) τ13=r3−r1cτ_{13}=\frac{r_3-r_1}{c}τ13=cr3−r1 (7)

其中，ccc 为声速。标准大气压、15°15°15° 条件下，声速为 340m/s340m/s340m/s。

根据麦克风阵列的的几何关系，由余弦定理，可以得到：

r22=r12+d2+2r1dcosθ1r_2^2=r_1^2+d^2+2r_1dcosθ_1r22=r12+d2+2r1dcosθ1 (8) r32=r12+4d2+4r1dcosθ1r_3^2=r_1^2+4d^2+4r_1dcosθ_1r32=r12+4d2+4r1dcosθ1 (9)

其中， τ12、τ13τ_{12}、τ_{13}τ12、τ13可以通过互相关GCC求得（下面会细说），且 c、dc、dc、d 为已知。结合(6) ~ (9) 即可求出未知量 r1、r2、r3、θr_1、r_2、r_3、θr1、r2、r3、θ ，结合坐标系可以求出s(k)s(k)s(k) 的坐标。

2. 远场模型

当 M<2L2λM<\frac{2L^2}{λ}M<λ2L2 时，符合 远场模型，此时声源到达麦克风阵列的波形视为 平面波。

对于两个麦克风的情况，只能计算出方位角，无法计算出方位距离。假设 τττ 为声波到达两个麦克风之间的时延，则有：

τ=dcosθcτ=\frac{dcosθ}{c}τ=cdcosθ (10) θ=arccoscτdθ=arccos\frac{cτ}{d}θ=arccosdcτ (11)

五. 三维空间阵列的声源定位系统实现

计算方法：几何求解、最小二乘法

上面所举的近场、远场的定位模型例子，分析都是最简模型下，平面声源的定位。下面以三维麦克风阵列为例，通过几何的方式求解三维空间下声源的定位。

1. 推导过程

一个三维麦克风阵列，麦克风分别为 mi(i=0、1⋅⋅⋅n)m_i(i=0、1···n)mi(i=0、1⋅⋅⋅n)，声源 SSS 符合近场模型。现以麦克风 m0m_0m0 为原点，如图建立直角坐标系。推导公式之前，需要先确定以下概念：

序号	概念	字符
1	麦克风的坐标	mim_imi i∈[0,n]i \in[0,n]i∈[0,n]
2	声源的估计坐标	SSS
3	声源 sss 到 mim_imi 与 m1m_1m1 的估计距离差	d^i\hat d_id^i
4	麦克风 mim_imi 到原点的的距离	RiR_iRi
5	声源 sss 到原点的的距离	RsR_sRs

如上图，根据三角形 m0mism_0 m_i sm0mis，由余弦定理有：

(Rs+di)2=Ri2−2miTs+Rs2(R_s+d_i)^2=R_i^2-2m_i^Ts+R_s^2(Rs+di)2=Ri2−2miTs+Rs2 (12)

公式 (12) 展开后得到：

Ri2−2miTs−2Rsdi−di2=0R_i^2-2m_i^Ts-2R_sd_i-d_i^2=0Ri2−2miTs−2Rsdi−di2=0 (13)

由于 did_idi 是通过估计的时延得到的，与实际值之间会有一个偏差，因此公式 (10) 不为零，可写为：

ε=(Ri2−di2)−2miTs−2Rsdiε=(R_i^2-d_i^2)-2m_i^Ts-2R_sd_iε=(Ri2−di2)−2miTs−2Rsdi (14)

此时已经得到目标值的误差函数，使用最小二乘法求解，问题可以转化为：估计声源坐标 s(x,y,z)s(x,y,z)s(x,y,z)，使得最终的误差平方和最小，即：

argmins(x,y,z)∑i=1n[(Ri2−di2)−2miTs−2Rsdi]2\underset {\ s(x,y,z)}{argmin}\displaystyle\sum_{i=1}^{n} [(R_i^2-d_i^2)-2m_i^Ts-2R_sd_i]^2 s(x,y,z)argmini=1∑n[(Ri2−di2)−2miTs−2Rsdi]2 (15)

将公式 (14) 写成矩阵形式：

ε=2RsD^+2Ms−δε=2R_s\hat D+2Ms-\deltaε=2RsD^+2Ms−δ (16)

其中，M=[m1m2m3....mn]=[x1y1z1x2y2z2x3y3z3.........xnynzn]M=\begin{bmatrix}m_1\\m_2\\m_3\\....\\m_n\end{bmatrix}=\begin{bmatrix}x_1&y_1&z_1\\x_2&y_2&z_2\\x_3&y_3&z_3\\...&...&...\\x_n&y_n&z_n\end{bmatrix}M=⎣⎡m1m2m3....mn⎦⎤=⎣⎡x1x2x3...xny1y2y3...ynz1z2z3...zn⎦⎤ D^=[d^1d^2d^3...d^n]\hat D=\begin{bmatrix}\hat d_1\\\hat d_2\\\hat d_3\\\ ... \\\hat d_n\end{bmatrix}D^=⎣⎡d^1d^2d^3 ...d^n⎦⎤ δ=[R12−d^12R22−d^22R32−d^32..........Rn2−d^n2]\delta=\begin{bmatrix}R_1^2-\hat d_1^2\\R_2^2-\hat d_2^2\\R_3^2-\hat d_3^2\\..........\\R_n^2-\hat d_n^2\end{bmatrix}δ=⎣⎡R12−d^12R22−d^22R32−d^32..........Rn2−d^n2⎦⎤

公式 (16) 可以简化为：

ε=Aμ−bε=A\mu-bε=Aμ−b (17)

其中，A=[MD^]A=\begin{bmatrix}M&\hat D\end{bmatrix}A=[MD^] μ=[sRs]\mu=\begin{bmatrix}s\\R_s\end{bmatrix}μ=[sRs] b=12δb=\frac{1}{2}\deltab=21δ

公式 (17) 的最小二乘解可以表示为：

μ^=(ATA)−1ATb\hat\mu=(A^TA)^{-1}A^Tbμ^=(ATA)−1ATb (18)

关于矩阵形式下最小二乘解推导过程可以参考下面的链接（点击跳转），此处不详细讲解。

公式 (18) 即为计算结果，下面对结果进行进一步化简：

定义沿 AAA 到 D^\hat DD^ 的投影矩阵为： PA=D^D^TD^TD^P_A=\frac{\hat D \hat D^T}{\hat D^T \hat D }PA=D^TD^D^D^T (19)
沿 D^\hat DD^ 到 AAA 的投影矩阵即为：PD=I−PA=I−D^D^TD^TD^P_D=I-P_A=I-\frac{\hat D \hat D^T}{\hat D^T \hat D }PD=I−PA=I−D^TD^D^D^T（ III 为单位矩阵） (20)
根据投影矩阵的性质，可以得到：A=PD[M0]A=P_D\begin{bmatrix}M&0\end{bmatrix}A=PD[M0] (21)
将公式 (21) 带入 (17) 中得到：ε=PDMs−bε=P_D M s-bε=PDMs−b (22)
公式（22）的最小二乘解即为最终简化结果。

最终简化结果：

s=(MTPDM)−1MTPDbs=( M^T P_D M )^{-1}M^TP_Dbs=(MTPDM)−1MTPDb

MMM只与麦克风阵列的个数与排列有关，PDP_DPD 与 bbb 在计算出时延估计的距离差 (did_idi) 后也可以求得。

六. 六元圆形麦克风阵列声源定位

硬件平台： Seed Studio 的六元麦克风阵列套件（官网)）

软件平台： Matlab 2019a

实际应用中，二维圆形阵列主要采用远场模型较多，因为可以在保证计算精度的同时，大大降低复杂性，减少运算量，且不需要较大的阵列间距，此处使用近场模型纯属因为要做毕设

基于时延法的麦克风阵列声源定位分析相关推荐

【声源定位】球面散乱数据插值方法/似然估计hybrid spherical interpolation/maximum likelihood (SI/ML) 麦克风阵列声源定位
1.软件版本 MATLAB2021a 2.本算法理论知识点球面散乱数据插值方法/似然估计SI/ML 麦克风阵列声源定位 3.算法具体理论这个部分的程序如下所示: 这个部分理论如下所示: 本文最后的 ...
麦克风阵列声源定位 GCC-PHAT
麦克风阵列声源定位 GCC-PHAT 麦克风阵列声源定位(一) 利用麦克风阵列可以实现声源到达方向估计(direction-of-arrival (DOA) estimation),DOA估计的其中一 ...
音视频开发（40）---麦克风阵列声源定位 GCC-PHAT
麦克风阵列声源定位 GCC-PHAT 版权声明:本文为博主原创文章,未经博主允许不得转载. https://blog.csdn.net/u010592995/article/details/79735 ...
matlab 声源定位csdn_麦克风阵列声源定位 GCC-PHAT（一）
麦克风阵列声源定位(一) 0 a" N0 Q" t t2 l$ t) F利用麦克风阵列可以实现声源到达方向估计(direction-of-arrival (DOA) estima ...
麦克风阵列声源定位四通道麦克风数据库及TDOA双曲交汇定位算法实验
麦克风阵列声源定位四通道麦克风数据库建立四通道麦克风数据库建立物理模型的建立,来源于文献:SLoClas: A DATABASE FOR JOINT SOUND LOCALIZATION AND C ...
麦克风阵列声源定位效果测试
下列图片如果不清楚可以直接访问淘宝链接,从链接中的网盘资料进行拉取.从此链接看到的购买可以跟客服说,提我可以便宜50块钱~~~ 店铺链接:首页-智能语音开发者联盟-淘宝网产品链接:https://i ...
麦克风阵列声源定位解决方案
其高科技: http://www.keygotech.com/cn/solution/ssl/array/noise-source-location-based-on-mic-array 一般来说,基 ...
麦克风阵列声源定位 SRP-PHAT
DOA 声源定位方法一般可分为三类,一种是基于TDOA的两步算法(two-stage algorithm),一种是基于空间谱估计如MUSIC等,还有就是基于beamforming的方法,也就是这里要介 ...
matlab双麦克风阵列声源定位,双麦克风阵列拾音算法可实现空间滤波和360全向自动声源定位降噪...
不仅如此,还发掘了它的隐藏技能: 送孩子去上辅导班,忙着做笔记往往就错过正在讲的知识点,自从让他带着这只录音笔去上课,回来温习的学习效果翻倍! 爸妈想要出国旅行又怕沟通障碍,自己常因为工作忙没法陪同. ...

基于时延法的麦克风阵列声源定位分析