2023年土耳其“双大震”InSAR同震形变观测及断层滑动分布反演

郝小叶 张永志,2 柳姣姣 程雅雯 吴文强 周 甜

1 长安大学地质工程与测绘学院,西安市雁塔路126号,710054

2 地理信息工程国家重点实验室,西安市雁塔路中段1号,710054

根据中国地震台网中心(CENC)测定,2023-02-06 09:17:37和18:24:50土耳其境内发生2次7.8级大地震,依次为帕扎尔西克地震和伊尔比斯坦地震,全球地震矩心矩张量(GCMT)和美国地质调查局(USGS)也相继公布“双大震”的震中位置及震源机制解(表1)。两次大震震中位置相距约96 km,初步判定均为走滑型破裂机制,且都发生于东安纳托利亚断裂系中[1]。据记载,土耳其境内近 25 a共发生8次7级以上地震,但在过去一个世纪北安纳托利亚断裂带超强的地壳运动背景下,此次地震序列急转于东安纳托利亚断裂系,且近乎同时发生2次大地震;此外,地震震中处于人口密集区,事后还多次发生震级较大的余震,人民生命财产安全和经济发展都遭受重创。因此,研究本次土耳其“双大震”的发震机理,有助于深入了解东安纳托利亚断裂系的孕震机制,对该区域地震危险性分析具有重要的参考价值。

表1 2023年土耳其“双大震”地震信息

近几十年来,InSAR技术在地球物理学领域得到广泛应用,因其具有高时空分辨率、全天时、全天候等优势,使得高效获取近场地表形变信息成为可能,特别是差分干涉测量(D-InSAR)与像素偏移追踪(POT)技术的联合分析可提取到更加完整、丰富的同震形变信息,现已多次应用于强震形变监测分析和震源参数反演研究中[2]。本文采用欧空局(ESA)提供的2023年土耳其“双大震”破裂区的完整Sentinel-1 A影像数据,解算此次事件的近断层同震形变场,并基于“四叉树”降采样后的像素偏移量进行发震断层的精细化滑动分布反演,进而分析和探讨此次地震序列的震源机制情况。

1 地质构造背景

2023-02-06土耳其“双大震”发生于阿拉伯板块和安纳托利亚块体的交界带——东安纳托利亚断裂系中,是世界上最活跃的大型转换断裂带之一,接纳了安纳托利亚块体的长期西向挤压运动[3-4]。安纳托利亚板块是在欧亚板块与阿拉伯板块的长期碰撞,以及非洲板块东北向的不断俯冲作用下分离产生的一个微板块,在各板块汇聚处形成以东安纳托利亚断层和北安纳托利亚断层2条大型走滑断裂为主的超强活动构造结[5]。地质资料表明,“双大震”破裂断层均处于以东安纳托利亚断层为主的断裂系中(图1)。

图(a)中蓝色方框表示T14升轨覆盖范围,棕色方框表示T21降轨覆盖范围,白色方框表示图(b)覆盖范围;红色实线表示安纳托利亚板块边界的主要断层结构(东安纳托利亚断层和北安纳托利亚断层数据取自文献[6]);红色五角星及震源沙滩球表示由GCMT提供的“双大震”震中位置及震源机制解,黑色和蓝色沙滩球分别表示自1960年以来6级以上历史地震

东安纳托利亚断裂带(EAFZ)是地中海区域的主要陆内转换断裂带之一,全长约550 km,为左旋走滑型断层,与右旋走滑性质的北安纳托利亚大型断裂带(NAFS)交会于东北向的卡尔勒奥瓦三联点处,且在西南方向延伸至安塔基亚附近与死海断裂带毗邻[4,7]。现有研究结果表明,东安纳托利亚断层滑动速率由东北端向西南端递减,形变速率值介于4～10 mm/a,远低于北安纳托利亚断层的平均滑动速率(21～25 mm/a)[8]。研究区内1960年以来共发生6次5级以上地震: 1964年宾格尔MS6.9地震、1971年宾格尔MS6.7地震、2003年宾格尔MW6.3地震、2010年埃拉泽MW6.1地震和2020年埃拉泽MW6.8地震[9],均与东安纳托利亚断层地质构造活动有关。东安纳托利亚断层的几何形状较为复杂,断层南端存在多条次生断层,其中索尔古和恰尔达克断层均呈现北倾的左旋走滑特性[4]。该地震序列的发生表明,东安纳托利亚断层近年来的地震活动仍十分活跃,因此解算此次“双大震”序列的同震形变场,分析发震断层的破裂机制,对东安纳托利亚断层近年来的运动特征及断层固有特性的进一步研究具有重要意义。

