明鋒， 楊元喜， 曾安敏

地理信息工程國家重點(diǎn)實(shí)驗(yàn)室，西安測(cè)繪研究所，西安 710054

0 引言

高精度海底大地基準(zhǔn)是研究海底地殼運(yùn)動(dòng)、水下導(dǎo)航定位、油氣勘探和環(huán)境監(jiān)測(cè)等海洋科學(xué)和海洋工程的基礎(chǔ)(楊元喜等，2017).全球?qū)Ш叫l(wèi)星系統(tǒng)-聲學(xué)測(cè)距組合定位系統(tǒng)(Global Navigation Satellite System-Acoustic ranging，GNSS-A)可對(duì)海底參考點(diǎn)精確測(cè)量，能夠?qū)⒋蟮鼗鶞?zhǔn)從陸地?cái)U(kuò)展至海底.GNSS-A技術(shù)主要包含水面GNSS動(dòng)態(tài)定位和水下聲學(xué)測(cè)距，最早由斯克里普斯海洋研究所的Spiess于20世紀(jì)80年代所提出(Spiess，1985)，之后其研究團(tuán)隊(duì)進(jìn)一步提高了定位精度和觀測(cè)效率，并給出了一些相關(guān)試驗(yàn)結(jié)果(Spiess et al.，1998；Kussat et al.，2005).

從20世紀(jì)90年代后期開始，日本的名古屋大學(xué)(Kimura et al.，2019)、東北大學(xué)(Kido et al.，2006；Honsho and Kido，2017；Honsho et al.，2019)以及日本東京海岸警衛(wèi)隊(duì)和東京大學(xué)聯(lián)合體(Fujita et al.，2006；Ikuta et al.，2008；Sato et al.，2011，2013；Watanabe et al.，2014；Yokota et al.，2016)等研究團(tuán)隊(duì)持續(xù)改進(jìn)GNSS-A技術(shù)，建立并維護(hù)較大規(guī)模的海底觀測(cè)網(wǎng)，以監(jiān)測(cè)環(huán)太平洋地震帶的板塊構(gòu)造運(yùn)動(dòng).法國科學(xué)家也于2014年在馬爾馬拉海中部的北安納托利亞斷層布設(shè)了由10個(gè)應(yīng)答器組成的聲學(xué)測(cè)距網(wǎng)絡(luò)以持續(xù)監(jiān)測(cè)斷層活動(dòng)(Lange et al.，2019).

近年來，國內(nèi)也有相關(guān)學(xué)者對(duì)GNSS-A開展了理論和實(shí)驗(yàn)研究，得出了一些有益結(jié)論(韓云峰等，2017；趙爽等，2018；Chen et al.，2020).Yang等(2020)介紹了海底大地控制網(wǎng)建設(shè)關(guān)鍵技術(shù)，在3 km深海建立了長期大地基準(zhǔn)站，并進(jìn)行了GNSS-A觀測(cè)，初步結(jié)果表明定位精度優(yōu)于5 cm.之后，進(jìn)一步的精化處理表明坐標(biāo)分量內(nèi)符精度優(yōu)于0.4 cm，單程測(cè)距的均方誤差為11 cm(Yang and Qin，2021).我國臺(tái)灣省中山大學(xué)也致力于在沖繩海槽弧后盆地西端建設(shè)GNSS-A觀測(cè)網(wǎng)(Chen et al.，2018).

本文簡(jiǎn)要介紹基于貝葉斯框架的GNSS-A數(shù)據(jù)處理模型，然后對(duì)某海域2019年深海GNSS-A觀測(cè)數(shù)據(jù)進(jìn)行重處理，重點(diǎn)分析航跡幾何結(jié)構(gòu)、參數(shù)初值以及參數(shù)估計(jì)策略對(duì)水下應(yīng)答器三維定位精度的影響.

1 貝葉斯解算模型

1.1 觀測(cè)方程

