改進(jìn)粒子濾波在海底地聲參數(shù)反演中的應(yīng)用

吳偉文，任群言，鹿力成，馬 力

(1.中國科學(xué)院聲學(xué)研究所中國科學(xué)院水聲環(huán)境特性重點(diǎn)實(shí)驗(yàn)室，北京 100190；2.中國科學(xué)院大學(xué)，北京 100049)

0 引言

聲信號在水中具有良好的傳播性能，因此常用來進(jìn)行水下目標(biāo)探測或信息傳遞。淺海環(huán)境中，聲傳播受海底的影響較大，準(zhǔn)確的地聲參數(shù)有助于分析淺海環(huán)境中聲傳播特性，因此地聲參數(shù)估計是海洋聲學(xué)中的一個基礎(chǔ)研究問題。

地聲參數(shù)可以通過直接測量法[1]獲得，但此方法時間成本和經(jīng)濟(jì)代價都很高。反演算法如匹配場處理[2]是比直接測量法更有效率的方式。匹配場處理將地聲參數(shù)反演問題作為優(yōu)化問題，其目的是找到一個模型，此模型包含一組地聲參數(shù)(假設(shè)模型的真實(shí)地聲參數(shù)向量包含在參數(shù)集合中)，使前向模型產(chǎn)生的數(shù)據(jù)與觀測數(shù)據(jù)達(dá)到最優(yōu)匹配。

匹配場處理可分為兩步。(1)選擇目標(biāo)函數(shù)，對觀測數(shù)據(jù)與前向模型生成數(shù)據(jù)之間的誤差進(jìn)行衡量；(2)基于最小化目標(biāo)函數(shù)的準(zhǔn)則搜索最優(yōu)的參數(shù)向量。目標(biāo)函數(shù)中待匹配的數(shù)據(jù)，一般可選擇為對地聲參數(shù)變化較敏感的物理量。例如，海底反射系數(shù)隨沉積層分層特性、聲速、密度的變化而明顯變化，因此海底反射系數(shù)[3]可以作為目標(biāo)函數(shù)中待匹配的物理量。本文將海底反射系數(shù)作為待匹配的物理量，利用優(yōu)化算法搜索地聲參數(shù)向量。由于地聲反演問題具有維度高，局部最優(yōu)多，前向模型非線性等特點(diǎn)，所以參數(shù)搜索方式一般選遺傳算法[4]、模擬退火[4]、MCMC[5]采樣等全局算法。上述優(yōu)化算法一般稱之為批處理反演算法，不僅搜索過程費(fèi)時而且一批數(shù)據(jù)只能生成一個反演結(jié)果，反演結(jié)果之間互相獨(dú)立，無法高效處理環(huán)境隨距離變化時的參數(shù)估計問題。

為了利用相鄰兩段距離之間地聲參數(shù)的聯(lián)系，減少運(yùn)算量，可以將序貫尋優(yōu)算法如粒子濾波應(yīng)用到地聲參數(shù)的反演中。在地聲參數(shù)隨距離變化緩慢的情況下粒子濾波算法已經(jīng)取得了良好的效果[6]。但是在地聲參數(shù)隨距離變化劇烈，比如沉積層分層個數(shù)改變的情況下，傳統(tǒng)粒子濾波可能會失效。為了解決此問題，Dettmer等[7]將帶橋接的粒子濾波應(yīng)用到地聲參數(shù)反演中，通過橋接過渡劇烈變化的目標(biāo)分布，以實(shí)現(xiàn)地聲參數(shù)隨距離變化劇烈時的反演，并通過仿真數(shù)據(jù)證明了此方法的有效性。本文在粒子濾波中增加橋接的同時改進(jìn)了粒子的采樣方式，由傳統(tǒng)的固定方差采樣改進(jìn)為隨機(jī)方差采樣，同時對沉積層厚度的采樣采用高維與低維混合采樣的方式。將此方法應(yīng)用到地聲參數(shù)反演中，通過仿真分析得出，此方法在地聲參數(shù)隨距離變化劇烈的情況下仍可得到較穩(wěn)定的反演結(jié)果。