2 土耳其“双大震”同震形变场

2023-02-06土耳其“双大震”发生后,欧空局先后发布覆盖该震区的Sentinel-1 A卫星影像(表2),本文采用卫星宽幅干涉模式时间基线较短的升、降轨干涉像对进行两轨法差分干涉,以获取地表形变信息。据美国地质调查局(USGS)和全球矩心矩张量(GCMT)机构测定的土耳其“双大震”震源参数信息,在差分干涉处理前,分别对Track14(升轨)和Track21(降轨)的主、从影像对进行拼接处理。获得相应的卫星影像后,基于GAMMA软件分别对升、降轨影像对作两轨法差分干涉,并采用美国宇航局(NASA)公布的SRTM 30 m数字高程模型消除地形相位的影响[10];设置影像距离向和方位向视数比为10∶2进行多视处理,抑制干涉相位的噪声误差,并使用Goldstein方法进行滤波,从而消除噪声因素的影响;采用最小费用流法(MCF)对影像作相位解缠处理,对相干系数大于0.4的像元进行解缠,并基于GACOS大气改正系统消除大气延迟误差;通过地理编码实现影像数据SAR卫星坐标系到WGS84坐标系的投影变换,得到土耳其“双大震”的升、降轨LOS向形变场(图2)。

图2 2023年土耳其“双大震”同震形变场 (D-InSAR LOS向形变)Fig.2 Coseismic deformation field of the Türkiye conjugate earthquakes in 2023 (D-InSAR LOS direction)

表2 Sentinel-1 A数据参数

由图2可知,升、降轨 LOS 向形变场呈现出2个明显的同震形变区域,表明此次“双大震”的同震形变分别集中在东安纳托利亚主断层和索尔古分支断层附近,最大LOS向形变值为1.76 m,出现在MW7.5地震形变区。但土耳其“双大震”的破裂断层附近存在大量失相干缺值现象,未完整表现出此次“双大震”的同震形变,会降低发震断层模型反演分析的准确性。为呈现出完整的近断层地表位移,本文基于GAMMA软件进行像素偏移追踪处理,该技术采用影像的幅度特征提取位移信息,不受影像对相干性的影响。研究表明,POT技术获取地表形变的准确性取决于像素尺寸大小及影像对的配准精度,精度可达到像素大小的1/10～1/100[11],但由于该技术基于像素集合获取位移,相较于LOS向形变场空间分辨率显著降低。本文利用影像强度的互相关技术获取SAR影像对的最佳匹配同名点,选取10∶2的距离向与方位向多视比及320×64的搜索窗口进行影像对的偏移追踪处理,使用一次多项式曲线拟合去除像素偏移量的轨道趋势误差,再选用9×9像素的中值滤波窗口进行去噪,最终得到地震前后的距离向像素偏移(图3)和方位向像素偏移(图4)。从图3、4可以看出,升、降轨数据获取的距离向偏移量效果明显优于方位向偏移量,这是因为,像素偏移追踪技术的精度还与影像分辨率有关,而Sentinel-1 A数据的距离向分辨率远大于方位向分辨率。此外,升、降轨同震形变场位移反向性显著,与东安纳托利亚断裂系的左旋走滑性质一致。经过POT技术获取的形变场表现出非常明显的变形不连续性,提取近断层位移不连续面来识别地震破裂迹线,这些位移不连续面很可能代表了MW7.8帕扎尔西克地震和MW7.5伊尔比斯坦地震的发震断层,初步估计2次大地震的发震断层破裂长度分别约为260 km和110 km。

图3 2023年土耳其“双大震”同震形变场 (距离向偏移量)Fig.3 Coseismic deformation field of the Türkiye conjugate earthquakes in 2023 (range offset)

图4 2023年土耳其“双大震”同震形变场 (方位向偏移量)Fig.4 Coseismic deformation field of the Türkiye conjugate earthquakes in 2023 (azimuth offset)

分析此次地震序列的同震形变场(图2、3)可知,升轨影像仅覆盖此次“双大震”的部分破裂区,但“双大震”破裂断层两侧的地表运动变化表现出良好的对称性,地表形变范围为-3～3 m,量级较大,其中正值表示地面向靠近卫星方向运动,负值表示向远离卫星方向运动。由于POT技术获取的形变场覆盖近断层区域,最大像素偏移量(2.94 m)远大于失相干的最大LOS向形变值(1.76 m)。在升轨形变场中,形变量正值出现在发震断层北端,负值则在发震断层南端;相反,降轨形变场中正、负形变值分别对应于发震断层的南、北两侧,表明此次“双大震”的破裂机制均呈现出左旋走滑特性,与USGS和GCMT所发布的震源机制解一致。

