山东CORS坐标时间序列周期与噪声特征分析

席换，曲国庆

(山东理工大学建筑工程学院，山东淄博 255049)

0 引言

GPS坐标时间序列不仅反映坐标分量的整体趋势，还反映出基准站的非线性变化.相关研究表明，由于全球卫星导航系统（GNSS）基准站坐标时间序列受到地球物理因素及系统误差的影响，不仅存在线性变化，还存在周期变化，尤其在高程方向表现得最为明显[1-6].分析坐标时间序列的非线性变化特征，有利于研究基准站的稳定性和基于非线性运动模型获得基准站的瞬时位置.坐标时间序列的噪声信息会影响坐标解算精度，造成测站非线性运动模型有一定的偏差.GNSS坐标时间序列不仅存在白噪声（WN）还存在有色噪声（CN），目前最优随机模型广泛地被描述为WN+闪烁噪声（FN）[7-14].就山东部分连续运行参考站（CORS）而言，唐江森等[15]研究了两年坐标时间序列的噪声特征，表明其最佳噪声模型也为WN+FN.通常认为时间序列跨度大于2.5年，才能准确估计时间序列的线性项、周期项及其精度[16].

本文在已有研究的基础上，以山东CORS 2015-01-01—2018-12-31连续四年的坐标时间序列为研究对象，深入分析CORS站坐标时间序列的周期与噪声特征，以确定山东CORS站点的最佳噪声模型.

1 数据处理策略

选取了较为均匀分布在山东地区的CORS站，如图1所示.

图1 站点分布

受到外界条件、多路径效应或传输信号干扰等因素的影响，观测时间序列存在粗差，不仅影响数据质量，而且影响站点速度的计算及噪声特性的分析，有必要将其剔除.采用拉依达准则法（3σ准则）探测并剔除观测序列的粗差.外界观测条件、站点受损或网络故障等都会影响时间序列的连续性和完整性，加之剔除粗差的历元，坐标时间序列会出现在某天或某几天的观测数据缺失的现象，经过对多种插值方法反复实验对比，采用三次多项式插值补全数据.

坐标时间序列在较长时间内受某些因素的影响，形成一种整体变动的趋势，该趋势使得坐标时间序列为非平稳序列.对坐标时间序列进行零均值化处理，实际上就是平稳化的过程，其拟合模型如下：

式中:yt为坐标时间序列；10−8分别为坐标时间序列拟合的常数值和一次项系数值；10−9为去除趋势项的残差序列.

2 CORS站坐标时间序列分析

2.1 周期特征分析

小波分析方法是一种同时在时、频两域内分析时间序列的方法，实现坐标时间序列周期项的探测与提取.其基本原理如下：

式(2)为小波分解公式.式中：AK，DK分别表示分解后第K层的低频和高频部分；H、G分别表示小波低通和高通滤波器；N表示信号长度.

式(3)为小波重构公式.式中：H*、G*为H、G的共轭转置.

通过小波分析提取出QUFU站的年、半年和季节周期项，如图2所示.

由图2可知，QUFU站N、E、U三个坐标方向均表现出一定的周期特征，具有明显的年周期项和半年周期项变化，整体周期性变化较为稳定.对27个CORS站的周期特征进行分析可得：1）U方向的半年周期变化特征最为明显，其次是年周期特征，最后是3～4个月的季节周期特征；2）N、E方向具有明显的年周期或半年周期，不同站点，差异较大；3）坐标时间序列的周期信号并不是严格的年周期、半年周期等周期项，而是接近于年周期、半年周期等周期项.不同测站所表现出来的周期特征存在着差异，同一测站的不同方向坐标分量的周期规律也有所不同，这是因为测站所处的地理位置不同，环境因素成为影响它们周期规律不同的主要因素.

图2 QUFU站周期项

2.2 噪声特征分析

2.2.1 噪声序列的获取

GPS基准站的拟合模型[17]可用下式表示

