引言

2016年9月10日18时9分在河北唐山市开平区（39.69°N，118.32°E）发生M_L4.7地震，震源深度10 km。本次地震发生在1976年唐山M7.8地震的余震区。唐山地区近百年来地震频繁，这与地质构造结构有关，其位于阴山-燕山隆起带东段，由NEE向丰台-野鸡坨断裂、宁河-昌黎断裂和NNW向滦县-乐亭断裂、蓟运河断裂相互切割成NE向的菱形块体。NNW向蓟运河断裂与NNE向沧东断裂相交，块体内部发育着地震活动性较强的NE向唐山断裂和近EW向丰台-丰南断裂。唐山断裂带由陡河断裂、唐山-巍山-长山南坡断裂和唐山-古冶断裂组成，走向NE-NNE，主要为唐山断裂附近，唐山-古冶断裂全长约30 km，南西段走向N30°E，东北段走向N50°E，断裂南端被近EW向的丰南断裂切断，如图1所示（虢顺民等，1977；国家地震局《一九七六年唐山地震》编辑组，1982；王同利等，2019）。

图 1  地震周边构造背景及观测站（点）分布

Figure 1.  Geological structure surrounding the epicenter and distribution of observation point

下载: 全尺寸图片 幻灯片

2016年8月至2016年10月，河北省唐山地区出现频繁的小震活动，地震均位于1976年唐山地震余震区内，在选定的研究范围内（39°N～41°N，117°E～1 120°E），中国地震台网共记录到621次M_L0以上地震，其中M_L0～0.9地震334次，M_L1.0～1.9地震233次，M_L2.0～2.9地震47次，M_L3.0～3.9地震6次，M_L4.0以上地震1次，最大地震为2016年9月10日发生在唐山开平的M_L4.7地震。地震活动日频次较大的时间集中在2016年8月21日至2016年9月24日，地震活动空间分布如图2所示。

图 2  研究区域2016年8月至2016年10月台网地震目录序列与空间分布

Figure 2.  Sequence and spatial distribution of network catalog from August to October 2016

下载: 全尺寸图片 幻灯片

1.   研究数据与方法

本研究范围内的37个固定台站连续观测数据来源于中国地震台网中心，台站分布如图1所示。研究方法采用基于图形处理器（GPU）的模板匹配定位识别技术（简称“GPU M＆L”），对开平M_L4.7地震前后2016年8月1日至2016年10月31日连续波形进行微震检测。为本研究提供可靠数据的定点观测台站和研究区域地质构造背景如图1所示，平均台站间距为20～50 km。

GPU M＆L技术是基于高性能GPU快速计算能力和图形识别能力，采用已知的地震事件波形作为模板事件，通过计算其与连续波形的互相关探测地震计未记录到的地震事件的技术，同时可进行地震事件检测与定位。

进行数据预处理和检测分析时，首先对连续波形进行预处理，该过程需先进行去均值和去线性趋势处理，然后采用2～8 Hz（该频段中P波和S波具有更高的信噪比）带通滤波进行预处理，同时为节约计算成本，将连续波形100 Hz的采样率降至20 Hz。经数据预处理后再进行模板事件选取，即将中国地震台网中心提供的地震目录中的地震事件从连续波形中提取出来（选取发震时刻前10 s至发震时刻后200 s）。在提取的目录波形中选取P波、S波震相清晰完整且至少有9个（信噪比＞5）台站分量同时记录到的地震事件作为模板事件，同时噪声窗口定义为P波到时前7 s至前1 s。根据以上挑选原则，从中国地震台网中心提供的原始目录中最终确定283个具有高信噪比、清晰震相的地震事件作为模板事件（王同利等，2019）。第三步通过Taup软件（Crotwell等，1999）为每个模板地震的所有台站记录的垂直分量和水平分量波形标记P波和S波到时，并通过人工方式读取波形进行拾取和校正。然后设定P波和S波到时的前1 s至后5 s为模板窗口长度，用于与连续波形计算滑动互相关系数。最后以模板事件所在的位置作为中心，设置固定尺寸网格进行搜索，假设搜索到的每个网格点均为检测到的新地震震源位置，则计算每个网格点与模板事件到每个台站的到时差，先将各台站每个分量的滑动互相关序列的发震时间减去其模板S波的走时，再减去检测与模板事件到每个台站的到时差进行发震位置校正，最后对所有台站分量的滑动互相关系数计算平均值，即为该检测事件位置的滑动互相关系数。每个检测事件重复以上操作可最终获得叠加的滑动互相关系数（王同利等，2019）。