3 土耳其“双大震”滑动分布反演

基于POT技术获取的地表距离向像素偏移量可以清晰展现此次地震序列的发震断层迹线(图3),与地震序列的余震分布也高度一致(图1(b))。本文通过获取的升、降轨距离向偏移量,使用基于约束条件下的最小二乘原理及最速下降法SDM反演程序[12]进行发震断层的滑动反演研究,其核心思想是基于负梯度方向的分析调整每次迭代计算的搜索方向,使得待优化的目标函数随着迭代次数的增加而逐渐减少。

为深入探讨土耳其“双大震”的震源机制,了解此次地震的震源参数及断层滑动分布趋势,以震中附近的距离向位移量为约束进行反演分析。首先选用“四叉树”方法对形变数据进行降采样,高密度选取影响力占比权重较大的近场形变区域,最终保留1 632个升轨数据点和1 714个降轨数据点。基于断层滑动的非均匀性构建分布式滑动模型,采用最速下降法进行发震断层的精细化滑动反演。在反演过程中将地壳视为均匀介质模型,发震断层细分为4 km×5 km小单元片;根据已发表的震源机制解(表1)和前人对土耳其活动断裂的研究成果[6,13]初步设定发震断层F1和F2(图3、4)的几何倾角初始值分别为88°和78°,搜索范围分别为(70°, 90°)和(40°, 90°)。通过不断改变拉普拉斯滑动因子以获取最优解,帕扎尔西克地震和伊尔比斯坦地震的断层滑动分布结果如图5、6所示,划分的子断层块滑动方向以箭头表示。

图5 2023年土耳其“双大震”断层滑动分布(二维)Fig.5 The fault slip distribution of the Türkiye conjugate earthquakes in 2023 (2D)

图6 2023年土耳其“双大震”断层滑动分布(三维)Fig.6 The fault slip distribution of the Türkiye conjugate earthquakes in 2023 (3D)

分析图5、6可知,此次地震序列的地表破裂主要集中在东安纳托利亚主断裂带及索尔古断层上。其中,东安纳托利亚主断层的破裂由帕扎尔西克地震触发,地表破裂长度约260 km,断层破裂深度在0～25 km范围内,最大滑动量约为10.5 m,位于东安纳托利亚断层的S5子断层地表区域,平均滑动角为3.74°,属于左旋走滑型兼少量逆冲分量地震;伊尔比斯坦地震的破裂断层主要为索尔古断层和恰尔达克断层,破裂深度范围为0～30 km,最大滑动量为11.8 km,平均滑动角为1.05°,以左旋走滑分量为主。此外,分析“双大震”序列后的余震分布和同震形变场可知,MW7.8地震还触发了努尔达吉-帕扎尔吉克分支断层上微弱的地表破裂,但USGS、GCMT解算的震源机制解及Li等[14]的研究结果均表明,MW7.8地震可能在该断层成核,推测该断层的微运动会显著影响主断裂带的地壳运动。两次地震发震断层的破裂机制均与主断裂带长期研究中的左旋走滑性质保持一致,且破裂断层机制与“双大震”后的余震分布结构高度匹配,表明反演结果具备合理性。

为进一步验证反演结果的可信度,对模型拟合的形变值和观测形变值进行对比分析,模拟残差效果如图7(图中两个黄色五角形,上为MW7.5地震,下为MW7.8地震;图(a)、(b)、(c)为升轨数据,图(d)、(e)、(f)为降轨数据)所示。滑动分布反演的模拟值与观测值之间的相关度达到92.74%,说明模型拟合效果良好,可信度高;拟合值和观测值之间的残差较小,基本分布在0值附近,断层线附近出现少许高残差值(约为0.61 m)是由于真实的发震断层线较为复杂,存在多个分支断层结构,且各断层的倾角存在多样性,与模拟选取的断层结构之间存在较大差异。另外,在形变值远场还存在部分较高残差值,且研究得出的地震矩与USGS、GCMT等机构存在较大差异,分析可能与MW6.7余震的地质活动有关。鉴于主震与余震之间的时间间隔及空间距离的紧密性,获取的形变值无法剔除此次强余震的影响。

