李赫 郭新毅 馬力

1) (中國科學(xué)院聲學(xué)研究所,中國科學(xué)院水聲環(huán)境特性重點實驗室,北京 100190)

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

海洋環(huán)境噪聲場中包含了海洋中的諸多信息,海底地聲參數(shù)是影響海洋環(huán)境噪聲場空間分布的主要因素之一.對于不同的海底分層結(jié)構(gòu),海底反射損失會根據(jù)沉積層厚度和各層聲速呈現(xiàn)出不同的臨界角和干涉條紋結(jié)構(gòu).本文利用Harrison能流理論,從理想反射系數(shù)出發(fā),分別考慮了聲速、密度、衰減系數(shù)、沉積層厚度等幾種參數(shù)對無沉積層和單層沉積層中反射系數(shù)的影響,并對單層沉積層海底的反射系數(shù)進(jìn)行了化簡,結(jié)合互易原理解釋了反射損失條紋結(jié)構(gòu)的形成機理.中國黃海某海區(qū)試驗結(jié)果表明,利用海洋環(huán)境噪聲空間方向譜獲得的海底反射損失,可以提取海底反射臨界角和干涉條紋信息,由此可估計出海底分層結(jié)構(gòu)、聲速和沉積層厚度等海底參數(shù)信息.

1 引 言

海洋環(huán)境噪聲是海洋中永恒存在的聲場,主要由風(fēng)浪、降雨、艦船、海洋生物、工業(yè)等因素形成,其中風(fēng)成噪聲在各個海域普遍存在,其頻段覆蓋幾百至幾十千赫茲[1].一方面海洋環(huán)境噪聲是影響聲吶工作性能的主要因素之一,只有對其進(jìn)行充分的了解研究才能充分提高水下設(shè)備的工作性能,所以國內(nèi)外學(xué)者紛紛對海洋環(huán)境噪聲進(jìn)行建模研究[2?5].另一方面,海洋環(huán)境噪聲場中也包含了諸多水體、海底、海面等環(huán)境信息[6],因此利用海洋環(huán)境噪聲進(jìn)行參數(shù)反演,獲取海水、海底等環(huán)境信息也成為時下的熱點問題之一[7,8].

海底地聲參數(shù)反演主要分為主動反演和被動反演兩部分.Zeng等[9]利用爆炸聲源對中國黃海海域海底衰減系數(shù)進(jìn)行反演;李夢竹等[10]提出了一種適用于低聲速沉積層的海底參數(shù)聲學(xué)反演方法.相比于主動反演,利用海洋環(huán)境噪聲進(jìn)行被動反演無需主動發(fā)射聲源信號,在實驗中僅需數(shù)據(jù)接收裝置,大幅度節(jié)省了實驗消耗以及工作任務(wù),并且對海洋生物沒有任何影響.周建波等[11]利用矢量水聽器通過海洋環(huán)境噪聲提取了聲場格林函數(shù);江鵬飛等[12]根據(jù)海洋環(huán)境噪聲的空間指向性在不同掠射角范圍內(nèi)對海底參數(shù)的敏感度不同,提出了分步反演的方法;駱文于[13]利用中國東海環(huán)境噪聲數(shù)據(jù)反演了液態(tài)半空間海底聲速、密度、衰減系數(shù);早在1982年,布列霍夫斯基[14]就曾利用海洋環(huán)境噪聲場的各向異性,來獲取海底反射信息;Harrison在2002年進(jìn)一步系統(tǒng)地研究了海洋環(huán)境噪聲中的能流理論,利用上行波和下行波的比值提取了海底反射損失,并對海底進(jìn)行了反演[15],并利用希爾伯特變換補全反射系數(shù)的虛部,通過計算垂直方向的沖擊響應(yīng)利用漂流浮標(biāo)對海底進(jìn)行成像[16];Siderius等[17]利用海洋環(huán)境噪聲提取了格林函數(shù)用來表征兩點間的沖擊響應(yīng),對海底分層進(jìn)行成像;Muzi等[18,19]在Harrison的基礎(chǔ)上利用互譜密度矩陣的特性對陣元進(jìn)行了理論擴展處理,提高了利用海洋環(huán)境噪聲提取的海底損失曲線的分辨率.

