采用子帶二值加權(quán)累積的經(jīng)驗(yàn)格林函數(shù)提取?

郭啟超 李風(fēng)華 楊習(xí)山

(1 中國(guó)科學(xué)院聲學(xué)研究所 聲場(chǎng)聲信息國(guó)家重點(diǎn)實(shí)驗(yàn)室 北京 100190)

(2 中國(guó)科學(xué)院大學(xué) 北京 100049)

0 引言

理論與實(shí)驗(yàn)表明，兩點(diǎn)之間的經(jīng)驗(yàn)格林函數(shù)(Empirical Green’s function, EGF)可以通過(guò)此兩點(diǎn)之間的環(huán)境噪聲互相關(guān)函數(shù)(Noise crosscorrelation function, NCF)累積后得到[1?3]。相比于主動(dòng)發(fā)射信號(hào)估計(jì)格林函數(shù)，NCF累積方法具有成本低、隱蔽性強(qiáng)與環(huán)境友好的特點(diǎn)，已被成功應(yīng)用于被動(dòng)聲層析[4?5]、陣列陣型校準(zhǔn)[6]與目標(biāo)無(wú)源成像[7?8]等方面。然而，由于實(shí)際海洋噪聲場(chǎng)的非均一性，往往需要累積較長(zhǎng)時(shí)間的NCF才能獲得具有一定信噪比的EGF，限制了NCF 理論在環(huán)境參數(shù)變化較快的場(chǎng)景中的應(yīng)用。

目前，改善短時(shí)間內(nèi)EGF提取質(zhì)量的方法主要分為空域?yàn)V波以及時(shí)頻域?yàn)V波兩類(lèi)。空域?yàn)V波方法主要是通過(guò)陣列波束形成[9?11]、矢量水聽(tīng)器[12?13]或矩陣分解[14?15]濾除兩點(diǎn)非端射方向的非擴(kuò)散噪聲干擾，可以顯著加快NCF 的收斂速度。時(shí)頻域?yàn)V波主要包括頻域白化等預(yù)處理手段[16]與一系列加權(quán)累積算法[17?20]等，在不同場(chǎng)景下獲得了一定的優(yōu)化效果。在加權(quán)累積算法中，時(shí)域二值加權(quán)累積算法[17](Time-domain binary weighted stacking,TBWS)通過(guò)選取可以使EGF 信噪比增大的NCF快拍在時(shí)域進(jìn)行累積，并結(jié)合全局最優(yōu)化算法，獲得了比均方根加權(quán)累積類(lèi)算法[18?19]更好的地震波EGF 提取效果。然而，對(duì)于帶寬較寬的環(huán)境噪聲提取EGF 場(chǎng)景，直接將NCF 在時(shí)域進(jìn)行選擇性累積，會(huì)丟失(或保留)每一快拍NCF 中有用(或無(wú)用)的頻率成分，降低EGF的提取質(zhì)量。

本文首先通過(guò)實(shí)測(cè)海洋環(huán)境噪聲數(shù)據(jù)驗(yàn)證了TBWS 算法相比于原始累積方法能夠有效提高EGF 信噪比，同時(shí)指出TBWS 算法在時(shí)域累積的不足，并提出一種基于子帶二值加權(quán)累積(Subband binary weighted stacking, SBWS)的改進(jìn)算法。新算法基于各頻域子帶可以產(chǎn)生具有一定時(shí)間分辨率的EGF 的假設(shè)，在每一個(gè)子帶內(nèi)根據(jù)“累積后提取的EGF 信噪比增加”準(zhǔn)則對(duì)每一快拍NCF進(jìn)行選擇性累積，將各子帶累積結(jié)果譜白化后在頻域拼接，再反傅里葉變換得到最終的EGF 提取結(jié)果。相比于TBWS 算法，改進(jìn)算法將二值加權(quán)思想從時(shí)域拓展到各頻域子帶，兼顧了EGF對(duì)噪聲快拍與頻率的選擇，可以有效保留(或?yàn)V除)每一快拍數(shù)據(jù)中的有用(或無(wú)用)的子帶頻率成分，較為明顯地提高了從海洋環(huán)境噪聲中提取EGF 的信噪比。將SBWS 算法與空域?yàn)V波方法結(jié)合，成功實(shí)現(xiàn)了利用海洋環(huán)境噪聲進(jìn)行2 h一次的海水平均聲速反演。

1 理論與方法

1.1 EGF提取基本理論

本小節(jié)介紹經(jīng)驗(yàn)格林函數(shù)(EGF)提取的基本理論。兩點(diǎn)之間的歸一化噪聲互相關(guān)函數(shù)(NCF)在頻域可以表示為

