海底控制點(diǎn)定位的半?yún)?shù)平差模型法

孫文舟,殷曉冬,暴景陽,曾安敏

1. 海軍大連艦艇學(xué)院軍事海洋與測繪系，遼寧 大連 116018； 2. 武漢大學(xué)測繪學(xué)院，湖北 武漢 430079； 3. 地理信息工程國家重點(diǎn)實(shí)驗室，陜西 西安 710054

高精度確定海底控制點(diǎn)的三維坐標(biāo)是建立海底大地控制網(wǎng)的關(guān)鍵環(huán)節(jié)[1]。船載GNSS結(jié)合水下聲學(xué)測距是實(shí)施水下聲標(biāo)(控制點(diǎn))交會定位的有效手段[2-10]，影響這種方法定位精度最主要的因素是聲速在時間和空間上復(fù)雜變化所引起的測距誤差。采用圓走航的方式可以保證水面觀測點(diǎn)相對于海底控制點(diǎn)對稱分布[11]，對提高控制點(diǎn)水平坐標(biāo)的精確度有重要幫助，然而垂直方向的解依然較差且不穩(wěn)定。

為了提高海底控制點(diǎn)垂直解的精確度，國內(nèi)外學(xué)者進(jìn)行了不同的嘗試。文獻(xiàn)[12]提出了利用觀測歷元之間差分的方法，這種方法可以有效地消除聲速誤差長周期項的影響，但是短周期項依然存在。文獻(xiàn)[13]利用傳播時間建立了聲速誤差模型，提出了在500 m水深情況下，不依賴聲速剖面的水下靜態(tài)目標(biāo)三維坐標(biāo)快速求解方法。文獻(xiàn)[14]認(rèn)為聲速誤差項為常數(shù)，通過南海的試驗驗證了該方法相比于傳統(tǒng)方法的有效性。近年來，壓力傳感器的測深精度不斷提高，又有學(xué)者提出了深度約束下的海底控制點(diǎn)坐標(biāo)確定方法[15-17]，進(jìn)一步提高了三維坐標(biāo)解的穩(wěn)健性和精度。

上述方法都是對聲速變化引起的測距誤差中可以參數(shù)化的部分進(jìn)行建模改進(jìn)。根據(jù)文獻(xiàn)[2,18]對水下聲學(xué)測距誤差的研究，由內(nèi)波引起的測距誤差形態(tài)復(fù)雜，難以用少量的參數(shù)表達(dá)，而引入?yún)?shù)過多又往往會導(dǎo)致法方程的病態(tài)甚至秩虧。因此，本文嘗試?yán)冒雲(yún)?shù)平差模型進(jìn)行控制點(diǎn)坐標(biāo)的解算，以期提高垂直解的精確度。

1 海底控制點(diǎn)坐標(biāo)半?yún)?shù)平差模型解算方法

1.1 聲速變化引起的測距誤差分析

影響聲速變化的主要因素是太陽輻射和海水運(yùn)動[19]。聲速剖面的變化包括日變化、季節(jié)變化以及附加之上的隨機(jī)擾動，聲速剖面的季節(jié)變化主要與太陽直射點(diǎn)在南北回歸線之間周期性的運(yùn)動有關(guān)，而聲速剖面日變化主要與太陽高度角的日變化有關(guān)，由外力擾動引起的海洋內(nèi)波造成水體微團(tuán)垂向運(yùn)動，使得上下層海水混合引起溫度的變化從而導(dǎo)致海水聲速變化[20]。文獻(xiàn)[2]對聲學(xué)測距誤差的研究表明，由聲速變化引起的測距系統(tǒng)誤差項包含兩部分，長周期項與短周期項，其中長周期誤差項與潮汐的影響有關(guān)，而短周期誤差項則主要與海洋內(nèi)波現(xiàn)象有關(guān)。之后文獻(xiàn)[18]在夏威夷島附近進(jìn)行了相關(guān)試驗，再一次驗證了這個結(jié)論，并且認(rèn)為聲速主要受溫度變化影響，鹽度變化對其影響微乎其微[18]。

