孫 廣, 朱翔宇, 李雙欽, 郭美軍, 翟 偉, 洪英杰

(西安航天天繪數(shù)據(jù)技術(shù)有限公司，西安 710100)

0 引言

目前，衛(wèi)星鐘差以及軌道的精度，在一定程度上決定了全球?qū)Ш叫l(wèi)星系統(tǒng)(Global Navigation Satellite System，GNSS)的服務(wù)精度，當(dāng)軌道平穩(wěn)變化時(shí)，軌道預(yù)報(bào)精度可支持實(shí)時(shí)精密定位，但衛(wèi)星的鐘差精度較低。此外，星載原子鐘作為星上的時(shí)間基準(zhǔn)，其性能很大程度決定了導(dǎo)航、定位和授時(shí)的精度[1-2]。大量時(shí)間累積的衛(wèi)星鐘差數(shù)據(jù)，可以分析評(píng)估星載原子鐘的健康狀態(tài)。GNSS的完好性[3]、系統(tǒng)性能評(píng)估[4]和GNSS衛(wèi)星鐘差計(jì)算[5]及預(yù)報(bào)[6-7]等領(lǐng)域都是當(dāng)前研究的熱點(diǎn)。同時(shí)，很多專業(yè)人員對(duì)GNSS星載原子鐘周期項(xiàng)計(jì)算以及鐘差改正評(píng)估都有一定的研究[2,7]。針對(duì)中國北斗衛(wèi)星導(dǎo)航系統(tǒng)(BeiDou Navigation Satellite System，BDS)原子鐘周期的研究成果較少。

我國的BDS，將在2020年實(shí)現(xiàn)全球覆蓋；而俄羅斯的GLONASS已經(jīng)是全球覆蓋，歐洲的Galileo系統(tǒng)也正在逐步發(fā)展。針對(duì)BDS星載原子鐘的性能分析[8]采用修正阿倫方差計(jì)算，對(duì)BDS的傾斜地球同步軌道(Inclined Geosynchronous Satellite Orbit，IGSO)和地球靜止軌道(Geostationary Earth Orbit，GEO)星載原子鐘的穩(wěn)定性進(jìn)行計(jì)算分析[9]；有研究人員對(duì)BDS衛(wèi)星鐘差數(shù)據(jù)的周期進(jìn)行了分析[10]。通常鐘差序列包含線性、二次曲線以及周期性變化規(guī)律。研究人員基本利用太陽軌道和衛(wèi)星軌道2個(gè)軌道面間的夾角，加入地影期相關(guān)系數(shù)的變化對(duì)鐘差的周期性進(jìn)行改正分析[11]。大量的周期性研究主要集中在全球定位系統(tǒng)(Global Positioning System，GPS)，其他系統(tǒng)鐘性能的研究較少?；诖?，本文通過2016年的全年多星定軌聯(lián)合解算的衛(wèi)星鐘差數(shù)據(jù)，研究分析了BDS、GPS鐘差周期特性指標(biāo)，為后續(xù)鐘差模型構(gòu)建提供參考。

1 算法原理

1.1 鐘差數(shù)據(jù)預(yù)處理

星載原子鐘長(zhǎng)時(shí)間在軌運(yùn)行，累計(jì)產(chǎn)生誤差，很多不確定的因素會(huì)對(duì)鐘差數(shù)據(jù)產(chǎn)生影響，所以鐘差數(shù)據(jù)中經(jīng)常會(huì)有大的粗差數(shù)據(jù)、缺失以及異常數(shù)據(jù)。因此，需要對(duì)衛(wèi)星鐘差數(shù)據(jù)進(jìn)行去粗差及異常處理。對(duì)于大量時(shí)間累積的衛(wèi)星鐘差數(shù)據(jù)相應(yīng)的頻率數(shù)據(jù)，本文采用抗差較好的中位數(shù)法(Median Absolute Deviation，MAD)[12-13]進(jìn)行粗差去除，如式(1)所示。

|yi|>(m+n·MAD)

(1)

利用式(1)進(jìn)行粗差判斷，并提出，式中，m=Median(yi),MAD= Median{|yi-m|/0.6745}。MAD算法的數(shù)據(jù)源是頻率數(shù)據(jù)，用鐘差數(shù)據(jù)做差分計(jì)算頻率數(shù)據(jù)，然后利用該方法進(jìn)行粗差去除。對(duì)于去除數(shù)據(jù)的點(diǎn)，一般進(jìn)行內(nèi)插替換該點(diǎn)數(shù)據(jù)或取0。上述處理方法會(huì)在原始數(shù)據(jù)中加入新數(shù)據(jù)。為了避免引入新數(shù)據(jù)，本文對(duì)異常鐘差數(shù)據(jù)設(shè)為空，表示該歷元鐘差數(shù)據(jù)缺失。

1.2 星載原子鐘鐘差模型

