文章目錄

• 一. 關(guān)于麥克風(fēng)陣列
• 二. 關(guān)于聲源定位
• 三. 基于廣義互相關(guān)(GCC)計(jì)算時(shí)延
• 四. 基于時(shí)延差的聲源定位法
• • 1. 近場(chǎng)模型
  • 2. 遠(yuǎn)場(chǎng)模型
• 五. 三維空間陣列的聲源定位系統(tǒng)實(shí)現(xiàn)
• • 1. 推導(dǎo)過(guò)程
• 六. 六元圓形麥克風(fēng)陣列聲源定位
• 七. 相關(guān)鏈接

一. 關(guān)于麥克風(fēng)陣列

麥克風(fēng)陣列： 麥克風(fēng)陣列是由一定數(shù)目的聲學(xué)傳感器（麥克風(fēng)）按照一定規(guī)則排列的多麥克風(fēng)系統(tǒng)，而基于麥克風(fēng)陣列的聲源定位是指用麥克風(fēng)拾取聲音信號(hào)，通過(guò)對(duì)麥克風(fēng)陣列的各路輸出信號(hào)進(jìn)行分析和處理，得到一個(gè)或者多個(gè)聲源的位置信息。

麥克風(fēng)陣列系統(tǒng)的聲源定位技術(shù)研究意義在于： 輸入的信息只有兩個(gè)方向難以確定聲源的位置，人類(lèi)的聽(tīng)覺(jué)系統(tǒng)主要取決于頭和外耳氣壓差聲波實(shí)現(xiàn)聲源定位。假使沒(méi)有這個(gè)壓力差，只能定位在平面上聲源的位置，但就無(wú)法知道聲音是從前面，或從后面?zhèn)鱽?lái)的。因此，由人的聽(tīng)覺(jué)系統(tǒng)，科技研發(fā)人員得到了靈感，使用多個(gè)麥克風(fēng)系統(tǒng)可以實(shí)現(xiàn)在三維空間中的聲源位置的定位，麥克風(fēng)的數(shù)量越多，所接收到的信息量也越多。聲源的聲源定位和聲源增強(qiáng)是實(shí)現(xiàn)智能處理的兩個(gè)關(guān)鍵問(wèn)題，而聲源定位是實(shí)現(xiàn)語(yǔ)音增強(qiáng)的前提和基礎(chǔ)。一個(gè)麥克風(fēng)的信息量較少，使得聲源定位所需的信息缺乏，而麥克風(fēng)陣列克服了上述缺點(diǎn)，充分利用每個(gè)麥克風(fēng)信號(hào)之間的數(shù)據(jù)相關(guān)性，并加以融合，可以實(shí)現(xiàn)聲源定位。

麥克風(fēng)陣列聲源定位技術(shù)的應(yīng)用： 廣泛應(yīng)用于國(guó)防、智能機(jī)器人、視頻會(huì)議及語(yǔ)音增強(qiáng)等眾多領(lǐng)域,尤其在當(dāng)下以智能辦公和智能家居為主要室內(nèi)場(chǎng)景的遠(yuǎn)場(chǎng)語(yǔ)音交互系統(tǒng)中。

二. 關(guān)于聲源定位

目前基于麥克風(fēng)陣列的聲源定位方法主要有三種：基于最大輸出功率的可控波束成形的定位方法、基于高分辨譜估計(jì)的定位方法、基于到達(dá)時(shí)延差估計(jì)的定位方法(Time Difference of Arrival，TDOA)。

基于最大輸出功率的可控波束成形定位法。 波束形成法的原理是將麥克風(fēng)接收到的信號(hào)進(jìn)行濾波加權(quán)求和來(lái)形成波束，按照一定的規(guī)律對(duì)聲源位置進(jìn)行搜索，當(dāng)麥克風(fēng)達(dá)到最大輸出功率時(shí)，為時(shí)搜索到的聲源位置即為真實(shí)的聲源方位。波束形成可分為常規(guī)的波束形成CBF（Conventional Beam Forming）、CBF+Adaptive Filter和自適應(yīng)波束形成ABF（Adaptive Beam Forming）。

