袁 兵，葉世榕，鮑立桂

（1.武漢大學衛(wèi)星導航定位技術研究中心，武漢 430072;2.福建新大陸電腦有限公司，福州 350015）

0 引言

傳統(tǒng)單點定位只需要一臺接收機就可以進行實時導航定位，條件要求少，因此廣泛應用于船舶導航、車輛定位、地質勘探等對精度要求不高的領域［1］。并且，多衛(wèi)星導航系統(tǒng)的發(fā)展使得單點定位有可能為人們提供更好的服務。GNSS多系統(tǒng)的集成無論是從連續(xù)性、可用性、可靠性、精度以及效率等各方面都更具優(yōu)勢［2］。國內外學者對偽距單點定位和多系統(tǒng)導航定位都有較多深入的研究:李春華等［3］提出一種北斗偽距單點定位加權最小二乘 （WLS）快速算法，使得偽距單點定位計算的復雜度和精度都有明顯改善;范磊等［4］分析了COMPASS和GPS兩個系統(tǒng)偽距單點定位之間的差異;陸亞峰等［5］分析得出北斗/GPS組合偽距單點定位精度明顯優(yōu)于單系統(tǒng);李鶴峰等［6-7］給出建立在時空統(tǒng)一上的 BDS/GPS/GLONASS多模融合定位模型;文獻 ［8-9］給出了GPS、BDS、GLONASS組合單點定位中各系統(tǒng)偽距之間的定權方法;張月超等［1］將Kalman濾波應用于GNSS偽距單點定位中，得到了精度更高，穩(wěn)定性更好的結果。

本文將Kalman濾波應用于GNSS偽距單點定位中，提出一種融合多模導航定位數(shù)據(jù)處理算法，在一定條件下將多個系統(tǒng)接收機鐘差融合為一個系統(tǒng)接收機鐘差，減少了多模導航定位估計參數(shù)，從而降低多模導航定位的最低條件，充分發(fā)揮多模GNSS較單系統(tǒng)具有的優(yōu)勢。

1 多模GNSS偽距單點定位方程

以GPS、COMPASS和GLONASS三系統(tǒng)為例，在傳統(tǒng)的誤差方程中，每個系統(tǒng)的接收機鐘差都作為一個獨立的參數(shù)來估計，除去三維坐標參數(shù)外，三系統(tǒng)的導航定位數(shù)據(jù)處理還需估計3個接收機鐘差參數(shù)，總共6個估計參數(shù)，誤差方程如式（1）所示:

事實上，對于同一接收機采集的多系統(tǒng)數(shù)據(jù)，各系統(tǒng)共同使用同一接收機鐘，則它們接收機鐘差應該是一樣的。但是，不同系統(tǒng)的授時系統(tǒng)不一，導致各系統(tǒng)時間系統(tǒng)間存在一定偏差，對固定兩個系統(tǒng)而言，兩個系統(tǒng)時間基準間偏差為一波動量較小的常數(shù)［10］。這些偏差表現(xiàn)為接收機鐘差系統(tǒng)間偏差，如式 （2）中的m、n所示。假設式 （2）中 m、n為常數(shù) （在第2節(jié)用實驗驗證m、n為波動量較小的數(shù)值，可以當做常數(shù)處理），得到:

由此可得，若得知m、n的值之后，便可以將三個系統(tǒng)的接收機鐘差參數(shù)融合為1個參數(shù)，在多模導航定位中只估計4個參數(shù)即可。

2 Kalman濾波

Kalman濾波是一種線性最小方差估計方法，其計算過程是不斷預測和修正的過程，無需存儲大量不同時刻的數(shù)據(jù)便可以解決大部分問題。

在初始歷元中，并不能確切知道m(xù)、n值。因此在衛(wèi)星數(shù)充足的情況下，將m、n同三維坐標以及GPS接收機鐘差一起作為估計參數(shù)，狀態(tài)方程如下

式中，Xk，k－1為 k 時刻的狀態(tài)預測值，Xk－1為k－1時刻的狀態(tài)估計值，Φk，k－1為狀態(tài)轉移矩陣。

k=0時，狀態(tài)方程初始值X0為第一個歷元用最小二乘算法迭代計算得出的偽距單點定位結果。

狀態(tài)轉移矩陣可表示為

式中，Ip、IΔt分別為 3×3、2×2 單位矩陣。

第k歷元狀態(tài)估計方程為

式中，Zk為第k歷元觀測序列，為第k歷元各衛(wèi)星偽距觀測量加上偽距觀測量初值及相關誤差改正 （對流層、電離層、地球自轉改正等）之后的結果;Hk為第k歷元設計矩陣，由測站至衛(wèi)星在三個方向的方向余弦組成;Kk濾波增益矩陣。

濾波增益矩陣Kk解算公式為

式中，Rk為系統(tǒng)觀測噪聲方差矩陣，它是根據(jù)第k歷元各衛(wèi)星高度角計算得到，高度角越大，Rk矩陣中相應的值越小;Pk，k－1為k時刻誤差方程矩陣預測值。

