2016年8月20日唐山震群小震震源机制

王 宁，蔡 璐，王莉婵，王亚玲，郭垚嘉，刘 丽

（1. 河北省地震局，石家庄 050021；2. 菏泽市地震局，山东 菏泽 274000）

0 引言

震源机制可以用来反映震源区域在地震发生前后的动力学运动过程，准确求解的震源机制可以用来判断震源类型，揭示地震区域应力场运动特征，并判定地震断层破裂面性质及其孕震构造。因此，震源机制的求解极为重要。震源机制反演方法很多，从最早1980年Helmberger等[1]发现Pnl波可以用于确定震源机制开始，国内外学者开展了很多研究工作，并取得了丰富的成果。关于中强震的震源机制的确定，现阶段运用较多并被接受的是改进的CAP方法[2]。而小震震源机制的确定，一直是学术界的难点。小震发生能量小，高频信号衰减快速，波形触发台站较少，其本身的局限性决定了其震源机制求解困难。即便如此，在小震求解方面的研究一直没有停止过。如利用P波初动和振幅比方法[3-5]确定小震震源机制等，虽然取得了一些成果，但现阶段并没有被广泛的认可。针对小震研究现状，本文采用了近几年被成功运用[6-8]的广义极性振幅技术（GPAT）对2016年唐山震群小震震源机制反演进行尝试，以期为唐山地区小震震源机制的研究提供一条途径。

唐山地区，地震构造复杂，活动断裂发育，是华北地震高发区，尤其是1976年唐山大地震之后，它一直是河北省重点地震监测区域。近几年来唐山地区地震发生尤其活跃，并以2016年8月唐山震群最为典型。本文以此震群ML2.0以上地震事件波形资料为基础，利用广义极性振幅技术（GPAT）反演这些事件的震源机制解，通过分析其相关特征，对判定震群的孕震机理、区域应力场相关特征等具有重要意义。

1 研究方法

广义极性振幅法（GPAT）是严川提出来的一种求解小震震源机制的非线性反演方法，是现有初动-振幅类方法的一般化[6]。此方法以震相初动极性和振幅为研究对象，并运用网格搜索技术来进行震源机制反演。

一般地，地震产生的地表位移可以表示为：

式中：j＝1，2，3，分别代表EW、NS和UD 三分向；P、S、F、O分别表示P波、S波、面波和其他震相。假设2个行矢量：

式中：i表示不同的台站，j为不同的分量，k分别代表直达P波、S波和面波；行矢量v1的元素a是带有极性的振幅绝对值最大值，行矢量v2的元素b代表P波初动极性。

按照式（3），对于给定的震源产生的波场可以写出类似的矢量：

式中：ɑm为观测矢量元素；βm为合成矢量元素；T为行矢量向列矢量的转置；w为振幅信息与P波初动信息的相对权重值。

假定震源的震源机制标量、地震矩、震源位置与实际震源的完全相同，则有=；而如果假定只有标量地震矩不同但其他参数与实际震源的相同，那么和线性相关。设定ρ为理论震源与实际震源参数的相关系数，则

另外，关于w相对权重的取值问题定义2个列向量X、Y，

两者的相关系数为：

如果w=0，则k=1，表明极性不参与计算，只有振幅参与计算；如果w→+∞，则k→-1，表明振幅对k值的影响很小。由此，根据w的取值不同，振幅和初动极性发挥不同的作用。

根据以上公式可知，βm的确定受参与计算的事件震源机制和震源位置的影响。如果震中位置确定，则受震源机制与震源深度影响，此时ρ是关于震源机制、震源深度的函数：

式中：ψ，θ，λ，d表示断层走向、倾角、滑动角、震源深度。此时震源机制反演变成非线性求解ρ1的问题。

2 资料选取与格林函数计算

2.1 资料选取

通过首都圈和“十五”测震台网的建设和改造，现河北省测震台网由1个台网中心、71个区域数字地震台站组成，台网中心汇集了省内所有台站数据，以及国家测震台网中心接收周边98个台站（北京29个、天津31个、山东14个、山西9个、内蒙古5个、河南5个、辽宁5个）的数据，共有168个台站的实时波形数据。基于台站布设分布情况，全省范围内地震监测能力达到ML2.0，有些地区已高达ML1.0。

本文使用了河北台网记录的2016年8月唐山震群资料（图1），并对每一次参与计算的地震事件进行了仔细的筛选。反演震源机制时按照震中距＜100 km、波形清晰、较高信噪比、P波震相记录应多于6个台站、排除重叠地震的条件进行筛选，选出59个（图2）符合条件的天然地震事件，其中ML2.0～2.9地震49个，ML3.0以上地震10个。

图 1 震群主震震中与台站分布图

