精密卫星钟差内插的三种方法及精度分析＊

王 俊，方书山

（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．