張小紅，陳興漢，郭 斐

武漢大學(xué)測(cè)繪學(xué)院，湖北 武漢430079

1 引 言

近年來，精密單點(diǎn)定位（precise point positioning，PPP）技術(shù)的快速發(fā)展使其成為獲取高精度測(cè)站坐標(biāo)的一種重要技術(shù)手段，已被廣泛應(yīng)用于地球動(dòng)力學(xué)、電離層和對(duì)流層延遲估計(jì)等諸多地學(xué)研究。但是，受衛(wèi)星定位幾何條件的限制，GNSS精密定位獲得的測(cè)站高程精度明顯次于平面位置精度（1.5～2倍）。因此，在利用GNSS跟蹤站的觀測(cè)資料計(jì)算地殼垂直運(yùn)動(dòng)量時(shí)，一些微弱的地學(xué)信號(hào)往往被淹沒在測(cè)站坐標(biāo)時(shí)間系列的噪聲中而無法準(zhǔn)確提取［1－2］。造成 GNSS定位高程精度偏低的主要原因是地球遮擋影響，接收機(jī)捕獲的衛(wèi)星信號(hào)僅來自地平線以上的可視衛(wèi)星，這種幾何上非對(duì)稱性的觀測(cè)條件使得接收機(jī)鐘差、天頂對(duì)流層延遲以及測(cè)站高程參數(shù)之間存在顯著的數(shù)學(xué)相關(guān)性（瞬時(shí)相關(guān)性高達(dá)80%～90%），這種強(qiáng)相關(guān)性使得高程方向的定位精度偏低［3－4］。盡管引入低高度角衛(wèi)星的觀測(cè)數(shù)據(jù)能夠在一定程度上削弱接收機(jī)鐘差、測(cè)站高程以及天頂對(duì)流層延遲三者之間的相關(guān)性，但是由于低高度角衛(wèi)星受到的系統(tǒng)誤差（大氣殘余誤差、多路徑效應(yīng)）和觀測(cè)噪聲（低信噪比）明顯較大，反而可能降低水平分量的定位精度。解決這一問題的一種有效途徑是引入外部約束信息。假設(shè)能夠使用一些高頻穩(wěn)度的接收機(jī)鐘，如原子鐘或超穩(wěn)石英鐘，充分利用其短期頻穩(wěn)的約束條件，對(duì)接收機(jī)鐘差進(jìn)行建模與預(yù)報(bào)，則有望顯著減少接收機(jī)鐘差參數(shù)個(gè)數(shù)，降低上述3種參數(shù)的相關(guān)性，從而更加準(zhǔn)確地分離（估計(jì)）出測(cè)站高程與天頂對(duì)流層延遲。

長(zhǎng)期以來，國(guó)內(nèi)外學(xué)者對(duì)時(shí)鐘性能評(píng)估與鐘差建模方面的研究主要是針對(duì)衛(wèi)星端高性能的原子鐘，并取得了豐碩的研究成果［5－10］。針對(duì)地面或星載GNSS接收機(jī)端的鐘差建模及其應(yīng)用研究相對(duì)較少，尚處于起步階段。近年來，隨著IGS跟蹤站網(wǎng)硬件設(shè)備的不斷升級(jí)，目前已有超過130個(gè)跟蹤站（接收機(jī)端）配置了高準(zhǔn)確度和高頻穩(wěn)度的氫原子鐘、銣原子鐘或銫原子鐘［11］；此外，一些科學(xué)試驗(yàn)衛(wèi)星如GRACE重力衛(wèi)星上也搭載了超穩(wěn)振蕩器［12］，其短期（1～1000s）頻穩(wěn)度高達(dá)1×10－13～3×10－13（文獻(xiàn)［13］）。這些高穩(wěn)定度的接收機(jī)鐘使得地面或星載接收機(jī)鐘差建模成為可能。

因此，本文在分析評(píng)價(jià)當(dāng)前IGS跟蹤站的幾類原子鐘性能的基礎(chǔ)上，擬重點(diǎn)研究并分析鐘差建模方法在精密單點(diǎn)定位中的應(yīng)用。本文第2節(jié)將利用Allan方差法簡(jiǎn)要評(píng)估現(xiàn)有IGS幾類原子鐘的穩(wěn)定性能，第3節(jié)給出了本文鐘差建模的二次多項(xiàng)式模型，第4節(jié)重點(diǎn)討論了鐘差建模在精密度單點(diǎn)定位中的應(yīng)用及效果。

