宋會(huì)杰,董紹武,2,王翔，王燕平,張繼海,屈俐俐，趙書紅,張首剛,2

(1.中國(guó)科學(xué)院 國(guó)家授時(shí)中心，西安 710600；2.中國(guó)科學(xué)院大學(xué) 天文與空間科學(xué)學(xué)院，北京 100049)

0 引言

原子鐘的噪聲及守時(shí)性能分析是守時(shí)工作中非常重要的一項(xiàng)研究工作。原子鐘噪聲是影響時(shí)間穩(wěn)定性的重要因素，分析原子鐘的噪聲特性，研究降低原子鐘噪聲的算法可以提升守時(shí)性能[1-4]。守時(shí)性能主要表現(xiàn)為產(chǎn)生標(biāo)準(zhǔn)時(shí)間的準(zhǔn)確性、穩(wěn)定性和實(shí)時(shí)性。利用原子鐘的比對(duì)測(cè)量數(shù)據(jù)分析原子鐘的性能，研究守時(shí)算法提升標(biāo)準(zhǔn)時(shí)間長(zhǎng)期性能[5-9]。當(dāng)前國(guó)內(nèi)外研究人員對(duì)守時(shí)算法的研究不斷取得新的進(jìn)展，比如2014年后，國(guó)際權(quán)度局(BIPM)對(duì)原子鐘進(jìn)行二次模型預(yù)報(bào)并根據(jù)可預(yù)測(cè)性方法取權(quán)，提升了標(biāo)準(zhǔn)時(shí)間的整體性能[10-11]。目前我國(guó)采用的守時(shí)算法也是基于可預(yù)測(cè)性方法取權(quán)。但是上述研究和算法都是基于國(guó)外原子鐘，比如美國(guó)的銫原子鐘、氫原子鐘和俄羅斯的氫原子鐘等，由于國(guó)產(chǎn)原子鐘的工作原理、材料、制造工藝與國(guó)外原子鐘不同，導(dǎo)致噪聲特性存在一定差異，因此研究國(guó)產(chǎn)鐘的守時(shí)性能，提出符合國(guó)產(chǎn)鐘特性的守時(shí)算法，解決我國(guó)標(biāo)準(zhǔn)時(shí)間自主產(chǎn)生和時(shí)間基準(zhǔn)保持的核心技術(shù)是當(dāng)前面臨的迫切問題。文中基于小波分解研究了國(guó)產(chǎn)原子鐘不同尺度下的噪聲情況，并與國(guó)外原子鐘進(jìn)行了比較，計(jì)算出不同尺度下的噪聲強(qiáng)度。通過經(jīng)典加權(quán)算法與Kalman濾波算法分別計(jì)算基于國(guó)產(chǎn)鐘的時(shí)間尺度，并比較了時(shí)間尺度的準(zhǔn)確度、穩(wěn)定度和實(shí)時(shí)性，整體上分析了基于國(guó)產(chǎn)鐘的守時(shí)性能。

1 國(guó)產(chǎn)原子鐘的噪聲分析

小波變換作為一種時(shí)頻局部化方法，其窗口是可變的，即在低頻部分具有較低的時(shí)間分辨率和較高的頻率分辨率，具有對(duì)信號(hào)的自適應(yīng)性。應(yīng)用于國(guó)產(chǎn)原子鐘的測(cè)量數(shù)據(jù)，分析不同尺度的噪聲強(qiáng)度。

1.1 小波變換的原理

(1)

‖Ψjk(t)‖=‖Ψ(t)‖j,k∈Z，

(2)

通過對(duì)國(guó)產(chǎn)鐘測(cè)量數(shù)據(jù)進(jìn)行離散小波變換，可以對(duì)數(shù)據(jù)進(jìn)行時(shí)頻分析，即測(cè)量數(shù)據(jù)在某一時(shí)間上，對(duì)應(yīng)某一尺度a0的成分。相應(yīng)的離散小波變換可表示為

(3)

對(duì)于離散小波變換，一個(gè)3尺度的分解如圖1所示。從圖1可以看出小波變換分離信號(hào)的不同成分，其中不同尺度的低頻系數(shù)表示去噪后的數(shù)據(jù)，不同尺度的高頻系數(shù)表示噪聲，對(duì)應(yīng)國(guó)產(chǎn)原子鐘不同頻率段的噪聲。