圖1為GNSS-A基本原理和本次試驗(yàn)應(yīng)答器配置示意圖.在水上單元，沿測(cè)線航行的測(cè)量船的精確位置通過GNSS動(dòng)態(tài)定位來解算，結(jié)合GNSS天線與換能器之間的臂長向量和測(cè)量船姿態(tài)向量，可得到換能器的絕對(duì)位置.在水下單元，換能器與應(yīng)答器之間的幾何距離使用聲學(xué)測(cè)距單元測(cè)量.將換能器位置、傳播時(shí)延和聲速剖面數(shù)據(jù)相結(jié)合，可較為精確地確定虛擬參考點(diǎn)(應(yīng)答器的質(zhì)心)在外部參考框架下(如國際地球參考框架(International Terrestrial Reference Frame，ITRF))的坐標(biāo)，進(jìn)而實(shí)現(xiàn)了海陸基準(zhǔn)統(tǒng)一.

圖1 GNSS-A基本原理示意圖Fig.1 Schematic of the basic principle of GNSS-A

在水下聲學(xué)測(cè)距中，若已知第j個(gè)海底應(yīng)答器的位置Xj，換能器位置P(t)、四維聲速場(chǎng)V(e,n,u,t)，在聲速水平分層假設(shè)下，設(shè)第i個(gè)聲脈沖與第j個(gè)應(yīng)答器的往返傳播時(shí)間為Ti，則Ti可依據(jù)Snell定律計(jì)算(Chadwell and Sweeney，2010)：

(1)

其中，ti+和ti-分別是第i個(gè)聲脈沖的發(fā)送和接收時(shí)刻，P(ti+)、P(ti-)分別是換能器在ti+、ti-時(shí)刻的位置.

P(t)=Q(t)+R(Θ(t))·M,

(2)

式中，M定義在船載坐標(biāo)系下，船艏為X，右舷為Y，垂直向上為Z，可事先通過經(jīng)典測(cè)量手段獲得.R(Θ)是將船載坐標(biāo)系轉(zhuǎn)換為全球坐標(biāo)系的轉(zhuǎn)換矩陣，即

海水中的聲速是海水的溫度、鹽度和靜壓力的函數(shù)(Del Grosso，1974).本文假設(shè)聲速隨時(shí)間平滑變化并呈水平分層結(jié)構(gòu).這意味著在某一時(shí)刻朝不同水下應(yīng)答器傳播的聲信號(hào)具有相同的聲速剖面.然而，由于海洋洋流、潮汐等因素，實(shí)際海水的溫度和鹽度在水平方向上的分布并不均勻，導(dǎo)致聲速剖面具有垂向和水平方向的時(shí)空變化.由于大部分聲速變化都局限在沿聲射線路徑的淺層部分，因此四維聲速場(chǎng)空間變化可以近似為時(shí)間變化，故可將V(e,n,u,t)分解為一個(gè)水平分層的靜態(tài)聲速剖面和一個(gè)隨時(shí)間變化的擾動(dòng)，即(Ikuta et al.，2008)：

exp(-γi)·τi(P(ti+),P(ti-),Xj,V0(u)),

(3)

式中，τi表示第i個(gè)聲脈沖依據(jù)參考聲速剖面V0(u)利用聲線跟蹤計(jì)算得到的模型值，V0(u)是深度的分段線性函數(shù)，γi和Xj為待估參數(shù).當(dāng)姿態(tài)值的變化足夠大或臂長M測(cè)量精度不高時(shí)，也可以指定M為待估參數(shù).exp(-γi)為改正系數(shù)，表征聲速場(chǎng)的空間和時(shí)間擾動(dòng)引起的實(shí)際時(shí)延與理論時(shí)延的差異率.γi定義為γi≡0.5×∑l=i+,i-Γ(tl,P(tl),Xj)，即聲脈沖發(fā)射和接收時(shí)刻聲速擾動(dòng)的平均，函數(shù)Γ(t,P,X)稱為聲速攝動(dòng)模型，其計(jì)算公式為：

