东祁连山北缘断裂带基于深度学习的密集台阵地震事件快速检测与定位研究1
2022-06-01杨少博王炳文张海江
杨少博 王炳文 高 级 张海江
1)中国科学技术大学地球和空间科学学院, 合肥 230026
2)安徽蒙城地球物理国家野外科学观测研究站, 合肥 230026
引言
地震事件检测与定位是很多地震学研究的基础。对于活动断层区域,可利用完备、可靠的地震事件检测与定位结果更好地了解区域地震活动性、刻画断层几何形态(Ross 等,2019a,2020)。近年来,密集台阵的大规模布设有利于更多的微地震活动监测(Meng 等,2018;Li 等,2018),使需处理的地震数量迅速增加。人工智能在地震数据处理方面的成功应用提供了一系列快速、可靠的自动化地震数据处理方法,包括地震事件检测(Perol 等,2018;Ross 等,2018;Yang 等,2021)、到时拾取(Zhu 等,2019;Mousavi 等,2020;Xiao 等,2021)、震相关联(McBrearty 等,2019;Ross 等,2019b)等,这些自动化方法能够达到甚至超越人工处理的精度,大幅度提高了数据处理效率。
Wang 等(2020)在Raton Basin(位于美国新墨西哥州和科罗拉多州的边界)布设了密集台阵,研究该区域诱发地震活动性。该研究将连续数据划分为多个12 s 的数据段,使用PhaseNet(Zhu 等,2019)拾取每个数据段中的P、S 波到时,使用REAL 算法(Zhang 等,2019)进行震相关联并初定位,使用VELEST 算法(Kissling 等,1994)优化定位结果,使用相对定位算法hypoDD(Waldhauser 等,2000)进行重定位。相比该地区之前的地震目录,该流程定位了更多的地震,更好地刻画了断层的几何形态。Liu 等(2020)基于PhaseNet 到时构建了2019 年加州Ridgecrest 地震序列目录。赵明等(2021b)对PhaseNet 做迁移学习,进而用于构建长宁地震的前震目录。Ross 等(2018)使用基于深度神经网络的Generalized Phase Detection 算法进行地震事件检测及拾取,使用PhaseLink(Ross 等,2019b)进行震相关联,使用传统定位方法获得该区域的地震目录。上述研究结果均表明将人工智能地震数据处理方法与传统定位方法相结合,能够快速、可靠地获得某地区的地震目录。
将上述地震事件检测与定位流程应用于东祁连山断裂带的密集台阵数据处理。由于密集台阵数据量大,在缺少GPU 加速的情况下,使用PhaseNet 拾取连续数据需耗费大量时间。因此,首先使用一个轻量的神经网络CNNDetector(Yang 等,2021)进行地震事件检测,然后使用PhaseNet 对检测到的事件进行到时拾取,最后使用传统方法进行地震定位和震级计算。相比中国地震台网正式地震目录,获得了更多的地震数据。
1 研究区域
祁连山位于青藏高原块体、塔里木块体和阿拉善块体的衔接带上,其构造形变受青藏高原挤压作用及周缘活动断裂的控制。始新世末-渐新世早期,伴随着印度板块和欧亚板块的碰撞,祁连山地区发生构造变形与隆升;渐新世晚期-中新世早期,祁连山地区构造稳定,地形起伏变小,形成大面积夷平面;中新世中晚期(10~17 Ma)约为断裂带起始活动时间,断裂带以 NNE-NE 向的挤压逆冲和地壳缩短增厚为变形特征,造成断裂带沿线山脉的快速隆升;中更新世以来,断裂带转变为左旋走滑为主兼具逆冲运动,调整着青藏高原东北缘地壳物质向东运动,最终被断裂带尾端的六盘山逆冲褶皱带吸收(张培震等,2004;戚帮申,2014;郭鹏,2019)。
晚第四纪以来,在印度板块向欧亚板块NNE 向持续的构造推挤过程中,青藏高原东北缘受北侧阿拉善地块和东侧鄂尔多斯块体的阻挡,在NE 向或NEE 向挤压作用下,构造变形十分强烈。祁连山吸收了青藏高原NE 向的挤压变形,张培震等(2004)计算得到青藏高原北部边界的阿尔金山、祁连山共吸收了青藏15%~17% 的地壳缩短;祁连-海原断裂带作为青藏高原东北缘主要的边界左旋走滑断层,在调节高原东北缘地壳物质相对于阿拉善地块的向东运动中起重要作用。
本研究区位于北祁连山东部,主要包含皇城-双塔断裂、冷龙岭断裂。冷龙岭断裂与海原断裂、老虎山断裂、毛毛山断裂、金强河断裂和托莱山断裂构成长约 1 000 km 的祁连-海原断裂带。冷龙岭断裂早期活动以挤压逆冲为主,自中更新世以来,断裂活动性质以左旋走滑为主,兼有正断层分量,晚第四纪时期主要表现为左旋走滑运动,现今滑动速率约为4~6 mm/yr (郭鹏等,2017;郭鹏,2019)。
从几何分段上可将皇城-双塔断裂带划分为西、中、东段,其中西段为皇城盆地断裂段,在晚更新世以来活动性弱,中、东段分别为上寺断裂段和冬青顶断裂段,活动性均较强。上寺断裂段以叠瓦状逆冲断构造为特征,冬青顶断裂段分为南、北枝,南枝为正断层,北枝为逆断层(侯康明,1998),1927 年古浪地震发生于东段。
2 研究数据
为确定研究区域精细的断裂结构,并监测祁连山断裂带附近的微地震活动性,于2019-07-20-2019-08-21 在祁连山断裂带布设了面状密集台阵。该台阵包含6 条测线,每条测线由40 个短周期仪器EPS-2-M6Q组成(图1),仪器频带范围0.2~150 Hz,采样率200 Hz。测线之间的距离约10 km,测线内台间距从测线两端至中间逐渐由数公里缩小至数百米,整个台阵覆盖范围约为3 000 km2(50 km×60 km)。对采集数据进行预处理,包括去均值、去线性化及重采样至100 Hz。图1 中蓝色三角表示台站,四个白色三角形为编号1 001、1 002、1 030、2016 的台站,黑色圆点表示历史地震事件(莘海亮,2020),黑色虚线表示断裂带,HSF 表示皇城-双塔断裂,LLLF 表示冷龙岭断裂,两个红色五角星为编号317 和编号623 的地震事件。
图1 研究区域台站及历史地震事件分布Fig. 1 Distribution of stations and historical earthquakes in the study area
3 数据处理流程和结果
参考Wang 等(2020)的数据处理流程,若直接使用PhaseNet 进行事件检测与到时拾取,需将24 h 连续数据切割成5 760 个15 s 的数据片段,在没有GPU 的服务器上处理240 个台站24 h 的数据需耗费约36 h。相比PhaseNet,CNNDetector 仅用于检测地震事件,其网络结构较简单,检测效率更高。在同台服务器中,处理240 个台站24 h 的数据仅需耗费约4 h,因此首先使用CNNDetector 检测地震事件,然后使用PhaseNet 拾取到时,再使用REAL 初定位,最后使用hypoDD 进行重定位。
3.1 地震事件检测
Yang 等(2021)提出基于卷积神经网络的地震事件检测算法(CNNDetector),此算法综合多个台站的三分量地震记录判断1 个时窗内是否存在地震事件,其检测精度高、误检率小,并在多个区域测试中表现出较好的泛化能力。进行祁连山密集台阵数据处理时,为检测到更多的微小地震,将240 个台站按区域分为12 组,对数据进行3~30 Hz 的带通滤波后,使用CNNDetector 对每组台站分别进行检测。设置阈值时,首先需满足能够检测到大多数已知地震事件的要求,然后尝试降低阈值,直至出现虚假检测结果,由此将阈值设置为0.5。在同一时窗内,若有2 组以上的CNNDetector 检测值超过阈值,则认为其为地震事件。在此标准下,共检测到1 036 个地震事件,其中2 个低阈值的地震事件如图2 所示,其地震位置由图1 红色五角星表示。
图2 CNNDetector 检测出的2 个地震事件波形Fig. 2 Waveforms of two seismic events detected by CNNDetector
3.2 到时拾取
由于台站数量多,手动拾取检测到的1 036 个地震事件需耗费大量时间,因此使用自动到时拾取算法进行拾取。PhaseNet 是目前应用较多的基于深度学习的到时自动拾取算法(Zhu 等,2019),其训练数据是北加州地震数据中心超过30 年记录到的约700 000 个地震波形,其对低信噪比的数据也具有一定拾取精度,且在不同地区具有较好的泛化能力(赵明等,2021a)。将P 波和S 波拾取阈值均设为0.3,共拾取到74 585个P 波到时和75 191 个S 波到时。地震事件317 的部分到时拾取结果如图3 所示(图1 中白色三角形为对应的接收台站位置),由图3 可知,部分信噪比较低的信号也能被准确拾取。
图3 PhaseNet 到时拾取结果和不同台站接收到的事件317 的波形图Fig. 3 Examples of PhaseNet arrival time picking results and waveforms of event 317 received at different stations
3.3 震相关联及初定位
手动拾取到时时会考虑台站间信号的关联性,从而判断待拾取的是地震信号还是干扰信号。但自动拾取不会考虑,因此在震相关联前的时距曲线中有很多偏离P、S 波时距曲线的点。因此对于自动拾取的震相到时,做进一步的震相关联是十分必要的。这里我们使用基于网格搜索的REAL 算法(Zhang 等,2019)进行震相关联及初定位。
选取台阵中心点处的USTClitho 2.0(Han 等,2022)速度模型作为震相关联及初定位使用的速度模型(图4)。REAL 的搜索范围为0.6°×0.6°×35 km,搜索网格尺寸为0.05°×0.05°×2 km,最终保留了550 个满足以下条件的地震事件(图5(a)):超过10 个P 波到时、5 个S 波到时、20 个P 波加S 波到时且台站间隙角<240°,其中震相走时残差需<0.5 s。震相关联后剩余43 152 个P 波到时及42 414 个S 波到时。震相关联后的时距曲线如图5(b)所示,异常点已剔除,其中灰色圆点为关联前的震相散点,红点为关联后剩余的震相散点。初定位后的残差分布如图5(c)所示。
图4 地震定位使用的一维速度模型Fig. 4 1D velocity models used for earthquake location
图5 REAL 震相关联及初定位结果Fig. 5 Seismic phase association and earthquake location results by the REAL algorithm
3.4 重定位
为进一步提高地震定位精度,使用相对定位算法hypoDD(Waldhauser 等,2000)进行重定位。hypo-DD 算法使用地震对的相对到时数据进行地震定位,由于震中距远大于地震对的距离,2 个地震近似拥有共同路径,因此可减小速度模型不准确带来的相对定位误差,且可更好地约束地震的相对位置。重定位时使用USTClitho2.0 速度模型,选取的地震对距离<10 km,相互关联的地震数目超过8 个。重定位后剩余517 个地震事件,地震分布如图6 中黑色圆点所示,其中灰色圆点为2013-2016 年重定位的历史地震分布(莘海亮,2020),灰色虚线为断裂带。重定位的迭代过程及残差分布如图7 所示。
图6 HypoDD 定位结果Fig. 6 HypoDD location results
图7 HypoDD 迭代过程及走时残差分布Fig. 7 HypoDD iterative process and residual distribution
3.5 震级计算
参考Wang 等(2020)研究方法计算震级:
式中,ML为 震级,A为振幅,d为震中距。
对事件波形数据去仪器响应、滤波至1~50 Hz 后,对于每个地震,在所有存在到时的台站上计算震级,取所有台站计算震级的中位数为最终的震级,如图8(a)所示,2 个地震目录中共同地震事件的震级对比如图8(b)所示,由此可知,2 种震级差异较小,说明式(1)计算的震级是可靠的。
图8 震级分布与对比Fig. 8 Magnitude distribution and comparison
4 讨论与结论
在祁连山密集台阵观测期间,中国地震台网正式地震目录中共有39 个位于台阵内的地震事件,而密集台阵共监测到550 个地震事件,约是地震台网地震数量的15 倍,尤其监测到了更多小于0 级的地震事件。通过对比发现,中国地震台网正式地震目录中的39 个地震事件全部被密集台阵监测。因此通过密集台阵监测,可大大提高微地震监测能力,可监测到-2 级的微震事件。
由于数据量大,为减少人工拾取工作量,提高数据处理效率,使用最新发展的基于人工智能的地震数据处理方法。在缺少GPU 加速的情况下,首先使用CNNDetector 进行地震事件检测,然后使用PhaseNet 拾取P 波和S 波到时,大大节约了计算时间。在人工智能拾取到时的基础上,结合传统定位方法REAL 及hypo-DD,精确地定位出517 个地震,这表明基于深度学习进行地震事件检测与到时拾取,并在此基础上进行地震定位是可靠、高效的。密集台阵重定位的地震位置与研究区域历史地震活动具有较好的一致性,均由西至东呈现出4 个条带状分布。最西侧的地震簇和冷龙岭断裂相关,走向与断裂走向类似,呈近NW-SE 向。与地表断层迹线相比,最西侧的地震簇位于断层东北侧,垂向剖面显示该断层向EN 向倾斜。与最西侧的地震簇相比,东侧的3 个地震条带走向和皇城-双塔断裂、冷龙岭断裂走向不一致,而是呈NNW-NNE 向,与最东侧的断层走向一致。地震定位结果表明,研究区域除地表显现的NW-SE/NEE-SWW 向的断层外,还存在活跃的NNW-NNE 向隐伏活跃断层。这些隐伏活跃断层与地表显现的断层呈一定斜交关系,且更活跃。与中国地震台网地震目录给出的地震定位相比,基于密集台阵给出的地震定位表现出了更明显的线性结构,体现了地震定位的高精度,可刻画出断层面的三维几何形态。
致谢 感谢国家重点研发计划(2018YFC1504102)的资助。感谢中国地震台网中心国家地震科学数据中心(http://data.earthquake.cn)提供的数据支撑,感谢莘海亮博士提供的历史地震目录。感谢审稿人提供的宝贵建议。