閆循鵬，卜憲海，劉洪霞，辛明真，陽凡林，2

(1. 山東科技大學(xué) 測繪科學(xué)與工程學(xué)院，山東 青島266590；2. 海島(礁)測繪技術(shù)國家測繪地理信息局重點(diǎn)實(shí)驗(yàn)室，山東 青島266590)

聲速剖面插值方法在多波束測深中的應(yīng)用

閆循鵬1，卜憲海1，劉洪霞1，辛明真1，陽凡林1，2

(1. 山東科技大學(xué) 測繪科學(xué)與工程學(xué)院，山東 青島266590；2. 海島(礁)測繪技術(shù)國家測繪地理信息局重點(diǎn)實(shí)驗(yàn)室，山東 青島266590)

聲速剖面通常采用站位式測量得到，分布較為稀疏，即使是走航式測量也難以得到足夠密集的聲速剖面。為此應(yīng)在現(xiàn)有的聲速剖面基礎(chǔ)上進(jìn)行插值以滿足各方面的需要。反距離加權(quán)插值法、線性插值三角網(wǎng)法和基于經(jīng)驗(yàn)正交函數(shù)的聲速剖面反演法是聲速剖面插值常用的三種方法，文章討論了其在多波束測深聲速剖面加密中的適用性和插值精度。結(jié)果表明，三種聲速剖面插值方法均可提高聲速改正精度，其中基于經(jīng)驗(yàn)正交函數(shù)的聲速剖面反演法插值精度高，適用性好；線性插值三角網(wǎng)法內(nèi)插精度較高，但外推能力較差；反距離加權(quán)插值法算法簡單，便于實(shí)現(xiàn)。

聲速剖面插值；反距離加權(quán)插值法；線性插值三角網(wǎng)法；經(jīng)驗(yàn)正交函數(shù)

聲速改正是多波束測深數(shù)據(jù)后處理的關(guān)鍵步驟之一，其所采用的聲速剖面的準(zhǔn)確性直接影響著波束點(diǎn)歸位計算的精度[1-3]。目前國內(nèi)有關(guān)多波束測量的規(guī)范對于聲速剖面的測量、區(qū)域分布和站位密度沒有特殊要求[4-5]，通常情況下進(jìn)行聲速改正時采用“最近原則”法，即采用距離測深位置最近的聲速剖面。若聲速剖面采集站位布設(shè)不合理，或在聲速變化復(fù)雜地區(qū)，如大陸架淺水區(qū)、河口沖淡區(qū)，聲速剖面采集具有滯后性，將會引入較大的聲速剖面代表性誤差，從而影響多波束測深最終成果質(zhì)量，造成假地形[6]。為此，Cartwright[7]在聲速結(jié)構(gòu)變化較快地區(qū)采用走航式聲速剖面儀來獲取較為密集的聲速剖面，但此方法通常無法采集到表層和底層聲速值。關(guān)永賢[8]提出了在多波束數(shù)據(jù)后處理中利用現(xiàn)有聲速剖面采用反距離加權(quán)插值法(inverse distance weighting，IDW)進(jìn)行插值，得到了更為密集的聲速剖面。LeBlanc[9]等的研究表明，經(jīng)驗(yàn)正交函數(shù)(empirical orthogonal function，EOF)是描述聲速剖面有效的基函數(shù)，并利用此函數(shù)構(gòu)建了深海聲速場模型。我國學(xué)者[10-13]在基于EOF的聲速剖面反演法上也取得了一些研究成果，其中部分研究成果已應(yīng)用于多波束測深聲速改正中。聲速剖面插值雖在多波束測深數(shù)據(jù)處理中得到一定的應(yīng)用，但目前對各種方法的插值精度和適用性分析較少。本文以實(shí)測和模擬的聲速剖面數(shù)據(jù)為例，比較分析了反距離加權(quán)插值法、線性插值三角網(wǎng)法(triangulation with linear interpolation，TLI)和基于EOF的聲速剖面反演法在多波束測深聲速剖面加密中的插值精度和適用性，對多波束測深中聲速站的布設(shè)和聲速改正有一定的參考價值。

1 聲速誤差對多波束測深精度的影響

1.1 表層聲速誤差對多波束測深的影響

海水表層溫度和鹽度易受日、風(fēng)、海流等因素影響，表層聲速變化大，若未利用表層聲速儀實(shí)時采集表層聲速，會帶來表層聲速誤差，其對多波束測深的影響主要體現(xiàn)在波束指向角的準(zhǔn)確性上，從而影響多波束測深精度。

