论文信息:Zhang, Y., Xu, Y., Mooney, W. D., & Chen, C. (2021). Local separation of potential field anomalies using equivalent sources: application for the 3-D structure of mantle uplift beneath Von Kármán crater, the Moon. Geophysical Journal International, 227(3), 1612–1623. https://doi.org/10.1093/gji/ggab307

位场(重力与磁力)观测是研究地下密度与磁性结构的重要途径,但在实际资料中,浅部小尺度目标的异常总是叠加在深部大尺度场的背景之上,如何把它们可靠地分离开来一直是位场资料处理与解释中的经典难题。本文以一作张壹为代表的研究,提出了一种基于等效源技术 + 迭代反演的位场异常局部分离新方法,并将其应用于月球 Von Kármán 撞击坑(VKC,嫦娥四号着陆区)下方地幔隆起的 3D 结构重建。下面围绕这篇 Geophysical Journal International 论文,梳理其动机、方法要点与主要结论。

为什么要做"局部分离":区域—剩余异常分离的局限

观测位场数据是多个互相叠加的地下源体所产生的重力或磁效应之和,各源体的异常往往难以单独提取。传统上,我们把反映深部大尺度背景的异常称为区域异常,把浅部小尺度目标所对应的异常称为剩余异常,二者的分离(regional–residual separation)是位场解释中极有价值的一步。几十年来发展出了大量方法,大体可分两类:

  • 基于数据的方法:利用区域与剩余异常在图形特征或频率上的差异来分离,例如多项式拟合(Beltrao et al., 1991)、最小曲率法(Mickus et al., 1991)、有限元(Mallick & Sharma, 1999)、三维主成分与纹理分析(Zhang et al., 2009)、双维经验模态分解(Hou et al., 2012)、维纳滤波(Pawlowski & Hansen, 1990)、独立成分分析(Forootan & Kusche, 2012)以及各类小波分析(Fedi & Quarta, 1998;Xu et al., 2009)等。
  • 基于模型的方法:利用位场的固有特性来分离,例如向上延拓作为标准分离滤波(Jacobsen, 1987;Zeng et al., 2007)、优化滤波参数(Pilkington & Cowan, 2006)、以及基于**等效层(equivalent layer)**概念的带通滤波等(Pawlowski, 1994;Guo et al., 2013)。

这些方法虽多,却普遍存在几个不足:其一,对被分离异常的质量缺乏可量化的评价指标,选择何种分离结果仍带有较强的主观性;其二,无论是小波分析还是球谐分析本质上都是信号处理手段,所生成异常的空间物理合理性(physical plausibility)并不总能得到保证。更重要的是,正如 Li & Oldenburg (1998) 指出的,以往绝大多数分离方法针对的是埋深明显不同的源体所产生、波长差异较大的异常;而当我们想要从观测中分离出波长相近的局部异常时(例如几处相邻源体的异常叠加),基于频率差异的方法便难以奏效,这正构成了本文要解决的"局部分离"问题。

方法:用等效源与迭代反演实现局部波长分离

受 Dampney (1969) 启发,观测位场异常可以用一组等效源来合成。这些等效源通常是以一层(或多层)离散模型单元呈现,通过反演确定单元的物理性质,使正演出的数据能复现观测值。由于等效源技术只依赖位场固有的正演—反演关系,它既能在笛卡尔坐标下用直角棱柱单元实现,也能在球坐标下用球面三角棱柱等单元实现(von Frese et al., 1981; Li et al., 2019)。

本文的分离思路相当直观(见其论文 Figure 1):给定观测异常 d 及其等效源 m,把等效源分成两组——局部源 mloc(大致置于目标异常 dloc 之下)与剩余源 mrem(其余部分),于是 d = dloc + drem。直接的做法是先反演出等效源,再分别正演局部源与剩余源所对应的异常,即可得到分离结果。

然而这里有一个关键难点:等效源的物理性质与实际位场源体并不相同(例如 Figure 1 中的两个多面体源),因此所选定的局部源 mloc 往往无法恰好只包含正确局部异常所需的信息,分离出的"局部"与"剩余"异常通常会各自混入对方的一部分,直接分离精度不足。为此,作者引入了一个迭代精化流程(见其 Figure 2):把上一步分离出的 dloc、drem 分别作为新的观测,再次进行分离,得到各自的局部与剩余分量(dloc_loc、dloc_rem、drem_loc、drem_rem),再按照"局部总=局部之局部+剩余之局部"的方式重组,如此反复迭代;当 dloc_rem 与 drem_loc 相等(即二者之差满足给定的数据不确定度判据)时停止,此时分离不再能被改进。

