本發(fā)明涉及一種封閉空間聲源定位方法，特別涉及一種基于封閉空間幾何信息建模的單傳聲器聲源定位方法。

背景技術(shù)：

封閉空間內(nèi)的聲源定位問(wèn)題是聲學(xué)研究中的一項(xiàng)重要內(nèi)容，在人們生活中的很多方面都有較廣的應(yīng)用。例如飛機(jī)艙室的噪聲源定位問(wèn)題，會(huì)議室的智能話筒傳音控制系統(tǒng)，以及地下車(chē)庫(kù)的車(chē)位識(shí)別系統(tǒng)等。這類(lèi)空間的主要特點(diǎn)是環(huán)境復(fù)雜且混響程度較高，某些小尺度的空間還有嚴(yán)重的聲波干涉、衍射效應(yīng)，對(duì)聲源定位技術(shù)的精度、穩(wěn)定性要求較高。

文獻(xiàn)“基于傳聲器陣列的聲源定位，信息技術(shù)，2016(6):136-138”采用了基于聲達(dá)時(shí)間差的定位技術(shù)，這種聲源定位方法一般分為2個(gè)步驟進(jìn)行：先進(jìn)行聲達(dá)時(shí)間差估計(jì)，并從中獲取傳聲器陣列中陣元間的聲延遲；再利用獲取的聲達(dá)時(shí)間差，結(jié)合已知的傳聲器陣列的空間位置進(jìn)一步定出聲源的方向。這種方法相對(duì)于可控波束形成技術(shù)和基于高分辨率的譜估計(jì)技術(shù)的精度較好，計(jì)算量較小。

從算法和實(shí)際應(yīng)用來(lái)看，這種聲源定位技術(shù)還存在著一些缺陷。首先，要獲得較高的定位精度，需要大規(guī)模的傳聲器陣列，這在經(jīng)濟(jì)上會(huì)造成很大的負(fù)擔(dān)；對(duì)于基于聲達(dá)時(shí)間差的定位技術(shù)，聲源定位實(shí)際上只是定位出了聲源的方向，具體的距離則無(wú)法給出，另外，這種方法受環(huán)境混響影響很大，當(dāng)混響程度較高時(shí)，精度很低。同時(shí)該方法在定位階段用的參數(shù)已經(jīng)是對(duì)過(guò)去時(shí)間的估計(jì)，因此不是最優(yōu)估計(jì)。

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

為了克服現(xiàn)有封閉空間聲源定位方法實(shí)用性差的不足，本發(fā)明提供一種基于封閉空間幾何信息建模的單傳聲器聲源定位方法。該方法在封閉空間邊界所包圍的流體區(qū)域內(nèi)布置n個(gè)節(jié)點(diǎn)，將此n個(gè)節(jié)點(diǎn)從1到n進(jìn)行編號(hào)。利用此n個(gè)節(jié)點(diǎn)構(gòu)建出封閉空間及內(nèi)部物體的幾何形狀，并根據(jù)此n個(gè)節(jié)點(diǎn)坐標(biāo)信息建立描述封閉空間數(shù)學(xué)性質(zhì)的系統(tǒng)矩陣；在實(shí)際聲源定位時(shí)，利用單傳聲器測(cè)得多個(gè)頻率下的聲信號(hào)，將此聲信號(hào)與系統(tǒng)矩陣共同運(yùn)算，得到聲信號(hào)的位置信息，從而實(shí)現(xiàn)聲源定位，實(shí)用性好。

本發(fā)明解決其技術(shù)問(wèn)題所采用的技術(shù)方案：一種基于封閉空間幾何信息建模的單傳聲器聲源定位方法，其特點(diǎn)是包括以下步驟：

步驟一、在封閉空間邊界所包圍的流體區(qū)域內(nèi)布置n個(gè)節(jié)點(diǎn)，將此n個(gè)節(jié)點(diǎn)從1到n進(jìn)行編號(hào)。n個(gè)節(jié)點(diǎn)一方面作為封閉空間的幾何信息在后期進(jìn)行建模計(jì)算，一方面用來(lái)表示未來(lái)聲源定位時(shí)的位置。

