基于海底地形-重力联合提高海洋重力空间分辨率的方法

摘要

本发明公开了一种基于海底地形‑重力联合提高海洋重力空间分辨率的方法，包括：对大地水准面的正常重力场进行高度改正和中间层改正后，求得海底控制点处的理论重力值；确定海底控制点处的自由空间重力异常、由水深变化引起的重力异常和大地水准面上的自由空间重力异常三者之间的关系；得到海底地形‑重力联合法的基本公式模型；求解得到离散海底控制点处的绝对重力值；求解得到格网化海底绝对重力值；通过已有格网化高精度海底地形数据及重力校正，得到海底地形‑重力联合法模型，并求解得到大地水准面上格网化的自由空间重力异常。本发明旨在保证高精度的前提下，提高海洋重力异常模型的空间分辨率。

说明书

技术领域

本发明属于海洋重力学、水下导航学等交叉技术领域，尤其涉及一种基于海底地形-重力联合提高海洋重力空间分辨率的方法。

背景技术

水下重力匹配导航技术旨在利用外部重力信息修正惯性导航系统随时间累积的误差，是真正的无源导航技术，有利于实现水下潜器自主、连续、精确和长航时的定位。基本原理：载体在航行过程中经过重力特征较丰富区域时利用海洋重力仪/重力梯度仪采集周围重力场信息，并与预先存储在导航系统中的重力基准图进行对比，依据相关准则判断二者之间的拟合度，从而确定最佳匹配位置。因此，获取高精度和高空间分辨率海洋重力基准图是提高水下重力匹配导航效果的关键因素之一。

随着重力场测量技术的发展，由于卫星测高技术反演海洋重力场高效便捷且能够有效地解决传统获取手段的数据稀疏、重复周期性差、舰船无法直接到达偏远海洋区域等问题，是目前获取高空间分辨率、高精度及全球海洋重力场的较优手段。(1)目前许多学者主要围绕卫星测高/重力及船测重力数据开展海洋重力场反演研究。1979年，Rapp首次利用Geo-3卫星得出若干个1°×1°和5°×5°的海洋重力场值，其精度分别为7.8mGal和2.7mGal。1989年，Hwang在Rapp的研究基础上，联合Geo-3和Seasat卫星所有数据，并加入海面地形、潮汐、卫星偏差等影响因素，计算结果与Rapp的重合部分计算结果均方根差为5.43mGal。1998年，Hwang联合多种测高卫星数据，利用简洁、严密、更加容易计算的Vening-Meinesz公式得到精度为9.9mGal的全球海洋重力场模型。1999年，许厚泽等利用垂线偏差法联合四年的T/P卫星和1年的ERS-1数据恢复出中国近海重力异常分辨率为30′×30′且精度为3.5mGal。2003年，李建成等利用Topex/Poseidon、ERS-1/2、Geosat/GM和ERM测高资料结合垂线偏差法反演出空间分辨率2.5′×2.5′且精度为9.35mGal的中国近海重力异常。2009年，李娜通过卫星测高数据和船测重力数据，在垂线偏差法计算的逆Vening-Meinesz公式中，证实了积分半径在30′～60′最佳。2009年，Sandwell等引入EGM2008全球海洋重力场模型和平均动态海面模型作为参考模型并结合双调和插值法，利用垂线偏差法反演出分辨率1′×1′和精度优于4mGal的全球海洋重力场模型V18.1。2014年，Sandwell联合现存的测高卫星数据和Cryosat-2卫星数据，把1′×1′的全球海洋重力场模型的精度从上一代的3-5mGal提高到2-4mGal。2014年，鲍李峰等提出双星伴飞的测高模式，成功获取分辨率为1′×1′区域高精度海洋重力场。(2)目前许多学者主要利用海洋重力场和船测水深数据反演海底地形研究。2011年，聂琳娟等利用卫星测高得到重力异常，基于非补偿板块模型，反演了南海海底地形，反演结果与LDEO船测深数据相比较，两者之差的RMS为253.3m，较ETOPO5海深模型精度和分辨率均有明显提高。2015年，Hu等提出了基于垂直重力梯度异常(VGG)，重力异常和船舶测深推导的全球海平面地形模型。船舶测深仪用于限制波长大于200km的海底地形，VGG比率可用于预测100-200km波长带的海底地形，重力异常用于恢复波长小于100km的海底地形。结果表明，在所讨论的区域中，该模型的精度优于ETOPO1、GEBCO和DTU10，并接近V15.1模型。2018年，范雕等应用自适应赋权技术，采用导纳函数方法构建菲律宾海域1'×1'海底地形模型。试验发现，当重力异常垂直梯度反演海深结果与重力异常反演海深结果的权比为2:3时，所构建的海深模型检核精度最高，而相较于ETOPO1海深模型和DTU10海深模型检核精度分别提高了约27.17％和39.02％。综上所述，目前未发现有公开文献显示基于高精度和高空间分辨率海底地形数据反演海洋重力场的研究。

