基于粒子濾波的波疊加算法

孫玲莉 楊超 郭輝 胡定玉

摘 要： 與傳統波疊加法通過逆問題求解獲得源強不同，提出一種基于粒子濾波的波疊加算法來實現正向求解源強。該方法基于波疊加理論建立狀態空間模型，根據初始化粒子的先驗信息預測并更新狀態向量。通過權重計算和重采樣，估計等效源的強度及其位置，進而重建三維輻射聲場。仿真分析與實驗驗證表明，該算法可避免逆問題求解中的不適定性難題，實現聲場的高精度重建。

關鍵詞： 粒子濾波; 波疊加; 聲源定位; 聲場重建

文章編號： 2095-2163（2021）03-0080-05 中圖分類號：TP391.9; TB52 文獻標志碼：A

【Abstract】Different from the traditional wave superposition method which is utilized to obtain the source intensity by solving the inverse problem， an innovative wave superposition algorithm based on particle filtering is proposed to realize the forward solving of source intensity. This method presents a state space model based on the wave superposition theory， then updates the state vector according to the prior information prediction of initializing particle. Through weight calculation and resampling， the intensity and the position of the equivalent source can be estimated. Thus， the three-dimensional sound field can be reconstructed. Simulation analysis and experimental verification show that the algorithm can avoid the ill-posed problem of solving the inverse problem and achieve high-precision reconstruction of the sound field.

【Key words】 particle filter; wave superposition; sound source localization; sound field reconstruction

0 引 言

波疊加法由Koopmann等人[1]在1989年提出，是一種求解聲輻射問題的積分方程法，其原理為：任何物體輻射的聲場均可以等效為該輻射體內部若干個不同源強的等效源所產生聲場的疊加。利用波疊加法近場聲全息進行聲場重建，既能克服基于空間Fourier變換法的近場聲全息[2]只能計算規則形狀聲源的缺點，又能避免基于邊界元法的近場聲全息[3]所帶來的復雜插值運算和奇異積分處理困擾，國內外學者對該方法進行了深入研究[4-8]。

利用波疊加法進行聲場重建時，其關鍵是獲取等效源強度和配置等效源位置。通常采用Tikhonov正則化方法[9-10]求解聲學逆問題中的病態方程，獲取準確的等效源強度，再合理配置等效源位置，進行聲場重建。Antoni[11]提出了貝葉斯正則化準則的基本理論，Pereira等人[12]進一步通過理論和實驗對比了等效源法和貝葉斯法。Chardon等人[13]針對傳統近場聲全息逆問題正則化方法的復雜性，提出了一種新的基于壓縮采樣的正則化方法。Yang等人[14]提出了一種基于卡爾曼濾波的混合局部近場聲全息技術，旨在解決等效源配置不確定性問題。

以上所提方法均是從求解聲學逆問題的角度出發獲得源強，為避免求解逆問題的不適定性，Bai等人[15-16]提出基于等效源模型的自適應濾波方法，對麥克風陣列信號進行處理，該方法利用遞歸方式估計未知源振幅，屬于聲學正問題。本文提出一種基于粒子濾波的波疊加方法（Particle Filter based WSM，PF-WSM），該方法以波疊加法為理論基礎，建立狀態空間模型，利用粒子濾波算法從正向估計等效源強度和位置，根據估計得到參數計算重建面聲壓值。本文安排如下：首先對PF-WSM算法進行理論描述，再仿真分析聲源位置已知和未知兩種情況的聲場重建結果，最后通過實驗驗證，證明了該方法的可行性。

1 理論背景

1.1 PF-WSM算法理論

傳統波疊加法通過Tikhonov正則化方法求解聲學逆問題中的病態方程，獲得準確的等效源強度，故被稱為Tikhonov正則化波疊加（TR-WSM）算法。與TR-WSM算法反解源強相比，PF-WSM算法正向估計源強，避免了求解逆問題的不適定性。

粒子濾波（Particle Filter，PF）是一種基于蒙特卡洛方法進行遞推貝葉斯估計的濾波算法[17-18]。通過采用一組具有相應權值的隨機粒子估計狀態的前驗分布，并通過觀測值逼近狀態的后驗分布，再求出這組粒子的加權求和作為狀態的估計值。