本文在Harrison能流理論基礎(chǔ)上,從理想環(huán)境反射系數(shù)出發(fā),解釋了不同海底參數(shù)對海底反射損失(bottom?loss,BL)曲線的影響,并利用實測海洋環(huán)境噪聲數(shù)據(jù)對海底分層結(jié)構(gòu)、聲速及沉積層厚度進(jìn)行估計.本文首先引入Harrison能流理論估計海底反射損失的方法,分析了不同海底分層情況下的曲線特點;然后利用簡正波理論,討論了淺海波導(dǎo)中近場連續(xù)譜和遠(yuǎn)場離散譜在噪聲估計海底反射損失中的貢獻(xiàn);并從理想反射系數(shù)出發(fā),分別討論了不同海底參數(shù)對無沉積層、單一海底沉積層兩種情況的影響,結(jié)合互易原理解釋了反射損失條紋結(jié)構(gòu)的產(chǎn)生機理,并根據(jù)噪聲獲得的海底反射損失曲線給出海底聲速、沉積層厚度的估計方法;最后介紹了 2016年冬中國黃海某海區(qū)海上試驗,對實測海洋環(huán)境噪聲數(shù)據(jù)進(jìn)行處理,估計其海底分層結(jié)構(gòu)、海底聲速和沉積層厚度.

2 利用海洋環(huán)境噪聲垂直方向譜估計海底反射損失

從淺海理想波導(dǎo)出發(fā),假設(shè)噪聲源均勻分布在海面,利用垂直陣接收環(huán)境噪聲信號,其空間結(jié)構(gòu)如圖1所示.波導(dǎo)中的聲波可分為兩部分： 接收陣上方的來波稱為下行波,從噪聲源直達(dá)接收陣;下方的稱為上行波,經(jīng)過海底一次反射后到達(dá)接收陣.在此近場情況下不考慮多次反射,后文將給出解釋.

圖1 淺海中利用垂直陣接收噪聲空間結(jié)構(gòu)Fig.1.Spatial structure of receiving ocean ambient noise using vertical line array in shallow water.

聲速剖面為等聲速的理想情況下,表面聲源的出射角θs、海底掠射角θb與接收角θr均相等,在此均用掠射角θb表示.對于某一角頻率ω下,海底反射損失BL有如下定義：

其中B為波束形成得到的不同角度和角頻率下的接收能量,即垂直方向能量譜,

其中E代表期望;H為共軛轉(zhuǎn)置;w是波束形成的加權(quán)向量;p為對應(yīng)角頻率ω的接收數(shù)據(jù),在M陣元接收時,為M維向量;Cω為噪聲互譜密度矩陣.

為說明此方法計算得到BL的有效性,先考慮極端情況,假設(shè)陣元間隔為0.2 m的垂直陣遍布在聲速為1500 m/s、水深為40 m的全海深中,忽略海底橫波的影響,海底參數(shù)如表1所列.

表1 海水及海底聲學(xué)參數(shù)Table 1.Parameters of ocean and sub?bottom.

互譜密度矩陣直接利用OASES軟件中的OASN模塊計算獲得[20],w設(shè)置為常規(guī)波束形成加權(quán)向量.圖2(a)為根據(jù)理想情況下真實反射系數(shù)公式求得的BL[21],后文將對該公式進(jìn)行討論;圖2(b)是根據(jù)噪聲垂直空間指向性利用OASN模塊仿真得到的;圖2(c)提取出了1800 Hz下兩種方法的掠射角BL曲線.圖2(a)和圖2(b)中的條紋與圖2(c)中的曲線峰值成對應(yīng)關(guān)系,陣元遍布全海深的理想狀態(tài)下,圖2(c)中兩條曲線的峰值一一對應(yīng),僅在幅度上存在一定偏差,因此利用Harrison能流理論估計海底反射損失是可行的.

圖2 (a) 真實反射損失BL;(b) 根據(jù)噪聲垂直空間指向性利用OASN模塊仿真得到的 ;(c) 1800 Hz下兩種方法的比較,實線為BL,虛線為Fig.2.(a) True BL;(b) computed by vertical direc?tionality of ocean ambient noise using OASN;(c) two meth?ods compare under 1800 Hz,the full line is BL,the imagin?ary line is.

