基于波形相似性的单台地震目录校正＊——以仙游地震序列为例

张红才廖诗荣蔡辉腾陈智勇张燕明黄玲珠巫升华陈思辉

1） 中国福州 350003 福建省地震局

2） 中国福建厦门 361000 中国地震局厦门海洋地震研究所

3） 中国福州 350001 至永建设集团有限公司

4） 中国福州 350001 福建卓越建设集团有限公司

引言

地震目录是地震台网最主要的产品之一，也是进行区域活动构造、地震预报和工程地震研究的重要基础资料.对于震群型事件，地震目录能够直观地展示出震群区域内地震能量的释放规律、区域应力状态和序列活动特点等，这对于序列趋势判定具有重要参考价值.近年来，随着我国地震监测台网建设的持续推进，台网规模不断增大，布局逐渐优化，台网的监测能力得到了逐步提升.目前，全国大部分地区地震监测能力已达到ML2.5，其中首都圈为ML1.0，华东地区为ML2.0.正在实施的“国家地震烈度速报与预警工程”项目建设完成后，我国的地震监测能力将得到进一步提升，为地震目录的获取奠定了重要的台网基础.

为了获得更为完整的区域地震事件目录，近年来也发展了多种具有较高实用性的自动化算法，如震源扫描算法（source-scanning algorithm，缩写为 SSA）（李文军，陈棋福，2006；李文军等， 2008；Janskyet al，2012；Liaoet al，2012；Zahradníket al，2015；Yuet al，2018），模板匹配滤波技术（matched filter technique，缩写为 MFT）（Peng ，Zhao，2009；Menget al， 2013；Yoonet al，2015；Zhang ，Wen， 2015；李璐等， 2017；Walteret al，2017），以及机器学习算法（Rosset al，2018；赵明等，2019a，b；Konget al，2019；Zhu，Beroza，2019； Zhuet al，2019；Mousaviet al，2020）等.其中震源扫描算法根据一系列潜在相位信号特征上的堆积能量识别事件；模板匹配滤波技术通过扫描连续波形找出模板事件附近发生的、与模板事件相关性较好的遗漏地震事件；机器学习算法则从大量的历史波形数据中提取、挖掘高维度特征信息，可实现高效的全空间搜索和函数的泛化，具备从低信噪比记录中可靠拾取震相到时的能力.

为保证产出地震目录的可靠性，地震台网分析人员通常在一定数量的台站（例如至少3个台站记录到5个以上清晰震相）观测到震相信息后才产出相关事件目录.由于地震事件的复杂性和观测台站分布的不均匀性，对于一些震级较小的事件，经常仅震中周边少数台站有触发记录.现有产出模式下，台网分析人员仍然需要对这些事件进行处理，即产出“单台”地震目录.对于震群型地震事件而言，这种情况更为普遍.受制于有限的震相信息，单台事件往往不能采用常规定位算法给出确切的震中位置，通常采用奇异值分解（singular valuedecomposition，缩写为SVD）得到震中方位角信息，然后根据单台震相到时差S－P获得震源距，最终给出参考震中位置.尽管单台地震目录中事件的震相大多数具有较高的信噪比，但由于无法提供较可靠的震中位置，因此单台地震目录在震后趋势判定等应用中一般未被考虑.若能够将这些单台事件合理地加以利用，那么便能够使之有效地服务于震后各项应急响应及科学决策工作，进一步提升地震台网的科技效益.

本文拟以福建地震台网观测报告中仙游地震序列的单台地震目录为研究对象，以多台目录震源位置为参考，考虑记录波形之间的相似性，将单台事件目录进行校正，给出较为合理的序列事件目录，以期实现单台事件目录的合理利用.

1 数据

仙游地震序列发生于福建金钟水库库区.该水库于2010年5月开始蓄水，同年8月库区开始出现地震序列活动，2012年后库区地震活动水平逐渐提高，其中2013——2014年和2018——2019年库区经历了两次较为明显的震群序列活动（图1）.截止至2020年12月，库区共记录到ML0.0以上地震4300余次，其中最大震级事件为2013年9月4日MS4.8（ML5.0）地震.目前库区地震活动依然较为活跃，主要与水库的蓄水过程密切相关.前期研究表明，仙游地震序列发生在少震、弱震区，其发震构造主要为沙县——南日岛断裂带，但无法确定具体是哪个次级断裂发震（李强等， 2015）.