在PF-WSM算法中，用等效源強度和位置描述系統狀態，源強與傳遞矩陣計算的聲壓描述觀測狀態，由此，狀態空間模型表示為：

其中，qn表示源強的系統狀態矢量;sn表示源位置的系統狀態矢量且不隨遞歸變化;yn表示觀測的聲壓矢量;n表示估計狀態遞歸次數。源強的系統噪聲和觀測噪聲服從均值為零、方差分別為σ21和σ23的正態分布，即v1～N0，σ21和v3～N0，σ23;源位置的系統噪聲服從均勻分布，即v2～U-σ2，σ2。G表示等效源到全息面的傳遞矩陣。

第i個粒子的權值更新方程表示為：

當經過一定次數的估計狀態遞歸后，會存在粒子退化問題，可通過重采樣[19]方法來解決。該方法是對后驗概率密度的離散近似表示再進行一次采樣，繁殖權值較高的樣本，淘汰權值較低的樣本，重新生成一個新樣本集合，以克服樣本權值退化問題。

1.2 PF-WSM算法流程

粒子濾波的核心是利用優選隨機粒子來實現狀態分布的估計，為更好地進行源強分布估計，給定粒子初始狀態分布為高斯分布，這是因為根據中心極限定理，隨著粒子數量的增加，大多數粒子均值變得相似并越來越接近總體均值。均值為μ和方差為σ2的高斯分布的概率密度函數為：

對于復數域的源強概率密度函數，式（7）可以看作是一個聯合概率密度函數的實部和虛部，即：

利用PF-WSM算法在2種情況下進行聲場重建對比，根據公式（1）～公式（3）建立狀態空間模型。對這2種情況可簡述為：

（1）聲源位置已知，其觀測方程中聲壓值的變化僅與粒子的源強有關。

（2）聲源位置未知，聲壓值的變化依賴于粒子的源強和源位置。

文中給出了PF-WSM算法的流程如圖1所示。具體過程如下：

（1）粒子初始化：第一次遞歸時，從高斯分布Nμ，σ2中隨機產生K個粒子qini=1，2，…，K，并平均分配權重。

（2）系統狀態矢量預測：根據式（1）～式（2）所示的狀態轉移方程，計算粒子qin在n+1次遞歸的預測粒子qin+1，服從一階馬爾可夫過程。

（3）粒子權重更新：利用式（3）計算n+1次粒子預測值對應的觀測聲壓，同時計算觀測聲壓與n+1次的實際測量聲壓之間的差值，根據差值更新每個粒子的權重（式（4）～式（5）），再將權重進行歸一化。

（4）源強分布估計：利用當前遞歸次數中每個粒子的權重win+1以及取值qin+1，估計系統狀態矢量（源強分布），計算公式為：

（5）重采樣：為了快速有效地選擇接近真實源強分布的粒子，在進行下一次系統狀態估計之前，需要按照粒子的權重對粒子進行重采樣。重復上述步驟，即可基于粒子濾波算法構建系統狀態估計器，實現源強分布估計。

（6）聲場重建：根據上述構建的系統狀態估計器，得到合理的估計源強分布，利用源強與重建面之間的傳遞矩陣求得重建聲壓：

其中，H表示等效源面到重建面的傳遞矩陣，Q表示PF-WSM算法估計源強分布。

2 數值仿真

2.1 仿真設置

為驗證PF-WSM算法在聲場重建中的可行性，將其與Tikhonov正則化波疊加 （TR-WSM）算法進行比較。單極子源的半徑為0.01 m，表面脈動速度為0.2 m/s，由25個單極子聲源組成的輻射聲場布置在xoy面，以聲源面中心為坐標原點o，其法面方向為z軸，建立空間直角坐標系。全息面平行于xoy面，距坐標原點0.2 m，在0.8×0.8 m2的平面上均勻布置25個傳聲器，間距為0.2 m。重建面尺寸為1×1 m2，網格間距為0.1 m，距坐標原點0.05 m。對全息面的測量數據加上隨機高斯白噪聲，假設聲場中的信噪比為20 dB。設置遞歸次數N為90次，粒子數K為100個。

仿真過程中，相對誤差的計算公式如下：

2.2 仿真分析