星載原子鐘在長(zhǎng)時(shí)間運(yùn)行會(huì)發(fā)生一定的頻率漂移，所以，構(gòu)建鐘差模型的參數(shù)包含衛(wèi)星鐘時(shí)頻特性的相位、頻率、頻率漂移的二次多項(xiàng)式模型。構(gòu)建的模型如式(2)。

Δti=a0+a1(ti-t0)+a2(ti-t0)2+Δi

(2)

式中，i=1,2,…,n，Δti表示第i歷元的衛(wèi)星鐘差；t0表示星鐘的參考時(shí)刻；ti是對(duì)應(yīng)的歷元時(shí)刻；需要估計(jì)的參數(shù)a0、a1和a2，對(duì)應(yīng)為t0的相位、頻率及頻率漂移率；Δi表示原子鐘觀測(cè)誤差。需要鐘差數(shù)據(jù)大于3個(gè)，可利用式(2)進(jìn)行估計(jì)參數(shù)計(jì)算。

1.3 星載原子鐘的周期項(xiàng)計(jì)算

本文利用傅里葉變換對(duì)擬合殘差序列進(jìn)行頻譜分析。在進(jìn)行衛(wèi)星鐘差周期項(xiàng)計(jì)算時(shí)，需要利用二次多項(xiàng)式對(duì)衛(wèi)星鐘差進(jìn)行擬合，可消除鐘差趨勢(shì)項(xiàng)，即可獲得擬合殘差。然后利用擬合殘差進(jìn)行周期項(xiàng)的計(jì)算。離散傅里葉變換(Discrete Fourier Transform，DFT)，其表達(dá)式為

(3)

式中，X(k)表示k時(shí)段對(duì)應(yīng)的頻譜值；X(n)表示相應(yīng)的擬合殘差序列；n表示殘差序列中的序號(hào)；N是殘差的個(gè)數(shù)，利用該式計(jì)算時(shí)，需要滿足N=2L，L=0,1,2,…,n的條件，若N不滿足要求，則進(jìn)行補(bǔ)0；然后即可求出殘差序列中各點(diǎn)的頻譜值。在實(shí)際仿真計(jì)算過程中，采用快速傅里葉變換 (Fast Fourier Transform，F(xiàn)FT)進(jìn)行求解。然后，根據(jù)殘差序列的頻譜圖即可進(jìn)行周期項(xiàng)提取。

1.4 周期項(xiàng)特性分析

首先，利用多星聯(lián)合解算的GNSS精密衛(wèi)星鐘差數(shù)據(jù)進(jìn)行預(yù)處理消除粗差，然后進(jìn)行二次多項(xiàng)式擬合，再通過擬合數(shù)據(jù)和原始數(shù)據(jù)做差得到擬合殘差數(shù)據(jù)，對(duì)該數(shù)據(jù)進(jìn)行頻譜分析可以獲得GNSS星載原子鐘的周期項(xiàng)。圖1給出了本文分析GNSS星載原子鐘周期特性的整個(gè)流程。

圖1 GNSS原子鐘周期性分析流程Fig.1 Flow of GNSS atomic clock periodic analysis

2 實(shí)驗(yàn)結(jié)果分析

本文采用2016年全年的多星定軌聯(lián)合解算的BDS、GPS衛(wèi)星鐘差數(shù)據(jù)，進(jìn)行數(shù)據(jù)預(yù)處理粗差剔除，采用鐘差二次多項(xiàng)式模型對(duì)衛(wèi)星鐘差數(shù)據(jù)進(jìn)行逐天擬合，最終得到相應(yīng)的擬合殘差。利用該擬合殘差，對(duì)BDS、GPS的周期特性進(jìn)行分析。實(shí)驗(yàn)統(tǒng)計(jì)結(jié)果中C代表BDS，G代表GPS。

2.1 BDS原子鐘周期項(xiàng)頻譜分析

利用BDS的擬合殘差數(shù)據(jù)對(duì)BDS原子鐘周期特性進(jìn)行了實(shí)驗(yàn)分析。對(duì)BDS所有的衛(wèi)星進(jìn)行了頻譜分析，本文以中地球軌道(Medium Earth Orbit，MEO)、GEO、IGSO三種衛(wèi)星鐘的一種為例，進(jìn)行頻譜分析。頻譜分析結(jié)果如圖2～圖4所示。

從圖2～圖4分析可知，BDS的GEO衛(wèi)星C02星的主周期項(xiàng)在24h、12h以及8h附近時(shí)，周期性明顯；圖3中，IGSO衛(wèi)星C07的周期項(xiàng)也分布明顯，基本與C02主周期項(xiàng)相同，24h的周期與GEO、IGSO衛(wèi)星軌道周期23h56min基本符合；而MEO類型的C11星的主周期項(xiàng)在12.8h左右，也與MEO衛(wèi)星軌道周期12h53min基本一致。BDS所有衛(wèi)星鐘主周期項(xiàng)分析統(tǒng)計(jì)結(jié)果如表1所示。