发明内容

本发明的技术解决问题：克服现有技术的不足，提供一种基于海底地形-重力联合提高海洋重力空间分辨率的方法，旨在保证高精度的前提下，提高海洋重力异常模型的空间分辨率。

为了解决上述技术问题，本发明公开了一种基于海底地形-重力联合提高海洋重力空间分辨率的方法，包括：

步骤1，对大地水准面的正常重力场进行高度改正和中间层改正后，求得海底控制点处的理论重力值

步骤2，根据重力地质法，确定海底控制点处的自由空间重力异常Δg

步骤3，基于步骤1和步骤2，得到海底地形-重力联合法的基本公式模型；

步骤4，基于海底地形-重力联合法的基本公式模型，求解得到离散海底控制点处的绝对重力值

步骤5，根据求解得到的离散海底控制点处的绝对重力值

步骤6，通过已有格网化高精度海底地形数据及重力校正，得到海底地形-重力联合法模型，并根据海底地形-重力联合法模型，求解得到大地水准面上格网化的自由空间重力异常Δg

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，对大地水准面的正常重力场进行高度改正和中间层改正后，求得的海底控制点处的理论重力值

其中，g(x,y)表示大地水准面处的正常重力，g

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，通过如下公式，确定大地水准面处的正常重力g(x,y)：

g(x,y)＝g

其中，g

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，通过如下公式，确定大地水准面到海底控制点的高度改正g

其中，h表示大地水准面距海底控制点深度，R表示地球平均半径。

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，其特征在于，通过如下公式，确定大地水准面与海底控制点之间的中间层改正g

g

其中，ρ表示中间层物质密度，G表示万有引力常数。

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，Δg

Δg

其中，Δg

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，基于步骤1和步骤2，得到海底地形-重力联合法的基本公式模型，包括：

确定海底控制点处的绝对重力值

联立公式(1)、(2)、(5)和(6)，推导得到海底地形-重力联合法的基本公式模型，表示如下：

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，基于海底地形-重力联合法的基本公式模型，求解得到离散海底控制点处的绝对重力值

根据式(7)，确定离散海底控制点处的绝对重力值的解算表达式：

其中，Δg

根据公式(8)，解算得到离散海底控制点处的绝对重力值

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，海底地形-重力联合法模型的表达式如下：

其中，g

在上述基于海底地形-重力联合提高海洋重力空间分辨率的方法中，

Δg

其中，h

本发明具有以下优点：

鉴于当前公布的全球海底地形空间分辨率(0.25′×0.25′)高于全球海洋重力异常空间分辨率(1′×1′)，这是由于当前丰富的海底地形数据并未参与海洋重力场的反演研究。通过卫星测高/重力及船测重力数据等主要反演手段获取更高分辨率的海洋重力场技术较难在短期内得到突破。因此，本发明综合考虑了海洋重力异常与海底地形之间的相关性，利用当前丰富的高精度和高空间分辨率海底地形数据，结合高精度船测重力数据，并联合重力校正和重力地质法，提出了新型海底地形-重力联合法，旨在保证高精度的前提下，提高海洋重力异常模型的空间分辨率，这是对目前主流卫星测高/重力及船测重力数据反演海洋重力场的补充探索研究。