而在實際海上試驗中,陣元個數(shù)是有限的,很難達(dá)到上述遍布全海深的理想狀態(tài).有限的陣元個數(shù)會帶來不同掠射角方向分辨率的下降,這樣會導(dǎo)致圖2(b)中條紋信息一定程度上的缺失,同時也使圖2(c)中峰值的模糊,但仍然有很多特征可以利用.圖3給出了兩種不同海底分層結(jié)構(gòu)的,保留表1中的參數(shù),圖中實線為間隔0.2 m的水聽器遍布全海深的結(jié)果,虛線為間隔0.2 m的42元水聽器陣,第一個陣元設(shè)置在水深30 m處,頻率均為1800 Hz.海底模型如圖4所示,其中圖4(a)為無限大液態(tài)聲學(xué)半空間海底,圖4(b)為帶有一層沉積層的海底.

圖3 1800 Hz下不同海底分層結(jié)構(gòu)下的 (實線為陣元遍布全海深,點劃線為42陣元,間隔均為0.2 m) (a) 海底為無限大液體聲學(xué)半空間;(b) 海底為單層沉積層Fig.3.The of different structure of sub?bottom under 1800 Hz： (a) Infinite acoustic half space;(b) 1 layer of sedi?ment.The full line corresponds to the condition that the elements set across the sea.The imaginary line corresponds to the condition that 42 elements set at the depth of 30 m.

圖4 海底分層模型 (a)液體無限大聲學(xué)半空間海底;(b) 存在一層沉積層海底Fig.4.Model of sub?bottom stratification： (a) Infinite acoustic half space;(b) 1 layer of sediment.

3 不同方向噪聲對 曲線影響分析

在淺海波導(dǎo)中,點聲源在大掠射角范圍內(nèi)發(fā)出的聲波很難進(jìn)行遠(yuǎn)距離傳播,因為這一部分經(jīng)過海底反射存在能量損失,反射波攜帶著一定海底信息,這一特性對于均勻分布在表面的噪聲源同樣適用.在分層介質(zhì)中,點源聲場的波函數(shù)的積分表達(dá)式為

為了描述準(zhǔn)確,這一部分在K/I模型下進(jìn)行討論[2].Kuperman給出海面噪聲源均勻情況下噪聲互譜密度的波數(shù)積分表達(dá)式

圖5 Pekeris分支割線、極點和積分的復(fù)波數(shù)平面Fig.5.Pekeris branch cut map.

和簡正波表達(dá)式

式中g(shù)為點源格林函數(shù),Ψm為模式函數(shù).由于接收器為垂直陣,式中r=0.對于近場連續(xù)譜部分采用(5)式波數(shù)積分表達(dá)式進(jìn)行計算,而遠(yuǎn)場離散譜部分則采用(6)式簡正波表達(dá)式,將其分別代入(3)式中求得不同譜域內(nèi)的.圖6給出了無限大液態(tài)聲學(xué)半空間海底情況下,不同譜域內(nèi)噪聲場的垂直指向性和對應(yīng)的,圖6(a)和圖6 (b)是利用(5)式計算近場連續(xù)譜部分,即 0＜kr＜k2,這一部分的貢獻(xiàn)主要來自大掠射角,其對應(yīng)的是大于臨界角的部分,而這部分經(jīng)過海底后存在能量泄漏,返回到接收器陣的聲波可以表達(dá)出海底反射的信息;圖6(c)和圖6(d)是(6)式計算得到的遠(yuǎn)場離散譜部分,此時k2＜kr＜k1,其貢獻(xiàn)主要來自于小掠射角,對應(yīng)小于臨界角的部分,這一角度內(nèi)發(fā)生全內(nèi)反射,海底沒有聲能吸收,對的大掠射角條紋部分沒有貢獻(xiàn);圖6(e)和圖6 (f)為全譜域結(jié)果,與僅考慮近場條件下的基本一致.在近場條件下大掠射角部分反射損失較高,所以本文不考慮多次海底反射的情況.

圖6 不同譜域內(nèi)噪聲場的指向性和 (a),(b)連續(xù)譜部分噪聲場;(c),(d) 離散譜部分噪聲場;(e),(f)全譜域噪聲場Fig.6.Vertical directivity and of ocean ambient noise in different spectral domain： (a),(b) Near field continuous spectrum;(c),(d) far field discrete spectrum;(e),(f) full spectrum domain.

