李變，屈俐俐，高玉平，3

1．中國科學(xué)院 國家授時(shí)中心，西安 710600 2．中國科學(xué)院 時(shí)間頻率基準(zhǔn)重點(diǎn)實(shí)驗(yàn)室，西安 710600 3．中國科學(xué)院大學(xué)，北京 100049

目前，精度最高的時(shí)間尺度國際原子時(shí)（TAI）、協(xié)調(diào)世界時(shí)（UTC）和地球時(shí)的另一個(gè)實(shí)現(xiàn)TT（BIPMXX）都是由國際權(quán)度局（BIPM）基于原子鐘建立并保持的。BIPM首先收集分布在全球80余個(gè)守時(shí)實(shí)驗(yàn)室約500臺(tái)原子鐘的數(shù)據(jù)資料，采用ALGOS算法加權(quán)平均算法得到自由原子時(shí)（EAL），然后利用10余臺(tái)基準(zhǔn)頻標(biāo)對(duì)EAL進(jìn)行頻率校準(zhǔn)后獲得TAI。UTC是TAI與世界時(shí)相互協(xié)調(diào)的產(chǎn)物，與TAI只差整數(shù)秒，截止目前二者相差37 s［1］。

參與EAL計(jì)算的原子鐘主要有銫鐘和氫鐘2種類型，分別約占總鐘數(shù)的59%和40%，TT（BIPMXX）是BIPM每年年初利用所有可用的基準(zhǔn)頻標(biāo)和原子鐘的數(shù)據(jù)資料，采用事后處理的方式獲得的滯后一年的時(shí)間尺度，也是目前全球最準(zhǔn)確、最穩(wěn)定的原子時(shí)尺度［2］。然而TT（BIPMXX）、TAI和UTC是滯后的紙面時(shí)間，無法滿足時(shí)間/頻率作為參考實(shí)時(shí)性的要求。因此，各個(gè)國家或地區(qū)的守時(shí)實(shí)驗(yàn)室都建立與保持著具有實(shí)時(shí)、連續(xù)、穩(wěn)定的物理信號(hào)UTC（k），作為UTC的物理實(shí)現(xiàn)用于實(shí)時(shí)測(cè)量比對(duì)，k為實(shí)驗(yàn)室名稱縮寫。

毫秒脈沖星通常被認(rèn)為是由雙星系統(tǒng)演化而來的，雙星系統(tǒng)中的正常脈沖星通過吸積伴星質(zhì)量來增加角動(dòng)量，將脈沖星自轉(zhuǎn)周期加速到毫秒量級(jí)。與年輕正常脈沖星相比，毫秒脈沖星具有輻射流量弱、自轉(zhuǎn)穩(wěn)定度高、脈沖形狀陡、計(jì)時(shí)觀測(cè)精度高、且很少發(fā)生自轉(zhuǎn)頻率的突變（glitch）等特征，有些自轉(zhuǎn)周期變化率可達(dá)到10?19～10?21，被譽(yù)為自然界最穩(wěn)定的“時(shí)鐘”。1982年Backer等［3］發(fā) 現(xiàn) 第 一 顆 毫 秒 脈 沖 星PSR B1937+21之后，人們便開始了毫秒脈沖星的計(jì)時(shí)觀測(cè)，關(guān)于脈沖星計(jì)時(shí)資料分析及脈沖星時(shí)算法的研究也隨之興起［4-7］，2012年澳大利亞Parkes天文臺(tái)的Hobbs等得到了可以與地球時(shí)TT（BIPMXX）相媲美的脈沖星時(shí)［8-9］。

受觀測(cè)設(shè)備、觀測(cè)技術(shù)、星際介質(zhì)傳播延遲等因素的影響，目前毫秒脈沖星最好的計(jì)時(shí)精度僅為100 ns，遠(yuǎn)低于原子鐘測(cè)量比對(duì)精度。但原子時(shí)的長期頻率穩(wěn)定度受原子鐘性能、季節(jié)性變化、周年效應(yīng)、以及原子時(shí)算法等因素的影響，長期穩(wěn)定度不如脈沖星時(shí)。由于脈沖星時(shí)沒有獨(dú)立定義時(shí)刻起點(diǎn)，因此脈沖星只能被視為長期穩(wěn)定度極高的頻率源，用于其他時(shí)頻設(shè)備（例如：原子鐘）輸出頻率長期性能的驗(yàn)證或校準(zhǔn)。