2 接收機(jī)原子鐘的性能評(píng)估

鐘差性能（頻穩(wěn)度）評(píng)估是接收機(jī)鐘差建模與預(yù)報(bào)的前提，目前常用的GNSS接收機(jī)鐘類型有石英鐘、銣原子鐘、銫原子鐘、氫原子鐘等［14］，但并非所有類型的接收機(jī)鐘都適合建模，特別是對(duì)于一些頻穩(wěn)度較差的石英鐘，即便采用復(fù)雜的鐘差模型也難以準(zhǔn)確描述時(shí)鐘的運(yùn)行特性。此外，溫度變化、空間環(huán)境（地面和星載）差異等因素也會(huì)對(duì)接收機(jī)鐘的穩(wěn)定性造成影響。表1給出了當(dāng)前IGS跟蹤站原子鐘配置統(tǒng)計(jì)［15］，圖1為3類原子鐘的全球分布圖。

表1 全球IGS跟蹤站配備高穩(wěn)原子鐘數(shù)目統(tǒng)計(jì)情況Tab.1 Number of atomic clocks in IGS tracking stations

圖1 配備高穩(wěn)原子鐘的IGS跟蹤站全球分布圖Fig.1 IGS tracking stations with high stability atomic clocks

Allan方差（或標(biāo)準(zhǔn)差）是目前最常用的時(shí)域頻率穩(wěn)定性分析方法［16］，它不僅可以用于計(jì)算接收機(jī)鐘的穩(wěn)定度，還可用于識(shí)別接收機(jī)鐘的噪聲類型、計(jì)算噪聲水平系數(shù)。為了克服傳統(tǒng)的Allan方差無法識(shí)別調(diào)相白噪聲和調(diào)相閃爍噪聲這一缺陷［17］，本文采用修正Allan方差來表征接收機(jī)鐘的時(shí)域穩(wěn)定度［18］，其計(jì)算公式為

式中，N為鐘差采樣數(shù)；τ為鐘差取樣（平滑）間隔；xi為歷元i對(duì)應(yīng)的鐘差值；n為鐘差平滑因子，一般取為。

噪聲類型識(shí)別方面，通過修正Allan方差雙對(duì)數(shù)圖的斜率來區(qū)分噪聲過程，相應(yīng)給出了振蕩器噪聲的修正Allan方差在雙對(duì)數(shù)圖上的表現(xiàn)形式，如圖2所示。從圖中可以得到，調(diào)相白噪聲（WPM）和調(diào)相閃爍噪聲（FPM）的斜率是－1，調(diào)頻白噪聲（WFM）的斜率是－0.5，調(diào)頻閃爍噪聲（FFM）的斜率是0，調(diào)頻隨機(jī)游走噪聲（RWFM）的斜率是0.5。

圖2 振蕩器噪聲類型的修正Allan方差Fig.2 The modified Allan variance of oscillator noise

為了分析接收機(jī)鐘的穩(wěn)定性，首先采用精密單點(diǎn)定位獲得了部分IGS跟蹤站30s間隔的接收機(jī)鐘差，然后對(duì)接收機(jī)鐘差序列采用修正Allan方差計(jì)算得到部分IGS跟蹤站接收機(jī)鐘（H：氫原子鐘；Cs：銫原子鐘；Rb：銣原子鐘；QUARTZ：石英鐘）的半日穩(wěn)定度。如圖3所示，其中橫軸代表時(shí)間間隔；縱軸代表對(duì)應(yīng)的Allan標(biāo)準(zhǔn)方差數(shù)值（數(shù)值越小，穩(wěn)定度越高）；兩條黑色虛線（GPS code clock，GPS phase clock）分別代表無電離層組合偽距和載波相位觀測(cè)值的等效噪聲下界［19］。