多波束平面陣換能器預(yù)形成波束時需要預(yù)知表層聲速，例如，當(dāng)預(yù)形成第k個波束時，波束指向角θk與表層聲速c0的關(guān)系為[14]：

(1)

式中，N為水聽器基元個數(shù)，k∈(-N/2，N/2)；f為水聽器頻率；l為基元間隔，通常取l=λ/2；λ為聲波波長，且c0=λf。從上式可以看出，當(dāng)換能器發(fā)射頻率一定時，預(yù)形成的波束指向角只與表層聲速c0相關(guān)?，F(xiàn)取實(shí)際表層聲速c0=1 500m/s，f=100 kHz，N=57，水深d=100 m的條件下預(yù)生成波束指向角，在對表層聲速c0分別添加±1、±2、±3m/s的誤差后，對實(shí)際波束指向角的影響如圖1(a)所示，其對多波束測深的影響如圖1(b)所示。

圖1 表層聲速誤差對多波束測深的影響

從圖1中可以看出，入射角為0°的中央波束指向角不受表層聲速誤差的影響，隨著波束入射角和表層聲速誤差的增大，波束指向角誤差也逐漸增大，在波束入射角大于60°時，波束指向角誤差增大速度加快。當(dāng)設(shè)定的表層聲速大于實(shí)際表層聲速時，波束指向角較原來增大，多波束條帶覆蓋寬度增加，測得的海底地形整體呈現(xiàn)出“哭臉”狀失真，即中央波束附近水深小于實(shí)際水深，邊緣波束水深大于實(shí)際水深；當(dāng)表層聲速小于實(shí)際表層聲速時，結(jié)果相反。此外，在表層聲速誤差為2m/s、指向角為75°時，波束指向角誤差達(dá)0.26°，水深誤差達(dá)1.39m，超過了1%水深的限差，故表層聲速誤差應(yīng)限定在1m/s內(nèi)。

1.2 聲速剖面代表性誤差對多波束測深的影響

相鄰的聲速剖面結(jié)構(gòu)具有相似性，但又無明顯的規(guī)律性，為了研究聲速剖面代表性誤差對多波束測深的影響，假設(shè)相鄰聲速剖面存在整體偏差，即聲速剖面各聲速節(jié)點(diǎn)加減同樣的值，聲速梯度保持不變，且不考慮表層聲速變化對波束指向角的影響。取一實(shí)測的聲速剖面，并對聲速剖面的各聲速節(jié)點(diǎn)分別添加±2、±4、±6m/s的誤差形成6個新的聲速剖面，預(yù)設(shè)水深d=54 m，結(jié)合常梯度聲線跟蹤算法，對比各聲速剖面對同一海底的聲速改正結(jié)果，如圖2所示。

圖2 聲速剖面代表性誤差對多波束測深的影響

由圖2可知，當(dāng)聲速剖面整體偏小時，聲速改正后的水深減小，邊緣波束略向下彎曲；當(dāng)聲速剖面整體偏大時，聲速改正后的水深增大，邊緣波束略向上彎曲。聲速剖面整體偏差越大，水深誤差越大，整體偏差6 m/s，水深誤差達(dá)0.24 m，已接近1%水深限差的一半。

可見，表層聲速誤差和聲速剖面代表性誤差均會影響多波束測深精度，且對邊緣波束影響更大，故聲速剖面?zhèn)€數(shù)采集不足或測區(qū)海洋環(huán)境復(fù)雜時應(yīng)采用合理的插值方法加密聲速剖面。

2 聲速剖面插值方法

2.1 反距離加權(quán)插值法

反距離加權(quán)插值法又稱為距離反比法，認(rèn)為被估未知點(diǎn)屬性與其周圍一定距離內(nèi)離散的已知點(diǎn)的屬性有關(guān)，利用這些已知點(diǎn)的數(shù)值進(jìn)行加權(quán)運(yùn)算即可以求取未知點(diǎn)數(shù)值，各已知點(diǎn)所給予的權(quán)重依據(jù)已知點(diǎn)與未知點(diǎn)的距離遠(yuǎn)近來確定，距離未知點(diǎn)越近的權(quán)重越大，即

圖3 線性插值三角網(wǎng)法示意圖

(2)

2.2 線性插值三角網(wǎng)法

線性插值三角網(wǎng)法需要將所有站位中相鄰的三個點(diǎn)連接成三角形，并且所有三角形的邊都不能與其余任意一個三角形相交叉，從而構(gòu)成一個覆蓋整個測區(qū)的三網(wǎng)形平面網(wǎng)。三角網(wǎng)中每個三角形內(nèi)的任意一點(diǎn)P(x,y)均位于該三角形三個頂點(diǎn)確定的平面上，其數(shù)值由三個頂點(diǎn)加權(quán)運(yùn)算求取，權(quán)系數(shù)由點(diǎn)P與三個頂點(diǎn)之間的線性關(guān)系確定。如圖3所示：