n個(gè)節(jié)點(diǎn)近似均勻分布于整個(gè)流體區(qū)域內(nèi)，封閉空間的邊界處需布有節(jié)點(diǎn)。當(dāng)封閉空間內(nèi)部存在物體時(shí)，物體內(nèi)部不布置節(jié)點(diǎn)，且物體邊界上應(yīng)布置有節(jié)點(diǎn)，即節(jié)點(diǎn)應(yīng)勾勒出物體形狀。節(jié)點(diǎn)布置完畢后，封閉空間及其內(nèi)部物體的邊界用Γ表示，邊界所包圍的流體區(qū)域用Ω表示。

步驟二、假定在封閉空間內(nèi)聲源的位置為r，它在單位時(shí)間內(nèi)向單位體積內(nèi)的空間提供了ρ₀q(r,t)的媒質(zhì)質(zhì)量。根據(jù)質(zhì)量守恒定律，媒質(zhì)中聲波的連續(xù)方程寫(xiě)為：

式中，ρ'為媒質(zhì)密度增量，ρ₀表示媒質(zhì)靜態(tài)密度，q為q(r,t)的簡(jiǎn)寫(xiě)，v為媒質(zhì)質(zhì)點(diǎn)速度，t表示時(shí)間，div為散度算子，在三維空間笛卡爾坐標(biāo)系中，

除了連續(xù)性方程之外，用來(lái)描述媒質(zhì)聲波的基本方程還有兩個(gè)，它們不受聲源的影響，分別為運(yùn)動(dòng)方程：

和物態(tài)方程：

p＝c₀²ρ' (3)以上兩式中，p代表聲壓，c₀代表聲速，grad為梯度算子，在三維空間笛卡爾坐標(biāo)系中，

由媒質(zhì)中聲波的三個(gè)基本方程得到有源情況下封閉空間中有關(guān)聲壓p的波動(dòng)方程：

式中，為拉普拉斯算子，在三維空間笛卡爾坐標(biāo)系中，

在頻域，聲源強(qiáng)度q(r,t)表示為

q(r,t)＝q_ω(r)e^jωt (5)

式中，ω為諧振頻率，q_ω(r)為在位置r處頻域內(nèi)的聲源強(qiáng)度。

封閉空間聲場(chǎng)各點(diǎn)的聲壓的頻率與聲源相同，聲壓表示為：

p(r,t)＝p_ω(r)e^jωt (6)

式中，p_ω(r)為在位置r處頻域內(nèi)的聲壓。

將式(6)及式(7)代入式(5)中，得到簡(jiǎn)諧聲源激勵(lì)下的聲波波動(dòng)方程為：

令式中k＝ω/c₀，稱其為波數(shù)，并且消去e^jωt，得到只依賴于空間坐標(biāo)的那部分方程，即室內(nèi)有源Helmholtz方程：

這樣就將聲壓的時(shí)域問(wèn)題轉(zhuǎn)換為頻域問(wèn)題，式(9)即為封閉空間聲場(chǎng)的控制方程。

在封閉空間中，邊界具有吸聲能力，其聲壓梯度表示為：

式中，n為封閉空間壁面外法線方向，ζ稱為比聲阻抗，滿足下式：

式中，Z為界面聲阻抗。

根據(jù)Galerkin型加權(quán)殘量法，為了求解式(8)，首先設(shè)一試函數(shù)為代入有源Helmholtz方程及其邊界條件，試函數(shù)產(chǎn)生殘量R和

根據(jù)伽遼金法確定權(quán)函數(shù)，有

由格林第一公式

式(13)簡(jiǎn)化為

在聲場(chǎng)中任意一點(diǎn)的聲壓用各節(jié)點(diǎn)聲壓來(lái)表示，即

式中，N_i為節(jié)點(diǎn)i處的形函數(shù)，p_i為節(jié)點(diǎn)i處的聲壓。

將式(16)代入式(15)，得

式中，為形函數(shù)的導(dǎo)數(shù)矩陣，其表達(dá)式為：

整理式(17)，得到

令

其中，K稱為剛度矩陣，M稱為質(zhì)量矩陣，C稱為阻尼矩陣，G稱為載荷矩陣。當(dāng)聲源位于位置r₀(x₀,y₀,z₀)處時(shí)，頻域內(nèi)的聲源強(qiáng)度表示為：

