IGS精密星歷內插及外推方法比較分析＊

劉 翔，時振偉

（山東科技大學 測繪科學與工程學院，山東 青島266590）

0 引 言

國際GNSS服務（IGS）發布的GPS精密星歷采樣間隔為15min，在實際定位中需要間隔更小的數據，這就需要高精度、快速地對精密星歷進行內插或外推。目前，利用IGS精密星歷求任意時刻衛星坐標的方法可分為插值法和擬合法兩類方法，而且已有大量文獻闡述了這兩類方法。只要選擇合適的階數，插值點的位置等條件就能到得到很好的內插精度，但外推精度有限。由于IGS精密星歷僅給出每天00∶00∶00－23∶45∶00的星歷數據，為了得到23∶45∶00－24∶00∶00的星歷數據，有兩種方法：一是將當天的精密星歷與下一天的“拼接”，這樣就可以沿用內插方法；二是僅用當天的精密星歷進行外推［1]。前一種方法需要收集處理更多的精密星歷數據，特別是當兩天的精密星歷中出現某顆衛星的數據中斷的情況，則該方法失效，必須依賴第二種方法。而要使用第二種方法，則必須提高星歷外推的精度。

簡要介紹了三角插值的基本模型，及衛星在慣性坐標系和地固系的實際運動特點，Lagrange插值和Chebyshev擬合的數學模型在很多文獻中都有涉及［2－5]，這里不再闡述。用 Matlab編制程序實現了這些插值方法，通過精密星歷內插及外推衛星位置的實際算例，對精度進行分析，得到了一些有意義的結論，證明了由于顧及了衛星的實際運動特點，三角插值在精密星歷外推能力方面優于Lagrange插值和Chebyshev擬合。

1 GPS衛星運動特點及三角函數插值算法

1.1 衛星運動特點

IGS給出的精密星歷是以地固系為空間基準的，是在ITRF參考框架內的。以2012年8月6日PRN01衛星為例，圖1和圖3分別示出了其在地固系和慣性坐標系下的運動特點，圖2和圖4分別示出了衛星的運動軌跡，其他GPS衛星的運動特點與此類似。圖1中，一天內除了Z坐標分量表現較好的12h周期性（GPS衛星軌道周期為11 h 58min）外，X和Y都是較復雜的曲線；圖2中衛星在地固坐標系下的三維運行軌跡不太規則，主要是因為在地心地固坐標系中，衛星的坐標還要受坐標旋轉框架的影響，而慣性坐標系的坐標框架是非旋轉的［6]，所以圖3中，X、Y、Z 三坐標分量都表現出約12h周期性。圖4中衛星在慣性坐標系下的三維運行軌跡較為規則，近似為橢圓形。

GPS衛星在地固系和慣性坐標系下的基本運動特點不同的根本原因，是由于在地固系下給出的衛星星歷引入了地球本身的運動，包括自轉、歲差、章動、極移等，由此影響了衛星自身運動的規律。

而在慣性坐標系下描述衛星運動，則更能反映衛星本身的運動特點，更有規律可循［6]。

1.2 三角函數多項式插值算法

基于軌道坐標的這個特征，可以利用基于三角函數多項式來進行插值，多項式的形式為［7]：

式中：f為衛星的坐標分量x，y，z；ω為頻率因子，ω＝2π／T；T 為GPS衛星運行周期；a0，a1，a2…an為待求系數；n為三角函數的項數（插值階數）；t為對應歷元時刻。根據衛星的運動周期，對于慣性坐標系下的軌道，運用基于三角函數多項式的插值方法反映了衛星位置在坐標框架下的周期特性，可能會使插值精度更高。當然先把精密星歷中的坐標轉換成慣性坐標系下的坐標，轉換公式為［6，8]

式中：E為旋轉真春分點時角；M為極移旋轉矩陣；N為章動旋轉矩陣；P為歲差旋轉矩陣；（x，y，z）TCTS為以協議地極方向CTP為指向的協議地球坐標；（x，y，z）TCTS協議天球坐標系；R1，R2，R3分別為旋轉矩陣；GAST為真春分點的格林尼治時角；xp，yp為極移分量；ε，Δε，Δψ 分別為黃赤交角、交角章動和黃經章動；ZA，θA，ζA稱為歲差參數。

根據上述公式得到慣性坐標系坐標，求出待求系數，基于所要插值的時間，得到實時慣性系坐標，完成插值后，再轉換到地固系即可。

2 計算實例及結果比較

2.1 內插精度分析

