基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波在北斗單向授時中的應(yīng)用

柯騰倫,孫光輝,吳榮剛,丁凱生,陳睿容

(1.北京遙感設(shè)備研究所,北京 100086; 2.西安電子科技大學(xué)雷達(dá)信號處理全國重點實驗室,西安 710071)

0 引言

當(dāng)今世界時間信號的精度在工程技術(shù)、基礎(chǔ)研究及日常生活等領(lǐng)域中的作用愈發(fā)凸顯,對授時精度的要求越來越高[1]。傳統(tǒng)授時方法如網(wǎng)絡(luò)時間協(xié)議授時法受存儲空間及計算能力等限制,授時范圍和授時精度不足[2],而通過全球衛(wèi)星導(dǎo)航系統(tǒng)(global navigation satellite system, GNSS)進(jìn)行授時具有低成本、高精度、廣域性、開放性和可對移動用戶服務(wù)等優(yōu)點,成為了當(dāng)前最主要的授時工具,并得到了廣泛的研究[3]。

GNSS單向授時是指接收機(jī)通過觀測量自主解算接收機(jī)鐘差,并利用其對接收機(jī)秒信號產(chǎn)生電路進(jìn)行修正,輸出標(biāo)準(zhǔn)時間信號以實現(xiàn)授時的方法[4]。GNSS單向授時又可以分為基于偽距測量的GNSS單向授時和GNSS精密單點定位 (precise point positioning, PPP)授時[5]。

基于偽距測量的GNSS單向授時使用廣播星歷和廣播星鐘,授時精度一般在10～20 ns。劉利等[6]對BDS-3提供的各種授時服務(wù)進(jìn)行了比較分析,結(jié)果表明RNSS(radio navigation satellite service)授時精度和RDSS(radio determination satellite service)雙向授時精度相當(dāng),介于SBAS(satellite-based augmentation service)授時精度和RDSS單向授時精度之間,在9ns左右;辛潔等[7]研究了BDS不同服務(wù)體制下的授時差異,從時空基準(zhǔn)、授時精度檢核、設(shè)備延時及衛(wèi)星健康狀態(tài)等方面,分析了基于已知位置的RNSS單/雙頻授時等模型之間的差異,結(jié)果表明RNSS單/雙頻授時精度分別可以達(dá)到8.56 ns和6.89 ns。張大眾等[8]利用iGMAS站數(shù)據(jù)分析了不同軌道衛(wèi)星的BDS單星授時精度,結(jié)果表明GEO、IGSO、MEO衛(wèi)星的授時精度分別為27.39、18.37和18.62 ns。

GNSS-PPP授時主要使用載波相位、偽距觀測數(shù)據(jù)以及高精度軌道和星鐘差等產(chǎn)品,授時精度相對較高。易卿武等[9]提出一種基于北斗三號B2b信號的精密單點授時方法,根據(jù)PPP算法實時估算接收機(jī)鐘差,結(jié)果表明與傳統(tǒng)的精密單點授時方法相比, 該方法授時精度更高,并且具有成本低、不依賴地面通信網(wǎng)絡(luò)和分析中心的優(yōu)勢。馬祥泰等[10]針對消電離層PPP模型組合噪聲大的問題,采用非組合PPP模型對BDS精密授時精度開展研究,結(jié)果表明在靜態(tài)和實時動態(tài)條件下,非組合PPP模型精密授時精度均優(yōu)于消電離層組合PPP模型。韓金陽等分析了不同機(jī)構(gòu)差分碼偏差(differential code bias, DCB)產(chǎn)品對選定的2個測站PPP授時精度的影響,實驗結(jié)果表明所選的2個測站使用不同機(jī)構(gòu)的DCB產(chǎn)品估計鐘差的均方差和時間偏差均優(yōu)于0.4 ns[11]。劉根友等[12]提出了基于PPP的云平臺高精度授時方法,采用PPP技術(shù)實時解算授時終端鐘差,通過馴服恒溫晶振輸出亞納秒精度的1 PPS(pulse per second),實現(xiàn)了長時間高精度的授時能力,精密單點授時精度優(yōu)于1 ns。

