2013年9月24日在巴基斯坦西南部的俾路支省发生了一次M_W7.7的浅源强烈地震.该地区处于欧亚大陆板块南缘的突出部位(图 1)，阿拉伯板块以约3 cm·a^-1的速率北向俯冲至欧亚大陆板块之下(Reilinger et al., 2006)，强烈的构造运动形成从伊朗东南部到巴基斯坦西南部的莫克兰(Makran)增生楔；东部以左旋走滑速率约3 cm·a^-1的查曼断裂系统(包括Chaman, Ghazaband和Ornach-Nal断层)与印度板块交汇(Szeliga et al., 2012).地震活动性显示该地区地震主要分布在板块汇聚边界，2013年之前有仪器记录的7级以上的地震包括南部阿拉伯板块和欧亚大陆板块碰撞俯冲带上1945年发生的M_W8.1地震(Byrne et al., 1992)和东部查曼断裂带Ghazaband断层上1935年发生的M_W7.7地震(Lawrence et al., 1981).这次地震发生在莫克兰增生楔东端，查曼断裂带的西侧，沿着已知的推覆构造Hoshab断层形成一个约220 km的弧形地表破裂带，但这次地震并不是断层的逆冲破裂，而是一次在北倾的断层面上(USGS W-phase震源机制解显示断层倾角45°，Global CMT显示断层倾角39°)以左旋走滑为主的地震事件.这次具有罕见运动模式的大地震引起了地震学界广泛关注，虽然缺乏近场台站观测资料，但得益于高分辨的卫星遥感数据获取的同震水平位移场，许多科学家对这次地震的运动学特征进行了详细研究.Avouac等(2014)，Barnhart等(2014)分别利用Landsat-8卫星光学影像获取的水平位移结合远场波形构建了由多个不同倾角的平面断层组成的地震断层模型反演获得了这次地震的破裂过程；Jolivet等(2014)和冯光财等(2015)则构建了铲状断层模型利用卫星光学影像资料获取的同震水平位移场研究这次地震的滑动分布及运动学特征.这些研究工作揭示了这次地震的主要运动学特征并探讨了这种逆冲型形态的弧形断层上为何会发生走滑破裂为主的强烈地震，提出了多种可能来解释这次地震特殊的断层运动模式，如动态应力的作用促进破裂的弧形扩展(Avouac et al., 2014; Barnhart et al., 2015); 刚性“球-座(Ball-and-Socket)”的顺时针旋转(Barnhart et al., 2014)；同一断层在不同应力状态下可能表现为走滑或逆冲的不同运动形态(Barnhart et al., 2015)；莫克兰增生楔东部一系列与发震断层平行的弧形构造一直以稳定的速率和相同的运动模式活动(Zhou et al., 2015).总之，认识这一特殊地震的动力学成因将可能推动地震学研究的发展，而这次地震的详细运动学特征是开展这些工作的基础.本文构建了一个铲型断层模型，通过远场体波和光学影像获取的同震水平位移场的联合反演，获取了2013年巴基斯坦俾路支M_W7.7地震的详细震源破裂过程.

图 1

Fig. 1

图 1 俾路支地震破裂轨迹及构造背景 图中白色细线为查曼断裂带的三条断层，红色实线为俾路支地震的地表断层轨迹，震源机制解为USGS W-phase矩张量解.粉色圆表示了余震分布.蓝色框为获取同震水平位移的Landsat-8卫星光学影像覆盖范围.左上图给出震源区构造背景和历史地震(M>7.0)的分布. Fig. 1 Tectonic setting and surface rupture trace of Balochistan earthquake The white thin lines show the three faults of Chaman fault system. The red solid line indicates the surface rupture from the event. Black focal mechanism corresponds the source solution derived from USGS W-phase. Pink dots indicate aftershocks, sized by their magnitudes. Blue rectangle indicates range of surface co-seismic displacements derived from Landsat-8 optical image. Top left inset: Regional tectonic map and historical earthquakes (M>7.0).