附图说明

图1是本发明实施例中一种基于海底地形-重力联合提高海洋重力空间分辨率的方法的步骤流程图；

图2是本发明实施例中一种研究区域数据的示意图；2(a)为1′×1′海底地形图，2(b)为船测控制点及检核点离散分布；

图3是本发明实施例中一种经验密度参数选取示意图；3(a)为相关系数曲线，3(b)为标准差曲线；

图4是本发明实施例中一种1′×1′海洋重力异常示意图；4(a)为模型V28.1，4(b)为海底地形-重力联合法反演模型；

图5是本发明实施例中一种1′×1′海洋重力异常模型V28.1与检核点误差散点图；5(a)为三维，5(b)为二维；

图6是本发明实施例中一种海底地形-重力联合法反演模型与船测检核点误差散点图；6(a)为三维，6(b)为二维；

图7是本发明实施例中一种重力异常误差绝对值统计直方图；7(a)为模型V28.1，7(b)为海底地形-重力联合法反演模型；

图8是本发明实施例中又一种研究区域数据的示意图；8(a)为0.25′×0.25′海底地形图，8(b)为1′×1′海洋重力异常模型V28.1控制点及检核点分布；

图9是本发明实施例中又一种经验密度参数选取示意图；9(a)为相关系数曲线，9(b)为标准差曲线；

图10是本发明实施例中一种海洋重力异常示意图；10(a)为1′×1′模型V28.1，10(b)为0.25′×0.25′海底地形-重力联合法反演模型；

图11是本发明实施例中一种0.25′×0.25′海底地形-重力联合法反演模型与V28.1检核点的重力异常误差散点图；11(a)为三维，11(b)为二维；

图12是本发明实施例中一种基于海底地形-重力联合法的0.25′×0.25′海洋重力异常加密模型V28.1示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明公开的实施方式作进一步详细描述。

鉴于当前公布的全球海底地形空间分辨率(0.25′×0.25′)高于全球海洋重力异常空间分辨率(1′×1′)，这是由于当前丰富的海底地形数据并未参与海洋重力场的反演研究。通过卫星测高/重力及船测重力数据等主要反演手段获取更高分辨率的海洋重力场技术较难在短期内得到突破。因此，本发明综合考虑了海洋重力异常与海底地形之间的相关性，利用当前丰富的高精度和高空间分辨率海底地形数据，结合高精度船测重力数据，并联合重力校正和重力地质法，提出了新型海底地形-重力联合法，旨在保证高精度的前提下，提高海洋重力异常模型的空间分辨率，这是对目前主流卫星测高/重力及船测重力数据反演海洋重力场的补充探索研究。

本发明的核心思想之一在于：公开了一种基于海底地形-重力联合提高海洋重力空间分辨率的方法，开展了联合高精度和高空间分辨率海底地形与高精度重力异常数据提高海洋重力异常基准图空间分辨率的探索性研究。第一，综合考虑了海洋重力异常与海底地形之间的相关性，以高精度和高空间分辨率海底地形与高精度船测重力异常为数据元，联合重力校正和重力地质法，提出新型海底地形-重力联合法；第二，以研究区域(经度131°～133°E，纬度36°～38°N)为例，联合1′×1′全球海底地形与高精度船测重力异常数据，基于新型海底地形-重力联合法反演的1′×1′海洋重力异常模型与船测检核点的标准差为4.91mGal，而1′×1′海洋重力异常模型V28.1与船测检核点的标准差为4.60mGal；且与船测检核点差值绝对值小于10mGal的点数，1′×1′新型海底地形-重力联合法反演模型占比约95.72％，1′×1′海洋重力异常模型V28.1占比约95.83％，表明海底地形-重力联合法反演模型精度接近于海洋重力异常模型V28.1。第三，以0.25′×0.25′全球海底地形和1′×1′全球海洋重力异常模型V28.1为基础数据，基于新型海底地形-重力联合法，反演了0.25′×0.25′全球海洋重力异常，其与V28.1检核点的标准差为0.36mGal，进而对1′×1′海洋重力异常模型V28.1实现了加密补充，旨在保证精度前提下构建全球更高空间分辨率的海洋重力基准图提供了技术参考。新型海底地形-重力联合法的优点为全球海洋重力异常图空间分辨率高、计算速度快。

