楊文可，占建偉，龔航，朱祥維，孫廣富

(國防科學(xué)技術(shù)大學(xué) 電子科學(xué)與工程學(xué)院，湖南 長沙，410073)

采用靜地軌道(GEO)衛(wèi)星的衛(wèi)星雙向時間頻率傳遞(TWSTFT)是國際上最重要的實現(xiàn)遠(yuǎn)程時間比對的方法之一，而包括地球自轉(zhuǎn)效應(yīng)(Sagnac效應(yīng))修正、電離層時延差、對流層時延差等在內(nèi)的雙向空間傳播時延差是造成 TWSTFT所測量的鐘差存在偏差的主要因素[1]。實際上，由于各種攝動力的存在[2]，GEO衛(wèi)星相對于地面站存在相對運動。國際權(quán)度局(BIPM)給出的TWSTFT數(shù)據(jù)處理指導(dǎo)[1]以理想GEO衛(wèi)星為基礎(chǔ)，給出了計算Sagnac效應(yīng)修正量、電離層時延差的方法，該方法是國際上一般處理TWSTFT數(shù)據(jù)的方法，本文簡稱其為BIPM法。雖然其中提到了衛(wèi)星運動對雙向Sagnac效應(yīng)修正量、雙向幾何距離時延差的影響量級，但未給出估算方法。研究衛(wèi)星運動對TWSTFT影響的文獻(xiàn)中，一般只分析了衛(wèi)星運動對 Sagnac效應(yīng)修正量的影響[3?6]及雙向幾何距離差[3?4,7?8]，而忽略了衛(wèi)星運動對電離層時延差、對流層時延差的影響。其中，Sagnac效應(yīng)修正量的分析一般以BIPM法中的Sagnac效應(yīng)修正公式為基礎(chǔ)[3,5?6]；分析雙向幾何距離差時，一般是以兩地面站與衛(wèi)星之間的距離差除以光速估算得到地面站信號到達(dá)衛(wèi)星的時差，再結(jié)合衛(wèi)星速度進(jìn)行估算[7?8]。以上分析方法未考慮信號真實傳播路徑，未綜合考慮信號傳播方向上經(jīng)歷的Sagnac效應(yīng)、幾何路徑時延、電離層時延、對流層時延等，并且一般以模型[3]或精度為幾十千米的軌道數(shù)據(jù)[4,8]來描述衛(wèi)星運動，存在誤差。在北斗導(dǎo)航系統(tǒng)中，基于GEO導(dǎo)航衛(wèi)星的TWSTFT實現(xiàn)站間時間同步[9]，且該GEO衛(wèi)星具有2 h擬合精度好于0.1 m、徑向誤差為厘米量級的廣播星歷[10]。不同于前述文獻(xiàn)所提方法，本文作者以廣播星歷數(shù)據(jù)為依托，對衛(wèi)星或地面站的坐標(biāo)進(jìn)行旋轉(zhuǎn)來修正Sagnac效應(yīng)，獲知真實信號傳播路徑，并綜合考慮幾何路徑時延、電離層時延修正、對流層時延修正等對衛(wèi)星與地面站之間的單向空間傳播時延建模，提出一種利用迭代法分別計算上行、下行信號時延，并精確計算TWSTFT中雙向空間傳播時延差的方法。

1 衛(wèi)星雙向時延模型

當(dāng)考慮衛(wèi)星相對于地面站的相對運動時，衛(wèi)星雙向時間頻率傳遞的時延模型如圖1所示。

圖1中，τUi為地面站i(i=1, 2)到衛(wèi)星的上行信號空間傳播時延，包含4部分：由地面站i到衛(wèi)星的偽幾何距離時延τUi|SatP、地球自轉(zhuǎn)效應(yīng)(Sagnac效應(yīng))時延修正τSCUi、電離層時延修正τUi|Ion、對流層時延修正τUi|Trop；τDi為衛(wèi)星到地面站i的下行信號空間傳播時延，也包含4部分：由衛(wèi)星到地面站i的偽幾何距離時延τDi|SatP、Sagnac效應(yīng)時延修正τSCDi、電離層時延修正τDi|Ion、對流層時延修正τDi|Trop；τTi和τRi分別為地面站i發(fā)射信號時延和接收信號時延；τSij(i=1, 2,j=1, 2,i≠j)為衛(wèi)星由地面站i向地面站j轉(zhuǎn)發(fā)信號時延。

