基于IGS精密星歷的衛星位置內插方法比較

馬俊，唐詩華，黃鷹，王文杰

(桂林理工大學土木與建筑工程學院，廣西桂林 541004)

基于IGS精密星歷的衛星位置內插方法比較

馬俊?，唐詩華，黃鷹，王文杰

(桂林理工大學土木與建筑工程學院，廣西桂林 541004)

在高精度的GPS測量數據處理過程中，獲取高精度衛星軌道是重要的環節，內插是獲取任意歷元的精密軌道信息的重要手段，本文利用IGS提供的精密星歷用MATLAB作為工具編程計算，用實例分析和比較了拉格朗日插值法、牛頓插值法、內維爾逐次插值法以及最小二乘擬合在求衛星坐標當中的應用效果。前三種方法在16階時內插精度最高，內維爾插值和牛頓插值有很多相似的特性，方便靈活，龍格現象沒有拉格朗日插值那么劇烈，將會得到廣泛應用。

GPS精密星歷；拉格朗日插值；牛頓插值；內維爾插值；最小二乘擬合

1 引 言

在GPS定位中，將GPS衛星的位置作為已知值，通過計算衛星至接收機的距離來計算測站點的坐標。因此，衛星的軌道信息是定位的基本要素，而衛星的軌道信息是通過星歷給出的，星歷誤差將成為一種起算數據誤差。這種誤差直接影響著單點定位的精度，對絕對定位也有一定的影響。為了降低星歷誤差對單點定位和相對定位的影響，提高目標在地心坐標系的絕對定位精度，IGS(International Geodetic Service)免費給出了若干精度等級和不同時間延遲的精密星歷。精密星歷是按一定時間間隔(通常為15 min)給出相應歷元時刻GPS衛星在地心坐標系的三維坐標、三維速度及衛星鐘改正數等信息。

GPS精密數據處理需要得到間隔1 s的衛星坐標。GPS接收機的采樣率一般為30 s或者15 s甚至更密，因此，要想得到某瞬時的衛星坐標就必須對精密星歷和進行高精度、快速的內插或擬合。其中較常用的插值方法是拉格朗日多項式插值法，和牛頓多項式插值法。近來也有一些學者和論述中提到內維爾逐次線性插值和切比雪夫多項式擬合法。本文著重分析了拉格朗日多項式插值、牛頓多項式插值、內維爾逐次插值以及最小二乘擬合算法，將其在精密星歷衛星位置內插應用方面做了詳細的比較，得出了一些有益的結論。

2 插值數學模型

2.1 拉格朗日多項式插值

一般情況，對于給定的n＋1個插值節點，x0

于是所求n次拉格朗日插值多項式為:

拉格朗日多項式函數模型簡單，插值效率高，收斂速度快，是經典的插值方法，但是要增加或者刪除一個插值節點需要重新構造多項式。

2.2 牛頓多項式插值

為f(x)在k＋1點x0，x1，x2，…xk的k階差商。

得到多項式牛頓插值多項式:

實際計算中我們可以把計算過程中的各階差商保存為一個差商表，如表1所示。這個差商表是個下三角陣，當增加節點時只需要在差商表最下邊增加一行就行了，原先的計算全部有效，差商表計算完成后，只需將差商表中對角線上的元素帶入牛頓多項式，就可計算出插值點結果。

牛頓多項式插值差商表 表1

2.3 內維爾插值

內維爾插值是一種線性逐次插值，主要是通過低一次多項式的組合來獲得高一次插值多項式。具體做法是，首先求出若干個一次插值多項式的值，利用它們之間的兩兩組合得到若干個二次插值多項式，進一步利用這些二次多項式間的兩兩組合來得到若干個三次插值多項式，依此類推。給定n＋1個已知點(xi，f(xi))令0階內插結果為，Ti，0＝f(xi)，下一階插值結果是上一階線性內插得到的，第j階插值的公式為:

可列出各階內維爾的插值表，如表2所示。它是一個下三角矩陣。事實上后一項就是前一項和左上角一項的線性內差，內插表中的任何一個元素都是一次內插的結果。內維爾插值同樣解決了增加節點的問題，其數學模型簡單，當第i階插值結果與第i－1階插值結果小于限差，終止插值。

內維爾插值流程表 表2

2.4 最小二乘擬合

對于給定的數據(xi，yi)(i＝0，1，2，…m)，要求在給定的函數:

使s?(x)滿足:

取得極小值。