4 BL曲線條紋結(jié)構(gòu)產(chǎn)生機理及海底參數(shù)估計

4.1 不同海底參數(shù)對無沉積層BL的影響

影響B(tài)L的海底地聲參數(shù)有3個,分別是海底衰減系數(shù)、聲速c和密度.圖7別為不同海底參數(shù)下的BL.由圖7(a)可見,僅在掠射角小于臨界角時,對BL有影響,在建立淺海遠(yuǎn)程傳播模型中起到了至關(guān)重要的作用[21].但第3節(jié)中,通過對近場連續(xù)譜和遠(yuǎn)場離散譜的分解已經(jīng)證明： 噪聲估計的主要來自于近場連續(xù)譜的貢獻(xiàn),在近場條件下,衰減系數(shù)又是一個不敏感參數(shù),因此可以忽略衰減系數(shù)帶來的影響.在分界面上掠射角θ1,θ2利用Snell定律建立關(guān)系

圖7 不同海底參數(shù)下的BL (a)衰減系數(shù) ;(b)聲速c;(c)密度Fig.7.The BL under different parameters of sub?bottom：(a) Attenuation coefficient ;(b) sound speed c;(c) dens?ity.

透射波水平射出,即θ2=0°時存在臨界角效應(yīng),此時,cosθ1=c1/c2,在海水聲速c1已知的情況下,反射系數(shù)的臨界角θ1完全由c2決定.除此之外,海底聲速對于反射系數(shù)的幅度也有一定影響.圖7(b)為不同海底聲速下的BL.海底密度ρ僅對BL的幅值起到影響,與臨界角不存在任何關(guān)系,如圖7(c).

4.2 單層沉積層BL條紋結(jié)構(gòu)產(chǎn)生機理

下面考慮存在一層沉積層時的反射情況,反射系數(shù)的表達(dá)式為

其中,R12,R23分別為兩個分界面上的反射系數(shù),φ=k2hsinθ2表示聲波穿過沉積層的垂直相移.(8)式的分子和分母中均含有 exp(2jφ) 一項,難以討論其各個參數(shù)對R的影響,特別是圖2反射系數(shù)中影響條紋的因素,為了讓各個參數(shù)在兩式中的作用明顯表達(dá),需要對其進(jìn)行一定程度上的化簡,下面從單一沉積層的反射模型開始考慮,如圖8.

圖8 單層沉積層海底反射模型Fig.8.Reflection model of sub?bottom with 1 layer of sedi?ment.

聲波在沉積層中傳播,會形成多次返回海水中的反射波,反射系數(shù)的原始表達(dá)為各個反射波的疊加

式中,T12表示海水進(jìn)入沉積層的透射系數(shù),T21代表沉積層進(jìn)入海水中的透射系數(shù).近場條件中,海面風(fēng)成噪聲源源級較低,可忽略海底多次反射后返回海水中反射波帶來的影響,僅考慮圖8中前兩條返回聲線,這種情況下反射系數(shù)化簡為

與(8)式相比,(10)式非常簡潔,便于討論不同參數(shù)與R的對應(yīng)關(guān)系.圖9分別為(8)和(10)兩式計算的海底反射損失BL和1800 Hz下的BL曲線,依然沿用表1中參數(shù).從圖9可以看出,兩式計算的結(jié)果基本一致,僅在幅值上存在細(xì)微差異.從級數(shù)的角度解釋,大掠射角近場條件下,R為一純小數(shù),(9)式中第三項之后均為第二項的高階小量,因此僅在幅值上存在一定影響.綜上,可以利用(10)式化簡后的反射系數(shù)R代替(8)式,對單層沉積層海底反射進(jìn)行討論.

用θc1和θc2表示兩個臨界角,當(dāng)掠射角θ＜θc1時,沒有透射波,(10)式僅存在R12一項,此時R為模近似1的復(fù)數(shù);當(dāng)θc1＜θ＜θc2時,(10)式中R12,T12,T21均為實數(shù),R23為模近似1的復(fù)數(shù),R12為隨θ單調(diào)遞增的負(fù)數(shù),其模單調(diào)遞減,T12=1+R12為正且單調(diào)遞增,T21同理,因此對于此部分振蕩完全由 exp(2jφ) 一項引起;當(dāng)θ＞θc2時,(10)式中除 exp(2jφ) 一項以外均為隨θ單調(diào)改變的實數(shù),其本身不存在極值點,因此這一部分與θc1＜θ＜θc2時一樣,BL曲線中的振蕩均由exp(2jφ)一項引起.在計算BL過程中,取R的模,所以exp(2jφ)項中僅有實部對BL曲線有貢獻(xiàn),其實部通過Euler公式可知