基于高分辨譜估計(jì)的定位法。 基于高分辨率譜估計(jì)的定位方法通過(guò)分解協(xié)方差矩陣估計(jì)聲源方位。適合多個(gè)聲源的情況，且聲源的分辨率與陣列尺寸無(wú)關(guān)，突破了物理限制。該方法的優(yōu)點(diǎn)是不受采樣頻率限制，且在一定程度下可以實(shí)現(xiàn)任意程度的定位，但是該方法計(jì)算復(fù)雜度較高，抗噪和抗混響性能較差，因此該方法適合在一些特定的環(huán)境下使用。這類(lèi)方法可以拓展到寬帶處理，但是對(duì)誤差十分敏感（如麥克風(fēng)單體誤差，通道誤差），適合遠(yuǎn)場(chǎng)模型，且矩陣運(yùn)算量巨大。

基于到達(dá)時(shí)延差估計(jì)的定位法。 TDOA（time difference of arrival）是先后估計(jì)聲源到達(dá)不同麥克風(fēng)的時(shí)延差，通過(guò)時(shí)延來(lái)計(jì)算距離差，再利用距離差和麥克風(fēng)陣列的空間幾何位置來(lái)確定聲源的位置。可分為T(mén)DOA估計(jì)（估計(jì)信號(hào)到達(dá)各麥克風(fēng)的時(shí)間差）和TDOA定位（運(yùn)用幾何關(guān)系確定聲源位置）兩步。

三. 基于廣義互相關(guān)(GCC)計(jì)算時(shí)延

時(shí)延估計(jì)有很多種，比較經(jīng)典就是廣義互相關(guān)函數(shù) (Generalized Cross Correlation, GCC) 估計(jì)時(shí)延，這里簡(jiǎn)單介紹基于廣義互相關(guān)函數(shù)估計(jì)時(shí)延的方法。

在噪聲存在情況下，一個(gè)由遠(yuǎn)處聲源發(fā)出的，并且被兩個(gè)不同空間中的麥克風(fēng)監(jiān)聽(tīng)的信號(hào)可以數(shù)學(xué)建模為：

$x_1=s_1(t)+n_1(t)$?????? ???????????????????(1)

$x_2=\alpha s_1(t?D)+n_2(t)$???????? ??????(2)

其中，$s(t)$ 是聲音信號(hào)，$n_1(t)、n_2(t)$是兩個(gè)聲音傳感器檢測(cè)噪聲。 三者是穩(wěn)定的隨機(jī)過(guò)程，且互不相關(guān)。

計(jì)算 $x_1$ 與 $x_2$ 的互相關(guān)函數(shù)：

$R_{x_1{x}_2}(\tau )=E[x_1(t)?x_2(t?\tau )]$???? ? ? ? ??? (3)

${\hat{R}}_{x_1{x}_2}(\tau )=\frac{1}{T?\tau }{\int }_{\tau }^T{x}_1(t)x_2(t?\tau )\mathrm{d}t$???? ? (4)

其中估計(jì)的時(shí)延 $D$ 為互相關(guān)函數(shù)值達(dá)到最大值時(shí)取得的 $\tau$ 值，即：

$\hat{D}=\underset{\tau }{\mathrm{argmax}}{R}_{x_1{x}_2}(\tau )$???? ? ? ? ? ? ? ? ? ? ? ? ? ? ? (5)

MATLAB 例程：