与之配套,作者还给出了用于确定变尺寸等效源的反演方法。反演目标函数为数据拟合项 ‖Wd(Gm − dobs)‖²,并用预条件共轭梯度法求解;等效源的构建上,根据观测异常横向梯度的缓急,在小梯度处用较大的模型单元、大梯度处用较小的模型单元(可借四分树算法实现),既保证了高精度又控制了反演规模。同时,作者引入了一个定量评价量 Δd(=分离过程中未被确定部分异常的 RMS),用以量化分离的完整度:Δd 越小,说明分离越完整;不同等效源几何与局部分割方式的组合都可以用 Δd 来客观比较,从而最大程度地减少了主观性。

合成算例:对磁法异常与重力异常的分割验证

在应用于实际数据之前,作者先用两组合成算例检验了方法。

  • 磁法算例:采用两个矩形块体 M1、M2(Figure 3)正演出合成 ΔT(总磁力异常)。由于两个异常在正方向几乎被叠加混在一起,无法从观测上直接分辨(其 Figure 5a)。作者把 M1 的异常视为局部异常、M2 的异常视为剩余异常,通过多组等效层中心深度与厚度的组合(Table 1)筛选,以最小的 Δd = 2.93 nT 完成分离。分离出的局部与剩余异常都能较好地复原合成数据,差异主要出现在局部异常负值部分及局部/剩余源的边界附近,二者的 RMS 差异约 7 nT。
  • 重力算例:采用五个块体模型(Figure 6)构造更复杂的合成重力异常。把 1 号和 4 号块体的异常视为局部异常、其余三块视为剩余异常,通过鞍部区域划定局部源范围,以最小 Δd = 5.94×10⁻³ mGal 完成分离。分离异常与合成异常的均值差约 2×10⁻² mGal、RMS 差约 3.1×10⁻² mGal,同样主要误差集中在局部/剩余源的边界。

两组算例都表明,该方法能够在波长相近的情况下以较低误差把任意形状的局部异常从整体中提取出来,Δd 也提供了客观的定量评价标准。

应用:Von Kármán 撞击坑下地幔隆起的 3D 结构

方法建立并验证之后,作者将其应用于月球 Von Kármán 撞击坑(VKC) 区域的局部重力异常分离与地幔隆起反演。VKC(176.3°E、44.45°S,直径约 186 km)是位于月球背面 南极—艾特肯(SPA)盆地西北部的尼普敦纪(约 3.9 Ga)撞击坑,也是**嫦娥四号(CE-4)**于 2019 年 1 月 3 日在月球背面软着陆的着陆区;其南部叠加了更老的 Von Kármán M 撞击坑(VKMC,176.2°E、47.2°S,直径约 225 km)。月球地幔密度高于上覆地壳,巨型撞击坑通常表现为正重力异常并对应撞击引起的地幔隆起,因此这些异常为研究地幔隆起的空间结构提供了宝贵窗口。

数据处理流程为:利用月球重力模型 JGGRX 1500E(Lemoine et al., 2013)提取研究区重力扰动,用多分辨率月表高程模型 LRO LTM05(Smith et al., 2010)正演地形成分(取平均地壳密度 2856 kg·m⁻³)并从观测中扣除,得到布格重力异常(图 10);再通过二阶多项式拟合扣除长波长区域趋势,得到剩余布格重力异常(图 11)。在此基础上,应用本文提出的等效源局部分离方法,从剩余布格异常中把与 VKC、VKMC 下方抬升地幔相关的局部重力异常(最大约 165 mGal)提取出来(图 13a),而 Δd 约为 0.06 mGal,远小于局部异常本身,可忽略不计。

位场局部分离结果
图1 剩余布格重力异常的局部分离结果(引自论文 Figure 13)

(a) 分离出的与抬升地幔相关的局部重力异常(最大约 165 mGal);(b) 分离出的区域(剩余)重力异常;(c) 分离过程的 Δd 异常(量值微小,约 ±1.0 mGal)。

