梅斌,孫立成,史國友,馬麟,李偉峰

(1.大連海事大學 a.航海學院； b.遼寧省航海安全保障重點實驗室，遼寧 大連 116026； 2.中國船級社，北京 100007)

0 引 言

在船舶航行中，由于海洋環(huán)境干擾的存在或者避讓操作需要，船舶會偏離計劃航線一定距離，形成偏航距離(cross track distance)。[1-2]偏航距離在船舶運動控制、靠泊、進出港操縱以及航道規(guī)劃和航線設計等領域都有應用。在船舶航跡保持控制中，船舶實際位置與計劃航線間的偏航距離是航跡保持控制的誤差[3]，該誤差一般作為控制調(diào)節(jié)量。對于船舶進出港操縱，船舶與岸線間的偏航距離可用于數(shù)字化精確引航。[4]在航道規(guī)劃和航線設計中，可根據(jù)偏航距離，在鉆井平臺附近、淺水區(qū)域設置安全緩沖區(qū)域。因此，研究偏航距離的計算方法具有一定的實用意義。

關于偏航距離計算方法，已有許多研究成果。LEE等[5]采用CCD攝像機影像辨別岸基導標顏色，計算船位與導標延伸線之間的偏航距離，實現(xiàn)了視覺導航，用于船舶自動靠泊?；贓CDIS的導航系統(tǒng)一般利用GPS提供的位置信息，計算船舶與岸線之間的偏航距離。[6-7]在汽車自動行駛、飛機航行和水下航行器控制中，需要計算偏航距離[8-12]：汽車的偏航距離定位儀器處于歐氏空間，使用歐氏距離；在飛機航行中，飛行員計算飛機偏航距離時將地球作為球體近似體，精度不高；水下航行器作業(yè)范圍較小，不需要進行大范圍的航線設計，因此水下航行器的偏航距離也是歐氏距離。以上計算方法不能用于船舶偏航距離的計算，這是因為它們采用的大地球體近似體不同。

在船舶航跡保持控制中:FOSSEN[3]采用平面坐標系上的歐氏距離作為控制量，并未論述大地坐標轉換，因此不能直接用地理坐標進行航跡偏差度量；張顯庫等[13]和鄭云峰等[14]采用墨卡托投影，將船舶的大地坐標投影到平面直角坐標系內(nèi)，以歐氏距離作為航跡保持控制的航跡偏差(簡稱基于墨卡托投影的歐氏距離算法，算法1)；在衛(wèi)星定位與導航中，使用地心地固坐標系(ECEF)作為慣性坐標系、北-東-地坐標系(NED)作為附體坐標系時，并沒有考慮球面上船舶軌跡的彎曲性質(zhì)[3]。在船舶航行中，船員習慣使用墨卡托投影進行導航，但是墨卡托投影存在長度變形特性，這個特征在高緯度區(qū)域更加顯著，故當緯度大于86°時不適合使用墨卡托投影。為研究高緯度區(qū)域船舶運動的偏航距離計算誤差，本文進行偏航距離概念探討和真值求解算法設計，分析多種算法的誤差大小和規(guī)律；通過分析地圖投影特點，選取適合高緯度區(qū)域偏航距離求解的投影算法；在此基礎上，根據(jù)航海領域的精度需求，推薦偏航距離計算方法。

1 偏航距離計算模型

船舶在地球上航行，偏航距離與地球球體形狀和球面線形有關。以下分別對偏航距離的內(nèi)涵、橢球面算法和后續(xù)使用的優(yōu)化算法進行分析和概述。

1.1 偏航距離內(nèi)涵

綜合IMO、IHO和IALA有關描述，偏航距離是垂直于計劃航線(預期航跡)的偏移量，當船舶與該段計劃航線的間距超出設置的偏移量時，電子海圖系統(tǒng)會發(fā)出警報[1-2，7，15]。在航行方法中，大洋航行可選用大圓航線(great circle route)、恒向線航線(rhumb line route)、等緯圈航線(parallel route)和混合航線(composite route)，然后再使用分點將選用的航線分成恒向線。在航用海圖上，分航段折線是恒向線航線；大地球體通常采用地球橢球體(Earth ellipsoid)。因此,高精度的偏航距離計算是地球橢球體上船位與恒向線間的距離計算。

如圖1所示，航線(恒向線)AiAi+1的延長線與赤道相交于點A0，AiAi+1的航程為S。航行中船舶位于點O2處，與計劃航線的偏航距離為O1O2。