本文利用反射系數(shù)構(gòu)造似然函數(shù)，分析其對沉積層厚度、聲速、密度的敏感性；給出了改進(jìn)粒子濾波的算法流程；然后進(jìn)行仿真數(shù)據(jù)處理；得到改進(jìn)粒子濾波在地聲參數(shù)隨距離變化劇烈的情況下仍可得到較穩(wěn)定的反演結(jié)果的結(jié)論。

1 海底反射系數(shù)

1.1 反射系數(shù)計算

海底反射系數(shù)對沉積層厚度以及聲速都較為敏感，因此常作為地聲參數(shù)反演中待匹配的物理量。反射系數(shù)R為反射到海水的所有反射波與入射波的比值，基于射線聲學(xué)理論，平面波反射系數(shù)計算原理如圖1所示[8]，R的計劃公式為

圖1 沉積層數(shù)為1時的海底反射圖Fig.1 Bottom reflection diagram when the number of sediment layer is 1

其中：ρ1、ρ2、ρ3分別為海水、沉積層和基底的密度；c1、c2、c3分別為海水、沉積層和基底的聲速(當(dāng)考慮海底聲衰減系數(shù)α2、α3時，海底聲速為復(fù)數(shù))；h2為沉積層厚度；R12為聲波由海水入射到沉積層的反射系數(shù)，T12為聲波由海水入射到沉積層的透射系數(shù)，R21為聲波由沉積層入射到海水的反射系數(shù)，T21為聲波由沉積層入射到海水的透射系數(shù)，R23為聲波由沉積層入射到基底的反射系數(shù)，T23為聲波由沉積層入射到基底的透射系數(shù)，R32為聲波由基底入射到沉積層的反射系數(shù)，T23為聲波由基底入射到沉積層的透射系數(shù)。相位延遲為

其中：k2=2πf/c2，f為聲源頻率,θ2為折射波方向(已知入射波方向θ1，聲速c1和c2，可根據(jù)斯涅爾定律計算θ2)。根據(jù)R12=-R21和T12T21=1-R212，反射系數(shù)R進(jìn)一步可表示為

多層沉積層的反射系數(shù)可以依公式(3)遞推獲得，也可通過OASR工具箱[9]計算獲得。

1.2 敏感性分析

反射系數(shù)隨沉積層層厚，聲速，密度的變化而變化，本文用反射系數(shù)構(gòu)造如式(4)所示的似然函數(shù)：

其中：xt為第t幀的地聲參數(shù)向量，包含沉積層深度、聲速、密度、衰減系數(shù)；yt為第t時刻的觀測值，對應(yīng)于反演問題中的第t幀的反射系數(shù)觀測向量(包含N個掠射角)；g為測量函數(shù)。研究似然函數(shù)對沉積層聲速、密度、衰減系數(shù)、厚度的敏感性，對反演的可行性以及反演結(jié)果分析具有重要作用。下面在四層沉積層的情況下，對似然函數(shù)進(jìn)行敏感性分析。

令聲源頻率分別為500 Hz和1 000 Hz，仿真參數(shù)如表1所示。

表1 仿真信道參數(shù)Table 1 Parameters of the simulation channel

圖2 Φ隨ci、ρi、αi、hi各參數(shù)變化的曲線圖Fig.2Variation curves of Φ with parameters ci、ρi、αi、hi

可以清楚地看到，Φ對聲速和沉積層深度敏感，對密度、衰減系數(shù)不敏感。同時可以看到Φ隨h變化時，局部最優(yōu)明顯。為了防止算法在對沉積層深度的搜索時陷入局部最優(yōu)，本文提出了沉積層厚度混合搜索方式，此方法將在第2.2節(jié)介紹。

2 改進(jìn)粒子濾波

2.1 帶橋接的粒子濾波