從圖3中可以看出，當(dāng)采樣間隔在30～2000s左右時(shí)，4款接收機(jī)鐘的噪聲類型主要表現(xiàn)為調(diào)頻白噪聲，當(dāng)采樣間隔取至2000s甚至更長(zhǎng)時(shí)，調(diào)頻閃爍噪聲占主導(dǎo)地位。通過比較4種類型的接收機(jī)鐘差在雙對(duì)數(shù)圖的縱坐標(biāo)值，可以明顯地看出氫原子鐘的穩(wěn)定度比銫原子鐘要好，銫原子鐘比銣原子鐘的穩(wěn)定度高，石英鐘的穩(wěn)定度最差。分析圖3中無電離層組合偽距和載波相位觀測(cè)值的噪聲下界不難發(fā)現(xiàn)，石英鐘的頻穩(wěn)度較差，即使在非常短的時(shí)間內(nèi)也無法滿足載波相位水平的鐘差建模要求（即模型預(yù)報(bào)誤差大于觀測(cè)噪聲）；銫原子鐘具有較好的短期穩(wěn)定度，但是當(dāng)平滑間隔在超過300s時(shí)，其預(yù)報(bào)誤差也將超出載波相位觀測(cè)值的噪聲水平；而氫原子鐘的穩(wěn)定度最高，可以滿足1～2h以內(nèi)的鐘差建模精度要求。

圖3 4款接收機(jī)鐘差的修正Allan方差Fig.3 The modified Allan variance of four types of receiver clock offset

根據(jù)Allan方差原理，鐘差隨機(jī)項(xiàng)中不同噪聲分量引起的模型預(yù)報(bào)誤差RMSx（τ）存在所示的對(duì)應(yīng)關(guān)系，見表2［16］。給定預(yù)報(bào)誤差的閾值RMSmax，利用表2中的模型預(yù)報(bào)誤差計(jì)算公式（RMSx（τ）≤RMSmax），即可以確定鐘差建模（模型參數(shù)）的有效時(shí)長(zhǎng)。以圖3中氫原子鐘差建模為例，設(shè)定預(yù)報(bào)誤差的閾值RMSmax為0.006m，平滑間隔2000s以內(nèi)的噪聲類型表現(xiàn)為調(diào)頻白噪聲，其預(yù)報(bào)誤差的計(jì)算公式為

式中，RMS（τ）為預(yù)報(bào)誤差；τ為預(yù)報(bào)的時(shí)間間隔；σy（τ）為平滑間隔為τ時(shí)的修正Allan標(biāo)準(zhǔn)差；c為光速。

表2 不同噪聲類型引起的模型預(yù)報(bào)誤差Tab.2 The prediction error of different types of noise

當(dāng)τ取最大值2000s時(shí)，由圖3知σy（τ）為10－14，根據(jù)式（2）計(jì)算得到的預(yù)報(bào)誤差為0.006m，滿足RMS（τ）≤RMSmax的條件，可以在30～2000 s內(nèi)任取采樣間隔進(jìn)行鐘差建模。當(dāng)采樣間隔超過2000s時(shí)，調(diào)頻閃爍噪聲占主導(dǎo)地位，調(diào)頻閃爍噪聲的預(yù)報(bào)誤差計(jì)算公式為

由圖3知，采樣間隔為2000～8000s以內(nèi)的σy（τ）基本在10－14左右，根據(jù)式（3）算得的τ為1665s。綜合以上分析，得到跟蹤站mac2鐘差建模的有效時(shí)長(zhǎng)取為2000s。

3 接收機(jī)鐘差模型

接收機(jī)鐘差可以用確定性變化分量和隨機(jī)性變化分量來描述［8］，即

式中，右邊前3項(xiàng)為鐘的確定性時(shí)間分量；a0、a1和a2依次代表接收機(jī)的鐘差、鐘速與鐘漂；εx（t）為接收機(jī)鐘差的隨機(jī)變化分量。接收機(jī)鐘差的系統(tǒng)性變化部分（確定性分量）可采用線性方程或二階多項(xiàng)式等確定性函數(shù)模型來表達(dá)［20］，而其隨機(jī)性變化分量只能從統(tǒng)計(jì)意義上來分析，其統(tǒng)計(jì)特性由噪聲冪律譜模型（式5）確定［21］。

瞬時(shí)時(shí)間偏差x（t）的功率譜密度可表示為