圖1 3尺度分解的組織形式

1.2 基于小波變換的國(guó)產(chǎn)鐘噪聲分析

基于小波變換理論對(duì)國(guó)產(chǎn)鐘的測(cè)量數(shù)據(jù)進(jìn)行噪聲分析，并與國(guó)外原子鐘的測(cè)量噪聲進(jìn)行比較。國(guó)產(chǎn)原子鐘編號(hào)選為Cs2025、Cs3050和Cs3059，國(guó)外原子鐘編號(hào)選為Cs3089。測(cè)量數(shù)據(jù)時(shí)間為2020年8月至2021年6月，測(cè)量參考為UTC(NTSC)。利用db5小波對(duì)國(guó)產(chǎn)鐘測(cè)量數(shù)據(jù)進(jìn)行6尺度分解，得到不同尺度的原子鐘噪聲，其結(jié)果如圖2，圖3，圖4和圖5所示。為了清晰，圖中只表示了3尺度分解。通過小波多尺度分解得出不同尺度噪聲的標(biāo)準(zhǔn)差，其結(jié)果如表1和表2所示。結(jié)合圖與表可以看出尺度越低，噪聲越小，隨著尺度的增大，噪聲變大。小波分解中尺度越小，對(duì)應(yīng)的頻率越高，尺度越大，對(duì)應(yīng)的頻率越低。隨著頻率的降低，噪聲變大，說明國(guó)產(chǎn)鐘在低頻段噪聲明顯，低頻段的噪聲對(duì)國(guó)產(chǎn)鐘性能的影響起主要作用。守時(shí)算法需重點(diǎn)考慮降低國(guó)產(chǎn)鐘低頻噪聲的影響。比較國(guó)產(chǎn)鐘Cs2025，Cs3050和Cs3059不同尺度的噪聲大小，在尺度1至尺度3，Cs2025和Cs3050的噪聲明顯低于Cs3059的噪聲，說明國(guó)產(chǎn)鐘Cs3059的高頻噪聲相對(duì)偏大。分析尺度4至6，Cs2025的噪聲明顯高于Cs3050和Cs3059，說明國(guó)產(chǎn)鐘Cs2025的低頻噪聲相對(duì)偏大。3臺(tái)原子鐘，Cs3050的低頻噪聲與高頻噪聲都相對(duì)偏小，表示Cs3050的長(zhǎng)期穩(wěn)定度和短期穩(wěn)定度最高。通過上述分析，得出不同原子鐘的噪聲情況不同，部分國(guó)產(chǎn)鐘低頻噪聲相對(duì)偏大，部分國(guó)產(chǎn)鐘高頻噪聲相對(duì)偏大。守時(shí)過程中，根據(jù)國(guó)產(chǎn)原子鐘不同尺度的噪聲情況進(jìn)行加權(quán)處理，能夠同時(shí)降低不同尺度噪聲的影響，提升守時(shí)的整體性能。

圖2 Cs2025不同尺度的噪聲

圖3 Cs3050不同尺度的噪聲

圖4 Cs3059不同尺度的噪聲

圖5 Cs3089不同尺度的噪聲

表1 原子鐘尺度(1～3)的噪聲標(biāo)準(zhǔn)差 單位：ns

表2 原子鐘尺度(4～6)的噪聲標(biāo)準(zhǔn)差 單位：ns

2 國(guó)產(chǎn)原子鐘的守時(shí)性能研究

守時(shí)性能通常表現(xiàn)為多臺(tái)原子鐘產(chǎn)生時(shí)間尺度的準(zhǔn)確性、穩(wěn)定性和實(shí)時(shí)性。在時(shí)間保持工作中，時(shí)間尺度相當(dāng)于一臺(tái)虛擬原子鐘，這臺(tái)虛擬原子鐘是通過測(cè)量鐘差實(shí)現(xiàn)的物理鐘的加權(quán)平均。虛擬鐘通過與某一臺(tái)鐘的偏差來表示。一般認(rèn)為虛擬鐘的性能優(yōu)于單臺(tái)原子鐘的性能[5]。下面主要通過時(shí)間尺度算法研究國(guó)產(chǎn)原子鐘的守時(shí)性能，重點(diǎn)研究一種Kalman濾波時(shí)間尺度算法用于國(guó)產(chǎn)鐘的守時(shí)性能，并與一種經(jīng)典的加權(quán)平均算法的守時(shí)性能進(jìn)行比較。

