本發(fā)明涉及一種聲源定位方法，特別涉及一種強混響環(huán)境下基于參數(shù)化貝葉斯字典學(xué)習(xí)的聲源定位方法。

背景技術(shù)：

文獻“Structured sparsity models for reverberant speech separation,Audio,Speech,and Language Processing,IEEE/ACM Transactions on,vol.22,pp.620-633,2014”借助房間內(nèi)部采集的多個傳聲器信號，在房間內(nèi)外設(shè)置多個網(wǎng)格點估計房間形狀以及多級虛聲源，其次利用房間的虛聲源模型估計房間的吸聲系數(shù)，最后使用凸優(yōu)化的數(shù)據(jù)算法來進行室內(nèi)多源分離的研究。然而，當(dāng)室內(nèi)混響較強時，空間內(nèi)設(shè)定網(wǎng)格點數(shù)目會極大地增加，且吸聲參數(shù)未知，該方法所涉及得估計變量較多，計算效率低下，求解算法精度不高，因而無法準(zhǔn)確計算獲得室內(nèi)聲源的位置。因此，對于室內(nèi)聲源定位，對未知房間的聲學(xué)參數(shù)進行準(zhǔn)確估計是合理估計房間聲源位置信息的前提條件。

技術(shù)實現(xiàn)要素：

為了克服現(xiàn)有聲源定位方法精確度差的不足，本發(fā)明提供一種強混響環(huán)境下基于參數(shù)化貝葉斯字典學(xué)習(xí)的聲源定位方法。該方法基于傳聲器鏡像簡化字典，在稀疏貝葉斯的框架下，實現(xiàn)聲源定位和房間反射系數(shù)的同步估計求解。利用本發(fā)明方法將字典預(yù)先計算出來，完成聲源定位和反射系數(shù)估計的效率可以得到提升，這樣在沒有顯著提升計算量的情況下，提高了聲源定位方法的精度。本發(fā)明對復(fù)雜室內(nèi)環(huán)境下的聲源定位有很大的潛在實用價值，也可作為封閉空間中聲學(xué)參數(shù)估計的方法。

本發(fā)明解決其技術(shù)問題所采用的技術(shù)方案：一種強混響環(huán)境下基于參數(shù)化貝葉斯字典學(xué)習(xí)的聲源定位方法，其特點是包括以下步驟：

步驟一、定位環(huán)境選取能夠用于矩形或根據(jù)形狀直接獲取虛源的房間。選定房間內(nèi)某一維或二維平面進行室內(nèi)網(wǎng)格劃分，假設(shè)室內(nèi)網(wǎng)格點數(shù)目為G。

步驟二、構(gòu)建信號模型：聲源定位時需要獲得的目標(biāo)信息是聲源在房間中的位置，這一信息的量化離散化為分布在室內(nèi)空間以及室外鏡像空間的有限個網(wǎng)格點,每個網(wǎng)格點的聲源幅度為x_g。假設(shè)聲源存在于任意一個網(wǎng)格點內(nèi)，任意位置的接收信號y(r,f)表示為式(1)

式中，f為當(dāng)前頻率，c為空氣中的聲速，G為室內(nèi)網(wǎng)格點數(shù)目，β_g為對應(yīng)虛源的反射系數(shù)，x_g為反映聲源位置的向量。在所有的G個網(wǎng)格點上，x_g只有很少的一些非零元數(shù)，對應(yīng)于房間空間的實際聲源位置和鏡像空間的虛聲源位置。在構(gòu)造字典時，只在房間內(nèi)部設(shè)定網(wǎng)格點，每個字典元素為聲源存在的內(nèi)部網(wǎng)格點與接收點鏡像之間的通道疊加。由此得到

y_f＝A_fx_f (4)