式中，fα（α＝－4，－3，－2，－1，0）代表5種噪聲類型的傅氏頻率；hβ（β＝－2，－1，0，1，2）為噪聲強(qiáng)度系數(shù)［22］。對(duì)于給定的原子鐘，隨著取樣時(shí)間的變化將表現(xiàn)出不同的隨機(jī)噪聲分量，但在較短的時(shí)間內(nèi)，一般只有1～2種噪聲起主導(dǎo)作用。

4 鐘差建模在精密單點(diǎn)中的應(yīng)用

4.1 附有鐘差約束的精密單點(diǎn)定位方法

在傳統(tǒng)的精密單點(diǎn)定位［23］中，采用Kalman濾波作為參數(shù)估計(jì)器，一般將接收機(jī)鐘差視作獨(dú)立的白噪聲（過程噪聲非常大），忽視了鐘差參數(shù)之間可能存在的短期相關(guān)性。因此，本文采用鐘差－鐘速二維狀態(tài)模型描述接收機(jī)鐘的動(dòng)態(tài)過程［24］，其狀態(tài)一步預(yù)測(cè)方程可表示為

式中，xp和xf分別代表接收機(jī)鐘差、鐘速；Δt為歷元間隔；wp和wf為狀態(tài)的過程噪聲。

根據(jù)誤差傳播定律［25］，得到狀態(tài)一步預(yù)測(cè)的協(xié)因數(shù)矩陣Qx，k為

式中，Φ為狀態(tài)轉(zhuǎn)移矩陣；Qw，k為系統(tǒng)過程噪聲，取決于鐘的穩(wěn)定度，可根據(jù)Allan方差或譜密度系數(shù)確定［26］

式中，Sp、Sf分別代表引起時(shí)差和頻差的噪聲譜振幅。

采用Kalman濾波算法即可獲得所有歷元的遞推解，為了確保所有歷元濾波解的精度和可靠性，本文采用雙向平滑濾波算法［27］對(duì)其結(jié)果進(jìn)行優(yōu)化。設(shè)前向?yàn)V波和后向?yàn)V波的解分別為Xk，f、Xk，b，對(duì)應(yīng)的估值協(xié)因數(shù)矩陣為Qk，f、Qk，b，則平滑濾波的估值及其協(xié)因數(shù)矩陣為

4.2 試驗(yàn)結(jié)果與分析

利用2010年4月18日WTRZ跟蹤站（配備氫原子鐘）的觀測(cè)數(shù)據(jù)（30s采樣率）和CODE分析中心提供的精密星歷和精密鐘差產(chǎn)品，依次采用以下兩種方案進(jìn)行動(dòng)態(tài)PPP解算：

方案1：采用傳統(tǒng)的逐歷元估計(jì)一維鐘差參數(shù)的方法，過程噪聲設(shè)置為3×105；

方案2：采用本文的鐘差建模方法估計(jì)鐘差－鐘速二維狀態(tài)［28］，過程噪聲由譜密度系數(shù)（本文根據(jù)氫原子鐘的特性設(shè)置了h0＝1×10－24；h－1＝4×10－29）確定。

基于上述兩種方案獲得的測(cè)站高程分量偏差、接收機(jī)鐘差參數(shù)以及測(cè)站高程分量與接收機(jī)鐘差參數(shù)之間的相關(guān)系數(shù)時(shí)序，如圖4—6所示。利用PPP逐歷元解算的方差協(xié)方差陣中的相關(guān)數(shù)據(jù)，計(jì)算得到各個(gè)歷元的相關(guān)系數(shù)ρ

式中，cov（xh，xt）代表高程分量與接收機(jī)鐘差參數(shù)的協(xié)方差；σxh代表高程分量的方差；σxt代表接收機(jī)鐘差參數(shù)的方差。

類似的，根據(jù)2010年4月18日WAB2跟蹤站（配備氫原子鐘）的觀測(cè)數(shù)據(jù)（30s采樣率）和CODE分析中心提供的精密星歷和精密鐘差產(chǎn)品采用上述兩種方案進(jìn)行動(dòng)態(tài)PPP解算，得到測(cè)站的高程分量偏差、高程分量與接收機(jī)鐘差參數(shù)之間相關(guān)系數(shù)如圖7和圖8所示。