q_ω(r)＝q_ωδ(r-r₀) (24)

其中

將式(24)代入式(19)中，得

最后，將式(20)、(21)、(22)、(23)代入式(19)并整理得到

式中，F(xiàn)＝j(luò)ρ₀N^Tq_ω。K、C、M均為n×n階的系數(shù)矩陣，各自的表達(dá)式分別為：M＝∫_ΩN^TNdΩ/c₀²，C＝∫_ΓN^TNdΓ/c₀ζ，各式中的N為形函數(shù)，在實(shí)際求解中，K、C、M中的積分運(yùn)算用求和運(yùn)算代替其中m為積分點(diǎn)的數(shù)量，m_b為邊界上積分點(diǎn)的數(shù)量，ξ_i為積分系數(shù)；ω為圓頻率；c₀為空氣中的聲速；ζ稱為比聲阻抗，滿足ζ＝Z/ρ₀c₀，ρ₀為空氣密度，Z為界面材料的聲阻抗；p為封閉空間內(nèi)任意位置處的聲壓，實(shí)際定位時(shí)為單傳聲器所測(cè)得的聲信號(hào)；F為n×1階的列向量，表示聲源相關(guān)信息，其表達(dá)式為F＝j(luò)ρ₀N^Tq_ω，其中q_ω表示聲源強(qiáng)度。

步驟三、根據(jù)節(jié)點(diǎn)坐標(biāo)，利用移動(dòng)最小二乘法，獲得步驟二中所涉及的形函數(shù)N。利用移動(dòng)最小二乘法構(gòu)建形函數(shù)。一個(gè)場(chǎng)函數(shù)u(x)在一點(diǎn)的近似值表示為：

其中是計(jì)算點(diǎn)x的鄰域范圍內(nèi)各節(jié)點(diǎn)的坐標(biāo)，為基函數(shù)向量，m為基函數(shù)的個(gè)數(shù)，a(x)＝[a₁(x),a₂(x),…a_m(x)]為待定系數(shù)向量。使用單項(xiàng)式基函數(shù)做運(yùn)算，在三維空間中常用的線性及二次單項(xiàng)式基函數(shù)分別為：

將求解域用節(jié)點(diǎn)離散后，在每個(gè)節(jié)點(diǎn)處定義一個(gè)權(quán)函數(shù)該函數(shù)只在支撐域內(nèi)不為零，在支撐域之外為零，在三維情況下，權(quán)函數(shù)的支撐域?yàn)榍蛐巍＿x定權(quán)函數(shù)后，就求得近似函數(shù)在節(jié)點(diǎn)處的誤差加權(quán)平方和：

令J取最小值，即

經(jīng)過(guò)整理后，得到下式：

A(x)a(x)＝B(x)u (32)

式中，A(x)，B(x)的含義為：

由公式(32)得到a(x)，將其代入式(28)得：

步驟四、在封閉空間內(nèi)部任意位置處設(shè)置一傳聲器。當(dāng)聲源有聲音發(fā)出時(shí)，傳聲器拾取到一段音頻信號(hào)f(t)，將其進(jìn)行傅里葉變換，得到此音頻信號(hào)的頻域信號(hào)F(ω)。

步驟五、將步驟四中得到的頻域信號(hào)F(ω)作為式(29)中的p值代入式(29)，求解后得到n×1階的列向量F，由于F中的j、ρ₀、q_ω均為常數(shù)，因此，求得的列向量中值最大的元素所代表的節(jié)點(diǎn)的位置即為聲源的位置。