在隨距離變化的地聲反演問題中，可以將地聲參數(shù)分為很多幀。本文的幀是指一小段范圍，在此范圍內(nèi)地聲參數(shù)可以近似認(rèn)為是不變的。由于每一幀之間地聲參數(shù)具有相關(guān)性，所以可以將粒子濾波等序貫尋優(yōu)方式應(yīng)用到地聲反演中。根據(jù)粒子濾波原理建立狀態(tài)方程和觀測方程分別為

式中：xt是第t時刻的狀態(tài)，在地聲反演問題中，即為第t幀的地聲參數(shù)向量，包含沉積層厚度、聲速、密度、衰減系數(shù)；xt-1為上一時刻的狀態(tài)，對應(yīng)于反演問題中上一幀的地聲參數(shù)向量；yt為第t時刻的觀測值，對應(yīng)于反演問題中的第t幀的反射系數(shù)觀測向量(包含N個掠射角)；f為狀態(tài)轉(zhuǎn)移函數(shù)，g為測量函數(shù)；g(xt)為以xt為地聲參數(shù)向量，代入反射系數(shù)計算公式，所求得的反射系數(shù)向量；nt和vt分別是狀態(tài)方程與觀測方程的噪聲，vt一般為高斯噪聲，則yt即為以g(xt)為均值的高斯分布。假設(shè)各掠射角的高斯噪聲獨(dú)立，似然函數(shù)L可以表示為

粒子濾波可能會由于似然變化較大而出現(xiàn)粒子枯竭的情況。為了應(yīng)對粒子枯竭，可以采用重采樣的方式，重采樣是指維持總粒子數(shù)量不變的情況下，復(fù)制權(quán)重大的粒子，同時舍去權(quán)重非常小的粒子，是一種有效防止粒子枯竭的方法。

將粒子濾波應(yīng)用到地聲參數(shù)反演中，在地聲參數(shù)隨距離緩慢變化時，傳統(tǒng)粒子濾波具有良好的效果。但當(dāng)?shù)芈晠?shù)隨距離劇烈變化時，似然函數(shù)變化劇烈，即使使用了重采樣，傳統(tǒng)粒子濾波仍可能會出現(xiàn)粒子枯竭[11]的情況，無法正確收斂到目標(biāo)分布。為了使粒子濾波可以在地聲參數(shù)劇烈變化的情況下收斂到目標(biāo)分布，有文獻(xiàn)使用了帶橋接的粒子濾波的方法。使用橋接可以使粒子濾波結(jié)果緩慢收斂到目標(biāo)分布，實(shí)現(xiàn)目標(biāo)分布之間的平滑過渡。帶橋接的粒子濾波原理為[12]

其中，1＞α1＞α2…＞αM＞0，M為橋接層數(shù)。橋接的作用是使目標(biāo)分布從逐漸收斂到，即為重要性采樣概率密度，當(dāng)m=1時，αM+1=1，。當(dāng)m=0時，，即為最終的目標(biāo)分布。通過橋接，目標(biāo)分布從q(x1:t)逐漸收斂到，即最終收斂到參數(shù)后驗(yàn)概率密度。橋接原理如圖3所示，每一步橋接中都使用重采樣從q(x1:t)平緩過渡到，以防止粒子枯竭。

圖3 橋接原理圖Fig.3 Bridging priciple diagram

2.2 方差隨機(jī)、單多層層厚混合采樣的粒子濾波

粒子濾波解決地聲反演問題時，狀態(tài)方程式(5)可以表示為式(15)形式：