圖2 BDS-C02星的鐘差數(shù)據(jù)頻譜分析結(jié)果 Fig.2 Spectral analysis based on BDS-C02 satellite clock bias

圖3 BDS-C07星的鐘差數(shù)據(jù)頻譜分析結(jié)果Fig.3 Spectral analysis based on BDS-C07 satellite clock bias

圖4 BDS-C11星的鐘差數(shù)據(jù)頻譜分析結(jié)果Fig.4 Spectral analysis based on BDS-C11 satellite clock bias

h

從表1中可以看出，GEO、IGSO衛(wèi)星主周期項(xiàng)基本為12h或24h；MEO主周期項(xiàng)為12.8h或6.4h，所有衛(wèi)星的主周期項(xiàng)基本為軌道周期的1倍或1/2倍。

表1中，GEO、IGSO鐘差均含有24h主周期項(xiàng)，而MEO衛(wèi)星鐘差的主周期項(xiàng)則包含相對(duì)較長(zhǎng)的周期項(xiàng)20.2h、28.1h，這些周期項(xiàng)基本接近24h，可以推斷，可能是由于光壓模型攝動(dòng)力對(duì)MEO衛(wèi)星的周期項(xiàng)產(chǎn)生的影響造成的。

2.2 GPS原子鐘周期項(xiàng)頻譜分析

利用GPS的擬合殘差數(shù)據(jù)對(duì)GPS原子鐘周期特性進(jìn)行了實(shí)驗(yàn)分析。以鐘類型為例，分別給出了G06、G08、G30這3顆衛(wèi)星鐘的頻譜分析圖，其中，G08衛(wèi)星鐘為銫鐘，其他2顆衛(wèi)星鐘為銣鐘。頻譜分析結(jié)果如圖5～圖7所示。

從圖5～圖7中可以看出，GPS衛(wèi)星的主周期項(xiàng)為12h、6h、24h，與GPS衛(wèi)星軌道周期11h58min保持一致，衛(wèi)星鐘差數(shù)據(jù)的主周期基本為衛(wèi)星軌道周期的1倍、1/2或2倍。銣鐘主周期項(xiàng)分布較為明顯，而銫鐘各周期項(xiàng)分布規(guī)律不清，主要周期項(xiàng)基本在12h附近。GPS衛(wèi)星的主周期項(xiàng)統(tǒng)計(jì)分析結(jié)果如表2所示。本文只給出了前6個(gè)主周期項(xiàng)，結(jié)果表明，GPS衛(wèi)星鐘的周期項(xiàng)與鐘的類型有關(guān)。

圖5 GPS-G06星的鐘差數(shù)據(jù)頻譜分析結(jié)果Fig.5 Spectral analysis results based on GPS-G06 satellite clock bias

圖6 GPS-G30星的鐘差數(shù)據(jù)頻譜分析結(jié)果Fig.6 Spectral analysis results based on GPS-G30 satellite clock bias

圖7 GPS-G08星的鐘差數(shù)據(jù)頻譜分析結(jié)果Fig.7 Spectral analysis results based on GPS-G08 satellite clock bias

h

從表2中可以看出，G08為銫鐘，周期性不好。G06和G30對(duì)應(yīng)的主周期略有不同，這主要是由于不同時(shí)期的不同衛(wèi)星軌道周期略有差異，星載鐘的主周期項(xiàng)與衛(wèi)星軌道周期基本為1倍、1/2或2倍關(guān)系。在解算過程中，不同算法對(duì)周期也有一定影響，但都在12h附近。由于鐘差本身存在一定的誤差，擬合后得到的鐘差數(shù)據(jù)也會(huì)存在一定的誤差，頻譜分析得到的鐘差周期與GPS衛(wèi)星的軌道周期近似一致。

3 結(jié)論

本文介紹了星載原子鐘周期項(xiàng)的提取過程，分析了2016年全年的BDS、GPS主周期項(xiàng)，最終的統(tǒng)計(jì)結(jié)果表明，BDS、GPS衛(wèi)星鐘除了銣鐘存在明顯的周期項(xiàng)，銫鐘周期項(xiàng)不明顯。BDS、GPS的主周期項(xiàng)基本在12h、24h、6h左右。通過本文的試驗(yàn)分析，對(duì)BDS、GPS的周期項(xiàng)可得如下結(jié)論：

1)軌道類型不同會(huì)影響衛(wèi)星鐘差周期項(xiàng)，GEO、IGSO原子鐘周期性較為明顯，而相同軌道類型的衛(wèi)星其鐘差周期項(xiàng)也有輕微差異；

2)衛(wèi)星的鐘差主周期項(xiàng)，分別為對(duì)應(yīng)衛(wèi)星軌道周期的1/2或1倍左右；

3)不同類型原子鐘的衛(wèi)星鐘差周期項(xiàng)也不同。