圖4 兩種方案解算的高程方向偏差（WTRZ）Fig.4 Time series of height errors（WTRZ）

圖5 兩種方案解算的接收機(jī)鐘差參數(shù)（WTRZ）Fig.5 Time series of receiver clock errors（WTRZ）

圖6 高程分量與接收機(jī)鐘差參數(shù)的相關(guān)系數(shù)時(shí)序（WTRZ）Fig.6 Correlation coefficients between height component and receiver clock offset（WTRZ）

圖7 兩種方案解算的高程方向偏差（WAB2）Fig.7 Time series of height errors（WAB2）

圖8 高程分量與接收機(jī)鐘差參數(shù)的相關(guān)系數(shù)時(shí)序（WAB2）Fig.8 Correlation coefficient between height component and receiver clock offset（WAB2）

本文采用正反向平滑濾波算法（后處理），而不是單向?yàn)V波處理，即后處理PPP正反向平滑濾波算法消除了PPP單向?yàn)V波的收斂過程。分析上述結(jié)果不難發(fā)現(xiàn)，采用將接收機(jī)鐘差視為白噪聲的鐘差逐歷元估計(jì)方案，由其得到的測(cè)站高程分量、天頂對(duì)流層延遲與接收機(jī)鐘差參數(shù)之間存在顯著的相關(guān)性（特別是測(cè)站高程分量與接收機(jī)鐘差參數(shù)之間的相關(guān)性高達(dá)80%～90%），由此導(dǎo)致動(dòng)態(tài)PPP高程方向的定位精度偏低，且波動(dòng)較大，這就使得一些微弱的地學(xué)信號(hào)往往被淹沒在測(cè)站坐標(biāo)時(shí)間系列的噪聲中而無法準(zhǔn)確拾取。方案二則充分利用了接收機(jī)鐘差參數(shù)之間的短期穩(wěn)定性，通過鐘差建模在一定程度上削弱了測(cè)站高程分量、天頂對(duì)流層延遲與接收機(jī)鐘差參數(shù)之間數(shù)學(xué)相關(guān)性，進(jìn)而改善了動(dòng)態(tài)PPP高程方向的定位精度，相應(yīng)的RMS提高了50%左右。

此外，本文還對(duì)比分析了上述兩種方案估計(jì)的天頂對(duì)流層延遲精度。利用2010年4月18日AMC2、TWTF、WAB 2、WDC3、WTZR 5個(gè)測(cè)站（配備氫原子鐘）的GPS觀測(cè)數(shù)據(jù)進(jìn)行動(dòng)態(tài)PPP解算獲得了各觀測(cè)站的天頂對(duì)流層延遲參數(shù)（ZPD），并與IGS分析中心提供的ZPD參考值進(jìn)行比較，統(tǒng)計(jì)兩種方案解算的天頂對(duì)流層延遲參數(shù)的外符合精度，見表3。結(jié)果表明，基于接收機(jī)鐘差建模的精密單點(diǎn)定位對(duì)各測(cè)站的ZPD參數(shù)估計(jì)均有不同程度的改善，尤其是WTZR和TWTF站的ZPD精度提高了近20%。

表3 兩種方案解算的天頂對(duì)流層延遲參數(shù)的精度比較Tab.3 Comparisons of the accuracy of tropospheric delay

對(duì)比前文高程分量的改善程度，接收機(jī)鐘差建模的方法對(duì)ZPD參數(shù)估值精度的改善幅度不及測(cè)站高程分量。這主要是因?yàn)閆PD參數(shù)與接收機(jī)鐘差參數(shù)之間的相關(guān)性次于測(cè)站高程分量與接收機(jī)鐘差參數(shù)之間的相關(guān)性。圖9給出了接收機(jī)鐘差參數(shù)與測(cè)站高程分量、ZPD參數(shù)的相關(guān)系數(shù)，虛線代表接收機(jī)鐘差參數(shù)與測(cè)站高程分量之間的相關(guān)系數(shù)，實(shí)線代表接收機(jī)鐘差參數(shù)與ZPD參數(shù)之間的相關(guān)系數(shù)，4條黑色虛線是判斷相關(guān)性強(qiáng)弱程度的臨界，接收機(jī)鐘差參數(shù)與ZPD參數(shù)的相關(guān)系數(shù)的絕對(duì)值基本維持在0.3～0.5之間；接收機(jī)鐘差參數(shù)與測(cè)站高程分量之間的相關(guān)系數(shù)的絕對(duì)值卻大于0.5。