把分离出的局部异常进一步用球坐标下的密度界面反演方法(Zhang et al., 2019)进行反演,即可重建地幔隆起的 3D 结构。反演时假设平均 Moho 深度 43 km、平均地幔密度 3356 kg·m⁻³,反演出的地幔隆起不确定度约 16 m。如图 2 的 3D 视图所示,VKC 南缘下方(原 VKMC 的中心位置)存在一个明显的主地幔隆起,其平均抬升深度约 20 km、边缘最大可达 10 km;隆起周围还观察到可能与撞击事件相关的、最大深度约 50 km 的地壳隆起(crustal bulge),符合典型的峰环盆地(peak-ring basin)地壳结构特征。最值得注意的是,在主隆起北侧、靠近 VKC 中心的位置,还存在一个次级地幔隆起分支

地幔隆起3D结构
图2 反演得到的 Von Kármán 撞击坑下方地幔隆起的 3D 结构(引自论文 Figure 14,为论文核心结果)

(a) 由重力数据反演得到的抬升地幔的 3D 视图,虚线为剖面 AA′ 与 BB′ 位置;(b)(c)(d) 分别为从南、北、西三个方向的侧视图。可见主隆起(原 VKMC 中心)及北侧靠近 VKC 中心的次级隆起分支。

这一"主隆起 + 次级隆起"的双重结构,为 VKC 地区经历了两期地幔隆起事件提供了直接证据:主隆起对应形成 VKMC 的撞击事件,而次级隆起对应其后形成 VKC 的撞击事件。尽管两个盆地直径相近,次级隆起却显著弱于主隆起,作者把这归因于两次撞击的**前存条件(pre-existing conditions)**不同——在 VKMC 撞击事件后,地壳结构得到加强、热状态演化改变(可能阻碍后续撞击过程中的地幔上涌),从而削弱了 VKC 撞击引起的隆起幅度。主隆起中心的地幔高程反而较低、边缘较高,说明其结构较为复杂,其形成可能涉及 VKMC 撞击对原有地壳结构的破坏与重组、VKC 撞击的短期效应以及盆地的长期撞击后均衡调整等多重过程。

结论与讨论

本文的主要贡献有三:其一,提出了一种基于等效源技术的局部分离方法,通过迭代反演优化分离的局部异常,并引入定量评价量 Δd 对分离质量进行客观评估,相较于小波分析或局部球谐分析等纯信号处理手段,其分离结果因由密度变化正演而来而始终物理上合理;其二,给出变尺寸等效源的求解算法,兼顾精度与计算效率;其三,成功将该方法应用于月球 Von Kármán 撞击坑的局部重力异常分离与地幔隆起反演,重建出对应 VKMC 与 VKC 两期撞击事件的主、次级地幔隆起 3D 结构。

作者也坦诚地指出了方法的局限:方法并非完全客观——仍需由使用者确定目标等效源的大致范围(不过位场异常的水平边界相对容易识别,且可用 Δd 客观地选择优选结果);此外,考虑到重力源的非唯一性,特征相似的异常可能由不同几何与物理性质的源体产生,因此在解释时还需结合其他地质与地球物理信息加以验证。对于地幔隆起的成因机制,作者亦指出需要撞击动力学模拟来更全面地理解撞击后的地幔结构演化。作为嫦娥四号着陆区,VKC 的深部地壳结构不仅反映了该区域的演化历史,也为未来 CE-4 数据的联合解释提供了重要背景。

延伸阅读与引用

  • Li, Y., & Oldenburg, D. W. (1998). Separation of regional and residual magnetic field data. Geophysics, 63(2), 431–439.
  • Zhang, Y., Mooney, W. D., Chen, C., & Du, J. (2019). Interface inversion of gravitational data using spherical triangular tessellation: an application for the estimation of the Moon’s crustal thickness. GJI, 217(1), 703–713.
  • Zhang, Y., Mooney, W. D., & Chen, C. (2018). Forward calculation of gravitational fields with variable resolution 3D density models using spherical triangular tessellation: theory and applications. GJI, 215(1), 363–374.
  • Dampney, C. N. G. (1969). The equivalent source technique. Geophysics, 34(1), 39–53.
  • Li, D., Liang, Q., Du, J., Sun, S., Zhang, Y., & Chen, C. (2019). Transforming total-field magnetic anomalies into three components using dual-layer equivalent sources. GRL, 46. https://doi.org/10.1029/2019GL084607

本文为对该论文的中文解读,图片均引自论文原文。