在IGS官方網站上下載SP3精密星歷文件，選取2012年8月6日的SP3精密星歷文件作為實驗數據，歷元時刻為00∶00∶00－23∶45∶00.采用對稱內插的方法，計算時取30min間隔的歷元作為內插點計算各個內插時間段中間時刻的衛星位置，以便將內插得到的衛星位置同精密星歷給出的衛星位置進行比較。在Matlab7.1環境下實現得到1號衛星7階～25階的內插結果，IGS向全球GPS用戶提供的GPS精密星歷，其誤差一般小于5cm，所以這幾種內插精度至少要小于5cm.

表1示出了Lagrange插值、Chebyshev擬合、三角函數插值的7階～25階的標準差，圖5～7示出了三種方法的標準差折線，圖8示出了這三種方法在X方向不同階次插值效果比較。

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

從表1可以看到，用Lagrange插值時，當階數取13時，內插精度最高，X，Y，Z3個坐標的最大標準差為3.23mm；用Chebyshev擬合時，階數為15時，內插精度最高，X，Y，Z3個坐標的最大標準差為3.45mm；用三角函數多項式插值時，當階數取15時，內插精度最好，X、Y、Z3個坐標的最大標準差為2.13mm.對比兩表可以發現，三角函數多項式插值，從7階到25階，各項誤差指標變化不大，但是Lagrange插值及Chebyshev擬合的標準差隨著階數的不同有變化相對較大。

在插值時只選擇了在中間時刻的插值，如選擇的插值位置不同，則會有不同的精度，用Lagrange插值及Chebyshev擬合衛星位置，無論選用多少階，內插的衛星位置圖呈現U型結構，即出現了最小值區域。在插值過程中，當插值點位于節點的兩端時，插值精度更低。

2.2 外推精度分析

選取PRN為1號的衛星在2012年8月6日23h45min 00s的衛星坐標為外推起點，推15 min，到24h00min 00s，與8月7日00h00min 00s的精密星歷坐標比較，表2分別示出了三種方法不同階次的外推誤差。

表2 三種不同方法的不同階次預測誤差

從表2可以看出，Lagrange插值的外推精度在13階時插值誤差最小，在24h00min 00s時的最大誤差達到19cm；Chebyshev的外推精度在13階插值誤差最小，在24h00min 00s最大誤差達到－18cm；在三角函數多項式插值的外推精度在15階時插值誤差最小，在24h00min 00s時的最大誤差達到－10cm.

通過以上分析，Lagrange多項式插值及Chebyshev擬合由于未考慮衛星的實際運動特點，當待插值點位于節點數據兩端時，內插精度會下降；而三角函數多項式插值由于考慮了衛星的實際運動特點，可以從一定程度上克服這一問題，外推精度要明顯優于前兩種方法，這是因為后者有效地利用了衛星在空間的運動特點，使得其插值函數更符合衛星的實際運動。

3 結 論

三種插值方法在對精密星歷進行內插時，只要選擇的插值階數合適以及插值點的位置合適，精度相當，都可以達到毫米級，因此無論采用哪種方法進行內插，都可以取得滿意的結果，但是三角函數多項式插值比Lagrange插值及Chebyshev擬合更穩定。在精密星歷外推時，三角函數多項式插值的精度要高于Lagrange插值及Chebyshev擬合，因此應避免用Lagrange插值及Chebyshev擬合進行外推。通過上文的誤差分析可知，僅采用觀測數據當天的精密星歷，利用三角函數多項式插值獲得觀測數據對應的精密衛星位置，其插值精度足以滿足一般的動態單點定位要求。

［1]劉偉平，郝金明，李作虎.由廣播星歷解算衛星位置、速度及精度分析［J].大地測量學地球動力學，2010，30（2）：144－147.

［2]劉大杰，陶本藻.實用測量數據處理方法［M].北京：測繪出版社，2004.

［3]孟大志，劉 偉.現代科學與工程計算［M].北京：高等教育出版社，2009.

［4]魏二虎，柴 華.GPS精密星歷插值方法比較研究［J].全球定位系統，2006，31（5）：13－15.

［5]孔巧麗.用Chebyshev多項式擬合GPS衛星精密坐標［J].測繪通報，2006，31（8）：1－3.

［6]孔祥元，郭際明，劉宗泉.大地測量學［M].武漢：武漢大學出版社，2008.

［7]張守建，李建成，邢樂林，等.兩種IGS精密星歷插值方法的比較分析［J].大地測量與地球動力學，2007，27（2）：80－83.

［8]張 勤，李家權.GPS測量原理及應用［M].北京：科學出版社，2009.