BDS目前已經(jīng)建成,其授時應(yīng)用研究已成為重要的發(fā)展方向,雖然與PPP授時相比,基于偽距測量的單向授時受星歷誤差、信號傳播延時誤差和天線位置誤差的影響授時精度相對較低,但是其不需要使用實時的精密衛(wèi)星軌道和鐘差改正數(shù)產(chǎn)品,也無需較長的時間等待PPP結(jié)果收斂[13],因此對基于偽距測量的BDS單向授時進(jìn)行研究以提高其授時精度具有重要的實際意義。本文通過基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波對經(jīng)由偽距解算獲得的鐘差進(jìn)行降噪處理以達(dá)到提高授時精度的效果。

1 基于偽距測量的BDS單向授時原理

在BDS單向授時中用戶時鐘與系統(tǒng)時間的精確同步以及系統(tǒng)時間與UTC時間的精確同步是實現(xiàn)高精度授時的關(guān)鍵所在[14]。單向授時法通過修正接收機(jī)本地時間tu來得到協(xié)調(diào)世界時tUTC,其計算公式如式(1)

tUTC=tu-δtu-δtUTC

(1)

其中,tUTC表示協(xié)調(diào)世界時;tu表示接收機(jī)本地時間;δtu=tu-tsys表示接收機(jī)本地時間相對于系統(tǒng)時間的偏移量,即接收機(jī)鐘差;δtUTC=tsys-tUTC表示系統(tǒng)時間相對于協(xié)調(diào)世界時的偏移量;tsys表示系統(tǒng)時間。

δtUTC通過導(dǎo)航電文中的模型參數(shù)計算獲得。如果用戶接收機(jī)靜止且天線位置已知,接收機(jī)只要通過跟蹤單顆可見衛(wèi)星,建立單個觀測方程就可以求解出δtu;如果用戶接收機(jī)位置未知,接收機(jī)需要通過跟蹤4顆或更多可見衛(wèi)星,通過最小二乘法求解出δtu。下文討論的都是用戶接收機(jī)位置未知情況下的BDS單向授時方法。

雖然單向授時在工程上被廣泛應(yīng)用,但是由偽距求解所得接收機(jī)鐘差δtu容易受噪聲影響變得粗糙、雜亂,所以直接使用它來校正接收機(jī)本地時間tu所得系統(tǒng)時間tsys的誤差較大,授時精度不高,時鐘抖動明顯。因此本文擬通過將接收機(jī)鐘差建模為鐘差(接收機(jī)時鐘與參考時鐘的初始相位偏差)、鐘速(接收機(jī)時鐘與參考時鐘的頻率偏差)和頻率漂移(接收機(jī)時鐘頻率的變化率),采用基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波對其進(jìn)行處理以降低噪聲引起的抖動,達(dá)到提高授時精度的效果。

2 基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波方程

晶振頻率在老化和溫度的影響下基本呈現(xiàn)二次型曲線變化,經(jīng)典鐘差模型的表達(dá)式如式(2)[15-16]

(2)

其中,δtu(t)表示鐘差;a表示初始相位偏差;b表示頻率偏差;c表示線性頻率漂移;εx(t)表示鐘差中的隨機(jī)分量。

對式(2)求關(guān)于時間t的二階導(dǎo)數(shù),可以得到頻率漂移關(guān)于時間的表達(dá)式

(3)

其中,εz(t)表示頻率漂移中的隨機(jī)分量。

圖1所示為通過本地時鐘解算鐘速前后歷元差分計算獲得的頻率漂移,說明本地時鐘頻率漂移在零值附近上下波動,該波動呈現(xiàn)出隨機(jī)噪聲的特性。根據(jù)式(3)和圖1的結(jié)果,接收機(jī)時鐘的頻率漂移可以建模為常數(shù)加隨機(jī)分量的形式。

圖1 本地時鐘解算鐘速前后歷元差分計算獲得的頻率漂移Fig.1 The frequency drift calculated by epoch difference of the local clock speed before and after resolution

“當(dāng)前”統(tǒng)計模型原本是作為一種機(jī)動目標(biāo)跟蹤算法在Singer模型的基礎(chǔ)上提出的。它將目標(biāo)加速度建模為加速度均值加零均值有色噪聲的形式,因此可以參照“當(dāng)前”統(tǒng)計模型對于目標(biāo)加速度的建模,將本地時鐘頻率漂移建模如式(4)[17]

(4)

(1)狀態(tài)方程

將接收機(jī)鐘差建模為鐘差、鐘速和頻率漂移并對其進(jìn)行基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波,其狀態(tài)方程如式(5)

(5)

“當(dāng)前”統(tǒng)計模型中狀態(tài)轉(zhuǎn)移矩陣A的計算公式如式(6)[17]