其中：nt一般選擇為高斯分布，對于地聲反演這種局部最優(yōu)多的優(yōu)化問題，選擇高斯分布可能會因隨機(jī)性的不足而導(dǎo)致反演結(jié)果易陷入局部最優(yōu)中。尤其是對于地聲參數(shù)隨距離劇烈變化的情況，即xt與xt-1相距甚遠(yuǎn)的情況，這時候采用高斯分布噪聲，更可能由于采樣隨機(jī)性的不足而導(dǎo)致反演陷入局部最優(yōu)。容易聯(lián)想到增大方差來增加粒子隨機(jī)性的方式，但如此一來由于搜索精度不足，也容易使反演問題陷入局部最優(yōu)。有文獻(xiàn)使用在粒子濾波后增加馬爾可夫鏈蒙特卡洛采樣[9]的方法，利用馬爾科夫鏈蒙特卡洛方法可以在目標(biāo)分布附近采樣，從而增加粒子多樣性。但是馬爾可夫鏈蒙特卡洛方法計算量巨大，代價高昂。

為解決此問題，在每一步中令nt的方差隨機(jī)，使用隨機(jī)方差采樣。令z=0+nt，用方差固定和方差高斯這兩種高斯采樣方式，分別進(jìn)行100 000個點(diǎn)的采樣。方差固定采樣是在方差固定為1的零均值高斯分布中采樣，隨機(jī)方差采樣相對于固定方差采樣，其方差由固定的1變?yōu)殡S機(jī)。隨機(jī)的方差可以在方差為1的零均值高斯分布中采樣獲得。如圖4所示為兩種采樣時，變量z的樣本分布情況。

圖4 兩種采樣方式(方差固定和方差高斯)時，變量z的樣本分布Fig.4 Sample distributions of two sampling methods(fixed variance and Gaussian variance)

從圖4中可以看到，在方差為高斯分布的情況下，樣本在中心更集中，同時，其樣本又有較好的兩端擴(kuò)展能力，即在距離中心較遠(yuǎn)的地方仍然具有相當(dāng)數(shù)量的樣本點(diǎn)。在距離相關(guān)的地聲參數(shù)反演問題中，每幀之間的地聲參數(shù)既可以緩慢變化，也可以劇烈變化，所以方差隨機(jī)的高斯采樣方式十分適合距離相關(guān)的地聲參數(shù)反演問題。本文中的方差隨機(jī)采樣方法，同樣達(dá)到了增加粒子多樣性的效果，進(jìn)一步防止了粒子枯竭，有助于目標(biāo)分布收斂。相較于固定方差的采樣，3.3節(jié)中的模擬結(jié)果驗(yàn)證了方差隨機(jī)采樣可以使反演結(jié)果具有更好的穩(wěn)定性。

對于采樣的維度，傳統(tǒng)粒子濾波中一般使用多維高斯分布采樣，采樣維度與待反演地聲參數(shù)的維度相同。但是根據(jù)敏感度曲線，沉積層厚度有局部最優(yōu)，使用高維高斯分布采樣雖然保證了全局性能，但由于搜索精度不夠，不容易找到最優(yōu)解。如果在搜索中混合低維采樣，相當(dāng)于提高了搜索精度，但同時全局性能不足，容易陷入局部最優(yōu)。本文提出了單層層厚搜索與多層層厚搜索混合的方式，層厚混合搜索的粒子濾波采樣時使用了沉積層厚度的單層多層混合改變的搜索方式，單層搜索指所有沉積層的層厚只改變?nèi)我庖粋€，多層搜索指所有沉積層的層厚全部改變，采用50%概率所有層層厚都變化，50%概率只改變?nèi)我粚訉雍竦姆绞健４朔椒ㄔ谠黾恿巳炙阉髂芰Φ耐瑫r，增加了搜索結(jié)果的精度。3.2節(jié)的仿真結(jié)果驗(yàn)證了此方法的有效性。

基于地聲參數(shù)反演的特點(diǎn)，本文使用了帶橋接、方差隨機(jī)、單層層厚與多層層厚混合搜索的粒子濾波方法，在提高搜索精度的同時保證了解的全局性，更好地處理了地聲參數(shù)隨環(huán)境劇烈變化的問題。

3 仿真結(jié)果分析

3.1 仿真環(huán)境及傳統(tǒng)粒子濾波反演結(jié)果