本發(fā)明的有益效果是：該方法在封閉空間邊界所包圍的流體區(qū)域內(nèi)布置n個(gè)節(jié)點(diǎn)，將此n個(gè)節(jié)點(diǎn)從1到n進(jìn)行編號(hào)。利用此n個(gè)節(jié)點(diǎn)構(gòu)建出封閉空間及內(nèi)部物體的幾何形狀，并根據(jù)此n個(gè)節(jié)點(diǎn)坐標(biāo)信息建立描述封閉空間數(shù)學(xué)性質(zhì)的系統(tǒng)矩陣；在實(shí)際聲源定位時(shí)，利用單傳聲器測(cè)得多個(gè)頻率下的聲信號(hào)，將此聲信號(hào)與系統(tǒng)矩陣共同運(yùn)算，得到聲信號(hào)的位置信息，從而實(shí)現(xiàn)聲源定位，實(shí)用性好。本發(fā)明將幾何信息建模引入到封閉空間的聲源定位之中，僅需單傳聲器即可完成精度較高的定位。解決了背景技術(shù)聲源定位方法需要使用傳聲器陣列進(jìn)行聲源定位的技術(shù)問(wèn)題，僅僅利用單傳聲器測(cè)得多個(gè)頻率下的聲信號(hào)，將此聲信號(hào)與系統(tǒng)矩陣共同運(yùn)算，得到聲信號(hào)的位置信息，從而實(shí)現(xiàn)快速實(shí)時(shí)聲源定位。

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

附圖說(shuō)明

圖1是本發(fā)明基于封閉空間幾何信息建模的單傳聲器聲源定位方法的流程圖。

圖2是本發(fā)明方法所涉及的封閉空間及傳聲器、聲源示意圖。

圖3是本發(fā)明方法所涉及的聲源點(diǎn)及接收點(diǎn)設(shè)置示意圖。

具體實(shí)施方式

參照?qǐng)D1-3。本發(fā)明基于封閉空間幾何信息建模的單傳聲器聲源定位方法具體步驟如下：

實(shí)際封閉空間為一矩形空間，其長(zhǎng)寬高分別為0.9m、1.0m、1.2m，其邊界為亞克力玻璃結(jié)構(gòu)。坐標(biāo)原點(diǎn)設(shè)置在封閉空間的底角處，傳聲器吊裝在封閉空間頂部，其坐標(biāo)為(0.5m，0.5m，1m)，聲源位置坐標(biāo)為(0.7m，0.8m，0.3m)。

步驟1、建立包含封閉空間幾何信息的節(jié)點(diǎn)模型，由于此空間形狀規(guī)則，因此節(jié)點(diǎn)為均勻分布，每個(gè)坐標(biāo)軸向上有9個(gè)節(jié)點(diǎn)，共有729個(gè)節(jié)點(diǎn)。

步驟2、獲取計(jì)算聲源位置的系統(tǒng)方程。

步驟3、計(jì)算形函數(shù)N。

步驟4、在(0.5m，0.5m，1m)處吊裝單傳聲器。當(dāng)聲源有聲音發(fā)出時(shí)，傳聲器拾取到一段音頻信號(hào)f(t)，將其進(jìn)行傅里葉變換，得到此音頻信號(hào)的頻域轉(zhuǎn)換信號(hào)，為F(ω)。

步驟5、將步驟4中得到的頻域信號(hào)F(ω)作為式(29)中的p值代入式(29)，求解后得到n×1階的列向量F，由于F中的j、ρ₀、q_ω均為常數(shù)，列向量F實(shí)際由N控制。因此，求得的列向量F中值最大的元素所代表的節(jié)點(diǎn)的位置即為聲源的位置。

在本實(shí)施例中，對(duì)10個(gè)聲源位置分別利用不同的聲音作為測(cè)試信號(hào)進(jìn)行了定位測(cè)試，每個(gè)聲源位置測(cè)試10次。最終測(cè)試正確率為96％，有效地證明了本發(fā)明的有效性。

更具體的方法步驟如下：

步驟一、設(shè)有一封閉空間，在其內(nèi)部要實(shí)現(xiàn)聲源定位。在封閉空間邊界所包圍的流體區(qū)域內(nèi)布置n個(gè)節(jié)點(diǎn)，將此n個(gè)節(jié)點(diǎn)從1到n進(jìn)行編號(hào)。n個(gè)節(jié)點(diǎn)一方面作為封閉空間的幾何信息在后期進(jìn)行建模計(jì)算，一方面用來(lái)表示未來(lái)聲源定位時(shí)的位置。

