北斗三號衛星共視時間比對性能分析

王威雄，董紹武，武文俊，王 翔，郭 棟，廣 偉

(1. 中國科學院國家授時中心, 西安 710600； 2. 中國科學院時間頻率基準重點實驗室, 西安 710600；3. 中國科學院大學天文與空間科學學院, 北京 101048)

0 引 言

時間是一個國家的重要戰略參數資源，在經濟建設，軍事國防和各類科學研究都有廣泛應用。國際標準時間協調世界時(Coordinated universal time, UTC)是在國際地球自轉服務組織(International earth rotation service, IERS)協助下由國際權度局(Bureau international des poids et mesures, BIPM)利用分布在全球80多個守時實驗室的500多臺原子鐘統一歸算產生[1]，這些原子鐘通過國際時間比對網參與UTC的計算。因此，遠距離時間比對是UTC產生中的重要環節。現國際上主要采用基于全球導航衛星系統(Global navigation satellite system, GNSS)的時間比對和地球同步軌道衛星(Geosy-nchronous earth orbit, GEO)的雙向時間比對(Two-way satellite time and frequency transfer, TWSTFT)兩種手段，但由于TWSTFT設備昂貴且受基線長度限制，只有少數十幾個實驗室才有TWSTFT鏈路，而GNSS時間比對仍是國際上廣泛采用的技術[2]。GPS共視法(Common view, CV)從20世紀80年代就被BIPM用于UTC的計算。為提高時間比對的可靠性、穩定性及準確性，BIPM鼓勵采用多手段冗余時間比對技術[3]。2009年俄羅斯GLONASS時間比對第一次被引入UTC的計算，2017年國際時間頻率咨詢委員會(Consultative committee for time and frequency, CCTF)國際原子時實驗室工作組會議也提出應積極推動北斗參與UTC的計算[4]。

北斗衛星導航系統是我國獨立自主開發的衛星導航系統。近年來，國內外學者針對北斗二號衛星開展了諸多時間比對研究。例如，Huang等[5]在歐洲內部對北斗GEO、傾斜地球同步軌道(Inclined Geosynchronous orbit, IGSO)以及中圓地球軌道(Medium earth orbit, MEO)不同星座的共視時間比對進行了分析[5]；Guang等[6]進一步開展了亞歐長基線北斗共視時間比對試驗；而Andreas等[7]也基于德國物理技術研究院(Physikalisch-technische bundesanstalt, PTB)的國家時間基準系統對北斗信號以及其亞歐共視比對結果進行了有效評估。上述諸多研究表明：北斗二號衛星時間比對性能穩定，但由于歐洲可視衛星數較少，比對性能仍有提升空間。文獻[6-7]中也都提到期望對北斗三號衛星進行更深入的研究。

目前北斗三號衛星正處于全球系統的組網階段，現已發射了20顆組網衛星，能夠向“一帶一路”國家和地區提供基本導航服務。全球系統預計在2020年底建成，屆時整個星座將由3顆GEO衛星、3顆IGSO衛星和24顆MEO衛星組成，為全球提供導航、測速及授時服務。北斗三號衛星配備了高性能的銣原子鐘和氫原子鐘，銣原子鐘天穩定度為E-14量級，氫原子鐘天穩定度為E-15量級，比北斗二號衛星星載鐘的穩定度提高了一個數量級[8-9]。我國標準時間由中科院國家授時中心(National time service center, NTSC)負責產生并對外發播，而北斗時(BDS time, BDT)通過UTC(NTSC)與UTC取得聯系[10]。本文依據CCTF發布的全球衛星導航時間比對標準(Common GNSS generic time transfer standard, CGGTTS)，利用目前NTSC和捷克光電研究院(Institute of photonics and electronics, TP)可觀測的北斗三號衛星對亞歐長基線共視時間比對性能進行了初步計算評估。

1 北斗共視比對算法實現