其中，之所以稱τUi|SatP和τDi|SatP為偽幾何距離時延，是由于該幾何距離時延未考慮 Sagnac效應(yīng)的影響，不同于將偽幾何距離時延經(jīng)過Sagnac效應(yīng)修正得到的地面站到衛(wèi)星(或衛(wèi)星到地面站)真實傳播幾何路徑時延，記為τUi|Sat(或τDi|Sat)。

圖1 衛(wèi)星雙向時間頻率傳遞的時延模型示意圖Fig. 1 Schematic diagram of two-way satellite time and frequency transfer with satellite in motion

當(dāng)記以鐘i為參考測量所得的地面站j的信號時延為PRij(i=1, 2，j=1, 2，i≠j)時，鐘1與鐘2之間的時差，記為T1?T2，可表示為[1]：

式中：0.5((τT1?τR1)?(τT2?τR2))為雙向地面站發(fā)射時延與接收時延差；0.5(τS12?τS21)為雙向衛(wèi)星轉(zhuǎn)發(fā)時延差，這兩者可通過標(biāo)定獲得[1,11]。因此，為精確計算雙向時延差，即式(1)等式右邊除雙向測量時差0.5(PR21?PR12)以外，還需確定雙向空間傳播時延差 0.5((τU1?τD1)?(τU2?τD2))，記為τUD，且可分解為

2 算法原理

對下行或上行信號傳播時延的估算是一個基于信號路徑估計的迭代計算的過程，基本流程如圖2所示。

圖2 單向信號傳播時延迭代計算流程Fig. 2 One-way signal propagation delay iterative calculation flow

2.1 下行信號傳播時延模型

如圖1所示，以衛(wèi)星到地面站1的下行鏈路為例，下行信號時延τDn1，可表示為

其中，地面站1信號接收時延可通過標(biāo)校獲得[11]。以下分析除τR1以外的時延計算模型，即

2.1.1 Sagnac效應(yīng)修正

在式(2)中對下行信號時延的 Sagnac效應(yīng)的修正是以時延修正量τSCD1的方式表示的。然而，本質(zhì)上，下行信號時延的 Sagnac效應(yīng)修正是為了獲得衛(wèi)星到地面站方向的信號傳播路徑。為此，本文對衛(wèi)星位置坐標(biāo)進(jìn)行修正，獲得信號實際傳播路徑，具體過程如下。

(1) 假設(shè)下行信號時延為τDn1，地面站1接收信號時刻為tr1，則可知衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21為

(2) 對衛(wèi)星廣播星歷進(jìn)行插值，獲得在ttr21時刻的地心地固(ECEF)坐標(biāo)系下衛(wèi)星位置{xs,ys,zs}，記為S；

(4) 地面站1在tr1的ECEF坐標(biāo)系下的坐標(biāo){x1,y1,z1}，記為GS1，與衛(wèi)星位置S′之間的連線方向是實際信號傳播路徑。其中，GS1為已知量。另外，據(jù)式(6)，S′由tr1，ttr21和S共同決定。已知tr1，同時，當(dāng)星歷已知時，S實際上由ttr21決定。由此，S′是ttr21的函數(shù)，且τR1已知，則據(jù)式(5)和式(4)可知S′實際上是τD1的函數(shù)，記為

2.1.2 幾何距離時延

當(dāng)經(jīng)過 Sagnac效應(yīng)修正后，偽幾何距離時延τD1|SatP修正為幾何距離時延τD1|Sat，即

其中：c為光在真空中傳播的速度，為299 792 458 m/s。

可見：幾何距離時延τD1|Sat是地面站1坐標(biāo)GS1與經(jīng) Sagnac效應(yīng)修正的轉(zhuǎn)發(fā)信號時刻衛(wèi)星坐標(biāo)S′的函數(shù)。若隱去已知量，結(jié)合式(7)，則τD1|Sat也可簡記為τD1的函數(shù)：

2.1.3 電離層時延修正

由于電離層折射引起的偽碼相位時延修正量τD1|Ion可寫為[1]

其中：s為傾斜因子，由信號傳播路徑?jīng)Q定，是經(jīng)Sagnac效應(yīng)修正的轉(zhuǎn)發(fā)信號時刻ttr21衛(wèi)星坐標(biāo)S′和地面站 1坐標(biāo)GS1的函數(shù)；TEC為地面站1頂端垂直方向上總電子數(shù)，可采用IGS提供的電離層網(wǎng)格參數(shù)[13]進(jìn)行時間和空間上的線性插值計算得到，是衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21和地面站1坐標(biāo)GS1的函數(shù)；fd為下行信號載波頻率，為已知量。

若隱去已知量，結(jié)合式(6)，則τD1|Ion同樣可簡記為τD1的函數(shù)：