2.1 鐘的模型

Kalman濾波時(shí)間尺度算法是一種實(shí)時(shí)的原子鐘狀態(tài)估計(jì)方法，估計(jì)參考鐘和理想鐘之間的差，并將此值作為修正量，計(jì)算時(shí)間尺度。算法原理如下,原子鐘的狀態(tài)方程可表示為：

(4)

y(tk+1)=y(tk)+δω(tk)+η(tk),

(5)

ω(tk+1)=ω(tk)+α(tk),

(6)

式(4)至(6)中,x(tk)，y(tk)和ω(tk)分別表示鐘在tk時(shí)刻的相位偏差，頻率偏差和頻率漂移偏差。tk和tk+1是相鄰取樣時(shí)間，具有如下關(guān)系，tk+1=tk+δ。ε(tk)，η(tk)和α(tk)分別表示tk時(shí)刻原子鐘的調(diào)頻白噪聲，調(diào)頻隨機(jī)游走噪聲和頻漂隨機(jī)游走噪聲。原子在tk時(shí)刻的測(cè)量模型表示為

zij(tk)=xi(tk)-xj(tk)+v(tk),

(7)

式(7)中，zij表示第ith臺(tái)鐘和第jth臺(tái)鐘之間的相位差測(cè)量，v表示測(cè)量噪聲。

Kalman濾波時(shí)間尺度問題是估計(jì)每臺(tái)原子鐘的狀態(tài)，但是實(shí)際測(cè)量數(shù)據(jù)是鐘差數(shù)據(jù)，因此建立鐘差的狀態(tài)方程估計(jì)原子鐘的狀態(tài)。通過方程(4)，(5)和(6)對(duì)鐘i和鐘j做差的鐘差方程表示為：

(8)

yij(tk+1)=yij(tk)+δωij(tk)+ηij(tk),

(9)

ωij(tk+1)=ωij(tk)+αij(tk)。

(10)

根據(jù)N-1個(gè)鐘差測(cè)量估計(jì)N臺(tái)鐘的狀態(tài)，這種情況的解是不穩(wěn)定的。根據(jù)中心極限定理，當(dāng)原子鐘數(shù)量很多時(shí)，噪聲的加權(quán)和為0。即有下列輔助方程：

(11)

(12)

(13)

式(11)至(13)中,ai(tk)，bi(tk)和ci(tk)分別表示tk時(shí)刻鐘i噪聲的權(quán)重。根據(jù)鐘差方程和輔助方程估計(jì)單臺(tái)原子鐘的狀態(tài)。

2.2 Kalman濾波時(shí)間尺度算法原理

首先定義tk鐘差狀態(tài)變量和噪聲向量為：

(14)

鐘差方程的向量形式表示為

(15)

(16)

式(16)中，H表示測(cè)量矩陣，R是測(cè)量標(biāo)準(zhǔn)偏差。Kalman濾波的遞推方程表示為：

(17)

(18)

(19)

(20)

Kij表示Kalman增益，表示為一個(gè)3-元素列向量：

(21)

式(20)結(jié)合式(14)的估計(jì)噪聲向量為

(22)

2.3 國(guó)產(chǎn)原子鐘的守時(shí)性能研究