近年來隨著脈沖星觀測(cè)設(shè)備性能提升，觀測(cè)到的毫秒脈沖星的數(shù)量有所增加、觀測(cè)精度不斷提高。尤其是我國FAST的投入運(yùn)行，將積累大量的、優(yōu)質(zhì)的毫秒脈沖星觀測(cè)資料［10］，同時(shí)中國SKA先導(dǎo)項(xiàng)目已順利開展，這些為改變目前時(shí)頻體系單一依賴原子鐘的格局，進(jìn)而提高當(dāng)前原子時(shí)體系的可靠性和穩(wěn)定性提供有力的數(shù)據(jù)及技術(shù)支持。

脈沖星時(shí)與原子時(shí)具有優(yōu)勢(shì)互補(bǔ)的特點(diǎn)，將脈沖星時(shí)用于原子時(shí)的檢測(cè)、校準(zhǔn)與駕馭，利用二者的優(yōu)勢(shì)，可以建立一種能夠在保持原子時(shí)短期穩(wěn)定度優(yōu)勢(shì)的基礎(chǔ)上，又具有脈沖星時(shí)長期穩(wěn)定度的優(yōu)勢(shì)綜合時(shí)間尺度。

1 駕馭銫鐘的參考毫秒脈沖星選取

國際脈沖星計(jì)時(shí)陣（IPTA）2019年發(fā)布了第二批脈沖星計(jì)時(shí)觀測(cè)數(shù)據(jù)（IPTA dr2）。這些觀測(cè)數(shù)據(jù)分別來自歐洲脈沖星計(jì)時(shí)陣（EPTA）、北美脈沖星計(jì)時(shí)陣（NANOGrav）和澳大利亞Parkes脈沖星計(jì)時(shí)陣（PPTA）。與IPTA發(fā)布的第一批數(shù)據(jù)（IPTA dr1）相比，IPTA dr2中包含65顆毫秒脈沖星，比IPTA dr1多了16顆［11-12］。

1.1 IPTA dr2概況

IPTA dr2中共有65顆毫秒脈沖星，其中EPTA包含42顆，觀測(cè)時(shí)間跨度7～18 a，NANO?Grav包含37顆，觀測(cè)時(shí)間跨度0.6～9.2 a，PPTA包含20顆，觀測(cè)時(shí)間跨度5.6～29.4 a。

與IPTA dr1不同，IPTA dr2分別采用脈沖星計(jì)時(shí)軟件TEMPO2和TEMPOnest擬合計(jì)時(shí)模型得到A、B兩個(gè)版本的脈沖星計(jì)時(shí)殘差。IPTA dr2具有以下特點(diǎn)：① 對(duì)于同一顆脈沖星多個(gè)觀測(cè)帶寬的情況，根據(jù)脈沖到達(dá)時(shí)間（TOA）的不確定度計(jì)算權(quán)重，并將權(quán)重最大的數(shù)據(jù)序列作為計(jì)時(shí)模型的參考，進(jìn)而得到其他時(shí)間序列相對(duì)于該參考的系統(tǒng)差；② 處理方法更全面，計(jì)時(shí)殘差精度更高。例如：對(duì)低偏心率脈沖雙星，采用ELL1模型擬合拉普拉斯-拉格朗日參數(shù)。對(duì)低偏心率中-高傾角的脈沖雙星，采用DDH模型對(duì)軌道周期的3次諧波振幅和連續(xù)諧波的振幅比進(jìn)行擬合等；③ 計(jì)時(shí)觀測(cè)的時(shí)間跨度相差大。毫秒脈沖星J1939+2134的觀測(cè)間隔最長，約29.4 a，J0931-1902的觀測(cè)間隔最短，只有0.6 a；④ 計(jì)時(shí)殘差數(shù)據(jù)質(zhì)量差異大，不僅表現(xiàn)在不同的脈沖星之間，而且同一顆脈沖星不同觀測(cè)帶寬之間的數(shù)據(jù)質(zhì)量差異也非常明顯。