預測誤差方程矩陣Pk，k－1計算公式為

一般情況下不考慮過程噪聲 Qk－1，因此式（8）可以簡化為

k=0時，P0為單位矩陣。

估計誤差方程矩陣計算公式為

由以上方程可知，只給定初值X0和P0，就可以根據(jù)k時刻觀測值Zk遞推求得k時刻狀態(tài)估計值Xk。

待經(jīng)過多次歷元解算得到一個較為穩(wěn)定的m、n值之后，可以將m、n作為常數(shù)參與到往后的歷元解算中，此時只需估計4個參數(shù) （3個坐標參數(shù)，1個接收機鐘差參數(shù)）。狀態(tài)方程為

相應狀態(tài)轉移矩陣變化為

式中，Ip為3×3單位矩陣。

3 實驗及結果分析

本文實驗數(shù)據(jù)為某個IGS站一天24小時的觀測文件 （包含GPS、北斗和GLONASS觀測數(shù)據(jù)），觀測時間為2014年1月16日，設置30s采樣率及15°截止高度角，總共2880個觀測歷元數(shù)據(jù)。

3.1 實驗

為了驗證對于同一臺接收機接收到的多系統(tǒng)數(shù)據(jù)而言，兩兩系統(tǒng)間接收機鐘差偏差值m、n為一在單天內波動較小的數(shù)值，將m、n值連同三維坐標以及GPS系統(tǒng)接收機鐘差作為式 （4）狀態(tài)方程參數(shù)，一起參與式 （4）～式 （10）的Kalman濾波估計。m、n值估計結果如圖1、圖2所示。為了直觀，接收機鐘差以及m、n值均用米作為單位（即它們的值為鐘差乘以光速c之后的結果）。圖1和圖2分別為Kalman濾波方程下各歷元得到的m值和n值，其橫坐標單位為個，每一個單位間隔代表一個歷元。

圖1 Kalman濾波方程下單天時段內各歷元m值Fig.1 Value of m of one day by Kalman filter algorithm

圖2 Kalman濾波方程下單天時段內方向下各歷元n值Fig.2 Value of n of one day by Kalman filter algorithm

表1 Kalman濾波方程解出的單天時段內m、n值均值與均方差列表Tab.1 The list of mean and variance of m and n of one day computed by Kalman filter algorithm

分析圖1、圖2及表1可知，m、n值在初始一段歷元中在一定較小范圍內波動，后逐漸趨于平緩，它們的均方差均小于0.4m，說明單天時段內m、n值波動范圍不大。同時從表1可以看出，單天時段內m、n值最大值與最小值之差均小于2m，這個誤差對偽距單點定位10m精度要求來說可以接受。但為了求得更精確的定位結果，選第50個歷元的Kalman濾波值作為m、n值的已知值參與到后續(xù)Kalman濾波迭代運算中。如圖3、圖4和圖5中紅色線和青色線所示，在第50個歷元固定m、n值之后，多模偽距單點定位結果并沒有受到太大影響。

圖3 單天Kalman濾波解各情況下X坐標值比較Fig.3 Comparation of X’s coordinate values in the three cases

圖4 單天Kalman濾波解各情況下Y坐標值比較Fig.4 Comparation of Y’s coordinate values in the three cases

為了對比在衛(wèi)星數(shù)充足情況下固定和不固定m、n值對定位結果的影響，以及驗證在三系統(tǒng)衛(wèi)星總數(shù)只有 （3+1）或者 （3+2）情況下固定m、n值之后依然可以進行多模單點定位，本文做了三個實驗，如圖3、圖4和圖5所示，縱坐標單位為米，橫坐標一個單位代表一個歷元，一個歷元30s，時間跨度從2014年1月16日0時0分0秒開始到2014年1月16日23時59分30秒結束。圖3、圖4和圖5中紅色線條所示為在衛(wèi)星充足情況下不固定m、n，每個系統(tǒng)的接收機鐘差均被獨立估計情況下得出的測站三維坐標與測站坐標真值之差;青色線條為衛(wèi)星充足情況下在第50個歷元之后，固定m、n得到的測站三維坐標與測站坐標真值之差;藍色線條為在第50個歷元之后，人為地將每個歷元的GLONASS可見衛(wèi)星減少為1顆，北斗可見衛(wèi)星減少為1～2顆，GPS可見衛(wèi)星減少為1～2顆，保證三系統(tǒng)可見衛(wèi)星總數(shù)為 （3+1）顆或者 （3+2）顆，然后開始固定m、n值而得到的測站三維坐標值與測站坐標真值之差。

由于采用了Kalman濾波，測站坐標沒有出現(xiàn)過大的跳變，結果較平滑，最后逐漸收斂到一個較接近真值的值域內。

圖5 單天Kalman濾波解各情況下Z坐標值比較Fig.5 Comparation of Z’s coordinate values in the three cases