式(2)右端表示傳聲器及其鏡像聲能量貢獻的總疊加，式(3)為字典原子計算公式，N₀為接收點鏡像級數(shù)，為n_s,i階鏡像對應(yīng)的聲壓反射系數(shù)。整理式(3)得到式(4)。式中，y_f∈C^M×1,x_f∈C^G×1,A_f∈C^M×G。

步驟三、接收點鏡像求解：對于尺寸為lx×ly×lz的矩形房間，假設(shè)接收點位置為(x₀,y₀,z₀)，某級接收點鏡像的位置為(x_i,y_i,z_i)，l、m、n是界于-∞和+∞之間的整數(shù)，則任意接收點鏡像的空間坐標(biāo)表示為

x_i＝2l·lx±x₀

y_i＝2m·ly±y₀ (5)

z_i＝2n·lz±z₀

鏡像的反射級數(shù)由下式給出

N＝|2l-Δ_l|+|2m-Δ_m|+|2n-Δ_n| (6)

Δl＝1表示求x_i式中x₀前取負(fù)；Δl＝0表示求x_i式中x₀前取正；Δm＝1表示求y_i式中y₀前取負(fù)；Δm＝0表示求y_i式中y₀前取正；Δn＝1表示求z_i式中z₀前取負(fù)；Δn＝0表示求z_i式中z₀前取正。

步驟四、構(gòu)造求解字典：按式(3)構(gòu)造反射系數(shù)β的不同量化數(shù)值下的字典A_f及其對β的一階導(dǎo)數(shù)A'_f和二階導(dǎo)數(shù)矩陣A″_f。

步驟五、位置向量x_f求解：在進行室內(nèi)聲定位時，房間內(nèi)的網(wǎng)格點數(shù)目G不僅遠遠大于傳聲器個數(shù)M，而且大于聲源個數(shù)。即x_f具有空間稀疏性，通過稀疏貝葉斯方法估計含有聲源位置信息的x。觀察信號由式(7)采樣得到。

I_M為M維的單位陣，x_f由下式得到：

x_f|α～CN(0,Λ^-1I) (8)

其中上式中CN表示復(fù)高斯分布，Gamma表示伽馬分布。上述概率模型通過變分貝葉斯-期望最大化算法求解。其中x_f的后驗概率通過式(9)得到：

表示對應(yīng)變量的估計值。

下面直接給出各參數(shù)估計的方法，其中服從均值為的Gamma分布。α₀服從均值為的高斯分布，Tr表示矩陣的跡。

步驟六、參數(shù)β的估計：上述各字典中包含反射系數(shù)β參數(shù)，估計該參數(shù)等價于使式(10)最小化：

該最優(yōu)化問題通過Newton法求解，迭代公式如下：

式中，A'_f表示對對應(yīng)的參數(shù)β求導(dǎo)，A″_f表示二階導(dǎo)數(shù)。通過對稀疏貝葉斯框架中的字典引入反射系數(shù)作為參數(shù)，并通過Newton法迭代解決該問題，通過分步求解，實現(xiàn)在估計聲源位置x_f的同時，估計出該反射系數(shù)β。改參數(shù)的估計一方面為房間聲學(xué)參數(shù)估計，同時也得到更準(zhǔn)確的房間聲學(xué)信息，達到更好的聲源定位效果。

步驟七、參數(shù)β量化：在房間信息完全未知的情況下，β的初值設(shè)置為0.5，通過帶參數(shù)的貝葉斯稀疏方法，分別利用式(9)和式(10)對聲源位置x_f和反射系數(shù)β進行估計，反射系數(shù)估計與聲源位置的估計是同時進行的，即在式(9)和式(10)之間反復(fù)進行，每步迭代都需要估計反射系數(shù)并通過量化尋找最接近的字典并進行下一步迭代，β不同將導(dǎo)致字典多次重新計算，帶來巨大的計算量。

為了避免重復(fù)計算字典，首先將β在其值域范圍(0～1)內(nèi)量化為若干離散數(shù)據(jù)點，然后將這些離散點的值對應(yīng)的字典在迭代前計算出來，將計算得到的β量化到對應(yīng)的離散值上并調(diào)用對應(yīng)的字典可顯著減少計算量。