為進(jìn)一步分析文獻(xiàn)[2,18]的試驗結(jié)論，研究聲速對測距誤差的影響規(guī)律，本文進(jìn)行了如下試驗，首先利用文獻(xiàn)[21]的標(biāo)準(zhǔn)方程構(gòu)建了背景聲速剖面[21-22]

(1)

式中，C1是聲道軸處的聲速；z1是聲道軸的深度；ε=7.4×10-3是擾動系數(shù)；B=1300為聲道厚度尺度；z是水層深度。

利用表層溫度的變化規(guī)律以及熱傳導(dǎo)方程構(gòu)建聲速時空變化的部分[23]

Δc=Δc0e-βzsinωt-βz-φ

(2)

式(2)可表示聲速剖面在不同時間不同深度的聲速日變化值，通過式(1)和式(2)可構(gòu)造一個穩(wěn)定變化的聲速剖面。

在海底控制點(diǎn)坐標(biāo)解算的過程中往往采用某一固定的聲速剖面，用上述方法構(gòu)建的聲速剖面簇的平均值作為固定聲速剖面，計算在水深為3000 m，聲速采樣間隔為5 m的條件下，聲速變化所引起的測距誤差。為了不引入入射角偏差而導(dǎo)致的誤差，仿真試驗采用垂直波束，設(shè)表層溫度日變化的振幅為0.2℃，每小時采樣一次，共采樣24次，表層海水聲速變化與測距誤差的計算結(jié)果如圖1所示。

圖1 表面聲速變化與測距誤差曲線Fig.1 Surface sound velocity variation and ranging error curve

溫度的日變化只影響表層海水聲速，圖1(a)是對表層海水等間隔取樣得到的海水聲速平均值在時間上的變化，圖1(b)是各采樣時間所對應(yīng)的測距誤差的變化。從圖中可以看出，測距誤差的變化規(guī)律與聲速的變化規(guī)律基本相似，周期相同但相位存在差異。

文獻(xiàn)[2,18]的試驗結(jié)果顯示，測距誤差長周期項的變化與潮汐的變化規(guī)律相近。這種現(xiàn)象是由于潮波的振動產(chǎn)生了日潮內(nèi)波或半日潮內(nèi)波，在潮內(nèi)波的影響下海水產(chǎn)生垂向運(yùn)動，從而導(dǎo)致了上層溫度較高的海水與下層溫度較低的海水進(jìn)行熱量交換，引起了海水溫度的變化，溫度變化是引起海水聲速變化最主要的因素，其周期與潮內(nèi)波的周期相同。測距誤差短周期項的變化與海洋內(nèi)波有關(guān)，與潮內(nèi)波的作用機(jī)理相同，海洋內(nèi)波造成了海水垂向運(yùn)動，引起了海水溫度的變化進(jìn)而導(dǎo)致聲速的變化。相比于聲速長周期項變化，短周期項變化較為復(fù)雜，因此與其對應(yīng)的測距誤差短周期項隨時間變化同樣具有復(fù)雜的形式，難以進(jìn)行參數(shù)化建模。

1.2 圓走航距離交會定位模型

水下控制點(diǎn)定位的觀測方程通常表示為

ρi=fXi,Xo+δρdi+δρvi+εi

(3)

式中，ρi是第i時刻換能器到海底應(yīng)答器的距離觀測值；f(Xi,Xo)是兩者之間的直線距離，Xi是第i時刻換能器的位置，Xo是海底應(yīng)答器的位置；δρdi是由應(yīng)答器電路延遲引起的系統(tǒng)性誤差，是可計算的改正量；δρvi是聲速變化及聲線彎曲引起的系統(tǒng)性誤差，為主要誤差；εi是第i時刻測量的隨機(jī)噪聲。

若δρdi已通過外部設(shè)備進(jìn)行改正，則式(3)線性化展開得[14]