因此BL曲線中的振蕩由聲速和沉積層厚度h共同決定.在海底聲速確定的情況下,對于某一掠射角θ0,BL曲線對于頻率f的振蕩周期可由(11)式確定,把頻率f看作變量,(11)式表示為

其振蕩周期為c2/2hsinθ2,在某一頻率fc范圍內(nèi)振蕩次數(shù)為 2hsinθ2/c2·fc,振蕩次數(shù)即BL的極值點數(shù),這與圖9(a)和圖9(b)中BL的條紋個數(shù)是對應(yīng)的.

根據(jù)互易原理,可以看作不同距離的表面噪聲源到達(dá)同一接收點的路徑,圖8中前兩條反射聲線看作入射聲線,原來的入射聲線看作匯聚后的出射聲線.兩條入射聲線的相位差等于聲線在沉積層中走過路線的相位變化,即垂直相位變化 2φ.在相位相差 2π 的整數(shù)倍時發(fā)生相長干涉,

圖9 (a) (8)式計算的真實BL;(b) (10)式計算化簡后的BL;(c) 1800 Hz下兩條BL曲線Fig.9.(a) True BL computed by equation (8);(b) BL after simplified computed by equation (10);(c) the curves of BL of two methods under 1800 Hz.

(13)式f0即為(12)式中f的峰值點,與其振蕩周期相對應(yīng),同樣也解釋了圖9中BL條紋的形成機理.

而在某一頻率確定的情況下,(12)式中sinθ2可用(7)式確定,把 sinθ2作為變量,(12)式的振蕩周期為c2/(2fh) ,從圖10可以看出BL在大掠射角范圍內(nèi)關(guān)于 sinθ2均勻振蕩.

圖10 (10)式計算化簡后的BL (橫軸利用(7)式轉(zhuǎn)為 sinθ2)Fig.10.Curve of BL after simplified computed by equation(10) with the horizontal axis converted to sinθ2 by equa?tion (7).

4.3 利用海洋環(huán)境噪聲獲得的估計海底地聲參數(shù)

BL曲線的臨界角效應(yīng)和周期振蕩這兩種特征可以分別用來確定單沉積層海底的分層聲速和沉積層深度.在利用噪聲獲得時,由于陣元的個數(shù)限制,會引起條紋的模糊,圖11(a)仿真了42陣元,第一陣元設(shè)置在30 m處的,其他參數(shù)依舊沿用表1,與圖2陣元遍布全海深的相比,條紋顯得模糊很多.圖11(b)單獨提取出了1800 Hz下的曲線,該曲線中依然可以看到較為明顯的臨界角信息,即幅值發(fā)生躍變的點,通過兩處臨界角θc1,θc2可以利用上文所述的Snell定律,在折射角為0時確定沉積層聲速

和基底聲速

而在θc2后,存在一個明顯“凸點”,這個“凸點”不一定是這條曲線本身的極值點,但實際上對應(yīng)的是理想情況下反射損失BL的一個極值點,由于有限陣元帶來的模糊,使得這個極值點變?yōu)榱恕巴裹c”.假設(shè)“凸點”的位置對應(yīng)為θt,那么利用θt與θc2便可以得到(12)式在某一固定頻率下關(guān)于sinθ2的半周期

圖11 (a) 42陣元設(shè)置在水下30 m處 ;(b) 1800 Hz下 曲線Fig.11.(a) with 42 elements at the depth of 30 m;(b) the curve of under 1800 Hz.

從圖11(b)可以看出,θt約為46°,θc2約為41°,利用(17)式估計出沉積層厚度約為2.84 m,與實際仿真設(shè)置的3 m相比偏小,這是因為有限陣元帶來的模糊,即空間采樣點不足導(dǎo)致曲線中的極值點之間的分辨率下降.也可以利用理論陣元合成的方法對噪聲的互譜密度矩陣進(jìn)行擴展,以提高曲線的分辨率[17,18].