对每个台站各分量的信噪比进行计算，当互相关系数及检测的地震事件波形信噪比均大于探测阈值时，则该位置被确定为检测地震的震源位置（Zhang等，2015）。检测地震的震级是其与对应的模板事件在各台站的S波水平分量到时前2 s至后2 s内最大振幅比的平均值。如果1个检测时间窗口中出现1个模板检测到多个事件或多个模板检测到1个事件，挑选相关系数最高的事件为最终检测事件。本研究设定地震检测阈值为9倍中位差，检测事件的时间窗口为2 s，如果2 s内检测到了多个事件，则相关系数最高的事件被确定（王同利等，2019；Liu等，2020）。

利用双差定位方法（HypoDD）对检测事件进行相对定位校正，从而获得精度较高的地震事件相对位置。王同利等（2019）采用该方法对2015年阿拉善M_S5.8地震进行了检测，廖诗荣等（2021）采用该方法对2021年云南双柏M_S5.1地震进行了检测，张雅楠等（2023）采用该方法对上海及邻区10年微震进行了检测，有效提升了地震目录的完备性。本研究将该方法运用在中级地震前后微震检测中，经GPU M&L技术探测和定位后，探测到2016年唐山开平M_L4.7地震共847个事件，完备性从M_L=0.2提升到了M_L=−0.1，如图3所示。

图 3  检测地震目录数量-震级分析

Figure 3.  The magnitude-number of catalogs after detection

下载: 全尺寸图片 幻灯片

2.   结果分析

震群的研究对认识区域地震的孕育和发展过程具有重要意义，尤其是对触发震群活动的物理过程及持续活动的动力学研究。Vidale等（2006）通过分析美国南加州地区的震群活动，提出级联触发模型是触发震群活动最常见的物理模型之一，其发震物理过程与主震触发余震的机制相同，主要建立在地震目录完备性和统计特征上（Helmstetter等，2003a，2003b；谭毅培等，2016）。因此，通过对唐山地区中强地震前后的地震事件进行检测和定位，以对其完备性和地震活动性进行探讨，为唐山地区震群深入研究和发震机理分析等奠定基础。

2.1   检测地震事件检验

在本次唐山开平地震震中周边150 km范围内共选取283个地震事件作为模板事件，通过GPU M&L技术检测微小地震，得到878个地震事件，检测到的地震事件个数约为模板事件的3倍，约为台网目录的1.4倍。其中，检测到主震前地震事件198个，检测到主震后地震事件680个，共847个有着清晰P波震相的检测事件进行双差重定位后为检测地震目录。采用最小二乘法对台网目录和检测地震目录进行b值拟合分析，结果如图4所示。由图4可知，检测地震事件的b值拟合效果明显优于台网目录，尤其是微震部分，且检测到的遗漏事件震级多在M_L2.0以下。

图 4  台网目录和检测地震目录b值拟合分析

Figure 4.  Statistical catalogs analysis of b value between the network and the detection

下载: 全尺寸图片 幻灯片

已有研究表明b值是检验地震目录个数和震级之间线性相关的最佳方法之一，Gutenberg等（1954）研究认为地震目录的震级与频次较好地服从G-R关系，表明地震频度对数与震级之间存在线性关系，其关系曲线斜率为b值。b值能够反映震源区域的应力状态，即某区域b值基本稳定则反映该区域处于正常地震活动状态，如果b值降低或升高则反映该区域的介质构造应力较强，国内外诸多震例均证实了该现象（马鸿庆，1978；黄德瑜等，1981；Felzer等，2004）。任雪梅（2011）通过运用不同完整性的地震目录分别计算了b值。黄玮琼等（1989）对b值法的定量分析表明b值法是检测地震目录完备性最好的方法之一。台网目录和检测地震事件相同震级的数量对比如图4所示，统计结果显示经过GPU M&L技术检测后的地震目录相比台网目录b值线性更好，这反映出地震目录中小震遗漏较多是造成b值线性的原因之一。

对于区域地震台网监测能力的评估标准，国际上普遍采用完备震级M_C标准，M_C表示在一个时空范围内，该台网100%能够监测到的最小地震震级（Rydelek等，1989）。但中强地震发生前后的某些小震难以用常规方法识别，这是由于该台网的信噪比低或该事件震相被其他震相干扰，造成地震目录不完整，从而影响地震活动性特性研究的准确性。台网目录和重定位后847个有着清晰P波震相的检测地震震级和数据变化显示，台网目录的完备震级为M_L=0.2，检测地震目录的完备震级提升为M_L=−0.1（图3），有效提升了地震目录在小震级范围的完备震级。

2.2   检测地震事件空间分析

