CN106066901B - 一种gnss自动化变形监测的基准点稳定性分析方法 - Google Patents

Abstract

本发明公开了一种GNSS自动化变形监测的基准点稳定性分析方法，该方法主要包含以下步骤：对GNSS观测数据进行解算，得到基准点的坐标残差时间序列；对基准点坐标残差时间序列中的粗差进行探测和剔除；提取并剔除序列中的线性趋势项和周期项；若剔除趋势及周期项后的坐标残差序列有缺失数据，则对时间序列进行插值；提取出残差时间序列的共模误差，并从各残差序列中剔除该共模误差；最后对所剩残差序列用量化计算法判定该序列是否服从正态分布，若是则判定基准点稳定，否则不稳定。本发明可以有效地对站点共性的误差进行剔除，提高坐标时间序列的信噪比，最终使得站点的坐标时间序列主要误差仅为观测随机噪声，从而判断出站点是否稳定。

Description

一种GNSS自动化变形监测的基准点稳定性分析方法

技术领域

本发明涉及一种GNSS自动化变形监测的基准点稳定性分析方法。

背景技术

GNSS具有覆盖范围广、不受气候条件限制、无需通视条件、自动化程度高等优势，因此，近年来GNSS已经广泛应用在变形监测领域，而且获得了令人满意的结果。

GNSS自动化变形监测系统中通常选择1～2个稳定点作为基准点，但如何判断基准点是否稳定，目前仍没有一套具体的方法。在实际GNSS监测过程中，由于板块运动、固体潮、土壤水等质量负荷等因素的影响，使站点坐标会产生相应误差，这种误差往往有几毫米甚至十几毫米，严重影响了对测站本身位置变化情况的判断，现有技术常采用复杂的建模方法减小误差，计算量大，且不能有效地对误差进行剔除。

发明内容

本发明的目的在于提出一种GNSS自动化变形监测的基准点稳定性分析方法。

为解决上述技术问题，本发明所采用的技术方案是：一种GNSS自动化变形监测的基准点稳定性分析方法，包括以下步骤：

1)对站点GNSS观测数据进行解算，得到站点坐标残差时间序列；

2)对站点坐标残差时间序列中的粗差进行探测和剔除；

3)提取并剔除站点坐标残差时间序列的线性趋势项和周期项；若线性趋势项与所在板块地壳运动速率不一致，则判定当前基准点不稳定；否则，进入步骤4)；

4)若剔除线性趋势项及周期项后的站点坐标残差时间序列有缺失数据，则对剔除线性趋势项和周期项后的站点坐标残差时间序列进行插值；

5)对插值后的站点坐标残差时间序列进行时空滤波，提取并剔除插值后的站点坐标残差时间序列的共模误差；

6)判定剔除共模误差后的站点坐标残差时间序列是否服从正态分布，若是，则判定基准点稳定，否则不稳定。

步骤2)中，粗差探测和剔除的公式为：

其中，v_i表示站点坐标残差时间序列中i天的残差值(v_i是步骤1得到的残差值，如果v_i满足此不等式，则认为这一天的结果是粗差，删除这一天的结果)，w表示样本窗口长度，表示以i为样本中心，取站点坐标残差时间序列中左右两侧各天的值作为样本，IQR为四分位数间距，即将一段序列按照由小到大的顺序排列，站点坐标残差时间序列75％处的值与站点坐标残差时间序列25％的值之差。

根据拟合函数，利用最小二乘法对线性趋势项和周期项进行剔除，拟合函数y(t_i)公式如下：

其中，a表示站点位置，b表示站点线速度，c、d表示年周期运动系数，e、f表示半年周期运动系数，t_i表示以年为单位的GPS历元，g表示发生在历元T_g处的n_g个阶跃式偏移量；H为阶梯函数，若在时刻T_gj处发生突变，则在此之前H＝0，在此之后H＝1。

与现有技术相比，本发明的有益效果为：本发明可以有效地对这些误差进行剔除，提高坐标时间序列的信噪比，使得站点的坐标时间序列主要误差仅为观测随机噪声，从而判断出站点是否稳定。

附图说明

图1为坐标残差时间序列；

图2为剔除线性趋势项和周期项后的残差序列；

图3为剔除共模误差后的残差时间序列；

图4为剩余残差统计直方图；

图5为本发明方法流程图。

具体实施方式

本发明方法流程如图1所示，其包括以下具体方法步骤：

1)获取基准点RINEX格式的观测数据文件，同时下载相应的IGS数据产品；观测数据文件和IGS数据产品是坐标解算的输入文件；

2)利用获得的数据文件及数据产品，采用精密单点定位方式对基准点坐标进行解算，得到基准点单日解的坐标时间序列；

3)采用四分位数间距判别准则进行粗差的探测和剔除。其公式如下：

其中，v_i表示残差序列中i天的残差值，w表示样本窗口长度，表示以i为样本中心，取序列中左右两侧各天值作为样本，IQR为四分位数间距，即将一段序列按照由小到大的顺序排列，序列75％处的值与序列25％的值之差；这种方法剔除粗差不会受到个别极大值或极小值的影响；