對(duì)比Yang和Qin(2021)的觀測(cè)方程(公式2)，式(5)本質(zhì)上仍是關(guān)于傳播時(shí)延的觀測(cè)方程，觀測(cè)量物理意義明確，避免了由于幾何距離計(jì)算不準(zhǔn)確而產(chǎn)生的模型誤差.如前所述，本文假設(shè)聲速水平分層結(jié)構(gòu)，然而如果聲速結(jié)構(gòu)中存在水平梯度，則梯度較慢(較快)一側(cè)的應(yīng)答器的位置會(huì)被錯(cuò)誤地確定為較遠(yuǎn)(較近)的位置(Tomita et al.，2019).因此，式(4)中的水平梯度是GNSS-A水下應(yīng)答器的水平位置估計(jì)中一個(gè)具有海洋學(xué)意義的系統(tǒng)偏差，解算出的水平梯度可以與物理海洋模型給出的結(jié)果相互驗(yàn)證.故利用式(5)既能夠提取淺層海水的聲速水平梯度信息，又能夠解算深層海底應(yīng)答器的幾何位置.

1.2 參數(shù)估計(jì)

式(5)可重新改寫為：

y=f(x)+e,

(6)

(7)

(8)

(9)

式中，m是模型參數(shù)的個(gè)數(shù).

由于μt,μMT、ρ2,λ2分別為觀測(cè)數(shù)據(jù)誤差方差協(xié)方差陣或參數(shù)先驗(yàn)分布中包含的未知參數(shù)，稱為超參數(shù)，可利用赤池貝葉斯信息量準(zhǔn)則(Akaike′s Bayesian Information Criterion，ABIC)(Akaike，1980)確定.本文采用網(wǎng)格搜索法取ABIC最小時(shí)的結(jié)果為最優(yōu)參數(shù)估值.

關(guān)于超參數(shù)確定和貝葉斯估計(jì)詳細(xì)推導(dǎo)過程，可參考Watanabe等(2020)，此處不再贅敘.

2 數(shù)據(jù)處理方案

2.1 試驗(yàn)基本情況

試驗(yàn)區(qū)域位于水深約3 km的某海域.2019年7月試驗(yàn)共布設(shè)5個(gè)應(yīng)答器，如圖2所示.其中4個(gè)應(yīng)答器(M01-M04)為錨系結(jié)構(gòu)，分別位于正方形4個(gè)頂點(diǎn)，正方形邊長約2 km；另1個(gè)應(yīng)答器(M05)為固定站，近似位于正方形的幾何中心.

圖2 5個(gè)應(yīng)答器(M01—M05)點(diǎn)位概略分布及走航構(gòu)型橫軸為E分量，縱軸為N分量(單位：m).Fig.2 The general distribution of 5 transponders (M01—M05) and surveying geometry configuration The horizontal axis is the E component, and the vertical axis is the N component (unit: m).

此次試驗(yàn)采用“繞行觀測(cè)”，測(cè)量航跡近似為圓形(大圓半徑約3 km，小圓半徑約1.5 km)和井字形(間距約1 km)(見圖2).為減弱多普勒效應(yīng)，觀測(cè)過程中船速不超過4 kn.船舶配備Veripos公司的LD7接收機(jī)，GNSS天線位置采用Veripos高精度定位系統(tǒng)Apex結(jié)果，該系統(tǒng)標(biāo)稱的定位精度水平分量≤5 cm，垂向分量≤12 cm.采用IXSEA Octans光纖羅經(jīng)測(cè)量船體姿態(tài)角，測(cè)量精度航向角0.1°，橫搖和俯仰角0.01°.

走航過程中采用AML Minos X設(shè)備約每10小時(shí)進(jìn)行一次往返聲速剖面測(cè)量，AML Minos X的標(biāo)稱精度為±0.006 m·s-1，精確度為±0.025 m·s-1(Ocean Scientific International Ltd，2015).整個(gè)觀測(cè)過程中共進(jìn)行4次往返聲速剖面測(cè)量，如圖3所示.從圖3可以看出，雖然4次聲速剖面整體趨勢(shì)一致，但在淺層海水(深度小于1 km)中單次往返測(cè)量聲速變化較大，同一深度最大差異達(dá)2 m·s-1，這已遠(yuǎn)超標(biāo)稱精度.

圖3 不同時(shí)段觀測(cè)的垂向聲速剖面(深度截止至約1000 m)其中紅色為前半段觀測(cè)，藍(lán)色為后半段觀測(cè).Fig.3 The vertical sound speed profile observed at different time intervals (the depth cuts off to about 1000 m) The red line denotes the first half of the observation, and the blue line denotes the second half of the observation.

2.2 參數(shù)估計(jì)策略