图7 2023年土耳其“双大震”滑动分布拟合结果Fig.7 The fitting results of distributed slip model of the Türkiye conjugate earthquakes in 2023

4 讨 论

2023年土耳其“双大震”发生在安纳托利亚板块与阿拉伯板块的交界带——东安纳托利亚断裂系中,Sentinel-1 A卫星成功覆盖此次地震序列的形变区域,震区植被覆盖稀疏,但矩震级较大,地震矩释放及破裂程度明显,导致影像数据的相干性较差,故联合D-InSAR技术与POT技术进行同震形变的获取分析。土耳其“双大震”序列的地震矩总释放量约为1.24×1021N·m,略大于Li等[14]联合POT技术获取的距离向偏移量与GPS速率反演得到的9.62×1020N·m地震矩量,可能和反演时引入的GPS数据以及选取的断层参数差异有关。该结果与USGS、GCMT发布信息之间存在差异,可能是受Sentinel-1A影像重访周期的限制,将2次事件的发震断层滑动反演同步进行所引起。

本文研究结果表明,MW7.8帕扎尔西克地震震源位于东安纳托利亚断层系的努尔达吉-帕扎尔吉克断裂上,该分支断层破裂程度较小,但会显著影响东安纳托利亚主断裂带的地壳运动,地表破裂长度约260 km,主断裂东北方向上S3、S5部分的滑动破裂明显。MW7.5伊尔比斯坦地震发生在东安纳托利亚断层系的分支断裂——索尔古断层和恰尔达克断层上,断层破裂深度在20～30 km之间,沿断层走向无明显横向变化,破裂长度约110 km,远小于MW7.8地震所引起的破裂长度,但此次“双大震”序列的最大滑动量出现在索尔古分支断裂的地表位置。研究表明,MW7.8帕扎尔西克地震从努尔达吉-帕扎尔吉克分支断裂沿北北东向延伸到东安纳托利亚主断裂带上后,沿主断裂带继续向东北方向扩展,随着能量消散破裂方向发生逆转,向西南方向沿主断裂带继续扩展。Zhao等[6]联合InSAR同震形变场及GPS数据得到的发震断层精细滑动分布结果也表明,帕扎尔西克地震所触发的同震破裂在东安纳托利亚主断裂东北部分较为明显,西南端的滑动程度较低,地震序列的最大滑动分布于第2次大地震的形变区域。Li等[14]以GPS数据和距离向偏移量为约束,反演发震断层的精细滑动分布认为,MW7.8地震引发的地表破裂由东北向西南方向逐渐递减,且在努尔达吉-帕扎尔吉克分支断裂上未发生明显滑动,而MW7.5地震的发震断裂出现3个明显的滑动区。由此可知,本文基于同震形变反演的发震断层滑动分布与已有研究成果基本一致,此次土耳其“双大震”地震事件均属于左旋走滑兼少量逆向倾滑分量的地震,该地震序列可解释为相对稳定的东安纳托利亚断裂系的一次超强能量活动。

5 结 语

本文选取欧空局公布的Sentinel-1A影像数据作为研究数据源,基于GAMMA软件使用两轨法和POT技术进行影像干涉,获取土耳其“双大震”同震形变场,基于分布式滑动模型进行断层精细化滑动分布反演,得到以下结论:

1)土耳其“双大震”升、降轨LOS向同震形变场在形态和形变量上呈现出明显的对称性,均属于典型的左旋走滑型地震,地震破裂的最大距离向偏移量约2.94 m,远大于LOS向最大形变量。

2)土耳其“双大震”事件均以左旋走滑兼部分倾滑运动为主,与GCMT、USGS等机构发布的震源机制解一致,与该区域左旋走滑的地质构造背景也一致。其中,MW7.8帕扎尔西克地震发生于东安纳托利亚断层系的努尔达吉-帕扎尔吉克断裂上,未见明显破裂,但会显著影响主断裂带的地壳运动,地表破裂超过260 km;MW7.5伊尔比斯坦地震发生在东安纳托利亚断层系的分支断裂——索尔古断层和恰尔达克断层上,破裂长度约110 km。

3)土耳其“双大震”断层滑动反演结果与原始观测值的拟合度达92.74%,最大滑动量为11.8 m,出现在MW7.5伊尔比斯坦地震的形变区,地震矩总释放量为1.24×1021N·m。

致谢：欧洲航天局(ESA)提供Sentinel-1 A卫星影像数据及GAMMA软件的技术支持,文中大部分图件使用GMT软件绘制,在此一并表示感谢。