1.2 參考毫秒脈沖星選取

用于銫鐘駕馭的脈沖星參考必須滿足測(cè)量精度高、觀測(cè)間隔盡可能長、長期穩(wěn)定度高，并且計(jì)時(shí)觀測(cè)數(shù)據(jù)的缺失盡可能少等條件。為此，首先在IPTA dr2中選取滿足觀測(cè)間隔τ>15 a、TOA測(cè)量精度σTOA<1.5 μs條件，并且計(jì)時(shí)觀測(cè)數(shù)據(jù)缺失盡可能少的6顆毫秒脈沖星作為候選參考（圖1）。

圖1是對(duì)初步選取的編號(hào)分別為J0437-4715、J1640+2224、J1713+0747、J0613-0200、J1744-1134和J1939+2134的6顆毫秒脈沖星不同觀測(cè)波段計(jì)時(shí)觀測(cè)數(shù)據(jù)，采用TEMPO2軟件處理后的計(jì)時(shí)殘差。圖中橫坐標(biāo)是約化儒略日（MJD）表示的時(shí)間序列，縱坐標(biāo)計(jì)時(shí)殘差。

圖1 6顆毫秒脈沖星的計(jì)時(shí)殘差Fig. 1 Timing residuals of 6 millisecond pulsars

利用毫秒脈沖星駕馭銫鐘的宗旨是在保持原子鐘短期穩(wěn)定度的基礎(chǔ)上，發(fā)揮脈沖星長期穩(wěn)定度的優(yōu)勢(shì)。因此，穩(wěn)定度是必須考慮的關(guān)鍵因素。脈沖星數(shù)據(jù)采樣具有不規(guī)則性，其穩(wěn)定度估計(jì)采用與Allan方差類似的σz(τ)估計(jì)方法［13］：

式中：τ為采樣間隔；c為τ間隔上多項(xiàng)式擬合的三次項(xiàng)系數(shù)； 表示統(tǒng)計(jì)平均。

對(duì)選取的6顆候選毫秒脈沖星，采用σz(τ)分別計(jì)算穩(wěn)定度，結(jié)果如圖2所示。

圖2 6顆毫秒脈沖星的穩(wěn)定度Fig. 2 Stabilities of 6 millisecond pulsars

計(jì)時(shí)噪聲也是決定毫秒脈沖星參考選取的一個(gè)重要考慮因素。計(jì)時(shí)噪聲由低頻信號(hào)組成，是計(jì)時(shí)殘差中不可預(yù)報(bào)的長期變化趨勢(shì)［14-17］。表1是6顆毫秒脈沖星的觀測(cè)間隔、TOA測(cè)量精度、計(jì)時(shí)殘差精度及是否含有計(jì)時(shí)噪聲等情況。

表1 6顆毫秒脈沖星比較Table 1 Comparison of 6 millisecond pulsars

根據(jù)表1、圖1及圖2中6顆候選毫秒脈沖星的穩(wěn)定度，選取J0613-0200作為銫鐘駕馭的參考，開展毫秒脈沖星駕馭銫鐘的實(shí)驗(yàn)研究。

2 J0613-0200駕馭銫鐘實(shí)驗(yàn)

2.1 實(shí)驗(yàn)系統(tǒng)

利用中國科學(xué)院國家授時(shí)中心（NTSC）的一臺(tái)5071A銫鐘（編號(hào)：Cs2928）、相位微調(diào)儀、計(jì)數(shù)器和相位頻率分配放大器等相關(guān)時(shí)頻設(shè)備以及IPTA dr2的計(jì)時(shí)觀測(cè)資料，建立單顆毫秒脈沖星駕馭銫鐘實(shí)驗(yàn)系統(tǒng)（圖3）。