本發(fā)明的有益效果是：該方法基于傳聲器鏡像簡化字典，在稀疏貝葉斯的框架下，實現(xiàn)聲源定位和房間反射系數(shù)的同步估計求解。利用本發(fā)明方法將字典預(yù)先計算出來，完成聲源定位和反射系數(shù)估計的效率可以得到提升，這樣在沒有顯著提升計算量的情況下，提高了聲源定位方法的精度。本發(fā)明對復(fù)雜室內(nèi)環(huán)境下的聲源定位有很大的潛在實用價值，也可作為封閉空間中聲學(xué)參數(shù)估計的方法。

分別對得到的采集信號做不同的信噪比處理并進行定位，選定頻率范圍為1kHz～2kHz。分別統(tǒng)計本發(fā)明方法與文獻方法下的計算時間，表1顯示了相應(yīng)方法在同一臺計算上完成一次定位任務(wù)的耗時情況，本發(fā)明方法提供的兩種耗時，一個是未量化反射系數(shù)條件下計算得到(左)，一種是量化反射系數(shù)條件下計算得到(右)。

表1不同方法耗時對比

表1表明在準(zhǔn)確估計反射系數(shù)的同時，計算耗時比文獻算法耗時少。本發(fā)明方法的計算時間不受房間混響程度的影響，它只于房間劃分的網(wǎng)格點數(shù)有直接的關(guān)系，網(wǎng)格點數(shù)目越多，求解向量維度越高，耗費時間越多。

下面結(jié)合附圖和具體實施方式對本發(fā)明作詳細說明。

附圖說明

圖1是本發(fā)明強混響環(huán)境下基于參數(shù)化貝葉斯字典學(xué)習(xí)的聲源定位方法的流程圖。

圖2是本發(fā)明方法用測試系統(tǒng)連接框圖。

圖3是10*10網(wǎng)格下的聲源信息矢量x迭代收斂過程。

圖4是20*20網(wǎng)格下的聲源信息矢量x迭代收斂過程。

圖5是不同信噪比下傳感器數(shù)目為10*10時的反射系數(shù)迭代求解均方誤差(MSE)曲線。

圖6是不同信噪比下傳感器數(shù)目為20*20時的反射系數(shù)迭代求解均方誤差(MSE)曲線。

圖中，MSE_x，MSE_beta，MSE_x0分別對應(yīng)x、β和文獻方法所得x₀。

具體實施方式

參照圖1-6。本發(fā)明強混響環(huán)境下基于參數(shù)化貝葉斯字典學(xué)習(xí)的聲源定位方法具體步驟如下：

步驟一、定位環(huán)境選取：本方法只可用于矩形或其他根據(jù)形狀可直接獲取虛源的房間。選定房間內(nèi)某一維或二維平面內(nèi)的進行聲源定位的室內(nèi)網(wǎng)格劃分，這里假設(shè)室內(nèi)網(wǎng)格點數(shù)目為G。

步驟二、構(gòu)建信號模型：聲源定位時需要獲得的目標(biāo)信息是聲源在房間中的位置，這一信息的量化離散化為分布在室內(nèi)空間以及室外鏡像空間的有限個網(wǎng)格點,每個網(wǎng)格點可能的聲源幅度為x_g。假設(shè)聲源可能存在于任意一個網(wǎng)格點內(nèi)，任意位置的接收信號y(r,f)可表示為式(1)

式中f為當(dāng)前頻率，c為空氣中的聲速，G為網(wǎng)格點的總數(shù)目，β_g為對應(yīng)虛源的反射系數(shù)，x_g為反映聲源位置的向量。在所有可能的G個網(wǎng)格點上，x_g只有很少的一些非零元數(shù)，對應(yīng)于房間空間的實際聲源位置和鏡像空間的虛聲源位置。在構(gòu)造字典時，只在房間內(nèi)部設(shè)定網(wǎng)格點，每個字典元素為聲源可能存在的內(nèi)部網(wǎng)格點與接收點鏡像之間的通道疊加。因此可得

