精密衛星鐘差內插的三種方法及精度分析＊

王 俊，方書山

（1．武漢大學 資源與環境學院，湖北 武漢430079；2．中國測繪科學研究院，北京100038；3．四川省國土勘測規劃研究院，四川 成都610031）

0 引 言

目前IGS網絡上能下載各個GPS分析中心提供的精密衛星鐘差CLK文件，文件里的數據內容包括各個GPS跟蹤站接收機鐘差及鐘速和衛星鐘差及其鐘速，大部分分析中心都提供有5min、30s甚至5s間隔的精密衛星鐘差產品［1］。在實際應用中，如高采樣率的精密單點定位，高速度的航攝飛機事后差分定位等，IGS網絡提供的精密衛星鐘差文件顯然不能滿足這些要求。為了由目前IGS網絡提供的鐘差文件得到更高采樣率（1s甚至是0．1s）的鐘差值，需要利用適當的插值方法來計算。衛星鐘非常敏感，極易受到周圍環境的影響，在短時間內會發生抖動，但是長期看來卻呈現一定的規律性。根據IGS網絡資源提供的已知節點的衛星鐘差或鐘差與鐘速，可以以一定精度內插出更高分辨率的衛星鐘差，滿足導航與定位的要求。

1 精密衛星鐘差內插方法

IGS提供的鐘差文件的文件名命名規則為:前三位是分析中心的名稱代碼，中間四位是GPS周，最后一位是 GPS日，后綴的“clk”或者“clk＿30s”表示是鐘差產品文件。如選取IGS精密鐘差星歷文件“igs16176．clk＿30s”為例子，“igs”表示機構名，“1617”表示 GPS周，“6”表示一周的第六天，“clk”表示是鐘差文件，30s為采樣率，單位s），它對應的時間為2011年1月8日，采用率30s．如圖1所示。

圖1 5號衛星鐘差文件

選中內容指的是GPS 5號衛星在2011年1月8日0點0分0秒的衛星鐘差是－9．534 979 027 717e－05（單位是 s），鐘速是 1．287 027 035 820e－11（單位是秒的平方），中間的標識符 “2”指的是鐘差數據類型，“2”表示包含鐘差和鐘速，“1”表示僅僅包含鐘差。分析其衛星鐘差的變化情況。選其中5號星一整天的衛星鐘差數據，事后精密鐘差變化如圖2所示。

從圖2中可以看出，衛星鐘差變化值（一天范圍內）有一定的變化趨勢。根據鐘差文件，利用一定的方法進行鐘差插值，可得到任意時刻點的鐘差值。

圖2 5號衛星一天內鐘差變化

精密鐘差插值的常見方法有很多，有線性內插法，二次內插法，拉格朗日內插法等［2－4］。利用這些方法進行插值并分析其插值的精度。

1．1 線性插值

選取精密鐘差星歷文件igs16176．clk＿30s（采樣率30s），時間為2011年1月8日，采樣率30s．選其中5號星一整天的衛星鐘差數據進行分析。

從圖2中可以看出，衛星鐘差變化值（一天范圍內）大致呈下降趨勢（有時上升趨勢），并呈線性趨勢分布，變化振幅也不大。可以利用線性內插的方法進行鐘差估算，內插公式為y＝ax＋b［5］．選取igs16176．clk文件（采樣率5min）5號衛星，根據它相鄰5min采樣率間隔的點進行線性內插，最終得到30s采樣率的數據，得到的結果與igs16176．clk＿30s（采樣率30s）文件中的同樣5號星的結果作比較。

此時igs16176．clk＿30s文件中5號星在內插點的值視為真值。得出內插值與真值的差值如圖3所示。

圖3 5號衛星24h線性內插值與真值的差值

從圖3中可以看出，通過線性內插后，24h內插值與真值的差異絕大多數在0．2ns以內。

1．2 Lagrange多項式插值

雖然衛星鐘差具有不規律的跳變性，但是在短期內鐘差的變化趨勢可以認為是近似平滑的，可利用 Lagrange（拉 格 朗 日 ）插 值 模 型 來 內 插［3－4，6］。Lagrange（拉格朗日）插值多項式Ln（x）滿足

若n次多項式lj（x）（j＝0，1，…，n）在n＋1個節點x0＜x1＜…＜xn上滿足條件

稱這n＋1個n 次多項式l0（x），l1（x），…，ln（x）為節點x0，x1，…，xn上的n 次插值基函數。可得到n次插值基函數為

顯然它滿足式（1），滿足式（3）的插值多項式Ln（x）可以表示為