解方程組，可得a0，a1，a2…an，帶入上式中，即可得所求函數。

3 實例分析和比較

3.1 精度分析

在IGS官方網站上下載SP3精密星歷文件，本文選取起始歷元在2011年2月23日0時0分0秒，終止歷元時刻在2011年2月23日23時45分0秒的SP3精密星歷文件作為實驗數據。本文選用衛星序號為14的歷元值作為實驗數據，取擬合的時間間隔為30 min，采用對稱內插的方法，計算時取30 min間隔的歷元作為內插點計算各個內插時間段中間時刻的衛星位置，以便將內插得到的衛星位置同精密星歷給出的衛星位置進行比較。在Matlab 7.1下進行實驗，得到14號衛星6階～25階的內插結果，IGS向全球GPS用戶提供的GPS精密星歷，其誤差一般小于5 cm，所以這幾種方法內插的精度至少要小于5 cm。

表3給出了拉格朗日插值法，牛頓多項式插值法，內維爾插值法的6階～25階的標準差，表4給出最小二乘擬合6階～50階的標準差。為了能更好地比較插值效果，筆者根據表3做出了圖1～圖4，圖1～圖3是前三種方法不同階次插值的標準差曲線圖，圖4是這前三種方法在x方向不同階次插值效果比較圖，以下將根據這些數據圖做詳細分析。

三種方法不同階次插值標準差 表3

最小二乘擬合各階次的標準差 表4

圖1 拉格朗日不同階次插值標準差

圖2 牛頓插值不同階次插值標

圖4 三種方法不同階次標準差

對于拉格朗日插值，從圖1可以看出，當階數過低或過高時，插值精度都較差，當階數取12階～21階之間時，x、y、z方向的標準差均小于50 mm，完全可以滿足各種應用。其中插值階數取14階～18階之間時，三個方向的標準差均小于10 mm，并且三個方向其中取16階時最小，由表1可知三方向的標準差分別是1.17 mm，0.78 mm及0.37 mm。

對于牛頓插值法，從圖2可以看出，和拉格朗日插值法一樣，當階數過低或過高時，插值精度都較差，當階數取12階～21階之間時，x、y、z方向的標準差均小于50 mm，完全可以滿足各種應用。其中插值階數取14階～19階之間時，三個方向的標準差均小于10 mm，取16階時均達到最小值，由表1可知三方向的標準差分別是1.88 mm，0.78 mm及0.37 mm，幾乎和拉格朗日插值在16階時的精度一樣，從圖2還可以看出從14階～17階之間，三個方向的曲線幾乎重合。

對于內維爾逐次插值，從圖3可以看出，從10階～19階x、y、z三個方向的標準差曲線圖與牛頓插值的標準差曲線圖走勢相同，當插值取12階～21階時，內插精度完全滿足各種應用需要。由表4可以看出從10階～19階，牛頓插值與內維爾插值在三個方向的標準差相同，所以維爾插值在16階時三個方向標準差均達到最小值，分別是1.88 mm，0.78 mm及0.37 mm。

由圖1至圖3中還可以看出從10階～20階之間，z方向插值的精度要比x、y方向的精度高。

對于最小二乘擬合，由表4可以看出，從16階～41階最小標準差都小于100 m，但是精度最小也只能達到分米級，小于16階或者大于41階的，精度都超過100 m，6階和50階都超過萬米，因此最小二乘擬合法不能作為衛星位置內插的選擇，因此沒有再為其作圖與其他方法進行詳細比較，下文就不再將其詳述

圖4給出了拉格朗日插值，牛頓插值及內維爾插值在x方向各階次的標準差曲線，從圖4中我們可以看出從9階～20階，牛頓插值和內維爾插值重合，從12階～21階三種方法標準差都小于50 mm，其中16階，17階時三條曲線重合，標準差同時達到最小值，但是從13階～18階拉格朗日插值標準差波動較大，19階～21階之間，拉格朗日插值法與其他兩種方法相比，標準差較小。

3.2 龍格現象

為了詳細研究比較龍格現象對拉格朗日插值、牛頓插值和內維爾插值的影響，筆者分別作出其在標準差較好的插值階數下部分歷元的誤差曲線圖，如圖5是拉格朗日插值法部分歷元的誤差曲線圖，圖6，圖7分別是牛頓插值法、內維爾插值法的部分歷元的插值曲線圖。

圖5 拉格朗日插值法16階～18階誤差曲線

圖6 牛頓插值法15階～17階誤差曲線