插值函數(shù)為：

(3)

式中，fi為頂點(diǎn)P(xi,yi)(i= 1,2,3)的已知數(shù)值，ωi為相應(yīng)的權(quán)系數(shù)：

(4)

由于確定權(quán)系數(shù)ωi的前提是待估點(diǎn)P位于三角形內(nèi)部，因此對于三角網(wǎng)外部的點(diǎn)線性插值三角網(wǎng)法不再適用。

2.3 基于EOF的聲速剖面反演法

基于EOF的聲速剖面反演法是一種對同一測區(qū)內(nèi)具有空間相關(guān)性的聲速剖面進(jìn)行特征向量分解，并結(jié)合采樣數(shù)據(jù)重構(gòu)該測區(qū)任意一點(diǎn)聲速剖面的反演方法。

原始的m個聲速剖面數(shù)據(jù)節(jié)點(diǎn)通常是非等深采樣，需對每一條聲速剖面進(jìn)行樣條插值得到n個垂直分布的等深節(jié)點(diǎn)，則內(nèi)插后的聲速剖面矩陣為：

(5)

C的協(xié)方差矩陣R：

(6)

對協(xié)方差矩陣R進(jìn)行特征分解得到：

RF=DF，

(7)

測區(qū)內(nèi)任意一點(diǎn)的聲速剖面c(z)可用前k階EOF表示為：

(8)

3 實(shí)例計算與分析

3.1 三種插值方法插值精度評定

2015年9月在黃海某區(qū)域共采集了14個站位的聲速剖面，站位分布及三角網(wǎng)法構(gòu)建的三角網(wǎng)如圖4所示，聲速剖面結(jié)構(gòu)如圖5所示?，F(xiàn)選取其中12條聲速剖面作為已知數(shù)據(jù)，2條聲速剖面a和b(站位如圖4所示)分別用于檢核插值算法的內(nèi)插與外推能力。

從圖4和圖5中可以看出，此次野外采集的聲速剖面站位分布得不是很均勻，聲速剖面結(jié)構(gòu)總體相似，但在水深約5 m處存在著聲速躍層，呈表層聲速變化較大、底層聲速穩(wěn)定變化、個別聲速剖面結(jié)構(gòu)與其他聲速剖面結(jié)構(gòu)差異略大的態(tài)勢，有明顯的時空分布差異。首先，對原始聲速剖面進(jìn)行等深節(jié)點(diǎn)插值處理，再分別利用上述三種方法對檢核剖面a和b所在位置進(jìn)行插值，并將插值結(jié)果與實(shí)測和最鄰近聲速剖面對比，如圖6所示。

圖4 聲速剖面站位分布

圖5 聲速剖面樣本結(jié)構(gòu)圖

圖6 實(shí)測聲速剖面與內(nèi)插聲速剖面對比

數(shù)值插值方法的好壞通常通過交叉驗(yàn)證進(jìn)行評價[15]。交叉驗(yàn)證的主要評定指標(biāo)有平均誤差、均方根誤差和標(biāo)準(zhǔn)殘差等，其中均方根誤差越小，待插點(diǎn)的預(yù)測值就越接近于其真實(shí)值，即插值方法越好。a和b兩站聲速剖面插值結(jié)果殘差統(tǒng)計如表1所示。

表1 殘差統(tǒng)計結(jié)果

通過上述圖6和表1可以發(fā)現(xiàn)：

1) 在內(nèi)插時,即a站位處，三種插值方法的內(nèi)插精度較采用“最近原則”均有所提高，其中基于EOF的聲速剖面反演法各指標(biāo)數(shù)值較小，即插值精度最高，反距離加權(quán)插值法和線性插值三角網(wǎng)法次之，插值精度大致相似。

2) 檢測剖面b位于構(gòu)建的三角網(wǎng)外，線性插值三角網(wǎng)法不再適用，退化為“最近原則”法，外推精度大為降低；反距離加權(quán)插值法均方根誤差最小，精度最高；基于EOF的聲速剖面反演法精度雖有所下降，但適用性仍然較好。

3) 將檢測剖面a、b插值結(jié)果對比可知，若插值處站位距離已知聲速站位較遠(yuǎn)，各插值方法的插值精度均會下降，采用“最近原則”法選取聲速剖面進(jìn)行聲速改正時，誤差會迅速增大，最大達(dá)2.05 m/s。因此，在多波束外業(yè)測量中，應(yīng)根據(jù)已有資料盡量滿足聲速站采樣的密度要求。