節(jié)點(diǎn)布置方式?jīng)]有嚴(yán)格規(guī)定，但是應(yīng)近似均勻分布于整個(gè)區(qū)域內(nèi)，封閉空間的邊界處需布有節(jié)點(diǎn)。當(dāng)封閉空間內(nèi)部存在物體時(shí)，物體內(nèi)部不布置節(jié)點(diǎn)，且物體邊界上應(yīng)布置有節(jié)點(diǎn)，即節(jié)點(diǎn)應(yīng)勾勒出物體形狀。節(jié)點(diǎn)布置完畢后，封閉空間及其內(nèi)部物體的邊界用Γ表示，邊界所包圍的流體區(qū)域用Ω表示。

步驟二、獲取計(jì)算聲源位置的系統(tǒng)方程。假定在封閉空間內(nèi)聲源的位置為r，它在單位時(shí)間內(nèi)向單位體積內(nèi)的空間提供了ρ₀q(r,t)的媒質(zhì)質(zhì)量。根據(jù)質(zhì)量守恒定律，媒質(zhì)中聲波的連續(xù)方程寫(xiě)為：

式中ρ'為媒質(zhì)密度增量，ρ₀表示媒質(zhì)靜態(tài)密度，q為q(r,t)的簡(jiǎn)寫(xiě)，v為媒質(zhì)質(zhì)點(diǎn)速度，t表示時(shí)間，div為散度算子，在三維空間笛卡爾坐標(biāo)系中，

除了連續(xù)性方程之外，用來(lái)描述媒質(zhì)聲波的基本方程還有兩個(gè)，它們不受聲源的影響，分別為運(yùn)動(dòng)方程：

和物態(tài)方程：

p＝c₀²ρ' (3)

上兩式中，p代表聲壓，c₀代表聲速，grad為梯度算子，在三維空間笛卡爾坐標(biāo)系中，

利用與推導(dǎo)無(wú)源波動(dòng)方程相類(lèi)似的方法，由媒質(zhì)中聲波的三個(gè)基本方程得到有源情況下封閉空間中有關(guān)聲壓p的波動(dòng)方程:

式中為拉普拉斯算子，在三維空間笛卡爾坐標(biāo)系中，

在頻域，聲源強(qiáng)度q(r,t)表示為

q(r,t)＝q_ω(r)e^jωt (5)

式中ω為諧振頻率，q_ω(r)為在位置r處頻域內(nèi)的聲源強(qiáng)度。

由于通常將封閉空間聲場(chǎng)作為線性系統(tǒng)考慮，因此空間內(nèi)部各點(diǎn)的聲壓的頻率與聲源相同，聲壓表示為：

p(r,t)＝p_ω(r)e^jωt (6)

式中p_ω(r)為在位置r處頻域內(nèi)的聲壓。

將式(6)及式(7)代入式(5)中，得到簡(jiǎn)諧聲源激勵(lì)下的聲波波動(dòng)方程為：

令式中k＝ω/c₀，稱其為波數(shù)，并且消去e^jωt，得到只依賴于空間坐標(biāo)的那部分方程，即室內(nèi)有源Helmholtz方程：

這樣就將聲壓的時(shí)域問(wèn)題轉(zhuǎn)換為頻域問(wèn)題，式(9)即為封閉空間聲場(chǎng)的控制方程。

在封閉空間中，邊界通常具有一定的吸聲能力，其聲壓梯度表示為：

式中n為封閉空間壁面外法線方向，ζ稱為比聲阻抗，滿足下式：

式中Z為界面聲阻抗，查表得到。

根據(jù)Galerkin型加權(quán)殘量法，為了求解式(8)，首先設(shè)一試函數(shù)為代入有源Helmholtz方程及其邊界條件，由于試函數(shù)通常不是精確解，因此將產(chǎn)生殘量R和

根據(jù)伽遼金法確定權(quán)函數(shù)，有

由格林第一公式

式(13)簡(jiǎn)化為

在聲場(chǎng)中任意一點(diǎn)的聲壓用各節(jié)點(diǎn)聲壓來(lái)表示，即

式中N_i為節(jié)點(diǎn)i處的形函數(shù)，p_i為節(jié)點(diǎn)i處的聲壓。

將式(16)代入式(15)，得

式中為形函數(shù)的導(dǎo)數(shù)矩陣，其表達(dá)式為：

整理式(17)，得到

令