聲源位置已知時，假設25個單極子源均勻分布在聲源面上，如圖2（a）所示，其源位置與全息面上傳聲器位置在xoy平面的俯視圖上一一對應。根據1.2節的算法流程計算出重建聲壓，圖2（b）比較了100個粒子的PF-WSM算法和TR-WSM算法聲壓重建誤差。源位置均勻分布在網格點上，PF-WSM算法基于源強為高斯分布的先驗信息僅估計源強，與TR-WSM算法相比，PF-WSM算法在聲壓重建中誤差較大。

聲源準確位置已知的情況下，TR-WSM算法的性能要優于PF-WSM算法。但在實際場景中，通常無法預先知道聲源的準確位置，對源位置假設不正確，導致聲場重建值偏離實際值，因此假設聲源位置已知是一種理想化的情況。由此，假設源位置是未知的，并且存在位置擾動，更符合真實情況。

聲源位置未知時，假設25個單極子源分布在聲源面上存在擾動，源位置相對于傳聲器分布位置如圖3（a）所示。圖3（b）是PF-WSM算法和TR-WSM算法的重建誤差對比。在源位置未知且存在擾動的情況下，PF-WSM算法同時估計源強和源位置，基于源強服從高斯分布和源位置服從均勻分布的先驗信息，在不斷遞歸之后，源位置越接近聲源真實位置，源強越接近真實源強，PF-WSM算法的聲場重建結果越接近真實輻射聲場。而TR-WSM算法在源位置未知的情況下，由于等效源配置不準確，重建誤差增大。圖3（b）中PF-WSM算法的重建誤差小于TR-WSM算法，說明PF-WSM算法在源強和源位置均未知的情況下，估計效果更好，而TR-WSM算法在未知源位置時的重建效果較差，其應用有一定局限性。圖4給出了2種算法在100 Hz時的聲場重建效果對比，圖4（a）為理論聲壓，圖4（b）為PF-WSM算法的重建聲壓，圖4（c）為TR-WSM算法的重建聲壓。由圖4可以看出，PF-WSM算法的聲場重建效果更好，其熱點最大值與理論聲壓相同，說明粒子濾波可以較準確地估計出源強和源位置，而TR-WSM算法的熱點最大值小于理論值，這是因為不知道聲源位置時，等效源配置不準確，影響聲場重建效果。

3 實驗驗證

為進一步檢驗算法的可行性，在半消聲室進行實驗驗證。實驗設置如圖5（a）所示。半消聲室尺寸為9.8 m×8.6 m×3.5 m，背景噪聲為18 dB（A），截止頻率為125 Hz。采用音箱模擬噪聲源，頻率為300 Hz，音箱紙盆中心位置為（0， 0， 0）。以音箱紙盆中心為坐標原點，音箱平面法向為z軸建立空間直角坐標系。傳聲器陣列平面平行于聲源面，陣列中心距坐標原點0.6 m，離地面高度為1.2 m。傳聲器陣列如圖5（b）所示，采用36個傳聲器組成3個同心圓的面陣列，3個同心圓半徑分別為：0.07 m、0.18 m、0.25 m。重建面大小與傳聲器陣列面相同，重建面距坐標原點0.05 m。PF-WSM算法中，粒子的先驗信息是基于高斯分布，36個虛源均勻分布在0.25 m×0.25 m的正方形平面上，間隔為0.05 m，距坐標原點-0.01 m。將PF-WSM算法和TR-WSM算法重建的聲壓與傳聲器測得的聲壓進行對比。圖6為實驗驗證結果。圖6（a）為傳聲器陣列測得的聲壓分布，圖6（b）和圖6（c）分別為PF-WSM算法和TR-WSM算法重建的聲壓分布。由式（11）計算PF-WSM算法和TR-WSM算法重建誤差，分別為5.89%和13.43%，可以得出結果：PF-WSM算法的重建效果優于TR-WSM算法的重建效果。

4 結束語