2.1.4 對流層時延修正

對流層時延修正模型有多種，本文使用 Collins模型，則τD1|Trop可寫為[14]

其中：h為地面站1高度；ε為地面站1對衛(wèi)星的仰角。

由于地面站1高度h可由其ECEF坐標(biāo)GS1轉(zhuǎn)換得到，地面站1對衛(wèi)星仰角ε是GS1和S′的函數(shù)。則隱去已知量，結(jié)合式(7)，τD1|Trop可簡寫為

2.2 上行信號空間傳播時延

仍以地面站2到衛(wèi)星到地面站1的鏈路為例，說明當(dāng)計算得到衛(wèi)星到地面站1的下行信號空間傳播時延τD1后，如何計算地面站 2到衛(wèi)星的上行信號空間傳播時延，記為τU2。

上行信號時延τUp2為

其中，地面站2信號發(fā)送時延可通過標(biāo)校獲得[11]。以下，分析除τT2以外的時延計算模型，即

當(dāng)忽略衛(wèi)星轉(zhuǎn)發(fā)信號時延，及期間的衛(wèi)星位置改變，認(rèn)為衛(wèi)星轉(zhuǎn)發(fā)信號時刻的位置即為衛(wèi)星接收信號時刻的位置，則已知衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21和該時刻ECEF坐標(biāo)系下的坐標(biāo)為S，則可按照類似下行單向時延的模型對上行單向時延進(jìn)行建模。

區(qū)別僅在于，修正Sagnac效應(yīng)時，需修正地面站2(信號發(fā)射方)的位置，獲知信號實際傳播路徑，其過程如下。

(1) 已知衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21，則可知地面站2發(fā)射信號時刻tt2為：

(2) 地球旋轉(zhuǎn)速率設(shè)為，將tt2時刻地面站2在ECEF坐標(biāo)系下的坐標(biāo){x2,y2,z2}，記為GS2，修正至衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21的ECEF坐標(biāo)系下坐標(biāo)}，記為：

(3) 衛(wèi)星轉(zhuǎn)發(fā)信號時刻ttr21的ECEF坐標(biāo)系下地面站2的坐標(biāo)到衛(wèi)星位置坐標(biāo)S之間的連線方向為實際信號傳播路徑。其中，據(jù)式(21)，由ttr21，tt2和GS2共同決定。而GS2為已知量，ttr21已知，tt2由τUp2和ttr21共同決定。隱去已知量，GS2 ′實際上是τU2的函數(shù)。

按照類似于下行單向時延的建模和推導(dǎo)過程，可得到：

其中：

即上行信號幾何距離時延τU2|Sat，電離層時延修正τU2|Ion，對流層時延修正τU2|Trop都是τU2的函數(shù)。

2.3 迭代法計算

將式(9)，(11)，(13)代入式(4)，得到：

另外，式(23)與式(21)具有類似的形式，可通過迭代法[15]求等式(21)和(23)的解τD1和τU2。以求解式(21)的解τD1的過程進(jìn)行說明：

(1) 設(shè)定初始解τD1(0)=130 ms：

(2) 第j(j≥1)次迭代計算過程如下：

(3) 設(shè)置門限λ(λ＞0)，收斂條件為：

其中，門限選取為時延有效數(shù)據(jù)位量級的1/10。若有效數(shù)據(jù)位量級為1 ps，則λ=0.1 ps。

(4) 當(dāng)滿足式(16)的收斂條件時，停止迭代計算，則

當(dāng)計算得到地面站2到衛(wèi)星再到地面站1以及地面站1到衛(wèi)星再到地面站2的雙向鏈路的下行單向時延τDn1和τDn2以及上行單向時延τUp1和τUp2時，可求得雙向空間傳播時延差τUD。

3 實驗結(jié)果

本文以北斗導(dǎo)航系統(tǒng)中相距3 376 km的地面站1、地面站2之間通過位于東經(jīng)140°的GEO衛(wèi)星建立的C波段TWSTFT鏈路的2 d實測數(shù)據(jù)為例，采用本文方法和BIPM法計算雙向傳播時延差，并分析兩者的差異。

圖3所示為GEO衛(wèi)星的星下點軌跡。由圖3可見：GEO衛(wèi)星的星下點軌跡具有以d為單位的周期波動(以下簡稱日波動)，經(jīng)度變化幅度為0.09°，緯度變化幅度為3.18°。

圖3 GEO衛(wèi)星的星下點軌跡Fig. 3 Sub-satellite point track of GEO satellite