圖1 偏航距離示意圖

1.2 橢球面算法

在偏航距離求解中,采用以下幾種橢球面算法進行長度計算和坐標轉換。

1.2.1 大地線主題正反解

大地線有3種不同定義，雖然這3種大地線微分方程不同，但是由它們導出的大地線方程相同。本文采用的大地線是連接曲面上兩點的最短線。假設已知圖1中點O1(φ1,λ1)，則偏航距離O1O2長度即為大地線長度。目前關于地球圓球體和地球橢球體上的大地線主題正反解(direct and inverse solution of geodesic)問題的求解方法很多，文森特公式精度高、適用廣，最適合于計算機計算[15]。

1.2.2 恒向線主題正反解

針對恒向線主題正反解(direct and inverse solution of rhumb)問題已有高精度的求解算法。史國友等[16]對恒向線主題正反解問題的求解算法進行了優(yōu)化。優(yōu)化后的算法具有不需要迭代和象限判定的優(yōu)點，在船舶航跡計算和電子海圖航線設計中通用性強,適合于計算機計算。本文使用該算法對圖1中航段AiAi+1進行恒向線正反解計算。

1.2.3 墨卡托投影正反解

墨卡托投影正反解(direct and inverse solution of Mercator projection)能夠?qū)崿F(xiàn)大地坐標系與墨卡托坐標系間的坐標轉換。張顯庫等[13]和鄭云峰等[14]在墨卡托投影坐標上以歐氏距離進行了偏航距離的計算。SNYDER[17]使用等量坐標和級數(shù)展開計算了墨卡托投影的正反解，其中等量坐標使用等量緯度正反解實現(xiàn)。本文使用該算法進行墨卡托正反解計算。

1.3 Nelder-Mead simplex算法

Nelder-Mead simplex算法是一種優(yōu)化算法，它通過迭代求解目標函數(shù)的最大值或最小值。這種優(yōu)化方法不同于COTTLE等[18]提出的單純形法，被稱為下山單純形法[19]。該算法的低維適用性和不需要求導的特性符合本文后續(xù)目標函數(shù)的求解要求，因此本文用其進行偏航距離的迭代求解。

2 偏航距離求解的算法設計

在船舶航跡保持控制、靠泊或進出港操縱等應用計劃航線進行導航的情境中，需要求解船位與一系列折線的最近距離，這個最近距離就是偏航距離。本文以用計劃航線進行導航為例進行算法設計。

2.1 直接求解

直接求解法(以下簡稱算法2)源于筆者研究團隊電子海圖模擬器，綜合運用了橢球面上大地線主題正反解和墨卡托投影正反解算法，直接計算偏航距離，計算流程見圖2。值得注意的是，在船舶定位和航線設計中，點坐標(計劃航線轉向點和船位點)使用大地坐標系，因此轉向點坐標用經(jīng)緯度表示。

2.2 迭代求解

圖2 偏航距離直接求解流程

圖3 偏航距離迭代求解示意圖

如圖4所示，迭代求解算法設計思路是：在確定計劃航線(預期航跡)航段后，迭代求解船舶與計劃航線的最近距離。迭代求解算法采用Nelder-Mead simplex算法。目標函數(shù)D連續(xù)且一階微分存在唯一零點，用優(yōu)化算法取得D的極小值。

圖4 偏航距離迭代求解

算法偽代碼如下：

輸入：船位點O2，計劃航線轉向點(A1，A2，…，An)

輸出：偏航距離

中間變量：航向CAiAi+1、航程SAiAi+1

% 確定計劃航線航段，計算航段航向、航程

FOR EACHA∈{A1,A2,…,An-1}

IF球面∠O2AiAi+1<=90°&∠O2Ai+1Ai<90°

[CAiAi+1,SAiAi+1]=InverseRhumb(Ai,Ai+1);

BREAK FOR;

END IF

END FOR

% 航段上任意點與船位的相對距離

相對航段起點Ai真方位C=CAiAi+1

% Nelder-Mead simplex算法

O1= DirectGeodesic(Ai,C,SDmin);%垂足

偏航距離 = InverseGeodesic(O1,O2);

偽代碼中:DirectRhumb、InverseRhumb分別為恒向線正、反解計算，DirectGeodesic、InverseGeodesic分別為大地線正、反解計算，Nelder_Mead_simplex為下山單純形法優(yōu)化計算。

3 算法誤差分析