研究國(guó)產(chǎn)原子鐘的守時(shí)性能，基于不同時(shí)間尺度算法分析國(guó)產(chǎn)原子鐘的準(zhǔn)確性、穩(wěn)定性和實(shí)時(shí)性。選取國(guó)產(chǎn)原子鐘，編號(hào)為Cs3050，Cs2025和Cs3059計(jì)算時(shí)間尺度，計(jì)算數(shù)據(jù)時(shí)間段為：2020年8月至2021年6月，數(shù)據(jù)采樣周期為1 h，國(guó)產(chǎn)鐘參考UTC(NTSC)的相對(duì)頻率偏差數(shù)據(jù)如圖6所示，為了表示清楚，Cs3059的相對(duì)頻率偏差都減去1×10-12。從圖6可以看出Cs2025的相對(duì)頻率偏差相比Cs3050與Cs3059明顯偏大，Cs3050與Cs3059的相對(duì)頻率偏差相當(dāng)。文中采用上述給出的Kalman濾波算法和經(jīng)典加權(quán)平均(類ALGOS)算法計(jì)算國(guó)產(chǎn)原子鐘的時(shí)間尺度，并分析時(shí)間尺度的準(zhǔn)確度、穩(wěn)定度和實(shí)時(shí)性?；趦煞N算法計(jì)算的時(shí)間尺度如圖7所示，得出本文采用的類ALGOS算法計(jì)算的時(shí)間尺度長(zhǎng)期波動(dòng)較大，基于上述Kalman濾波算法的時(shí)間尺度長(zhǎng)期波動(dòng)相比本文采用的類ALGOS算法明顯減小。分析原因，本文采用的類ALGOS算法采用前一個(gè)月的速率作為本月的預(yù)報(bào)值，如果本月速率相比前一個(gè)月的速率有一個(gè)較大的變化，會(huì)引起時(shí)間尺度較大的波動(dòng)。Kalman濾波算法是一種實(shí)時(shí)算法，給定初始狀態(tài)下根據(jù)測(cè)量鐘差數(shù)據(jù)實(shí)時(shí)修正原子鐘的狀態(tài)估計(jì)，提高了時(shí)間尺度估計(jì)的準(zhǔn)確性。表3給出了基于兩種算法的時(shí)間尺度準(zhǔn)確度指標(biāo)，本文采用的類ALGOS算法的時(shí)間尺度與參考UTC(NTSC)的最大偏差為216.42 ns,明顯大于基于Kalman算法的最大偏差128.42 ns；另外，類ALGOS算法的時(shí)間尺度與參考UTC(NTSC)偏差的平均值為33.19 ns，基于Kalman算法的時(shí)間尺度與參考UTC(NTSC)偏差的平均值為39.40 ns，兩者相差不大；最后，基于本文類ALGOS算法的時(shí)間尺度的RMS為89.17 ns，明顯大于基于Kalman濾波算法的RMS值52.63 ns，說明基于本文的類ALGOS算法的時(shí)間尺度長(zhǎng)期波動(dòng)大?；谝陨戏治?，本文采用的Kalman濾波算法優(yōu)于類ALGOS算法計(jì)算的時(shí)間尺度性能。

圖6 國(guó)產(chǎn)原子鐘的相對(duì)頻率偏差數(shù)據(jù)

圖7 基于兩種算法的國(guó)產(chǎn)鐘時(shí)間尺度

表3 基于兩種算法的時(shí)間尺度的準(zhǔn)確度 單位：ns