3.2 實驗分析

分析圖3、圖4和圖5中的紅色線條和青色線條，即衛(wèi)星數(shù)充足情況下不固定m、n（紅色線條所示）和固定m、n（青色線條所示）兩種模式得出的測站坐標結果在X、Y、Z方向上基本一致，由此可知將在第50個歷元計算得出的m、n值從第51個歷元開始，作為m、n的已知值參與后續(xù)Kalman濾波迭代運算之后，m、n值在各歷元的波動量并沒有引起多模偽距單點定位太大的變化，在固定m、n值和不固定m、n值情況下偽距單點定位結果幾乎一致。因此，在單天時段內，取一定數(shù)量歷元計算所得的m、n平均值作為m、n已知值代入式 （3），代替GLONASS系統(tǒng)和北斗系統(tǒng)接收機鐘差，參與到后續(xù)Kalman濾波計算中是可行的。

圖3、圖4和圖5的藍色線條表示了在三系統(tǒng)可見衛(wèi)星總數(shù)只有4顆或者5顆情況下，固定m、n所得測站坐標值與真值之差，其精度要低于在衛(wèi)星數(shù)充足情況下固定m、n所得測站坐標的精度，但是仍在偽距單點定位所能接受的誤差范圍之內。其精度較差原因可能是由于在極少數(shù)可見衛(wèi)星數(shù)情況下，衛(wèi)星空間結構較差，以及多余觀測值較少導致的。

4 結論

本文提出了一種在多模偽距單點定位中將多系統(tǒng)接收機鐘差融合為一個系統(tǒng)接收機鐘差的算法，并推導了該算法的遞推公式。該算法建立了Kalman濾波方程，將GLONASS系統(tǒng)和北斗系統(tǒng)接收機鐘差與GPS系統(tǒng)接收機鐘差之差 （如式 （2）所示）同位置參數(shù)和GPS系統(tǒng)接收機鐘差一起作為估計參數(shù)參與到Kalman濾波迭代運算中，并選取經(jīng)過一定歷元計算得到的Kalman濾波值作為m、n已知值參與到后續(xù)Kalman濾波迭代運算中，從而將三系統(tǒng)狀態(tài)方程估計參數(shù)從6個減少到4個。該算法降低了多模定位的最低條件，充分發(fā)揮了多模GNSS相對于單系統(tǒng)的可見衛(wèi)星總數(shù)較多等優(yōu)勢，即使在可見衛(wèi)星總數(shù)只有 （3+1）或者（3+2）顆情況下依然可以進行三系統(tǒng)的偽距單點定位。實驗結論表明，融合多模導航定位數(shù)據(jù)處理算法可以達到10m的偽距單點定位精度，符合標準單點定位的精度要求。本文做了單天時段內三系統(tǒng)導航定位的實驗，在以后的工作中可以將該算法推廣到雙系統(tǒng)、四系統(tǒng)或者是更多系統(tǒng)的導航定位數(shù)據(jù)處理中，其理論和算法同三系統(tǒng)是一樣的。同時試驗時間跨度可以從單天時段延伸到數(shù)日或數(shù)月情況下。

［1］張月超，陳義，胡川．Klman濾波在GNSS偽距單點定位中的應用［J］．全球定位系統(tǒng)，2013，38（6）:31-57．

［2］郭斐，張小紅，王明華．GNSS多系統(tǒng)集成的兼容性問題［J］．測繪信息與工程，2012，40（3）:13-15．

［3］李春華，蔡成林，鄧克群．一種北斗偽距單點定位的加權最小二乘 （WLS）快速算法［J］．重慶郵電大學學報 （自然科學版），2014，26（4）:466-472．

［4］范磊，鐘世明，歐吉坤．COMPASS與GPS偽距單點定位精度分析［A］．第四屆中國衛(wèi)星導航學術年會，中國湖北武漢，2013．

［5］陸亞峰，樓立志，馬緒瀛，等．北斗與GPS組合偽距單點定位精度分析［J］．全球定位系統(tǒng)，2013，38（6）:1-6．

［6］李鶴峰，黨亞民，秘金鐘，等．BDS/GPS/GLONASS融合定位模型及性能分析［J］．測繪通報，2014（9）:1-5．

［7］李鶴峰，黨亞民，秘金鐘，等．BDS與GPS、GLONASS多模融合導航定位時空統(tǒng)一［J］．大地測量與地球動力學，2013，33（4）:73-78．

［8］何俊，袁小玲，曾琪．GPS/BDS/GLONASS組合單點定位研究［J］．測繪科學，2014，39（8）:124-128．

［9］尚夢云，高暉，常青，等．GPS/GLONASS/BDS三系統(tǒng)組合定位的定權方法［J］．太赫茲科學與電子信息學報，2014，12（3）:374-378．

［10］王黨衛(wèi)．時間統(tǒng)一系統(tǒng)研究［J］．現(xiàn)代導航，2012，6:450-455．