图1 仙游地震序列 M-t图Fig.1 M-t plot for the Xianyou earthquake sequence

本研究收集了福建台网产出的仙游序列事件观测报告，并从数据库中截取了相应的波形记录.图2分别展示了台网编目产出的多台目录（图2a）及单台目录（图2b）的震中空间分布.由图可见，多台目录中的震中分布与库区周边的构造之间具有较好的相关性，而单台目录中的震中位置则在一定空间内呈随机状分布，基本位于以最近的XYSC台为圆心，22 km为半径的圆形区域内.两个目录中事件震级的统计结果表明（图3），多台目录中的事件以中小地震事件为主，震级范围为ML−0.6——5.0，而单台目录中的事件则以微震事件为主，事件震级集中于ML−1.0——0.4.需要说明的是，本研究所选用的多台事件均采用HypoSAT方法（Schweitzer，1997）进行绝对定位，定位过程中选用华南一维速度模型（范玉兰等，1990）（表1）.

图2 仙游地震序列编目震中分布（a） 多台地震目录；（b） 单台地震目录Fig.2 Epicentral distribution of the Xianyou earthquake sequence catolog（a） The earthquake events recorded by multiple stations；（b） The earthquake events recorded by a single station

图3 仙游地震序列编目震级情况统计直方图Fig.3 Histogram of cataloged earthquake magnitudes for the Xianyou earthquake sequence

表1 本文 HypoSAT 定位所采用一维华南速度模型（范玉兰等，1990）Table 1 The 1D Southern-China velocity model that used during HypoSAT earthquake location in this paper （Fan et al，1990）

2 方法

2.1 多台编目事件的精定位

为了使单台事件目录得到更合理的校正，本文首先采用Waldhauser和Ellsworth（2000）提出的丛集地震相对定位方法——双差地震定位法（HypoDD）将台网产出的多台目录重新进行精定位，以提供更为准确的参考震源位置（Bergen，Beroza，2018）.HypoDD方法认为如果两次地震震源之间的距离远小于震源到观测台站及介质速度变化的尺度，即可忽略不同震源到观测台站之间的射线路径差异，进而可以认为同一台站记录的两个事件之间相同震相的走时差是由震源空间位置的差异造成的.该方法利用信号的走时差反演震源位置，能够有效地消除震源至台站的共同传播路径效应，对地壳速度模型的依赖性相对较小（王未来等，2014；房立华等，2018）.

2.2 根据波形相似性创建事件丛集

由于震群序列通常集中发生于一个较小的空间区域内，以本文研究的仙游地震序列为例，多台及单台事件编目的震中均位于金钟水库库区周边（图2），且研究表明该序列在不同时间段内的发震断层参数也不尽一致（袁丽文，李强，2016），因此呈现出较明显的事件丛集特征，相同事件丛集的记录波形高度相似.本文研究采用公式

对所有的多台事件波形求取互相关系数r，式中，xi和yi分别为同一台站两个多台事件记录波形，¯ 和¯分别为两个多台事件的地震波记录期望.为方便计算，本文仅选取Pg震相到时后3 s窗长内的记录波形.根据事件记录的主要频率构成，计算相关系数前首先采用2——10 Hz无相移带通滤波器对记录波形进行处理.设定波形相关最小阈值为0.8，当互相关系数r超过阈值0.8时即认为其属于同一事件丛集.将所有事件丛集整合创建事件模板库，这样将大幅提高单台事件与多台事件匹配的计算效率，定位结果较为可靠.需特别说明的是，由于本文研究的仙游序列中单台事件大多数为XYSC台单台事件，因此本研究中事件丛集挑选、事件模板库创建及单台目录校正均基于该台站记录波形进行.

2.3 单台事件目录校正

多台事件丛集模板库创建后，采用相同的滤波频带对单台事件记录波形进行处理，然后将其与多台事件丛集模板库中的记录波形进行匹配，从中寻找相关程度最高的事件丛集，而后在该事件丛集中再次搜索最佳匹配多台事件波形，即认为单台事件发生于最佳匹配的多台事件的位置，进而实现对以多台事件目录为参考的单台事件目录的校正.本文所提出的单台目录校正方法流程示意如图4所示.