% 導(dǎo)入兩個(gè)麥克風(fēng)的音頻數(shù)據(jù) [y_0,Fs] = audioread('音軌-0.wav'); [y_1] = audioread('音軌-1.wav');fprintf('采樣頻率：%d\n ······\n', Fs);% 取出兩段音頻前2048個(gè)采樣點(diǎn)，并作互相關(guān)處理，繪制曲線。 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個(gè)點(diǎn)為一幀，計(jì)算互相關(guān)得到的時(shí)延，繪制出兩段音頻的時(shí)延變化。 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的樹(shù)莓派六麥克風(fēng)套件，拿手機(jī)放歌作為聲源圍繞麥克風(fēng)陣列移動(dòng)，錄制音頻，并保存為6個(gè)wav音頻文件，采用率為44100Hz。對(duì)其中兩個(gè)麥克風(fēng)音頻文件做處理，如果取初始的2048個(gè)采樣點(diǎn)做互相關(guān)處理就能得到下面圖一，圖一的峰值的橫坐標(biāo)即為估計(jì)時(shí)延。如果以2048個(gè)采樣點(diǎn)作為一幀，對(duì)兩段音頻的每一幀進(jìn)行互相關(guān)，并找到每幀的估計(jì)時(shí)延，繪制出來(lái)就能得到時(shí)延的變化曲線，如下面圖二。

注意： 這里的時(shí)延并不是指具體的時(shí)間，而是采樣點(diǎn)數(shù)。想要得到具體時(shí)間還需要除以采樣頻率。以2048采樣點(diǎn)為一幀，就是以0.464s為一幀（$\frac{2048}{44100}\approx 0.464$s）。

四. 基于時(shí)延差的聲源定位法

在對(duì)麥克風(fēng)陣列進(jìn)行建模之前，我們需要分清楚什么近場(chǎng)與遠(yuǎn)場(chǎng)。顧名思義，離麥克風(fēng)近則符合近場(chǎng)模型 ，離得遠(yuǎn)則符合遠(yuǎn)場(chǎng)模型 。

假設(shè)$L$ 為陣列間距，$\lambda$ 為聲波波長(zhǎng)，$M$為聲源與麥克風(fēng)的距離，我們定義 $\frac{2L^2}{\lambda }$ 為遠(yuǎn)近場(chǎng)臨界值。

當(dāng) $M<\frac{2L^2}{\lambda }$ 時(shí)，符合近場(chǎng)模型，此時(shí)聲源到達(dá)麥克風(fēng)陣列的波形視為球面波。

當(dāng) $M<\frac{2L^2}{\lambda }$ 時(shí)，符合遠(yuǎn)場(chǎng)模型，此時(shí)聲源到達(dá)麥克風(fēng)陣列的波形視為平面波。

可聽(tīng)聲的機(jī)械波頻帶為20Hz ~20000Hz，機(jī)械波波長(zhǎng)大約在1.7cm ~ 17m（聲速取340m/s）。然而在現(xiàn)實(shí)生活中，過(guò)高頻率或過(guò)低頻率的聲波都是非常少量的，以人的聲音為例，人語(yǔ)音頻帶的范圍大概為300Hz至3400Hz，波長(zhǎng)范圍為0.1m~ 1.12m，當(dāng)陣列間隔取4cm時(shí)，遠(yuǎn)近場(chǎng)臨界值范圍為 0m~3.2m。若取中間值，則可以認(rèn)為1.6m內(nèi)符合近場(chǎng)模型，1.6m外符合遠(yuǎn)場(chǎng)模型。

1. 近場(chǎng)模型

當(dāng) $M<\frac{2L^2}{\lambda }$ 時(shí)，符合 近場(chǎng)模型，此時(shí)聲源到達(dá)麥克風(fēng)陣列的波形視為 球面波。

近場(chǎng)模型至少需要3個(gè)麥克風(fēng)，以最簡(jiǎn)單的3麥克風(fēng)模型為例（如圖：）。假設(shè) 分別表示第一個(gè)麥克風(fēng)與第二和第三個(gè)麥克風(fēng)之間的時(shí)延，那么有：