3.2 三種插值方法對多波束聲速改正的影響

為了分析三種插值方法對多波束聲速改正的影響，利用同一組多波束回波時間觀測序列，分別采用實(shí)測和插值得到的聲速剖面進(jìn)行聲線跟蹤，并以實(shí)測聲速剖面的聲線跟蹤結(jié)果為基準(zhǔn)，計算a、b兩站位內(nèi)插的聲速剖面引起的水平位移誤差和深度誤差，如圖7和圖8所示。

圖7 a站位三種聲速剖面插值方法對多波束聲速改正的影響

圖8 b站位三種聲速剖面插值方法對多波束聲速改正的影響

從圖7和圖8中可以看出，若采用傳統(tǒng)的“最近原則”法選取聲速剖面進(jìn)行聲速改正時，隨著波束入射角的增大，引起的多波束測深誤差也逐漸增加，造成條帶邊緣處地形上翹或下彎。在波束入射角大于65°時，a、b兩站最大測深誤差分別達(dá)到0.25和0.45 m，后者已接近《海道測量規(guī)范》限定的1%水深誤差。但在采用聲速剖面插值后，聲速改正精度得到提高，在波束入射角小于70°時水平位移誤差小于0.1 m，水深誤差小于0.15 m，滿足多波束測量精度要求。

4 結(jié)論

本研究通過對常用三種聲速剖面插值方法的實(shí)現(xiàn)難易程度、插值精度、適用性以及插值后的聲速剖面對多波束聲速改正的影響進(jìn)行詳細(xì)對比分析，得到以下結(jié)論：

1) 三種插值方法在一定程度上均可提高多波束聲速改正精度；

2) 基于經(jīng)驗(yàn)正交函數(shù)的聲速剖面反演法插值精度高，適用性好，使用較少的參數(shù)即可表示聲速與水深的復(fù)雜函數(shù)關(guān)系，便于聲速場的構(gòu)建；

3) 線性插值三角網(wǎng)法和反距離加權(quán)插值法算法簡單，便于實(shí)現(xiàn)；線性插值三角網(wǎng)法對于三角網(wǎng)內(nèi)的位置插值精度較反距離加權(quán)插值法高，但外推能力較差。

在實(shí)際工作中可根據(jù)具體情況選用某種插值方法，或綜合使用。但值得注意的是，聲速剖面插值只是對外業(yè)測量不足的一種補(bǔ)救措施，不應(yīng)過于依賴。

[1]李家彪,鄭玉龍,王小波，等.多波束測深及影響精度的主要因素[J]海洋測繪，2001(1)：26-32.

[2]陽凡林,韓李濤,王瑞富,等.多波束聲納水柱影像探測中底層水域目標(biāo)的研究進(jìn)展[J].山東科技大學(xué)學(xué)報(自然科學(xué)版),2013,32(6):75-83. YANG Fanlin,HAN Litao,WANG Ruifu,et al.Progress in object detection in middle and bottom-water based on multibeam water column image[J].Journal of Shandong University of Science and Technology(Natural Science),2013,32(6):75-83.

[3]周堅,周青,呂良,等.關(guān)于多波束聲速剖面改正問題的探討[J].海洋測繪,2014,34(4):62-65. ZHOU Jian,ZHOU Qing,Lü Liang,et al.Discussion of multibeam sound velocity profile correction[J].Hydrographic Surveying and Charting,2014,34(4):62-65.

[4]GB 12327-1998.海道測量規(guī)范[S].

[5]S-44. IHO standards for hydrographic surveys[S].

[6]GENG X,ZIELINSKI A. New methods for precise acoustic bathymetry[J]. Processing of the Canadian Hydrographic Conference,1998:255-256.

[7]CARTWRIGHT D. Multibeam bathymetric surveys in the Fraser river delta,managing severe acoustic refraction issues[D]. Fredericton:University of New Brunswick,2003.

[8]關(guān)永賢,屈小娟.多波束測深中聲速剖面的橫向加密方法[J].海洋測繪,2009,29(5):54-56. GUAN Yongxian,QU Xiaojuan.The interpolation of sound velocity profile in horizontal plane in multibeam echosounding[J]. Hydrographic Surveying and Charting,2009,29(5):54-56.

[9]LEBLANC L R,MIDDLETON F H. An underwater acoustic sound velocity data model[J].Journal of the Acoustical Society of America,1980,67(6):2055-2062.

