0 引言

岩石压力实验表明，岩石的电阻率随压力的改变有显著变化。一般认为，在地震发生之前存在一个应力积累的过程，在积累过程中，地下物质的电导率逐渐发生变化。冯志生等(2004、2009)基于磁测深原理提出了地磁谐波振幅比方法，总结了中强震、强震前地磁谐波振幅比的异常特征(冯志生等，2004、2009、2013)。研究表明，地震发生前后地磁谐波振幅比的趋势性异常变化特征与直流地电阻率的特征相似，表现为下降变化，并有转折—恢复的变化过程特征，地震基本发生在趋势性下降变化后的转折或恢复过程中。5级左右的地震异常持续总时间为1~2年，震前异常为0.5~2年；而6级左右地震的异常持续总时间为2~3年，震前异常为1~2年；7级以上地震的地磁谐波振幅比异常震前持续时间均在4年以上。异常最关键的特征为同一周期2个方向出现不同步变化，并至少有2个周期出现不同步变化。同时，异常低点时间具有随周期迁移的特征。

地磁谐波振幅比方法对实际地球介质的假定和约束较少，目前对于地震孕育过程中深部电阻率异常变化的机理并不十分清楚，在深部电阻率异常变化所引起的变化磁场异常可能很微弱的情况下，该方法相对电磁感应的各种正反演方法和地磁转换函数来说可能具有一定意义(冯志生等，2004)。目前，国内越来越多的学者参与到地磁谐波振幅比的研究工作中，并总结积累了更多的震例(李琪等，2016；蒋延林等，2016；戴勇等，2016；倪晓寅等，2016；何畅等，2017；艾萨等，2017；刘长生等，2017)。然而，这些研究工作基本均为单台异常特征分析总结，对震中附近多个台站的异常时空变化特征综合分析工作则开展较少。

2017年11月23日重庆武隆发生5.0级地震，震中附近多个台站出现异常。本文分析了这些台站异常的时空变化综合特征，研究结果作为谐波振幅比方法的重要进展，对今后应用该方法开展地震分析预报工作具有一定借鉴意义。

1 资料、基本原理及计算方法

武隆地震500km范围内的地磁台站分别为涪陵江东台、武隆仙女山台、石柱黄水台、万州天星台、奉节荆竹台、巫山建坪台、四川成都台、巴中台、江油台和湖北恩施台(图 1)。本文采用上述台站的磁通门磁力仪(仪器型号GM4，采样率1s)分钟值数据进行计算。

地磁谐波振幅比可表示为

