多項式插值算法處理GNSS精密鐘差數據適應性研究*

曹可勁 馬恒超 朱銀兵 李 豹

（海軍工程大學電氣工程學院 武漢 430033）

1 引言

定位數據的后處理中，對精密鐘差產品進行插值計算是初始步驟，最常采用的是拉格朗日、牛頓和切比雪夫多項式插值法。目前IGS分析中心提供的最終產品，鐘差精度為0.1ns，優于廣播星歷約2個數量級［1，11]，折算成等效衛星測距誤差SISRE已下降到厘米級，鑒于多項式插值固有的特性，在此精度上進行插值計算應仔細處理。例如文獻［2]對北斗三種衛星軌道坐標采用了拉格朗日和牛頓插值法分別進行了計算，認為應該根據不同情況選擇不同階數；文獻［3]對拉格朗日、牛頓、滑動拉格朗日三種插值方式進行了對比分析，認為對于鐘差插值線性插值與多項式插值的精度差別不大；文獻［4]認為切比雪夫插值法精度比拉格朗日插值法高，但階數高了受影響比較大，精度變化比較劇烈。以上分析基于各種多項式插值的特點，但都沒有給出IGS數據插值處理的一般原則，本文選擇最具代表性的滑動拉格朗日插值法進行分析，找出基本規律和選擇的策略。

2 問題描述

下載IGS某天的30s、5min間隔的精密鐘差產品，以及5min和15min的精密星歷產品，分別做兩點間一階線性插值和多階滑動拉格朗日插值［4]：將5min間隔的精密鐘差數據內插生成30s間隔的數據，將15min間隔的數據內插生成5min間隔的數據。將內插得到的值與原始30s間隔的精密鐘差數據進行對比，分析內插精度。其中拉格朗日內插時，5min間隔內插成30s間隔采用通常的10階運算，15min間隔內插成5min間隔采用8階運算，取30個點進行評估。內插精度用誤差絕對值的平均值表示，結果如圖1、圖2所示。

圖1 5min內插到30秒精度比較

圖2 15min內插到5min精度比較

從以上計算發現，兩種插值過程中除個別衛星，絕大部分拉格朗日插值的精度都比直接線性插值精度低，我們對其它日期的數據進行計算，發現該現象普遍存在。

3 問題分析

3.1 時鐘漂移特征

目前導航衛星搭載有三種類型的原子鐘：銣鐘、氫鐘和銫鐘。從長期穩定性而言，三種類型的鐘是逐步提高的，但銣原子鐘在質量、體積、功耗等方面占有優勢，氫鐘在短期和中期穩定度指標方面占有優勢，銫鐘的準確度和漂移率指標綜合指標最好［5，12]。因此美國GPS采用了銫原子鐘和銣原子鐘結合的方式，歐盟的伽利略、俄羅斯的三代格洛納斯以及我國正在建設的北斗三號，采用銣原子鐘和被動型氫原子鐘相結合的方式［6]。

相對于軌道誤差，采用不同的原子鐘和搭配方式，衛星的時鐘誤差變化規律要復雜得多。為了普通用戶使用方便，廣播星歷中將時鐘漂移簡化為二次曲線規律，時鐘模型采用了一個2小時更新的簡單二階線性函數表示［7]：

其中系數af0、af1、af2以及參考時間toc隨著衛星導航電文更新。由于該模型主要考慮的衛星鐘差十分鐘級到小時級的變化特性，因此對時間參數t的敏感性很弱，造成時鐘預報誤差只能到納秒以上的精度。

3.2 精密鐘差數據特征

IGS提供的精密鐘差產品分為兩種格式：一種是以5min和30s間隔由單獨的時鐘文件給出，另一種以5min和15min間隔包含在精密星歷文件中給出。對比廣播星歷的時鐘漂移模型，將各精密時鐘數據畫成曲線，直觀考察變化規律。任意取2016年9月3日凌晨開始的15min間隔、5min間隔和30s間隔的數據分別進行分類，可以得到三種情況：

類型一：三種數據間隔的線性都比較好，且趨勢相同，如圖3所示。

類型二：隨著數據間隔減少，線性度下降，尤其30s間隔數據情況惡劣，如圖4所示。

類型三：所有間隔的數據都不呈線性，且趨勢不同，如圖5所示。

從分析的數據來看，隨著間隔縮小，鐘差線性度變差，相對抖動加大。從分析的數據觀察，絕大部分鐘差數據都屬于后兩種情況。