由1.2節(jié)可知，參數(shù)估計(jì)采用迭代法進(jìn)行數(shù)值求解，當(dāng)模型的非線性較強(qiáng)時(shí)，若初值選擇不當(dāng)，則最終的收斂結(jié)果可能并不是全局最優(yōu).另外，Watanabe等(2020)對(duì)應(yīng)答器陣列情形進(jìn)行了分析，對(duì)單個(gè)應(yīng)答器定位的適用性和穩(wěn)定性則沒有涉及.

本文結(jié)合本次試驗(yàn)實(shí)際，分單站(固定站)和多站(固定站+錨系站)兩種模式對(duì)上述算法的性能進(jìn)行分析.對(duì)于單站模式，主要分析：(1)參數(shù)初始值對(duì)結(jié)果的影響；(2)走航圖形結(jié)構(gòu)對(duì)結(jié)果的影響. 單站模式下的參數(shù)估計(jì)策略分為：

策略1：不估計(jì)臂長參數(shù)M，參數(shù)初值不加偏差；

策略2：不估計(jì)臂長參數(shù)M，參數(shù)初值加0.5 m偏差；

策略3：估計(jì)臂長參數(shù)M，且參數(shù)初值不加偏差；

策略4：估計(jì)臂長參數(shù)M，且參數(shù)初值加0.5 m偏差.

需要指出的是，由于垂向分量與垂向聲速剖面強(qiáng)相關(guān)，因此上述估計(jì)策略施加0.5 m偏移量僅針對(duì)水平分量.

對(duì)于多站模式，主要分析錨系站的定位精度，估計(jì)策略區(qū)分是否估計(jì)參數(shù)M.

表1 應(yīng)答器初始位置Table 1 The initial position of the transponder

表2 不同走航圖形結(jié)構(gòu)Table 2 Different surveying geometry configurations

3 結(jié)果分析

3.1 單站模式的結(jié)果

表3 單站模式的結(jié)果(M05)(單位：m)Table 3 The result of single station mode (M05) (unit: m)

3.2 多站模式的結(jié)果

多站模式下各應(yīng)答器ENU三個(gè)分量的改正量列于表4.

表4 多站模式的結(jié)果(單位：m)Table 4 The result of multi-station mode (unit: m)

3.3 分析與比較

3.3.1 走航幾何構(gòu)型對(duì)結(jié)果的影響

從表3可以看出，對(duì)于單站模式：

(1) 大圓(方案1)和小圓(方案2)的定位精度相當(dāng)，但是大圓+小圓(方案3)的精度略高；井字形測(cè)線三個(gè)分量定位精度最低(方案4)，加入圓形測(cè)線后精度有所改善(方案5)；方案6的水平分量定位精度最高，但垂向分量誤差較大；方案7的定位精度僅次于方案6，但其垂向精度與方案1-3相當(dāng).

(2) 整體上來看，圓形測(cè)線定位精度要優(yōu)于井字測(cè)線，特別是對(duì)于垂向分量；而圓形和井字測(cè)線組合有助于增強(qiáng)水平分量定位的魯棒性.對(duì)于單一入射角的圓形測(cè)線，由于對(duì)垂向分量的觀測(cè)孔徑過小，無法分辨垂向分量改正量和聲速偏差，從而造成兩者強(qiáng)相關(guān)；而不同入射角組合的大圓+小圓測(cè)線(方案3)則能有效分辨應(yīng)答器垂向分量改正量和聲速偏差，因此垂向分量精度較高(李昭，2016).這從另一方面也說明，高精度垂向聲速剖面是精確確定垂向分量的前提條件(趙建虎等，2016).

(3) 井字形測(cè)線定位精度最低，原因之一可能是走航時(shí)間相對(duì)較長(約20 h)，聲速剖面采樣頻率偏低，聲速代表性誤差較大(趙建虎和梁文彪，2019)；另一個(gè)可能性較大的原因是，由于內(nèi)波或其他海洋學(xué)因素導(dǎo)致淺層海水聲速時(shí)空變化較大，而走航時(shí)間已遠(yuǎn)大于淺層海水的水平聲速結(jié)構(gòu)持續(xù)時(shí)間，導(dǎo)致同時(shí)求解應(yīng)答器位置和聲速水平梯度的方法已不再適用(Iinuma et al.，2021).