分析基于兩種算法的時(shí)間尺度的穩(wěn)定度，穩(wěn)定度通過重疊Allan偏差表示，重疊Allan偏差越小穩(wěn)定度越高,其結(jié)果如圖8所示。從圖8可以得出，不同國(guó)產(chǎn)原子鐘穩(wěn)定度差別明顯，Cs2025和Cs3050穩(wěn)定度性能10天內(nèi)相差不大。10～83 d天內(nèi)，Cs3050的穩(wěn)定度明顯優(yōu)于Cs2025,Cs3059穩(wěn)定度在相應(yīng)取樣時(shí)間差于Cs2025和Cs3050。說明相比較于Cs3050，國(guó)產(chǎn)鐘Cs2025具有較差的長(zhǎng)期穩(wěn)定度，Cs3059短期穩(wěn)定度和長(zhǎng)期穩(wěn)定度都明顯低于其余兩臺(tái)國(guó)產(chǎn)原子鐘?；谖闹蠯alman濾波算法計(jì)算國(guó)產(chǎn)鐘的時(shí)間尺度，算法建立了鐘差的相位偏差，頻率偏差和頻率漂移偏差的狀態(tài)方程，并且狀態(tài)方程中包括鐘差的調(diào)頻白噪聲，調(diào)頻隨機(jī)游走噪聲和頻漂隨機(jī)游走噪聲。Kalman濾波算法通過使線性方差取極小值，同時(shí)降低了鐘差狀態(tài)方程中的多種噪聲，提升了時(shí)間尺度的穩(wěn)定度。圖8給出了本文采用的類ALGOS算法計(jì)算的時(shí)間尺度的重疊Allan偏差和文中Kalman濾波算法計(jì)算的時(shí)間尺度的重疊Allan偏差。本文采用的類ALGOS算法計(jì)算時(shí)間尺度的重疊Allan偏差在相應(yīng)取樣時(shí)間均小于單臺(tái)原子鐘的重疊Allan偏差，30 d的重疊Allan偏差為2.9×10-14，83 d的重疊Allan偏差為1.1×10-14?；贙alman濾波算法計(jì)算時(shí)間尺度在取樣時(shí)間超過5 d后，重疊Allan偏差明顯減小，穩(wěn)定度顯著提高，30 d的重疊Allan偏差為1.0×10-14，83 d的重疊Allan偏差為5.5×10-15。對(duì)于取樣時(shí)間小于5 d，基于Kalman濾波算法計(jì)算時(shí)間尺度的重疊Allan偏差高于類ALGOS算法計(jì)算的重疊Allan偏差，甚至高于單臺(tái)原子鐘的重疊Allan偏差。分析原因，采用文中Kalman濾波時(shí)間尺度算法，首先需要給定鐘差數(shù)據(jù)的初始狀態(tài)，基于Kalman濾波算法估計(jì)預(yù)測(cè)狀態(tài)，根據(jù)量測(cè)數(shù)據(jù)實(shí)時(shí)修正預(yù)測(cè)狀態(tài)，給出鐘差的估計(jì)，這是一個(gè)遞推的過程。當(dāng)原子鐘可預(yù)測(cè)性變差時(shí)，Kalman濾波中量測(cè)數(shù)據(jù)不斷的修正預(yù)測(cè)值，用于提高估計(jì)的準(zhǔn)確性，這樣就降低了時(shí)間尺度的短期穩(wěn)定度。本文采用的類ALGOS算法計(jì)算時(shí)間尺度，根據(jù)前一個(gè)月的速率作為當(dāng)前月的速率預(yù)報(bào)值，速率是恒定的，然后根據(jù)月速率的變化計(jì)算權(quán)重，這樣降低了對(duì)時(shí)間尺度短期穩(wěn)定度的影響。

圖8 國(guó)產(chǎn)原子鐘與時(shí)間尺度的重疊Allan偏差

對(duì)于實(shí)時(shí)性方面，本文采用的類ALGOS算法需要利用前一個(gè)月的速率作為當(dāng)月速率的預(yù)報(bào)值，是一種滯后的時(shí)間尺度算法。Kalman濾波算法根據(jù)預(yù)測(cè)方程和量測(cè)方程遞推的估計(jì)原子鐘的狀態(tài)，是一種實(shí)時(shí)的時(shí)間尺度算法。

3 結(jié)語

本文基于小波分解得到國(guó)產(chǎn)原子鐘和國(guó)外原子鐘的相應(yīng)分解尺度的噪聲，比較不同分解尺度的噪聲大小，分析國(guó)產(chǎn)鐘的性能。守時(shí)過程中，通過對(duì)多臺(tái)國(guó)產(chǎn)鐘在不同分解尺度加權(quán)處理，能夠降低不同分解尺度的噪聲，提升守時(shí)的性能。時(shí)間尺度的準(zhǔn)確度、穩(wěn)定度和實(shí)時(shí)性是守時(shí)性能的一個(gè)重要表現(xiàn)，時(shí)間尺度是由多臺(tái)原子鐘通過算法產(chǎn)生的一個(gè)紙面時(shí)間，國(guó)產(chǎn)鐘守時(shí)的性能與時(shí)間尺度的性能有直接的關(guān)系。本文研究了基于國(guó)產(chǎn)鐘的Kalman濾波時(shí)間尺度算法和本文采用的類ALGOS時(shí)間尺度算法，通過分析得出，Kalman濾波算法計(jì)算得到的時(shí)間尺度的準(zhǔn)確度優(yōu)于類ALGOS算法的準(zhǔn)確度。5 d內(nèi)，類ALGOS算法計(jì)算的時(shí)間尺度的穩(wěn)定度優(yōu)于文中Kalman濾波算法計(jì)算的穩(wěn)定度；5～80 d，Kalman濾波算法計(jì)算的時(shí)間尺度的穩(wěn)定度優(yōu)于類ALGOS算法計(jì)算的穩(wěn)定度。在具體實(shí)施方面，Kalman濾波算法具有較強(qiáng)的實(shí)時(shí)性。