如图1，在本实施例中，该基于海底地形-重力联合提高海洋重力空间分辨率的方法，包括：

步骤101，对大地水准面的正常重力场进行高度改正和中间层改正后，求得海底控制点处的理论重力值

在本实施例中，基于地球物理勘探领域的重力校正知识，可以对大地水准面的正常重力场进行高度改正和中间层改正后，求得海底控制点处的理论重力值

其中，g(x,y)表示大地水准面处的正常重力，g

优选的，可以通过如下公式，确定大地水准面处的正常重力g(x,y)：

g(x,y)＝g

其中，g

优选的，可以通过如下公式，确定大地水准面到海底控制点的高度改正g

其中，h表示大地水准面距海底控制点深度，R表示地球平均半径。

优选的，可以通过如下公式，确定大地水准面与海底控制点之间的中间层改正g

g

其中，ρ表示中间层物质密度，如，中间层物质为海水，海水密度取理论值1.03g/cm

步骤102，根据重力地质法，确定海底控制点处的自由空间重力异常Δg

在本实施例中，Δg

Δg

其中，Δg

步骤103，基于步骤101和步骤102，得到海底地形-重力联合法的基本公式模型。

在本实施例中，自由空气异常(free air anomaly)指地球自然表面上的绝对重力值，经过高度改正归算到大地水准面后，再减去正常重力值得到的重力差。即，控制点处的重力值为控制点处理论值和控制点所在高度上的自由空间重力异常之和，即可以确定海底控制点处的绝对重力值

进一步的，联立公式(1)、(2)、(5)和(6)，可以推导得到海底地形-重力联合法的基本公式模型，表示如下：

步骤104，基于海底地形-重力联合法的基本公式模型，求解得到离散海底控制点处的绝对重力值

在本实施例中，基于公式(7)，可以改写得到离散海底控制点处的绝对重力值的解算表达式：

其中，Δg

进一步的，根据公式(8)，可以解算得到离散海底控制点处的绝对重力值

步骤105，根据求解得到的离散海底控制点处的绝对重力值

步骤106，通过已有格网化高精度海底地形数据及重力校正，得到海底地形-重力联合法模型，并根据海底地形-重力联合法模型，求解得到大地水准面上格网化的自由空间重力异常Δg

在本实施例中，基于公式(7)，可以改写得到海底地形-重力联合法模型的表达式如下：

其中，g

进一步的，格网化的水深变化引起的重力异常Δg

Δg

其中，h

最后，可以根据公式(9)，可以求解得到大地水准面上格网化的自由空间重力异常Δg

在上述实施例的基础上，下面对本发明实施例所述的基于海底地形-重力联合提高海洋重力空间分辨率的方法进行验证。

1、新型海底地形-重力联合法验证

