袁 兵，葉世榕，鮑立桂

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

0 引言

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

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

1 多模GNSS偽距單點(diǎn)定位方程

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

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

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

2 Kalman濾波

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

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

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

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

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

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

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

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

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

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

預(yù)測誤差方程矩陣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個坐標(biāo)參數(shù)，1個接收機(jī)鐘差參數(shù)）。狀態(tài)方程為

相應(yīng)狀態(tài)轉(zhuǎn)移矩陣變化為

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

3 實驗及結(jié)果分析

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

3.1 實驗

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

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

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

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

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

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

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

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

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

3.2 實驗分析

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

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

4 結(jié)論

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

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

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

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

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

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

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

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

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

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

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