用算法1、算法2和迭代求解法都能求解偏航距離。用迭代求解法得到的值逼近真值、精度高，將其分別與用算法1和算法2得到的計算結果比較，驗證算法1和算法2的精度。

3.1 誤差的分析過程設計

如圖5所示，為驗證算法精度，分別以緯度φ、船舶相對O點的真方位Ψ和恒向線OO2長度作為自變量，設置取值范圍分別為0°～90°、0°～360°和1～10 n mile(與NMEA0183中RMB語句中偏航距離范圍(0～9.9 n mile)對應)[20]，全局驗證橢球體上的偏航距離。

用算法1、算法2和迭代求解法計算船位點O2與航線OAi+1間的偏航距離。將用算法1和算法2求解的值分別減去用迭代求解法得到的值作為誤差。如果誤差大于0，則用算法1和算法2求解的值大于迭代求解值。

圖5 偏航距離計算誤差的分析過程設計

3.2 算法1誤差分析

將用算法1與迭代求解法計算出的偏航距離的差(即算法1的誤差)作為驗證指標，檢驗算法1的可行性和計算精度，結果見圖6。

a)OO2=1 n mile時的誤差b)真方位為40°時的誤差

圖6算法1的偏航距離計算誤差與恒向線OO2長度、O點緯度和船舶相對O點的真方位的關系

如圖6a所示:隨著船位相對O點的真方位(0°～360°)的變化，誤差幾乎不變；隨著緯度的增加，誤差出現(xiàn)較大變化，這種特征與墨卡托投影特性一致。

如圖6b所示：當恒向線OO2長度為1 n mile、緯度為31°時，誤差接近0.179 n mile(331 m)；當恒向線OO2長度為5 n mile、緯度為60°時，誤差接近5 n mile。在船舶避碰或大洋航行的應用中，誤差的影響較為明顯。出現(xiàn)誤差的原因是，墨卡托投影保證了等角和航海視覺習慣[21]，但無法保證等積、等距離。

3.3 算法2誤差分析

墨卡托投影導致子午線(弧長)發(fā)生形變[16]，算法1誤差較大。算法2中使用了墨卡托投影，因此有必要驗證算法2在全局緯度上的計算精度。將用算法2與用迭代求解法計算出的偏航距離的差(即算法2的誤差)作為驗證指標，檢驗算法2的可行性和計算精度，結果見圖7。

a)OO2=1 n mile時的誤差b)O點緯度為89°時誤差與真方位的關系

c)真方位為0～360°時誤差最大值與緯度的關系d)O點緯度為89°時直接與迭代求解的垂足間距離與真方位的關系

圖7算法2的偏航距離計算誤差與恒向線OO2長度、O點緯度和船舶相對O點的真方位的關系

圖7展示了算法2的偏航距離計算精度。由圖7a可知，O點緯度在0～80°范圍內(nèi)時誤差較小，在80～90°范圍內(nèi)時誤差較大；船舶相對O點的真方位在0～360°范圍內(nèi)時誤差呈現(xiàn)雙峰形狀。由圖7a～7d可知，OO2=4，7，10 n mile時的誤差大于OO2=1 n mile時的誤差。

圖7b解釋了船舶相對O點的真方位對誤差的影響。誤差最大值出現(xiàn)在船舶相對O點的真方位為90°或270°左右時，即船舶航向為0°或180°時，與迭代求解值呈非線性增函數(shù)關系。

圖7c展示了高緯度區(qū)域算法2誤差的數(shù)量級與緯度的關系。當O點緯度在89.6°～89.7°范圍內(nèi)且OO2=10 n mile時，誤差為400～900 m；在高緯度區(qū)域，誤差非線性增長。

圖7d展示了當O點緯度為89°時，用算法2和迭代求解法求出的垂足不在同一位置，即算法2求解的空間物標跟蹤定向航行船舶的最近距離點精度不準。

綜上可知：低緯度區(qū)域(0～80°)算法2計算精度較高；當船舶在高緯度(80°～90°)區(qū)域航行時，算法2近距離計算精度較高，遠距離計算誤差較大。

4 改進的直接求解算法及其誤差分析

由于墨卡托投影的性質(zhì)，高緯度(80°～90°)區(qū)域算法2的誤差較大?，F(xiàn)有投影變換都不能兼顧等角、等積和等距離條件，各種投影具有各自的特點，在極端情況下都會出現(xiàn)投影變形。

