2014版EAL算法中原子鐘的權重分析

屈俐俐,李變

(1.中國科學院 國家授時中心,西安 710600；2.中國科學院 時間頻率基準重點實驗室,西安 710600)

0 引言

國際權度局(BIPM)基于全球70多個守時實驗室自由運轉的原子鐘數據資料,采用ALGOS算法計算得到自由原子時(EAL),然后利用基準頻標(PFS)對EAL進行頻率修正后導出國際原子時(TAI)。

BIPM從1973年起采用ALGOS算法用于EAL的計算,計算時主要考慮EAL的長期穩定度。四十多年來,為了適用原子鐘和比對鏈路性能的提升,BIPM對EAL的計算方法做過多次修改,其中針對ALGOS算法中的加權方法改進的次數最多。最近20年來BIPM對用于EAL計算的ALGOS算法(以下簡稱“EAL算法”)做過幾次大的修改,第一次是在1996年將計算周期由每兩月改為每月計算一次[1],數據點的間隔由每10 d改為每5 d；第二次是在1998年,對算法中的權系統做了調整,由最大絕對權重改為最大相對權重[2]；第三次是2001年最大權采用ωmax=A/N,其中N是每月參加EAL計算的鐘的數量,A是經驗值[3]；最近一次BIPM對EAL算法進行修改是2014年[4],本文主要針對2001年和2014年權重算法修改前后,參與EAL計算的原子鐘權重結構進行分析研究。

1 EAL算法

EAL算法自1973年以來,期間經過幾次修改,但每次修改只是少許調整(計算間隔改變、相對權重替代絕對權重和最大權原則等),一直沿用至2013年。該算法對于沒有頻率漂移,并且主要表現為隨機游走頻率調制噪聲的原子鐘(如：銫原子鐘)是非常適合的。然而,目前參與EAL計算的原子鐘除了5071A銫原子鐘以外,還有大量的氫原子鐘。這種以微波激射振蕩器為基礎的自激型頻標,因具有極純的信號頻譜表現出非常優秀的短期穩定度。但是因物理結構及自身工作特性,氫原子鐘存在頻率漂移,且頻率漂移量相對穩定。為了充分發揮氫原子鐘的優勢,G.Panfilo等人針對2001版EAL算法中的取權方法與頻率預報方法提出了相應的改善方法。

1.1 2001版的EAL算法

為了提高時間尺度的可靠性,2001版EAL算法是在原算法的基礎上,增加了最大權限制原則,解決了時間尺度過多依賴某一臺性能優秀的原子鐘的問題。之后該算法也被大多數時間實驗室用于地方原子時的建立。

1.1.1 2001版的時間尺度計算方法

2001版的時間尺度計算方法基本公式如下[5]：

(1)

xi(t)=EAL(t)-hi(t),

(2)

hi′(t)=ai(t0)+Bip(t-t0),

(3)

式(1)～(3)中:hi(t)為鐘i讀數；xi(t)是EAL(t)和hi(t)的差；hi′(t)是hi(t)預測值；ai(t0)是鐘i相對于EAL(t)在時刻t0的相位差；Bip(t)是鐘i相對于EAL(t)在[t0,t]的預報頻率；ωi(t)表示鐘i的權重。

EAL每個月計算一次,采用的數據是通過遠程時間比對得到的各守時實驗室的鐘與UTC(PTB)的比對結果,數據間隔5 d,數據點對應的時刻是約化儒略日(MJD)尾數為4和9的UTC 00:00的標準歷元。

計算所用的方程為：

Xij(t)=xj(t)-xi(t),i=1,2,…,N,i≠j,

(4)

(5)

t=t0+mT/6,m=0,1,…,6,T=30 d,

(6)[注]在實際計算中,一年中某個月會出現T=35 d的情況,這時t=t0+mT/7,m=0,1,…,7。

式(4)～(6)中,t0是上個時間段的最后一個歸算歷元,也是本次計算的第一個歸算歷元,t是計算EAL的時刻,T是當前計算的時間段。N是參加計算的原子鐘數,ωi(t)表示鐘i的權重。當前時間段[t0,t0+T]的預報速率Bip(t)等于前一時間段[t0-T,t0]上的速率。

1.1.2 2001版的權系統計算方法

(7)

式(7)中,k為時間段索引,Bip,Ik(t)是鐘i在時間段Ik內的速率值,< >表示統計平均。

權值計算公式為:

(8)

雖然性能優秀的鐘應該對EAL的貢獻更大,但是會導致EAL對性能優秀的極少數原子鐘的依賴性過大,降低EAL的可靠性,所以對單鐘的最大權進行限制。ωmax表示最大權,當ωi(t)≥ωmax,則ωi(t)=ωmax。

ωmax=A/N,

(9)

式(9)中,A為經驗常數,在2014版算法之前,A值采用2.5,N是參加計算的鐘數。

1.2 2014版的EAL算法

2001版的EAL算法,速率預報只考慮鐘的相位和頻率影響,即一次預報模型,不考慮頻率漂移。而2014版的EAL算法中速率預報算法加入二次項來描述鐘的頻率漂移,也就是鐘性能采用二次預報模型。