雙向幾何距離時延差τUD|Sat包含偽幾何距離時延差τUD|SatP和Sagnac效應(yīng)修正量τUD|SC。如圖4所示，使用BIPM法，即認(rèn)為GEO衛(wèi)星相對地面站無相對運動時，τUD|Sat為定值?84.149 ns；使用本文方法，即考慮衛(wèi)星運動以及信號傳播路徑時，發(fā)現(xiàn)τUD|Sat呈現(xiàn)幅度為0.412 ns的日波動，均值為?84.266 ns，與BIPM法結(jié)果偏差77 ps。

圖4 雙向幾何距離時延差(包含雙向Sagnac效應(yīng)修正)τUD|SatFig. 4 Two-way geometry signal path difference (including two-way Sagnac effect correction) τUD|Sat

圖5 雙向電離層時延的雙差(本文方法與BIPM法計算所得雙向電離層時延差τUD|Ion 的差值)Fig. 5 Double difference of two-way ionospheric corrections(the difference of two-way ionospheric correction differenceτUD| Ion calculated by method introduced in this work and BIPM method)

圖 5所示為本文方法計算得雙向電離層時延差τUD|Ion與BIPM法計算得雙向電離層時延差τUD|Ion的差值。由于衛(wèi)星運動造成地面站到衛(wèi)星，以及衛(wèi)星到地面站的雙向信號路徑改變，使得雙向電離層時延修正量較衛(wèi)星相對地面站靜止的情況有所改變。并且，由于電離層電子濃度在白天較大，在夜晚較小，使得該差值呈現(xiàn)日波動，且幅值白天較大，夜晚較小，最大幅度為0.120 ns，均值為31 ps。

雙向?qū)α鞑顣r延差如圖6所示。從圖6可見：當(dāng)不考慮衛(wèi)星相對地面站運動時(BIPM法)，雙向?qū)α鲗訒r延差τUD|Trop幾乎為0；使用本文方法計算，發(fā)現(xiàn)D|Trop也呈現(xiàn)日波動，幅度為0.3 ps，均值為?0.9 ps。

同時考慮τUD|Sat，τUD|Ion和τUD|Trop，得到衛(wèi)星雙向傳播時延差τUD。衛(wèi)星雙向傳播時延的雙差(本文方法與 BIPM 法計算所得衛(wèi)星雙向傳播時延差τUD的差值)如圖 7所示。從圖 7可見：對比本文方法結(jié)果與BIPM 法結(jié)果，衛(wèi)星星下點軌跡的日波動引入了τUD幅度為0.406 ns，均值為?0.087 ns的日波動。

圖6 雙向?qū)α鲗訒r延差τUD|TropFig. 6 Two-way tropospheric correction difference τUD| Trop

圖7 衛(wèi)星雙向傳播時延的雙差(本文方法與BIPM法計算所得衛(wèi)星雙向傳播時延差τUD的差值)Fig. 7 Double difference of two-way propagation delay (the difference of non-reciprocity of satellite two-way propagation delay τUD calculated by method introduced in this work and BIPM method)

4 結(jié)論

(1) 提出了一種估計信號傳播路徑，結(jié)合衛(wèi)星星歷和迭代法精確計算衛(wèi)星雙向時間頻率傳遞(TWSTFT)中雙向傳播時延差的新方法。該方法的特點在于：通過修正衛(wèi)星位置或者地面站位置來進(jìn)行Sagnac效應(yīng)修正，并由此獲得信號真實傳播路徑，并以此為基礎(chǔ)，綜合考慮地面站到衛(wèi)星的幾何距離時延、Sagnac效應(yīng)、電離層時延修正、對流層時延修正，計算雙向空間傳播時延差。該方法只依賴于衛(wèi)星星歷、兩地面站坐標(biāo)和接收信號時刻，可推廣應(yīng)用于采用IGSO衛(wèi)星等非靜地軌道衛(wèi)星的TWSTFT中。另外，還可繼續(xù)研究兩地面站時間與衛(wèi)星星歷的參考時間之差造成的影響。

(2) 本文方法能夠精確估計衛(wèi)星雙向傳播時延差，反映衛(wèi)星運動對TWSTFT的影響；并發(fā)現(xiàn)衛(wèi)星運動不但引起雙向幾何距離時延差(包含 Sagnac效應(yīng)修正)幅度達(dá)0.412 ns的日波動和均值達(dá)77 ps的偏差，對于實驗中的C波段TWSTFT鏈路而言，還造成了雙向電離層時延差幅度達(dá)0.120 ns的日波動和均值達(dá)31 ps的偏差，必須從雙向測量結(jié)果中扣除，否則將在鐘差中疊加具有偏差的日波動，影響對站間鐘差的評估。