(4)

(5)

若進(jìn)行了n次觀測，可將式(5)寫成矩陣形式

V=BX+δρv-L

(6)

由上文分析可知，由聲速變化引起的測距誤差變化復(fù)雜，尤其是測距誤差的短周期項，難以對其進(jìn)行參數(shù)化的建模，若將每次觀測的系統(tǒng)誤差δρv都作為未知參數(shù)，則此時未知參數(shù)的個數(shù)為n+3，而觀測方程的個數(shù)為n，因此式(6)無法得到未知參數(shù)的唯一解。本文采用半?yún)?shù)模型進(jìn)行控制點(diǎn)三維坐標(biāo)的解算。

1.3 控制點(diǎn)坐標(biāo)解算的半?yún)?shù)模型

根據(jù)式(4)，可寫出如下矩陣形式的觀測方程

L=Bx+S+Δ

(7)

式中，S=[δpv1δpv2…δpvn]T，是由聲速變化引起的測距系統(tǒng)誤差的n維未知向量。式(7)即為海底控制點(diǎn)坐標(biāo)解算的半?yún)?shù)模型[24-25]。

式(7)的誤差方程是

(8)

由于未知參數(shù)多于方程的個數(shù)，因此無法直接利用最小二乘原理進(jìn)行求解，需要選擇新的平差準(zhǔn)則[22]

VTPV+αSTRS=min

(9)

式中，α是光滑因子，對兩參數(shù)V和S起平滑作用；R是一個給定的正定矩陣，本文中選定R=GTG，其中

同時，若令所有歷元聲速變化引起的測距誤差和為0，則又可設(shè)置一個限制條件

eTS=0

(10)

利用式(8)—式(10)可構(gòu)造拉格朗日函數(shù)

2γeTS

(11)

根據(jù)所構(gòu)造的拉格朗日函數(shù)即可推導(dǎo)出式(8)—(9)的解，這里直接給出推導(dǎo)結(jié)果

(12)

將計算得到的坐標(biāo)改正數(shù)加上所設(shè)定坐標(biāo)初始值即可獲得海底控制點(diǎn)的三維絕對坐標(biāo)。

(13)

式中，δpvl是測距誤差的長周期項；δpv0是長周期項測距誤差的振幅；T是測距誤差的長周期項的周期；φ是初始相位。對長周期項參數(shù)化處理后，短周期項仍然認(rèn)為不可參數(shù)化，同樣根據(jù)式(12)計算坐標(biāo)改正數(shù)，加上坐標(biāo)初始值獲得控制點(diǎn)三維絕對坐標(biāo)。

2 仿真試驗分析

為了驗證半?yún)?shù)平差模型解算海底控制點(diǎn)坐標(biāo)的有效性，本文進(jìn)行了如下仿真試驗：水深設(shè)置為3000 m，應(yīng)答器布放于海底，測量船以應(yīng)答器為中心，半徑為3000 m的圓航跡航行，速度為6節(jié)(約3 m/s)，每10 s實(shí)施一次水聲測距。解算的坐標(biāo)系采用笛卡兒坐標(biāo)系，應(yīng)答器為原點(diǎn)，x軸指向東，y軸指向北，z軸指向天頂。測距誤差的仿真采用與文獻(xiàn)[12]相同的方法，表示為

(14)

式中，誤差共包含4項，第1項和第2項分別是聲速變化引起的短周期項和長周期項，第3項是測區(qū)相關(guān)性誤差，第4項為隨機(jī)性誤差。參照文獻(xiàn)[2]在北太平洋的試驗和文獻(xiàn)[18]在夏威夷島附近的試驗，測距誤差中的長周期項與潮汐變化規(guī)律近似，周期為12 h，振幅約20 cm，由內(nèi)波引起的測距誤差短周期項的周期從幾十分鐘到幾個小時不等，振幅約為12 cm，因此令c1=12 cm、c2=20 cm、c3=2 cm，Tw=20 min，隨機(jī)誤差滿足方差為5 cm的高斯分布，t0是初始時刻，t是任意時刻，x是t時刻測量船的三維坐標(biāo)，x′是海底應(yīng)答器的三維坐標(biāo)。除此之外，設(shè)置潮汐的周期為12 h，振幅為5 m，海面波浪的周期為20 s，振幅為2 m，測量船水平方向的定位精度為5 cm，垂直方向定位精度為10 cm。