3.3.2 應(yīng)答器位置初值對(duì)結(jié)果的影響

從表3可以看出，對(duì)于單站模式：

(1)整體來上說參數(shù)初值對(duì)結(jié)果影響較大，方案6(大圓+小圓+井字)水平分量對(duì)參數(shù)初值依賴性均較低；方案4(井字)影響最大，如前所述，這可能是由于同時(shí)求解應(yīng)答器位置和聲速水平梯度的方法已不再適用，此時(shí)位置參數(shù)初值施加偏差對(duì)模型求解更加不利.

(2)對(duì)于圓形航跡(方案1—3)，當(dāng)初值較精確時(shí)，結(jié)果與初值偏差不大，但當(dāng)初值顯著偏離時(shí)，模型解算不能正確收斂到合理的估計(jì)值附近，導(dǎo)致時(shí)延殘差增大.這表明，雖然圓形航跡水平對(duì)稱，但參數(shù)初始值誤差仍可導(dǎo)致應(yīng)答器水平分量系統(tǒng)性偏差.系統(tǒng)性偏差與初值偏差的定量關(guān)系可通過模擬計(jì)算分析得到，限于篇幅，本文不做進(jìn)一步的討論.

(3)雖然本文僅對(duì)水平分量施加偏移量，但從策略2、4可以看出，在包含井字測(cè)線的方案中其對(duì)垂向分量仍有較大影響(見表3，方案4—6)，其原因還需進(jìn)一步分析.

3.3.3 多站模式定位精度

從表4可以看出，相對(duì)于單站模式，多站模式的定位精度均較差，這主要是由于錨系站自身運(yùn)動(dòng)，而定位時(shí)則又假定其位置不變所導(dǎo)致.實(shí)際計(jì)算中方案4、5最終迭代已不能收斂，導(dǎo)致估計(jì)結(jié)果偏差較大.同時(shí)注意到，解算過程中由于假設(shè)固定站與錨系站存在相關(guān)性，因而導(dǎo)致固定站(M05)的定位結(jié)果相對(duì)單站模式偏差也較大.

對(duì)比錨系站(M01—M04)和固定站(M05)，可以發(fā)現(xiàn)不同方案中固定站的定位結(jié)果差異相對(duì)較小，水平分量差異最大為21.90 m(方案4)，最小為0.27 m(方案1)；垂向分量最大為7.76 m(方案4)，最小為0.05 m(方案1).這也說明對(duì)于高精度海底基準(zhǔn)站定位應(yīng)采取固定站模式，盡量減少因應(yīng)答器自身運(yùn)動(dòng)而引起的誤差.

3.3.4 估計(jì)臂長參數(shù)對(duì)結(jié)果的影響

對(duì)于多站模式，從表4可以看出，幾乎在所有方案中，估計(jì)臂長參數(shù)對(duì)位置分量改善并不明顯.這表明此種模式下，臂長參數(shù)已無法完全補(bǔ)償因錨系站運(yùn)動(dòng)、聲速剖面不準(zhǔn)確等引起的系統(tǒng)誤差，因而必須對(duì)錨系站運(yùn)動(dòng)建模或?qū)r(shí)延殘差進(jìn)行函數(shù)擬合，以盡可能削弱系統(tǒng)誤差的影響(Yang and Qin，2021).

4 討論

4.1 定位精度估計(jì)

Yamada等(2002)采用模擬數(shù)據(jù)分析了單海底應(yīng)答器的定位精度，結(jié)果表明在深海(～3 km)環(huán)境下，如果能夠獲得每條測(cè)線的精確聲速，可以將與聲學(xué)定位相關(guān)的誤差抑制在大約18 cm以內(nèi)(1倍中誤差)；采用實(shí)測(cè)數(shù)據(jù)分析表明單應(yīng)答器定位精度約為28 cm.考慮到本文貝葉斯估計(jì)中并沒有納入GNSS定位誤差和姿態(tài)傳感器誤差，在施加3 m先驗(yàn)約束的情形下，應(yīng)答器位置后驗(yàn)精度約為0.3 m，該精度與Yamada等(2002)的結(jié)果相當(dāng).