${\tau }_{12}=\frac{r_2?r_1}{c}$?????? ??(6)

${\tau }_{13}=\frac{r_3?r_1}{c}$?????? ??(7)

其中，$c$ 為聲速。標(biāo)準(zhǔn)大氣壓、$15°$ 條件下，聲速為 $340m/s$。

根據(jù)麥克風(fēng)陣列的的幾何關(guān)系，由余弦定理，可以得到：

$r_2^2=r_1^2+d^2+2r_{1}dcos{\theta }_1$?????(8)

$r_3^2=r_1^2+4d^2+4r_{1}dcos{\theta }_1$????(9)

其中， 可以通過(guò)互相關(guān)GCC求得（下面會(huì)細(xì)說(shuō)），且 為已知。結(jié)合(6) ~ (9) 即可求出未知量 ，結(jié)合坐標(biāo)系可以求出$s(k)$ 的坐標(biāo)。

2. 遠(yuǎn)場(chǎng)模型

當(dāng) $M<\frac{2L^2}{\lambda }$ 時(shí)，符合 遠(yuǎn)場(chǎng)模型，此時(shí)聲源到達(dá)麥克風(fēng)陣列的波形視為 平面波。

對(duì)于兩個(gè)麥克風(fēng)的情況，只能計(jì)算出方位角，無(wú)法計(jì)算出方位距離。假設(shè) $\tau$ 為聲波到達(dá)兩個(gè)麥克風(fēng)之間的時(shí)延，則有：

$\tau =\frac{dcos\theta }{c}$??????????????????????(10)

$\theta =arccos\frac{c\tau }{d}$??????????????(11)

五. 三維空間陣列的聲源定位系統(tǒng)實(shí)現(xiàn)

計(jì)算方法：幾何求解、最小二乘法

上面所舉的近場(chǎng)、遠(yuǎn)場(chǎng)的定位模型例子，分析都是最簡(jiǎn)模型下，平面聲源的定位。下面以三維麥克風(fēng)陣列為例，通過(guò)幾何的方式求解三維空間下聲源的定位。

1. 推導(dǎo)過(guò)程

一個(gè)三維麥克風(fēng)陣列，麥克風(fēng)分別為 ，聲源 $S$ 符合近場(chǎng)模型。現(xiàn)以麥克風(fēng) $m_0$ 為原點(diǎn)，如圖建立直角坐標(biāo)系。推導(dǎo)公式之前，需要先確定以下概念：

序號(hào)概念字符

1	麥克風(fēng)的坐標(biāo)	$m_i$ ? $i\in [0,n]$
2	聲源的估計(jì)坐標(biāo)	$S$
3	聲源 $s$ 到 $m_i$ 與 $m_1$ 的估計(jì)距離差	${\hat{d}}_i$
4	麥克風(fēng) $m_i$ 到原點(diǎn)的的距離	$R_i$
5	聲源 $s$ 到原點(diǎn)的的距離	$R_s$

如上圖，根據(jù)三角形 $m_0{m}_{i}s$，由余弦定理有：

$(R_s+d_i{)}^2=R_i^2?2m_i^{T}s+R_s^2$??????????????????????(12)

公式 (12) 展開(kāi)后得到：

$R_i^2?2m_i^{T}s?2R_s{d}_i?d_i^2=0$?????????????????????????(13)

由于 $d_i$ 是通過(guò)估計(jì)的時(shí)延得到的，與實(shí)際值之間會(huì)有一個(gè)偏差，因此公式 (10) 不為零，可寫(xiě)為：

$\epsilon =(R_i^2?d_i^2)?2m_i^{T}s?2R_s{d}_i$??????????????????????(14)

此時(shí)已經(jīng)得到目標(biāo)值的誤差函數(shù)，使用最小二乘法求解，問(wèn)題可以轉(zhuǎn)化為：估計(jì)聲源坐標(biāo) $s(x,y,z)$，使得最終的誤差平方和最小，即：

