多面函数和移动法综合模型支持下的区域垂直形变场拟合

2018-03-07郑洪艳苏广利王家庆

测绘通报 2018年2期

田晓，郑洪艳，苏广利，王家庆

(中国地震局第一监测中心，天津 300180)

地壳垂直运动是地壳运动的一个重要分量，由于空间大地测量技术测定的高程分量精度低于水平分量，因此利用重复高精度水准测量仍然是目前研究地壳垂直运动的主要方法。但由于水准重合点是离散的，且水准网内部没有数据，因此动态平差得到的垂直运动速率在时间域和空间域都是不连续的，实际上地壳垂直运动是一种连续的形态变化。为了实现在时空域内客观准确地描述形变分布特征，往往需要基于离散的测点数据构建合适的数学模型拟合(最佳逼近)整个地区的地壳垂直形变速率面[1]。

在我国，从20世纪80年代末期开始，黄立人、赵承坤、杨国华和陶本藻等研究了多面函数拟合法在地壳垂直运动中的应用，几位学者针对多面函数拟合中核函数的选择、平滑因子的选择及核函数中心点的选择作了细致的研究，并提出了一系列独特的见解[2-6]。刘经南等[7]研究了多面函数拟合法在建立我国地壳平面运动速度场模型中的应用。王文利等[8]利用多面函数拟合了我国陆地垂直运动速率图。刘同文等[9]应用多面函数模型对西秦岭地区的垂直形变速率作了研究和探讨。秦姗兰等[10]利用水准资料研究了西秦岭地区的垂直形变，应用多面函数对区域整体形变特征进行拟合。

综上所述，研究地壳垂直运动中，多面函数法是一种非常传统并得到广泛应用的方法，而其他方法应用并不多。本文应用多面函数和移动法综合拟合模型，并以首都圈和晋冀蒙两个地区的垂直运动为例，与其他几种模型作对比研究。

1 模型方法数学原理

1.1 多面函数模型

多面函数模型是由美国Hardy教授于1977年提出的，并建议用于地壳垂直形变分析。该方法的基本思想为：任何数学表面和任何不规则的圆滑表面，总可以用一系列有规则的数学表面的和以任意精度逼近[9]。根据这一思想，假设地球表面上有一点(x,y)，其垂直速率值ξ(x,y)可以用下式来表示

(1)

式中，aj为待定系数；Qx,y,xj,yj是多面函数的核函数。核函数有多种类型，常用的有距离型核函数和锥面函数，其中距离型核函数表达式为

(2)

式中，δ为平滑因子，通常可取一小正数或0；μ一般取1/2、-1/2、3/2。当μ=1/2时，核函数称为正双曲面函数；当μ=-1/2时，核函数称为倒双曲面函数；当μ=3/2时，核函数称为三次曲面函数。

锥面函数表达式为

(3)

式中，c为待定参数。

设有m个已知水准点xi,yi，选取其中n(n≤m)个点作为核函数的中心点xj,yj，令Qij=Q(xi,yi,xj,yj)，则各数据点应满足

(4)

由此可列出误差方程

V=QX-ξ

(5)

根据最小二乘原理可求得待定系数X，即

(6)

待定系数求出后，根据式(4)可计算测区各未知点的垂直形变速率值。

需要注意的是，应用多面函数拟合垂直形变速率时，要根据测区大小和地形情况确定合适的核函数、平滑因子和中心点个数。

1.2 移动曲面模型

移动曲面模型的基本原理是以每一个待定点为中心，选取待定点周围适量的已知点，拟合出一个多项式曲面，进而内插出待定点上的函数值，是一种特殊的局部函数逼近法。

为了选取邻近的数据点，以待定点Pxp,yp为圆心，以R为半径作圆，凡落在圆内的数据点即被选用。应用时要根据实际情况选择合适的半径R。

定权时要根据已知点与待定点间的距离进行定权，常采用的权有以下几种形式

在进行拟合前，可先将参与拟合的各点坐标进行中心化处理，即将各点坐标减去待定点P的坐标，这样转换后就变成以P为原点的坐标

(7)

(8)

最后求待定点P的垂直形变速率时，直接代入坐标0,0即可，从而在很大程度上方便了计算，而且有利于改善方程的数值稳定性和提高解算精度[11]。

1.3 多面函数和移动法加权综合模型

该方法是基于综合预报理论提出的，它并不是一种全新的方法，而是首先采用多面函数和移动法进行单独拟合，然后对其拟合值进行加权综合而得到最终拟合值[12]。用数学模型可表示为

(9)

文献[12]中并没有给出该模型解的具体解法和表达式，本文给出W解的具体解法如下：

先构造拉格朗日极值函数