(6)

“當(dāng)前”統(tǒng)計模型中輸入關(guān)系矩陣U的計算公式如式(7)[17]

(7)

在本模型中,系統(tǒng)輸入量為接收機(jī)鐘差和鐘速,系統(tǒng)狀態(tài)量為接收機(jī)鐘差,鐘速和頻率漂移,因此輸入關(guān)系矩陣B的計算公式如式(8)

(8)

式(6)～(8)中的T表示采樣周期;e表示自然常數(shù)。

(2)觀測方程

基于“當(dāng)前”統(tǒng)計模型卡爾曼濾波的觀測方程如式(9)

yk=Cxk+vk

(9)

3 基于“當(dāng)前”統(tǒng)計模型卡爾曼濾波的BDS單向授時法

基于“當(dāng)前”統(tǒng)計模型卡爾曼濾波的BDS單向授時法的主要流程如圖2所示。

圖2 基于“當(dāng)前”統(tǒng)計模型卡爾曼濾波的BDS單向授時法流程Fig.2 BDS one-way timing service process based on the current statistical model Kalman filter

1)通過偽距、偽距變化率觀測方程分別解算獲得接收機(jī)鐘差和鐘速。

3)利用濾波后的接收機(jī)鐘差和鐘速對本地時鐘進(jìn)行校正。

3.1 濾波過程

(1)時間更新過程

狀態(tài)矢量先驗估計值的計算公式如式(10)

(10)

狀態(tài)矢量先驗估計值均方誤差陣的計算公式如式(11)

(11)

“當(dāng)前”統(tǒng)計模型中過程噪聲矢量協(xié)方差矩陣Qk的計算公式如式(12)[17]

(12)

(13)

其中,amax表示給定的最大頻率漂移。

(2)測量更新過程

增益矩陣的計算公式如式(14)

(14)

其中,Kk表示增益矩陣;Rk表示測量噪聲矢量的協(xié)方差矩陣。

測量噪聲矢量的協(xié)方差矩陣Rk通過采樣統(tǒng)計的方法進(jìn)行計算

(15)

狀態(tài)矢量后驗估計值的計算公式如式(16)

(16)

狀態(tài)矢量后驗估計值的均方誤差陣的計算公式如式(17)

(17)

3.2 本地時鐘修正

圖3所示為接收機(jī)本地時鐘修正示意圖,包含本地時鐘修正和偽距變化率修正。

(1)本地時鐘修正

首歷元僅使用鐘差修正相位值,往后使用鐘差修正相位值的同時使用鐘速修正頻率控制字。首次修正后的每個歷元都要進(jìn)行相位修正,其計算方法如式(18)

(18)

從第2個歷元開始,間隔一拍使用鐘速修正頻率控制字,其計算方法如式(19)

(19)

(2)偽距變化率的修正

由于在進(jìn)行積分多普勒計算時包含中頻,所以在計算偽距變化率時需要將中頻對應(yīng)的載波相位積分值減去

(20)

上文對本地時鐘的修正只改變觀測量的鎖存時間,但是作為載波環(huán)和碼環(huán)頻率基準(zhǔn)的板上晶振或外部時鐘并沒有被改變,所以需要使用鐘速對偽距變化率進(jìn)行額外修正

(21)

此外對本地時鐘的相位修正會使得積分時長改變,導(dǎo)致載波相位積分值發(fā)生變化從而影響偽距變化率,對此同樣需要加以修正,修正方式如式(22)

(22)

4 實驗結(jié)果與授時性能分析

為了更直觀清晰地體現(xiàn)基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波對于BDS單向授時精度的改善效果,設(shè)計了相關(guān)實驗并進(jìn)行了數(shù)據(jù)采集與分析。

“當(dāng)前”統(tǒng)計模型卡爾曼濾波的效果與機(jī)動頻率和最大頻率漂移的取值有關(guān),所以需要研究不同機(jī)動頻率和最大頻率漂移下的授時精度。為了方便對實驗結(jié)果進(jìn)行統(tǒng)計分析并排除時間因素的影響,所以選擇2022年9月10日至30日每天的同一時段(10∶00至11∶00)采集1小時的數(shù)據(jù)。天線架設(shè)于北京的某大樓樓頂,由于自行研制的導(dǎo)航接收機(jī)只能接收北斗B1I頻點,所以實驗結(jié)論都是基于北斗B1I信號。

(23)