Lagrange（拉格朗日）插值模型簡單，結構緊湊，是經典的插值方法。但Lagrange的插值多項式和每一個節點都有關，當節點個數改變時，需要重新計算，且當增大插值階數時容易出現龍格現象。根據試驗，采用8階的Lagrange插值，不僅計算量少，而且精度也足夠高。具體例子為:選取igs16176．clk文件（采樣率5min）中一次40min區間的5號衛星（試驗中鐘差文件選取對應時間為2011年1月8日，跨度為0點0分0s至0點40分0s共40min），根據5min采樣率取得9個節點x0，x1，…，x8，以及對應9個節點的鐘差值y0，y1，…，y8，可以得到如表格1所示的結果。

表1 Lagrange插值節點數據表

根據公式

通過上式內插30s采樣率的數值。算出的結果與igs16176．clk＿30s文件中5號星在對應的點（當作真值）作差，差的大小如圖4所示。

圖4 5號衛星拉格朗日8階內插值與真值的差值

從圖4中可以看出，通過8階lagrange內插后，非邊緣部分，40min的內插值與真值的差異少于0．2ns．插值區間邊緣部分，精度不太好，避免此缺點的方法可以采用滑動區間進行插值，始終取區間中間部分值。

1．3 埃爾米特插值

根據精密星歷或者鐘差文件的格式，可以看到描述鐘差（不管是接收機還是衛星鐘）的形式為:時刻、鐘差、鐘速。在數學上可以定義為節點x上的函數f（x）及其導數f′（x），這種情況的數值分析可以歸結為埃爾米特插值［7］。

對于給定的函數表，如表2所示。

表2 埃爾米特插值形式表

其中xi∈ ［a，b］，且xi互異，尋求一個2n＋1次多項式為使H2n＋1（）x滿足插值條件:

其幾何意義是曲線y＝H2n＋1（x）與曲線y＝f（x）不但在xi處重合，而且在xi處有公切線。為求得H2n＋1（）x多項式，構造兩組2n＋1次多項式，使αj（x）與βj（x），i＝0，1，…，n，滿足條件:

作為重要的特例，當n＝1時，可以得到滿足插值條件 H3（x0）＝y0、H3（x1）＝y1、H′3（x0）＝m0、H′3（x1）＝m1的兩點三次埃爾米特插值多項式:

相應的余項為

對于精密衛星鐘差而言，可以采用分段兩點三次埃爾米特插值來對原始的5min采樣率數據進行30s的插值。得到30s采樣率的精密鐘差。具體過程為:

選取igs16176．clk文件（采樣率5min）中40 min區間的5號衛星（試驗中鐘差文件選取對應時間為2011年1月8日，跨度為0點0分0秒至0點40分0秒共40min），根據5min采樣率取得對于9個節點x0，x1，…，x9，以及對應9個節點的鐘差值y0，y1，…，y9，鐘速值m0，m1，…，m9，可以得到如表3所示的結果。

根據表3，每5min區間分成一組，共8組。利用公式（10）內插30s采樣率的數值。算出的結果與igs16176．clk＿30s文件中5號星在對應的點（當作真值）作差，差的大小如圖5所示。

試驗例子結果顯示，30s的內插鐘差在與igs16176．clk＿30s文件中5號星在對應的點衛星鐘差的差值均小于0．2ns．

表3 兩點三次埃爾米特內插數據表

圖5 5號衛星兩點三次埃爾米特插值與真值的差值

2 結 論

高分辨率的衛星鐘差可以通過所介紹的一次插值、拉格朗日插值、埃爾米特二次插值三種插值方法來實現。使用IGS網絡資源的5min采樣率的衛星鐘差CLK文件，利用介紹的三種方法對其分別進行插值，都可以得到取值區間內30s或更高的采樣率的衛星鐘差。由于埃爾米特插值顧及到了鐘速，其插值的精度要比一次插值與拉格朗日插值都高。從插值精度而言，三種插值方法的精度優于0．2ns，能夠滿足一般導航與定位的要求。

［1］ KAUBA I．A guide to using intern－ational GPS service（IGS）products［R］．IGS Central Bureau，IGS，Pasadena，2009．

［2］ 孫志忠，袁慰平，聞震初．數值分析［M］．南京:東南大學出版社，2006．

［3］ KRESS R．Numerical analysis［M］．New York:Springer－Verlag，1998．

［4］ 王沫然．MATLAB與科學計算 ［M］．北京:電子工業出版社，2006．

［5］ 李明峰，江國焰，張 凱．IGS精密星歷內插與擬合法精度的比較［J］．大地測量與地球動力學，2008，28（2）:76－80．

［6］ 洪 櫻，歐吉坤，彭碧波．GPS衛星精密星歷和鐘差三種內插方法的比較［J］．武漢大學學報·信息科學版，2006，31（6）:516－519．

［7］ DING Yu，CHEN Yangquan．Advanced applied mathematical problem solutions with MATLAB［M］．New York:Springer－Verlag，2008．