試驗中將本文的方法與其他兩種方法進(jìn)行對比，方法1是直接對觀測值利用最小二乘進(jìn)行解算，方法2是采用文獻(xiàn)[12]所提出的差分算法進(jìn)行解算，本文方法為方法3。

第1次仿真試驗的采樣數(shù)為4320，采樣時長剛好為測距誤差長周期項的一個周期，試驗中分別考慮存在內(nèi)波和不存在內(nèi)波兩種情況，表1給出了不同方法解算所得控制點(diǎn)坐標(biāo)偏差的統(tǒng)計參數(shù)。

表1 不同方法解算所得控制點(diǎn)坐標(biāo)偏差的統(tǒng)計參數(shù)

對比水平方向和垂直方向坐標(biāo)偏差的精確度(MSE)可以發(fā)現(xiàn)，聲速誤差對控制點(diǎn)垂直解的精確度影響較大，而對水平解的精確度影響小，這與圓走航測距定位消除水平方向誤差的設(shè)計相吻合，也被國內(nèi)外研究所證實(shí)。對比兩個試驗采用差分算法計算的垂直坐標(biāo)偏差的精確度，可驗證文獻(xiàn)[12]的結(jié)論，即差分可以消除測距誤差中長周期項的影響，但是短周期項的影響依舊存在。當(dāng)存在內(nèi)波時，水平方向依舊可以達(dá)到較高精確度，但是垂直方向解的精確度會明顯變差。將3種方法進(jìn)行橫向比較，當(dāng)不存在內(nèi)波時，差分算法是3種方法中垂直解精確度最高的，高于半?yún)?shù)模型和最小二乘法，而在實(shí)際的海洋環(huán)境中，更普遍的情況是存在海洋內(nèi)波，此時，差分算法垂直解的精確度降低，而最小二乘法與半?yún)?shù)模型解的精確度較高。

第2次試驗的采樣數(shù)為720，采樣時長為測距誤差長周期項周期的1/6，試驗中同樣分別考慮存在內(nèi)波與不存在內(nèi)波兩種情況，不同方法解算控制點(diǎn)坐標(biāo)偏差的統(tǒng)計參數(shù)見表2。

表2 不同方法解算所得控制點(diǎn)坐標(biāo)偏差的統(tǒng)計參數(shù)

