0 引言

据中国地震台网测定，2021年5月22日青海省果洛州玛多县发生7.4级地震，震中位于34.59°N，98.34°E，震源深度17km。野外调查和余震精定位结果显示，此次玛多M_S7.4地震发生在巴颜喀拉块体北部边界以南的走向为285°的江错断裂带(高程约4500m)上(李智敏等，2021；王未来等，2021)，该发震断裂是与巴颜喀拉块体北边界东昆仑断裂近平行的次级断裂，距离东昆仑断裂约70km，并接近东昆仑断裂带玛沁—玛曲段处的几何大拐弯。青藏高原的大型走滑断裂将其划分为一系列活动块体，其中，巴颜喀拉块体位于青藏高原东北部，强震频发，是青藏高原内部拥有复杂的变形模式以及近几十年来地震活动最频繁的次级块体(邓起东等，2002、2014；Gan et al，2007)，该块体狭长，EW向约为2000km、SN向约为200km(马玉虎等，2013)，以东昆仑走滑断裂为北边界，甘孜—玉树走滑断裂为南边界，东边界为龙门山逆冲断裂，块体边界断裂活动性强，强震频发，复杂的变形模式和构造应变主要集中在块体边界断裂和逆冲断层上，使得该构造活动区域具有地震发生频度极高、强度极大的特点。巴颜喀拉块体内部发生的强震相对较少，但发育有许多与东昆仑断裂和甘孜—玉树—鲜水河断裂近平行的次级断裂，如玛多—甘德断裂、昆仑山口—江错断裂、巴颜喀拉山主峰断裂、达日断裂、清水河断裂等。从区域上看，昆仑山口—江错断裂可认为是东昆仑断裂带南侧的次级断裂，该断裂带具有最大剪应变率高的特征(张培震等，2022)，潘家伟等(2021)分析认为2021年玛多地震是2001年昆仑山口西M_S8.1地震应变能沿巴颜喀拉块体北部边界向东南方向进一步迁移的结果。对玛多地震资料的初步分析表明，本次地震为左旋走滑断裂，并有微小的正断层分量。然而，该地震没有发生在东昆仑断裂上的原因以及孕震区的潜在地震危险性仍然不清楚，研究这次地震为认识发震机制和区域地震危险性提供了难得的契机。

由强震造成的地表同震形变可通过GNSS和InSAR等现代大地测量技术获取，在地震破裂过程、震源机制和断层滑动分布等方面的研究中得到了广泛的应用(屠泓为等，2016；屈春燕等，2017；单新建等，2017；王阅兵等，2022；Zheng et al，2023；伍吉仓等，2020)。玛多地震发生后，许多研究机构利用GNSS和InSAR等大地测量数据对玛多地震进行研究。王迪晋等(2022)生成的近场GNSS观测结果较为明显地表现出发震断裂两侧的运动特征，研究发现，近场形变在WN-ES方向具有清晰的拉张变形特点，EN-WS方向具有挤压变形特点，近场形变特点与远场的四象限变形特点相吻合，左旋走滑特征在发震断裂两盘的近场同震形变表现与震源机制解相符。余鹏飞等(2022)基于GNSS连续观测、InSAR观测等大地测量数据得到了此次玛多地震的同震形变场，并反演了发震断层的破裂分布模型。不同机构在震后短时间内做出了相关研究(华俊等，2021；李志才等，2021；余鹏飞等，2022；Jin et al，2021；He et al，2021)，但结果有所差异，尤其是最大位错及形变区域差别最为明显。华俊等(2021)基于InSAR升降轨数据模拟了地下位错模型，最大滑动量约6m，破裂带东西两侧滑动量有较大差异，但由于InSAR数据失相干并缺少断层近场约束，近场反演相对残差较大。李志才等(2021)基于21个连续运行卫星定位基准站反演了此次地震的发震断层的几何参数破裂分布，由于数据量小，地下位错模型相对粗糙。Jin等(2021)基于Sentinel-1A/B合成孔径雷达数据，反演获得同震位错模型，认为最大滑动位于震源以东，且最大值位于深度3~4km处。He等(2021)使用干涉合成孔径雷达干涉图和SAR像素偏移反演估计断层的几何形状和同震断层滑动分布，并利用累积实现变形进行震后滑动反演，认为最大同震滑动约5m，早期震后滑移主要发生在主要同震滑动凹凸体的下倾方向。基于此，本文拟在前人的研究基础上收集更多的资料，勾勒更加吻合于InSAR相干性、地表观测和余震分布的破裂迹线，结合GNSS、Sentine-1 InSAR观测资料，运用ENVI SARscape、SDM等软件进行对比计算，获得玛多地震更为精细的震源参数及同震位错分布并对同震场进行计算，以此对本次玛多地震进行深入研究。