通过采用HypoDD方法对GPU检测结果进行重定位（杨智娴等，2003；Hauksson等，2005），如图5所示，图中红圈为台网目录，蓝圈为检测目录，圆圈大小随震级的增大而增大。在定位过程中，考虑P波震相走时精度高于S波震相走时，因此赋予P波震相权重为1.0，S波震相权重为0.6，还设定事件对与台站之间的最大距离为300 km、事件可以组成的最大事件对为10个、事件对之间的间距＜40 km、地震对的最大距离为25 km等。对本次唐山开平地震前后847个有着清晰P波震相的检测地震事件进行了重定位，重定位后的地震目录相对台网目录如图5（a）所示，在主震周边明显增强，地震事件的长轴方向呈NEE向，与唐山-古冶断裂走向基本一致。主震发生前后唐山-古冶断裂的地震活动性明显增强，尤其是主震前，断裂附近震群发育，主震后的地震活动性逐渐减弱。重定位后的地震事件随震源深度的变化如图5（b）所示，重定位后的地震事件弥补了震源深度为5～20 km的微小地震，使地震事件总体分布较均匀。分别对沿唐山断层走向（图5（a）中AA′）和垂直于断层走向（图5（a）中BB′）进行了剖面地震事件分析，结果如图5（c）和图5（d）所示，大多数的检测地震事件震源较深，且多在断层附近分布。台网目录和检测地震目录的空间分布同样显示出了上述特征（图6），检测事件主要分布在震源深度为5～15 km的部位，且基本集中在主震附近。因此，认为本次唐山开平地震的发震断裂为唐山-古冶断裂的可能性较大。

图 5  重定位后地震事件分布

Figure 5.  Distribution of detection earthquake

下载: 全尺寸图片 幻灯片

图 6  检测地震事件震源空间分布

Figure 6.  Source spatial distribution of detection seismic

下载: 全尺寸图片 幻灯片

对水平方向和震源深度方向重定位误差进行了计算。重定位后所有事件均方差（RMS）平均值为0.157 s，平均定位误差在EW、NS、UD向分别为0.045、0.054、0.056 km。已有研究表明重定位方法在深度方向的误差大于水平方向的误差（房立华等，2013；王未来等，2014；梁姗姗等，2017；朱艾斓等，2005），这不是方法定位的误差，而是由于观测台站分布在地表、在震中附近台站较稀疏、观测精度低等原因导致的误差。本次唐山开平地震重定位后主震前40 d和主震后50 d地震序列分析结果如图7所示。由图7可知，重定位后的地震目录有效补充了本次地震序列，其中主震后24 h检测到的地震事件最多，约为台网目录的3倍，主震后45 d检测到的微小地震是台网目录的2.5倍左右，且检测到的余震事件多发生在主震发震断裂附近。GPU M&L技术检测后的地震事件分布更均匀，尤其是补充了主震后的地震事件。从断层走向来看，地震事件主要集中在主震周边15 km范围，垂直于断层走向的地震事件主要分布在主震周边10 km范围，这可能与发震断裂的构造结构有关。

图 7  检测地震事件时序分布

Figure 7.  Time serials distribution of detection seismic

下载: 全尺寸图片 幻灯片

3.   结论

地震目录的最小完备震级是地震台网监测能力的重要参数之一，通常将其定义为检测能力覆盖范围内的所有地震事件中检测到的最小震级，因此通过GPU M＆L技术检测的地震目录增强了地震台网的检测能力。利用GPU M＆L技术检测时，在台网目录中选取波形清晰的地震事件作为模板事件，计算模板事件波形和连续地震信号之间波形的相关性，根据信噪比大小给不同台站分量赋予不同的权重系数进行互相关滑动检测缺失的地震事件。该技术显著提高了微震事件探测能力，同时利用GPU内核进行加速，提高了微震事件检测速度。对于与模板事件间距越小且震级越大的检测事件，对应的互相关系数越大，即事件波形相似度越高。

本研究开展了中级地震的微震检测，首先从研究区域内地震台网地震事件目录中挑选了283个地震事件作为模板事件，然后对唐山开平M_L4.7地震周边37个台站2016年8月1日至2016年10月31日地震波形进行扫描，共识别出878个地震事件。检测到的地震事件数量约为模板地震目录的3倍，约为台网目录的1.4倍，使地震目录的完备震级由原来的M_L=0.2提升为M_L=−0.1，取得了显著效果。

GPU M&L技术进行微震检测后，利用HypoDD方法对847个地震事件进行了重定位。重定位后的结果显示，本次唐山开平地震主震发生前后区域地震活动性较强，地震事件空间分布较集中，主要集中在震中附近，且以小震为主。地震事件主要沿唐山—古冶断裂呈线性分布，分析认为本次唐山开平地震发震及小震活跃的断层为唐山-古冶断层。