采用的仿真環(huán)境設(shè)置如下：地聲參數(shù)隨距離變化的100幀海底環(huán)境，在第10幀和第11幀之間沉積層分層數(shù)目改變，同時沉積層聲速和密度也有較大變化。在第20幀和第21幀之間以及第80幀和第81幀之間沉積層分層數(shù)目也發(fā)生改變。第21幀沉積層由4層變?yōu)?層，第81幀沉積層變?yōu)?層。第1幀共4層沉積層，地聲參數(shù)參見表1中。在所有幀中，沉積層總厚度設(shè)置為固定的10 m。粒子數(shù)量設(shè)置為10 000個。水聽器數(shù)量設(shè)置為13個，收集在36°～79°之間共13個掠射角的反射系數(shù)。發(fā)射聲源發(fā)射兩個單頻信號，頻率分別為1 000 Hz和500 Hz。實(shí)驗(yàn)裝置如圖5所示。

圖5 實(shí)驗(yàn)布設(shè)示意圖Fig.5 Layout diagram of experiments

模擬的真實(shí)地聲參數(shù)分布如圖6～8所示。利用傳統(tǒng)粒子濾波算法反演結(jié)果如圖9～11所示。圖9中，每幀聲速用10 000個粒子的均值表示。圖10中，每幀密度用10 000個粒子的均值表示。圖11中，每幀衰減系數(shù)用10 000個粒子的均值表示。

圖6 模擬真實(shí)環(huán)境聲速Fig.6 Simulated sound speed in real environment

圖7 模擬真實(shí)環(huán)境密度Fig.7 Simulated density in real environment

圖8 模擬真實(shí)沉積層衰減系數(shù)Fig.8 Simulated attenuation in real environment

圖9 傳統(tǒng)粒子濾波沉積層聲速c反演結(jié)果Fig.9 Acoustic velocity c inversion results of sedimentary layer obtained by inversion with traditional particle filtering

圖10 傳統(tǒng)粒子濾波沉積層密度ρ反演結(jié)果Fig.10 Density ρ of sediment layer obtained by inversion with traditional particle filtering

圖11 傳統(tǒng)粒子濾波沉積層衰減系數(shù)α反演結(jié)果Fig.11 Attenuation coefficient α of sediment layer obtained by inversion with traditional particle filtering

可見傳統(tǒng)粒子濾波在地聲參數(shù)劇烈變化時，無法正確估計地聲參數(shù)，獲得的地聲參數(shù)有較大誤差。

3.2 層數(shù)混合采樣

使用帶橋接的粒子濾波，一定程度上能夠處理地聲參數(shù)隨距離劇烈變化的反演問題。但是僅僅使用橋接并不夠穩(wěn)定，粒子濾波會出現(xiàn)無法收斂或收斂到局部最優(yōu)的情況。基于地聲反演的特點(diǎn)，本節(jié)將2.2節(jié)提出的混合層厚搜索的方式，與只有橋接但沒有混合層厚搜索的粒子濾波做比較。

由于仿真環(huán)境設(shè)置為第11幀變化最劇烈，第11幀將是粒子濾波最具挑戰(zhàn)性的1幀。接下來對比這三種采樣方式在第11幀時的性能差異。隨機(jī)生成30個不同的仿真環(huán)境分別用單層層厚改變、多層層厚均改變、層厚混合改變的采樣方式進(jìn)行帶橋接的粒子濾波，記錄通過這三種粒子濾波后，第11幀的反射系數(shù)最小均方誤差(反射系數(shù)最小均方誤差即指通過多層橋接后得到的最小均方誤差)。為了使均方誤差不至于太小，這里的均方誤差沒有除以樣本數(shù)量(26)，而是直接用|R-R?|2表示，其中R為觀測反射系數(shù)，R?為用所有粒子均值求得的反射系數(shù)。各掠射角的反射系數(shù)的標(biāo)準(zhǔn)差設(shè)置為0.01。將這三種方式所得的30次實(shí)驗(yàn)的最小均方誤差按照從大到小排列，得到如圖12所示曲線。