圖7 內維爾值法14階～17階誤差曲線

從圖5～圖7中可以看出龍格現象對三種插值方法的影響主要集中在前兩個內插節點和后兩個內插節點，三種方法在16階，17階時插值精度最好，其他階次在前兩個和后兩個插值節點震蕩都較明顯。對于拉格朗日插值法16階精度最好的內插節點是從第9個～35個歷元，精度都在10 mm以內；牛頓插值法16階精度最好的內插節點也是從第9個～35個歷元，精度同拉格朗日；對于內維爾插值法，精度小于10 mm的內插節點是從第7個歷元～第39個歷元。總的來說，內維爾逐次插值法在16階次精度小于10 mm的內插節點最多，兩端龍格震蕩最為平緩。

4 結 論

(1)最小二乘擬合精度最好的時候只能達到分米級，因此不能作為GPS衛星位置內插的方法。而將拉格朗日多項式插值、牛頓多項式插值和內維爾逐次線性插值用于GPS精密星歷插值處理，無論使用哪一種方法，插值階數都不應少于12階，最高不能超過21階，更高階數對精度提高的貢獻非常有限，甚至會產生更大的差值，并且插值過程中盡量將插值點位置選在節點的中間，并盡可能選取偶數個節點。在進行16階，17階的插值時，三種插值方法精度最高，三個方向的標準差最大不超過3 mm。隨著節點數的增加，插值精度反而下降。

(2)拉格朗日插值模型簡單、結構緊湊，利用基函數就很容易得到插值函數，是經典的插值方法。但是當插值節點增加、減少或其位置變化時，整個插值多項式的結構都會改變，需要重新計算增加了算法的復雜度，且當增大插值階數時容易出現龍格現象。

(3)內維爾插值和牛頓插值具有十分相似的特性，在14階～19階，三個方向的精度均達到毫米級。模型簡單、計算量小，在內插節點兩端龍格震蕩現象比拉格朗日插值平緩；程序方便靈活，當加入新的節點時，可保留原有計算結果而不降低程序的效率，對于短時間內插衛星位置無疑是一個好的選擇，是一種較實用的內插方法。

[1] 陳兆林，朱彥武，張書畢.一種GPS精密星歷的多項式插值方法[J].工程勘察，2007(4)

[2] 卓寧.GPS數據處理中衛星精密星歷插值方法分析[J].戰術導彈技術，2010(4)

[3] 陳遠鴻.基于精密星歷的計算衛星位置的方法[J].測繪信息與工程，2008(2)

[4] 李明峰，江國焰，張凱.IGS精密星歷內插與擬合法精度的比較[J].大地測量與地球動力學，2008，28(2)

[5] 李勇軍，丁士軍.不同星歷誤差對靜態單點定位精度的影響與分析[J].測繪工程，2010，6(19)

[6] 王磊，羅年學.GPS精密星歷的軌道內插方法比較[J].地理空間信息，2009，7(5)

[7] 馬昌鳳，林偉川.現代數值計算方法(MATLAB版)[M].北京:科學出版社，2008

[8] 張德豐.MATLAB數值分析與應用[M].北京:國防工業出版社，2007

[9] Jan Kouba，Pierre Héroux.Precise Point Positioning Using IGS Orbit and Clock Products[J].GPS Solutions，2001，5(2)

IGS Precise Ephemeris Based on the Comparison of the Ways for the Satellite Position Interpolation

Ma Jun，Tang Shihua，Huang Ying，Wang Wenjie
(Guilin University of Technology，Construction Engineering College，Guilin 541004，China)

In the high－precision GPS measurement data processing，to obtain high－precision satellite orbit is an important link.Interpolation is an important means for obtaining any precise orbit information in any epoch.This article provided by IGS precise ephemeris using MATLAB as a tool for programming comparatives and analyses the ways which contain Lagrange interpolation，Newton interpolation，Neville interpolation and least square fitting with the example that these ways is used to obtain the position of satellite.The first three in the 16－order time interpolation have the highest accuracy，Neville interpolation and Newton Interpolation has many similar features，it is convenient and flexible.Its runge phenomenon is less severe than Lagrange interpolation.So it will be widely used.

GPS Precise Ephemeris；Lagrange Interpolation；Newton interpolation；Neville interpolation；leastsquares approximation

2011—04—01

馬俊(1985—)，男，碩士研究生，主要研究方向:GPS單點定位技術。

1672－8262(2011)05－89－05

P228

A