本文提出一種基于粒子濾波的波疊加算法進行聲場重建。該方法在源位置未知且隨機的情況下，能夠有效利用參數的先驗知識，通過粒子濾波估計器，達到求解源強的效果。特別是在聲源位置未知且存在擾動的情況下，源位置是非高斯分布的，傳統波疊加法由于等效源配置不準確而聲場重建效果不好，PF-WSM算法則利用這一特點，可以更好地重建聲場。PF-WSM算法不僅可以估計源強，還可以估計源位置，甚至在源位置未知的情況下，估計性能更好，即對源位置的估計越準確，重建的聲壓就越準確。另一方面，粒子數的增加直接影響估計狀態向量的精度，粒子數越多，估計狀態越接近真實狀態。本文在權衡運算時間和估計精度之后，選擇100個粒子進行估計，若增加粒子數（如1 000個粒子），估計值將越接近真實值，但運算時間相應增加，后續可采用改進粒子濾波算法進行優化，降低采樣復雜度，提升運算速度。

參考文獻

[1] ??KOOPMANN G H， SONG L， FAHNLINE J B. A method for computing acoustic fields based on the principle of wave superposition[J]. The Journal of the Acoustical Society of America， 1989， 86（6）： 2433-2438.

[2] HALD J. Patch near-field holography using a new statistically optimal method[J]. Proceedings of Inter-noise， 2003（975）： 2203-2210.

[3] UCHIDA H，SAIJYOU K. Data extrapolation method for boundary element method-based near-field acoustical holography[J].The Journal of the Acoustical Society of America， 2004， 115（2）： 785-796.

[4] 楊殿閣， 李兵， 王子騰， 等. 運動聲源識別的動態波疊加方法研究[J]. 物理學報， 2012， 61（5）：054306.

[5] FERNANDEZ-GRANDE E， XENAKI A， GERATOFT P. A sparse equivalent source method for near-field acoustic holography[J]. The Journal of the Acoustical Society of America， 2017， 141（1）： 532-542.

[6] 胡定玉， 李再幃， 方宇. 非自由聲場中目標聲場還原與重建的等效源方法[J]. 聲學學報， 2017，42（4）： 465-475.

[7] XU Zhongming， WANG Qinghua， HE Yangsong， et al. A monotonic two-step iterative shrinkage/thresholding algorithm for sound source identification based on equivalent source method[J]. Applied Acoustics， 2018， 129： 386-396.

[8] BI Chuanxing， LIU Yuan， ZHANG Yongbin， et al. Extension of sound field separation technique based on the equivalent source method in a sparsity framework[J]. Journal of Sound and Vibration， 2018， 442： 125-137.

[9] TIKHONOV A N. Solution of incorrectly formulated problems and the regularization method [J]. Soviet Math Dokl， 1963， 4： 1035-1038.

[10]LECLERE Q. Acoustic imaging using under-determined inverse approaches： Frequency limitations and optimal regularization[J]. Journal of Sound and Vibration， 2009， 321（3-5）： 605-619.

[11]ANTONIO J. Bayesian focusing：A unified approach to inverse acoustic radiation[C]//Proceedings of ISMA 2010. ?Leuven， Belgium： Katholieke Universiteit Leuven， 2010： 35-45.

[12]PEREIRA A， LECLRE Q， ANTONI J. A theoretical and experimental comparison of the equivalent source method and a Bayesian approach to noise source identification[C]//BeBec-2012，2012： 4.

[13]CHARDON G， DAUDET L， PEILLOT A， et al. Nearfield Acoustic Holography using sparsity and compressive sampling principles[J]. Journal of the Acoustical Society of America， 2012， 132（3）： 1521-1534.

[14]YANG C， WANG Y S， Guo H. Hybrid patch near-field acoustic holography based on Kalman filter[J]. Applied Acoustics， 2019， 148： 23-33.

[15]BAI M R， CHEN C C. Kalman filter-based microphone array signal processing using the equivalent source model[J]. Journal of Sound and Vibration， 2012， 331（22）： 4940-4955.

[16]BAI M R， AGARWAL A， CHEN C C， et al. Bayesian approach of nearfield acoustic reconstruction with particle filters[J].The Journal of the Acoustical Society of America， 2013， 133（6）： 4032-4043.

[17]王法勝， 魯明羽， 趙清杰， 等. 粒子濾波算法[J]. 計算機學報， 2014， 37（8）： 1679-1694.

[18]昝孟恩， 周航， 韓丹， 等. 粒子濾波目標跟蹤算法綜述[J]. 計算機工程與應用， 2019， 55（5）： 8-17.

[19]孔紅山， 李小鵬， 郁濱. SIR粒子濾波的改進算法[J]. 計算機工程與設計， 2020，41（7）： 1899-1904.