这次地震的震源区人烟稀少，缺乏强震动、GPS等地震监测设施，没有可用于震源研究的近场观测资料.监测地震形变的有力工具InSAR观测卫星当时正处于更新换代期，没能获得完整的震源区形变数据.震区所处的干旱山区极利于卫星光学影像成像，通过美国2013年发射的Landsat-8卫星获取的高精度多波段影像资料，Avouac等(2014), Jolivet等(2014)，冯光财等(2015)利用Cross-correlation技术通过计算震前和震后高精度光学影像的像素偏移量，获得了不同分辨率的同震水平位移数据.本文使用Avouac等(2014)获得的240 m×240 m分辨率的同震水平形变场资料用于反演研究，数据的误差评估为东西向分量0.47 m，南北向分量0.6 m.本文对数据进行了重采样，选取了覆盖主要地震形变区均匀分布的128×128个点的矩形形变场数据进行实际反演计算，图 2b和2c分别为以图像形式给出的反演计算用东西向和南北向水平同震位移资料，图中可以清晰直观地看到这次地震同震水平位移场特征和地震破裂的地表轨迹.

图 2

Fig. 2

图 2 用于反演计算的俾路支地震同震位移场和远场地震台站分布 (a)远场地震台站分布，黑色三角为选用P波记录台站，蓝色圆为选用SH波记录台站；(b)同震位移东西向分量；(c)同震位移南北向分量. Fig. 2 The co-seismic surface displacements and global distribution of stations used for source model inversion of Balochistan earthquake (a) Distribution of stations selected for inversion of teleseismic waveforms, the black triangles and blue circles represent the stations used for P and SH wave respectively. (b) E-W offsets and (c) N-S offsets of co-seismic horizontal displacements.

高分辨的同震位移资料可以很好地约束断层滑动分布，结合地震波形资料则能获得详细准确的地震断层滑动的时间过程.本文从IRIS(地震学研究联合会)下载地震波形数据，选取其中信噪比较高并且沿方位角分布比较均匀的33个远场P波波形(震中距位于30° < Δ < 90°范围之内)和33个SH波波形资料用于震源破裂过程反演研究，台站分布见图 2a.

根据光学影像获取的同震水平形变场显示的地震断层地表破裂轨迹，本文构建了一个由10段具有不同走向的子断层组成的几何模型来模拟整个地震断层的形态，考虑可能的隐伏破裂，在断层模型两端进行扩展，即断层模型长度大于实际地表破裂轨迹.USGS确定的震源位置并不在断层面上，基于震源机制解显示的断层走向、倾角和深度信息，将震源点设置于第5子断层上.根据前人平直断层震源模型给出的不同分段的断层倾角和莫克兰增生楔基底滑脱面的深度(Avouac et al., 2014; Jolivet et al., 2014)，初定一个倾角由浅至深从75°变化至15°的铲型断层面，修改断层倾角产生不同的断层几何模型，根据反演试算资料的拟合误差不断调整，最终确定如表 1所示的断层模型几何参数.断层模型由71(走向方向)×9(倾向方向)个断层元组成，断层元长宽尺度为4 km×3 km.

表 1 Table 1

表 1 断层模型几何参数 Table 1 Geometrical parameters of fault model 子断层 1 2 3 4 5 6 7 8 9 10 走向(°) 203 206.5 220 202 208 222.2 229.5 238 240 250 沿倾向方向倾角(°) 75 75 65 55 55 60 75 75 75 75 75 75 60 55 55 60 75 75 75 75 65 65 55 50 50 55 65 65 65 65 55 55 50 45 45 50 55 55 55 55 45 45 45 45 45 45 45 45 45 45 35 35 35 35 35 35 35 35 35 35 25 25 25 25 25 25 25 25 25 25 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 长度(km) 36 16 20 36 12 24 16 16 44 64

表 1 断层模型几何参数 Table 1 Geometrical parameters of fault model

有限断层模型反演的待求解参数为每个断层元的滑动角、滑动量、破裂传播速度以及震源时间函数的上升时间(姚振兴和纪晨，1997; 王卫民等，2008).同震位移响应的计算采用基于分层介质模型的广义反射透射系数矩阵方法(谢小碧和姚振兴，1989)，震源区地壳结构模型获取自Crust2.0(Bassin et al., 2000).反演使用非线性最优化方法—模拟退火技术，采用波形相关系数与同震位移矢量差的归一化加权叠加来定义目标函数，并增加模型光滑约束，可表示为

