小区域GNSS高程异常拟合方法研究

谢萌丽,姜永涛,刘国仕

(1.南阳师范学院,地理科学与旅游学院 测绘工程，河南 南阳 473000;2.昆明理工大学,国土资源工程学院 资源勘查工程与技术，云南 昆明 650031;3.邵阳学院,湖南 邵阳 422000)

0 引 言

我国工程建设普遍采用正常高系统,而正常高一般采用等级水准测量方式获取,耗时、繁琐且效率低．全球卫星导航系统(GNSS)定位技术可以快速、精确地获取大地高[1],但是,由于未能获得同等精度的高程异常值,导致转换的正常高达不到高程精度要求[2],制约了GNSS高程定位技术的应用．

区域高程异常值可通过对已知GNSS水准点的高程异常拟合获取,常用方法有等值线图法[3]、解析内插法[4]、曲面拟合法[5-6]、加权平均法[7]、非格网GNSS散点数据考虑地形改正法[8]和BP神经网络法[9]等．研究显示,在高程异常变化比较平缓大区域,可首选二次多项式或双线性内插高程异常拟合方法[10-11];若已知高程异常点呈直线排列可采用三次样条曲线、多项式曲线拟合法[12]．现有文献主要针对大面积区域高程异常拟合研究,缺乏对小区域高程异常拟合方法和精度的探讨．为此,本文利用湖南省新化县四等水准测量和E级GNSS测量成果,研究县城区域尺度GNSS高程异常拟合方法及其精度,对于小区域高程拟合方法的选取有借鉴意义．

1 数据与方法

1.1 数据来源

对湖南省新化县分别进行四等水准测量及E级GNSS测量,获取控制点平面坐标、大地高及正常高数据,数据及控制点分布如图1和表1所示．

图1 GNSS水准点点位分布

表1 GNSS水准起算点和检核点

1.2 空间插值原理

1.2.1 反距离加权插值法

反距离加权插值法假定每个已知数据点都会对待定点产生一种局部影响,这种影响随着距离的增大而减小．反距离加权高程异常插值实质是由插值点周围已知点的高程异常值加权平均求算的[7]．

设在插值点周围选m个已知点,高程异常为ζi(i=1,2,…，m),则待定点的反距离加权高程异常值为

(1)

1.2.2 多项式插值法

多项式插值法是由已知高程异常数据拟合出一个由多项式函数定义的平滑曲面,并根据该函数求算待定点高程异常值．通常多项式曲面无法通过每个已知点,一般利用最小二乘方法获得最佳拟合结果．

常用的二次曲面拟合法,公式如下:

ζ(x,y)=a0+a1x+a2y+a3x2+

a4y2+a5xy,

(2)

式中：ai(i=0,1,…，5)为待定系数；(x,y)为高程异常点的平面坐标．

由已知数据可根据最小二乘原理求算公式(2)中的拟合系数ai,而后可以求得待定点处的高程异常值．

全局多项式插值法是利用所有数据进行多项式拟合,阶数视数据表面包含的弯曲个数而定;而局部多项式插值法是对位于指定重叠邻域内的数据分别进行多项式拟合,通过样本属性对搜索邻域进行定义,且邻域相互重叠．全局多项式插值法适用于分析数据长波变化,而局部多项式插值法适用于分析短波信息．

1.2.3 径向基函数插值法

径向基函数(RBF)插值法是在多维空间中寻找一个能够最佳匹配已知高程异常数据的曲面,然后对待定点进行插值,获取高程异常值．RBF为一系列精确插值法的组合,其拟合曲面通过每一个已知点,且曲面总曲率最小[13],拟合函数为

(3)

式中：φ(di)(i=1,2,…,n)为n个径向基函数；wi为权．

RBF中权wi的计算,是根据采用的RBF,由高程异常拟合数据集(ζi,di)(i=1,2,3，…，n)求算．其中,di为待定点到已知数据点ζi的距离．

1.2.4 核平滑障碍插值法

核平滑障碍插值法是将每个已知点处的局部趋势融合形成预测面,再与待定点相匹配进行插值[14]．核平滑障碍插值法可看做局部多项式插值法的变型,用于解决局部不稳定问题．其优点在于通过使用岭参数将少量偏差引入方程,解决大预测标准误差和失真预测的问题．但岭参数应尽可能小,以维持拟合模型的稳定性[14]．

1.2.5 克里金插值法

克里金插值是根据待定点周围有限邻域内的若干已知点数据,在考虑已知点与待定点的空间位置关系以及变异函数提供的结构信息,对待定点进行的一种线性无偏最优估计[15-16]．克里金插值公式如下:

(4)

(5)

克里金插值的适用条件是区域变量存在空间相关性,在使用时须检查数据是否存在正态分布、是否平稳、是否有趋势项,若存在趋势,需要对趋势进行多项式移除后,才能进行空间插值．

普通克里金插值法适用于满足二阶平稳假设的区域变化量,而泛克里格插值法适用于非平稳或具有漂移存在的区域变化量．小区域高程异常满足二阶平稳变化条件,且数据存在两种趋势,西北高东南低,因此在使用泛克里金插值法时应进行一阶和二阶移除比较．