$\left\{ \begin{array}{l} {Y_{ZHx}} = \left| {\frac{{Z\left(\omega \right)}}{{{H_x}\left(\omega \right)}}} \right|\\ {Y_{ZHy}} = \left| {\frac{{Z\left(\omega \right)}}{{{H_y}\left(\omega \right)}}} \right| \end{array} \right.$

(1)

$\left\{ \begin{array}{l} \frac{{Z\left(\omega \right)}}{{{H_x}\left(\omega \right)}} = \frac{{Z\left(\omega \right)}}{{{H_y}\left(\omega \right)}} = i\frac{\lambda }{\theta }\\ {\theta ^2} = \sigma \mu \omega i + {\lambda ^2} \end{array} \right.$

(2)

其中，Y_ZHx为NS向谐波振幅比值，Y_ZHy为EW向谐波振幅比值；Z(ω)、H_x(ω)、H_y(ω)分别为地磁垂直分量、SN向水平分量和EW向水平分量的观测谱值；μ为磁导率，σ为电导率，ω为圆频率，λ为变化磁场的波数(常数)，θ为与λ、ω、σ有关的常数，i为虚数单位。

具体计算步骤为：①将每天的分钟采样数据分成4段，计算每段5~65min周期所有频点的傅里叶谱(FFT)；②分别计算各段各频点的谐波振幅比，对其进行平均和365天滑动平均，消除噪声和年变。

2 计算结果及初步分析

根据上述计算原理及步骤，计算结果如图 2所示，其中涪陵江东台、武隆仙女山台、石柱黄水台、万州天星台、奉节荆竹台和湖北恩施台存在震前异常，巫山建坪台数据干扰较大，故在文中未做分析。距离较远的四川巴中台、江油台和成都台异常特征不明显。

对NS向谐波振幅比值(Y_ZHx)和EW向谐波振幅比值(Y_ZHy)异常特征更为详细的定量描述和统计见表 1、表 2。

3 异常特征分析 3.1 异常时空同步性变化特征

由图 2、表 1和表 2可以发现：

(1) 各台站的NS向谐波振幅比值(Y_ZHx)10min、20min和30min周期自2016年开始出现同步下降变化，2017年开始部分台站部分周期恢复上升；而40min、50min和60min周期的变化与此相反，自2016年开始出现同步上升变化，2017年开始部分台站部分周期变为下降变化。10min、20min和30min周期的变化与40min、50min和60min周期的变化为反相位，在30~40min之间似乎存在分层现象，短周期(浅部)的变化与长周期(深部)变化反相位，依据电磁学原理，认为在30~40min之间存在一个EW向的面电流。

(2) 除万州天星台外，其他5个同台站EW向谐波振幅比值(Y_ZHy)各周期自2017年开始几乎均为上升变化，万州天星台与此相反。

(3) 自2015年开始，距离震中最近的武隆仙女山台(20km)、涪陵江东台(58km)和石柱黄水台(105km)EW向谐波振幅比值(Y_ZHy)低点时间存在由长周期(深部)向短周期(浅部)迁移的现象，距离震中较远的万州天星台(160km)、湖北恩施台(180km)和奉节荆竹台(223km)未出现这种迁移现象。

另外，从目前的分析来看，谐波振幅比值幅度与地震三要素关系不十分密切，但作为地震异常分析惯例，一般需要给出异常幅度，其结果或可用于后续的研究。

3.2 异常空间分布与震源机制

从台站及震中分布(图 1)可以看到，异常台站在空间上分布较为集中，主要分布在震中以北及北东区域。距离震中最近的异常台站为武隆仙女山台，震中距约20km；最远的为奉节台，震中距约223km；震中距280km及以上的巴中、江油和成都台无异常；震中距223~280km之间无台站，说明此次5.0级地震可能引起的异常范围不超过280km。

李翠平等(2018)利用CAP方法最终反演得到的武隆5.0级地震的震源机制解(图 3)：节面Ⅰ走向24°、倾角16°、滑动角-108°；节面Ⅱ走向223°、倾角75°、滑动角-85°；P轴方位140°、仰角60°；T轴方位309°、仰角30°。结果显示此次地震主要是因NW-SE向的拉张作用而发生的。

我们计算了地磁台站的方位角(极坐标下的极角)和离源角(极坐标下的极径)，将台站与震源机制在同一坐标系叠加，结果如图 3所示，可以看出，异常台站几乎均位于震源机制解的拉张区(P波向下的区域)，其中4个台站在拉张区域内，2个台站在拉张区域边界，奉节荆竹台和湖北恩施台的投影位置重叠。

4 结论与讨论

(1) 异常台站在空间上分布集中，几乎均位于地震震中的东北侧，且从台站在震源机制解图上的投影分布看，异常台站几乎均位于震源机制解的拉张区(P波向下的区域)，有象限分布的特征。钱复业等(1982)在对1976年唐山7.8级地震前电阻率异常特征分析时发现，电阻率下降的异常台站位于P波初动向下的区域，异常台站呈象限分布；冯志生等(2013)研究认为，谐波振幅比异常变化特征与直流电阻率相似。本文研究结果也存在象限分布的特征，与前人研究较为一致。

(2) 各台站NS向地磁谐波振幅比值(Y_ZHx)变化存在分层现象，10~30min短周期(浅部)的变化与40~60min长周期(深部)变化反相位，从电磁学机理看，在30~40min周期之间存在一个EW向的面电流。

(3) 近震中台站的异常低点时间有由长周期(深部)向短周期(浅部)迁移的现象，该现象若有普遍性，则对今后判定震中有一定意义。

综上所述，本文分析发现的异常空间分布与震源机制解关系、异常在周期上的分层现象以及近震中台站异常低点时间由长周期向短周期的迁移现象，作为谐波振幅比方法的重要进展，首次利用地磁台阵分析了不同震中距台站的异常特征。上述研究是否具有普遍性，需要积累更多的震例进一步分析。

致谢: 本文在撰写过程中得到冯志生研究员、李鸿宇高级工程师、龙锋高级工程师的帮助，在此表示衷心感谢。