其中，e^w和e^d分别表示波形拟合差和同震位移拟合差，p^w和p^d分别表示波形拟合和同震位移拟合的权重，O_k(t)和f_k(t)分别为第k个台站的观测和理论地震图，N_w为地震台站总数；D_k^o和D_k^s分别为第k个水平同震位移观测值和理论值矢量，N_d为水平位移场资料总数；S和λ分别代表模型光滑度参数和其权重(Ji et al., 2003).在实际计算中，以资料的总体拟合度最优为准则, 不同资料和光滑因子的权重通过多次试算调整来确定，权重系数p^w，p^d和λ的取值分别为0.8，0.2和0.2.

联合同震水平位移场和远场体波波形反演获得了2013年巴基斯坦俾路支M_W7.7地震的震源破裂过程模型(图 3).结果显示这次地震事件在一个弧形铲状断层上发生以左旋走滑为主具少量逆冲分量的断层破裂滑动，破裂由东侧查曼断裂的西端分支开始，沿弧形的Hoshab断层向莫克兰增生楔内部扩展，同时东侧的查曼断裂只发生了少量的滑动，破裂持续时间约65 s.主要滑动分布在15 km以上的浅层地壳内，破裂上升时间较短，平均约4 s，破裂扩展速度平均约2.7 km·s^-1，最大滑动量达到1320 cm.反演得到的标量地震矩M₀为6.1×10²⁰ Nm, 换算矩震级为M_W7.8，接近GlobalCMT的结果.图 4给出这次事件的地震矩释放率(震源时间函数)，较高的破裂传播速度和较低的破裂上升时间揭示了这次地震能量释放集中、迅速，表明断层面上的高应力降，使得断层西段的高滑动区很少发生余震，余震集中发生在断层东北段未发生破裂或破裂滑动很小的地段(图 1)，包括主震后第四天发生的最大一次M_W6.8余震.

图 3

Fig. 3

图 3 断层面上滑动分布及破裂扩展时间 震源位置(初始破裂点)以五角星标示，箭头表示破裂的滑动方向，滑动量大小以色标给出的不同颜色表示.等值线代表破裂前锋扩展时间.左上图为断层滑动分布的地表投影；右上图为有限断层模型三维示意图. Fig. 3 Slip distribution and times of rupture front on the fault The hypocenter is denoted by a star. Slip amplitude is shown in color and the slip direction (rake angle) is indicated with arrow. Contours show the times of rupture front in seconds. Surface projection of the fault model is shown in the upper left. The upper right is a 3D view of the finite fault model.

图 4

Fig. 4

图 4 标量地震矩与时间关系(震源时间函数) Fig. 4 Source time function for the Balochistan earthquake

本文反演给出的震源模型较好地拟合了水平向同震位移数据和远场体波波形资料(图 5，6).水平位移场的平均拟合差为48 cm，略小于观测资料平均误差，从拟合残差图 5可看到，除去零星分布的高误差点，主要的误差分布在地表破裂带附近，这是由于近断层同震位移变化剧烈，实际破裂轨迹局部变化与断层几何模型不完全吻合，同时本文没有人为去除近断层的错区数据点造成近断层局部数据点的高残差.在图 5中还可以更清晰地看到同震位移场数据存在的条带状误差，这表明利用完整的高空间分辨率同震位移场观测资料参与反演是可以克服这些资料预处理阶段很难去除的观测误差，得到真实的地震滑动分布.GlobalCMT给出的机制解显示这次地震存在很大的非双力偶分量，特殊的弧形断裂几何形态可能是造成非双力偶分量的原因，同时也不能排除存在某种垂直于断层面的运动对非双力偶分量的贡献，这可能也是本文远场波形某些方位的资料拟合不佳的原因.

图 5

Fig. 5