y_f＝A_fx_f (4)

式(2)右端表示為傳聲器及其鏡像的聲能量貢獻的總疊加，式(3)為字典原子計算公式，N0為接收點鏡像級數(shù)，為n_s,i階鏡像對應(yīng)的聲壓反射系數(shù)。整理式(3)可得式(4)。式中，y_f∈C^M×1,x_f∈C^G×1,A_f∈C^M×G。

步驟三、接收點鏡像求解：對于尺寸為lx×ly×lz的矩形房間，假設(shè)接收點位置為(x₀,y₀,z₀)，某級接收點鏡像的位置為(x_i,y_i,z_i)，l、m、n是界于-∞和+∞之間的整數(shù)，則任意接收點鏡像的空間坐標(biāo)可表示為

x_i＝2l·lx±x₀

y_i＝2m·ly±y₀ (5)

z_i＝2n·lz±z₀

鏡像的反射級數(shù)可以由下式給出

N＝|2l-Δ_l|+|2m-Δ_m|+|2n-Δ_n| (6)

Δl＝1表示求x_i式中x₀前取負(fù)；Δl＝0表示求x_i式中x₀前取正；Δm＝1表示求y_i式中y₀前取負(fù)；Δm＝0表示求y_i式中y₀前取正；Δn＝1表示求z_i式中z₀前取負(fù)；Δn＝0表示求z_i式中z₀前取正。

步驟四、構(gòu)造求解字典：按式(3)構(gòu)造反射系數(shù)β的不同量化數(shù)值下的字典A_f及其對β的一階導(dǎo)數(shù)A'_f和二階導(dǎo)數(shù)矩陣A″_f。

步驟五、位置向量x_f求解：在進行室內(nèi)聲定位時，房間內(nèi)的網(wǎng)格點數(shù)目G不僅遠遠大于傳聲器個數(shù)M，而且大于聲源個數(shù)。即x_f具有空間稀疏性，通過稀疏貝葉斯方法估計含有聲源位置信息的x。觀察信號由式(7)采樣得到。

y_f|x_f,α₀～CN(A_fx_f,α₀^-1I_M) (7)

I_M為M維的單位陣，x_f由下式得到：

x_f|α～CN(0,Λ^-1I) (8)

其中上式中CN表示復(fù)高斯分布，Gamma表示伽馬分布。上述概率模型可通過變分貝葉斯-期望最大化算法求解。其中x_f的后驗概率可通過式(9)得到：

表示對應(yīng)變量的估計值。

下面直接給出各參數(shù)估計的方法，其中服從均值為的Gamma分布。α₀服從均值為的高斯分布，Tr表示矩陣的跡。

步驟六、參數(shù)β的估計：上述各字典中包含反射系數(shù)β參數(shù)，估計該參數(shù)可以等價于使式(10)最小化：

該最優(yōu)化問題可通過Newton法求解，迭代公式如下：

式中A'_f表示對對應(yīng)的參數(shù)β求導(dǎo)，A″_f表示二階導(dǎo)數(shù)。通過對稀疏貝葉斯框架中的字典引入反射系數(shù)作為參數(shù)，并通過Newton法迭代解決該問題，通過分步求解，可實現(xiàn)在估計聲源位置x_f的同時，估計出該反射系數(shù)β。改參數(shù)的估計一方面為房間聲學(xué)參數(shù)估計，同時也可得到更準(zhǔn)確的房間聲學(xué)信息，可達到更好的聲源定位效果。

