海洋重力不规则测线网平差模型对比分析

马 龙,郑彦鹏,华清峰,郭延良,赵 强,邢 健

(1.自然资源部 第一海洋研究所,山东 青岛266061;2.中国科学院 海洋地质与环境重点实验室,山东 青岛266071;3.自然资源部 海底矿产资源重点实验室,广东 广州510075;4.青岛海洋科学与技术国家试点实验室 海洋地质过程与环境功能实验室,山东 青岛266061;5.国家海洋局 东海海洋环境调查勘察中心,上海200137)

海洋调查过程中,受仪器设备、海况、定位精度等多种因素影响,重力观测值内包含一定大小的偶然误差和系统误差。海洋重力测网平差作为内业数据处理、质量控制的最后一步,其关键性不言而喻。而传统的测线网平差方法主要针对规则测线网下的系统差检验及测线误差调整[1-2]。部分学者针对定位误差等对重力观测值的影响[3-5],提出航迹拟合、提高厄特渥斯改正的测线网平差方法。随着定位精度的不断提高,研究逐渐注重于主测线和联络线本身的交点平差改正,基于最小二乘法建立方程组系数矩阵求解平差值的方法[6-14]被广泛应用于国内外的海洋重力测线网数据的平差处理中。

目前,海洋重力资料测网平差主要基于半系统差调整和最小二乘两种平差模型。在规则测线网中,2种平差模型的平差过程近似等效,结果一致。然而,传统的测线网仅依据主测线和联络线的交点个数进行划分,以m条主测线、n条联络线组成的测线网为例,具有m×n个交点时为规则测线网,因走向、长度差异造成交点数少于m×n个交点的测线网则为不规则测线网。该分类方法并未考虑同一区域内多个航次的往复测量、国际共享资料的逐年累加,乃至航次内的补充测线等问题。本文将结合测线分布和航次实例,重新划分复杂测线网内的测线类型,针对目前已经被广泛运用的半系统差调整和最小二乘两种模型进行讨论,提出一种适用于不规则复杂测线网下重力测线数据的平差方法,最终确定最佳的平差处理方式。

1 数据资料和平差方法

1.1 数据资料

在实际调查航次中,随着航次内补充测线及更多历史资料的加入,测区测线的分布更为复杂,此时分别采用最小二乘和半系统差调整两种平差模型计算所获得的结果则不尽相同。本文以国际共享MW9006航次的船测重力资料为例[15],开展数据平差模型选择和结果的研讨工作,以期建立一种适用于不规则测线网下重力数据的平差方法。MW9006航次为美国Moana Wave综合科考船执行的一次基金科考航次,自关岛阿加尼亚湾启航,菲律宾马尼拉港返航(该航次航迹线见图1a)。按照测线航迹分布、航速等因素,剔除弯折段、停航段及无交点等无效数据,将测线分为50条主测线(M),28条联络线(L)和3条补充测线(S)三类(图1b)。

图1 MW9006航次测线分割Fig.1 Segmentation of survey line No.MW9006

1.2 平差方法

海洋重力数据是在观测平台不断运动的情况下测量获得的,除包含系统误差外,还含有设备、海况等[8]各种因素引起的偶然误差。海洋重力测线网交点处的主测线(g ij)、联络线(g ji)交点观测值可表示为

式中,i和j为测线编号;g t为交点真实值;εij、βji分别为第i条主测线和第j条联络线的系统误差;Δij为随机误差,交点差(d ij)为

数据处理过程中将各个测点视为等精度观测。采用半系统差调整方式计算时,分别对主测线和联络线进行一次半系统差调整后,交点差(d ij)则变为[1]

迭代计算至交点差代数和为0。最小二乘法平差时,计算研究区测线交点差平方和(D):

对交点差平方和(D)求导并构建线性方程组A X=B,利用最小二乘平差原理确定方程组系数。数列X即为每条测线上的平差值,其物理意义为每条主测线和联络线调差值的总和为0。

假设测区仅存在m条主测线、n条联络线,当交点数等于或小于m×n个交点时为规则测线网。然而在航次开展过程中往往存在与主测线和联络线均存在交点的补充测线,这种情况在添加历史航次测线时则会更加突出,定义该测线网为不规则测线网。从计算原理及公式推导分析,在规则重力测线网下,2种模型的平差目的和最终结果均为测区内所有测线的调差值总和为0[2],即在此情况下2种平差方式近似等效。

2 规则测线网数据平差分析

在不考虑补充测线的情况下,国际共享MW9006航次有效测线组成的测线网为交点数小于50×28点的规则测线网,分别采用半系统差调整、最小二乘两种平差模型对其进行计算和对比分析。平差前计算测线交点差见图2。

图2 平差前交点差分布(×10-5 m·s-2)Fig.2 Distribution of crossover before adjustment(×10-5 m·s-2)

由图1可见,测区共有50条主测线(M),28条联络线(L),构建了78阶×78阶系数矩阵(A):