1 区域构造背景与现今地壳形变

2021年玛多地震震中处于巴颜喀拉块体内部靠近北部边界的位置(图 1)，位于玛多县城以南38km的黄河乡附近，距块体北部主边界东昆仑断裂带的东西大滩—托索湖段约70km。青藏高原可由多个次级块体组成，巴颜喀拉块体是青藏高原上最具代表性的块体，自新生代以来一直受到侧向压缩，表现出较为均匀的左旋剪切带特性(甘卫军等，2004)，并被数个大断层所包围。巴颜喀拉块体周围边界断裂带的构造运动不仅造成了整个块体E向“逸出”，且加剧了强震的发生(邓起东等，2010；李海兵等，2021)。剪切应变率高值区及边缘通常会引发强震，本次玛多地震便发生在左旋剪切应变率的高值区，玛多及其邻近区域具备EW向断裂产生走滑破裂错动的构造区域背景(江在森等，2003)。

GNSS观测结果显示，青藏高原地壳受到印度板块的碰撞挤压作用，使得巴颜喀拉块体向N61.45°E方向运动，其运动速率约为21mm/a，同时该高原地壳绕东喜马拉雅构造结旋转(Wang et al，2001；张培震等，2003；Zhang et al，2004)。活动断裂研究表明，东昆仑断裂带中阿尼玛卿挤压阶区东侧的玛沁—玛曲段比西侧的托索湖段的左旋滑动速率大2mm/a(Li et al，2011)。与东昆仑断裂带和玉树—甘孜—鲜水河断裂带近平行的次级断裂，其运动速率要比巴颜喀拉块体南北边界的小得多(徐锡伟等，2007；潘家伟等，2021)。玉树—甘孜—鲜水河断裂带东侧的鲜水河断裂比西侧的甘孜玉树断裂的左旋滑动速率至少大2mm/a(Wen et al，2003；Chen et al，2016；Bai et al，2018)。

2 2021年玛多地震同震InSAR数据处理及分析

Sentinel-1系列卫星是由欧洲航天局(ESA)发射的C波段并能在IW(Interferometric Wide)观测模式下进行宽幅数据扫描(幅宽250km)的主动微波遥感卫星，其能有效地提高影像质量，也能增强影像的干涉性能。利用Sentinel-1卫星在玛多地震后开展了应急观测，基于此，本文得到了时间基线为6天的升、降轨InSAR同震形变场数据。

结合前人的研究(华俊等，2021；余鹏飞等，2022；颜丙囤等，2022)，本文使用ENVI SARscape软件对升、降轨InSAR数据进行处理。在干涉处理过程中，为消除地形相位和进行地理编码，采用美国宇航局公开发布的STRM 3sec数字高程模型，为降低噪声和提高信噪比，采用10︰2的多视处理和自适应滤波技术对干涉图进行处理(Goldstein et al，1998)，之后采用最小费用流法(Minimum Cost Flow，MCF)进行相位解缠，为消除残余轨道误差，采用多项式拟合的方法进行轨道精炼和重去平，最后经相位转形变得到2021年玛多M_S7.4地震的同震形变场(图 2)。