步驟七、參數(shù)β量化：在房間信息完全未知的情況下，β的初值一般設(shè)置為0.5，通過帶參數(shù)的貝葉斯稀疏方法分別利用式(9)和式(10)對聲源位置x_f和反射系數(shù)β進行估計，反射系數(shù)估計與聲源位置的估計是同時進行的，即在(9)和(10) 之間反復(fù)進行，每步迭代都需要估計反射系數(shù)并通過量化尋找最接近的字典并進行下一步迭代，β不同將導(dǎo)致字典多次重新計算，帶來巨大的計算量。

為了避免重復(fù)計算字典，首先將β在其值域范圍(0～1)內(nèi)量化為若干離散數(shù)據(jù)點(本發(fā)明取1000個點的均勻量化方法)，然后將這些離散點的值對應(yīng)的字典在迭代前計算出來，將計算得到的β量化到對應(yīng)的離散值上并調(diào)用對應(yīng)的字典可顯著減少計算量。

應(yīng)用實施例：

房間為一封閉矩形空間。其長寬高分別為lx＝4m，ly＝5m，lz＝3m。

步驟1：在封閉房間內(nèi)部選定同一高度，布置15陣元的線性傳聲器陣列，陣列位置為x＝1.44m，z＝1.5m，y方向間隔為0.2m，起始端傳聲器位置為y＝0.1m。

步驟2：根據(jù)房間尺寸，選定z＝1.5m的平面為主要參考區(qū)域，x、y方向各劃分10段均勻的進行網(wǎng)格劃分，即房間內(nèi)部的網(wǎng)格點為100個(或20段共400點)，網(wǎng)格大小為40cm×50cm(20cm×25cm)。網(wǎng)格坐標(biāo)取網(wǎng)格中心位置坐標(biāo)。

步驟3：根據(jù)房間壁面材料，通過查表獲取各壁面吸聲系數(shù)或直接測量室內(nèi)混響時間確定房間的平均聲壓反射系數(shù)β₀。

步驟4：隨機選取一個預(yù)劃分的網(wǎng)格點作為聲源位置，按照附圖2連接電腦，以及相應(yīng)的采集設(shè)備。

步驟5：開啟所有設(shè)備并驅(qū)動聲源，持續(xù)播放預(yù)定的白噪聲信號并設(shè)定所采集時域信號的采樣頻率并進行錄制，保存數(shù)據(jù)。對各保存數(shù)據(jù)做短時傅立葉變換，獲取各傳聲器信號的頻域表示y_f。

步驟6：選定坐標(biāo)參考點，測量各傳聲器在房間內(nèi)的相對位置，根據(jù)房間內(nèi)部網(wǎng)格點劃分，求解各傳聲器鏡像的位置信息r_i以及各鏡像的反射級數(shù)n_s,i。

步驟7：根據(jù)各傳聲器鏡位置r_i、鏡像反射級數(shù)n_s,i以及各壁面反射系數(shù)初值β₀進行量化求解，根據(jù)公式(5)求解房間內(nèi)各網(wǎng)格點對應(yīng)的字典原子，形成各頻率點下的字典A_f及其導(dǎo)數(shù)矩陣A'_f、A”_f。

步驟8：基于y_f和A_f，根據(jù)公式(9)獲取各單獨頻率下的位置向量x_f以及聯(lián)合各頻率點數(shù)據(jù)求解的位置向量x_t，根據(jù)公式(10)進行房間反射系數(shù)估計，該步驟一般重復(fù)多次直到位置向量x_t及反射系數(shù)β不在變化。

步驟9：根據(jù)求解結(jié)果x_f、x_t中最大值對應(yīng)的網(wǎng)格點即為該頻率點下使用本方法定位得到的聲源位置。

本方法中反射系數(shù)在不同初值及網(wǎng)格劃分情況下迭代收斂的情況參照圖3、圖4。其中粗虛線是真實值，可以看到不同初值下10*10、20*20網(wǎng)格劃分情況下，利用本發(fā)明方法最終都能收斂到反射系數(shù)的真實值。