采用Morre-Penrose广义逆求解由最小二乘法构建的线性方程组,获得每条测线的平差值数列X。半系统差调整方式平差精度与平差次数正相关,平差到一定次数后,交点均方根趋近于最小二乘平差结果,精度与平差次数相关性见图3。

经半系统差调整第一次平差后,交点差均方根由(±4.620×10-5)m/s2减小至(±1.385×10-5)m/s2(图3),并随着平差次数叠加了150次后趋于(±1.381×10-5)m/s2,且第二次与第一次平差后精度差值小于10-8m/s2,已经完全达到了外业调查的精度要求。但在实际计算中,为进一步保障数据精度,建议在运算能力允许的情况下尽可能多地增加运算次数。相较而言,最小二乘平差模型平差过程中,虽然构建系数矩阵的过程比较复杂,但后期运算量小,运算一步即可获得预期效果;半系统差调整方法虽然运算量大,但更易于编程实现。

图3 半系统差调整次数与交点差精度对比Fig.3 Comparison of half-system adjustment times with crossover difference accuracy

图4 平差后交点差分布(×10-5 m·s-2)Fig.4 Distribution of crossovers after adjustment(×10-5 m·s-2)

图5 平差后交点差统计Fig.5 Histogram of crossover errors after adjustment

规则测线网经2种模型平差处理后的最终结果见图4和图5。结果表明,海洋重力规则测线网下,经上述2种平差模型处理后结果一致。分析对比平差前后(图2和图4)的结果,交点差平均值由平差前-0.65×10-5m/s2降为平差后的0 m/s2,均方根由平差前的(±4.62×10-5)m/s2降为平差后的(±1.38×10-5)m/s2,最终主测线、联络线之间的误差值均被控制在合理的范围内,交点差分布也基本呈正态分布。

3 不规则测线网数据平差分析

在规则测网下,随着半系统差调整次数的逐渐增多,其结果与最小二乘平差趋近,两者方法近似于等效。然而,将测区内的补充测线(S)加入之后,针对如何合理有效地将补充测线数据整合到已经平差结束的测线网,两种平差模型必然需要选择完全不同的平差方式。以本测区测线为例,S79,S80和S81三条补充测线分别与主测线和联络线存在交点。

3.1 不规则测线网下最小二乘法平差

如果采用最小二乘法平差,针对该补充测线数据可采用4种方式进行分析。方式1:构建系数矩阵时,先暂时舍弃补充测线与主测线的交点,即将补充测线设定为主测线进行平差计算。方式2:构建系数矩阵时,先暂时舍弃补充测线与联络线的交点,即将补充测线设定为联络线进行平差计算。方式3:以与主测线、联络线的交点个数更少和交点误差更小为原则,逐条将补充测线设定为对应的测线类型并构建系数矩阵。方式4:将全部补充测线视为一个系统,经最小二乘平差后的规则测线网视为另一个高精度系统,将两个系统做差,将补充测线系统调差至已平差的高精度规则测线网系统中。交点差计算时,与主测线相交时,将补充测线视为联络线;与联络线相交时,将补充测线视为主测线。最终在计算数据精度时,以均方根作为评价参数,其计算公式为

式中,N为计算误差的点数。平差前研究区补充测线与主测线和联络线的交点差分布见图6。

补充测线与主测线交点多,但与联络线交点少且误差相对更小,故在采用方式3计算时,2条补充测线均被设定为联络线计算。因此,本测区方式2与方式3的计算结果相同,最终计算结果如图7所示。

由上述计算结果可知,采用方式1～4平差后,交点差均方根较平差前的(±4.62×10-5)m/s2均明显减小。由于采用方式1和方式4平差后,将补充测线设定为主测线或独立系统,对测线进行平差时未考虑部分交点(方式1)或引入补充测线系统(方式4),导致补充测线的部分交点与经最小二乘平差后的测线交点差进一步增大,从而引入了更多的误差。然而,依据测线分布的特点,若采用方式2和方式3则结果相同,且数据整体基本满足正态分布,比前两种方式要好。

3.2 不规则测线网下半系统差调整平差

半系统差调整方法基于全部测点等精度观测,可分2种方式平差:1)经方式4处理后,再对全部平差结果测线遍历进行半系统差调整平差(方式5);2)遍历全部原始测线完成数据平差(方式6)。经计算,方式5和方式6平差后交点差对比见图8,方式6半系统调整平差后交点差分布见图9。

图9 方式6半系统调整平差后交点差分布Fig.9 Distribution of crossover errors after semi-system adjustment by using method six

采用方式5和方式6经平差后,测线交点差相等,但不是每条测线的平差值都相同(图8)。规则测线网下,全部测线的调差值总和为0,然而引入补充测线经半系统差调整平差时,方式6测线平差值绝对值普遍比方式5大。经计算,方式5的调差值总和为18.86(绝对值总和为275.30),方式6的调差值总和为-71.42(绝对值总和为290.80)。分析认为,由于方式5判断原主测线、联络线组成的测线网系统精度更高,在完成补偿测线的独立系统平差时,先完成了整体调差(消除相对系统误差),所以盲目引入补充测线,参与全测线半系统调差,会造成一定的系统误差。