式中：y(ti)为CORS站原始坐标时间序列；a为观测序列的初始位置；b为线性速度；ti为时间；c、d、e、f分别为年、半年周期项系数；gi为tgj时刻的阶跃；H为阶梯函数；vi为噪声序列.若顾及更复杂的情况，如某时刻测站运动速率的改变或某一事件发生后测站运动速率呈指数衰减的情况等，可采用文献[8]提供的坐标分量每日解观测序列参数模型.在数据处理策略中可实现坐标时间序列线性趋势项、阶跃的去除，再利用小波分析提取坐标时间序列的周期项，由此获得坐标分量的噪声序列.

2.2.2 谱指数与最大似然估计

噪声的功率谱密度P(f)与噪声频率f之间存在着某种幂次关系，即

式中，α为谱指数.式(5)两边取对数为

明显地，α在双对数直角坐标系中为拟合直线的斜率.不同的谱指数对应不同的噪声类型，α=0对应WN，α=−1对应FN，α=−2对应RWN.

最大似然估计（MLE）是确定不同噪声模型下噪声分量的方法，使得在相应噪声模型下噪声序列与其协方差的概率密度最大[18]，即

式中:X为时间序列列向量；σw,σrw,σf分别为WN、RWN和FN；N为序列长度；QXX为X的协方差.

根据式(6)谱指数的定义，计算求得CORS站各坐标分量的谱指数，如表1所示.

由表1可知，CORS站水平方向的谱指数均在−1～0之间，垂直方向除ZAZH站外，谱指数均在−1～0之间，说明CORS站各坐标分量噪声均不具有纯WN的特性，且同一测站不同方向的噪声模型有所不同.为确定CORS站坐标分量的最佳噪声模型，根据最大似然估计准则，以ZAZH站为例，计算了五种噪声模型：WN、WN+FN、WN+RWN、FN+RWN和WN+FN+RWN的MLE值，结果如表2所示.

由表2可知，不同噪声模型对应的MLE值不同，WN模型MLE值最小，说明坐标分量中不仅含有WN，还含有CN.水平方向上，WN+FN模型和WN+FN+RWN模型的MLE值最大且相等，说明在没有RWN模型存在的情况下，WN+FN模型与WN+FN+RWN模型所估计出来的噪声分量是相同的.因此，ZAZH站水平方向的最佳噪声模型为WN+FN或WN+FN+RWN；垂直方向上，MLE值为−4 202.15所对应的最佳噪声模型为WN+FN+RWN.由此可以看出同一CORS站的不同方向可以有不同的最佳噪声模型.假设CORS站的噪声模型为WN+FN+RWN，在此模型下计算了山东地区27个CORS站的噪声分量，如表3所示.

表2 ZAZH噪声模型M LE 值

表3 WN+FN+RWN噪声模型下的噪声分量估计值

由表3可知，观测时间序列坐标分量具有不同的噪声特征.N方向上，有25个站可用模型WN+FN描述，2个站可用模型WN+FN+RWN描述；E方向上，20个站可用模型WN+FN描述，7个站可用模型WN+FN+RWN描述；U方向上，16个站可用模型WN+FN描述，11个站可用模型WN+FN+RWN描述.因此山东CORS站的最佳噪声模型为：N方向上采用WN+FN模型，E和U方向上采用WN+FN+RWN模型.同时，U方向的噪声分量估计值要明显高于其它两个方向，这与高程方向精度低于水平方向精度的结论是一致的.

3 结束语

本文以山东27个CORS站的坐标时间序列为研究对象，通过小波分析方法提取周期项进行分析，得出了坐标时间序列中并不严格地存在年周期项、半年周期项等周期信号，而是与之接近的周年、半周年项等周期信号；另外，在观测序列的坐标分量中，U方向的周期信号较水平方向明显.

CORS站的谱指数表明坐标时间序列中不仅存在WN，还存在CN.根据MLE方法，估计了CORS站各坐标分量在WN+FN+RWN模型下的噪声分量，得出了CORS站各坐标分量可以有不同的噪声模型，N方向的最佳噪声模型为WN+FN，E、U方向的最佳噪声模型为WN+FN+RWN的结论.