曹鳳志,方書山
(1.四川二灘國(guó)際工程咨詢有限責(zé)任公司,四川 成都610072;2.中國(guó)測(cè)繪科學(xué)研究院,北京100038)
目前IGS網(wǎng)站上能下載各個(gè)GPS分析中心提供的精密衛(wèi)星鐘差CLK文件,文件里的數(shù)據(jù)內(nèi)容包括各個(gè)GPS跟蹤站接收機(jī)鐘差及其鐘速和衛(wèi)星鐘差及其鐘速,采樣間隔有5min及30s兩種。衛(wèi)星的鐘差在短時(shí)間內(nèi)雖然會(huì)發(fā)生抖動(dòng),但是長(zhǎng)期看來卻呈現(xiàn)一定的規(guī)律性,根據(jù)已知節(jié)點(diǎn)的衛(wèi)星鐘差及其鐘速可以利用相應(yīng)的方法內(nèi)插出任意節(jié)點(diǎn)的衛(wèi)星鐘差。
IGS提供的鐘差文件名的命名規(guī)則為:前三位是分析中心的名稱代碼,中間四位是GPS周,最后一位是 GPS日,后綴的“clk”或者“clk_30s”表示是鐘差產(chǎn)品文件。如選取IGS精密鐘差星歷文件“igs16176.clk_30s”為例子,“igs”表示機(jī)構(gòu)名,“1617”表示 GPS周,“6”表示一周的第六天,“clk”表示是鐘差文件,30s為采樣率,單位是s),它對(duì)應(yīng)的時(shí)間為2011年1月8日,采用率30s.在文件里面,截取一部分內(nèi)容出來,如截取內(nèi)容“AS G16 2011 01 08 21 59 0.0000 2-.523 556 281 012E-04-.323 556 769 012E-10”,此內(nèi)容指的是 GPS 16號(hào)衛(wèi)星G16在2011年1月8日21點(diǎn)59分0秒的衛(wèi)星鐘差和鐘速分別是-0.523 556 281 012E-04和-0.323 556 769 012E-10,中間的標(biāo)識(shí)符 “2”指的是鐘差數(shù)據(jù)類型,“2”表示包含鐘差和鐘速,“1”表示僅僅包含鐘差。
分析其衛(wèi)星鐘差的變化情況。選其中16號(hào)星一整天的衛(wèi)星鐘差數(shù)據(jù),事后精密鐘差變化如圖1所示。
圖1 16號(hào)衛(wèi)星一天內(nèi)鐘差變化
從圖1中可以看出,衛(wèi)星鐘差變化值(一天范圍內(nèi))有一定的變化趨勢(shì)。同時(shí)根據(jù)鐘差文件,要想得到任意時(shí)刻點(diǎn)內(nèi)的鐘差值,可以想到利用一定的方法進(jìn)行鐘差插值。
在已知的數(shù)據(jù)點(diǎn)之間利用某種算法估算出新的數(shù)據(jù)點(diǎn)稱為插值。精密鐘差插值的常見方法有很多,有線性內(nèi)插法,二次內(nèi)插法,拉格朗日內(nèi)插法等,這些方法進(jìn)行插值時(shí)并沒有考慮鐘差的鐘速信息。根據(jù)前面引言介紹的,衛(wèi)星鐘差的變化成一定的趨勢(shì),根據(jù)鐘差文件的內(nèi)容,可知衛(wèi)星的鐘差信息由時(shí)刻,鐘差和鐘速組成,在數(shù)學(xué)上可以定義為在節(jié)點(diǎn)x上的函數(shù)f(x)及其導(dǎo)數(shù)f′(x),這種數(shù)值分析可以歸結(jié)為埃爾米特插值[1]。
對(duì)于給定的函數(shù)如表1所示。
表1 埃爾米特插值形式
其中,xi∈ a,[b]且xi互異,尋求一個(gè)2n+1次多項(xiàng)式
其幾何意義是曲線y=H2n+1(x)與曲線y=f(x)不但在xi處重合,而且在xi處有公切線。
為求得H2n+1(x)多項(xiàng)式,構(gòu)造兩組2n+1次多項(xiàng)式,使αj(x)與βj(x )(i=0,1,…,n) ,滿足條件
其中l(wèi)j(x)是n+1個(gè)節(jié)點(diǎn)x0、x1、…、xn上的n次基本插值多項(xiàng)式。
將αj(x)與βj(x)帶入公式(2),則有埃爾米特插值多項(xiàng)式
作為重要的特例,當(dāng)n=1時(shí),由上可以得到滿足插值條件
對(duì)精密衛(wèi)星鐘差,可以采用兩點(diǎn)三次埃爾米特插值來對(duì)原始的5min數(shù)據(jù)進(jìn)行30s的插值。
采用兩點(diǎn)三次埃爾米特插值可以從初始IGS鐘差CLK文件中的5min精密鐘差得到30s的精密鐘差,用于后續(xù)的精密鐘差應(yīng)用。以“igs16176.clk”文件中的PRN 16號(hào)星為例,說明內(nèi)插的結(jié)果及其精度[3]。選取“igs16176.clk”文件(采樣率5min)中一次5min區(qū)間的16號(hào)衛(wèi)星(試驗(yàn)中鐘差文件選取對(duì)應(yīng)時(shí)間為2011年1月8日,跨度為0點(diǎn)0分0秒至0點(diǎn)5min 0s共5min),根據(jù)5min采樣率取得2個(gè)節(jié)點(diǎn)x0,x1,以及對(duì)應(yīng)2個(gè)節(jié)點(diǎn)的鐘差值y0,y1,鐘速值m0,m1,可以得到如表2所示。
表2 兩點(diǎn)三次埃爾米特插值數(shù)據(jù)表
根據(jù)表2,利用公式(4)進(jìn)行兩點(diǎn)三次埃爾米特內(nèi)插30s的鐘差。為了顯示內(nèi)插的精度變化,和下載的鐘差文件“igs16176.clk_30s”(采樣30s,鐘差值當(dāng)做真值),進(jìn)行相同時(shí)間節(jié)點(diǎn)的比較[4],其差值見圖2.
圖2的結(jié)果顯示利用“igs16176.clk”文件進(jìn)行兩點(diǎn)三次埃爾米特30s采樣插值,內(nèi)插鐘差與“igs16176.clk_30s”文件在公共點(diǎn)衛(wèi)星鐘差的差值不超過0.3ns,內(nèi)插精度足夠高[5-6]。
圖2 16號(hào)星衛(wèi)星兩點(diǎn)三次埃爾米特插值與真值的差異
高分辨率的衛(wèi)星鐘差可以通過埃爾米特插值來實(shí)現(xiàn)[7]。它使用IGS網(wǎng)絡(luò)資源的5min采樣率的衛(wèi)星鐘差CLK文件,對(duì)其進(jìn)行兩點(diǎn)三次埃爾米特插值,可以得到取值區(qū)間內(nèi)30s或更高的采樣率的衛(wèi)星鐘差。從插值精度而言,精度已足夠,能夠滿足一般導(dǎo)航與定位的需求。
[1]孫志忠,袁慰平,聞?wù)鸪?數(shù)值分析[M].南京:東南大學(xué)出版社,2006.
[2]RAINER K.Numerical analysis[M].Springer-Verlag,1998.
[3]KAUBA I.A guide to using international GPS service(IGS)products[R].IGS Central Bureau,IGS,Pasadena,2009.
[4]王沫然.MATLAB與科學(xué)計(jì)算 [M].北京:電子工業(yè)出版社,2006.
[5]李明峰,江國(guó)焰,張 凱.IGS精密星歷內(nèi)插與擬合法精度的比較[J].大地測(cè)量與地球動(dòng)力學(xué),2008,28(2):76-80.
[6]洪 櫻,歐吉坤,彭碧波.GPS衛(wèi)星精密星歷和鐘差三種內(nèi)插方法的比較[J].武漢大學(xué)學(xué)報(bào)·信息科學(xué)版,2006,31(6):516-519.
[7]DING Yu,CHEN Yang-quan.Advanced applied mathematical problem solutions with MATLAB[M].Springer-Verlag,2008.