4)采用拟合函数模型，利用最小二乘法对各参数进行估计，剔除线性趋势项和周期项；若线性趋势项与所在板块地壳运动速率不一致，则即可判断该基准点不稳定，否则继续后面的步骤进一步判断；拟合函数公式如下：

其中，a表示站点位置，b表示站点线速度，c、d表示年周期运动系数，e、f表示半年周期运动系数，t_i表示以年为单位的GPS历元，g表示发生在历元T_g处的n_g个阶跃式偏移量，是由天线位置变化或地震影响等因素造成的站点位置突变。H为Heaviside阶梯函数，假设在时刻T_gj处发生突变，则在此之前H＝0，在此之后H＝1；这种估计方法是对参数的线性无偏估计，并且利于编程实现；

5)若剔除趋势及周期项后的坐标残差序列有缺失数据，则使用线性插值法和三次样条插值法对时间序列进行插值，其中对于序列中小于三天的空隙采用三次样条插值，大于三天的空隙使用线性插值；采用这种插值方法，这样既可以使插值后的坐标时间序列足够准确，也可以保持其整体的趋势不变，插值后的数据也方便了后续的主成分分析计算；

6)使用主成分分析法对基准点坐标时间序列进行时空滤波，提取出残差时间序列的共模误差，并从各残差序列中剔除该共模误差除；该方法可以将测站共同的误差进行剔除，剩下线性无关的误差；

7)对所剩残差序列用量化计算法判定其是否服从正态分布，若是则判定基准点稳定，否则不稳定。

广州南沙单双频GPS混合地面沉降监测系统在一片约100km²区域内布设了11个GPS观测点，分别命名为GD01至GD11。其中，四个为双频点，分别为GD03,GD05,GD08,GD10，其余7个点是单频点。GD05点位于基岩上，作为稳定的基准点，GD08点位于有深基桩的建筑物屋顶上，是一个相对稳定的点，其余点位于地面或无深基桩的屋顶上。因此，需要对GD08点的稳定性进行研究，判断其是否能和GD05点一起作为区域内的稳定基准点。

按照上述方法，对GD08点进行稳定性分析。获取了系统从2012年6月至2013年9月采集的GD05和GD08点的RINEX格式观测文件，并下载相应的IGS数据产品，利用GPS解算软件进行坐标解算。解算结果如图1所示。

对解算结果采用四分位数间距判别准则进行粗差的探测和剔除，之后利用拟合函数模型提取残差序列的线性趋势项和周期项。结果如下表1所示：

表1残差序列线性趋势项和周期项结果

由上表可以看出，GD05和GD08在三个方向上的位移速率都很接近，N、E、U方向分别只有1.21mm/yr,0.56mm/yr,0.94mm/yr的差距，因此可以初步判断两个点有着相同的位移趋势。

采用线性插值法和三次样条插值法对时间序列进行插值，补全残差时间序列中的空缺。之后，采用主成分分析法对GD05和GD08点坐标时间序列进行滤波，剔除两个序列的共模误差，剔除后的结果见图3。

对所剩残差进行量化统计，结果如图4：

通过量化法计算，当取置信水平等于0.95时，GD05和GD08点三个方向的剩余残差时间序列均满足正态分布。通过上述分析，可以看出，经过一系列的处理与分析，GD05和GD08点坐标残差时间序列最终呈正态分布，也就是说，GD05和GD08点均不存在地表位移，从而证明了GD08可以和GD05一样作为整个区域网的稳定基准点。

Claims (2)

1.一种GNSS自动化变形监测的基准点稳定性分析方法，其特征在于，包括以下步骤：

1)对站点GNSS观测数据进行解算，得到站点坐标残差时间序列；

2)对站点坐标残差时间序列中的粗差进行探测和剔除；粗差探测和剔除的公式为：

其中，v_i表示站点坐标残差时间序列中i天的残差值，w表示样本窗口长度，表示以i为样本中心，取站点坐标残差时间序列中左右两侧各天的值作为样本，IQR为四分位数间距，即将一段序列按照由小到大的顺序排列，站点坐标残差时间序列75％处的值与站点坐标残差时间序列25％的值之差；

3)提取并剔除站点坐标残差时间序列的线性趋势项和周期项；若线性趋势项与所在板块地壳运动速率不一致，则判定当前基准点不稳定；否则，进入步骤4)；

4)若剔除线性趋势项及周期项后的站点坐标残差时间序列有缺失数据，则对剔除线性趋势项和周期项后的站点坐标残差时间序列进行插值；

5)对插值后的站点坐标残差时间序列进行时空滤波，提取并剔除插值后的站点坐标残差时间序列的共模误差；

6)判定剔除共模误差后的站点坐标残差时间序列是否服从正态分布，若是，则判定基准点稳定，否则不稳定。

2.根据权利要求1所述的GNSS自动化变形监测的基准点稳定性分析方法，其特征在于，根据拟合函数，利用最小二乘法对线性趋势项和周期项进行剔除，拟合函数y(t_i)公式如下：

其中，a表示站点位置，b表示站点线速度，c、d表示年周期运动系数，e、f表示半年周期运动系数，t_i表示以年为单位的GPS历元，g表示发生在历元T_g处的n_g个阶跃式偏移量；H为阶梯函数，若在时刻T_gj处发生突变，则在此之前H＝0，在此之后H＝1。