$\underset{s(x,y,z)}{\mathrm{argmin}}\sum _{i=1}^n[(R_i^2?d_i^2)?2m_i^{T}s?2R_s{d}_i{]}^2$??????(15)

將公式 (14) 寫(xiě)成矩陣形式：

$\epsilon =2R_s\hat{D}+2Ms?\delta$??????????????????????????????????????????(16)

其中，$M=?\begin{array}{c}{m}_1\\ m_2\\ m_3\\ ....\\ m_n\end{array}?=?\begin{array}{ccc}{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{array}?$ ??? ?$\hat{D}=?\begin{array}{c}{\hat{d}}_1\\ {\hat{d}}_2\\ {\hat{d}}_3\\ ...\\ {\hat{d}}_n\end{array}?$ ??????? $\delta =?\begin{array}{c}{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{array}?$

公式 (16) 可以簡(jiǎn)化為：

$\epsilon =A\mu ?b$???????????????????????????????????(17)

其中，$A=\left[\begin{array}{cc}M& \hat{D}\end{array}\right]$ ????? $\mu =\left[\begin{array}{c}s\\ R_s\end{array}\right]$ ?????? $b=\frac{1}{2}\delta$

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

$\hat{\mu }=(A^{T}A{)}^{?1}{A}^{T}b$???????????????????????(18)

關(guān)于矩陣形式下最小二乘解推導(dǎo)過(guò)程可以參考下面的鏈接（點(diǎn)擊跳轉(zhuǎn)），此處不詳細(xì)講解。

公式 (18) 即為計(jì)算結(jié)果，下面對(duì)結(jié)果進(jìn)行進(jìn)一步化簡(jiǎn)：

定義沿 $A$ 到 $\hat{D}$ 的投影矩陣為： $P_A=\frac{\hat{D}{\hat{D}}^T}{{\hat{D}}^T\hat{D}}$ ????????????????????????????????????????????????????(19)

沿 $\hat{D}$ 到 $A$ 的投影矩陣即為：$P_D=I?P_A=I?\frac{\hat{D}{\hat{D}}^T}{{\hat{D}}^T\hat{D}}$（ $I$ 為單位矩陣）?????(20)

根據(jù)投影矩陣的性質(zhì)，可以得到：$A=P_D\left[\begin{array}{cc}M& 0\end{array}\right]$ ??????????????????????(21)

將公式 (21) 帶入 (17) 中得到：$\epsilon =P_{D}Ms?b$??????????????????????????????(22)

公式（22）的最小二乘解即為最終簡(jiǎn)化結(jié)果。

最終簡(jiǎn)化結(jié)果：

$s=(M^T{P}_{D}M{)}^{?1}{M}^T{P}_{D}b$

$M$只與麥克風(fēng)陣列的個(gè)數(shù)與排列有關(guān)，$P_D$ 與 $b$ 在計(jì)算出時(shí)延估計(jì)的距離差 ($d_i$) 后也可以求得。

六. 六元圓形麥克風(fēng)陣列聲源定位

硬件平臺(tái)： Seed Studio 的六元麥克風(fēng)陣列套件（官網(wǎng))）

軟件平臺(tái)： Matlab 2019a

實(shí)際應(yīng)用中，二維圓形陣列主要采用遠(yuǎn)場(chǎng)模型較多，因?yàn)榭梢栽诒ＷC計(jì)算精度的同時(shí)，大大降低復(fù)雜性，減少運(yùn)算量，且不需要較大的陣列間距，此處使用近場(chǎng)模型純屬因?yàn)橐霎呍O(shè)😝。公式都是根據(jù)三維陣列來(lái)推導(dǎo)的，之所以使用二維平面陣列是因?yàn)?#xff1a;沒(méi)錢(qián)搭三維陣列，但又需要點(diǎn)實(shí)驗(yàn)數(shù)據(jù)放進(jìn)畢業(yè)論文中。

