2022年门源MS6.9地震震区三维速度与发震机制研究

尹欣欣, 邱江涛, 李敏娟, 许康生, 张 波, 张胜霞

(1. 甘肃兰州地球物理国家野外科学观测研究站, 甘肃 兰州 730000;2. 甘肃省地震局, 甘肃 兰州 730000;3. 中国地震局第二监测中心, 陕西 西安 710054)

0 引言

据中国地震台网中心正式测定,2022年1月8日1时13分(北京时间)在青海海北州门源县(37.77°N,101.26°E)发生了MS6.9地震,震源深度10 km。自主震发生后截止到1月10日,此次地震共造成了青海海北州门源、祁连、刚察县共1 662户5 831人受灾,受伤9人,损坏房屋4 052间。据现场调查队报告,门源MS6.9地震形成了长约22 km长度的地表破裂带。虽然震中人口密度低、人员稀疏,但位于震中村落房屋受损较为严重:甘肃省永昌县距震中80 km,有149户受灾,受灾人数为502人,一般损坏房屋为296间,个别牲畜棚圈倒塌;在永昌县,地震造成的直接损失达232.3万元。地震造成多条铁路停运,其中,兰新高铁浩门至军马场路段,隧道和线路设施受到严重的损坏,导致兰新高铁长期中断运行。根据现场烈度评估报告结果,地震最高烈度为Ⅸ度(9度)(https://t.ynet.cn/baijia/32080013.html)。此次门源地震序列较为复杂,呈分段性质,西段沿托莱山断裂展布,而东段却与冷龙岭断裂很好的重合,主震位置位于托莱山断裂和冷龙岭断裂相接之处,属于应力集中点[1]。

强震的发展和发生涉及地球内部的动力学或构造运动,其产状与地壳和上地幔的深部结构、物理性质和动力环境密切相关。地壳的速度结构是反映地壳介质性质的重要参数,与构造运动密切相关的三维速度结构的特征提供了有关地震位置和震源介质的重要信息。结合地震活动的空间分布,速度结构是理解地震发生环境和机制的重要依据。此次门源地震的走滑性质与2016年1月21日的门源MS6.4地震的逆冲性质截然不同,序列拓展体现了一定的复杂性,为后续地震危险性分析造成了很大困难。部分学者认为2016年门源地震的发震断裂为冷龙岭断裂在青藏高原北东向拓展过程中产生的伴生断裂,表现出逆冲特征,是青藏高原北东向拓展过程中的一次地震事件[2-3]。为深入探讨门源地区的深部孕震环境,本文使用中国地震台网中心给出的地震观测报告,选取(36°～39°N,101°～104°E)范围的地震事件作为基础数据,利用双差成像方法[4]进行了精定位与层析成像处理,对2022年1月8日门源MS6.9地震进行案例研究。

1 区域背景构造简介

2022年1月8日门源MS6.9地震震中位置处于青藏高原一级地块和阿拉善一级地块过渡区,过渡区主要发育祁连山—河西走廊活动断裂系(图1)。活动断裂由南至北大体分为四个区块:一为高原内部断裂系,如走向NNW的日月山断裂;二为祁连山北缘断裂系,包括大型海原断裂(托莱山断裂段、冷龙岭断裂段、金强河断裂段)及其南北侧扩展断裂系(如门源盆地北缘断裂、民乐—大马营断裂、皇城—双塔断裂等);三为河西走廊内部断裂系,如民乐—永昌断裂、武威盆地南缘断裂、榆木山东缘断裂等;四为走廊北侧的阿拉善南缘断裂系,如龙首山南缘断裂、龙首山北缘断裂、河西堡断裂等。区域破坏性历史地震和现代中强震密集,主要发生在冷龙岭断裂沿线及其北扩断裂系、河西走廊内部断裂系和阿拉善南缘断裂系(图1)。区域内强震主要包括1927年古浪8级地震、1954年山丹7级地震、1986年门源6.5级地震、2003年民乐—山丹6.1级、2016年门源6.4级地震等。现代小震表现为丛集分布,在冷龙岭、龙首山东段尤为密集;现代强震也伴生较小的地震密集区,如2003年民乐—山丹6.1级地震、2016年门源6.4级地震等[5]。

F1: 日月山断裂;F2: 达板山断裂;F3: 木里—江仓断裂;F4: 大通山北缘断裂;F5: 门源盆地北缘断裂;F6: 马牙雪山断裂;F7: 金强河断裂;F8: 冷龙岭断裂;F9: 托莱山断裂;F10: 肃南—祁连断裂;F11: 祁连山北缘断裂;F12: 榆木山东缘断裂;F13: 民乐—大马营断裂;F14: 皇城—双塔断裂;F15: 民乐—永昌断裂;F16: 武威盆地南缘断裂;F17: 龙首山南缘断裂;F18: 河西堡断裂;F19: 龙首山北缘断裂;F20: 天桥沟—黄羊川断裂图1 门源地震的区域构造背景Fig.1 Tectonic setting of Menyuan earthquake

震中附近的活动断裂主要为海原断裂冷龙岭段(F8)、海原断裂托莱山段(F9)、肃南—祁连断裂(F10)、门源盆地北缘断裂(F5)、民乐—大马营断裂(F13)、皇城—双塔断裂(F14)等。门源地震的震中位于海原断裂冷龙岭段、托莱山段、肃南—祁连断裂的交点附近,海原断裂冷龙岭段、托莱山段均为全新世活动断裂,运动性质以左旋走滑为主,两段之间呈左阶拐弯;肃南—祁连断裂为晚更新世活动断裂,运动性质以逆断为主,断裂倾向南西;门源断裂可能为全新世活动断裂,活动性质兼具左旋走滑和逆冲分量;民乐—大马营断裂为全新世活动断裂,逆断分量显著[6]。

2 数据介绍

本研究收集了青海、甘肃、四川、陕西等省区域数字地震台网2009年1月至2022年2月的地震观测资料,这些震相报告所使用的台站主要为各省固定数字测震地震台站,另外有少许流动观测台站,台站仪器以宽频带为主,型号包括BBVS-60、CMG-3ESPC、KS-2000等,采样率为100 Hz。门源地震区及其周边地震台站分布密集,覆盖研究区方位角范围大。这确保了地震序列目录的完整性,以及射线交叉点的密度和成像数据反演的可靠性。截止到2022年2月8日,中国地震台网中心人工编目余震共计3 709次,单台事件2 690次,最少台站记录4个以上地震事件727次,其中有两次5级以上余震事件,分别是1月8日2时9分MS5.1以及1月12日18时20分MS5.2地震,MS≥3.0地震一共13次。为了保证足够的地震事件数,本文从中国地震台网中心下载了门源地震及周边(36°～39°N,100°～103°E)范围内2009年1月1日—2022年2月8日之间的地震台网观测报告。为了保证地震震中位置的可靠性,选取最小台站记录数为4,按照如上挑选条件共下载地震事件14 639次,震级范围M0.3～M6.9,原始震源深度分布范围0～30 km。另外由于研究区属于青海、甘肃交界处,按照中国地震台网中心工作规定,两省地震台网都需要进行地震编目分析,因此直接导出的地震目录中掺杂着重复地震目录,为确保数据唯一性,本文首先去除了不同台网间记录的边界上重复地震事件,在去重复地震事件的处理上以地震位置所在省份为准,只保留各台网给出的自己省份地震事件,地震震中分布以及时间分布如图2所示。

图2 地震台网目录震中时间空间分布Fig.2 Temporal and spatial distribution of epicenters in catalogue of seismic network

在震相数据制作过程中,由于用于反演的地震震相准确性直接关系到反演结果的好坏,本文为严格挑选地震震相,设定了震中距500 km为距离阈值,且走时不超过绿色实线限制内的震相,最终得到P波震相为97 504个,S波震相86 469个,具体分布如图3所示。

图3 本文所使用P、S震相数据Fig.3 P- and S- seismic phase data used in this paper

3 研究区地震重定位与三维速度反演

3.1 双差成像方法介绍

本研究获取地震相对定位结果所采用的TomoDD方法目前已在各个研究区域得到广泛应用并获得了众多研究成果[4,7]。该方法由于结合了双差地震定位(hypoDD)[8]和地震层析成像[4]两种理论,通过三维速度结构与地震定位的联合反演迭代,相比传统的绝对地震定位法可得到更精准的震源位置[7]。与其他层析反演方法不同,TomoDD使用P波和S波的绝对走时和相对走时,确保了震源参数确定的高精度。TomoDD的基本思想是一对事件之间的震源距离小于事件站距离。考虑到给定站点的每对事件,TomoDD最小化观察到的和预测的走时之间的残差,迭代求解震源之间的相对位置向量。TomoDD算法使用伪弯曲射线追踪理论来计算事件和站点之间的射线路径,反演过程使用可变加权过程。第一次反演是通过在绝对走时上放置一个大权重(1.0)和在相对走时上放置一个小权重(0.1)来执行的。获得的模型再用作第二次反演的输入,其中相对走时此时具有较大的权重(1.0),而绝对走时受到时间惩罚权重降低(0.1)。生成的模型将作为反演最后一步的输入。绝对走时和微分时间(CTDT)结合并同时使用TomoDD在迭代最小二乘程序中反演,该程序利用LSQR方法。

(1)

(2)

根据式(2)可以得到地震震中附近小尺度范围的速度结构以及相对震源位置。对地震事件的配对是使用ph2dt方法来确定的,该方法主要最用近邻方法来确定地震对[7]。在执行地震配对过程中,本研究设定震中间距50 km作为事件对的遴选阈值。此外,设定单一地震事件最多可与10个地震组成相应的地震对。对地震定位与成像研究区(100.00°～103.00°E,36.00°～39.00°N)的14 869个地震事件,经过配对后共得到14 665次地震事件,以及279 875对P波震相资料、265 983对S波震相数据。完成地震配对后本文再利用TomoDD方法进行层析成像以及地震重定位处理,最终获得14 024个地震事件的重定位结果。

3.2 速度模型以及参数选取

地震定位以及速度模型的反演对初始速度模型的准确度要求较高,有良好的初始速度模型,可以令反演更快地收敛并更容易找到准确结果。门源地区地震频繁,前人研究成果较多可供参考,因此在速度模型选取中,本文参考前人结果[7,9]并多次试错后,根据棋盘格测试速度模型恢复程度来作为判断依据,最终选取一维速度模型如表1所列。

表1 本文所用速度模型结构

在反演过程中,对于数据量较小的计算通常采用阻尼 LSQR 算法(最小二乘 QR 因子分解)来解决带有阻尼的最小二乘问题。总走时残差的L2范数被用作迭代和求解方程的目标函数。在成像过程中,慢度的变化受平滑因子的约束,而地震位置和慢度的变化受阻尼因子的约束。由于反演结果的稳定性很大程度上取决于平滑和阻尼因子的大小,因此在反演前检查阻尼和平滑因子的选择尤为重要。与速度模型一样,本文在参考前人研究结果后,经过多次尝试,最终确定阻尼和平滑因子分别为400和15。

门源地区的三维速度结构,已经有不少前人研究结果,但分辨率较低,本文为提高分辨率,得到更精确的三维速度结果,将主要反演区域设定为(37.5°～38.0°N,101.1°～101.8°E,以下简称研究区),节点设置为0.1°,即每隔0.1°设立一个反演网格。

3.3 地震重定位结果

地震台网在做地震编目时,采取了多种定位方法,如单纯型、hypoSAT定位方法等,最终给出的地震目录震源位置尤其是深度上存在着较大的误差[10],通常会采取一些地震精定位的方法来减小定位误差,本文为提高地震定位精度以及反演区域速度模型,使用了Zhang和Thurber[4]提出的双差层析成像(TomoDD)方法对2009年1月1日—2022年2月8日之间的门源及周边地区地震事件进行了精定位处理。最终重定位的地震震中与深度分布如图4所示。重定位后,最终获得14 024个地震事件,空间分布相对于原始地震目录分布地震更加集中,与断层相关性更高。2016年门源MS6.4地震(37.66°N,101.62°E)的深度为9.1 km,与前人研究结果一致[11],2022年门源MS6.9地震的震源深度为13.3 km,震中位置为(37.76°N,101.25°E);重定位地震目录的深度分布图4(c)比原始地震目录图4(b)更趋于均匀以及正态分布,主要集中在5～15 km,反映出门源及周边地区的地震震源深度较浅;10 km左右的地震分布最多,说明该深度可能为门源及周边地区的主要孕震区。

图4 反演前后震中分布与震源深度对比Fig.4 Comparison of epicenter distribution and focal depth before and after inversion

3.4 合成棋盘格测试

地震震相分析中,由于S波出现在P波之后,往往受尾波叠加干扰导致震相拾取准确率较低,因此在成像部分我们只讨论P波成像结果。我们根据之前TomoDD研究中常使用的合成棋盘测试的分辨率估计来评估最终模型的分辨率。在本文研究中,P初始模型的速度受到±5%的扰动,并在每个节点上交替。最终深度切片的合成棋盘格测试结果如图5所示。

图5 研究区棋盘格测试结果Fig.5 Chessboard test results of study area

在合成棋盘分辨率测试中,我们设定交替节点处的地震速度变化为初始输入速度的±5%。另外将高斯噪声添加到走时选择上—在P波上添加了0.08 s。分辨率最高的切片包括5～20 km深度范围,大约80%以上的棋盘格被恢复,25 km深度处的棋盘格测试恢复度大概在30%左右。在模型的较浅部分(0 km)以及较深部分(25 km以下)观察到大片“花状”,说明速度扰动没有很好地恢复,这主要是由于地表复杂的地形以及地震在地表和25 km以下分布较少,导致速度扰动恢复度差。根据棋盘格测试结果,本文只讨论5～25 km的成像结果。

3.5 三维速度结果分析与讨论

本文在反演过后得到了P、S波的成像结果,总体上P、S波成像结果比较一致,但在局部小范围内有所差异,S波的棋盘格恢复度相对于P波来说较差,因此,为充分认识门源主震区的三维速度结构与发震机理,本文只讨论震相识别准确率较高的P波成像结果,从P波水平成像结果(图6)上看,研究区具有明显的水平不均匀性以及垂向不均匀性,这充分反映了研究区强烈的速度结构不均匀性容易产生地震的构造现象。为查看研究区的连续深度速度结构变化,本文分别在2016年门源地震以及2022年门源地震主震位置做了两条沿经度和纬度的深度剖面(图7),剖面的水平位置为图6中红色虚线所示。

(黄色☆为2016年门源6.4地震,AA′、BB′、EE′为2016门源地震主震处的深度剖面;红色☆为2022年门源6.9地震,CC′、DD′、FF′为2022门源地震主震处的深度剖面)图7 P波剖面成像结果Fig.7 P-wave profile imaging results

如图6显示,2016年门源MS6.4以及2022年门源MS6.9地震都发生在长约40 km宽约20 km的高速体边缘,与前人研究结果相一致[7],其中2016年地震位于高速体的东端,2022年门源地震位于更深处的高速体西端末梢位置,在15 km深度,高速异常体面积大为减小,说明10～15 km中间是分界面,可能为滑脱层。10 km深度处的冷龙岭断裂以及门源盆地北缘断裂与该高速异常体近乎平行,可能为主导门源地区地震的高速异常体,受青藏高原地块整体构造向北东挤压的大应力环境,高速异常体向北东方向挤压,但受到研究区东北角的高速坚硬物质阻挡,导致高速异常体往深部插入[图7(a)],因而造成东北的断层上盘上升,产生了2016年逆冲性质的门源MS6.4地震。2022年门源地震余震最东端为101.53°E,离2016年门源主震(101.62°E)仍相隔约10 km,并没有引起2016年主震位置处的破裂,这也许是高速体已紧嵌入东北地层中,形成新的稳定地块。

图7为2016年门源6条剖面的成像结果,地震分布取剖面线两侧各0.05°之间的地震事件,从图7(d)、7(f)上都可以看出2022年门源MS6.9地震深度分布具有明显的拐弯现象,地震分布主要沿着高速体的边缘,EE′和FF′两条剖面上可以看出两次地震序列的倾角较陡。2022年门源MS6.9地震的余震展布上,在水平空间上主要沿冷龙岭、托莱山断裂分布,剖面上显示陡立、略南倾的地震条带。此次2022年门源地震为典型的走滑型事件,与祁连断裂性质接近,可以解释为青藏高原地块向阿拉善地块挤压过程中,遇阻而顺时针滑动旋转的结果。

4 结论

本文讨论了门源及周边地区三维速度结构与地震活动、深部孕震环境和青藏高原东北缘构造背景之间的关系。所得结果为了解门源地区地震发展的动力机制、预测该地区地震活动趋势和地震灾害风险提供了依据,并有助于最大限度地减少未来此类灾难造成的损害。经过本文的成像处理研究后,我们得到了此次地震以及2016年门源MS6.4地震的二维断层面几何形状、三维构造环境。根据本文研究结果主要得到以下结论:

(1) 研究区地震震源深度分布主要集中在5～25 km范围内,其中2016年门源MS6.4地震(37.66°N,101.62°E)的深度为9.1 km,与前人研究结果相一致[11],2022年门源MS6.9地震的震源深度为13.3 km,震中位置为(37.76°N,101.25°E),更靠近冷龙岭断裂带;

(2) 2022年门源MS6.9地震的余震分布上,平面上沿冷龙岭、托莱山分布,剖面上显示陡立、略南倾的地震条带;

(3) 成像结果显示,2016年门源MS6.4地震以及2022年门源MS6.9地震都发生在高速体的边缘,速度结构与断裂、地震序列吻合较好;

(4) 2022年门源MS6.9地震位于高速体的西端末梢位置,是该高速体受青藏高原东北缘顺时针旋转导致的滑动产生的走滑型地震事件。

致谢:感谢甘肃、青海测震站网提供的震相报告资料,另外本文做图使用了GMT绘图软件[12],两位匿名审稿人提出的建设性意见对本文质量提升很大,在此深表感谢。