圖3中：f為10 MHz頻率信號(hào)；p為1 PPS信號(hào)。時(shí)間基準(zhǔn)為NTSC建立和保持的中國標(biāo)準(zhǔn)時(shí)間UTC（NTSC），作為計(jì)數(shù)器的開門信號(hào)。計(jì)數(shù)器測(cè)量Cs2928及其被駕馭后輸出的實(shí)驗(yàn)系統(tǒng)時(shí)間（APTCs）與UTC（NTSC）的相位差，用于實(shí)驗(yàn)數(shù)據(jù)的分析處理。TEMPO2軟件用于分析處理IPTA dr2中的J0613-0200計(jì)時(shí)觀測(cè)資料，并基于J0613-0200的參數(shù)文件，定期模擬用于Cs2928駕馭的計(jì)時(shí)觀測(cè)資料。

圖3 駕馭實(shí)驗(yàn)系統(tǒng)Fig. 3 Experimental system for steering

2.2 J0613-0200模擬數(shù)據(jù)

IPTA dr2發(fā)布的J0613-0200計(jì)時(shí)觀測(cè)數(shù)據(jù)截止到2014年5月18日（MJD=56795），為了開展毫秒脈沖星駕馭銫鐘的實(shí)驗(yàn)研究，根據(jù)J0613-0200的自轉(zhuǎn)參數(shù)、天體測(cè)量參數(shù)及計(jì)時(shí)殘差的統(tǒng)計(jì)特性，利用TEMPO2軟件中的fake插件，每月定期模擬用于Cs2928駕馭的J0613-0200計(jì)時(shí)觀測(cè)資料。圖4是利用J0613-0200已有觀測(cè)資料和模擬數(shù)據(jù)計(jì)算得到的計(jì)時(shí)殘差。圖中橫坐標(biāo)是MJD表示的時(shí)間序列，縱坐標(biāo)是J0613-0200的計(jì)時(shí)殘差。圖中黑色曲線（Obs）是根據(jù)IPTA dr2中J0613-0200（MJD=50931～56795）的計(jì)時(shí)觀測(cè)數(shù)據(jù)，得到的計(jì)時(shí)殘差，灰色曲線（Sim）是根據(jù)J0613-0200模擬數(shù)據(jù)得到的計(jì)時(shí)殘差。

圖4 J0613-0200計(jì)時(shí)殘差Fig. 4 Timing residuals of J0613-0200

J0613-0200模擬數(shù)據(jù)中計(jì)時(shí)紅噪聲的功率譜密度P(f)為［18］

式中：P0為振幅；fc為截止頻率；α為譜指數(shù)。

2.3 實(shí)驗(yàn)系統(tǒng)銫鐘性能分析

守時(shí)原子鐘組通常由銫鐘和氫鐘組成，銫鐘擁有優(yōu)良的長期穩(wěn)定度（銫鐘的閃爍底噪1×10?14～1×10?15，典型值5×10?15），氫鐘 擁有出色的短期穩(wěn)定度（典型氫鐘的閃爍底噪1×10?15～1×10?16，典型值5×10?16），但因存在頻率漂移，其長期穩(wěn)定度取決于頻率漂移的變化與大小。本實(shí)驗(yàn)的目的是證明毫秒脈沖星駕馭原子鐘的可行性，因此，選擇沒有頻率漂移或頻率漂移很小的銫鐘作為被駕馭的原子鐘。

Cs2928是NTSC一臺(tái)連續(xù)工作近8年的銫鐘，該鐘的頻率穩(wěn)定性略高于同類銫鐘平均水平，并且擁有完整的歷史測(cè)量資料，噪聲特性長期持續(xù)穩(wěn)定。圖5是Cs2928的穩(wěn)定度曲線。由圖5可知：Cs2928的噪聲主要包括白色調(diào)頻噪聲（WFM）、頻率隨機(jī)游走調(diào)頻噪聲（RWFM）和閃爍調(diào)頻噪聲（FFM）。Cs2928的頻率漂移量值很小（<1×10?17），被淹沒在隨機(jī)游走頻率噪聲中。因此，Cs2928的駕馭實(shí)驗(yàn)不考慮頻率漂移的影響，只需根據(jù)歷史測(cè)量資料，結(jié)合其噪聲特點(diǎn)，準(zhǔn)確預(yù)報(bào)其頻率變化。