為了研究算法對于不同類型接收機(jī)參考晶振的適用性,實驗中分別采用溫補(bǔ)晶振(TCXO)和恒溫晶振(OCXO)作為接收機(jī)的時鐘源采集了相關(guān)的實驗數(shù)據(jù),實驗組成結(jié)構(gòu)圖如圖4所示。

圖4 實驗組成結(jié)構(gòu)圖Fig.4 Experimental composition structure diagram

表1和表2分別表示溫補(bǔ)晶振和恒溫晶振所對應(yīng)的15組使用經(jīng)濾波處理后的鐘差和鐘速對本地時鐘進(jìn)行修正得到的鐘差最大波動范圍的統(tǒng)計結(jié)果,圖5和圖6是其對應(yīng)曲線。

表1 本地時鐘鐘差最大波動范圍統(tǒng)計結(jié)果(溫補(bǔ)晶振)Tab.1 Statistical results of maximum fluctuation range of local clock error(TCXO) ns

表2 本地時鐘鐘差最大波動范圍統(tǒng)計結(jié)果(恒溫晶振)Tab.2 Statistical results of maximum fluctuation range of local clock error(OCXO) ns

圖5 不同機(jī)動頻率和最大頻率漂移對應(yīng)的本地時鐘鐘差最大波動范圍(溫補(bǔ)晶振)Fig.5 Maximum fluctuation range of local clock error corresponding to different maneuvering frequencies and maximum frequency drift(TCXO)

圖6 不同機(jī)動頻率和最大頻率漂移對應(yīng)的本地時鐘鐘差最大波動范圍(恒溫晶振)Fig.6 Maximum fluctuation range of local clock error corresponding to different maneuvering frequencies and maximum frequency drift(OCXO)

對于溫補(bǔ)晶振,使用未經(jīng)濾波的解算鐘差和鐘速直接對本地時鐘進(jìn)行修正得到的鐘差最大波動范圍為-31.23～35.56 ns,峰峰值為66.79 ns。濾波修正后所得鐘差的峰峰值與機(jī)動頻率和最大頻率漂移的取值有關(guān)。當(dāng)最大頻率漂移為1 ns/s2且機(jī)動頻率為10-2時最小為42.27 ns,上下界分別為21.89 ns和-20.38 ns,相比于未濾波修正情況減小了36.7%;當(dāng)最大頻率漂移為10 ns/s2且機(jī)動頻率為10-1時最大為70.69 ns,上下界分別為33.37 ns和-37.32 ns,相比于未濾波修正情況增加了5.8%。對于恒溫晶振,使用未經(jīng)濾波的解算鐘差和鐘速直接對本地時鐘進(jìn)行修正得到的鐘差最大波動范圍為-10.23～16.20 ns,峰峰值為26.43 ns。濾波修正后所得鐘差的峰峰值與機(jī)動頻率和最大頻率漂移的取值有關(guān)。當(dāng)最大頻率漂移為0.1 ns/s2且機(jī)動頻率為101時最小為7.452 ns,上下界分別為4.881 ns和-2.571 ns,相比于未濾波修正情況減小了71.8%;當(dāng)最大頻率漂移為10 ns/s2且機(jī)動頻率為101時最大為24.300 ns,上下界分別為12.72 ns和-11.58 ns,相比于未濾波修正情況減小了8.0%。

表3和表4分別表示溫補(bǔ)晶振和恒溫晶振所對應(yīng)的15組使用經(jīng)濾波處理后的鐘差和鐘速對本地時鐘進(jìn)行修正得到的鐘差標(biāo)準(zhǔn)差的統(tǒng)計結(jié)果,圖7和圖8是其對應(yīng)曲線。

表4 本地時鐘鐘差標(biāo)準(zhǔn)差統(tǒng)計結(jié)果(恒溫晶振)Tab.4 Statistical results of standard deviation of local clock error(OCXO) ns

圖7 不同機(jī)動頻率和最大頻率漂移對應(yīng)的本地時鐘鐘差標(biāo)準(zhǔn)差(溫補(bǔ)晶振)Fig.7 Standard deviation of local clock error corresponding to different maneuvering frequencies and maximum frequency drift(TCXO)

圖8 不同機(jī)動頻率和最大頻率漂移對應(yīng)的本地時鐘鐘差標(biāo)準(zhǔn)差(恒溫晶振)Fig.8 Standard deviation of local clock error corresponding to different maneuvering frequencies and maximum frequency drift(OCXO)

