分类：电子论文 时间：2022-01-20 热度：586
摘要：陆地重力观测相较于航空和卫星重力观测，距离场源更近，观测精度相对较高，其静态异常和时变数据已广泛应用于研究多种地球动力学问题. 21 世纪以来，绝对重力观测技术发展迅速，陆地观测网络日益完善，高精度陆地重力观测数据产品逐渐丰富，基于这些产品的大地测量和地球物理研究不断取得新进展. 本文总结了近十几年来高精度陆地重力观测数据在大地测量和地球物理领域的应用进展情况，包括基于重力异常数据构建重力场和大地水准面模型、建立地壳物性结构模型、反演 Moho 界面形态和估计岩石圈有效弹性厚度，以及利用时变重力数据构建时变重力场模型、探测微弱动力学信号、估计地壳构造变形速率和分析与火山、地震过程的可能关联，最后探讨分析了陆地重力测量的未来发展趋势，可为中国大陆重力观测系统建设与发展规划提供参考.
关键词：陆地重力观测;重力仪;重力异常;时变重力;研究进展
0 引 言
重力场是地球的基本物理场，其静态和动态(时变)特征反映了地球内部和表层物质的分布变化和运动状态. 重力场观测数据既可为地球内部结构和地球动力学等研究提供基础观测信息(胡明城， 2003;许厚泽等，2012;孙和平等，2017)，也在经济建设和国家安全等领域发挥重要作用，广泛应用于矿产资源勘探、水资源变化监测、地震活动性分析以及国防军事等领域(Hinze et al., 2013;祝意青, 2013a;陈石等, 2015;Van Camp et al., 2017;胡敏章等, 2019;姚宜斌等, 2020).
陆地重力观测在陆地表面固定测点设置重力仪进行单次、重复或连续的重力测量工作，是获取地球重力场数据的主要手段之一(刘经南等, 2000; Niebauer, 2015). 陆地重力观测采用的重力仪有两种，分别是相对重力仪和绝对重力仪，其中前者测定重力加速度的相对差值，后者测定的是绝对重力加速度值. 实用型高精度重力仪的发展经历了一个长期过程，目前可实现 μGal 量级(1 μGal = 1× 10−8 m/s2)陆地重力测量. 1970 年代，美国 LaCoste Romberg 公司推出的金属弹簧相对重力仪将陆地重力观测精度由 mGal 量级(1 mGal = 1×10−5 m/s2)提升到几十 μGal 量级(LaCoste Romberg G 型标称精度优于 20 μGal, Crossley et al., 2013)，并且便于移动，被广泛应用于大范围重力测量(贾民育， 2000). 随着电子机械技术的迅猛发展，高精度绝对重力仪和超导重力仪也逐渐实现了商业化和便携化(Niebauer, 2015)，如美国Micro-g 公司生产的FG5 型绝对重力仪，观测精度约为 2 μGal(Francis and van Dam, 2005; Niebauer et al., 2005; Niebauer, 2015);美国 GWR 公司生产的 OSG(observatory superconducting gravimeter)标准型超导重力仪的观测精度可达 0.02 μGal，可移动 iGrav 轻便型超导重力仪的观测精度可达 0.05 μGal(Crossley et al., 2013).
根据测量方式的不同，陆地重力观测得到的数据主要有两种，一种是静态的重力异常(以下称陆地重力异常数据);另一种是重力场随时间的变化(以下称陆地时变重力数据). 陆地重力观测具有距离场源近和测量精度高的特点，观测结果包含丰富的高频信息(Saleh et al., 2012; Bomfim et al., 2013; 陈晓东等, 2020). 其中，陆地重力异常数据对研究地球内部结构和密度分布有显著贡献(Torge, 1989; 马宗晋等, 2006; Hinze et al., 2013);同时，陆地时变重力数据的空间分辨率高、测点可重复性强，适于分析区域性、近地表的质量和应力变化特征(Abe et al., 2012; Kennedy et al., 2014; 许才军和尹智, 2014; Chen et al., 2016; Van Camp et al., 2017; 陈石等, 2018)，长期观测资料也为探测地球内部微弱动力学信号提供了新的数据支持(孙和平等，2017).
随着陆地重力观测技术的迅速发展、观测网络的逐渐完善以及观测数据的不断积累，近年来相关领域基于陆地重力观测的研究成果不断涌现. 本文将简要梳理和介绍陆地重力观测在大地测量学和地球物理学(地球内部结构、地球动力学、水资源监测以及地震预测)领域的研究工作.
1 陆地静态重力数据的应用
本文中的陆地静态重力数据，是指地表观测获得的重力异常数据，主要包括：自由空气重力异常(free-air gravity anomaly, FGA)和布格重力异常(Bouguer gravity anomaly, BGA). 地面重力数据经观测改正后，再做空间改正可获得自由空气重力异常 FGA. FGA 与地形相关性强，对 FGA 进一步做层间改正和地形改正，可获得适于反映地下密度不均匀性的完全布格重力异常(郭俊义, 2001; Hofmann-Wellenhof and Moritz, 2006; Fullea et al., 2008; Pašteka et al., 2017). FGA 和 BGA 在全球变化范围从负几百到近一千 mGal(Balmino et al., 2012). 下文如无特别说明，BGA 均指完全布格重力异常.
1.1 构建重力场和大地水准面模型
陆地实测重力异常在物理大地测量学领域的一个重要应用就是构建(精化)重力场及大地水准面模型(李建成等, 2003;Hofmann-Wellenhof and Moritz, 2006;韩建成, 2012). 目前得到广泛应用的 EGM2008 全球重力场模型(完整展开到 2 190 阶，2 159 次)，空间分辨率 5′×5′，在构建时采用了陆地重力、卫星重力和卫星测高等多源数据，其中全球地面测量 FGA 数据约 130 万(Pavlis et al., 2012, 2013). EGM2008 较其上一代全球重力场模型 EGM96(完整展开到 360 阶次)在空间分辨率和精度上都有了巨大提升，在我国大陆地区的整体精度约 10 mGal(章传银等，2009). 空间分辨率达 2′×2′的 XGM2019e 全球重力场模型(完整展开到 5 400 阶次)同样采用多源数据，其中地面测量 FGA 数据约 180 万(Zingerle et al., 2020). XGM2019e 较 EGM2008 在空间分辨率和精度上均有提升，且该模型未来会继续更新(Wan et al., 2020; Zingerle et al., 2020). 目前公开发布的全球地球重力场模型都可在 ICGEM(international centre for global Earth models)官方网站(http://icgem.gfz-potsdam.de/)上找到(Ince et al., 2019)，官网模型信息页面给出了各模型所用数据源，除仅利用卫星重力数据解算的模型外，基本都采用了陆地重力数据. 基于某一地区所测陆地重力数据可以精化区域(局部)重力场模型，如 Denker(2013)、Bucha 等(2016)和吴怿昊等(2016).
陆地重力异常数据对建立高精度、高分辨率大地水准面模型至关重要. 我国 2′×2′ 重力似大地水准面 CNGG2011 在计算中采用了超过 100 万点的陆地重力异常数据(李建成，2012). 2017 年由国际大地测量协会(International Association of Geodesy)的两个工作组联合发起并组织了“Colorado geoid experiment”( Claessens and Filmer, 2020; Jiang et al., 2020; Wang et al., 2020; Işık et al., 2021)，美国 NGS(National Geodetic Survey)为此实验提供了 Colorado 地区的多种数据，其中包括近 6 万点陆地重力观测数据. Hwang 等(2020)确定中国台湾区域 30″×30″大地水准面时，采用了 6 000 多点陆地重力异常数据. 在 2019～2020 年珠峰高程测量工作中，建立作为珠峰区域高程基准的重力似大地水准面模型使用了 8 000 多点陆地重力异常数据(党亚民等，2021).
陆地重力观测目前仍存在很多数据稀疏区和空白区，以全球来看，EGM2008 收集的陆地重力数据约有 12% 的数据空白区，主要集中在非洲、南美洲和南极洲(Pavlis et al., 2012);以我国来看，西部重力数据稀疏，青藏高原存在大面积空白区(李建成，2012). 在重力稀疏区和空白区，重力场和大地水准面模型精度受到明显限制. 以我国 CNGG2011 大地水准面为例，精度评估表明其平均精度为 0.13 m，在东部陆地重力数据较为密集的 18 省平均精度为 0.07 m，西部陆地重力较稀疏的 9 省精度为 0.14 m，而存在数据空白区的西藏，精度为 0.22 m(李建成，2012). 如何克服陆地重力稀疏区和空白区的影响，是进一步精化重力场和大地水准面模型要面临的主要问题.
1.2 建立地壳物性结构模型
陆地重力异常数据 BGA 反映的物质密度分布不均匀性是重力学关注的重要信息. 重力反演方法对横向密度变化敏感，是探测地壳密度构造的有效手段之一(Nabighian et al., 2005; Hinze et al., 2013). BGA 的精度受重力测量和各项改正误差的影响. 利用 BGA 研究区域尺度(比例尺大于 1:100 万)的地壳构造和密度分布等问题时，《区域重力调查规范(DZ/T 0082—2006)》(区域重力调查规范编写组，2006)建议陆地重力测量精度不低于 0.4 mGal，BGA 精度不低于 1 mGal. 国外相关研究并无对 BGA 精度的明确建议，例如 Hinze 等(2013)仅建议进行重力勘探等研究时陆地重力测量精度应满足 0.01～1 mGal.
由于 BGA 是地球内部密度不均匀和界面起伏的综合反映，需对 BGA 数据进行处理和分析，分离出岩石圈密度结构产生的重力效应，然后构建重力异常与密度模型的对应关系，最终通过一定算法估计最优的模型参数(王新胜，2012;石磊等， 2015;李红蕾等，2017;雷晓东等，2021). 需要指出的是，由于反演问题固有的多解性，数学意义上最优的反演结果可能与真实情形并不相符，此时加入先验约束条件可以改善重力反演的效果(Li and Oldenburg, 1998; 楼海, 2001; 陈石等, 2014; 申重阳等, 2015; 张晰等, 2016; 李红蕾等, 2021). 除单独利用重力异常进行反演外，还可引入重力梯度数据以 进 一 步 提 高 反 演 结 果 的 准 确 性 和 分 辨 率(Zhdanov et al., 2004; 秦朋波和黄大年, 2016; 李红蕾等, 2017; Tian et al., 2021).
重力反演方法还适合与其他地球物理学方法(如地震学和大地电磁等方法)进行联合反演，可有 效 提 高 反 演 精 度 并 降 低 多 解 性 . Maceira 和 Ammon(2009)发展了一种基于面波与 BGA 重力异常的联合反演方法，基于此方法获得了塔里木盆地和准噶尔盆地的地壳和上地幔三维 S 波速度模型. 该模型提高了 20 km 浅层地质结构的分辨率，进而发现塔里木盆地东部和西部的地壳和上地幔 S 速度存在差异. 王新胜等(2013)利用 BGA 和地震波数据反演了青藏高原东北缘三维密度结构，发现该区域中上地壳存在有利于地震孕育和发生的高低密度异常相间分布. 郭震等(2015)利用面波和 BGA 联合反演方法得到了山西断陷带地壳三维 S 波速度结构，结果显示该区域上地壳的速度结构与地 表 地 质 构 造 具 有 较 强 关 联 性 . Moorkamp 等(2011)提出了一个用于地震、大地电磁和重力数据的三维联合反演框架，并利用模拟实验分析了联合反演中各物理参数的不同耦合方法. 国内外其他学者也开展了引入重力的地壳密度结构联合反演，取得了一系列研究成果(Bailey et al., 2012; 索奎等, 2015; Syracuse et al., 2016; Zhao et al., 2018; Afonso et al., 2019; 李海龙等, 2020; 李红蕾等, 2021).
重力反演或联合反演所采用的 BGA 数据，其计算需要已知大地水准面以上的物质密度. BGA 由 FGA 数据扣去层间改正和地形改正得到，FGA 仅与点位高程和纬度相关，但层间改正和地形改正与大地水准面和地表之间的物质密度分布密切相关. 目前一般将此密度取为地壳的平均密度，即 2.67 g/cm3 . 当研究区域较大，尤其是跨越多个构造区时，近地表的岩石密度不能简单用地壳密度平均值代替. 许多学者提出了基于 FGA 估计近地表岩石密度的方法和技术(Murata, 1993; Nawa et al., 1997; 牛源源等, 2019). 其中，牛源源等(2019)提出的基于 FGA 和贝叶斯分析的剖面近地表岩石密度估计方法，在获得较可靠的近地表岩石密度的同时，给出光滑连续分布的 BGA. 在我国云南地区的应用表明，该方法所估计的地表浅层密度分布与该地区物性资料及地质特征较为吻合. 上述工作表明，通过地表重力测量能获得较可靠的近地表岩层参考密度信息，这对于构建浅层地壳密度模型以及精化地壳模型都具有重要意义. 例如在我国川滇地区，目前已积累大量 μGal 级精度的重力探测剖面数据，若能利用更可靠的浅层密度分布来计算 BGA，将会为地壳密度结构反演提供更准确的约束，从而获得更优的反演结果.
1.3 反演 Moho 界面形态
Moho 面地球内部最重要的分界面之一，其形态特征与地壳上地幔的密度分布和均衡状态关系密切(冯锐，1985). 应用 BGA 数据反演 Moho 面形态，迄今已有大量研究(黄建平等, 2006; Shin et al., 2007; Sjöberg and Bagherbandi, 2011; Tenzer and Chen, 2014; 吴招才等, 2017; 鲁宝亮等, 2018; 王鑫等, 2020; 张杰等, 2020). 此外，也有学者采用重力梯度数据反演 Moho 面形态(Ye et al., 2016; Chen, 2017). 与反演地壳密度结构相似，以地震学方法给出的 Moho 面起伏作为约束，可获得更优的反演结果(郭良辉等, 2012; 陈石等, 2015; 张杰等, 2020).
值得指出的是，将基于 BGA 反演得到的地壳厚度(Moho 面深度与地表高程之和)与基于均衡理论给出的均衡状态下的地壳厚度进行对比，二者差异可反映研究区域构造单元的均衡状态，为分析复杂的深部动力学过程提供重要参考. 王谦身等(2009)研究了四川盆地和龙门山造山带的 Airy 均衡重力异常，发现成都平原和川西山区均处于较均衡状态、龙门山地区处于不均衡状态，这一结果反映出龙门山断裂体系深部物质的强烈交换与转移，表明该区域存在发生构造运动的趋势背景. 对于龙门山区域的上述不均衡状态，张永谦等(2010)、 Fu 等(2014)和付广裕等(2018)分别进行了验证，得出了与王谦身等(2009)一致的结论.
1.4 估计岩石圈有效弹性厚度
岩石圈有效弹性厚度(effective elastic thickness, Te) 是 表 征 岩 石 圈 力 学 强 度 的 一 个 定 量 指 标(Forsyth, 1985). 定义一个覆盖在非黏性流体之上的弹性板，其在相同载荷下可产生与真实岩石圈相同的弯曲，则该弹性板的厚度即 Te(Burov and Diament, 1995; 郑勇等, 2012). Te 可反映岩石圈在地质时间尺度的长期载荷作用下抵抗变形的能力(Chen et al., 2013; 胡敏章等, 2015)，Te 越大意味着岩石圈抵抗变形的能力越强，即越难发生挠曲变形;反之，岩石圈则越易发生挠曲. 对于本文关注的陆地区域，通过估计不同构造和地质演化单元的 Te，可为认识大陆岩石圈的长期变形、理解构造运动(如克拉通破坏和青藏高原隆升等)的动力学机制研究提供参考(陈石等, 2011, 2013; 胡敏章等,
岩石圈在上部地形负载作用下发生挠曲变形的直接体现是引起重力异常，所以重力异常和地形数据成为估计 Te 的理想数据. 目前基于重力异常和地形数据估计陆地区域 Te 比较常用的是谱方法(Kirby, 2014; 胡敏章等, 2020)，如 BGA-地形相干性方法(Forsyth, 1985; Pérez-Gussinyé et al., 2004; Chen et al. 2013; 李永东等, 2013)、FGA-地形导纳法等(陈石等, 2011; 佘雅文等, 2018; 付广裕和王振宇, 2020)、Moho 地形导纳法(杨亭等, 2013)，分别采用导纳法和相干性方法(Pérez-Gussinyé et al., 2009; Chen B et al., 2018)以及导纳法和相干性方法联合反演(Audet, 2019; Lu et al., 2020; Shi et al., 2020). 针对谱方法中存在的频率泄露问题， Simons 等(2000)以及 Kirby 和 Swain(2011)分别提出了基于多窗分析和 Fan 小波分析的改进方法. 除谱方法外，部分学者通过直接求解岩石圈在地形载荷作用下发生挠曲形变的偏微分方程来计算 Te，如空间域的有限差分法(Jordan and Watts, 2005; 姜效典等, 2014; 胡敏章等, 2020).
目前对于中国大陆地区 Te 的计算，由于所用重力数据、处理过程和估计算法上的不同，不同研究给出的 Te 存在较大差异，并且所得结果的空间分辨率较低，普遍采用的谱方法也不利于对 Te 的横向变化特征展开研究(胡敏章等，2020). 还需指出的是，目前中国区域 Te 计算大多基于重力场模型数据计算，缺乏地表实测重力异常的约束;如果采用基于地表实测重力的 BGA 或 FGA 数据来估计 Te，预期可提高估计精度，为基于 Te 的地球物理解释和应用提供更加精确可靠的数据支持.
2 陆地时变重力数据的应用
构成地球系统的大气、海洋、地幔和地核(液态外核与固体内核)等圈层，其内部物质在日月等天体引力及各种构造活动的作用下不停地迁移运动，导致地球重力场产生潮汐和非潮汐变化(孙和平等， 2017). 图 1 给出了不同周期、振幅特点的潮汐和非潮汐重力场信号时频特征(Crossley et al., 2013). 其中，地壳长期形变和地震过程等与地壳内部变形相关的非潮汐时变信号多具有低频特征，量级约 10 μGal 甚至更小. 由于量级微小，重力场时变信号的观测依赖于高精度重力仪的革新.
时变重力数据获取的方法主要有两种：一是通过在固定台站设置重力仪连续观测获得连续重力变化;二是通过对野外固定点的定期复测获得重力变化信息. 前者具有观测精度高(μGal 级)和时变采样率高(可达秒采样)的优势，但建站成本高、台站数量少;后者观测精度在 10 μGal 级且时间采样间隔长(半年至多年)，但设站较为自由、观测密度高，有利于区域性重力变化信息的获取. 目前在连续重力观测方面，以 gPhone 弹簧型和 GWR 超导型为代表的两种相对重力仪是固定台站常用设备，其中 GWR(OSG 型和 iGrav 型)也是当前最高精度的重力观测设备，观测精度优于 0.05 μGal，采样间隔最高可达 1 s(Crossley et al., 2013);在固定点定期复测方面，常用设备包括 FG5 型和 A10 型绝对重力仪，以及 LaCoste Romberg(G)型、 CG5 型和 Burris 型弹簧型相对重力仪(赵云峰等， 2018). 其中，FG5 和 A10 的重力测量精度分别可以达到 2 μGal 和 10 μGal(Crossley et al., 2013).
2.1 构建时变重力场模型
我国陆地重力观测系统多以绝对和相对重力仪组网、定期复测的形式获得重力场随时间变化，例如“中国地壳运动监测网络”重力观测系统(许厚泽，2003;祝意青等，2012). 各期重力观测数据经各项改正和重力网平差后，通过不同差分方式(两期差分或以某期为基准差分)获得重力场时变信息(李辉等，2009;祝意青等, 2012, 2013b)，可为区域构造活动分析、地震风险性评估、地震前兆特征提取以及水资源变化监测等研究提供数据支持. 基于以上信息建立时变重力场模型时，国内学者目前主要采用空间域方法(申重阳等, 2009; Shen et al., 2015; Hao et al., 2021)，如最小二乘配置法(祝意青等，2013a，2017;徐伟民等，2021)、 Kriging 方法(祝意青等，2015，2016)和等效源方法(张贝等，2021).
除空间域方法外，谱方法也逐渐被应用到时变重力场建模研究中. 玄松柏等(2012)基于小波变换方法获得了 2000～2007 年中国大陆重力场时变信息的多尺度分解，给出了不同小波尺度下的重力场时变信息及其与构造活动的关联. 陈石等(2018)将 Slepian 局部谱分析方法引入到地面时变重力场建模中，利用“中国地壳运动监测网络”2005 年和 2008 年两期重力观测数据，建立了我国大陆区域 72 阶的重力场模型. 韩建成等(2021)基于 Slepian 方法和我国华北重力测网 2011～2013 年期间多期重力数据，获得了华北地区 120 阶(约 150 km 空间分辨率)的 2011～2013 年变化重力场模型，并与研究区域内基于 Slepian 方法和 GRACE (Gravity Recovery and Climate Experiment)卫星数据得到的重力变化进行了比较. Han 等(2021)将 Slepian 方法应用到我国川滇地区，获得了该地区 2015～2017 年 120 阶时变重力场模型，发现 2015～ 2017 年同期差分重力变化受西昌—玉溪一带的近南北向构造和北西向的红河断裂控制.
图 2 所示即为韩建成等(2021)基于 Slepian 局部谱分析方法获得的 120 阶华北地区 2011～ 2013 年尺度重力变化，所恢复重力场信号较好地都集中在研究区域(洋红色闭合曲线)内部.
重复重力测点位置空间分布不均匀是目前时变重力场建模面临的主要限制. 在点位分布稀疏的地区，不仅重力测网对异常场源的分辨能力受到限制，重力场建模的精度也受到明显影响，韩建成等(2021)报告了在华北地区基于 EGM2008 模型的模拟实验结果，以研究区域内 EGM2008 模型 120 阶 FGA(变化范围约 58 mGal)模拟重力变化，然后基于 Slepian 方法解算 120 阶模型;解算得到的 120 阶模型的重力变化范围约 55 mGal，解算误差为 2.50 mGal，在研究区域内观测点位密集的京津地区，解算误差为 2.05 mGal;在观测点位较为稀疏的河南、安徽和山东区域，解算误差显著增大至 2.98 mGal. 综上，在测点分布稀疏地区对重力测网进行加密，将大幅提高重力场建模的精度.——论文作者：韩建成1,2，陈 石1,2*，李红蕾1,2，张 贝1,2，卢红艳1,2，侍 文1,2，徐伟民1,2，贾路路1,2
文章名称：陆地高精度重力观测数据的应用研究进展