∫_Ωjρ₀ωN^Tq_ωdΩ＝G (23)

其中K稱為剛度矩陣，M稱為質(zhì)量矩陣，C稱為阻尼矩陣，G稱為載荷矩陣。當(dāng)聲源位于某一特定位置r₀(x₀,y₀,z₀)處時(shí)，頻域內(nèi)的聲源強(qiáng)度表示為：

q_ω(r)＝q_ωδ(r-r₀) (24)

其中

將式(24)代入式(19)中，得

G＝∫_Ω-jρ₀ωq_ωδ(r-r₀)N^Tdv＝-jρ₀ωq_ωN^T (26)

最后，將式(20)、(21)、(22)、(23)代入式(19)并整理得到

式中F＝j(luò)ρ₀N^Tq_ω。

式中K、C、M均為n×n階的系數(shù)矩陣，各自的表達(dá)式分別為：M＝∫_ΩN^TNdΩ/c₀²，C＝∫_ΓN^TNdΓ/c₀ζ，各式中的N為形函數(shù)，在實(shí)際求解中，K、C、M中的積分運(yùn)算用求和運(yùn)算代替其中m為積分點(diǎn)的數(shù)量，m_b為邊界上積分點(diǎn)的數(shù)量，ξ_i為積分系數(shù)；ω為圓頻率；c₀為空氣中的聲速；ζ稱為比聲阻抗，滿足ζ＝Z/ρ₀c₀，ρ₀為空氣密度，Z為界面材料的聲阻抗，查表得到；p為封閉空間內(nèi)任意位置處的聲壓，實(shí)際定位時(shí)為單傳聲器所測(cè)得的聲信號(hào)；F為n×1階的列向量，表示聲源相關(guān)信息，其表達(dá)式為F＝j(luò)ρ₀N^Tq_ω，式中q_ω表示聲源強(qiáng)度。

步驟三、計(jì)算形函數(shù)。根據(jù)節(jié)點(diǎn)坐標(biāo)，利用移動(dòng)最小二乘法，獲得步驟2中所涉及的形函數(shù)N。利用移動(dòng)最小二乘法構(gòu)建形函數(shù)。一個(gè)場(chǎng)函數(shù)u(x)在一點(diǎn)的近似值表示為：

其中是計(jì)算點(diǎn)x的鄰域范圍內(nèi)各節(jié)點(diǎn)的坐標(biāo)，為基函數(shù)向量，m為基函數(shù)的個(gè)數(shù)，a(x)＝[a₁(x),a₂(x),…a_m(x)]為待定系數(shù)向量。通常使用單項(xiàng)式基函數(shù)做運(yùn)算，在三維空間中常用的線性及二次單項(xiàng)式基函數(shù)分別為：

將求解域用節(jié)點(diǎn)離散后，在每個(gè)節(jié)點(diǎn)處定義一個(gè)權(quán)函數(shù)該函數(shù)只在一個(gè)有限區(qū)域(支撐域)內(nèi)不為零，在區(qū)域之外為零，在三維情況下，權(quán)函數(shù)的支撐域通常為球形。常用的權(quán)函數(shù)有高斯函數(shù)、樣條函數(shù)等。選定權(quán)函數(shù)后，就求得近似函數(shù)在節(jié)點(diǎn)處的誤差加權(quán)平方和：

令J取最小值，即

經(jīng)過(guò)整理后，得到下式：

A(x)a(x)＝B(x)u (32)

式中A(x)，B(x)的含義為：

由公式(32)得到a(x)，將其代入式(28)得：

步驟四、在封閉空間內(nèi)部任意位置處設(shè)置一傳聲器。當(dāng)聲源有聲音發(fā)出時(shí)，傳聲器拾取到一段音頻信號(hào)f(t)，將其進(jìn)行傅里葉變換，得到此音頻信號(hào)的頻域轉(zhuǎn)換信號(hào)，為F(ω)。

步驟五、將步驟四中得到的頻域信號(hào)F(ω)作為式(29)中的p值代入式(29)，求解后得到n×1階的列向量F，由于F中的j、ρ₀、q_ω均為常數(shù)，因此，求得的列向量中值最大的元素所代表的節(jié)點(diǎn)的位置即為聲源的位置。