5 利用噪聲能流理論獲取海底參數(shù)海上試驗驗證

實驗數(shù)據(jù)于2016年12月在中國黃海某海區(qū)采集,當(dāng)?shù)貙嶋H水深35 m,冬季海水近似看作等聲速,聲速為1487 m/s,實測聲速剖面如圖12所示.實驗采用陣元間隔為0.2 m的25元船掛陣采集噪聲數(shù)據(jù),將實測數(shù)據(jù)分為20組,每組時長5 min,對每組數(shù)據(jù)分別進(jìn)行處理.圖13和圖14分別給出了第9組和第11組數(shù)據(jù)處理結(jié)果,1600 Hz下兩圖的曲線中明顯存在兩次突變,因此存在兩個臨界角,符合單層沉積層海底模型,可以利用上文所述的方法對海底分層結(jié)構(gòu)及地聲參數(shù)進(jìn)行估計.圖中已標(biāo)注出θc1,θc2,θt的位置,θc1,θc2的選取依據(jù)主要是根據(jù)曲線中幅值突變的程度,選取其突變的初始角度,θt選取的是θc2后的第一個“凸點”,然后利用(14),(15),(17)三式對該海區(qū)海底聲速和沉積層厚度進(jìn)行估計.圖15為20組數(shù)據(jù)計算結(jié)果的統(tǒng)計分布直方圖,其中圖15(a)為沉積層聲速,圖15(b)為基底聲速,圖15(c)為沉積層厚度.沉積層聲速的估計平均值為1544 m/s,這與同海區(qū)的反演結(jié)果1531 m/s相近[22];基底聲速的估計平均值為1935 m/s,相比沉積層聲速其波動范圍略大,主要集中在1900－1960 m/s之間,目前暫無其他方法對這一結(jié)果進(jìn)行對比;沉積層厚度的估計平均值為2 m,考慮到有限陣元帶來的模糊影響,實際沉積層厚度可能略高于估計值.

圖12 實驗實測海水聲速剖面Fig.12.Sound speed profile of the experiment.

圖13 第9組數(shù)據(jù)處理結(jié)果 (a) ;(b) 1600 Hz下的 曲線Fig.13.Results of the 9th set of data： (a) ;(b) curve of under 1600 Hz.

圖14 第11組數(shù)據(jù)處理結(jié)果 (a) ;(b) 1600 Hz下 曲線Fig.14.Results of the 11th set of data： (a) ;(b) curve of under 1600 Hz.

圖15 數(shù)據(jù)計算結(jié)果統(tǒng)計分布直方圖 (a)沉積層聲速;(b)基底聲速;(c)沉積層厚度Fig.15.Statistical distribution histogram of 20 sets of data：(a) Sound speed of the sediment;(b) sound speed of the basement;(c) thickness of the sediment.

為驗證所估計沉積層厚度的可靠性,試驗船搭載了淺地層剖面儀(sub?bottom profiler),該地層剖面儀發(fā)射調(diào)頻信號,利用聲波在水中和水下沉積層內(nèi)傳播和反射的特性來探測海底地層結(jié)構(gòu).圖16給出了從試驗點向東北方向掃描的海底地形,35 m深處為海底表層,在其下方存在沉積層,圖16(b)中用黑色標(biāo)出,厚度約為2.2 m,噪聲估計出的沉積層厚度與此結(jié)果基本符合.

圖16 淺地層剖面儀探測試驗海區(qū)海底結(jié)構(gòu)(圖(b)中用黑點標(biāo)明海底分界面)Fig.16.Sub?bottom structure of the experimental area de?tected by sub?bottom profiler.The interface is marked by black spot in panel (b).

6 結(jié) 論

本文根據(jù)Harrison能流理論,假設(shè)噪聲源均勻分布在海面,利用垂直陣接收的海洋環(huán)境噪聲數(shù)據(jù)對海底反射損失進(jìn)行估計,討論了無限大液體聲學(xué)半空間海底和存在一層沉積層兩種情況的曲線特點,并將近場連續(xù)譜和遠(yuǎn)場離散譜兩部分的貢獻(xiàn)分開考慮,從理想反射系數(shù)出發(fā),推導(dǎo)了各個海底參數(shù)對BL的影響.