不规则测线网下最小二乘法平差中的方式2和方式3与半系统差调整平差后(方式5),交点差均方根绝对值仅相差0.03×10-5m/s2,且两者整体均成正态分布。然而,方式2和方式3交点差分布范围为-8×10-5～8×10-5m/s2比方式5(范围为-10×10-5～10×10-5m/s2)更小,所以前者的误差分布整体均一性更好,故针对本航次平差情况,建议采用方式2和3进行测线平差。

4 讨 论

部分学者通过数据分析和确定参数的方式选择平差模型或方法[6,9],处理过程具有一定的主观性,因而在一定程度上平差结果不具有很强的复制性。本文基于国际共享MW9006航次船载实测资料,采用目前已经被广泛运用的最小二乘和半系统差调整两种方式展开研究讨论,提出一种适用于不规则测线网下重力数据的平差方法(图10)。

通过程序设计与实现[16-17],在规则测线网平差时,采用最小二乘和半系统调差模型平差处理一致,均能实现明显降低交点差的目标。然而,一旦涉及与主测线、联络线均有交点的补充测线时,最小二乘平差方式可通过舍弃交点或增添补充测线系统的方式完成测线平差,最终导致损失部分交点数据而降低平差精度;半系统差调整可通过遍历测线的方式计算全部交点并参与平差,但存在交点差重复参与计算的问题。随着补充测线的逐步增多,需要依据补充测线与主测线和联络线的交点分布及数值大小情况,确定补充测线在参与计算时属于何种测线,并依次采用上述2种模型进行计算。

经过分析,在实际分类时共存在6种情况,具体分类见3.1、3.2节。由于半系统调差采用不规则测线网内所有测线遍历计算并平差的方式计算,其结果可以作为首要备选结果(方式6)。与此同时,按照补充测线与主、联络线交点个数及相对误差情况,依次将补充测线设定为主测线(方式1)、联络线(方式2)、部分主测线和部分联络线(方式3)、独立的测线系统(方式4)、方式4后再半系统差调整平差(方式5),加之首选结果共6种情况。通常情况下而言,方式3和方式5计算的结果相对更为合理。最终将最合理的结果与半系统差调整平差结果(方式6)进行分析对比,按照测线精度、交点差均方根值大小及平差结果整体的均一性程度选择最终的平差方式及结果。以国际共享MW9006航次为例,虽然平差后均方根绝对值仅相差0.03×10-5m/s2,但采用将补充测线识别为联络线的方式获得的结果交点差分布范围更小,且方式3与方式2结果一致,故最终采用方式2和方式3对整个数据进行平差处理。

图10 海洋重力不规则测线网数据平差流程Fig.10 Adjustment flow chart of marine gravity and magnetic irregular survey network

5 结 论

本文基于国际共享MW9006航次船载实测资料,针对目前已经被广泛应用的最小二乘和半系统差调整两种模型进行探讨,提出了一种普遍适用于不规则复杂测线网下重力测线数据的平差方法。主要结论如下:

1)规则测线网下的平差结果表明,随着半系统差调整平差次数的增加,其平差结果趋近于最小二乘平差,交点差平均值由平差前-0.65×10-5m/s2降为平差后的0 m/s2,均方根由平差前(±4.62×10-5)m/s2降为平差后(±1.38×10-5)m/s2,满足平差要求,结果基本呈正态分布。然而,一旦涉及与主测线、联络线均有交点的补充测线时,两种平常方法又均有局限。最小二乘平差只能损失部分交点数据,避免其参与平差而降低精度;半系统差调整存在交点差重复参与计算的问题。

2)基于补充测线与主测线、联络线交点的分布及大小等因素,将补充测线按4种类型参与最小二乘平差,或分两种类型将主测线、联络线、补充测线定为一个整体系统,利用半系统差调整平差,采用共6种类型的处理方式获得平差结果。其中方式5和方式6经平差后,测线交点差相等,但并非每条测线的平差值均相同,方式6测线平差值绝对值普遍比方式5大,经计算方式5的调差值总和为18.86(绝对值总和为275.30),方式6的调差值总和为-71.42(绝对值总和为290.80)。经本文分析,盲目引入补充测线参与全测线半系统调差,会造成一定的系统误差。方式5和方式6判断依据是主测线、联络线的规则测线网系统与补充测线独立系统的测线精度比较,如果前者精度高,采用方式5,反之则选择方式6。

3)在对6种方式进行评价时,选取了测线精度、交点差均方根及平差结果整体均一性三个要素进行分析评估,并最终确定本文实例中方式2和方式3为最佳的平差方式,最终航次内114个交点均方根降为(±1.84×10-5)m/s2。

本文所采用的计算实例,补充测线数据少且相对规律,随着补充测线及历史测线的逐渐增多,分析的复杂性也随之剧增,但可借鉴图10所示的平差流程图开展计算。同理,由于数据类型、误差来源不同,本文仅针对海洋重力资料展开讨论,但考虑方法、模型相同,该方法同样适用于磁力及水深测网数据。