對于溫補(bǔ)晶振,使用未經(jīng)濾波的解算鐘差和鐘速直接對本地時鐘進(jìn)行修正得到的鐘差標(biāo)準(zhǔn)差為9.552 ns。濾波修正后所得鐘差的標(biāo)準(zhǔn)差同樣與機(jī)動頻率和最大頻率漂移的取值有關(guān)。當(dāng)最大頻率漂移為1 ns/s2且機(jī)動頻率為101時最小為5.489 ns,相比于未濾波修正情況減小了42.5%;當(dāng)最大頻率漂移為10 ns/s2且機(jī)動頻率為10-1時最大為9.287 ns,相比于未濾波修正情況減小了2.8%。對于恒溫晶振,使用未經(jīng)濾波的解算鐘差和鐘速直接對本地時鐘進(jìn)行修正得到的鐘差標(biāo)準(zhǔn)差為4.075 ns。濾波修正后所得鐘差的標(biāo)準(zhǔn)差同樣與機(jī)動頻率和最大頻率漂移的取值有關(guān)。當(dāng)最大頻率漂移為0.1 ns/s2且機(jī)動頻率為101時最小為1.345 ns,相比于未濾波修正情況減小了67.0%;當(dāng)最大頻率漂移為10 ns/s2且機(jī)動頻率為100時最大為3.687 ns,相比于未濾波修正情況減小了9.5%。

由以上分析可知,相比于未濾波修正,在選取合適的機(jī)動頻率和最大頻率漂移的情況下,濾波修正得到的本地時鐘鐘差的峰峰值以及標(biāo)準(zhǔn)差都能夠得到明顯的改善。實驗結(jié)果表明,對于溫補(bǔ)晶振,本地時鐘鐘差的峰峰值最多可以減小36.7%,而標(biāo)準(zhǔn)差最多可以減小42.5%;對于恒溫晶振,本地時鐘鐘差的峰峰值最多可以減小71.8%,而標(biāo)準(zhǔn)差最多可以減小67.0%。另外,對于分別以溫補(bǔ)晶振和恒溫晶振作為時鐘源的接收機(jī)時鐘,在算法達(dá)到最佳改善效果時,對應(yīng)的最大頻率漂移和機(jī)動頻率的取值是不同的。由此可見濾波修正授時方案的授時性能相較于未濾波修正授時方案的授時性能有明顯的提高,基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波對于BDS單向授時精度的改善效果十分明顯。

以上結(jié)論是在使用溫補(bǔ)晶振和恒溫晶振作為接收機(jī)時鐘源、接收北斗B1I信號的條件下得出的。理論上本文所提出的方法不限定接收機(jī)、晶振與導(dǎo)航信號頻點,適用于所有單向授時方法,但是其他條件下對應(yīng)濾波算法的參數(shù)選取會有一定差異,算法的改善效果有待通過后續(xù)實驗進(jìn)行驗證。

5 結(jié)論

本文針對由偽距誤差與噪聲引起的接收機(jī)解算鐘差不夠平滑、準(zhǔn)確導(dǎo)致用戶接收機(jī)本地時間校準(zhǔn)所得系統(tǒng)時間的誤差較大、授時精度不高的問題,提出了一種基于“當(dāng)前”統(tǒng)計模型的卡爾曼濾波單向授時方法。對解算獲得的鐘差和鐘速進(jìn)行濾波并將濾波結(jié)果修正到本地時鐘以提高授時精度。使用以溫補(bǔ)晶振和恒溫晶振作為時鐘源的接收機(jī)進(jìn)行實驗,結(jié)果表明:

1)相比于未濾波修正,在選取合適的機(jī)動頻率和最大頻率漂移的情況下,溫補(bǔ)晶振時鐘的峰峰值最多可以減小36.7%,標(biāo)準(zhǔn)差最多可以減小42.5%;恒溫晶振時鐘的峰峰值最多可以減小71.8%,標(biāo)準(zhǔn)差最多可以減小67.0%。

2)對于以不同晶振作為時鐘源的接收機(jī),算法達(dá)到最佳改善效果時對應(yīng)的最大頻率漂移和機(jī)動頻率的取值有所不同。

本文提出的方法實現(xiàn)簡單,成本低廉,只需要在接收機(jī)終端添加相應(yīng)算法,就能夠提升基于偽距測量的GNSS單向授時的授時精度,具有工程實用價值。