二維陣列缺少$z$軸的基元，使用第六章中的公式求解出的聲源坐標(biāo)值沒(méi)有實(shí)際意義，但是可以通過(guò) $y$ 與 $x$ 的比值得到聲源的偏航角。

${\theta }_{yaw}=arctan\frac{y}{x}$????????(23)

MATLAB 代碼：

#代碼寫(xiě)得太爛，自行腦補(bǔ)

現(xiàn)以樹(shù)莓派Seed Studio麥克風(fēng)陣列套件為硬件平臺(tái)進(jìn)行環(huán)境聲音的獲取，將各個(gè)麥克風(fēng)所獲取到的音頻儲(chǔ)存為wav格式文件，采樣率為44100Hz，將音頻文件在PC機(jī)中使用Matlab進(jìn)行處理運(yùn)算。麥克風(fēng)陣列為六元均勻圓陣，屬于二維平面陣列，陣元之間的間距為4cm。

如上圖所示，為六個(gè)麥克風(fēng) 收集到的聲音曲線，其中橫坐標(biāo)表示時(shí)間，縱坐標(biāo)表示聲音幅值。麥克風(fēng)錄制到的信號(hào)是一段長(zhǎng)時(shí)間的音頻片段，因此需要進(jìn)行分幀處理，此處以2048個(gè)采樣點(diǎn)為一幀，由于時(shí)間=采樣點(diǎn)數(shù)/采樣頻率，因此以2048個(gè)采樣點(diǎn)為一幀即以0.0464秒為一幀。我們假設(shè)一幀時(shí)間內(nèi)，聲音信號(hào)是聯(lián)合平穩(wěn)的，聲源相對(duì)于參考坐標(biāo)系的相位運(yùn)動(dòng)足夠小可以忽略不計(jì)。

如上圖所示，左圖為該陣列對(duì)靜態(tài)聲源的估計(jì)曲線，右圖則為對(duì)動(dòng)態(tài)聲源的估計(jì)曲線。

左圖： 在空曠且安靜的室內(nèi)環(huán)境，使用手機(jī)播放歌曲作為聲源，在距離麥克風(fēng)陣列0.5m處保持靜止，同時(shí)記錄下六個(gè)麥克風(fēng)錄制的音頻信息，采樣率為44100Hz，以2048個(gè)采樣點(diǎn)為一幀，使用聲達(dá)時(shí)延法計(jì)算每一幀聲源的方位角。

右圖： 在空曠且安靜的室內(nèi)環(huán)境，使用手機(jī)播放歌曲作為聲源，在距離麥克風(fēng)陣列中心0.5m處繞麥克風(fēng)陣列中心勻速旋轉(zhuǎn)，同時(shí)記錄下六個(gè)麥克風(fēng)錄制的音頻信息，采樣率為44100Hz，以2048個(gè)采樣點(diǎn)為一幀，使用聲達(dá)時(shí)延法計(jì)算每一幀聲源的方位角。

七. 相關(guān)鏈接

相關(guān)鏈接：

1. 矩陣形式下最小二乘解推導(dǎo)：https://zhuanlan.zhihu.com/p/87582571

2. 關(guān)于投影矩陣與最小二乘：https://www.cnblogs.com/bigmonkey/p/9897047.html

3. 關(guān)于投影矩陣概念及其性質(zhì)：https://blog.csdn.net/niu_123ming/article/details/86371308

4. 參考文獻(xiàn)：基于六元空間陣列的聲源定位系統(tǒng)實(shí)現(xiàn)

5. 參考文獻(xiàn)：基于麥克風(fēng)陣列的室內(nèi)三維聲源定位優(yōu)化算法

總結(jié)

以上是生活随笔為你收集整理的基于时延法的麦克风阵列声源定位分析的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問(wèn)題。

如果覺(jué)得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。