1.2.1 2014版時間尺度計算方法

2014版的hi(t)預測公式如下[6-7]：

(10)

(11)

ai,Ik(tk)=EAL(tk)-hi(tk)=xi(tk)。

(12)

(13)

新增加的頻率漂移參數的求解如下：

(14)

1.2.2 2014版權系統計算方法

在原子時尺度算法中,每個鐘的權重與鐘的頻率穩定度成反比。2014版的權重算法更注重鐘頻率的“可預測性”,因此,只要鐘的頻率漂移足夠穩定,那么它的頻率預測值和實際測量值的偏差就很小。氫原子鐘通常存在明顯的頻率漂移,且頻率漂移的不確定性很小,所以在新算法中氫鐘會獲得較大的權重。

(15)

式(15)中,下標i表示第i臺鐘,Ik表示計算時間間隔。

2014版算法認為新的測量數據具有更可靠的統計特性,因此距離計算時刻越近的測量數據應該給予較大的權重,計算公式如下：

(16)

式(16)中,下標i是第i臺鐘,j代表計算間隔,Mi代表月份,BIPM時間部規定原子時的計算至少需要累積4個月的數據資料,因此4≤Mi≤12。

為了避免某個鐘權重過大,2001版與2014版均采用最大權限制原則。根據文獻[6-7],公式(9)中A值由原來的2.5改為4。

2 權系統修正前后原子鐘的取權結構分析

Panfilo.G.等人在2014版EAL算法中,針對氫原子鐘具有良好短期穩定度及頻率漂移相當穩定的特點,對2001版EAL算法中權重計算提出了相應的改善方法,以提高氫原子鐘在EAL計算中的權重。

2.1 EAL計算中的權分布

為了避免更換算法對EAL穩定度的影響,BIPM采用2006年至2013年的歷史數據進行試算,文獻[4-6]的結果表明：2014版采用平衡了氫鐘和銫鐘的權重分布的算法,提高了氫鐘在EAL計算中的權重。圖1和圖2是BIPM分別采用2001版和2014版EAL算法,計算EAL時不同類型原子鐘的權重分布情況。

圖1 2001版EAL算法中權重的構成 圖2 2014版EAL算法中權重的構成

圖1中基準頻標(PFS)的權重為3.8%,銫鐘(Cs)的權重高達86.4%,而氫鐘(HM)的權重只有9.8%。圖2是2017年4月BIPM采用2014版算法計算EAL時不同類型鐘的權重分布情況,其中基準頻標占4.6%,銫鐘的權重降低到25%,而氫鐘的權重提高到70.4%。不同版本算法中基準頻標所占權重份額較低、量值相對穩定,且與本文分析內容關系較小,所以后面的分析忽略與基準頻標相關的內容。

2.2 原子鐘類型與權重關系

由于每月參加EAL計算的原子鐘數量不完全相同,總鐘數呈緩慢增加趨勢,且不同類型原子鐘數量的比例關系也不完全確定,因此分析歷年來不同類型原子鐘的權重結構變化,只能基于不同類型原子鐘數量的變化比例和權重變化比例的增量來衡量。

圖3和圖4分別是氫鐘和銫鐘在2009年至2017年參加EAL計算時,氫鐘和銫鐘的數量及它們所取權重的比例關系圖。

圖3 2009-2017年參加EAL計算的氫鐘數量及其權重比例

圖4 2009-2017年參加EAL計算的銫鐘數量及其權重比例

由圖3和圖4可以看出：①2009年以來氫鐘占總鐘數數量比例增加不到10%,但氫鐘所占權重比例在采用2014版算法后,從2014年前的10%～15%,大幅度提高到54%,并且之后一路走高,在2017年氫鐘所占權重超過了75%；②2014年之前,銫鐘數量占總鐘數的75%左右,權重基本保持在總權重的85%～88%之間,采用2014版算法后,銫鐘的權重大幅下降,目前約占總權重的20%。

為了分析2014版算法使用前后氫鐘和銫鐘對EAL的貢獻,本文定義鐘的效率=權重/數量(只針對取到權的鐘),統計結果如下：2009年至2013年,氫鐘的效率是50%,銫鐘的效率是114%；2014年至2017年,氫鐘的效率是227%,銫鐘的效率是45%。即氫鐘的效率提高了3.54倍,銫鐘的效率降低了1.53倍。

2.3 不同取權區間內鐘的權重分布

在EAL計算中,利用公式(10)計算出的最大權每月并不相同,為了統計方便在本文中對權重進行歸一化計算,最大權為1。為了下面分析的需要,將權重區間按0.2間隔劃分,即ω=0(取0權的鐘),0<ω≤0.2,0.2<ω≤0.4,0.4<ω≤0.6,0.6<ω≤0.8,0.8<ω<1,ω=1(取最大權的鐘)。

圖5和圖6分別是2011年至2017年EAL計算中,不同權重區間內的原子鐘數量占總鐘數的百分比和原子鐘權重占總權重的百分比。

圖5 2011-2017年不同取權區間內原子鐘數量的分布