Sentinel-1 A/B升、降轨影像能完整覆盖此次地震震中及周边区域，同震形变场显示整体呈椭圆形，断层走向即形变场长轴方向整体上为105°。由于受到卫星进行升、降轨观测时的LOS向不同的影响，使升、降轨同震形变场上的发震断层南北两盘的形变量符号相反(图 2(b)、2(d))。此次地震产生了非常剧烈的地表位移，升轨最大LOS向位移达1.1m，降轨最大LOS向位移达0.9m，地表形变范围广，地表破裂长度约160km，清晰可见地表破裂附近由于形变梯度大而产生严重的失相关，导致干涉条纹缺失且不连续，断层远场的升、降轨干涉条纹整体上连续光滑(图 2(a)、2(c))。靠近断层近场的极震区位置有多个不连续的干涉条纹密集的同心圆，推测玛多地震在破裂过程中受到沿断层走向存在的凹凸体的影响，导致地表位错量沿断层走向呈分段分布。干涉条纹以发震断层为中心，断层南盘的干涉条纹较北盘相对稀疏，两盘呈对称分布。干涉条纹非常密集，体现出具有较大的形变梯度特点，可能出现同震破裂分布较浅的结果。

由于InSAR观测获得的同震形变场有数据量过大，在空间上连续且相关性高的特点，反演时如计算所有的InSAR形变数据会极大地增加计算时间并影响反演结果。因此，在反演断层滑动分布之前，通过四叉树降采样法(Jonsson et al，2002)对升、降轨InSAR同震形变场进行降采样处理，从而使InSAR形变数据点既能较完整体现同震形变场的位移特征，又能约束断层面的同震滑动分布。

3 断层滑动分布反演 3.1 断层参数

本文结合余震精定位结果，通过InSAR相干性图和历史断裂构造分布构建了一条走向约105°的发震断层地表破裂迹线。依照升、降轨干涉图、位移图展现的与之对应的形变特点与余震分布区域以及在主发震断裂东侧体现出马尾状分叉特点，在主断层东南侧对应区域构建了一条分支断层。根据InSAR同震形变场和余震分布区域，将地表破裂长度分别设置为180km和30km，断层面宽度设置为30km。综合震源机制解、余震分布区域、野外实地调查结果，将断层倾角限制在70°~90°范围内。通过多次尝试计算，最终发现主断层和分支断层倾角分别为85°和68°时残差最小，故确定85°和68°作为反演主断层和分支断层的倾角。

3.2 滑动分布反演

基于弹性半空间位错模型，以InSAR和GNSS同震形变场为约束(其中GNSS同震位移参考了王阅兵等(2022)的研究成果)，采用SDM(Steepest Descent Method)方法(Wang et al，2013)对玛多地震进行滑动分布反演。本文将泊松比设置为0.25，设置滑动角范围为-20°~20°，最大滑动量10m，剪切模量设为30GPa。为获取更好的发震断层同震滑动分布，对断层进行4km×4km的格网离散化。根据观测误差定权，将GNSS、InSAR升、降轨观测值相对权重设为1、0.25、0.25。通过分析位错模型的粗糙度和拟合残差之间的折中曲线来确定光滑因子(图 3)，确定为0.10，最终得到了玛多地震断层滑动分布(图 4)。

结果显示，主断层同震滑动以左旋走滑为主，断层走向为282°或102°，沿断层走向的破裂长度约为160km，倾角为85°，倾向北，平均滑动角为-11.25°。最大滑动出现在地下1.87km处，滑动4.34m，滑动区域集中在地下0~17km深度范围内(图 4)。分支断层最大滑动量3.26m，对应深度1.74km，平均滑动角-11.39°，反演得到的矩震级为M_W7.4。整体上滑动分布模型西侧的滑动量大于断层东侧，表现出断层滑动的分段性特点。在NW方向20km外，滑动分布逐渐向地壳深部发展，并在接近鄂陵湖段时，逐渐向地壳浅部发展，造成近3m的地表破裂。但在震中SW方向未出现明显的破裂向更深部发展的趋势，推测可能是震中SE方向破裂速度远大于NW方向(Lyu et al，2022)。同时，最大滑动量出现在破裂带近地表处(深度1.74km)，表现出同震滑动破裂到地表的特征，这与野外调查结果(李智敏等，2021；Pan et al，2022)相一致。

