中国大陆岩石圈震磁异常年变统计分析*

蔡苏苏，陈 斌

(中国地震局地球物理研究所，北京 100081)

0 引言

地震的孕育是一个复杂的过程，其影响因素繁多，其中地磁场异常变化与地震孕育过程具有关联性。我国正式开展广泛和系统的震磁观测与深入研究，始于20世纪50—60年代(丁鉴海等，2011)。许多科学家针对地震与震中附近磁异常的相关性进行了研究，如李家发等(1992)提出地震多沿地磁异常及梯度带、磁异常的端点转折交汇处分布。研究表明地震震中与单个要素的零值线存在相关性，若干与地震相关的地球内部磁场异常特征也被提出，如内源场的磁偏角零值周围区域往往易发地震(哈奇基扬，2008；顾春雷等，2012)；岩石圈磁场的要素与地震具有相关性：如魏荣强和于磊(2012)提出中国大陆及周边地区大多数≥5.0地震震中位于垂直分量的零值线上或附近，李琪等(2006)研究得出1998年张北6.2级地震震中基本上位于0～3 nT的分量等值线上。但Lei等(2018)分析认为-5～-3 nT才是地震高发区，该区域聚集了2004—2007年53.2%的≥5.0地震。

岩石圈磁场变化量与地震对应性的研究也呈现不同震例异常对应一个或多个不同地磁要素的现象。比如2013年芦山7.0地震对应的异常要素为地磁场总强度、垂直分量、地磁场水平分量以及磁偏角(倪喆等，2014a，b)；2015年九江—瑞昌5.7地震震中100 km附近存在明显的水平分量异常(顾左文等，2006)，震中对应的异常要素为磁偏角和磁倾角(顾春雷等，2010)；汤筱麒等(2019)研究认为甘川交界地区≥5.0地震对应的异常应为水平分量的磁弱变区域边界。

总结前述研究发现，各震例异常分析针对的地磁要素并不完全相同，即便是对同一地磁要素的分析，异常判定指标阈值也不尽相同。寻求与地震存在相关性较好的地磁异常对应的地磁要素及异常的某种定量描述，是地震地磁研究需要解决的问题，也是进一步提出异常机制解释的基础。

前人在重力、地震相关研究方面提出了一些统计和检验方法值得借鉴：如对研究区域的目标数据进行重力数据网格化(贾民育等，1995，1996，2000)，分析地震震中所在网格的数据特征量，得出≥5.0地震主要发生在重力变化正—负转换的高梯度带上(胡敏章等，2015)。在中长期地震危险性分析中，蒋长胜等(2011)将青藏高原东北缘地区进行了网格间距尺度为0.2°的划分，并运用算法和算法对各网格进行回溯性预测。

考虑到岩石圈磁场分布的复杂性、地震震中位置与各个地磁要素异常特征的相关性、统计方法的多样性，本文对2010—2020年中国大陆地区的岩石圈磁场年变数据进行区域网格化划分，并针对网格中心点对应的≥6.0地震相关地磁要素进行相关统计与计算分析。

1 岩石圈磁场数据处理

1.1 数据来源

中国地震局流动地磁团队每年进行一次地磁矢量地面测量，测量的地磁要素为地磁总强度、磁偏角和磁倾角，经过换算可以得到其他4个地磁要素的数值，即地磁场北向分量、东向分量、垂直分量和水平分量。测量得到的地磁场数据包括外源场、主磁场和岩石圈磁场，因此需对野外获得的磁场数据进行数据处理，以去除岩石圈磁场以外的外源变化磁场以及主磁场部分的影响。本文将2010—2020年相邻两年的中国大陆(19.5°～44.5°N，95.5°～126.5°E)岩石圈磁场数据进行时序差分获得这一时间段内的岩石圈磁场年变数据。笔者从中国地震台网中心下载2010—2020年≥6.0的地震目录，从中删除处于地面矢量地磁测网有效监测范围之外的地震事件及所有余震事件，最终共选取15次≥6.0地震事件进行研究。

1.2 地磁数据处理

原始观测数据是地球磁场内源场和外源场的总和，地球磁场内源场包括起源于核幔边界的主磁场以及起源于地壳居里点温度以上的岩石圈磁场；外源场包括电离层场、磁层场及感应场(徐文耀，2003)。要想获得岩石圈磁场变化数据，需要对观测数据进行处理，去除主磁场和外源场的影响。野外观测数据处理一共分为3步(董超等，2021；宋成科，2021)：

(1)日变通化改正。为消除流动地磁观测数据包含的地磁场日变化等外源场成分，利用邻近台站的连续观测地磁场分钟值对数据进行日变通化改正(陈斌等，2017)，测量时间段内通化零时为地方时00∶00—03∶00。

(2)长期变化改正。为消除观测数据中的地球主磁场长期变化成分，采用自然正交分量模型方法(陈斌，2011)，对每期日变通化数据集进行长期变化改正，改正日期为上一期的日变通化日。