地面上任意兩個守時實驗室進行北斗共視時間比對是由分別位于兩地的多頻多模接收機在同一時刻同時觀測一顆或多顆北斗衛星，通過偽距測量得到本地時與北斗系統時BDT之間的偏差，該偏差包括衛星、信號傳播路徑和接收機端引入的各項誤差，兩地通過數據交換扣除各項誤差后求差，就可以獲得兩實驗室之間的時差[11-13]，比對原理如圖1所示。

圖1 北斗共視時間比對原理Fig.1 The principle of common view time comparison by BeiDou

設A地的本地時間為tA，B地的本地時間為tB，ΔtABDT是A地與BDT的差，ΔtBBDT是B地與BDT的差，有：

ΔtABDT=tA-tBDT

(1)

ΔtBBDT=tB-tBDT

(2)

式(1)減去式(2)可得兩地的時間差：

ΔtABDT-ΔtBBDT=tA-tBDT-(tB-

tBDT)=tA-tB=ΔtAB

(3)

為標準化數據處理流程，國際時間頻率咨詢委員會在1993年規定了GPS時間比對的格式標準GGTTS V01，1998年該標準又加入了俄羅斯GLONASS衛星導航系統，并命名為CGGTTS V02，隨著歐洲Galileo和中國北斗等衛星導航系統的不斷發展，2015年CCTF的GNSS時間傳遞組再次對時間比對標準進行了擴展，并更名為CGGTTS V2E，該版本與舊版本數據格式兼容，可在不同版本間完成共視時間比對[13]，其中幾項重要的參數信息如表1所示。

表1 CGGTTS標準文件的參數信息Table 1 Parameter information of CGGTTS file

在CGGTTS時間比對標準框架下，其處理算法是連續跟蹤一個或多個衛星并采集13 min的偽碼觀測數據，觀測值1 s一個共780個，每15個測量值作一個二次多項式擬合，然后對擬合結果作各項誤差修正，取其擬合的中間值。13 min的跟蹤周期內有52個擬合中間值，最后將52個擬合中間值進行線性一次擬合，取中間值是本地時間尺度和衛星系統時的偏差[14]。在誤差項處理中，電離層時延和對流層時延分別采用雙頻無電離層組合模型和標準的NATO模型來修正，地球自轉效應使用Sagnac效應算法來修正，設備時延經校準后和測量的線纜時延和參考時延一起寫入文件表頭用于共視時間比對的計算[15-16]。

2 數據平滑濾波與異常值剔除

在共視比對中，由于衛星鐘的不確定性，以及環境干擾及測量噪聲因素的存在，比對結果中會出現一些異常值和隨機噪聲，本文采用3σ準則進行異常值剔除和Vondrak平滑法進行濾波處理。

2.1 Vondrak平滑濾波

Vondrak濾波是由捷克天文學家J.Vondrak在whittaker修勻的基礎上提出的一種濾波方法，其既適用于等間隔的觀測數據又適用于非等間隔的觀測數據，且可以在觀測數據變化規律未知的情況下對觀測數據進行有效的平滑，其基本準則是[17]：

(4)

式中：pi為觀測值的權，yi為觀測值，y′i為待求的平滑值，ε=1/λ是一個給定的無量綱正數，定義為平滑因子。F為平滑值與觀測值的擬合度，S為待求平滑曲線的平滑度。Vondrak濾波實質上就是在觀測值的絕對擬合和絕對平滑之間選擇一條折衷的曲線。ε作為調節折衷程度的參數，ε越小，曲線平滑程度越強，反之，平滑程度越弱[18]。

但在實際應用中，對于相同精度的觀測數據，若采用同一平滑因子對不同采樣間隔和不同長度的數據進行平滑，則得到的平滑曲線具有不同的平滑程度，必須選用不同的平滑因子才能得到同樣平滑程度的曲線，這就給Vondrak的實際運用帶來了不便。為此，Vondrak在1976年對其基本準則進行了改進，將式(4)中的擬合度和平滑度分別用平均值來代替：