由圖5可以看出：①兩種算法在ω=1,0.2<ω≤0.4和ω=0 3個權值區間內,原子鐘數量所占比例相當；②2014版算法在0.8<ω<1,0.6<ω≤0.8和0.4<ω≤0.6 3個權值區間內,原子鐘數量所占比例明顯低于2001版算法：2014版算法中落在這3個權值區間內鐘數量總和所占比例不足5%,而2001版算法中鐘數量總和所占比例約為25%；③在0<ω≤0.2權值區間,2014版算法中鐘數量所占比例超過50%,2001版算法中鐘數量所占比例約為35%。

圖6 2011-2017年不同權值區間內鐘權重的分布

如圖6所示：①2001版算法中取得最大權(ω=1)的鐘所占的權重比例約為40%,且該比例比較穩定；②在采用2014版算法后,取得最大權(ω=1)的鐘所占權重的比例從2014年的55%左右開始逐年增長,到2017年達到64%左右；③在0<ω<1區間,2001版算法中權值分布比2014版算法更合理。

2.4 不同取權區間內鐘類型的分布

圖7和圖8分別是2011年至2017年在EAL計算中,氫鐘和銫鐘在不同權值區間鐘數量的分布比例。

由圖7可知：①在ω=1權值,2001版算法中氫鐘數量約為10%,2014版算法中氫鐘數量在40%左右；②在0<ω≤0.2權值區間,2001版算法中氫鐘數量超過50%,2014版算法中氫鐘數量在25%左右；③在其他權值區間,兩種算法中氫鐘數量相當。

由圖8可知：①在ω=1權值,2001版算法中銫鐘數量約為15%,2014版算法中銫鐘數量小于1%；②在0<ω≤0.2權值區間,2001版算法中銫鐘數量小于35%,2014版算法中銫鐘數量在70%左右；③在0.4<ω≤0.6,0.6<ω≤0.8和0.8<ω<1權值區間,2001版算法中銫鐘數量總和在30%左右,2014版算法中銫鐘數量總和小于4%；④在0.2<ω≤0.4權值區間,兩種算法中銫鐘數量相當。

圖7 2011-2017年EAL計算中氫鐘在不同權值區間占氫鐘總數比例

圖8 2011-2017年EAL計算中銫鐘在不同權值區間占銫鐘總數比例

3 結論

參與EAL計算的原子鐘數據來源于全球70多個守時實驗室,這些實驗室配置的守時原子鐘主要是銫原子鐘和主動型氫原子鐘,這兩類鐘占參與EAL計算的總鐘數的98%,一般情況下,氫鐘的數量占總鐘數的1/3。

銫鐘頻率波動較大,幾乎沒有頻率漂移；氫鐘存在頻率漂移,但大多數氫鐘的頻率漂移量相對穩定。因此氫鐘在采用2014版算法計算EAL時,可以取得較高的權重。

BIPM提出2014版算法的目的是：①利用氫鐘良好的頻率預測特性,提高EAL的穩定性和準確度,這已在BIPM相關文獻中證明[6,8-11]；②平衡氫銫兩種原子鐘在EAL計算時的權重比例。

通過分析2014年以來3年多EAL的計算結果,該算法首先提高了氫鐘在計算EAL時的權重,氫鐘權重份額在這3年當中呈現不斷增長的趨勢,從2014年的54%提高到2017年的75%,權重在3年間增加了21%,況且根據圖3和圖4該比例仍將會繼續提高,這樣,未來幾年很有可能會出現EAL算法修改前的銫鐘權重完全占據主導地位的狀況,只不過這次是由銫鐘更換為氫鐘而已；其次,2014版算法中取最大權(ω=1)的鐘數量不超過15%,但是所占權重比例過高；權值在0<ω≤0.2區間的鐘數量偏大；在0.2<ω<1區間的鐘數量甚少。這是由于最大權門限和經驗常數A不匹配導致的,關于二者關系的分析將在另一篇文獻中給出。

通過以上的分析,對采用類ALGOS算法中權重的確定提出以下建議：

①依據守時鐘組的配置情況,適當平衡不同類型鐘的權重；

②確定取最大權鐘的數量及其所占權重份額,取最大權的鐘的數量不超過總鐘數的15%,所占權重份額約50%左右為宜[11]；

③分析落在不同權值區間鐘的數量、鐘類型的分布狀況,避免出現2014版算法中兩頭過大,中間區域幾乎沒有鐘的超級啞鈴狀的情況出現；

④算法投入使用前,應對算法使用后可能出現的狀況進行評估,避免權重份額出現趨勢性改變;

⑤經驗常數A的確定需經過大量試算。A值過大會使計算的時間尺度只依賴于少數性能好的鐘,從而降低時間尺度的可靠性；反之,又會使性能優秀的鐘對時間尺度的貢獻過小,降低時間尺度的準確性和穩定性。

綜上所述,守時實驗室在借鑒其他時間尺度算法時,不可直接照搬,需結合自身原子鐘資源的情況進行詳細分析并進行大量試算,最終找到滿足守時實驗室要求的時間尺度算法。