劉瑞華,丁其金
(中國(guó)民航大學(xué) 電子信息與自動(dòng)化學(xué)院,天津300300)
北斗衛(wèi)星導(dǎo)航系統(tǒng)(BDS)是我國(guó)自主建設(shè)、獨(dú)立運(yùn)行的衛(wèi)星導(dǎo)航系統(tǒng),正在逐步完善,2020年前后,將向全球提供服務(wù)[1]。其應(yīng)用范圍也越來越廣,在農(nóng)業(yè)、交通及航空等各個(gè)行業(yè)中發(fā)揮著越來越重要的作用。因此,對(duì)北斗接收機(jī)的需求也越來越大,各導(dǎo)航設(shè)備研制單位競(jìng)相展開了北斗接收機(jī)的研制工作。雖然不同的應(yīng)用環(huán)境,對(duì)北斗接收機(jī)要求的側(cè)重點(diǎn)不同,但接收機(jī)位置解算都是必須考慮的部分。本文對(duì)北斗接收機(jī)的位置解算過程進(jìn)行分析研究,并結(jié)合MATLAB軟件,分別使用北斗空間控制接口文件所提供的算法解算北斗衛(wèi)星位置和最小二乘法解算北斗接收機(jī)位置,并進(jìn)行精度分析,驗(yàn)證了算法的正確性和可行性。
北斗衛(wèi)星位置解算是接收機(jī)位置解算中至關(guān)重要的一部分,因?yàn)樾l(wèi)星位置是否準(zhǔn)確直接影響了接收機(jī)的定位精度。北斗衛(wèi)星導(dǎo)航系統(tǒng)采用三種軌道衛(wèi)星,即地球同步軌道 (GEO)、傾斜地球同步軌道(Inclined IGSO)和中地球軌道(MEO)組成的混合星座[1],因此,北斗衛(wèi)星位置的解算方法與GPS略有不同。其中IGSO、MEO衛(wèi)星計(jì)算方法相同,GEO衛(wèi)星因其軌道傾角和偏心率均很小,需采用坐標(biāo)系旋轉(zhuǎn)的方法增大軌道傾角以提高衛(wèi)星位置的計(jì)算精度[2-3]。北斗官方文件[4]給出了北斗衛(wèi)星位置解算的星歷參數(shù)用戶算法,這里不再贅述,圖1為算法流程框圖,表1示出了星歷參數(shù)及其定義。
表1 星歷參數(shù)及其定義
本節(jié)使用諾瓦泰FlexPark6北斗接收機(jī)2018年3月2日2點(diǎn)至23點(diǎn)采集到的實(shí)測(cè)星歷數(shù)據(jù),按照文獻(xiàn)[4]中的算法進(jìn)行北斗衛(wèi)星位置解算,并將解算出的衛(wèi)星位置分別和精密星歷以及諾瓦泰 FlexPark6北斗接收機(jī)所提供的衛(wèi)星位置進(jìn)行對(duì)比,分析所解算出的衛(wèi)星位置精度。本文使用武漢大學(xué)提供的代號(hào)為WUM的精密星歷,徑向誤差優(yōu)于10 cm,比廣播星歷精度高出1-2個(gè)數(shù)量級(jí)[5];根據(jù)諾瓦泰給出的FlexPark6技術(shù)指標(biāo),該接收機(jī)水平定位精度達(dá)到1.5 m(RMS),時(shí)間精度為20 ns,因此均可以作為真值來評(píng)估所解算衛(wèi)星位置的精度。
按照文獻(xiàn)[4]中的算法,進(jìn)行北斗衛(wèi)星位置解算,在解算過程中需注意以下幾個(gè)問題:
1) 在計(jì)算時(shí)間差tk時(shí),必須考慮周變換的開始或結(jié)束。
2) 按照上述算法所解算的衛(wèi)星位置是衛(wèi)星在信號(hào)發(fā)射時(shí)刻的位置,而北斗接收機(jī)位置解算過程中所使用的卻是信號(hào)接收時(shí)刻的衛(wèi)星位置。事實(shí)上,如果直接將求解的衛(wèi)星位置和接收機(jī)所提供的衛(wèi)星位置對(duì)比,會(huì)發(fā)現(xiàn)所有類型衛(wèi)星位置三維坐標(biāo)互差絕對(duì)值的RMS均在350 m左右。原因是信號(hào)接收時(shí)刻的地心地固坐標(biāo)系,是由在信號(hào)發(fā)射時(shí)刻的坐標(biāo)系繞地球自轉(zhuǎn)軸(即坐標(biāo)系的Z軸)自西向東旋轉(zhuǎn)而成的[6-7]。因此,當(dāng)以諾瓦泰北斗接收機(jī)所提供的衛(wèi)星位置作為真值進(jìn)行比較時(shí),要結(jié)合偽距觀測(cè)量,進(jìn)行地球自轉(zhuǎn)誤差的修正。
3) 統(tǒng)一時(shí)間系統(tǒng)。利用廣播星歷所計(jì)算的北斗衛(wèi)星位置采用的是北斗時(shí),而精密星歷以及北斗接收機(jī)所提供的北斗衛(wèi)星位置采用的卻是GPS時(shí),因此對(duì)比時(shí)需要統(tǒng)一時(shí)間系統(tǒng)。本文將時(shí)間統(tǒng)一為GPS時(shí)。
4) 精密星歷按照15 min的時(shí)間間隔給出衛(wèi)星在空間上的三維坐標(biāo),采樣間隔太大,不利于對(duì)比分析。本文采用滑動(dòng)拉格朗日插值法對(duì)所提供的衛(wèi)星位置進(jìn)行插值,以獲得較小的采樣間隔。文獻(xiàn)[8]證明了采用滑動(dòng)式拉格朗日插值法可以獲得mm級(jí)的精度,滿足研究需求。
在解決了上述問題后,便可以分別與精密星歷和北斗接收機(jī)所提供的衛(wèi)星位置進(jìn)行對(duì)比分析。將接收機(jī)所提供的衛(wèi)星位置作為真值,三個(gè)坐標(biāo)互差的絕對(duì)值分別設(shè)為|dX1|、|dY1|、|dZ1|;將精密星歷所提供的衛(wèi)星位置作為真值,三個(gè)坐標(biāo)互差的絕對(duì)值分別設(shè)為|dX2|、|dY2|、|dZ2|。數(shù)據(jù)處理結(jié)果如圖2至圖7及表2所示。
星座類型衛(wèi)星號(hào)|dX1|/|dX2|/m|dY1|/|dY2|/m|dZ1|/|dZ2|/m最大值RMS最大值RMS最大值RMS三維坐標(biāo) RMS1/m RMS2/mGEOC012.25/3.731.99/3.002.76/7.132.38/6.170.22/0.920.07/0.453.10/6.88C032.70/5.262.38/4.501.12/3.490.90/2.920.14/1.070.05/0.622.55/5.40IGSOC073.12/1.981.50/0.901.64/1.670.70/0.712.64/2.371.64/1.652.33/2.01C090.86/1.770.58/1.600.13/1.980.10/1.860.92/0.680.67/0.450.89/2.50MEOC114.10/1.261.83/0.694.19/3.031.93/1.612.66/1.511.06/0.662.86/1.88C124.00/1.972.12/0.802.34/1.971.11/1.123.18/2.511.63/1.322.89/1.91
對(duì)統(tǒng)計(jì)結(jié)果進(jìn)行分析可以看出:1) 將北斗接收機(jī)所提供的北斗衛(wèi)星位置作為真值,經(jīng)過對(duì)地球自轉(zhuǎn)誤差的修正,三維坐標(biāo)RMS均在5 m以內(nèi)。即在位置解算過程中,如果不考慮地球自轉(zhuǎn)的影響,將會(huì)產(chǎn)生較大的誤差;2) 將精密星歷所提供的北斗衛(wèi)星位置作為真值,三維坐標(biāo)RMS均在10 m以內(nèi);以C07、C09號(hào)衛(wèi)星為代表的IGSO衛(wèi)星和以C11、C12號(hào)衛(wèi)星為代表的MEO衛(wèi)星三維RMS小于以C01和C03號(hào)衛(wèi)星為代表的GEO衛(wèi)星,即IGSO和MEO衛(wèi)星精度好于GEO衛(wèi)星。以上結(jié)果證明了所解算的衛(wèi)星位置的正確性,同時(shí)驗(yàn)證了星歷參數(shù)用戶算法的正確性。
本文使用基于線性迭代的最小二乘法解算接收機(jī)位置,因?yàn)樵摲椒ǖ奈恢媒馑銉H使用當(dāng)前數(shù)據(jù),因此復(fù)雜度低,容易實(shí)現(xiàn),是北斗接收機(jī)中最常用的方法[9]。其主要思想是根據(jù)當(dāng)前參數(shù)矩陣迭代當(dāng)前的解,直到誤差小于用戶預(yù)定的門限值,其解算過程如下:
1) 偽距觀測(cè)量修正
從諾瓦泰接收機(jī)中所獲得的偽距ρi(i表示不同的衛(wèi)星)包含了衛(wèi)星鐘差、電離層延時(shí)及對(duì)流層延時(shí)等誤差[10],需對(duì)其進(jìn)行修正。即:
ρi=ri+cδtu-cδts+I+T+ερi,
(1)
(2)
2) 構(gòu)建偽距觀測(cè)量方程
假設(shè)可以得到n顆衛(wèi)星的偽距觀測(cè)量,則偽距觀測(cè)量方程為
(3)
3) 線性化偽距觀測(cè)方程
假設(shè)用戶位置和接收機(jī)鐘差的初始值為(x0,y0,z0,δtu0),基于這個(gè)初始值將式(3)進(jìn)行一階泰勒級(jí)數(shù)展開得到:
(4)
定義如下矢量:
hi
(5)
dx0[(xu-x0),(yu-y0),(zu-z0)]T,
(6)
式(5)中前三列由衛(wèi)星位置與用戶估計(jì)位置之間的方向余弦組成,第四列為常值1;dx0為線性化后初始點(diǎn)和真實(shí)位置之間的修正量。將式(5)、(6)帶入到式(4)并略去高階項(xiàng)整理得:
δρi=hidx0+ερi
(7)
(8)
將式(8)整理為矩陣的形式:
δρ=Hdx0+ερ.
(9)
4) 求解修正量
dx0=(HTH)-1HTδρ.
(10)
5) 更新接收機(jī)位置
x1=x0+dx0
(11)
6) 判斷牛頓迭代的收斂性
以更新后的接收機(jī)位置和鐘差為初始位置,重復(fù)步驟2)-5),對(duì)于第k次更新來說,
(12)
判斷‖dxk-1‖是否小于預(yù)定的門限值,預(yù)定的門限值取決于用戶的精度需求,本文將門限值設(shè)為0.001。如果已經(jīng)小于門限值,則可以將當(dāng)前這一次迭代計(jì)算后的更新值(即xk)作為接收機(jī)的位置和鐘差。即:
xk=xk-1+dxk-1.
(13)
為了驗(yàn)證上述算法的正確性,本文選取天津某已知基準(zhǔn)點(diǎn)2018年3月2日2點(diǎn)至23點(diǎn)諾瓦泰北斗接收機(jī)接收到的偽距信息結(jié)合第2節(jié)所解算的北斗衛(wèi)星位置進(jìn)行接收機(jī)位置解算,并以諾瓦泰接收機(jī)所提供的位置信息作為真值進(jìn)行精度分析。上文提到諾瓦泰FlexPark6接收機(jī)水平定位精度達(dá)到1.5 m(RMS),時(shí)間精度為20 ns,因此可以作為真值來評(píng)估所解算接收機(jī)位置的精度。數(shù)據(jù)處理結(jié)果如圖8至圖10及表3所示。
X軸誤差/mY軸誤差/mZ軸誤差/m 最大值1.813.472.15 平均值0.060.170.09 RMS0.100.270.20 三維RMS0.35
由圖8和圖9可以看出所計(jì)算的接收機(jī)位置分布較為密集,三維方向上波動(dòng)不大。圖10和表3為對(duì)比結(jié)果,可以看出: 1) 三個(gè)方向上誤差最大均不超過5 m,誤差平均值均小于0.5 m,三維RMS也小于0.5 m; 2) 通過圖10可以看出,三個(gè)方向上誤差都集中在1 m以下,個(gè)別時(shí)間段出現(xiàn)較大誤差,經(jīng)分析發(fā)現(xiàn),出現(xiàn)較大誤差的原因是衛(wèi)星個(gè)數(shù)減少,DOP值變大。結(jié)果表明,定位精度均優(yōu)于北斗官方性能規(guī)范中要求的精度[11],證明了所解算接收機(jī)位置的正確性,從而證明了最小二乘法用于北斗接收機(jī)位置解算的可行性。
本文利用實(shí)測(cè)數(shù)據(jù)驗(yàn)證了北斗接收機(jī)位置解算中所涉及到的算法的正確性,并指出衛(wèi)星位置及接收機(jī)位置解算過程中需要注意的問題。數(shù)據(jù)計(jì)算結(jié)果顯示,利用文獻(xiàn)[4]中星歷參數(shù)用戶算法解算出的北斗衛(wèi)星位置以及利用最小二乘法解算出的接收機(jī)位置精度均滿足北斗官方性能規(guī)范要求。接收機(jī)位置解算作為接收機(jī)基本功能之一,所得出的位置精度直接體現(xiàn)了接收機(jī)的性能好壞。因此,對(duì)位置解算過程中涉及到的相關(guān)算法進(jìn)行驗(yàn)證,為北斗接收機(jī)的設(shè)計(jì)和研發(fā)提供了更可靠的理論依據(jù),具有重要的現(xiàn)實(shí)意義。