將試驗2與試驗1垂直坐標(biāo)偏差的精確度進(jìn)行對比可得，無論是否存在內(nèi)波，當(dāng)觀測時長不是測距誤差長周期項的整數(shù)倍時，最小二乘法垂直解的精確度都會降低。差分方法的結(jié)論與前一試驗一致，整體精確度降低，這主要與觀測數(shù)據(jù)的數(shù)量減少有關(guān)。半?yún)?shù)模型的垂直解則依舊可以得到較高的精確度。圖2—圖5中測距誤差是指采用文獻(xiàn)[12]的方法仿真得到的真實(shí)測距系統(tǒng)誤差，而誤差估值是指半?yún)?shù)模型在不同試驗條件下對每次測距誤差的估值，可以看出半?yún)?shù)模型可以較準(zhǔn)確地描述測距誤差的變化規(guī)律，對比圖2和圖4、圖3和圖5，可以發(fā)現(xiàn)在內(nèi)波存在試驗中，半?yún)?shù)模型對測距誤差估值的結(jié)果明顯發(fā)散，精度明顯下降，表現(xiàn)為圖3中測距誤差與誤差估值的中誤差分別為0.17和0.22 cm；圖5中測距誤差與誤差估值的中誤差分別為0.10和0.18 cm。對比表2中半?yún)?shù)模型在存在和不存在內(nèi)波兩種條件下垂直坐標(biāo)偏差的精確度，可以得出內(nèi)波存在時垂直解誤差的精確度由6.43增加到25.53，證明內(nèi)波會影響半?yún)?shù)模型解算的精確度，但依舊小于最小二乘法的253.37和差分算法的389.20。從表2中還可以看出，兩種條件下，雖然最小二乘法垂直解的精度高于半?yún)?shù)模型，但由于其準(zhǔn)確度較低，所以其精確度遠(yuǎn)小于半?yún)?shù)模型。第2次試驗對測距誤差長周期項進(jìn)行了參數(shù)化處理，結(jié)果顯示，無內(nèi)波時對振幅的估計為20.02 cm，而當(dāng)內(nèi)波存在時，對振幅的估計為15.92 cm，相對于所設(shè)定的值20 cm，估值的準(zhǔn)確性明顯降低，對應(yīng)垂直方向坐標(biāo)解的誤差也增大，這個結(jié)果表明，當(dāng)測距誤差同時包含參數(shù)項和非參數(shù)項時，非參數(shù)項誤差會影響對參數(shù)項估值的準(zhǔn)確性。

圖2 半?yún)?shù)模型對測距誤差的估值(無內(nèi)波)Fig.2 Ranging error estimated by semi-parametric model (no internal waves)

圖3 半?yún)?shù)模型對測距誤差的估值(有內(nèi)波) Fig.3 Ranging error estimated by semi-parametric model (with internal waves)

圖4 半?yún)?shù)模型對測距誤差的估值(無內(nèi)波)Fig.4 Ranging error estimated by semi-parametric model (no internal waves)

圖5 半?yún)?shù)模型對測距誤差的估值(有內(nèi)波) Fig.5 Ranging error estimated by semi-parametric model (with internal waves)

在控制點(diǎn)三維坐標(biāo)解算的過程中，最小二乘法和差分算法要想達(dá)到最優(yōu)解，都需要較為苛刻的條件，最小二乘法要求觀測的時長為測距誤差長周期項的整數(shù)倍，而差分算法要求不存在短周期項的測距誤差，前者顯然實(shí)際工程應(yīng)用價值較低，因為測量船連續(xù)圓走航12 h甚至24 h難度大，而后者所需要不存在內(nèi)波現(xiàn)象的海洋環(huán)境條件更加難以滿足。因此半?yún)?shù)模型的優(yōu)勢則明顯突出出來，它不但可以取得較為精確的坐標(biāo)解算值，而且對客觀環(huán)境的依賴度相對較低。

3 結(jié) 論

針對傳統(tǒng)方法確定海底控制點(diǎn)三維坐標(biāo)垂直解精確度差的問題，本文提出了利用半?yún)?shù)平差模型處理測距系統(tǒng)誤差的方法，通過分析和試驗得出了以下結(jié)論：

(1) 測距誤差的周期性變化與聲速周期性變化有關(guān)。海水通過渦動熱傳導(dǎo)和內(nèi)波所引起對流熱傳導(dǎo)，將上層海水熱量向下層海水傳遞，引起了溫躍層聲速的變化，從而產(chǎn)生了測距誤差周期性的變化規(guī)律。

(2) 由海洋內(nèi)波等隨機(jī)海洋現(xiàn)象引起的測距誤差短周期項(非參數(shù)項)不會影響控制點(diǎn)坐標(biāo)水平解的精確度，但會影響控制點(diǎn)坐標(biāo)垂直解的精確度。

(3) 利用半?yún)?shù)平差模型可以有效處理由聲速變化引起的測距系統(tǒng)誤差，從而得出優(yōu)于傳統(tǒng)方法的海底控制點(diǎn)三維坐標(biāo)解算結(jié)果。