(3)主磁场剥离。以IGRF-13为观测区域地球主磁场参考场模型(Alken，2021)，对应各测点位置的日变通化结果减去主磁场部分，得到当期岩石圈磁场结果。将相邻两期岩石圈磁场相减，可获得岩石圈磁场空间变化模型。

1.3 地磁矢量场数据网格化

基于上述数据处理过程得到的岩石圈磁场数据，收集2010—2020年的岩石圈磁场年变数据并运用曲面样条插值方法(余志伟，1987)生成0.1°×0.1°的岩石圈磁场年变数值。将研究区域(18°～54°N，73°～135°E)划分为0.5°×0.5°的网格区域，以每个网格的中心点代表该网格，生成网格中心经纬度数据点，再从中筛选出有地震的网格中心经纬度数据点。计算岩石圈磁场年变数据与网格特征值的梯度，再挑选出≥6.0地震的震前1年的岩石圈磁场数据(表1)进行研究。

分别定义水平矢量模||表示为：

(1)

水平矢量模的梯度|grad|(周建兴，2008)表示为：

(2)

(3)

(4)

矢量模||表示为：

(5)

矢量模的梯度|grad|表示为：

(6)

(7)

(8)

式中：|grad|、|grad|分别为直角坐标系下两矢量模的梯度值；||为水平(东西)向水平矢量模的梯度值；||为垂直(南北)向水平矢量模上的梯度值；||为水平方向上矢量模的梯度值；||为垂直方向上矢量模的梯度值；Δ、Δ、Δ为北向、东向和垂直分量的年变数值；、为两个方向的间距，此处均为0.5°。||的首列数据为原矩阵第二列与第一列数据之差除以，最后一列则为最后两列之差除以；同理，可以得到||、||、||。

2 地震前一年岩石圈磁场磁异常变化特征分析

2.1 地磁矢量异常变化特征分析

2010—2020年≥6.0地震共15次，挑选网格中震前1年的岩石圈磁场年变数据(如2013年4月20日7.0芦山地震，挑选2011—2012年的岩石圈磁场年变数据)，绘制了所有年份对应的岩石圈磁场空间等变线分布图。本文以2011—2012年≥6.0地震的磁偏角、垂直分量、总强度、水平矢量模||、矢量模||及其梯度的等值线图(0.5°×0.5°)为例进行说明。图1、2显示≥5.0地震事件共8次，其中5.0～5.9地震事件6次、6.0～6.9地震事件1次(即岷县—漳县6.6地震)、7.0及以上地震事件1次(即芦山7.0地震)。

图1a、b分别为磁偏角年变率及其梯度的空间等值线分布，等值线间隔分别为0.1′、0.1′/°，图中的黑色实线分别为磁偏角的零值线与梯度值为0.5′/°的等值线。岷县—漳县6.6地震(图中紫红色点)与芦山7.0地震(图中红色点)的震中均位于磁偏角的零值线附近(图1a)及磁偏角梯度值为0.5′/°的等值线附近(图1b)。

图1c、d分别为2011—2012年总强度年度变化及其梯度的等值线图，等值线间隔分别为2 nT、2 nT/°，图中的黑色实线分别为总强度的零值线与梯度值为10 nT/°的等值线。岷县—漳县6.6地震与芦山7.0地震震中位于总强度的零变线附近(图1c)；图1d显示总强度的梯度值介于0～42 nT/°，芦山地震与岷县—漳县地震震中位于10 nT/°等值线附近，震中周围200 km附近区域总结度梯度值比较小，属于低值区；岷县—漳县地震震中东北部100 km附近的区域总结度梯度值为15～35 nT/°，属于中高值区。

图1e、f分别为震前垂直分量年度变化及其梯度的等值线图，等值线间隔分别为2 nT、2 nT/°，图中的黑色实线分别为垂直分量的零值线与梯度值为10 nT/°的等值线。岷县—漳县6.6地震与芦山7.0地震的西北方向为负值区，东南方向为正值区，两片区域大致以零变线为轴线呈对称分布，其震中位于垂直分量正负异常区交界附近(图1e)。图1f显示震中位于垂直分量梯度的10 nT/°附近，芦山7.0地震震中附近100 km仍为低梯度带，无异常梯度带出现。

图1 2011—2012年磁偏角(a)、磁偏角梯度(b)、总强度(c)、总强度梯度(d)、垂直分量(e)和垂直分量梯度(f)等值线图

2.2 水平矢量模与矢量模及梯度特征

图2a、b为2011—2012年水平矢量模||及其梯度的岩石圈磁场等变线图，等值线间隔分别为2 nT、2 nT/°，图中黑色实线分别为水平矢量模||及其梯度分别为10 nT、10 nT/°的等值线。≥6.0地震震中均位于图中黑色实线附近。水平矢量模||在芦山7.0地震震中周围100 km左右均大于10 nT，在小于100 km的区域小于10 nT，但震中附近并未出现高梯度带。岷县—漳县6.6地震震中东部存在分块的中高值区，其震中附近的水平梯度值与梯度值分别小于10 nT与10 nT/°。