式(1)中，yl(ω)表示第l快拍的頻域NCF，L為快拍總數(shù)，p1(ω)與p2(ω)分別表示兩個(gè)接收器記錄的噪聲信號(hào)的頻域形式，操作符(·)H與|·|分別表示對(duì)信號(hào)取共軛復(fù)數(shù)與絕對(duì)值，歸一化(即頻域白化)的作用是消除可能存在的強(qiáng)能量頻率干擾。利用空域?yàn)V波方法，將兩條陣列指向端射方向(如圖1(a)中的角度θ0)的波束輸出互相關(guān)代替原本的兩個(gè)單接收器聲壓互相關(guān)，可以顯著地提高每一快拍NCF 的質(zhì)量。

原始累積方法直接將多個(gè)快拍的yl(ω)相加，然后取反傅里葉變換，得到兩接收器位置之間的EGF：

式(2)中，x(t)與x(?t)分別表示因果EGF 與非因果EGF，它們關(guān)于原點(diǎn)對(duì)稱，jω為相移因子[21]。

1.2 子帶二值加權(quán)累積

本小節(jié)介紹基于子帶二值加權(quán)累積提取EGF的方法。首先將每一快拍頻域NCFyl(ω)按照一定的頻率間隔均勻劃分為K個(gè)子帶，在每個(gè)子帶中，賦予每一快拍NCF 一個(gè)二值加權(quán)系數(shù)，然后將多個(gè)快拍的NCF 加權(quán)累積，得到各子帶的累積結(jié)果Yk(ω)：

式(3)中，yk,l(ω)表示第k個(gè)子帶中的第l快拍頻域NCF，wk,l為對(duì)應(yīng)的二值加權(quán)系數(shù)，其取值為0 或1，L為快拍總數(shù)。wk,l的大小通過(guò)各子帶“累積后提取的EGF信噪比增加”準(zhǔn)則確定，即

式(4)中，SNRk,l與SNRk,l+1分別表示將yk,l(ω)累積前后反傅里葉變換得到的第k子帶EGF 的信噪比。本文沿用文獻(xiàn)[9]中對(duì)EGF信噪比的定義，即預(yù)期相干峰值到達(dá)時(shí)間段內(nèi)的包絡(luò)最大值與剩余時(shí)間殘余波動(dòng)的水平之比，其中剩余時(shí)間波動(dòng)的水平是從EGF在足夠大的非預(yù)期到達(dá)時(shí)間段內(nèi)(本文任意選擇為3～10 s)的標(biāo)準(zhǔn)差估計(jì)出來(lái)的。

式(5)中，2.3～2.4 s 為EGF 預(yù)期相干峰值出現(xiàn)的區(qū)間，3～10 s 為殘余擾動(dòng)存在的一段區(qū)間，操作符Hilbert(·)與std(·)分別表示對(duì)信號(hào)求希爾伯特變換以及標(biāo)準(zhǔn)差。根據(jù)干涉原理，如果兩個(gè)NCF 之間是相干的，那么兩者累積后信噪比必然增加[17]。因此EGF 信噪比增加準(zhǔn)則的目的是保證加入累積過(guò)程的NCF之間是相干的。

將式(3)中各子帶的累積結(jié)果譜白化后在頻域拼接，

式(6)中，譜白化的作用是防止某個(gè)頻段信號(hào)能量過(guò)強(qiáng)引起時(shí)域EGF 波形畸變。對(duì)Y(ω)進(jìn)行反傅里葉變換即可得到EGF的最終提取結(jié)果。

2 海試實(shí)驗(yàn)

2.1 實(shí)驗(yàn)數(shù)據(jù)介紹

實(shí)驗(yàn)于2018年5月在南中國(guó)海北部海域開(kāi)展，實(shí)驗(yàn)布放如圖1(a)所示。海深約100 m，兩套海底接收系統(tǒng)相距約3.5 km，每條陣列有15 個(gè)陣元，同步記錄海洋環(huán)境噪聲，采樣頻率5000 Hz。在實(shí)驗(yàn)過(guò)程中利用溫深儀(Temperature depth sensor, TD)每30 s 測(cè)量對(duì)應(yīng)深度的溫度，并利用聲速經(jīng)驗(yàn)公式獲得水體聲速，實(shí)測(cè)聲速剖面如圖1(b)所示。選取頻段為20～400 Hz[4]的噪聲數(shù)據(jù)進(jìn)行處理，首先將兩個(gè)陣列接收的噪聲數(shù)據(jù)的全部觀察時(shí)間分割成彼此互相連接的若干個(gè)時(shí)間段，每段時(shí)間長(zhǎng)度為10 s(如對(duì)于兩個(gè)小時(shí)噪聲數(shù)據(jù)，則每條陣列上每個(gè)陣元可以產(chǎn)生720 段的噪聲信號(hào))。接著根據(jù)第1.1節(jié)空域?yàn)V波方法，將兩個(gè)陣列的波束輸出做互相關(guān)，共得到720 段互相關(guān)輸出。將各段互相關(guān)輸出編號(hào)為l= 1,···,L，則每一段被稱為一個(gè)快拍的頻域NCF。算法處理中，每一快拍NCF 在頻域等間隔劃分為10個(gè)子帶，子帶帶寬為38 Hz。