2.4 駕馭算法

實(shí)驗(yàn)系統(tǒng)旨在產(chǎn)生和保持既具有Cs2928短期穩(wěn)定度，同時(shí)又具有J0613-0200長期穩(wěn)定度特性的時(shí)間/頻率信號(hào)。因此，銫鐘駕馭實(shí)驗(yàn)系統(tǒng)采用雙重駕馭的算法。首先，基于Cs2928的歷史數(shù)據(jù)特征及噪聲特點(diǎn)，以地方原子時(shí)TA（NTSC）為參考，計(jì)算用于每天駕馭的日駕馭量（offset1）。由2.3節(jié)內(nèi)容可知，Cs2928的頻率漂移量值很小，被隨機(jī)游走噪聲淹沒。因此，每天的日駕馭量僅考慮Cs2928的頻率估計(jì)；其次，每月根據(jù)fake插件模擬的J0613-0200計(jì)時(shí)觀測(cè)資料，分析計(jì)算用于每月駕馭的月駕馭量（offset2）。

對(duì)于月駕馭，需要評(píng)估J0613-0200和實(shí)驗(yàn)系統(tǒng)產(chǎn)生的時(shí)間APTCs之間的偏差。以TAI為共同參考，評(píng)估J0613-0200與APTCs之間的偏差：

采用卡爾曼濾波方法，估計(jì)濾波后的相位和頻率［19-20］，即

式中：B0和B1分別為濾波后的相位和頻率；H為測(cè)量矩陣。

設(shè)X為狀態(tài)向量；Φ為狀態(tài)轉(zhuǎn)移矩陣；Zk為測(cè)量輸入數(shù)據(jù)(J0613-0200)?APTCs；P為狀態(tài)向量的協(xié)方差矩陣；Q為濾波器的過程噪聲矩陣；R為測(cè)量噪聲矩陣，則

式中：

濾波器初始值由前一計(jì)算間隔的(J0613-0200)?APTCs計(jì)算得到，卡爾曼濾波估計(jì)的頻率即為用于每月駕馭的月駕馭量offset2：

這樣，每天的駕馭量是這2項(xiàng)之和：

3 實(shí)驗(yàn)結(jié)果

J0613-0200駕馭Cs2928的日駕馭量每天定時(shí)送入相位微調(diào)儀對(duì)Cs2928進(jìn)行初次駕馭；月駕馭量每月根據(jù)J0613-0200的模擬數(shù)據(jù)計(jì)算，對(duì)Cs2928進(jìn)行二次駕馭，輸出實(shí)驗(yàn)系統(tǒng)時(shí)間APTCs。利用時(shí)間間隔計(jì)數(shù)器將APTCs時(shí)間信號(hào)與NTSC保持的我國標(biāo)準(zhǔn)時(shí)間UTC（NTSC）進(jìn)行測(cè)量比對(duì)，用于實(shí)驗(yàn)系統(tǒng)時(shí)間APTCs的評(píng)估。

作為國際標(biāo)準(zhǔn)時(shí)間的UTC由TAI閏秒后得到，通常BIPM滯后15～45 d發(fā)布，采樣間隔為5 d?？焖賲f(xié)調(diào)世界時(shí)（UTCr）是BIPM于2012年開始，在40多個(gè)守時(shí)實(shí)驗(yàn)室和國際時(shí)間頻率咨詢委員會(huì)（CCTF）的大力支持下，正式啟動(dòng)的一項(xiàng)研究。近年來，UTCr和UTC的一致性很好，最大偏差優(yōu)于2 ns，每周三BIPM公布最新結(jié)果，并且UTCr的采樣間隔為1 d。因此，實(shí)驗(yàn)系統(tǒng)時(shí)間APTCs的評(píng)估以BIPM每周公布的快速協(xié)調(diào)世界時(shí)（UTCr）為參考。