图2c、d为2011—2012年矢量模||及其梯度的等值线图，等值线间隔分别为2 nT、2 nT/°，图中黑色实线分别为矢量模||及其梯度分别为10 nT、10 nT/°的等值线。由图2 c与图2a对比可知，水平矢量模和矢量模的≥6.0地震震中附近的值均为10 nT左右，但矢量模||加入垂直分量后其高值异常区更为明显。图2b与图2d显示≥6.0地震震中附近并未出现高梯度带。

图2 2011—2012年水平矢量模(a)、水平矢量模梯度(b)、矢量模(c)和矢量模梯度(d)等值线图

从图2与表1的数据可以看出，统计10年以来的水平矢量模和矢量模梯度均值略小于10 nT/°，与对应的2011—2012年水平矢量模和矢量模的梯度特征值有2～4 nT/°的误差，但若单独考虑水平矢量模和矢量模的异常梯度形态与地震震中的位置则较难统一量化，而图2震中附近两要素模大致在10 nT。

对比分析2011—2012年≥6.0地震地磁要素等值线图，发现震中主要位于磁偏角、总强度、垂直分量的零变线附近，震中附近水平矢量模和矢量模值在10 nT左右；对应的梯度值并未出现较大的高梯度带，显示高梯度带与地震震中的对应关系并不明显。

2.3 震例统计结果

表1为2010—2020年15次≥6.0地震的磁偏角、总强度、垂直分量、水平矢量模||与矢量模||及各要素梯度所对应的震前岩石圈磁场年度变化(0.5°×0.5°)。

表1 2010—2020年MS≥6.0地震地磁场要素及梯度网格中心值统计

磁偏角的值介于-2.26′～2.02′，磁偏角梯度值介于0.29′/°～1.49′/°，均值分别为-0.31′与0.80′/°，均方差值分别为1.05′与0.33′/°。40%的磁偏角的值分布在±0.5′之间，66.67%的梯度值分布在0.5′/°～1′/°。

总强度的值介于-15.16～22.46 nT，其梯度值介于0.53～10.33 nT/°；53.33%的总强度的值分布在±4 nT之间，13.33%的总强度的梯度值分布在8～10 nT/°。垂直分量的值介于-16.13～12.16 nT，其梯度值介于2.13～14.63 nT/°；53.33%的垂直分量的值分布在±6 nT之间，6.67%的垂直分量的梯度值分布于8～10 nT/°。

水平矢量模||的值介于1.75～22.40 nT，其梯度值介于1.58～12.22 nT/°；20.00%的水平矢量模||分布在10～12 nT，26.67%的梯度值分布于8～10 nT/°。矢量模||介于3.78～25.49 nT，其梯度值介于3.43～13.16 nT；26.67%的矢量模||分布在12～14 nT，20.00%的梯度值分布在8～12 nT/°。

3 结论

本文选取2011—2012年中国大陆(19.5°～44.5°N，95.5°～126.5°E)的岩石圈磁场年变数据，基于网格划分，统计分析各要素值的零值线以及各要素梯度值与≥6.0地震震中位置的相关性，发现地磁要素零值线与震中位置相关性各异，要素值的高梯度带仅存在于个别震例中，主要结论如下：

(1)2011—2012年岩石圈磁场等变线图中，≥6.0地震震中位于磁偏角，总强度，垂直分量零值线附近，水平矢量模||与矢量模||的值均为10 nT左右。

(2)15次≥6.0地震中，分布在磁偏角的值为±0.5′之间的震例占比为40%，分布在总强度的值为±4 nT之间的震例占比为53.33%，分布在垂直分量的值为±6 nT之间的震例占比为53.33%，分布在水平矢量模||为10～12 nT的震例占比为20.00%，分布在矢量模||为12～14 nT的震例占比为26.67%。

(3)66.67%的震例分布在磁偏角梯度值为0.5′/°～1′/°的地区，13.33%的震例分布在总强度梯度值为8～10 nT/°的地区，6.67%的震例分布在垂直分量梯度值为8～10 nT/°的地区，26.67%的震例分布在水平矢量||模梯度值为8～10 nT/°的地区，20.00%的震例分布在矢量模||梯度值为8～12 nT/°的地区。梯度数据分布太分散，导致高梯度带与地震震中对应并不准确。

统计2010—2020年研究区≥6.0地震震中位置与各要素之间的相关性发现，总强度、磁偏角、垂直分量单个地磁要素的零值附近指标较其他地磁要素对地震的指示更显著，但震例比例仅为50%左右，仍不太理想，这或许与仅用网格中心点作为整个网格特征值的近似处理过于简单有关。与过去的认识不同的是，高梯度带的震例比例远较预想低，一方面可能是网格特征值的简单处理引起，但另一方面可能是一个客观存在，还需要继续分析这一现象。下一步将研究网格特征值的选取方法，并综合多要素联合分析以期得到更有效的异常特征指标。

中国地震局流动地磁观测团队提供了岩石圈磁场数据以及开展了前期集中数据处理工作，中国地震台网中心提供了地震目录，审稿专家对本文提出了宝贵修改意见，在此一并表示感谢。