数据区域一：研究数据范围以日本海区域(经度131°～133°E，纬度36°～38°N)为例，GEBCO网站提供的1′×1′海底地形图(由0.25′×0.25′全球海底地形图稀疏处理提取，如图2(a)所示)；加州大学圣迭戈分校网站(https://topex.ucsd.edu/cgi-bin/get_data.cgi)提供的1′×1′海洋重力异常模型V28.1，用于与所反演的海洋重力场精度对比分析；源于美国国家地球物理数据中心(https://www.noaa.gov)提供的船测重力异常数据，船测点分布如图2(b)所示，共有船测重力异常数据点11626个，其中蓝色点作为控制点有8720个，用于参与海洋重力异常的反演计算过程；红色点作为检核点有2906个，用于对经验密度差的确定及重力异常反演结果的外部检核。

密度差Δρ为经验密度参数(无实际物理意义)，但却至关重要，在反演前需要选择较为合理的参数。本发明采用迭代法进行密度差确定，选取8720个重力异常控制点和2906个检核点，利用控制点计算不同密度差常数下的海洋重力场模型，采用最邻近法得到检核点位置的预测重力异常，通过比较预测值与检核点之间相关系数与差值标准差，分别绘制相关系数变化曲线图3(a)和差值标准差变化曲线图3(b)，将相关系数达到最大时的Δρ作为该区域的经验密度参数，可以看出相关系数最大时，标准差也达到最小值。

据图3可知，密度差参数在一定范围内，随着数值的增大相关系数随之增大，差值标准差逐渐减小；当密度差大于某一常数时，相关系数随密度差增大而降低，标准差逐渐增大。通过上述解算拟合的两条曲线，最终确定拟合曲线的相关系数达最大值0.96时，密度差参数取值约为5.3g/cm

基于新型海底地形-重力联合法，可求取8720个海底控制点处的实测重力值，并通过适当插值求得1′×1′海底实测重力值，并通过已有1′×1′海底地形数据，获得基于海底地形-重力联合法反演模型(图4(b)所示)。

本发明采用外部测点检核法，对图4中的模型V28.1和基于海底地形-重力联合法反演模型进行精度分析，是将模型利用最邻近法插值到2906个检核点上，求得预测值，通过计算预测值与检核点值之间的差异，对模型精度进行评价。

图5表示1′×1′海洋重力异常模型V28.1与检核点的误差散点图，它是将1′×1′海洋重力异常模型V28.1(图4(a))利用最邻近法插值到2906个检核点上，并与检核点重力异常数据之差的散点图，描述了图4(a)与船测检核点之间的误差。图6表示1′×1′海底地形-重力联合法反演模型与检核点的误差散点图，它是将基于海底地形-重力联合法所反演的1′×1′海洋重力异常模型(图4(b))利用最邻近法插值到2906个检核点上，并与检核点重力异常数据之差的散点图，描述了图4(b)与检核点之间的误差。图5(b)和图6(b)的横坐标代表检核点的序号，纵坐标代表重力异常差值大小。

据图5和图6可知，两种模型与检核点的重力异常误差均集中于0mGal附近，大致介于-10～10mGal之间，表明1′×1′海底地形-重力联合法反演模型精度接近于1′×1′海洋重力异常模型V28.1。此外，图5与图6中个别检核点重力异常误差值较大，主要有两点分析：①可能某些检核点存在较大数据粗差；②这是由于地球并非密度均衡的球体，而不论海洋重力异常模型，还是海底地形-重力联合法反演模型的构建过程均是基于大量控制点数据的数学反演，无法排除某些非控制点位置密度突变的不确定性，因此说明这些重力异常误差值较大的检核点位置应该存在密度的突变，即可能存在断崖高峰等地质构造或油气矿藏等自然资源。

图7表示1′×1′海洋重力异常模型V28.1、1′×1′海底地形-重力联合法反演模型分别与船测检核点的重力异常误差绝对值统计直方图。图中横坐标表示误差绝对值，纵轴表示所占的百分比。据图7可知，与船测检核点差值绝对值小于5mGal的点数，海底地形-重力联合法反演模型占比约82.35％，海洋重力异常模型V28.1占比约84.82％；小于10mGal的点数，反演模型占比约95.72％，模型V28.1占比约95.83％。1′×1′海底地形-重力联合法反演模型与检核点重力异常的标准差和均值分别为4.91mGal和0.28mGal；而1′×1′海洋重力异常模型V28.1与检核点重力异常的标准差和均值分别为4.60mGal和1.13mGal。结果表明：1′×1′海底地形-重力联合法反演模型精度接近于1′×1′海洋重力异常模型V28.1。

这是由于本发明反演过程中以高精度和高空间分辨率海底地形为基础数据，结合了高精度船测重力异常数据来构建框架，因而反演精度较高，接近于海洋重力异常模型V28.1，证明了新型海底地形-重力联合法的有效性和可靠性。但是，由于高精度船测重力异常数据属于局部分布，无法获取全球分布均匀的船测重力异常值。因此，依靠高精度局部船测重力异常数据反演全球、高精度和高空间分辨率海洋重力异常模型无法实现。鉴于此，本发明对于海底地形-重力联合法的全球应用思路进行了分析探讨。

2、新型海底地形-重力联合法的全球应用

目前公布的全球海洋重力异常模型的最高空间分辨率为1′×1′(加州大学圣迭戈分校网站)，全球海底地形最高为0.25′×0.25′(GEBCO，General Bathymetric Chart ofthe Oceans)。在精度损失较小前提下，为了获取更高分辨率全球海洋重力异常模型，本发明以0.25′×0.25′全球海底地形为基础数据，结合全球分布均匀的1′×1′海洋重力异常模型V28.1(用此模型替代高精度局部船测重力数据)，基于海底地形-重力联合法，反演可得0.25′×0.25′全球海洋重力异常模型，进而对1′×1′全球海洋重力异常模型V28.1进行加密补充。

数据区域二：研究数据范围以中国南海区域(经度112°～113°E，纬度10°～11°N)为例，GEBCO网站提供的0.25′×0.25′全球海底地形图(如图8(a)所示)，以及加州大学圣迭戈分校网站(https://topex.ucsd.edu/cgi-bin/get_data.cgi)提供的1′×1′全球海洋重力异常模型V28.1。在研究区域内1′×1′海洋重力异常模型V28.1共有原始格网点3782个，其中3242个蓝色点作为控制点，540个红色点作为检核点(如图8(b)所示)。

采用相同的迭代流程，通过比较预测值与检核点之间相关系数与差值标准差，分别绘制相关系数变化曲线图9(a)和差值标准差变化曲线图9(b)，将相关系数达到最大时的Δρ作为该区域的经验密度参数，由上述解算拟合的两条曲线可知，最终拟合曲线的相关系数达最大值0.99时，密度差参数取值约为6.3g/cm

基于海底地形-重力联合法，可反演获得0.25′×0.25′海洋重力异常模型，如图10(b)所示，并与图10(a)1′×1′海洋重力异常模型V28.1对比分析。

图11表示0.25′×0.25′海底地形-重力联合法反演模型与V28.1检核点的重力异常误差散点图，它是获取海底地形-重力联合法反演模型(图10(b))对应的540个检核点预测值，并与V28.1检核点重力异常数据之差的散点图，描述了图10(b)与V28.1检核点之间的误差。图10(b)的横坐标代表船测检核点的序号，纵坐标代表重力异常差值大小。据图10可知，海底地形-重力联合法反演模型与V28.1检核点的重力异常误差集中于0mGal附近。

据图11可知，0.25′×0.25′海底地形-重力联合法反演模型与V28.1检核点重力异常的标准差、均值、最大值和最小值分别为0.36mGal、-0.01mGal、1.21mGal和-1.29mGal。即表示反演模型与1′×1′海洋重力异常模型V28.1相比标准差约为0.36mGal，因此，新型海底地形-重力联合法为反演高空间分辨率海洋重力基准图提供了可行性和有效性。

利用图10(a)1′×1′海洋重力异常模型V28.1将0.25′×0.25′海底地形-重力联合法反演模型图9(b)中对应的格网点替换，获得图12基于海底地形-重力联合法的0.25′×0.25′海洋重力异常加密模型V28.1。即在不改变1′×1′海洋重力异常模型V28.1格网数据的前提下，利用0.25′×0.25′全球海底地形数据和全球分布均匀的1′×1′海洋重力异常模型V28.1，基于新型海底地形-重力联合法，反演得到0.25′×0.25′全球海洋重力异常模型，是对1′×1′海洋重力异常模型V28.1的加密补充。

本发明虽然已以较佳实施例公开如上，但其并不是用来限定本发明，任何本领域技术人员在不脱离本发明的精神和范围内，都可以利用上述揭示的方法和技术内容对本发明技术方案做出可能的变动和修改，因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化及修饰，均属于本发明技术方案的保护范围。

本发明说明书中未作详细描述的内容属于本领域专业技术人员的公知技术。