圖12 三種采樣方式得出的第11幀的最小均方誤差Fig.12 The minimum mean square errors of the 11th frame obtained by the three sampling methods

從圖12可以看出，層厚的混合改變采樣方式相較于層厚均改變的采樣方式有更小的均方誤差，前者更容易獲得全局最小的反演結(jié)果；相較于層厚單層改變的采樣方式反演效果也有小幅度提升。粒子濾波采用層厚混合采樣的方式時，第11幀地聲參數(shù)有更大的概率收斂到全局最優(yōu)。

3.3 方差隨機(jī)搜索

狀態(tài)方程中的噪聲一般設(shè)為方差較大的高斯擾動，設(shè)置此方差較大的原因是為了增加粒子濾波的搜索范圍，使地聲參數(shù)至少有從某一幀采樣到下一幀真值的可能性。為了對比隨機(jī)方差的粒子濾波與固定方差的粒子濾波對反演效果的影響，進(jìn)行了17次實(shí)驗(yàn)，分別使用隨機(jī)方差(均值為0，方差為1的高斯分布)采樣、固定方差0.6(狀態(tài)方程噪聲方差的0.6倍，下同)采樣、固定方差0.8采樣、固定方差1采樣、固定方差1.2采樣、固定方差1.4采樣。沉積層厚度采用層厚均改變的采樣方式。得到的最小均方誤差結(jié)果如圖13所示。

由圖13可以看出，隨機(jī)方差的采樣方式大大提高了第11幀地聲參數(shù)收斂到全局最優(yōu)的可能性。注意到方差過大(1.4)或過小(0.6)將會出現(xiàn)第11幀最小均方誤差大于1的情況，說明此時粒子濾波效果已十分不佳。同時還可以看出，所有固定方差的采樣方式，最小均方誤差曲線都不足17個點(diǎn)，這是由于某些實(shí)驗(yàn)中，固定方差的采樣方式已經(jīng)不能使粒子濾波穩(wěn)定運(yùn)行，即沒有搜索到一個最優(yōu)或局部最優(yōu)點(diǎn)，這將導(dǎo)致粒子濾波失效。而隨機(jī)方差采樣的方式則可大幅提高粒子濾波的穩(wěn)定性，幾乎每次粒子濾波都能正常收斂。

圖13 不同方差采樣得出的第11幀的最小均方誤差Fig.13 The minimum mean square errors of the 11th frame obtained by different variance sampling methods

為了減少運(yùn)算量，可以令橋接數(shù)量可變，橋接數(shù)量最大為9層，令反射系數(shù)均方誤差基本收斂時停止橋接。實(shí)驗(yàn)采用3.2中的沉積層層厚混合改變的采樣方式，統(tǒng)計所用橋接數(shù)量。圖14是方差固定和方差隨機(jī)時橋接數(shù)量隨幀數(shù)的變化，圖15是將這些橋接數(shù)量按照從大到小排列的結(jié)果。

圖14 方差固定和方差隨機(jī)時橋接數(shù)量隨幀數(shù)的變化Fig.14 Variation of the number of bridges with the number of frames when variance is fixed or random

圖15 當(dāng)幀數(shù)增加到100時橋的數(shù)量按降序排列圖Fig.15 Diagram of the number of bridges arranged in descending order when the number of frames is increased to 100

從圖15可以看出，方差隨機(jī)時平均橋接數(shù)量小于方差為固定值時的平均橋接數(shù)量。隨機(jī)方差采樣平均橋接數(shù)量為4.65層，固定方差采樣平均橋接數(shù)量為5.43層，方差隨機(jī)時算法收斂速度快。可以看到在地聲參數(shù)隨距離變化劇烈時需要的橋接數(shù)量較多，變化緩慢時需要的橋接數(shù)量較少，這是可以預(yù)見的。

3.4 改進(jìn)粒子濾波反演結(jié)果