圖1 實(shí)驗(yàn)布放示意圖及海水聲速剖面Fig.1 Experimental layout and sea sound speed profile

2.2 實(shí)驗(yàn)數(shù)據(jù)處理結(jié)果

時(shí)間結(jié)構(gòu)準(zhǔn)確性是應(yīng)用EGF 進(jìn)行聲層析或陣型校準(zhǔn)的基礎(chǔ)。本小節(jié)驗(yàn)證基于本文提出的子帶二值加權(quán)累積算法(SBWS)提取EGF 時(shí)間結(jié)構(gòu)的準(zhǔn)確性。圖2(a)為利用SBWS算法在兩小時(shí)海洋環(huán)境噪聲中提取的典型EGF，其中紅色散點(diǎn)為使用Bellhop 通過(guò)實(shí)測(cè)海水聲速剖面計(jì)算的本征聲線到達(dá)時(shí)間。EGF在2.3～2.4 s 時(shí)間段內(nèi)存在4 個(gè)明顯的波包，并且包絡(luò)的谷值或峰值位置與理論到達(dá)時(shí)間吻合較好。圖2(b)為72 h 內(nèi)利用EGF 第二個(gè)峰值到達(dá)時(shí)間反演的海水深度上的平均聲速與TD 實(shí)測(cè)值的對(duì)比，兩者的均方根誤差為0.32 m/s。海試數(shù)據(jù)驗(yàn)證了基于SBWS 算法能夠準(zhǔn)確地提取EGF的時(shí)間結(jié)構(gòu)，聲速反演結(jié)果能正確反映水體溫度的日變化周期。

圖2 提取的EGF 及反演的海水平均聲速Fig.2 The extracted EGF and the inverted average sea sound speed

頻率劃分間隔(即子帶帶寬)是SBWS 算法的重要參數(shù)，下面討論子帶帶寬對(duì)所提算法的影響。圖3為72 h 平均后EGF 的信噪比與聲速反演誤差隨子帶帶寬的變化情況?？梢钥吹剑珽GF 信噪比隨著子帶帶寬減小而逐漸增大，原因是子帶帶寬越窄，則算法對(duì)噪聲的頻率挑選越精細(xì)，從而盡可能地保留(或?yàn)V除)對(duì)EGF 信噪比有用(或無(wú)用)的頻率成分。當(dāng)子帶帶寬取最大帶寬380 Hz 時(shí)，算法退化為時(shí)域二值加權(quán)累積算法(TWBS)。另一方面，由于時(shí)間分辨率與頻率帶寬成反比，帶寬變窄會(huì)使各子帶EGF 包絡(luò)變寬，可能導(dǎo)致EGF 中錯(cuò)誤位置的峰值能量增加，影響到達(dá)時(shí)間估計(jì)的準(zhǔn)確性。如圖3(b)所示，聲速反演誤差隨著子帶帶寬變窄先平穩(wěn)變化，當(dāng)小于38 Hz 時(shí)明顯增大。為了盡可能提高EGF 信噪比，同時(shí)保證時(shí)間結(jié)構(gòu)的準(zhǔn)確性，將SBWS算法的子帶帶寬設(shè)置為38 Hz。

圖3 EGF 信噪比及聲速反演誤差隨子帶帶寬的變化情況Fig.3 The signal-to-noise ratio of EGF and the inversion error of sound speed with subband bandwidth

以EGF 信噪比為評(píng)價(jià)指標(biāo)，驗(yàn)證SBWS 算法相比于已有算法的優(yōu)越性。圖4展示了兩種典型情況下從2 h 海洋環(huán)境噪聲中提取的EGF 及其對(duì)應(yīng)的信噪比隨累積時(shí)間的變化曲線。由圖4(a)可知，相比于原始累積方法，TBWS 與SBWS 算法可以較為明顯地提高EGF 的信噪比。原因如圖4(b)的EGF 信噪比曲線所示，加權(quán)累積方法通過(guò)二值選擇過(guò)程，避免了干擾成分進(jìn)入到EGF 的累積過(guò)程中。另一方面，由圖4(c)可以看到，TBWS 算法存在提取質(zhì)量反而不如原始累積方法的情況。原因如圖4(d)所示，在2～40 min 內(nèi)，原始累積方法得到的EGF信噪比存在先下降后上升的情況，TBWS算法在信噪比下降初期就不再累積，雖然保證了局部時(shí)間段內(nèi)的最優(yōu)化，但也丟失掉了40 min 之后可以使EGF 信噪比提升的快拍。而本文提出的SBWS 算法能夠在TBWS 算法舍棄的NCF 中尋找有利的子帶頻率成分，相比于原始累積與TBWS算法可以一定程度上提升EGF的信噪比。