(5)

根據式(5)得到的平滑方法，對于相同精度的觀測數據，無論取樣間隔和數據長度如何，選用同一平滑因子就可得到具有同樣平滑程度的結果。

Vondrak濾波的關鍵是平滑因子的選取，目前常用的選取方法有觀測誤差法，頻率響應法和交叉證認法等[19]。本文根據原始共視比對數據的幅值頻譜圖，采用頻率響應法計算恰當的平滑因子作為最后確定的最佳平滑因子。頻率響應函數的解析表達式為：

(6)

式中：f為頻率，ε為平滑因子，A為在不同頻率和平滑因子下的頻率響應函數。

對于濾波結果的評價標準采用互相關系數(R)和濾波后噪聲的均方根誤差(NRMSE)來評估。

1)互相關系數R的計算公式為：

(7)

式中：cov(yi,y′i)是觀測值yi與濾波值y′i的協方差，σy為觀測值yi的標準差，σy′為濾波值y′i的標準差。R越大，說明濾波值與原始觀測值相關性越大，即濾波值與原始值越接近。

2)濾波后噪聲的均方根誤差NRMSE計算公式為：

(8)

式中：yi為觀測值，y′i為濾波值，n為觀測值個數。NRMSE越小，說明濾波值與原始觀測值擬合性越好，測量精度越高。

2.2 異常值剔除

3σ準則是隨機誤差服從正態分布情況下的異常值判別方法。設服從正態分布的一維隨機變量X的概率密度為：

(9)

式中：x∈(-∞,+∞)，μ為X的數學期望，σ2為X的方差。

正態隨機變量X出現在給定區間(μ-kσ，μ-kσ)內的概率(k為正數)為：

(10)

(11)

查概率積分表可得，當k=3時，正態隨機變量X出現在區間μ±3σ范圍內的概率約為99.7%。利用3σ準則剔除共視比對結果中異常值的流程圖如圖2所示。

圖2 3σ準則異常值剔除流程圖Fig.2 Flow chart of 3σ criterion outlier elimination

3 實驗與結果分析

選取2019年3月28日(MJD 58570)2019年4月27日(MJD 58600)間NTSC和TP的多頻多模接收機觀測獲得的BeiDou和GPS CGGTTS數據。該試驗中，每臺接收機的時間頻率參考都為國家時間基準UTC(k)(k指不同的守時實驗室)，參與此次時間比對實驗的接收機配置信息如表2所示。

表2 接收機配置信息Table 2 Receiver configuration information

3.1 零基線共鐘比對

為測試接收機的穩定性和對北斗三號衛星測距碼的噪聲水平進行評估，需要對接收機進行零基線共鐘時間比對(Common clock difference, CCD)實驗。利用共視法對捷克光電研究院的兩臺同類型接收機TP01和TP02開展了北斗三號衛星CCD實驗，并和GPS CCD的結果進行了比較。利用比對結果的標準差(Standard deviation, STD)對噪聲水平進行評估，比對結果如圖3和圖4所示。由于NTSC目前暫未有兩臺可接收北斗三號衛星信號的接收機，所以未進行CCD實驗，但參與本次比對的接收機與TP相同，故認為TP的CCD結果也適用于NTSC。

圖3 北斗三號衛星零基線共鐘比對結果Fig.3 The CCD results of BDS-3 satellite

圖4 GPS零基線共鐘比對結果Fig.4 The CCD results of GPS

由圖3和圖4可以看出，TP01和TP02接收機基于北斗三號衛星的零基線時間比對結果穩定，其STD為0.74 ns，與GPS零基線共鐘比對結果的STD處于同一量級，但由于北斗三號衛星正處于全球組網階段，觀測到的BDS-3衛星較GPS衛星數少，所以BDS-3衛星的測量噪聲水平略差于GPS，但可以用于遠距離的高精度時間比對。