[1]Bauch A, Piester D, Fujieda M, et al. Directive for operational use and data handling two-way satellite time and frequency transfer (TWSTFT)[R]. Paris: BIPM, 2011: 6?8.

[2]周善石. 基于區(qū)域監(jiān)測網(wǎng)的衛(wèi)星導(dǎo)航系統(tǒng)精密定軌方法研究[D]. 上海: 中國科學(xué)院研究生院, 2011: 31?36.ZHOU Shanshi. Studies on precise orbit determination theory and application for satellite navigation system with regional tracking network[D]. Shanghai: Graduate University of Chinese Academy of Sciences, 2011: 31?36.

[3]Hongwei S, Imae M, Gotoh T. Impact of satellite motion on two-way satellite time and frequency transfer[J]. Electronics Letters, 2003, 39(5): 482?483.

[4]Piester D, Bauch A, Fujieda M, et al. Studies on instabilities in long-baseline two-way satellite time and frequency transfer(TWSTFT) including a troposphere delay model[C]// The 39th Annual Precise Time and Time Interval Systems and Applications Meeting. Long Beach, 2007: 211?222.

[5]武文俊, 李志剛, 楊旭海, 等. 衛(wèi)星雙向時間頻率傳遞中的Sagnac效應(yīng)[J]. 宇航學(xué)報, 2012, 33(7): 936?941.WU Wenjun, LI Zhigang, YANG Xuhai, et al. Sagnac effect on two-way satellite time and frequency transfer[J]. Journal of Astronautics, 2012, 33(7): 936?941.

[6]Tseng W H, Feng K M, Lin S Y, et al. Sagnac effect and diurnal correction on two-way satellite time transfer[J]. IEEE Transactions on Instrumentation and Measurement, 2011, 60(7):2298?2303.

[7]Parker T E, Zhang V. Sources of instabilities in two-way satellite time transfer[C]// Proceedings of the 2005 IEEE International Frequency Control Symposium and Exposition. Vancouver, 2005:745?751.

[8]楊旭海, 李志剛, 李偉超, 等. 基于小傾角 IGSO衛(wèi)星的雙向時間傳遞試驗[C]// 2007全國時間頻率學(xué)術(shù)會議論文集. 宜昌,2007: 372?376.YANG Xuhai, LI Zhigang, LI Weichao, et al. Two-way time transfer experiment based on small obliquity IGSO satellite[C]//2007 National Time and Frequency Research Conference Proceeding. Yichang, 2007: 372?376.

[9]HAN Chunhao, YANG Yuanxi, CAI Zhiwu. Beidou navigation satellite system and its time scales[J]. Metrologia, 2011, 48(4):s213?s218.

[10]黃勇, 胡小工, 王小亞, 等. 中高軌衛(wèi)星廣播星歷精度分析[J].天文學(xué)進(jìn)展, 2006, 24(1): 81?87.HUANG Yong, HU Xiaogong, WANG Xiaoya, et al. Precision analysis of broadcast ephemeris for medium and high orbit satellites[J]. Progress in Astronomy, 2006, 24(1): 81?87.

[11]李星, 耿淑敏，李垣陵. 雙向時間同步系統(tǒng)的設(shè)備時延校準(zhǔn)技術(shù)研究[J]. 中國空間科學(xué)技術(shù), 2011, 31(4): 23?28.LI Xing, GENG Shumin, LI Yuanling. Equipment delay calibration technique for two-way time synchronization system[J]. Chinese Space Science and Technology, 2011, 31(4):23?28.

[12]Kaplan E D, Hegerty C J. GPS原理與應(yīng)用[M]. 2版. 寇艷紅,譯. 北京: 電子工業(yè)出版社, 2008: 229.Kaplan E D, Hegerty C J. Understanding GPS principles and applications[M]. 2nd ed. KOU Yanhong, trans. Beijing:Publishing House of Electronics Industry, 2008: 229.

[13]Kouba J. A guide to using International GNSS service (IGS)products[EB/OL]. [2009?05?13]. http://igscb.jpl.nasa.gov/igscb/resource/pubs/UsingIGSProductsVer21.pdf

[14]Farrell J, Barth M. The global positioning system and inertial navigation[M]. New York: McGraw Hill, 1998: 156.

[15]葉其孝, 沈永歡. 實用數(shù)學(xué)手冊[M]. 北京: 科學(xué)出版社, 2010:740.YE Qixiao, SHEN Yonghuan. Hand book of applied mathematics[M]. Beijing: Science Press, 2010: 740.