2020年土耳其MW5.7地震InSAR同震形变特征与破裂滑动分布

李承涛 李 琦 谭 凯 鲁小飞 左 啸

1 中国地震局地震研究所，武汉市洪山侧路40号，430071 2 引力与固体潮国家野外观测研究站，武汉市森林大道196号，430000 3 中国地震局地震大地测量重点实验室，武汉市洪山侧路40号，430071 4 湖北省地震局，武汉市洪山侧路48号，430071 5 武汉市测绘研究院，武汉市万松园路209号，430022

据中国地震台网测定，北京时间2020-06-14 22：24土耳其(40.75°E, 39.40°N)发生MW5.9地震，震源深度约10 km。根据土耳其灾害和紧急情况管理局(Turkey’s Disaster and Emergency Management Authority, AFAD)发布的相关研究报告(https://deprem.afad.gov.tr/depremdokumanlari/1836)，此次地震矩震级为MW5.7，持续时间约15.04 s。AFAD依据T50降轨数据计算出，主要的视线向形变量约7 cm。

如图1所示，研究区构造较为复杂，主要受安那托利亚板块、阿拉伯板块、欧亚大陆板块三大板块的共同作用，汇聚于卡尔勒奥瓦三联点(KTJ)，其结果使得安那托利亚板块向西不断运动，发育的主要断层为北安那托利亚断层与南安那托利亚断层[1]。其中，北安那托利亚断层是以右旋走滑为主的断层，且兼具少量的逆冲特性，是一条1 200 km长的走滑断层，从卡尔勒奥瓦向西延伸；南安那托利亚断层以左旋走滑为主[2-3]。GPS数据解算出北安那托利亚断层右旋走滑量为24 mm/a，南安那托利亚断层左旋走滑量为10 mm/a[4]。1976～2020年MW≥5地震目录( https://www.globalcmt.org/CMTsearch.html)统计结果显示，截至2020-07-01，该区域共发生49次5级以上地震，地震主要沿着这两条大断层分布。其中MW≥6.0地震共5次，最近的一次为2020-01-24发生在南安那托利亚断层附近的MW6.8地震。

图1中蓝色矩形框(40.55°～41.05°E, 39.25°～39.55°N)内为选取的重点研究区域。本次地震发生在该范围内，位于北安那托利亚断层附近，该区域历史地震活动性较强，聚集了大量5～6级地震，研究本次地震的相关特征，有利于分析所在区域的地震危险性。

本次地震附近GPS站点稀疏且未观测到明显的同震形变。另外，地震波数据一般受台站分布的影响，反演得到的发震断层的走向和倾角精度存在较大误差，且存在中强地震震中定位不精确的现象。如图1所示，GCMT、AFAD 、USGS等3个机构提供的震中位置存在差异，这给后续的同震破裂研究带来诸多不便。近场的大地测量数据，比如InSAR数据，能够直观地获取近场形变数据，可以很好地约束发震断层的几何形状参数。

红色矩形代表T50降轨范围，3个红色震源机制解分别由GCMT、AFAD、USGS提供，蓝色震源机制解表示2020-06-15发生的MW5.5余震图1 构造地震背景Fig.1 Tectonic and seismic background

综上，本文利用Sentinel-1A数据，通过D-InSAR技术提取同震形变场，采取贝叶斯自举优化法研究发震断层几何形状特征[5]，用有限断层方法反演断层破裂滑动分布[6]，为准确认识此次土耳其MW5.9地震的形变特征以及发震构造提供参考。

1 InSAR数据数据处理和形变分析

本文使用欧空局(https://sentinel.esa.int/)升降轨Sentinel-1A干涉宽幅数据，对于升轨T43数据，震前影像拍摄时间是2020-06-02，震后影像拍摄时间是2020-06-14(具体时间是15:18:35.017，距震后约1 h)，经过处理后的干涉图没有明显的同震信号，是噪声过大所致，所以舍弃该结果。考虑到2020-06-15 06:51:31(UTC)该区域发生MW5.5余震，如果影像时间跨度包含余震事件，则在后续的同震形变提取中难以剔除，所以选取的影像时间跨度上应尽量不包含MW5.5余震信号。综上，本文采用降轨T50数据，震前影像拍摄时间是2020-06-03，震后影像的结束拍摄的准确时间为2020-06-15 03:17:50.848(UTC)，影像相关的具体信息如表1所示。

表1 Sentinel-1A的详细参数

使用InSAR scientific computing environment (ISCE)[7]对降轨数据进行D-InSAR处理，采用欧空局发布的精确轨道文件来减小轨道误差的影响，采用30 m分辨率的SRTM数字高程模型消除地形的影响[8]，利用Goldstein滤波方法提高信噪比[9]，基于Snaphu(statistical-cost, network-flow phase-unwrapping algorithm)[10]算法进行相位解缠。然后，通过地理编码获取土耳其地震降轨的同震形变场(图2)。利用generic atmospheric correction online service(GACOS)模型[11](http://ceg-research.ncl.ac.uk/v2/gacos/)消除对流层大气延迟所引起的误差，最终得到经大气改正后的形变场(图3)。

图2 土耳其地震降轨的同震干涉图Fig.2 Coseismic interferogram of descending track of Turkey earthquake

图3 降轨LOS向形变场GACOS大气校正Fig.3 GACOS atmospheric correction for LOS deformation field of descending track

从图3(a)～(c)可以看出，通过GACOS模型改正，削弱了部分对流层的影响，大气校正的数值范围为 -0.95～0.8 cm，同震形变场得到改善，改正后的视线向最大沉降量为 -7.75 cm，最大抬升量为8.87 cm。结合干涉图可以发现，此次土耳其MW5.9地震的地表形变影响范围约15 km×10 km，降轨形变场长轴方向近EW向，LOS向形变场呈椭圆状(关于东西向基本对称)。

图4为降轨形变场横截面示意图，图4(a)中AB、CD表示降轨LOS向形变场中横截面的位置，图4(b)为对应横截面的LOS向同震位移变化。横截面AB(黑线表示)靠近CENC震中，从A点至B点，LOS向形变值由约-0.08 m变化至0.08 m，且呈现中心对称，显示出走滑地震特性；横截面CD变化相似，范围在-0.06～0.06 m。LOS向形变量绝对值大于0.02 m的区域范围约10 km2，表明此次地震影响范围较小。此外，图中还显示出不同机构确定的震中，黑色五角星表示CENC确定的震中，白色五角星代表本文利用InSAR数据确定的震中，位于LOS向形变场的形变中心附近。

图4 LOS向形变场的截面分析Fig.4 Cross section analysis of LOS deformation field

2 断层几何形状和滑动分布反演

为加快反演进程，采用四叉树方法[12]对降轨数据进行降采样，考虑到近场形变的重要性，在降采样的过程中，对近场形变进行加密采样，对远场区域进行稀疏采样，最终保留150个降轨数据点。

2.1 断层几何形状反演

利用贝叶斯自举优化法(Bayesian bootstrap optimization, BABO)估算发震断层相关的几何形状参数[5]。BABO基于直接搜索，从定义的模型空间中抽取随机模型参数，然后计算综合模型值，并将其与目标的观测数据进行比较。BABO优化的核心是对观测数据和预测数据之间的差异进行逐点计算，求取观测值与模型值最小化差异：

min ‖dobs-dsynth‖2=

(1)

式中，dobs为观测值，dsynth为模拟值，i(i=1,2,3,…)为观测点个数。

BABO将断层面作为一个均匀面来处理，优化得到断层面参数(断层长、宽、倾角、走向、滑动角等)。本文的观测数据是经过四叉树降采样处理后的T50降轨Sentinel1-A数据，所有优化参数序列如图5所示，可用于检查收敛性以及查看模型参数是否推动了给定的边界范围，较高的失配值为冷蓝色，较低的失配值为红色，失配值越低表明参数优化越好。对所有估计样本值进行排序，置信区间的下限取排序后的5个百分位数，置信区间的上限取排序后的95个百分位数，取两者之间的数值，得到各参数的90%置信区间。各参数迭代次数均为50 000，最终收敛状况都较好，且失配值随着迭代的不断进行而不断地降低。

图5 参数序列图Fig.5 Sequence plots for all parameters

图6表示模型参数的分布，x轴绘图范围表示给定的初始参数边界，y轴表示概率密度函数(probability density function，PDF)，本文取高斯核密度函数。图中显示模型参数分布十分集中，误差较小。最终各参数的具体数值见表2，以下仅列出最优解：东偏移、北偏移表示UTM坐标系下相对于参考点(40.68°E，39.35°N)(GCMT震中)的偏移量，物理含义上表示反演得到的理论震中位置，转化成经纬度即为(40.754°E，39.389°N)。此外，断层上边界中心点深度为1 308 m，断层长度为5 680 m，断层宽度为5 710 m，滑动量为0.223 m，断层走向、倾角、滑动角分别为257.48°、79.69°、154.2°。

红色曲形实线为高斯核密度; 红色虚线为参数的平均值; 90% 置信区间表示为中间的红色阴影区，阴影区最宽边界代表最小、最大值； 黑色实线表示GCMT参数图6 2020年MW5.7土耳其地震断层模型参数的分布Fig.6 Parameter distribution of the 2020 MW5.7 Turkey earthquake fault

表2 反演的断层几何形状参数

在反演断层几何参数时，得到的残差结果如图7所示。图中圆点表示采样点，矩形表示反演断层所在位置，矩形中的黑点表示震中，断层左上角点经纬度为(40.780 4° E，39.390 8° N)，对应深度为1.308 km。在靠近震中附近，LOS向视向线模型值都能够较好地拟合观测值，降轨数据残差范围是-0.018～0.012 m，残差的RMS为6 mm，总体来看模拟效果较好。

图7 断层几何参数反演获得的视线向形变残差Fig.7 LOS deformation residuals of fault geometry inversion

2.2 破裂滑动分布反演

基于上述BABO反演获得的几何形状参数，采用有限断层方法对发震断层滑动分布进行反演。为了获得较精细的地表变形，本文的震源矩形断层采用三角格网[13]构建，断层位错引起的地表变形通过弹性位错公式计算，其本质上属于非负最小二乘方法，即寻求观测值拟合度和滑动分布粗糙度最小化：

‖W(GS-d)‖2+β2‖LS‖2=min

(2)

式中，d取InSAR视向线数值；W为观测值的权矩阵；G为格林函数，可以通过基于半无限空间地壳分层的弹性位错模型[14]计算获得；β为子断裂滑动平滑因子；L为拉普拉斯二阶差分算子；S为待求参数，包含子断层滑动量以及InSAR轨道相关的改正数[15]。

断层参数设置：通过观测InSAR干涉图和LOS向形变场的空间特征，参考贝叶斯自举优化法得到的结果(表2)，走向取257.5°，倾角取79.7°，延伸反演断层的长宽，设置断层长25.5 km，宽9.0 km。坐标原点取断层左上点(40.780 4°E，39.390 8°N)，对应深度1.308 km。反演采用降轨相关数据。

最终获得的土耳其地震同震滑动分布如图8所示。该地震主要破裂区域长度约8 km，宽度约6 km，破裂没有出露地表，位于地表以下约0.8 km，最大埋深约8.9 km。主要破裂集中在深度3～6 km的区域，较深区域的子断层滑动角几乎垂直，滑动角随着深度的增大而增大，在接近地表的地方子断层几乎呈水平右旋滑动，主要滑动区域的平均滑动角为146.8°。主要破裂区左侧存在少量的滑动现象，可能是由于本文采用单一断层模型所致。同震破裂主要表现为右旋走滑，且兼具少量逆冲特性。断层面上最大滑动量约0.57 m，位于(40.769°E，39.394°N)处，埋深约4.278 km。假设地壳的剪切模量为3×105bar，此次地震释放的地震矩约4.54×1017Nm，对应矩震级MW5.7。滑动分布模型的形变残差如图9所示，模拟效果较好，降轨拟合残差的范围是为-3.9～4.0 cm，残差的RMS为0.7 cm。

图8 利用InSAR数据反演的土耳其地震滑动分布Fig.8 Coseismic slip distribution of Turkey earthquake inverted by InSAR

图9 滑动分布模型的形变残差Fig.9 The deformation residuals of the distributed slip model

3 讨 论

3.1 震中位置分析

根据表3，InSAR反演的震中与CENC测定的震中最接近，相差1.3 km。结合图4分析，InSAR得到的震中更符合实际震中特征，更为精确，且更接近最大形变值所处的位置，而与USGS、AFAD测定的震中在平面空间上分别相差约4.1 km、4.4 km，与GCMT测定的震中距离最远，相差约7.7 km，这可能主要受测震台网的分布、测震站台的位置以及地震震级(较小)等影响。近场数据测定的震中明显要比USGS、GCMT更加准确，这表明用近场InSAR形变场数据可以更好地确定中强地震的震中位置。

表3 震源参数比较

3.2 潜在危险性分析

Poyraz[16]研究认为，NAF最东端的平均滑移率为21 mm/a，闭锁深度为12.72 km。此次MW5.9地震发生在NAF的东段，造成震中局部地区近东西向的右旋剪切作用。为评估此次地震对土耳其地区带来的影响，假设摩擦系数为0.4，接收断层的走向为257.48°，倾角为79.69°，滑动角为154.2°(即上述反演结果，见表2)，得到最大滑动量所处深度4.278 km对应的库伦破裂应力(CFS)变化[17]，ΔCFS>0表示库伦破裂应力增加，ΔCFS<0表示库伦破裂应力减少，如图10所示。值得注意的是，当ΔCFS≥0.1 bar时，就有诱发后续地震的可能性[18]。由图可见，南安那托利亚断层所在区域库伦破裂应力略有增加，增量小于0.1 bar，表明此次地震对其影响较小。而对于震中的西边部分区域，库伦破裂应力变化较大，区域库伦破裂应力增量达到0.2 bar，2020-06-15 MW5.5地震就发生在这片区域，进一步表明此次土耳其MW5.9地震直接触发了MW5.5余震。图中其他库仑破裂应力增量超过0.1 bar的区域，尤其是北安那托利亚断层的部分区域，如果处于应力集中状态且接近临界值，那么未来其地震危险性都将大大提高。

绿色五角星表示本文反演的震中图10 深度位于4.278 km时的库伦破裂应力变化Fig.10 The Coulomb failure stress (CFS) change at a depth of 4.278 km

4 结 语

本文利用Sentinel-1A宽幅数据获取2020-06-14土耳其地震的同震形变场，并采用BABO和有限断层方法反演断层几何形状参数和断层破裂滑动分布，得到以下几点认识：

1)土耳其地震降轨LOS向同震形变场显示，断层北侧抬升，南侧沉降。降轨LOS向同震形变场最大沉降量约-7.75 cm，最大抬升形变量约8.87 cm，且关于东西向大致对称。

2)断层几何参数反演显示，破裂断层走向约257.48°±0.65°，倾角约79.69°±0.98°，滑动角约154.2°±3.8°，震中位置为(40.754°E，39.389°N)，通过与其他机构公布的震中比较得到，近场InSAR获取的LOS向形变场数据可以更好地确定中强地震的震中位置。

3)断层滑动分布反演表明，破裂区域长度约8 km，宽度约6 km，破裂没有出露地表，破裂的最浅埋深约0.8 km，最大埋深约8.9 km，平均滑动角约146.8°。此次土耳其地震是一次以右旋走滑为主兼具少量逆冲特性的地震。最大滑动量约0.57 m，对应深度约4.278 km。假设地壳的剪切模量为3×105bar，此次地震释放的地震矩约为4.54×1017Nm，对应矩震级MW5.7，与AFAD提供的震级一致。

4)对于ΔCFS≥0.1 bar的区域，未来危险性值得更多关注。

5)此次土耳其地震受近东西向的潜伏断层控制，由于靠近北安那托利亚断层，可能还受到该断裂的共同作用。

致谢:感谢欧洲空间局(ESA)提供Sentinel-1A卫星升降轨SAR数据。