图4 本文单台目录校正方法流程示意图Fig.4 Flow chart of the single-station event catalog correction process in this paper

3 结果

3.1 多台事件精定位

应用HypoDD算法对仙游地震序列（2010年8月至2020年12月）的多台事件重定位，结果如图5所示.重定位中仅使用人工拾取震相的走时差，区域速度模型采用金震等（2018）利用天然地震及18次人工爆破观测数据得到的中国台湾海峡西侧地壳一维速度结构（表2），其中，事件对与台站之间最大距离设定为200 km，事件对搜索半径设定为50 km，事件邻居上限设定为200个，震相最小连接数目及最小震相对数目均设置为1.定位选用多台事件共4 391次，HypoDD成功重定位4 339次事件.重定位后，多台事件的水平位置及深度测定残差均有所改善，走时残差显著降低（图5），对比图2a与图6可见，重定位震源沿库区周边断裂条带状分布特征更明显，事件丛集特征也更清晰，说明采用HypoDD方法获得了较好的序列分布位置，能够为后续的目录校正提供相对准确的参考震源位置信息.

图5 应用 HypoSAT （蓝色）和 HypoDD （红色）方法的仙游地震序列多台事件重定位结果残差对比（a） 东西向；（b） 南北向；（c） 深度方向；（d） 事件走时残差均方根值Fig.5 Earthquake relocation residuals comparison of multiple-station events for the Xianyou earthquake sequence byusing the HypoSAT algorithm （blue） and HypoDD algorithm （red） respectively（a） The east-west component；（b） The north-south component；（c） The depth component；（d） The root mean square of event travel time residual

图6 应用HypoDD方法得到的仙游地震序列多台事件重定位结果Fig.6 Earthquake relocation results for multi-station events of the Xianyou earthquake sequence by using the HypoDD algorithm

表2 HypoDD 相对定位所采用的一维速度模型Table 2 The 1D velocity model that used during the HypoDD relative location

3.2 多台事件丛集模板库

以记录波形互相关系数为参考，本文将仙游地震序列的多台事件共划分为161个事件丛集，图7为其中两个事件丛集包含的多台事件波形，它们分别由12个和20个多台事件构成.由图可见，各事件丛集内的多台事件记录波形间存在高度的一致性，互相关系数均超过0.80，震相到时差Sg−Pg也几乎相同（两次事件丛集分别为1.18±0.02 s和1.21±0.03 s），且不同事件丛集的记录波形也存在明显的差异，说明本文事件丛集挑选方法具有较好的实用性.随后，对各事件丛集中对应的多台事件重定位震源位置的分析表明，除个别事件外，在空间展布上同一事件丛集中的多台事件也是聚集的，这也说明本文所划分的161个事件丛集是合理的.此外需要说明的是，每个事件丛集中所包含的多台事件数量并不尽相同，为了保证事件丛集的可靠性，本文研究中设置一个事件丛集中至少包含3个多台事件，所创建的161个事件丛集中最多一个则包括248个多台事件.

图7 两个多台事件丛集及其多台波形的垂直向记录，每个丛集波形以Pg震相到时为参考截取Fig.7 Two multiple-station events clusters and UD channel waveform in these two clusters separately，and in each cluster waveforms are aliened based on Pg arrival

3.3 单台事件目录校正

如前文所述，本文以多台事件重定位目录为参考，应用波形互相关系数在161个事件丛集模板中搜索各个单台事件最佳匹配丛集（图8），然后进一步搜索最佳匹配事件，最终实现了单台目录的校正，结果如图9所示，共对2 638次单台事件中的2 596个事件进行校正.相较于图2b，校正后的单台目录明显更为合理，其震中不再呈以单台为中心的圆形随机状分布，匹配得到的单台目录也可提供较为可靠的震源深度，与库区周边主要发震断裂之间的相关程度也更高.单台事件的参考震源位置校正完成后，对事件震级值进行校正，结果如图10所示.大多数（约80%）事件的震级校正值为0.0，约11%的单台事件震级增加0.1个震级单位，5%的事件震级减小0.2个震级单位，其余4%的单台事件也进行了不同程度的校正.尽管由于记录台站数量较少，该震级量算结果可信度较低，但仍在一定程度上说明库区周边存在着较密集的微震活动.综上可以看出，校正后单台目录的实用性显著增强，为后续业务应用提供了必要支撑.