銫鐘駕馭實(shí)驗(yàn)從2021年2月下旬開始，圖6為MJD=59262～59323的Cs2928未經(jīng)駕馭的運(yùn)行結(jié)果相對(duì)于UTCr的相位偏差UTCr－Cs2928和Cs2928經(jīng)駕馭后的實(shí)驗(yàn)結(jié)果APTCs相對(duì)于UTCr的相位偏差UTCr－APTCs。

圖6 J0613-0200未經(jīng)駕馭和駕馭Cs2928實(shí)驗(yàn)結(jié)果Fig. 6 Experimental results of J0613-0200 Unsteering and steering Cs2928

由于APTCs的相位偏差是任意的，因此為了便于分析，UTCr?APTCs曲線在起點(diǎn)MJD=59262處被人為對(duì)齊到零。由圖6（b）可知，APTCs與UTCr之間存在1.39×10?15的頻偏。主要原因是：①J0613-0200與UTC/TAI之間存在2.82×10?15的頻率偏差；②J0613-0200的模擬數(shù)據(jù)與觀測(cè)數(shù)據(jù)存在1.27×10?15的頻偏；③Cs2928的一次駕馭采用歷史測(cè)量數(shù)據(jù)；④模型估計(jì)的不確定性也是一個(gè)影響因素。因此，APTCs與UTCr之間存在該量值的頻偏是合理的，結(jié)果符合毫秒脈沖星駕馭銫鐘的預(yù)期，證明脈沖星駕馭銫鐘實(shí)現(xiàn)脈沖星時(shí)間尺度頻率的可行性。

APTCs與UTCr之間固定的頻率偏差，并不影響2個(gè)時(shí)間尺度的相對(duì)穩(wěn)定性。截止目前，實(shí)驗(yàn)系統(tǒng)運(yùn)行7個(gè)多月，APTCs相對(duì)于UTCr的頻率穩(wěn)定度曲線如圖7所示。

圖7 UTCr與APTCs之間固定的頻率穩(wěn)定度Fig. 7 Frequency stabilities between UTCr and APTCs

4 結(jié)論

脈沖星時(shí)與原子時(shí)具有優(yōu)勢(shì)互補(bǔ)的特點(diǎn)，利用NTSC的原子鐘資源、IPTA dr2和相關(guān)時(shí)頻設(shè)備建立了毫秒脈沖星駕馭銫鐘實(shí)驗(yàn)系統(tǒng)，開展單顆毫秒脈沖星駕馭銫原子鐘的實(shí)驗(yàn)研究。

基于第3節(jié)的實(shí)驗(yàn)結(jié)果，可得出以下結(jié)論：

1）APTCs與UTCr之間的頻率偏差符合單顆毫秒脈沖星駕馭銫6s鐘的預(yù)期結(jié)果，證明了脈沖星駕馭原子鐘的可行性。如果利用性能優(yōu)秀的多顆毫秒脈沖星建立綜合脈沖星時(shí)，將有效克服因單顆毫秒脈沖星引起的頻率偏差，取得更好的結(jié)果。

2）本駕馭實(shí)驗(yàn)采用的數(shù)據(jù)資料是J0613-0200的模擬數(shù)據(jù)和NTSC的原子鐘測(cè)量數(shù)據(jù)，并未用BIPM每月發(fā)布的CirT公報(bào)。在消除APTCs相對(duì)于UTCr的固定頻偏后，可以達(dá)到國際電聯(lián)ITU對(duì)守時(shí)實(shí)驗(yàn)室±100 ns的要求，即|UTCr－APTCs|<100 ns。

3）實(shí)驗(yàn)結(jié)果出現(xiàn)較大的相位波動(dòng)，主要是因?yàn)槭蹸s2928自身原因或工作環(huán)境等因素的影響，使得用于一次駕馭的預(yù)報(bào)頻率與其實(shí)際頻率不完全符合，從而出現(xiàn)較大的頻率波動(dòng)，說明原子鐘是影響APTCs的關(guān)鍵因素之一。因此，需要從準(zhǔn)確度、穩(wěn)定度、頻率復(fù)現(xiàn)等方面綜合考慮，選擇性能優(yōu)秀的原子鐘作為被駕馭的頻率源。