針對高緯度區(qū)域，SKOPELITI等[21]推薦電子海圖系統(tǒng)使用的等距方位投影和正軸正射投影，以保持形狀不變。此外隨著極地規(guī)則[22]的生效和IHO下屬ACHR(Arctic Regional Hydrographic Commission)對極地區(qū)域的持續(xù)關注[23]，高緯度區(qū)域航行越來越受到重視。

本文選擇極點為投影中心的等距方位投影作為直接求解的投影變換，即圖2第2步使用正軸等距方位投影(azimuthal polar equidistant projection)直接求解(以下簡稱算法3)。這種投影的特點是：極點處無變形；距離極點越遠，面積變形越大，角度變形越大。[21]將算法3用于圖5偏航距離計算誤差的分析過程設計中，并進行與用迭代求解法得出的結果之差(即算法3的誤差)的計算。如圖8所示：當緯度在80°～90°范圍內(nèi)，恒向線OO2長度分別為1、4、7、10 n mile時，對應誤差均值分別為0.000 9、0.003 3、0.002 7、0.011 8 m，誤差很??；在O點緯度為0°時誤差達到最大；誤差峰均有4個。算法3有效提高了高緯度區(qū)域偏航距離的求解精度。

a)OO2=1 n mile時誤差b)OO2=4 n mile時誤差

c)OO2=7 n mile時誤差d)OO2=10 n mile時誤差

圖8算法3偏航距離的計算誤差與緯度和真方位的關系

5 實 例

5.1 中緯度區(qū)域

在水上水下通航安全的工程中，存在計算鉆井平臺與航道距離的案例。為評估鉆井平臺對航行的影響，需要計算鉆井平臺與航道(延長線)的距離。如圖9a所示，鉆井平臺處于黃驊港進出港航道邊緣。計算垂足O1和偏航距離并與電子海圖測量值比對，見圖9b。

a)鉆井平臺與航道的位置關系

b)迭代次數(shù)與偏航距離的關系

電子海圖活動距標圈、算法2、迭代求解法和算法1得出的偏航距離見表1。由于電子海圖活動距標圈是手動測量的，所以誤差較大。算法2與迭代求解法得出的結果相近。算法1的偏航距離的計算誤差較大。算法1是用墨卡托坐標計算的，計算的垂足用墨卡托坐標表示，因此沒有垂足的經(jīng)緯度。

表1 中緯度區(qū)域鉆井平臺與航道(延長線)間偏航距離計算結果比對

5.2 高緯度區(qū)域

假設高緯度區(qū)域存在轉向點A1(89°36′N，0°E)、轉向點A2(89°42′N，0°E)、船位點O2(89°36′N，24°30′E)。與中緯度區(qū)域?qū)嵗愃?，高緯度區(qū)域船位點O2與恒向線A1A2間的偏航距離計算結果見表2。由于電子海圖活動距標圈是手動測量的，所以誤差較大。算法3與迭代求解法得出的結果相近，并且二者與算法2得出的結果相差約432 m。與中緯度區(qū)域應用一樣，算法1的偏航距離計算誤差較大，并且是墨卡托坐標上的值，因此沒有垂足O1的經(jīng)緯度。

表2 高緯度區(qū)域船位點O2與恒向線A1A2間的偏航距離計算結果

6 結束語

適當?shù)剡x擇直接求解中的投影變換，可以達到減小偏航距離計算誤差、精確求解的目的。雖然墨卡托投影解決了大地坐標的球面與航用海圖的平面之間的矛盾，但是帶來了投影變形，從而導致高緯度區(qū)域偏航距離的直接求解存在誤差。

本文通過將迭代求解法分別與算法1(基于墨卡托投影的歐氏距離算法)和算法2(基于墨卡托投影的直接求解算法)比較，發(fā)現(xiàn)后2種方法均存在不同程度的計算誤差。為提高求解精度，提出算法3(基于正軸等距方位投影的直接求解算法)。

結論如下：用迭代求解法求得的偏航距離逼近真值，求解誤差最?。凰惴?的適用范圍較小，當緯度在80°～90°范圍內(nèi)時用算法2求解的偏航距離存在誤差；當緯度在80°～90°范圍內(nèi)時用算法3求解的偏航距離可提高2%的計算精度。

建議將算法2和算法3分別用于中低緯度和高緯度區(qū)域，以實現(xiàn)全局緯度上偏航距離的求解。

本文研究為北極航區(qū)航行、進出港操縱等情況下的船舶導航、運動控制、航路規(guī)劃提供參考。