5 結(jié) 論

針對(duì)當(dāng)前許多IGS跟蹤站均配置有高性能原子鐘的這一現(xiàn)狀，本文首先采用修正Allan方差評(píng)估了不同類型接收機(jī)鐘的短期頻穩(wěn)度及鐘差建模的可行性，結(jié)果表明，不同類型接收機(jī)鐘的穩(wěn)定性有所差異，在給定模型預(yù)報(bào)誤差限值條件下（預(yù)報(bào)誤差小于觀測(cè)噪聲水平），銫原子鐘僅能滿足較短時(shí)間（數(shù)分鐘）內(nèi)的建模精度要求，氫原子鐘具有良好的短期穩(wěn)定度，能夠滿足1～2h內(nèi)的鐘差建模精度要求。實(shí)例說明本文方法能有效改善動(dòng)態(tài)PPP高程方向的定位精度，對(duì)天頂對(duì)流層延遲參數(shù)的估計(jì)也有一定改善。

［1］CRYER J D，CHAN K S.Time Series Analysis：with Applications in R［M］.2nd ed.Berlin：Springer，2008.

［2］GU Guohua.Recent Progress in Researches on Crustal Movement through GNSS（GPS）Observations［J］.Recent Developments in World Seismology，2007，343（7）：9－15.（顧國(guó)華.GNSS（GPS）觀測(cè)研究地殼運(yùn)動(dòng)的新進(jìn)展［J］.國(guó)際地震動(dòng)態(tài)，2007（7）：9－15.）

［3］ROTHACHER M，BEUTLER G.The Role of GPS in the Study of Global Change［J］.Physics and Chemistry of the Earth，1998，23（9）：1029－1040.

［4］DACH R，BEUTLER G，HUGENTOBLER U.Time Transfer Using GPS Carrier Phase：Error Propagation and Results［J］.Journal of Geodesy，2003，77（1）：1－14.

［5］HESSELBARTH A，WANNINGER L.Short－term Stability of GNSS Satellite Clocks and Its Effects on Precise Point Positioning［C］∥Proceedings of ION GNSS 2008.Savannah，GA：［s.n.］：1885－1863.

［6］HUANG Guanwen，ZHANG Qin，WANG Jigang.Research on Estimation and Prediction of GPS Satellite Clock Error［J］.Journal of Geodesy and Geodynamics，2009，29（6）：118－122.（黃觀文，張勤，王繼剛.GPS衛(wèi)星鐘差的估計(jì)與預(yù)報(bào)研究［J］.大地測(cè)量與地球動(dòng)力學(xué)，2009，29（6）：118－122.）

［7］GUO Hairong，YANG Yuaxi，HE Haibo，et al.Determination of Covariance Matrix of Kalman Filter Used for Time Prediction of Atomic Clocks of Navigation Satellites［J］.Acta Geodaetica et Cartographica Sinica，2010，39（2）：147－149.（郭海榮，楊元喜，何海波，等.導(dǎo)航衛(wèi)星原子鐘Kalman濾波中噪聲方差－協(xié)方差的確定［J］.測(cè)繪學(xué)報(bào)，2010，39（2）：147－149.）

［8］DELPORTE J，BOULANGER C，MERCIER F.Short－term Stability of GNSS On－board Clocks Using the Polynomial Method［C］∥European Frequency and Time Forum（EFTF）.［S.l.］：IEEE，2012，117－121.

［9］GONG Hang，YANG Wenke，LIU Zengjun，et al.Estimation Method of BDS On－board Clock Short－term Stability Combining Satellite Two－way with One－way Carrier Ranging［J］.Journal of National University of Defense Technology，2013，35（3）：158－163.（龔航，楊文可，劉增軍，等.衛(wèi)星雙向與單向載波聯(lián)合的北斗星載鐘短穩(wěn)評(píng)估方法［J］.國(guó)防科技大學(xué)學(xué)報(bào)，2013，35（3）：158－163.）