图8 一个单台事件与多台事件记录匹配波形示例Fig.8 Waveforms of matched events for a single-station event and a multi-station event

图9 仙游地震序列事件目录校正结果Fig.9 Corrected earthquake catalog of the Xianyou earthquake sequence

图10 仙游地震序列单台事件目录震级偏差校正值统计Fig.10 Statistical of magnitude deviation correction value of single-station events for the Xianyou earthquake sequence catalog

4 讨论与结论

伴随着地震台网的不断加密和台站数量的增加，编目产出中待分析震相数量成倍增长，传统的依赖人工分析为主的地震目录产出模式将难以为继，迫切需要研发以计算机自动处理为主的新型地震目录产出模式.2018年开始，中国地震局即开展了“地震自动编目系统”的研发与示范应用，该系统由中国地震局地球物理研究所牵头，福建省地震局等单位参与.目前，地震自动编目原型系统已基本研发完成，并于2021年在首都圈、四川、云南、甘肃、福建等若干地震活动性较高的地震台网逐步推广使用.“地震自动编目系统”应用先进的人工智能算法，从实时地震观测数据中拾取震相到时信息，并采用稳健的震相关联算法和定位方法生成事件目录.离线测试结果表明，“地震自动编目系统”表现出优良的处理能力，其处理结果已基本达到甚至超过现有人工地震编目能力，相同情况下系统能够产出更完备的地震目录.更重要的是，漏报事件和误报事件发生率均维持在较低的水平.

本文针对“地震自动编目系统”产出大量单台事件目录的需求，研究相关解决方案.结果表明，应用本文方法能够有效地实现单台事件目录的校正，产出目录结果能够为震后应急响应工作提供一定参考.本文单台目录校正结果的准确性和实用性很大程度上依赖于多台事件丛集划分的科学性和合理性，这是本方法的局限性.已有研究表明（胡淑芳等，2018），仙游序列在空间上呈丛集分布，且具有明显的分区现象，事件发震时间受库区水位变化的调制作用较为明显.本文研究主要以实际观测波形之间的相似性为标准，将序列划分为161个丛集，各丛集内的事件在时间上又相互重叠.在对单台目录校正时，若能够合理增加考虑事件发生时序性特征，即在一定时间段内搜索相似波形，那么校正后目录的实用性将能够得到进一步提升.由于校正后的单台事件目录是以多台事件目录为参考，针对校正后单台目录的震源深度分析可见，仙游地震序列中的单台事件大都与浅源多台事件相关，震源深度普遍小于10 km.邱毅等（2014）采用矩张量反演方法，分析了仙游地震序列的发震构造及震区的区域应力场，结果表明该序列是由位于震区的金钟水库蓄水造成震区应力调整而引起的.因此，本文分析认为这些单台事件主要是与库区浅部的应力变化密切相关，与周边主要构造活动的相关程度较小.此外，本研究所分析单台事件以人工编目分析产出结果为基础，受到人工分辨能力限制，震源深度更深的事件可能由于信号过于微弱未能被人工拾取.应用近年来发展的一系列基于机器学习的震相拾取方法，如 PhaseNet （Zhu，Beroza，2019），GPD （Rosset al，2018），EQTransformer （Mousaviet al，2020）等，可能会从波形记录中识别出更多被遗漏的单台事件，这将有助于进一步完善地震目录.

本文基于记录波形相似性讨论了序列地震事件中单台事件目录的校正与应用，提出了一个实用的单台目录校正算法.本研究的主要目的是完善地震序列事件目录，提高单台目录的实用性.为了提供较可靠的参考位置信息，本文首先应用HyppDD算法对多台事件进行重定位，然后划分事件丛集再进行匹配搜索.张盛峰等（2014）分析认为，尽管重定位后地震目录的完备性较常规目录有一定的下降，但在地震丛集的内部下降程度较低，因此可以通过考虑地震丛集来克服.因此，本文认为以双差定位目录为参考的矫正后的单台目录在一定程度上仍可用于统计地震学参数的计算.

审稿专家对本文提出了有益的建议，福建地震台为本研究提供了台站观测记录及事件观测报告，文中图件使用Generic Mapping Tools和Matplotlib软件进行绘制，作者在此一并表示感谢.