圖3 9號衛星鐘差曲線

圖4 5號衛星鐘差曲線

圖5 8號衛星鐘差曲線

3.3 多項式插值模型的適用性

各類多項式插值方法要保證精度，需要測量值符合某個函數的變化規律［8]，插值中所構造的插值函數應該盡可能逼近該函數。以典型的拉格朗日插值法為例，在給定的k+1個樣本點：(x0，y0)，…，(xk，yk)基礎上構造一個階數為n的多項式lk(x)，其中xj對應著自變量的位置，而yj對應著函數在這個位置的取值。假設任意兩個不同的xj都互不相同，那么拉格朗日插值多項式為

其中每個lk(x)為拉格朗日基本多項式（或稱為插值基函數），其表達式為

該多項式的特點是在xj上取值為1，在其它點上取值為0。由式中可知，拉格朗日插值曲線可以精確無誤地通過已知樣本點，同時保證了插值點也是線性連續的。雖然通常來說，拉格朗日插值多項式的階數越高，插值函數越精確，但在實際中容易造成數值劇烈變化，即所謂的“龍格現象”［9]，反而降低了插值精度。

對比上一節中的時鐘漂移特性曲線，這種以多項式擬合時鐘誤差曲線的方式，在小尺度間隔情況下，并不能準確代表導航衛星時鐘漂移特性，即相對于分鐘以上時間粒度，秒級間隔呈現的時鐘漂移隨機誤差，不能用分鐘級間隔的樣本數據等效表示。

4 插值模型選擇策略

4.1 數據預分析

首先對用于插值的樣本數據進行線性度分析，采用最小二乘構建標準線性曲線，對于給定樣本點 (xi，yi)，i=1，2，…，N，存在Q≥ 0，求一次函數y=a+bx，使總誤差：

即參數a，b應滿足此時

以得到的Q值作為參考，本文以每顆星10個連續樣本值計算Q值。

4.2 仿真計算

1）取2016年9月4日00時00分至24時00分的5分鐘間隔的GPS精密鐘差數據計算分析。計算各衛星的Q值，如圖6所示，并分別用一階線性相鄰內插法、多階滑動拉格朗日插值法得到30s間隔的鐘差數據，以已知30s間隔精密鐘差數據為參考值，計算各星鐘差隨時間變化的誤差情況，其中1階線性插值誤差和9階滑動拉格朗日插值誤差如圖7、8所示。文中使用的滑動拉格朗日插值法，只采用原樣本點中部的插值結果，該區間相對其他樣本點區間精度較高［10]，如使用9階，需10個樣本值，采用第5和6樣本值間的插值。

圖6 各GPS衛星鐘差的Q值

圖7 1階線性相鄰內插后鐘差的誤差

圖8 9階滑動拉格朗日插值后鐘差的誤差

2）Q值與各插值法精度的關系

根據1）中圖示，每顆星Q值及誤差均隨時間變化，計算各星Q值的平均值，可將所有GPS衛星分為三組，如表1所示。

表1 GPS衛星鐘差的分類情況

選擇9號、5號、8號星為各類型代表，比較線性插值與不同階滑動拉格朗日插值情況，其插值精度如表2所示。

表2 三種分類的插值計算誤差對比

從計算結果比較，類型一對應的衛星插值精度最好，隨著真值數據的線性程度下降，拉格朗日插值模型的計算誤差相對于線性插值模型先增大后減小，線性插值精度整體優于拉格朗日插值；對表中各衛星而言，隨著滑動拉格朗日插值階數提高，誤差逐漸增大。對其他更多數據進行計算，總體上呈現了該特點。

5 結語

拉格朗日算法在實際應用于IGS衛星時鐘改正數插值計算時，只有在真值數據線性度較好（Q值較小）時，才可能提高插值的靈敏度，而且要選擇合適的階數，不宜取得過高。隨著導航衛星原子鐘技術的精度和穩定性提升，目前IGS測量得到的衛星時鐘漂移越來越呈現線性特性，尤其時間在分鐘以上間隔時更是如此。而且，目前IGS提供的衛星時鐘改正數精度已經在十分位納秒等級［7～8]，無論采用哪種插值算法，精度都可以達到千分位納秒等級，對定位精度的影響區別很小。鑒于拉格朗日算法的計算復雜性，以及存在誤差放大的風險，在進行高精度定位時，采用最簡單的線性插值反而更可靠。