图 5 同震水平位移场的拟合残差分布 (a)东西向分量拟合残差；(b)南北向分量拟合残差.红色线条表示地表破裂轨迹；蓝色线框表示断层模型的地表投影位置. Fig. 5 Model residuals for the co-seismic horizontal displacements (a) and (b) show the model residuals of E-W and N-S directions, respectively. The red line indicate the surface rupture trace, and blue rectangle indicate surface projection of the fault model geometry.

图 6

Fig. 6

图 6 有限断层模型的远场P和SH波形拟合 左列为P波垂向位移理论图(黑线)与记录(红线)的对比；右列为SH波切向位移理论图(黑线)与记录(红线)的对比；并给出每个记录的方位角(左上)、震中距(左下)、振幅比(右上)和台站名(右). Fig. 6 Comparison of the observed teleseismic P (left) and SH (right) records (red) with the synthetic waveforms (black) Numbers on the left indicate the azimuth and epicentral distance of the stations, the number above each waveform is the amplitude ratio and station names also indicated.

考察目前已发表的俾路支地震震源破裂模型(Avouac et al., 2014；Barnhart et al., 2014；Jolivet et al., 2014；冯光财等, 2015)，它们在主要滑动分布特征上存在一致性，但在更细节的地震运动学特征则存在较大差异.主要的差异表现在两点：(1)最大破裂滑动是否发生在地表；(2)最大滑动量是10 m还是15 m.本文反演获取的震源模型显示最大滑动并不是发生在地表断层(或并不只限于地表断层)，这与Avouac等(2014)和Barnhart等(2014)同样使用同震位移和远场波形联合反演的结果相似而不同于Jolivet等(2014)和冯光财等(2015)只使用同震位移资料反演的结果，这说明了波形资料在震源反演中的重要约束作用.本文反演获取的这次地震破裂最大滑动量13.2 m，大于Avouac等(2014)和冯光财等(2015)给出的10 m的最大滑动量结果，更接近Barnhart等(2014)和Jolivet等(2014)给出的15 m最大滑动量.造成这一差异的原因应为模型几何参数以及参与反演的同震位移数据点密集度不同所致.

通过震源模型计算得到的地表断层两侧的垂向变化，在西段可以达到断层上盘平均抬升约2.5 m，这一结果接近Zhou等(2015)利用高分辨卫星3D影像对近断层同震垂向位移的测量数据.考察滑动性质沿地震断层从东北段到西南段的变化，我们发现逆冲分量与走滑分量的比率存在较显著的增大(表 2).这种变化可能暗示连接莫克兰增生楔和查曼转换断层的这次地震发震断层Hoshab断层长期处于这种运动模式，而不是以纯逆冲或纯走滑的形式交替活动(Zinke et al., 2014; Barnhart et al., 2015).在斜向应力条件下，弧形铲状的Hoshab断层和一系列与其平行且具有类似形态的断层始终存在一个正应力分量，可能长期积累走滑和逆冲应变，从而在一个较长复发周期下形成具有类似运动模式的强烈地震，上述推断可能是形成具有特殊运动模式的2013年巴基斯坦俾路支M_W7.7地震的原因.

表 2 Table 2

表 2 子断层逆冲走滑分量比率 Table 2 The ratio of dip-slip to strike-slip for fault segments 子断层 1 2 3 4 5 6 7 8 9 10 逆冲/走滑 0.08 0.03 0.09 0.10 0.06 0.08 0.18 0.21 0.24 0.14

表 2 子断层逆冲走滑分量比率 Table 2 The ratio of dip-slip to strike-slip for fault segments

本文利用均匀分布的同震水平位移场数据联合远场波形资料反演获得了2013年巴基斯坦俾路支M_W7.7地震详细的震源过程模型，为深入认识这次具有罕见运动模式的强烈地震事件提供了科学依据.研究工作中采用大数据量均匀分布的同震位移资料参与反演计算，结果表明这种以完整图像的形式使用位移场资料的方式可以克服资料可能存在的系统误差获得可信的结果，为快速应用越来越丰富的空间对地观测资料研究地震震源过程提供了一种新的思路.

谨此祝贺姚振兴先生从事地球物理教学科研工作60周年.感谢J.P.Avouac提供2013巴基斯坦俾路支M_W7.7地震同震位移场资料.文中部分图件使用GMT软件绘制.