劉陽,彭競,龔航,胡旖旎,歐鋼
(1.國防科技大學 電子科學學院,湖南 長沙 410073; 2.北京跟蹤與通信技術(shù)研究所,北京100094)
守時實驗室通過建立溯源鏈路將本地時間和遠端高標準時間比對,依據(jù)比對結(jié)果調(diào)整本地時間,進而使本地時間同步到遠端高標準時間上,達到維持本地時間穩(wěn)定、準確的目的[1-2]. 然而,在部分守時實驗室沒有外部溯源鏈路或者溯源鏈路失效的情況下,需要本地鐘組自主守時,維持本地時間的穩(wěn)定和準確. 這對本地鐘組和自主守時算法提高了要求.
原子時算法分為兩大類:加權(quán)平均算法和卡爾曼濾波算法. 加權(quán)平均算法主要包括:國際權(quán)度局(BIPM)采用的ALGOS算法[3-6]、美國國家標準與技術(shù)研究院(NIST)采用的AT1算法[7-8];卡爾曼濾波算法,主要包括:NIST在20世紀80年代使用的標準Jones-Tyron 型卡爾曼濾波算法,和在此基礎(chǔ)上改進得到的一系列卡爾曼濾波算法[9]. AT1算法和卡爾曼濾波算法生成的時間尺度是實時的,適用于實驗室守時.
目前,守時實驗室用于守時的原子鐘主要包括銫鐘和氫鐘. 不同類型的原子鐘在時間保持中有各自的優(yōu)勢[10]. 銫鐘是發(fā)展較早的一種原子鐘,由于其優(yōu)良的長期穩(wěn)定度性能,過去一直作為主要的守時原子鐘. 然而,其不足在于短期穩(wěn)定度較差. 氫鐘是近年來新發(fā)展起來的一種原子鐘,由于其優(yōu)良的短期穩(wěn)定度,以及其各項指標的不斷提升,逐漸作為新型的原子鐘參與守時. 其不足在于長期穩(wěn)定度較差,而且頻率漂移較大[11].
針對上述現(xiàn)狀,研究氫鐘銫鐘聯(lián)合守時算法有一定的實際意義. 有論文對此進行過研究:文獻[10]利用氫鐘短期穩(wěn)定度好的特點,用作測量銫鐘噪聲的參考,對銫鐘噪聲濾波建立時間尺度. 文獻[12]對銫鐘扣除頻差,氫鐘扣除頻差和頻漂后,加權(quán)平均計算綜合原子時駕馭氫鐘. 文獻[13]通過小波變換對氫鐘濾除噪聲,用氫鐘和銫鐘生成時間尺度. 文獻[14]研究了銫原子噴泉鐘駕馭氫鐘的預估算法.
以上氫鐘銫鐘聯(lián)合守時研究中,在注重短期穩(wěn)定度的場合,由于加權(quán)平均算法中對于取權(quán)方式的選擇,往往導致氫鐘權(quán)重較大,而銫鐘權(quán)重較小,銫鐘的長期穩(wěn)定度性能沒有充分發(fā)揮.針對這一問題, 本文提出一種新的綜合利用銫鐘和氫鐘的方案:利用銫鐘數(shù)據(jù)通過卡爾曼濾波計算氫鐘頻漂并且對氫鐘頻漂校準,基于銫鐘數(shù)據(jù)和頻漂校準后的氫鐘數(shù)據(jù),利用單狀態(tài)變量卡爾曼濾波生成本地原子時,在保持原子時短期穩(wěn)定度的前提下,提高其長期穩(wěn)定度性能.這對于守時實驗室自主守時情況下,維持本地時間的長期穩(wěn)定度性能具有意義.
銫鐘模型[15]:
x(t)=x0+y0t+xr(t).
(1)
氫鐘模型:
x(t)=x0+y0t+0.5at2+xr(t).
(2)
式中:x0為初始時刻相位偏差;y0為初始時刻頻率偏差;a為頻率漂移系數(shù);xr(t)為隨機項,也就是原子鐘的噪聲.可見,銫鐘模型是一次多項式加隨機項形式,氫鐘模型是二次多項式加隨機項形式.
氫鐘和銫鐘模型中,多項式部分是確定性部分,隨機項部分是原子鐘噪聲的體現(xiàn). 公認的原子鐘噪聲模型是以下五種噪聲的線性疊加:頻率隨機游走噪聲、頻率閃爍噪聲、頻率白噪聲、相位閃爍噪聲、相位白噪聲[16].
大量實驗顯示,對于典型的地面氫鐘和銫鐘(地面鐘主要是和星載鐘相區(qū)別),其主要的噪聲為頻率白噪聲和頻率隨機游走噪聲[17].
引起時差波動的主要是頻率白噪聲,引起頻差波動的主要是頻率隨機游走噪聲. 銫鐘的頻率白噪聲相對較大,因而銫鐘的短期時差波動較大,短期穩(wěn)定度較差;氫鐘的頻率隨機游走噪聲相對較大,因而氫鐘的頻差波動較大,長期穩(wěn)定度較差.
頻漂是原子頻標的確定性變化分量,很難對其進行直接測量. 一般先直接測量原子鐘的相位或頻率,然后通過下列方法估計頻率漂移項[18-19].
表1 頻漂估計方法
基于相位數(shù)據(jù)的頻漂估計:
二次擬合:用二次多項式對相位差數(shù)據(jù)進行擬合,二次項系數(shù)乘以二得到頻漂估計值.
二次差分均值:對相位差數(shù)據(jù)做二次差分取均值即可得到頻漂估計值.
三點擬合:用相位數(shù)據(jù)的起點、中點和終點估計頻漂,如式(3)所示:
D=4[x(end)-2x(mid)+
x(start)]/(Nt)2.
(3)
基于頻率數(shù)據(jù)的頻漂估計:
線性擬合:用一次多項式對頻率數(shù)據(jù)進行擬合,一次項系數(shù)就是頻漂估計值.
兩段擬合:用頻率數(shù)據(jù)的前一半均值和后一半數(shù)據(jù)均值計算頻漂,如式(4)所示:
D=2[y(2ndhalf)-y(1sthalf)]/(Mt).
(4)
對數(shù)擬合:采用如下模型對頻率數(shù)據(jù)進行擬合:
y(t)=a·ln(bt+1).
(5)
式中:a,b是擬合系數(shù). 根據(jù)式(5),頻率漂移可表示為
(6)
2.1節(jié)提到的是已有的估計頻漂的算法,由于噪聲的存在,導致上述頻漂估計存在一定誤差. 本文提出利用卡爾曼濾波算法對原子鐘噪聲建模,估計原子鐘的頻漂,能夠在噪聲存在的情況下更準確地估計頻漂.
單臺原子鐘的卡爾曼濾波模型[20-25]:設(shè)x1(t),x2(t)和x3(t)分別是原子鐘的三個狀態(tài):相位、頻率和頻漂,x(t)是相位x1(t)的測量值. 單臺原子鐘的模型為
x(t)=x1(t)+n1(t),
(7)
(8)
(9)
(10)
(11)
(12)
(13)
單臺原子鐘的模型的離散形式為
(14)
x(k)=x1(k)+n1(k).
(15)
式中:
w3(kτ)]dt,
(16)
(17)
(18)
驅(qū)動噪聲協(xié)方差矩陣Q和測量噪聲協(xié)方差陣R分別是:
(19)
(20)
根據(jù)以上單臺鐘的卡爾曼濾波模型,估計Q矩陣和R矩陣,并且對狀態(tài)矢量和濾波誤差協(xié)方差矩陣賦初值,然后迭代計算就可以估計出原子鐘的頻漂.
本文提出兩級卡爾曼濾波算法計算原子時:第一級卡爾曼濾波用于估計氫鐘頻漂并且對氫鐘頻漂校準,模型及參數(shù)估計如2.2節(jié)所述;第二級卡爾曼濾波用于生成原子時. 算法流程如圖1所示.
圖1 兩級卡爾曼濾波算法流程圖
第二級卡爾曼濾波算法采用文獻[17]提出的單狀態(tài)變量卡爾曼濾波,即每臺鐘只有鐘差一個狀態(tài)變量.
鐘組的狀態(tài)方程表示為
Xk+1=Φk·Xk+uk.
(21)
式中:Xk=[x1(k),x2(k),…,xN(k)]T是鐘組在k時刻的鐘差;xi(k)(i=1,2,…,N)是鐘組中第i臺鐘的鐘差;Φk=IN×N是狀態(tài)轉(zhuǎn)移矩陣;uk是鐘組的過程噪聲,第i臺鐘的過程噪聲方差表示為
(22)
鐘組的過程噪聲協(xié)方差矩陣表示為
(23)
鐘組的觀測方程表示為
Z(k)=H·X(k)+V(k).
(24)
式中:H表示為
(25)
式中,V是測量噪聲,表示為Vk=[V1(k),V2(k),…,VN-1(k)]T,Vi(k)(i=1,2,…,N-1)是每一組鐘差的觀測噪聲,它們之間相互獨立. 測量噪聲協(xié)方差矩陣表示為
(26)
根據(jù)上述單狀態(tài)變量卡爾曼濾波模型和噪聲協(xié)方差矩陣,給定初始狀態(tài),按照卡爾曼濾波迭代計算就可以計算出每臺鐘的狀態(tài),即每臺鐘和綜合原子時的鐘差,用每臺鐘的鐘差減去卡爾曼濾波結(jié)果就得到本地原子時.
若采用實測數(shù)據(jù),氫鐘的實際頻漂未知. 為了驗證本文提出的卡爾曼濾波算法估計氫鐘頻漂的準確性,本文使用仿真數(shù)據(jù),通過人為給氫鐘加入頻漂的方式,使氫鐘的頻漂變成已知量.
本文采用如式(1)、(2)所示的原子鐘模型,使用文獻[26]中所描述的方法,仿真生成2臺氫鐘和1臺銫鐘的數(shù)據(jù). 仿真數(shù)據(jù)間隔1000 s,仿真生成86 400個數(shù)據(jù)點. 氫鐘和銫鐘的平方擴散系數(shù)σ2按照典型氫鐘和銫鐘的指標選取,如表2所示. 原子鐘確定性部分參數(shù)配置如表3所示. 為了方便分析,設(shè)置觀測噪聲為零.
表2 隨機部分參數(shù)配置
表3 確定性部分參數(shù)配置
將生成的氫鐘和銫鐘數(shù)據(jù)轉(zhuǎn)換成氫鐘相對銫鐘數(shù)據(jù)的相位差. 由于典型地面守時型氫鐘和銫鐘的主要噪聲是頻率白噪聲和頻率隨機游走噪聲,所以在已有的估計頻漂的方法中選用三點擬合法計算頻漂,用于和本文提出的卡爾曼濾波算法估計氫鐘的頻漂作對比分析. 頻漂計算結(jié)果如表4所示.
表4 頻漂真值和頻漂估計結(jié)果對比(天漂)
由表4可見,由于噪聲的存在,三點擬合法估計的氫鐘的頻漂存在一定誤差. 而本文提出的卡爾曼濾波算法通過對噪聲建模,減小了噪聲的影響,估計出的氫鐘的頻漂相對三點擬合法更接近氫鐘真實的頻漂. 表明本文提出的卡爾曼濾波方法估計氫鐘的頻漂更準確.
通過本文提出的卡爾曼濾波算法估計得到氫鐘頻漂后,再通過從鐘差數(shù)據(jù)中減去二次項的方法對氫鐘頻漂校準.
將原始鐘差數(shù)據(jù)分別用本文提出的兩級卡爾曼濾波算法和經(jīng)典的ALGOS算法以及AT1算法計算得到原子時,結(jié)果如圖2所示.
圖2 綜合原子時計算結(jié)果對比
由圖2可以看到,ALGOS算法和AT1算法計算的原子時存在明顯頻漂,這說明氫鐘的頻漂反映在了生成的原子時中. 而用本文提出的兩級卡爾曼濾波算法處理數(shù)據(jù),經(jīng)過對氫鐘頻漂校準后,計算的原子時幾乎不存在頻漂.
通過計算生成的原子時的阿倫方差,衡量生成的原子時的穩(wěn)定度性能,阿倫方差計算結(jié)果如表5所示.
表5 原子時穩(wěn)定度計算結(jié)果
圖3 綜合原子時穩(wěn)定度對比
根據(jù)表5和圖3的結(jié)果可以發(fā)現(xiàn),本文提出的兩級卡爾曼濾波算法計算的原子時相對AT1和ALGOS算法計算的原子時的穩(wěn)定度在千秒穩(wěn)、萬秒穩(wěn)、天穩(wěn)、十天穩(wěn)和月穩(wěn)性能都有提升,特別對于長穩(wěn),十天穩(wěn)的提升在10-15量級,月穩(wěn)的提升在10-14量級. 表明本文提出的兩級卡爾曼濾波算法在對氫鐘頻漂校準后,生成的原子時在保持短期穩(wěn)定度性能的基礎(chǔ)上,長期穩(wěn)定度性能得到了提升.
本文提出的兩級卡爾曼濾波算法,通過對氫鐘頻漂校準,計算得到的本地原子時長期穩(wěn)定度得到提升. 說明在計算本地原子時的時候需要考慮氫鐘頻漂,同時表明本文提出的兩級卡爾曼濾波算法使銫鐘的長期穩(wěn)定度性能得到更有效的發(fā)揮. 這對于自主守時情況下,保持本地原子時的長期穩(wěn)定有一定意義.
為了衡量頻漂校準算法計算出的氫鐘頻漂的準確性,需要知道氫鐘真實的頻漂,而實測數(shù)據(jù)中氫鐘頻漂未知,因而本文使用了仿真數(shù)據(jù).仿真數(shù)據(jù)相對實測數(shù)據(jù)的差異主要是:仿真數(shù)據(jù)只考慮原子鐘主要的兩種噪聲:頻率白噪聲和頻率隨機游走噪聲,不考慮測量噪聲.而實測數(shù)據(jù)中,除了原子鐘這兩種主要噪聲,還存在其他噪聲,以及測量噪聲.針對實測數(shù)據(jù)如何建模處理是進一步的研究工作.