1.3 实验方法

如图2所示，大地高H与正常高Hγ之差为高程异常ζ,三者之间关系为

H=Hγ+ζ.

(6)

图2 大地高与正常高关系

似大地水准面(ζ)为一个光滑的曲面,可采用空间插值方法(数学曲面)进行高程拟合研究．为研究县城区域尺度不同GNSS高程异常拟合方法及其精度,本文实施方案如下:

1) 确定数据点和检核点．为详细研究不同高程拟合方法的内插和外推精度,本文对GNSS水准点进行数据点和检核点的划分,如表1和图1所示,共选取16个数据点和8个检核点．检核点中,有4个在数据点范围内,用于评定内插精度;4个检核点位于已知数据点外,用于评定外推精度．

2) 利用ArcGIS地统计工具,进行GNSS水准高程拟合方法研究及精度评定．

ArcGIS地统计工具中的插值方法可分为空间确定性插值和地统计插值．前者包括反距离加权、多项式插值法、RBF插值以及核平滑障碍插值法,后者包括普通克里金、泛克里金等插值方法[14,17]．对GNSS水准高程拟合方法及其精度评定过程如下:

1)选用不同插值方法,利用GNSS数据点进行空间插值,获取插值拟合函数;

2)根据拟合函数,求算检核点的高程异常拟合值;

3)统计检核点高程异常拟合值与实测值的残差,采用中误差指标,分析和探讨不同插值方法内插和外推精度．

2 结果与分析

2.1 空间插值结果

对湖南省新化县GNSS水准数据,不同拟合方法的高程异常曲面结果如图3所示.

(a)反距离权重插值法 (b)全局线性多项式插值法

(c)全局二次多项式插值法 (d)全局三次多项式插值法

(e)径向基函数插值法 (f)局部多项式插值法

(g)核平滑障碍插值法 (h) 普通克里金插值法

(i)泛克里金插值法(一阶移除) (j) 泛克里金插值法(二阶移除)

图3 不同拟合方法的高程异常结果

图3可以看出,湖南省新化县高程异常整体呈现西北高,东南低的趋势．不同的空间插值方法得到的高程异常曲面(似大地水准面)不同,相同拟合方法不同参数的结果也具有较大差异．

2.2 拟合精度评定

本文对高程异常拟合精度评定中,1～4点为内插点,分别对应E009、E013、E015、E016,5～8为外推点,分别对应E006、E007、E024、E025．根据检核点高程异常拟合值与实测值的差值,得到的精度统计结果如表2和图4所示.

表2 不同插值法的内插与外推精度 m

如图4和表2,在县域尺度高程拟合方法中,内插精度最高的为RBF插值法,精度为±0.0158 m,其次为泛克里金二阶插值法、泛克里金一阶插值法和全局二阶多项式插值法等．径向基函数为一系列精确插值法的组合,其拟合曲面通过每一个已知点,且曲面总曲率最小,因此内插精度最高．

外推精度最高的为局部多项式插值法,精度为±0.0104 m,其次为全局二阶多项式、克里金插值和全局线性多项式插值．在多项式插值中,多项式阶数应依据数据趋势而定,对于县域尺度GNSS高程异常数据,由于局部多项式插值法更能体现曲面的短波变化,因此,插值精度较高．

其余插值方法,如反距离权插值法是按距离定权,建立函数关系进行拟合和预测的,该方法内插精度高于外推精度,但整体拟合精度不高．克里金插值法外推精度整体高于内插精度,原因是克里金插值可反映高程异常的整体曲面特征,但是在使用泛克里金二阶插值法时,点8(外推点)无法生成预测值,说明该方法有一定局限性．在克里金插值方法中,普通克里金插值精度高于泛克里金一阶插值,这是由于普通克里金插值法适用于满足二阶平稳假设的数据,而泛克里金插值法适用于非平稳或具有漂移存在的区域变化量[17]．

在本文研究插值方法中,RBF插值方法内插精度最高,其中误差为±0.0158 m,局部多项式插值法外推精度最高,其中误差为±0.0104 m;综合来说,局部多项式插值法整体精度最高．因此,在小区域高程异常拟合中,可以优先考虑局部多项式插值法．

3 结 论

本文采用湖南省新化县城区四等水准测量和E级GNSS测量数据,利用ArcGIS地统计工具,研究小区域GNSS水准高程拟合方法及精度,主要结论有: 1)在本文研究的插值方法中,径向基函数插值法内插精度最高,其中误差为±0.0158 m,局部多项式插值法外推精度最高,精度为±0.0104 m;综合来说,局部多项式插值法整体精度最高;2)在多项式插值中,局部多项式插值法更能体现数据的短波变化,在小区域GNSS高程异常拟合中内插和外推精度较高；3)普通克里金插值法适用于满足二阶平稳假设的数据,由于县域尺度高程异常变化属于平稳特征信号,因此普通克里金插值法在克里金插值中精度最高．本文研究结果对于小区域高程拟合方法的选取有一定的借鉴意义．