在粒子濾波中使用橋接的同時，同時使用狀態(tài)方程方差隨機(jī)以及沉積層層厚混合采樣的方式，得到如圖16～18的反演結(jié)果。圖16～18中，每一幀的聲速、粒子、密度均用粒子均值代替。

圖16 改進(jìn)粒子濾波聲速c反演結(jié)果Fig.16 Acoustic velocity c of sediment layer obtained by inversion with improved particle filtering

圖17 改進(jìn)粒子濾波密度ρ反演結(jié)果Fig.17 Density ρ of sediment layer obtained by inversion with improved particle filtering

圖18 改進(jìn)粒子濾波聲衰減系數(shù)α反演結(jié)果Fig.18 Attenuation coefficient α of sediment layer obtained by inversion with improved particle filtering

由反演結(jié)果可見，改進(jìn)粒子濾波在沉積層參數(shù)隨距離劇烈變化的情況下，仍可較準(zhǔn)確地反演地聲參數(shù)。真實(shí)環(huán)境中，沉積層層數(shù)會發(fā)生改變。但在改進(jìn)粒子濾波反演結(jié)果中，沉積層層數(shù)一直使用的都是4層。可以看到在20～80幀，以及80～100幀中，會出現(xiàn)某些層極薄且和上層或下層聲速、密度相近的情況，此時反演結(jié)果中這些層將不含明顯的分層顯示，所以一直使用4層沉積層仍能取得較好的反演結(jié)果。改進(jìn)粒子濾波在第11幀、第21幀、第81幀聲速和密度都有較好的反演結(jié)果。改進(jìn)粒子濾波在第11幀時反射系數(shù)擬合效果如圖19所示，反射系數(shù)曲線擬合效果較好。

圖19 第11幀反射系數(shù)R擬合效果Fig.19 Fitting results of the 11th frame reflection coefficient R

3.5 反演結(jié)果的邊緣概率密度

為了獲得某個地聲參數(shù)的分布情況，分析了反演結(jié)果的邊緣概率密度。圖20～22為第11幀聲速、密度和衰減系數(shù)反演結(jié)果的邊緣概率密度圖。圖20～22中粉色虛線為真值。

從圖20～22中可見聲速密度和衰減系數(shù)均在真值附近，其中聲速的不確定性較小，密度和衰減系數(shù)的不確定性大，這是由于似然函數(shù)對聲速的敏感性高，對密度和衰減系數(shù)的敏感性低。不確定性分析對應(yīng)于第1節(jié)敏感性的分析，證實(shí)了反演結(jié)果的可靠性。

圖20 聲速c邊緣概率密度圖Fig.20 Marginal probability density diagram of sound velocityc

圖21 密度ρ邊緣概率密度圖Fig.21 Marginal probability density diagram of density ρ

圖22 衰減系數(shù)α邊緣概率密度圖Fig.22 Marginal probability density diagram of attenuation coefficient α

4 結(jié)論

針對海底地聲參數(shù)隨距離變化的地聲參數(shù)反演問題，本文提出了一種改進(jìn)的粒子濾波序貫估計方法。在帶橋接的粒子濾波基礎(chǔ)上，使用了沉積層厚度混合變化采樣的方式，增加全局搜索能力的同時使搜索結(jié)果更加精確。在狀態(tài)方程中使用噪聲方差隨機(jī)采樣，增加了粒子的多樣性，同時減少了所需的橋接系數(shù)數(shù)量。相較于傳統(tǒng)粒子濾波，該方法在地聲參數(shù)隨距離劇烈變化的情況下，依舊能較穩(wěn)定地反演地聲參數(shù)。相較于只有橋接的粒子濾波，本文方法的反演結(jié)果更加穩(wěn)定，計算效率更高。

利用不同掠射角和不同頻率的反射系數(shù)作為待擬合的觀測數(shù)據(jù)，采用改進(jìn)粒子濾波的方法進(jìn)行地聲參數(shù)反演。仿真結(jié)果證明了在沉積層參數(shù)劇烈變化的環(huán)境中，本文方法仍然具有良好的效果。