3.2 單站觀測數據分析

利用CGGTTS標準共視文件在比對時間段內對NTSC和TP的觀測數據分別進行分析。兩個實驗室觀測到不同衛星數的觀測量占整個觀測量的百分比統計如圖5～6所示。

圖5 NTSC的可視衛星數統計Fig.5 The visual satellite number statistics of NTSC

圖6 TP的可視衛星數統計Fig.6 The visual satellite number statistics of TP

從圖5～6可以看出，由于北斗二號衛星是區域衛星導航系統，亞太地區覆蓋范圍比歐洲區域好，因此NTSC觀測到的北斗二號衛星數目比TP多，北斗三號衛星是全球衛星導航系統，目前在軌服務的衛星和GPS星座衛星相似都為MEO衛星，所以北斗三號衛星和GPS在亞太和歐洲地區都有相似的覆蓋情況，在NTSC和TP觀測到的北斗三號衛星數量大多都分布在2～5顆，GPS衛星數量大多為5～9顆。

相關研究發現測量信號的碼偽距多路徑噪聲與衛星高度角有關，且高度角越高，多路徑噪聲越小[6]。圖7和圖8分別為觀測期間(MJD 58570-58571)在NTSC和TP不同高度角區間的衛星觀測量統計情況，由于北斗二號衛星的GEO衛星(C01～C05)在NTSC的可觀測數多于TP，因此在頻率分布直方圖中未統計GEO衛星。

圖7 NTSC觀測衛星的高度角統計情況Fig.7 Altitude statistics of NTSC visual satellites

圖8 TP觀測衛星的高度角統計情況Fig.8 Altitude statistics of TP visual satellites

從圖7、8可以得出，北斗二號衛星在NTSC的觀測高度角約有70%都在50°以上，信號覆蓋情況優于TP，數據的觀測質量也更好。而NTSC和TP的北斗三號衛星和GPS衛星在不同高度角區間內的觀測數量在本地和異地都比較相近。

3.3 共視比對結果分析

在GNSS共視時間比對時，較多的可用衛星數可以平均出更好的時間比對結果。圖9給出了NTSC和TP在相同時刻能夠共視到的衛星數。

圖9 NTSC和TP共視衛星數Fig.9 The satellites number of CV between NTSC and TP

從圖9得知，NTSC和TP共視到的GPS衛星數大多為2～3顆，北斗二號衛星數大多為2～4顆，北斗三號衛星數最少，大多為1～2顆。

基于北斗三號衛星，依據GNSS標準共視數據處理規范并利用式(3)得到了幾何基線為7500 km上的NTSC-TP共視時間比對結果。由于衛星或者觀測條件的限制，GNSS時間比對還存在不同程度的粗差。因此，根據圖2中的流程先后對原始計算結果利用頻率響應法進行了Vondrak濾波平滑和3σ準則粗差剔除。由于Vondrak濾波實質是低通濾波，這里選用的頻率響應值A=0.2，原始共視比對結果頻譜圖低頻區域中較高頻率幅值最大的頻率為f，依據式(6)可得平滑因子ε=3493000。對原始比對結果進行濾波平滑和粗差剔除后的三種鏈路的共視比對結果如圖10所示。為便于分析，圖10中對GPS和北斗二號衛星共視結果進行了常數平移。圖11～13給出了三種鏈路時間比對結果濾波前后的幅值頻譜變化情況。

圖10 NTSC-TP鏈路共視比對結果Fig.10 The CV results of NTSC-TP link

表3采用原始結果和濾波結果的互相關系數和濾波后噪聲的均方根誤差對濾波結果進行了評估。

表3 三種共視鏈路濾波結果對比Table 3 Comparison of filtering results of three CV links