值得注意的是，沿走向约40km，深度15km处，即野马滩大桥西侧附近(图 5)，出现了约3m左右的滑动量，表现出地震破裂由浅部向深部(约15km)发展再向鄂陵湖段浅部破裂的趋势，体现出滑动面深部几何结构比浅部更加复杂的特点。结合玛多地震三维速度结构结果，吴萍萍等(2022)认为野马滩大桥到黄河乡段断地震发生次数较少，且断层破裂段存在向北倾斜的高速体，深部存在低速区；且该段断层迹线走向发生变化，因此认为玛多地震在破裂到野马滩大桥时破裂不充分，并在野马滩大桥西侧较深部发生约3m的滑动量。

从InSAR观测的模拟值和残差值分布图来看(图 6)，通过SDM方法获取的同震形变场与实际观测值一致性较高。近断层破裂迹线处的残差较高，可能与本次地震造成了严重的近场失相干有关。通过对比GNSS水平同震位移的观测值与模拟值(图 7)，发现区域内GNSS模拟值与观测值拟合结果较好。其中，断层南盘GNSS测站拟合残差相对北盘更小，观测值与联合反演模拟值的相关系数为97.71%。

4 模拟值与理论值比较

通过使用滑动分布反演结果，运用PSCMP/PSGRN程序(Wang et al，2006)来正演计算此次玛多地震水平方向理论同震位移，验证联合InSAR和GNSS滑动分布反演的合理性，并与基于InSAR和GNSS的大地测量数据通过SDM方法反演得到的反演模拟值进行对比(图 8、图 9)。InSAR和GNSS水平方向同震位移对比(图 8)结果显示，发震断裂周围地区的水平同震位移有呈四象限分布的特点，该结果符合左旋走滑地震特征，也与震源机制解结果一致。结果表明，通过SDM方法反演得到的位移值和PSCMP/PSGRN正演得到的位移值拟合得较好。其中，InSAR数据在分支断层处存在数个点位移值相差较大的现象，推测原因是由分支断层的存在导致的位移模拟值不同。图 8中基于GNSS数据得到的水平同震位移基本呈四象限对称分布，这与王阅兵等(2022)得出的结果相同，显示此次玛多地震是典型的左旋走滑地震。

5 讨论 5.1 同震滑动分布模型

同震滑动分布模型显示主断层滑动以左旋走滑为主，并造成了长度约160km的地表破裂，破裂区域主要集中在0~17km深度范围内，野马滩大桥西侧15km深度处也有小范围破裂，同震滑动分布模型与前人研究(华俊等，2021；余鹏飞等，2022；李志才等，2021；Jin et al，2021)有所不同。华俊等(2021)的破裂分布集中于鄂陵湖段，断层东部滑动分布较少，余鹏飞等(2022)的破裂分布集中于东草阿龙段，鄂陵湖段破裂规模较小，李志才等(2021)的断层滑移分布模型受GNSS测站数量限制，仅具有两个凹凸体。Jin等(2021)的破裂分布结果显示，同震滑动分布集中在深度3~4km的地壳最上层，同震滑动分布模型最大滑动量为4.34m，该结果与最大地震滑动量与震级关系式(Wells et al，1994)估计的最大滑动量基本相同。本文由SDM方法得到的矩震级为M_W7.4，这与各机构得到的结果相比偏大，推测是受到InSAR数据中的震后余滑信息的影响。