令

(10)

K=WTFTY-FW

(11)

另外，由式(10)可得

(12)

由式(11)和式(12)可采用迭代法解得W

(2) 将W0代入式(11)计算K。

(3) 将K代入式(12)计算W。

(4) 判断W-W0<10-6,若是，退出循环；若否，令W0=W继续重复步骤(2)—(4)。

经过以上迭代过程即可解得W为W*，则相应综合模型拟合结果为

(13)

1.4 多面函数和移动法综合拟合模型

多面函数法对整体性变化的拟合效果较好，但需要确定合适的模型参数，如果选取不当，容易出现拟合结果波动较大的现象，对局部变化拟合效果较差；而移动拟合法是以待定点为中心，以一定距离为半径，用一个多项式曲面拟合待定点函数值的局部拟合法，具有较好的灵活性[13]。为此，可以结合两者优点，采用先进行多面函数拟合，再对残差进行移动法拟合的综合模型拟合垂直形变速率。这很类似于组合法确定似大地水准面中的“移去—恢复”法，具体步骤如下：

(2) 求出m个拟合点的速率拟合残差v

v=ξN-ξ

(14)

(15)

2 算例与分析

2.1 拟合方案介绍

本文选用首都圈地区2001年和2005年，以及晋冀蒙地区2002年和2006年各两期水准数据进行研究，并对两个地区的数据采用自由网平差法进行动态平差，计算垂直运动速率。两段动态平差的单位权中误差分别为1.03和1.15 mm/a，说明这两个区域的数据整体实测精度比较好。剔除平差后速率中个别量值过大的突变点，首都圈地区共选用719个具有平差后速率的数据点(如图1所示)，晋冀蒙地区共有627个具有平差后速率的数据点(如图2所示)。由于对一个区域垂直形变场的拟合研究，主要为了得到未测点的垂直形变速率(即水准环内部区域的垂直形变速率)，因此为了更充分地得出结论，采用以下4种方案分别进行拟合计算：

方案1：均匀选取首都圈地区143个点作为外部检核点(约占总点数的1/5)，其余576个点作为拟合点参与函数模型的参数拟合(针对路线上数据)。

方案2：选取首都圈地区5条水准环内部路线全部点作为外部检核点(共80个点)，其余639个点作为拟合点参与函数模型的参数拟合(针对环内部数据)。

方案3：均匀选取晋冀蒙地区125个点作为外部检核点(约占总点数的1/5)，其余502个点作为拟合点参与函数模型的参数拟合(针对路线上数据)。

方案4：选取晋冀蒙地区5条水准环内部路线全部点作为外部检核点(共60个点)，其余567个点作为拟合点参与函数模型的参数拟合(针对环内部数据)。

本文采用外符合精度及检核点速率拟合值与实际观测值之差(残差)的绝对值大于2 mm/a的点的个数来衡量拟合的精度。外符合精度σ的表达式为

(16)

式中，Δi表示检核点速率拟合值与实际观测值之差；n为检核点个数，并称 Δi>2 mm/a的点为病态点。

图1 首都圈地区垂直形变速率(相对于北京)

图2 晋冀蒙地区垂直形变速率(相对于大同)

2.2 各种模型拟合结果比较

对首都圈地区的拟合计算有方案1和方案2，各种模型拟合结果见表1。

表1 首都圈地区各种模型拟合结果比较

为了直观地看出各种模型的拟合效果，给出不同模型在方案1和方案2中拟合检核点速率残差情况，如图3和图4所示。

图3 方案1不同模型检核点拟合残差

图4 方案2不同模型检核点拟合残差

对晋冀蒙地区的拟合计算有方案3和方案4，各种模型拟合结果见表2。

表2 晋冀蒙地区各种模型拟合结果比较

为了直观地看出各种模型的拟合效果，给出不同模型在方案3和方案4中拟合检核点速率残差情况，如图5和图6所示。

图5 方案3不同模型检核点拟合残差

图6 方案4不同模型检核点拟合残差

2.3 算例结果分析

(1) 从表1和表2中可以发现，4种模型中两种综合模型的拟合精度优于两种单一模型，但是特性不同。多面函数和移动法加权综合模型虽然精度上优于多面函数模型，但是病态点个数并不具有优势，这是因为加权综合会使得单一模型的各种误差传递到最终结果中。而多面函数和移动法综合拟合模型不管在精度上还是病态点个数上均优于多面函数模型。

(2) 多面函数和移动法综合拟合模型在两个地区4种不同情况下都是4种模型中拟合精度最高的，病态点数也是最少的，并且从图4—图7可以看出多面函数和移动法综合拟合模型的拟合残差是最稳定的。