從圖10和表3可以看出，GPS共視鏈路、北斗二號衛星共視鏈路和北斗三號衛星共視鏈路的Vondrak平滑值與時間比對計算值有很好的吻合性，可以反映出原始共視比對數據的具體特征。二者具有很強的相關性，濾波值和原始共視結果的相關系數可以達到0.985，0.957和0.968。從圖11～13可知，原始比對結果中含有很強的高頻噪聲，需對噪聲加以抑制或消除。從Vondrak濾波后的頻譜圖可見，經濾波后高頻噪聲得到了很好的抑制。因此，Vondrak濾波可以對共視比對數據進行有效平滑。濾波后三種共視比對鏈路噪聲的均方根誤差分別為0.839 ns，1.417 ns和1.161 ns。由圖9，圖10及表3可知，在目前北斗三號衛星比北斗二號衛星共視可視衛星少的情況下，其共視比對精度較北斗二號衛星有了明顯提高，提高幅度約19%。北斗三號衛星正處于全球組網階段，目前其共視比對精度較GPS略低。

圖14為原始共視結果和濾波值的殘差統計圖及正態擬合曲線。由圖可以看出，三種比對鏈路實測的共視比對結果相對于平滑結果的殘差分布服從正態分布。一方面證明本文通過3σ準則剔除異常值的恰當性，另一方面也說明利用Vondrak濾波主要濾掉的是比對鏈路中的隨機噪聲，平滑后的曲線能夠保留原始比對結果的固有信息，具有很高的保真度。

三種遠程共視比對鏈路的時間方差和修正阿倫方差如圖15、16所示。

從圖15、16可以看出，三種共視鏈路1 d的時間和頻率穩定度都在1 ns量級和10-14量級。當平均時間小于10000 s時，GPS共視比對鏈路的穩定度最高，北斗三號衛星優于北斗二號衛星。當平均時間大于10000 s時，三者的穩定度指標基本相當。

圖11 GPS共視比對結果的頻譜圖Fig.11 The spectrum of GPS CV results

圖12 北斗二號衛星共視比對結果的頻譜圖Fig.12 The spectrum of BDS-2 satellite CV results

圖13 北斗三號衛星共視比對結果的頻譜圖Fig.13 The spectrum of BDS-3 satellite CV results

圖14 原始結果與濾波值殘差分布Fig.14 Residual distribution of original results and filtering results

圖15 不同鏈路的時間方差Fig.15 Time variance of different links

圖16 不同鏈路的修正阿倫方差Fig.16 Modified Allan variance of different links

4 結 論

利用目前中捷兩地可接收的北斗三號衛星完成了NTSC和TP之間亞歐長基線共視時間比對初步試驗與評估，對于北斗系統早日正式納入UTC的歸算提供了技術參考。基于北斗三號衛星完成了零基線共鐘時間比對，利用CGGTTS數據分析了單站可視衛星數及高度角情況，通過Vondrak濾波對共視比對結果進行濾波平滑與粗差剔除，并和北斗二號系統及GPS系統做了比較驗證，并得出以下結論：

1)北斗三號衛星的零基線共鐘時間比對結果連續穩定，其標準差為0.74 ns，可以用于高精度時間比對。

2)Vondrak濾波可以對北斗共視結果進行有效平滑，濾波值具有很高的保真度。

3)在當前中捷之間北斗三號衛星共視可視衛星數比北斗二號衛星共視可視衛星數少一半的情況下，北斗三號衛星的共視比對精度為1.16 ns，北斗二號衛星的共視比對精度為1.42 ns，其比對精度較北斗二號衛星有明顯提升，提升幅度約19%。

4)北斗三號衛星，北斗二號衛星和GPS共視鏈路1 d的時間和頻率穩定度都在1 ns量級和10-14量級。10000 s以內的鏈路穩定度GPS共視最高，北斗三號衛星優于北斗二號衛星。當平均時間大于10000 s時，三者的穩定度指標基本相當。

致 謝

感謝捷克科學院光電研究院提供的北斗和GPS RINE以及CGGTTS觀測文件。