本次玛多地震是青藏高原地区数十年来唯一一次发生在块体内部的7级以上地震，发震断裂位于昆仑山口—江错断裂附近(王未来等，2021；徐志国等，2021；Zhu et al，2021)。李智敏等(2021)和Pan等(2022)的野外地质调查显示，此次地震破裂带长约160km，西起鄂陵湖南部，东至昌麻河以东，破裂带分段特点鲜明，并基于走向变化和几何结构将地表断裂带分为4段。野外地质调查结果显示，地表破裂带的西段、中西段和东段的最大位错分别约为2.9m、1.9m、1.8m(潘家伟等，2021)。结合断层滑动分布模型和InSAR干涉图，确定沿断层走向存在4个凹凸体，最大滑动量分别为3m、2.96m、4.34m、2.74m，这与野外实地调查确定的鄂陵湖段、黄河乡段、冬草阿龙湖段和昌麻河段相对应(图 5)。两者结合表明，玛多地震破裂过程复杂，有多处破裂到地表，且表现出震源东侧比西侧滑动量大和西侧出现较深部的滑动分布的特点，体现出本次玛多地震在破裂过程中有双侧破裂且破裂不对称的特点，这与发震断裂西侧由3次子地震事件组成的2001年昆仑山口西8.1级地震有相似之处(屠泓为等，2016；许力生等，2004)。

5.2 滑动分布模型的比较

在前人的研究中，利用不同的观测资料和方法得到了2021年玛多地震的几种滑动分布模型。静态GNSS位移和InSAR干涉图的大地测量数据可用于单独或联合反演断层滑动分布(华俊等，2021；李志才等，2021；王阅兵等，2022；余鹏飞等，2022)。本文的滑动分布模型与王阅兵等(2022)、余鹏飞等(2022)的模型均涉及发震断裂东端的分支断层，但由于缺少余震分布所显示的西段分支断层相关的大地测量观测，因此滑动分布模型中未体现最西侧的分支断层。华俊等(2022)的模型采用了更简单的断层几何形态，未考虑东侧分支断层。大多数研究模型的最大滑动量为4~6m，这些模型的滑动分布也有所不同。包括本文的滑动模型在内，大多数模型中存在至少3个以上的凹凸体，但李志才等(2021)的模型仅有2个凹凸体存在。通过对这些模型进行比较，很难直接区分哪一个模型更好，但这仍然有助于理解和分析2021年玛多地震的特征。如果用更多观测数据对滑动分布模型进行约束，将有助于得到更接近真实的模型结果。

6 结论

(1) 利用欧洲航天局Sentinel-1卫星升、降轨InSAR数据，得到了2021年玛多地震同震形变场，最大LOS向形变约1.1m。基于联合InSAR和GNSS大地测量数据反演的断层滑动分布模型和GNSS水平方向同震位移场，认为玛多地震是典型的左旋走滑地震。玛多地震造成了160km左右的地表破裂，最大滑动量为4.34m，对应深度较浅(1.74km)，反演矩震级为M_W7.4。

(2) InSAR同震干涉图和滑动分布模型表明沿发震断层走向方向存在4个凹凸体，发震断层的地表破裂存在分段性特点。其中，凹凸体滑动规模较大的东草阿龙段和鄂陵湖段最大滑动量分别为4.34m和3m，滑动范围较小的黄河乡段和昌麻河段最大滑动量分别为2.90m和2.74m。实地调查显示黄河乡段地表破裂程度较小且存在向北倾斜的高速体和深部低速区，与本文得到的较深部的滑动分布一致。

(3) InSAR和GNSS水平方向同震位移对比结果显示发震断裂周围地区的水平同震位移有呈四象限分布的特点，该结果符合左旋走滑地震特征，也与震源机制解结果一致，由此证明了联合InSAR和GNSS滑动分布反演的合理性。

致谢: 感谢王海涛研究员、江在森研究员、汪荣江研究员等专家提供的宝贵修改建议，文中图件使用GMT软件绘制，在此一并感谢。