2.2 格林函数计算

唐山震群波形资料震中距范围选择在0～100 km。此次唐山震群震级较小，波形记录的频率较高，故将震群的滤波范围设置在1～4 Hz之间，采样频率确定为滤波频率上限的5倍。格林函数库计算的方法为反射－折射率法[9]，深度范围为0～35 km，深度间隔为1 km。计算格林函数时选用的速度模型为河北台网现阶段在用的地壳速度模型—华北速度模型（表 1）。

图 2 震群震中分布图

表 1 华北速度模型

3 GPAT震源机制解结果分析

根据华北速度模型和近震台站记录，利用GPAT方法反演得到震群的最佳震源机制解（表2）。图3为地震2016—08—21 17:46:51.7实际观测数据与理论合成数据的对比结果，反演结果中观测数据初动符号与合成数据初动符号极性完全一致，观测数据振幅与合成数据振幅相关性很好，可见震源机制解的可信赖度很高。

表 2 震群震源机制解及断层类型

续表 2

图 3 2016-08-21 17:46:51.7GPAT震源机制反演结果

据震源机制解的3个应力轴倾角大小分类，Zoback将震源机制解分为正断型(NF)、正断为主兼走滑型(NS)、走滑型(SS)、逆断型(TF)、逆断为主兼走滑型(TS)、无法确定型(U)共6种[10]。据此标准，对本文所获取的59个震源机制解（表2）进行分类统计分析，不同类型的震源机制解用不同的颜色表示（图4），正断型和正断兼走滑分量类型为41次，占69.5%；走滑型地震为8次，占总数的13.6%；逆断型和逆断兼走滑分量类型为1次，占1.7%；无法确定型为9次，占15.3%。由此可见，本次震群以正断型为主，应力状态主要为剪切拉张。

GPAT震级计算采用的是矩震级MW[11-12]，台网观测报告采用近震震级ML[13]，二者之间定义不同，必然存在一些差异。图5表明：震群矩震级MW和近震震级ML有从小到大一致性越来越好的趋势，对于ML4.5以上地震是否一致性继续变好，受震群震级范围限制，没有具体数据支持，有待后续研究。另外对于ML3.6以下的地震事件，矩震级一般大于近震震级，且事件震级越小，二者差距越明显，最大差距可达到1级，这和先前研究结果基本一致[7]。

图 4 震源机制解分布图

图 5 矩震级MW与地方震级ML对比图

为了验证震源深度的正确性，把GPAT计算产出的所有事件的G深度与观测报告产出的震源深度（L深度）进行对比（图6a），并把相同事件的2种震源深度的平均值（M深度）假设为震源的正确深度，可以看出G深度和L深度有很好的相关性。为了做出更加精确的定量评价，计算了G深度与M深度的深度偏差（图6b），其中最小偏差为0 km，最大偏差为2.5 km，平均偏差为0.19 km，其中偏差小于2.0 km的占91.5%，整体偏差较小，而0.19 km的偏差可能是由于起始破裂点与“矩心”之间的差别产生的[6]，也可能是人工拾取震相误差或者选择台站不同导致的。总而言之，对于小震，GPAT深度与观测报告产出的震源深度具有相当好的一致性，它能够反映震源的真实深度。

图 6 深度对比图

4 结论

结合地震发震背景分析，震群发生于唐山断裂带中NE向的唐山-古冶断裂附近，主震为9月10日ML4.3地震。GPAT确定的震源机制节面Ⅰ走向为132°，倾向65°，滑动角−67°；节面Ⅱ走向为267°，倾向 33°，滑动角−130°；T轴方位角 205°，仰角17°；B轴方位角 302°，仰角 21°；P轴方位角 79°，仰角63°。节面Ⅰ与唐山-古冶断裂带东北段走向基本一致，判断主破裂面为节面Ⅰ，类型为正断型，兼右旋走滑活动，这与CAP和振幅比方法计算的结果总体一致[14]。另外，震群内部事件震源机制结果相似性很高，也为震源机制解的正确性提供了可靠依据。

以2016年8月河北台网产出的事件波形资料为研究对象，利用GPAT方法计算唐山震群震源机制，得出以下结论：

1）震群整体震源机制解结果相似性好，结果可信度很高；震源机制解结果以正断型为主，应力状态以剪切拉张为主。

2）矩震级一般大于近震震级，且事件震级越小，二者差距越明显，最大差距可达到1级；受震群震级范围限制，对于较大地震震级相关性暂时无法判断，有待后续研究。

3）GPAT深度与台网日常观测报告产出的震源深度具有相当好的一致性，GPAT深度能够反映震源的真实深度。

总之，GPAT方法可以作为目前台网反演小震震源机制、矩震级和震源深度的一种有效的、实用的手段。