袁 兵,葉世榕,鮑立桂
(1.武漢大學(xué)衛(wèi)星導(dǎo)航定位技術(shù)研究中心,武漢 430072;2.福建新大陸電腦有限公司,福州 350015)
傳統(tǒng)單點(diǎn)定位只需要一臺接收機(jī)就可以進(jìn)行實(shí)時(shí)導(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)定位計(jì)算的復(fù)雜度和精度都有明顯改善;范磊等[4]分析了COMPASS和GPS兩個(gè)系統(tǒng)偽距單點(diǎn)定位之間的差異;陸亞峰等[5]分析得出北斗/GPS組合偽距單點(diǎn)定位精度明顯優(yōu)于單系統(tǒng);李鶴峰等[6-7]給出建立在時(shí)空統(tǒng)一上的 BDS/GPS/GLONASS多模融合定位模型;文獻(xiàn) [8-9]給出了GPS、BDS、GLONASS組合單點(diǎn)定位中各系統(tǒng)偽距之間的定權(quán)方法;張?jiān)鲁龋?]將Kalman濾波應(yīng)用于GNSS偽距單點(diǎn)定位中,得到了精度更高,穩(wěn)定性更好的結(jié)果。
本文將Kalman濾波應(yīng)用于GNSS偽距單點(diǎn)定位中,提出一種融合多模導(dǎo)航定位數(shù)據(jù)處理算法,在一定條件下將多個(gè)系統(tǒng)接收機(jī)鐘差融合為一個(gè)系統(tǒng)接收機(jī)鐘差,減少了多模導(dǎo)航定位估計(jì)參數(shù),從而降低多模導(dǎo)航定位的最低條件,充分發(fā)揮多模GNSS較單系統(tǒng)具有的優(yōu)勢。
以GPS、COMPASS和GLONASS三系統(tǒng)為例,在傳統(tǒng)的誤差方程中,每個(gè)系統(tǒng)的接收機(jī)鐘差都作為一個(gè)獨(dú)立的參數(shù)來估計(jì),除去三維坐標(biāo)參數(shù)外,三系統(tǒng)的導(dǎo)航定位數(shù)據(jù)處理還需估計(jì)3個(gè)接收機(jī)鐘差參數(shù),總共6個(gè)估計(jì)參數(shù),誤差方程如式(1)所示:
事實(shí)上,對于同一接收機(jī)采集的多系統(tǒng)數(shù)據(jù),各系統(tǒng)共同使用同一接收機(jī)鐘,則它們接收機(jī)鐘差應(yīng)該是一樣的。但是,不同系統(tǒng)的授時(shí)系統(tǒng)不一,導(dǎo)致各系統(tǒng)時(shí)間系統(tǒng)間存在一定偏差,對固定兩個(gè)系統(tǒng)而言,兩個(gè)系統(tǒng)時(shí)間基準(zhǔn)間偏差為一波動量較小的常數(shù)[10]。這些偏差表現(xiàn)為接收機(jī)鐘差系統(tǒng)間偏差,如式 (2)中的m、n所示。假設(shè)式 (2)中 m、n為常數(shù) (在第2節(jié)用實(shí)驗(yàn)驗(yàn)證m、n為波動量較小的數(shù)值,可以當(dāng)做常數(shù)處理),得到:
由此可得,若得知m、n的值之后,便可以將三個(gè)系統(tǒng)的接收機(jī)鐘差參數(shù)融合為1個(gè)參數(shù),在多模導(dǎo)航定位中只估計(jì)4個(gè)參數(shù)即可。
Kalman濾波是一種線性最小方差估計(jì)方法,其計(jì)算過程是不斷預(yù)測和修正的過程,無需存儲大量不同時(shí)刻的數(shù)據(jù)便可以解決大部分問題。
在初始?xì)v元中,并不能確切知道m(xù)、n值。因此在衛(wèi)星數(shù)充足的情況下,將m、n同三維坐標(biāo)以及GPS接收機(jī)鐘差一起作為估計(jì)參數(shù),狀態(tài)方程如下
式中,Xk,k-1為 k 時(shí)刻的狀態(tài)預(yù)測值,Xk-1為k-1時(shí)刻的狀態(tài)估計(jì)值,Φk,k-1為狀態(tài)轉(zhuǎn)移矩陣。
k=0時(shí),狀態(tài)方程初始值X0為第一個(gè)歷元用最小二乘算法迭代計(jì)算得出的偽距單點(diǎn)定位結(jié)果。
狀態(tài)轉(zhuǎn)移矩陣可表示為
式中,Ip、IΔt分別為 3×3、2×2 單位矩陣。
第k歷元狀態(tài)估計(jì)方程為
式中,Zk為第k歷元觀測序列,為第k歷元各衛(wèi)星偽距觀測量加上偽距觀測量初值及相關(guān)誤差改正 (對流層、電離層、地球自轉(zhuǎn)改正等)之后的結(jié)果;Hk為第k歷元設(shè)計(jì)矩陣,由測站至衛(wèi)星在三個(gè)方向的方向余弦組成;Kk濾波增益矩陣。
濾波增益矩陣Kk解算公式為
式中,Rk為系統(tǒng)觀測噪聲方差矩陣,它是根據(jù)第k歷元各衛(wèi)星高度角計(jì)算得到,高度角越大,Rk矩陣中相應(yīng)的值越小;Pk,k-1為k時(shí)刻誤差方程矩陣預(yù)測值。
預(yù)測誤差方程矩陣Pk,k-1計(jì)算公式為
一般情況下不考慮過程噪聲 Qk-1,因此式(8)可以簡化為
k=0時(shí),P0為單位矩陣。
估計(jì)誤差方程矩陣計(jì)算公式為
由以上方程可知,只給定初值X0和P0,就可以根據(jù)k時(shí)刻觀測值Zk遞推求得k時(shí)刻狀態(tài)估計(jì)值Xk。
待經(jīng)過多次歷元解算得到一個(gè)較為穩(wěn)定的m、n值之后,可以將m、n作為常數(shù)參與到往后的歷元解算中,此時(shí)只需估計(jì)4個(gè)參數(shù) (3個(gè)坐標(biāo)參數(shù),1個(gè)接收機(jī)鐘差參數(shù))。狀態(tài)方程為
相應(yīng)狀態(tài)轉(zhuǎn)移矩陣變化為
式中,Ip為3×3單位矩陣。
本文實(shí)驗(yàn)數(shù)據(jù)為某個(gè)IGS站一天24小時(shí)的觀測文件 (包含GPS、北斗和GLONASS觀測數(shù)據(jù)),觀測時(shí)間為2014年1月16日,設(shè)置30s采樣率及15°截止高度角,總共2880個(gè)觀測歷元數(shù)據(jù)。
為了驗(yàn)證對于同一臺接收機(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濾波估計(jì)。m、n值估計(jì)結(jié)果如圖1、圖2所示。為了直觀,接收機(jī)鐘差以及m、n值均用米作為單位(即它們的值為鐘差乘以光速c之后的結(jié)果)。圖1和圖2分別為Kalman濾波方程下各歷元得到的m值和n值,其橫坐標(biāo)單位為個(gè),每一個(gè)單位間隔代表一個(gè)歷元。
圖1 Kalman濾波方程下單天時(shí)段內(nèi)各歷元m值Fig.1 Value of m of one day by Kalman filter algorithm
圖2 Kalman濾波方程下單天時(shí)段內(nèi)方向下各歷元n值Fig.2 Value of n of one day by Kalman filter algorithm
表1 Kalman濾波方程解出的單天時(shí)段內(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,說明單天時(shí)段內(nèi)m、n值波動范圍不大。同時(shí)從表1可以看出,單天時(shí)段內(nèi)m、n值最大值與最小值之差均小于2m,這個(gè)誤差對偽距單點(diǎn)定位10m精度要求來說可以接受。但為了求得更精確的定位結(jié)果,選第50個(gè)歷元的Kalman濾波值作為m、n值的已知值參與到后續(xù)Kalman濾波迭代運(yùn)算中。如圖3、圖4和圖5中紅色線和青色線所示,在第50個(gè)歷元固定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é)果的影響,以及驗(yàn)證在三系統(tǒng)衛(wèi)星總數(shù)只有 (3+1)或者 (3+2)情況下固定m、n值之后依然可以進(jìn)行多模單點(diǎn)定位,本文做了三個(gè)實(shí)驗(yàn),如圖3、圖4和圖5所示,縱坐標(biāo)單位為米,橫坐標(biāo)一個(gè)單位代表一個(gè)歷元,一個(gè)歷元30s,時(shí)間跨度從2014年1月16日0時(shí)0分0秒開始到2014年1月16日23時(shí)59分30秒結(jié)束。圖3、圖4和圖5中紅色線條所示為在衛(wèi)星充足情況下不固定m、n,每個(gè)系統(tǒng)的接收機(jī)鐘差均被獨(dú)立估計(jì)情況下得出的測站三維坐標(biāo)與測站坐標(biāo)真值之差;青色線條為衛(wèi)星充足情況下在第50個(gè)歷元之后,固定m、n得到的測站三維坐標(biāo)與測站坐標(biāo)真值之差;藍(lán)色線條為在第50個(gè)歷元之后,人為地將每個(gè)歷元的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é)果較平滑,最后逐漸收斂到一個(gè)較接近真值的值域內(nèi)。
圖5 單天Kalman濾波解各情況下Z坐標(biāo)值比較Fig.5 Comparation of Z’s coordinate values in the three cases
分析圖3、圖4和圖5中的紅色線條和青色線條,即衛(wèi)星數(shù)充足情況下不固定m、n(紅色線條所示)和固定m、n(青色線條所示)兩種模式得出的測站坐標(biāo)結(jié)果在X、Y、Z方向上基本一致,由此可知將在第50個(gè)歷元計(jì)算得出的m、n值從第51個(gè)歷元開始,作為m、n的已知值參與后續(xù)Kalman濾波迭代運(yùn)算之后,m、n值在各歷元的波動量并沒有引起多模偽距單點(diǎn)定位太大的變化,在固定m、n值和不固定m、n值情況下偽距單點(diǎn)定位結(jié)果幾乎一致。因此,在單天時(shí)段內(nèi),取一定數(shù)量歷元計(jì)算所得的m、n平均值作為m、n已知值代入式 (3),代替GLONASS系統(tǒng)和北斗系統(tǒng)接收機(jī)鐘差,參與到后續(xù)Kalman濾波計(jì)算中是可行的。
圖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)致的。
本文提出了一種在多模偽距單點(diǎn)定位中將多系統(tǒng)接收機(jī)鐘差融合為一個(gè)系統(tǒng)接收機(jī)鐘差的算法,并推導(dǎo)了該算法的遞推公式。該算法建立了Kalman濾波方程,將GLONASS系統(tǒng)和北斗系統(tǒng)接收機(jī)鐘差與GPS系統(tǒng)接收機(jī)鐘差之差 (如式 (2)所示)同位置參數(shù)和GPS系統(tǒng)接收機(jī)鐘差一起作為估計(jì)參數(shù)參與到Kalman濾波迭代運(yùn)算中,并選取經(jīng)過一定歷元計(jì)算得到的Kalman濾波值作為m、n已知值參與到后續(xù)Kalman濾波迭代運(yùn)算中,從而將三系統(tǒng)狀態(tài)方程估計(jì)參數(shù)從6個(gè)減少到4個(gè)。該算法降低了多模定位的最低條件,充分發(fā)揮了多模GNSS相對于單系統(tǒng)的可見衛(wèi)星總數(shù)較多等優(yōu)勢,即使在可見衛(wèi)星總數(shù)只有 (3+1)或者(3+2)顆情況下依然可以進(jìn)行三系統(tǒng)的偽距單點(diǎn)定位。實(shí)驗(yàn)結(jié)論表明,融合多模導(dǎo)航定位數(shù)據(jù)處理算法可以達(dá)到10m的偽距單點(diǎn)定位精度,符合標(biāo)準(zhǔn)單點(diǎn)定位的精度要求。本文做了單天時(shí)段內(nèi)三系統(tǒng)導(dǎo)航定位的實(shí)驗(yàn),在以后的工作中可以將該算法推廣到雙系統(tǒng)、四系統(tǒng)或者是更多系統(tǒng)的導(dǎo)航定位數(shù)據(jù)處理中,其理論和算法同三系統(tǒng)是一樣的。同時(shí)試驗(yàn)時(shí)間跨度可以從單天時(shí)段延伸到數(shù)日或數(shù)月情況下。
[1]張?jiān)鲁惲x,胡川.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é)報(bào) (自然科學(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].測繪通報(bào),2014(9):1-5.
[7]李鶴峰,黨亞民,秘金鐘,等.BDS與GPS、GLONASS多模融合導(dǎo)航定位時(shí)空統(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é)報(bào),2014,12(3):374-378.
[10]王黨衛(wèi).時(shí)間統(tǒng)一系統(tǒng)研究[J].現(xiàn)代導(dǎo)航,2012,6:450-455.