另一方面，觀測(cè)時(shí)延殘差與超參數(shù)的選擇密切相關(guān)，且僅能反映模型的內(nèi)符精度，而不能反映定位的準(zhǔn)確度.本文算例中，如圖4所示，單站模式的觀測(cè)時(shí)延殘差均方差(RMS)約為0.6 ms，等價(jià)于0.9 m(假設(shè)聲速為1500 m·s-1)；對(duì)于多站模式，觀測(cè)時(shí)延殘差RMS約為4.5 ms，等價(jià)于6 m.從圖4還可看出，對(duì)比固定站，多站模式各應(yīng)答器殘差均呈非高斯分布特征，表明仍存在潛在的系統(tǒng)誤差.Yang和Qin(2021)采用擬合殘差的方法減弱系統(tǒng)誤差的影響，因此其結(jié)果不確定度均遠(yuǎn)小于本文的結(jié)果.

圖4 單站模式和多站模式方案一的傳播時(shí)延殘差直方圖Fig.4 The histograms of delay residual of single station mode and scheme 1 in multi-station mode

4.2 算法的數(shù)值穩(wěn)定性和效率

在單站模式下，由于海底觀測(cè)網(wǎng)沒有網(wǎng)形約束，幾何構(gòu)型差，不合適的超參數(shù)易引起對(duì)角陣奇異，導(dǎo)致迭代計(jì)算失?。欢诙嗾灸Ｊ较?，應(yīng)答器觀測(cè)時(shí)延之間的相關(guān)性有助于抑制方差協(xié)方差陣奇異，提高算法的數(shù)值穩(wěn)定性，但同時(shí)也增加了一定的計(jì)算量.

在多站模式中，由于觀測(cè)量較多，矩陣求逆占用相當(dāng)大的計(jì)算資源，在保證幾何構(gòu)型對(duì)稱分布的前提下，可對(duì)原始觀測(cè)數(shù)據(jù)降低采樣，以提高計(jì)算速度.

5 結(jié)論

本文利用貝葉斯估計(jì)對(duì)3 km深海GNSS-A實(shí)測(cè)數(shù)據(jù)的定位性能進(jìn)行了分析，主要結(jié)論如下：

(1) 對(duì)于單個(gè)固定應(yīng)答器，走航時(shí)間較短的圓形測(cè)線的定位精度最高，而走航耗時(shí)較長的井字形定位精度最低；圓形測(cè)線對(duì)參數(shù)初值要求較高.

(2) 從本文的算例可以看出，對(duì)于圓形測(cè)線或井字測(cè)線，估計(jì)臂長參數(shù)，在一定程度上都有可能提高定位結(jié)果的可靠性.然而，臂長參數(shù)估計(jì)與航跡幾何構(gòu)型有關(guān)，對(duì)于有些圖形，估計(jì)臂長參數(shù)的能力不夠，即引入不適定性問題，將會(huì)導(dǎo)致結(jié)果不可靠.本文算例中，由于臂長參數(shù)的先驗(yàn)精度較高，在測(cè)線對(duì)稱分布條件下，估計(jì)該參數(shù)并不能顯著提升定位精度.

(3) 固定站的定位精度遠(yuǎn)高于錨系站，應(yīng)答器陣列能夠顯著增強(qiáng)解算的穩(wěn)定性，但會(huì)增加計(jì)算量.

(4) 同時(shí)求解應(yīng)答器位置和聲速水平梯度的方法僅適用于聲速水平梯度變化不大的情形，因此若采用此方法進(jìn)行海底應(yīng)答器定位，建議盡量減少走航時(shí)間，以滿足水平聲速結(jié)構(gòu)近似保持不變的假設(shè).

由于本文僅對(duì)圓形、井字形測(cè)線及其組合進(jìn)行了分析，限于觀測(cè)條件沒有對(duì)其他類型測(cè)線進(jìn)行分析.另外，由于GNSS-A走航觀測(cè)需耗費(fèi)大量船時(shí)，如何確定最優(yōu)的海底應(yīng)答器分布以及對(duì)應(yīng)的走航幾何構(gòu)型，既能提高觀測(cè)效率，又能保證定位精度仍需進(jìn)一步分析研究.