［10］HAUSCHILD A，MONTENBRUCK O，STEIGENBERGER P.Short－term Analysis of GNSS Clocks［J］.GPS Solutions，2013，17（3）：295－307.

［11］WANG K，ROTHACHER M.Stochastic Modeling of Highstability Ground Clocks in GPS Analysis［J］.Journal of Geodesy，2013，87（5）：427－437.

［12］DUNN C，BERTIGER W，BAR－SEVER Y，et al.Application Challenge－instrument of Grace－GPS Augments Gravity Measurements：Twin Satellites Trail Each Other in Earth Orbit.As They Pass over Contours in the Gravity Field，They First［J］.GPS World，2003，14（2）：16－29.

［13］WEINBACH U，SCHON S.GNSS Receiver Clock Modeling when Using High－precision Oscillators and Its Impact on PPP［J］.Advances in Space Research，2011，47（2）：229－238.

［14］ALLAN D W.Statistics of Atomic Frequency Standards［J］.Proceedings of the IEEE，1966，54（2）：221－230.

［15］CLKLOG.A Summary File of the Deployment History for GPS Receiver，Antenna，F(xiàn)requency Standards，and Other Equipment at IGS Stations［EB/OL］.2010－01－20［2014－05－30］ftp：∥igsws.unavco.org/igscb/station/general/loghist.txt.

［16］ALLAN D W.Time and Frequency（Time－domain）Characterization，Estimation，and Prediction of Precision Clocks and Oscillators［J］.IEEE Transactions on Ultrasonics，F(xiàn)erroelectrics and Frequency Control，1987，34（6）：647－654.

［17］ALLAN D W，WEISS M A.A Frequency－domain View of Time－domain Characterization of Clocks and Time and Frequency Distribution Systems［C］∥Proceedings of the 45th Annual Symposium on Frequency Control：667－678.

［18］ALLAN D W，HOWE D A，WALLS F L.Characterization of Clocks and Oscillators［M］.Gaithersburg：US Department of Commerce，National Institute of Standards and Technology，1990.

［19］WEINBACH U，SCH?N S.Improved GRACE Kinematic Orbit Determination Using GPS Receiver Clock Modeling［J］.GPS Solutions，2013，17（4）：511－520.

［20］GUO Hairong.Research on Theory and Method about Time－frequency Characteristic of Atomic Clocks of Navigation Satellites［D］.Zhengzhou：Information Engineering University，2006：30－33.（郭海榮.導(dǎo)航衛(wèi)星原子鐘時(shí)頻特性分析理論與方法研究［D］.鄭州：解放軍信息工程大學(xué)，2006：30－33.）

［21］BROWN R G，HWANG P Y.Introduction to Random Signals and Applied Kalman Filtering：with MATLAB Exercise and Solutions［J］.Wiley Global Education，2012（1）：1－10.

［22］BARNES J A，CHI A R，CUTLER L S，et al.Characterization of Frequency Stability［J］.IEEE Transactions on Instrumentation and Measurement，1971，1001（2）：105－120.

［23］KOUBA J.A Guide to Using International GNSS Service（IGS）Products［J］.International GNSS，2009，13（8）：1－10.

［24］STRANG G，BORRE K.Linear Algebra，Geodesy，and GPS［M］.Siam：Wellesley－Cambridge Press，1997.

［25］ZUMBERGE J F，HEFLIN M B，JEFFERSON D C，et al.Precise Point Positioning for the Efficient and Robust Analysis of GPS Data from Large Networks［J］.Journal of Geophysical Research，1997，102（B3）：5005－5017.

［26］HERRING T A，DAVIS J L，SHAPIRO I I.Geodesy by Radio Interferometry：The Application of Kalman Filtering to the Analysis of Very Long Baseline Interferometry Data［J］.Journal of Geophysical Research，1990，95（B8）：12561－12581.

［27］GELB A.Applied Optimal Estimation［M］.Cambridge：MIT Press，1974.

［28］SEEBER G.Satellite Geodesy：Foundations，Methods and Applications［M］.2nd ed.New York：Walter de Gruyter，2003.