[10]沈遠(yuǎn)海,馬遠(yuǎn)良,屠慶平. 聲速剖面的分層經(jīng)驗(yàn)正交函數(shù)表示[J].西北工業(yè)大學(xué)學(xué)報,2000,18(1):90-93. SHEN Yuanhai,MA Yuanliang,TU Qingping.On expression of ocean sound profile by layered empirical orthogonal function(EOF)[J].Journal of Northwestern Polytechnical University,2000,18(1):90-93.

[11]丁繼勝,周興華,唐秋華,等.多波束勘測聲速剖面場的EOF表示方法[J].武漢大學(xué)學(xué)報:信息科學(xué)版,2007,32(5):446-449. DING Jisheng,ZHOU Xinghua,TANG Qiuhua,et al.Expression of multibeam echo sounding sound velocity profile with empirical orthogonal functions[J].Geomatics and Information Science of Wuhan University,2007,32(5):446-449.

[12]孫文川,暴景陽,金紹華.聲速剖面EOF表示對多波束水深數(shù)據(jù)的影響研究[J].海洋測繪,2014,34(6):21-24. SUN Wenchuang,BAO Jingyang,JIN Shaohua.Effects of sound speed profiles represented by EOF to multibeam bathymetric data[J].Hydrographic Surveying and Charting,2014,34(6):21-24.

[13]成芳,陳爽,金紹華,等.陸緣深水區(qū)聲速剖面EOF延拓方法研究[J].海洋測繪,2016,36(2):26-29. CHENG Fang,CHEN Shuang,JIN Shaohua,et al.An expanding method of sound velocity profiles via EOF in marginal deepwater areas[J].Hydrographic Surveying and Charting,2016,36(2):26-29.

[14]朱小辰,肖付民,劉雁春,等.表層聲速對多波束測深影響的研究[J].海洋測繪,2007,27(2):23-25. ZHU Xiaochen,XIAO Fumin,LIU Yanchun,et al.Research on the influence of surface sound velocity in multibean echo sounding[J].Hydrographic Surveying and Charting,2007,27(2):23-25.

[15]顧春雷,楊漾,朱志春.幾種建立DEM模型插值方法精度的交叉驗(yàn)證[J].測繪與空間地理信息,2011,34(5):99-102. GU Chunlei,YANG Yang,ZHU Zhichung.Accuracy cross-validation of several interpolation methods of DEM[J]. Geomatics & Spatial Information Technology,2011,34(5):99-102.

(責(zé)任編輯：高麗華)

Application of Sound Velocity Profile Interpolation Methods in Multibeam Echosounding

YAN Xunpen1,BU Xianhai1,LIU Hongxia1,XIN Mingzhen1,YANG Fanlin1,2

(1.College of Geomatics, Shandong University of Science and Technology, Qingdao, Shandong 266590, China; 2. Key Laboratory of Surveying and Mapping Technology on Island and Reef, National Administration of Surveying， Mapping and Geoinformation, Qingdao, Shandong 266590, China)

Because of the few and scattered distribution of sound velocity profile which is usually obtained by station measurement, it is very difficult to obtain sound velosity profile with enough density even by underway measurement. To solve this problem, the existing sound velocity profile were interpolated to meet the needs of sound velocity correction in this paper. The applicability and accuracy of three commonly used methods in sound velocity profile encryption of multibeam echosounding were discussed, namely, the inverse distance weighted interpolation method, triangulation with linear interpolation, and the sound velocity profile inversion method based on the empirical orthogonal function. Results show that all the three methods can improve the sound velocity correction accuracy. The sound velocity profile inversion method based on the empirical orthogonal function has highest accuracy and better applicability. Triangulation with linear interpolation has higher interpolation accuracy but poor extrapolating ability. Inverse distance weighted interpolation method is simple and easy to implement.

sound velocity profile interpolation; inverse distance weighted interpolation method; triangulation with linear interpolation; empirical orthogonal function(EOF)

2016-07-28

國家自然科學(xué)基金項(xiàng)目(41376108)；測繪公益性行業(yè)科研專項(xiàng)經(jīng)費(fèi)資助項(xiàng)目(201512034)；海洋公益性行業(yè)科研專項(xiàng)經(jīng)費(fèi)資助項(xiàng)目(201305034)

閆循鵬(1990—)，男，山東濟(jì)寧人，碩士研究生，主要從事海洋測量方面的研究.E-mail:skdyxp@163.com 陽凡林(1974—)，男，湖北荊州人，教授，博士生導(dǎo)師，主要從事海洋測繪及GNSS應(yīng)用等方面的研究，本文通信作者.E-mail:yang723@163.com

P229

文章編號：1672-3767(2017)02-0027-07