圖4 兩種典型情況下提取的EGF 及其信噪比變化曲線Fig.4 The EGF and its signal-to-noise ratio curve of two typical cases

將3 種方法的提取結(jié)果經(jīng)過(guò)72 h 平均后，得到相對(duì)EGF信噪比曲線(初始值均置零)如圖5所示。可以看到3種方法的信噪比變化趨勢(shì)符合20lg，其中T為累積時(shí)間，這與NCF 提取EGF 的理論相符[12]。并且相比于原始累積方法，TBWS 與SBWS 算法提取的EGF 信噪比均有一定程度的提高，如120 min 時(shí)EGF信噪比分別提高了1.5 dB 與2.8 dB。海試實(shí)驗(yàn)數(shù)據(jù)驗(yàn)證了SBWS算法相比于原始累積法與TWBS 算法，能夠有效提高EGF 的信噪比。

圖5 三種方法提取的EGF 信噪比曲線Fig.5 The signal-to-noise ratio curve of EGF extracted by three methods

本小節(jié)繼續(xù)討論SBWS 算法能夠提升EGF 信噪比的原因。圖6為72 h 平均后從2 h 噪聲中提取的各子帶EGF 信噪比以及對(duì)應(yīng)的子帶信噪比累積曲線。從圖6(a)可知，較高頻率(>172 Hz)的子帶EGF 信噪比較低，并且從圖6(b)中各子帶信噪比累積曲線可知，在同一累積時(shí)間下，子帶頻率越高，則提取的EGF 信噪比越低。圖6說(shuō)明了影響EGF提取速度的主要干擾來(lái)自于較高頻率的噪聲成分。圖7為各子帶中被濾除的快拍數(shù)占快拍總數(shù)的比例。由圖7可知，各子帶中均濾除了50%左右的快拍并且不同子帶的濾除率不同，說(shuō)明有利于EGF提取的最佳頻段是隨時(shí)間變化的。SBWS算法通過(guò)對(duì)各子帶NCF 進(jìn)行選擇，實(shí)現(xiàn)了在累積之前對(duì)每一快拍NCF 中的干擾頻率(主要存在較高頻段)自動(dòng)濾除，即在每一時(shí)刻逼近最佳頻段的最優(yōu)解，因此能夠有效提高全頻帶EGF 的信噪比。值得注意的是，盡管只利用較低頻率的子帶就可以提取具有較高信噪比的EGF，但帶寬變窄會(huì)導(dǎo)致EGF 時(shí)間分辨率下降，而本文提出的SBWS算法不會(huì)產(chǎn)生此問(wèn)題。

圖6 各子帶EGF 信噪比及對(duì)應(yīng)的信噪比累積曲線Fig.6 The signal-to-noise ratio and the corresponding signal-to-noise ratio curve of EGF in each sub-band

圖7 各子帶快拍的濾除率Fig.7 The filtering rate of snapshot in each subband

3 結(jié)論與討論

海洋環(huán)境噪聲中有利于EGF 提取的最佳頻段是隨時(shí)間變化的。本文提出一種采用子帶二值加權(quán)累積的海洋環(huán)境噪聲互相關(guān)提取EGF 算法。該算法通過(guò)濾除每個(gè)子帶中令EGF 信噪比下降的快拍數(shù)據(jù)，實(shí)現(xiàn)了在累積之前對(duì)每一快拍噪聲的頻率自動(dòng)挑選，相比于原始累積方法與時(shí)域二值加權(quán)累積算法，可以有效提高從帶寬較寬的海洋環(huán)境噪聲中提取EGF的信噪比，算法性能通過(guò)了海試試驗(yàn)數(shù)據(jù)的驗(yàn)證。與時(shí)域二值加權(quán)累積算法相同，子帶二值加權(quán)累積算法是以“EGF 信噪比增加”作為加權(quán)系數(shù)的判定標(biāo)準(zhǔn)，該算法的不足是當(dāng)EGF峰值位置與真實(shí)的到達(dá)時(shí)間存在誤差時(shí)，“EGF信噪比增加”準(zhǔn)則會(huì)造成錯(cuò)誤的峰值能量增加。在獲取足夠的信道參數(shù)等先驗(yàn)知識(shí)條件下，通過(guò)縮小預(yù)期峰值的到達(dá)區(qū)間，可以一定程度上改善這個(gè)問(wèn)題。