-
地壳垂直运动是地壳运动的一种重要表现形式,长期观测成果可直接服务环境预测、防灾减灾,将极大推动大陆动力学、区域应力场以及地壳运动等方面的研究,为国家与地区的基础设施建设提供垂直运动信息。大地测量是定量描述地壳垂直运动的重要方法之一,在国内外得到广泛应用。
文献[1-5]利用全球导航卫星系统(global navigation satellite system,GNSS)、水准等技术获取垂直形变信息,分析了区域地壳垂直运动趋势。GNSS、水准观测数据源于不同观测仪器、误差模型、处理方法、参考基准,这导致GNSS、水准的垂直形变结果可能自身存在误差,且相互之间存在系统偏差,进而影响垂直形变分析的精度。另外水准网一般沿公路、铁路等路线布设,水准网中间存在较大空白区,导致其垂直形变量的可靠性不高。而GNSS控制网布设一般考虑网型结构,按照一定密度布设控制点,这些GNSS控制点能够进一步填补水准网的空白,有效提升垂直形变信息的可靠性。随着多期GNSS、水准观测数据的积累,GNSS水准点具有同时期的GNSS与水准观测数据,这有利于分析大地高变化与正常高变化之间的系统偏差,为实现两类高程变化融合统一提供了基础数据。为了充分发挥不同类型数据的优势,杨元喜等[6-7]深入研究了大地测量数据融合理论与方法,提出了函数模型、随机模型方面的融合方法,充分发挥大地测量数据对计算结果的贡献。对于GNSS、水准数据融合,文献[8-14]研究了联合平差的速率基准、观测值选取、权重调整、函数模型选取等问题,充分发挥不同类观测值对垂直运动速度的贡献,提高垂直运动速度的估计精度。GNSS与水准数据可以通过联合平差方法实现融合,也可以使用GNSS水准点多期成果,利用拟合逼近方法实现GNSS与水准数据融合。因此,基于拟合逼近原理,本文提出一种自适应最小二乘配置的GNSS与水准融合方法,该方法以水准点获取的正常高变化速度为基础,利用GNSS水准点建立大地高变化速度与正常高变化速度的转换关系,进一步修正GNSS点的大地高变化速度后,再根据融合后的高程变化速度结果,采用反距离加权法建立垂直运动格网模型。
-
精密水准测量是垂直运动监测的经典技术手段,国内外学者通常将GNSS、星载合成孔径雷达干涉测量(interferometric synthetic aperture radar,InSAR)垂直运动结果与水准垂直运动结果相比较,检验垂直运动成果的可靠性。若建立控制点与GNSS水准点两类高程变化速度之间差异的函数关系,就能实现GNSS与水准两类垂直运动速度的融合。对于函数关系建模,可以采用多面函数[15]、最小二乘配置[16]、神经网络[17]等拟合逼近方法。多面函数能够有效逼近曲面,无法有效抑制系统误差影响。神经网络理论能够削弱不确定因素等影响,但容易陷入极小值,且泛化能力差。最小二乘配置能够较好地解决两类高程变化量之间差异中含有倾向性、随机性因素的问题,控制系统误差等因素影响,提高两类高程变化量之间差异的拟合效果,但协方差矩阵一般根据经验确定。自适应最小二乘配置[16, 18]利用自适应因子调整观测噪声协方差矩阵与信号协方差矩阵的比例,使其更加接近真实关系。因此,将水准正常高变化量作为基础,不考虑垂线偏差造成的大地高与正常高变化速度差异,利用自适应最小二乘配置方法建立GNSS水准点的大地高与正常高速度之间的转换关系,进一步实现非GNSS水准点的垂直运动速度融合。
最小二乘配置的大地高与正常高的速度差异观测量分为固定参数效应和随机效应两部分,前者表示倾向性部分或系统性部分,后者为随机性影响部分。其公式可表示为:
$$\mathit{\boldsymbol{L}} = \mathit{\boldsymbol{AX}} + \mathit{\boldsymbol{BS}} + \mathit{\boldsymbol{\varepsilon }}$$ (1) 式中,L为大地高变化速度与正常高变化速度之间的差异;A、B分别表示系统性参数以及信号参数对应的系数矩阵,维数分别为n × tn、n × ts, n为观测值个数,tn为系统性参数个数,ts为信号参数个数;X为倾向性参数向量,维数为tn × 1;S为信号向量,维数为ts × 1。
最小二乘配置估计值$\mathit{\boldsymbol{\widehat X}}$、信号估计值$\mathit{\boldsymbol{\hat S}}$与未知点信号估计值$ \mathit{\boldsymbol{\hat S'}}$分别为:
$$ \left\{ {\begin{array}{*{20}{l}} {\mathit{\boldsymbol{\widehat X}} = {{\left( {{\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{L}}^{ - 1}\mathit{\boldsymbol{A}}} \right)}^{ - 1}}{\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{L}}^{ - 1}\mathit{\boldsymbol{L}}}\\ {\mathit{\boldsymbol{\widehat S}} = {\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{S}}}{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{L}}^{ - 1}(\mathit{\boldsymbol{L}} - \mathit{\boldsymbol{A\hat X}})}\\ {{{\mathit{\boldsymbol{\widehat S}}}^\prime } = {\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_{\mathit{\boldsymbol{S'S}}}}\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{S}}^{ - 1}\mathit{\boldsymbol{\hat S}}} \end{array}} \right. $$ (2) 式中,ΣL为观测值协方差矩阵;ΣS为信号的协方差矩阵;ΣS′S为已知点与未知点的互协方差矩阵。
$$ {\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{L}}} = \mathit{\boldsymbol{B}}{\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{S}}}{\mathit{\boldsymbol{B}}^{\rm{T}}} + {\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\mathit{\boldsymbol{\varepsilon }} } $$ (3) 式中,Σε为观测噪声的协方差矩阵。
$$ \alpha=\frac{\hat{\sigma}_{0 S}^{2}(i)}{\hat{\sigma}_{0 \varepsilon}^{2}(i)} $$ (4) 式中,$ \hat{\sigma}_{0 \varepsilon}^{2}(i) $为观测噪声的第i次验后单位权方差;$ \hat{\sigma}_{0 S}^{2}(i) $为信号的第i次验后单位权方差。自适应因子α具体计算方法可以参考文献[16, 18-19]。
观测噪声协方差矩阵$\mathit{\boldsymbol{ \boldsymbol{\widehat \varSigma} }}_\varepsilon ^{i + 1} $为:
$$ \mathit{\boldsymbol{ \boldsymbol{\widehat \varSigma} }}_\varepsilon ^{i + 1} = \frac{1}{\alpha }\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_\varepsilon ^i $$ (5) 如此反复迭代,直至α趋于1。
-
目前,我国局部区域出现了严重的非均匀沉降,这也给精确描绘区域垂直运动规律增加了难度。若采用某曲面函数整体拟合大区域垂直运动,可能造成垂直运动剧烈区域存在较大偏差。因此,为了更加可靠地反映垂直运动信息,本文采用水准点以及融合后的GNSS点垂直运动速度作为拟合量,充分利用垂直运动剧烈特征点变化信息,使用反距离加权法实现格网化,建立真实可靠的高分辨率垂直运动格网模型。
反距离加权法[20-21]是一种整体计算方法,根据待定格网点与控制点之间距离确定权重。该方法依据各点与待定点的距离平方大小确定权值,具体公式如下:
$$ P_{k}=\left\{\begin{array}{ll} {\frac{1}{R^{2}}, } & {S_{k} \leqslant R} \\ {\frac{1}{S_{k}^{2}}, } & {S_{k}>R} \end{array}\right. $$ (6) 式中,R为搜索半径;Sk为第k点与待定点之间的距离。
格网点垂直运动速度Vi,j为:
$$ {V_{i, j}} = \frac{{\sum\limits_{k = 1}^N {{P_k}} V_k^\prime }}{{\sum\limits_{k = 1}^N {{P_k}} }} $$ (7) 式中,N为控制点总个数;Vk′为第k点垂直运动速度;i、j分别为格网点的行、列号。
-
算例采用2006、2015年山东省似大地水准面精化GNSS和水准观测数据。水准网观测按照二等水准网观测纲要分别于2004—2006年、2015年进行观测,GNSS控制网按照C级网观测纲要也于2004—2006年、2015年同步进行数据采集。水准点速度是两期似大地水准面精化水准网成果之差除以时段长度,两期水准网平差都采用青岛原点为基准点,且假设青岛原点垂直运动速度为0,山东省区域两期水准网重合点共有490个。两期GNSS控制网重合153个GNSS点,且都为GNSS水准点,水准点与GNSS点分布如图 1所示。为了实现GNSS水准高程变化速度有效融合,采用如下步骤进行处理:(1)使用GAMIT/GLOBK软件处理,选取相同的起算基准点,获得GNSS点在国际地球参考框架ITRF2008下不同历元的成果,计算GNSS控制点的垂直运动速度。(2)采用水准网动态平差软件获取水准网重合点垂直运动速度。(3)选取大地高变化与正常高变化一致的GNSS水准点,建立正常高变化速度与大地高变化速度的函数关系,实现正常高变化速度与大地高变化速度的融合。(4)利用融合后的高程变化速度建立垂直运动模型。
-
为了验证算法可靠性,算例从153个GNSS点中选取60个点作为融合点(见图 2),选取原则如下:(1)大地高变化速度与正常高变化速度方向尽量一致。(2)大地高变化速度与正常高变化速度之间差异的绝对值不大于1.5 mm/a。对于GNSS、水准高程变化融合,从60个融合点中,选取40个点作为拟合点,20个点作为检核点,将大地高变化速度与正常高变化速度差异作为拟合量,采用3种方案进行计算,方案1:多面函数;方案2:最小二乘配置;方案3:自适应最小二乘配置。
拟合点与检核点误差统计见表 1,拟合点残差见图 3,检核点残差见图 4。
表 1 速度融合结果精度统计
Table 1. Accuracy Statistics of Velocity Fusion Results
方案 拟合点残差/(mm·a-1) 检核点残差/(mm·a-1) 最大值 最小值 均方根误差 最大值 最小值 均方根误差 1 0.71 -0.69 0.36 1.62 -0.84 0.72 2 0.77 -0.58 0.33 1.34 -0.49 0.69 3 0.80 -0.63 0.35 1.33 -0.49 0.68 由表 1、图 3、图 4的计算结果可以看出,3种方案都能有效减小GNSS高程变化速度与水准高程变化速度之间的差异,实现GNSS高程变化速度与水准高程变化速度之间的有效融合。自适应最小二乘配置比多面函数方法拟合精度与检核精度有所提高。另外,在数据准备过程中,要考虑GNSS与水准的起算基准、有色噪声等因素影响。
为了建立高精度区域垂直运动模型,本文实验利用了60个融合点作为拟合点,采用自适应最小二乘配置建立融合关系,获取剩余93个GNSS点的垂直运动速度。
-
利用490个水准点及转换后的93个GNSS点,即583个高程变化点,采用反距离加权方法获取垂直运动格网模型,其中搜索半径R为1 km,山东省区域垂直运动趋势如图 5所示。
图 5 山东省2006—2015年垂直运动等值线图
Figure 5. Contours of Vertical Movement Velocities in the Shandong Province (2006―2015)
从图 5中可以看出:
1)近十年来,山东垂直运动总体上呈西降东升、北降南升的趋势。西北、西南部地区沉降严重,整体呈下降趋势,垂直运动速度平均20.0 mm/a;德州、菏泽、广饶地区出现3个沉降中心,其中广饶地区沉降最严重,十年平均沉降速度70.0 mm/a左右,这主要是由于地下水开采造成的;而在胶东半岛呈缓慢的上升趋势,垂直运动速度约5.0 mm/a。
2)采用反距离加权方法构建垂直运动速度模型,最大限度地发挥距离较近GNSS点与水准点对格网点的贡献,提高垂直运动模型的真实性与精度。
3)水准网一般沿公路布设,其包围区域内为空白区。分布均匀的GNSS网点能够有效补充水准网空白区,提高垂直运动建模所需的控制点分辨率,提升垂直运动模型的可靠性。
-
GNSS与水准是地壳垂直运动速度研究的重要技术。利用自适应最小二乘配置方法建立GNSS水准点的大地高变化速度与正常高变化速度之间的融合模型,可实现GNSS控制点的大地高变化速度与正常高变化速度的统一,有效发挥GNSS高时间分辨率和水准的高精度优势。本文利用反距离加权法建立了山东省高精度、现势性强的垂直运动速度模型,可以为省内基础工程建设以及地面沉降灾害防治等提供决策依据。山东省垂直运动趋势整体上西部、北部下降,其平均下降速度达20.0 mm/a左右,另外在广饶、德州、菏泽区域存在几个沉降中心,特别是广饶地区沉降速度达到70.0 mm/a左右,而在东部、东南部垂直运动趋势为上升,这些不均匀地面沉降对该区域内的基础设施建设造成了严重威胁,大大降低了高程基准成果的可靠性和现势性。随着各省市连续运行参考站(continuously operating reference stations,CORS)系统的进一步完善,利用多期水准与连续GNSS观测数据,削弱各类数据之间的系统偏差,采用自适应融合方法实现多源数据有效融合,建立国家、区域地壳垂直运动速度模型,维护高程基准,是当前大地测量的重要研究方向。
Analysis of Vertical Deformation with the Adaptive Fusion of GNSS and Leveling Elevation Variation
-
摘要: 全球导航卫星系统(global navigation satellite system,GNSS)与水准是获取垂直形变量的重要手段,由于不同的仪器、误差模型、数据处理方法等因素,GNSS、水准两类高程变化结果之间存在系统偏差问题,由此提出了一种GNSS与水准高程变化的自适应融合方法。首先利用自适应最小二乘配置建立GNSS与水准重合点两类高程变化量之间的转换模型;然后利用转换模型修正其他GNSS控制点的高程变化量;最后利用反距离加权法建立区域高程变化速度格网模型。实验以山东省水准网与GNSS控制网数据为例,利用高程变化较一致的60个GNSS水准点建立高程变化转换模型,实现了GNSS与水准高程变化的有效融合,削弱了GNSS大地高变化与正常高变化之间的系统偏差,提高了山东省区域高程变化分析的可靠性。Abstract: Global navigation satellite system(GNSS) and leveling are two important approaches to obtain the numerical vertical deformation. Due to different instruments, error models, data processing methods and other factors, systematic deviations exist between the results of elevation variation derived from the two approaches. Therefore, an adaptive fusion method of GNSS and leveling elevation variation is proposed. Firstly, the adaptive least-square configuration is used to establish the transition model for elevation variation between GNSS and leveling coincidence points. Then, we use the model to correct the elevation variation of other GNSS control points. Finally, the inverse distance weighting algorithm is used to establish the regional elevation change velocity grid model. The experiment takes the data of Shandong Provincial leveling network and GNSS control network as an example. The elevation change transition model is established by using 60 GNSS leveling points with relatively elevation changes to achieve effective fusion of GNSS and leveling elevation changes, which weakening the systematic deviation between GNSS geodetic height variation and normal height variation, and improving the reliability of regional elevation change analysis in Shandong Province.
-
表 1 速度融合结果精度统计
Table 1. Accuracy Statistics of Velocity Fusion Results
方案 拟合点残差/(mm·a-1) 检核点残差/(mm·a-1) 最大值 最小值 均方根误差 最大值 最小值 均方根误差 1 0.71 -0.69 0.36 1.62 -0.84 0.72 2 0.77 -0.58 0.33 1.34 -0.49 0.69 3 0.80 -0.63 0.35 1.33 -0.49 0.68 -
[1] Holdahl S R, Morrison N L. Regional Investigations of Vertical Crustal Movements in the U.S. Using Precise Relevelings and Mareograph Data[J]. Tectonophysics, 1974, 23(4):373-390 doi: 10.1016/0040-1951(74)90073-0 [2] Kamil K, Jacek R.Adjustment of Vertical Crustal Movement Network on the Basis of Last Three Leveling Campaigns in Poland[J].Reports on Geodesy, 2012, 92(1):123-134 [3] Kooi H, Johnston H, Lambeck K, et al.Geological Causes of Recent(~100yr) Vertical Land Movement in the Netherlands[J].Tectonophysics, 1998, 299(4):297-316 doi: 10.1016/S0040-1951(98)00209-1 [4] Ching K, Hsieh M L, Johnson K M, et al. Modern Vertical Deformation Rates and Mountain Building in Taiwan from Precise Leveling and Continuous GPS Observations, 2000-2008[J]. Journal of Geophysical Research, 2011, 116:B08406 http://www.researchgate.net/publication/251434124_Modern_vertical_deformation_rates_and_mountain_building_in_Taiwan_from_precise_leveling_and_continuous_GPS_observations_20002008 [5] 顾国华.GPS观测得到的中国大陆地壳垂直运动[J].地震, 2005, 25(3):1-8 http://d.old.wanfangdata.com.cn/Periodical/diz200503001 Gu Guohua. Vertical Crustal Movement Obtained from GPS Observation in China's Mainland[J].Earthquake, 2005, 25(3):1-8 http://d.old.wanfangdata.com.cn/Periodical/diz200503001 [6] 杨元喜, 曾安敏.大地测量数据融合模式及其分析[J].武汉大学学报·信息科学版, 2008, 33(5):500-503 http://ch.whu.edu.cn/CN/abstract/abstract1674.shtml Yang Yuanxi, Zeng Anmin. Fusion Modes of Various Geodetic Observations and Their Analysis[J]. Geomatics and Information Science of Wuhan University, 2008, 33(5):500-503 http://ch.whu.edu.cn/CN/abstract/abstract1674.shtml [7] 杨元喜, 曾安敏.顾及几何观测信息和地球物理模型的形变参数自适应滤波解[J].中国科学:地球科学, 2009, 39(4):437-442 http://www.cnki.com.cn/Article/CJFDTotal-JDXK200904006.htm Yang Yuanxi, Zeng Anmin. Adaptive Filtering for Deformation Parameter Estimation in Consideration of Geometrical Measurements and Geophysical Models[J]. Scientia Sinica Terrae, 2009, 39(4):437-442 http://www.cnki.com.cn/Article/CJFDTotal-JDXK200904006.htm [8] 黄立人, 韩月萍, 高艳龙, 等.GNSS连续站坐标的高程分量时间序列在地壳垂直运动研究中应用的若干问题[J].大地测量与地球动力学, 2012, 32(4):10-14 http://d.old.wanfangdata.com.cn/Periodical/dkxbydz201204003 Huang Liren, Han Yueping, Gao Yanlong, et al. Several Issues in Application of Elevation Component Time Series of GNSS CORS in Vertical Crustal Movement Studying[J]. Journal of Geodesy and Geodynamics, 2012, 32(4):10-14 http://d.old.wanfangdata.com.cn/Periodical/dkxbydz201204003 [9] 高艳龙, 郑智江, 韩月萍, 等.GNSS连续站坐标在天津地面沉降监测中的应用[J].大地测量与地球动力学, 2012, 32(5):22-26 http://www.cnki.com.cn/Article/CJFDTotal-DKXB201205006.htm Gao Yanlong, Zheng Zhijiang, Han Yueping, et al. Application of CORS of GNSS in Land Subsidence Monitoring in Tianjin[J]. Journal of Geodesy and Geodynamics, 2012, 32(5):22-26 http://www.cnki.com.cn/Article/CJFDTotal-DKXB201205006.htm [10] 郝明.基于精密水准数据的青藏高原东缘先进地壳垂直运动与典型地震同震及震后垂直形变研究[D].北京: 中国地震局地质研究所, 2012 http://www.cnki.com.cn/Article/CJFDTotal-GJZT201307008.htm Hao Ming. Present Crustal Vertical Movement of Eastern Tibetan Plateau and Coseismic and Postseimic Vertical Deformation of Two Typical Earthquakes[D]. Beijing: Institute of Geology, China Earthquake Administration, 2012 http://www.cnki.com.cn/Article/CJFDTotal-GJZT201307008.htm [11] Baryla R, Paziewski J. The Concept of the GNSS Control Network Densification with Precise Leveling for Ground Deformation Monitoring[C]. The 9th International Conference Environmental Engineering, Vilnius, Lithuania, 2014 [12] 王庆良, 崔笃信, 王文萍, 等.川西地区现今垂直地壳运动研究[J].中国科学:地球科学, 2008, 38(5):598-610 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=zgkx-cd200805007 Wang Qingliang, Cui Duxin, Wang Wenping, et al. Present Vertical Crustal Displacements of Western Sichuan Region[J]. Scientia Sinica Terrae, 2008, 38(5):598-610 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=zgkx-cd200805007 [13] 许才军, 王乐洋.大地测量和地震数据联合反演地震震源破裂过程研究进展[J].武汉大学学报·信息科学版, 2010, 35(4):457-462 http://ch.whu.edu.cn/CN/abstract/abstract917.shtml Xu Caijun, Wang Leyang. Progress of Joint Inversion of Geodetic and Seismological Data for Seismic Source Rupture Process[J]. Geomatics and Information Science of Wuhan University, 2010, 35(4):457-462 http://ch.whu.edu.cn/CN/abstract/abstract917.shtml [14] 许才军, 尹智.利用大地测量资料反演构造应力应变场研究进展[J].武汉大学学报·信息科学版, 2014, 39(10):1135-1146 http://ch.whu.edu.cn/CN/abstract/abstract3084.shtml Xu Caijun, Yin Zhi. Progress in Inversion for Tectonic Stress-Strain Fields Using Geodetic Data[J]. Geomatics and Information Science of Wuhan University, 2014, 39(10):1135-1146 http://ch.whu.edu.cn/CN/abstract/abstract3084.shtml [15] 张菊清, 刘平芝.抗差趋势面与正交多面函数结合拟合DEM数据[J].测绘学报, 2008, 37(4):526-530 doi: 10.3321/j.issn:1001-1595.2008.04.021 Zhang Juqing, Liu Pingzhi. Combining Fitting Based on Robust Trend Surface and Orthogonal Multiquadrics with Application in DEM Fitting[J]. Acta Geodaetica et Cartographica Sinica, 2008, 37(4):526-530 doi: 10.3321/j.issn:1001-1595.2008.04.021 [16] 杨元喜, 张菊清, 张亮.基于方差分量估计的拟合推估及其在GIS误差纠正中的应用[J].测绘学报, 2008, 37(2):152-157 doi: 10.3321/j.issn:1001-1595.2008.02.004 Yang Yuanxi, Zhang Juqing, Zhang Liang. Variance Component Estimation Based Collocation and Its Application in GIS Error Fitting[J]. Acta Geodaetica et Cartographica Sinica, 2008, 37(2):152-157 doi: 10.3321/j.issn:1001-1595.2008.02.004 [17] 郭春喜, 聂建亮, 王斌, 等.区域似大地水准面拟合方法及适用性分析[J].大地测量与地球动力学, 2013, 33(1):103-107 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=dkxbydz201301023 Guo Chunxi, Nie Jianliang, Wang Bin, et al. Methods for Fitting Local Quasigeoid and Their Applicability Analysis[J].Journal of Geodesy and Geodynamics, 2013, 33(1):103-107 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=dkxbydz201301023 [18] Yang Y, Zeng A, Zhang J.Adaptive Collocation with Application in Height System Transformation[J].Journal of Geodesy, 2009, 83(5):403-410 doi: 10.1007/s00190-008-0226-9 [19] 郭春喜, 聂建亮, 王斌, 等.利用自适应最小二乘配置的GPS水准与重力似大地水准面的拟合[J].武汉大学学报·信息科学版, 2013, 38(1):44-47 http://ch.whu.edu.cn/CN/abstract/abstract6069.shtml Guo Chunxi, Nie Jianliang, Wang Bin, et al. The Fitting of GPS-Level and Gravity Quasigeoid Based on Adaptive Collocation[J]. Geomatics and Information Science of Wuhan University, 2013, 38(1):44-47 http://ch.whu.edu.cn/CN/abstract/abstract6069.shtml [20] 樊子德, 李佳霖, 邓敏.顾及多因素影响的自适应反距离加权插值方法[J].武汉大学学报·信息科学版, 2016, 41(6):842-847 http://ch.whu.edu.cn/CN/abstract/abstract5473.shtml Fan Zide, Li Jialin, Deng Min.An Adaptive Inverse-distance Weighting Spatial Interpolation Method with the Consideration of Multiple Factors[J]. Geomatics and Information Science of Wuhan University, 2016, 41(6):842-847 http://ch.whu.edu.cn/CN/abstract/abstract5473.shtml [21] 马晨旭, 许才军, 张琼, 等.引入变差函数的反距离加权法及应用研究[J].大地测量与地球动力学, 2010, 30(1):83-87 http://d.old.wanfangdata.com.cn/Periodical/dkxbydz201001018 Ma Chenxu, Xu Caijun, Zhang Qiong, et al. On New Form of Inverse Distance Weight Interpolation Method with Variogram Function and Its Application[J]. Journal of Geodesy and Geodynamics, 2010, 30(1):83-87 http://d.old.wanfangdata.com.cn/Periodical/dkxbydz201001018 -