北京时间2017年8月8日21时19分，我国四川省阿坝州九寨沟县发生了M_W6.5(M_S7.0)地震，震源位置位于33.20°N，103.82°E，深度20 km(源于中国地震台网中心(China Earthquake Networks Center, CENC)).根据中国地震局(2017)发布的烈度调查结果，此次地震的最大烈度为Ⅸ度，极震区震感强烈.截至2017年8月13日20时，九寨沟地震已造成25人死亡，525人受伤，7万多间房屋不同程度受损.同时，地震发生时间正值九寨沟旅游旺季，对九寨沟旅游业的冲击也颇为严重.主震之后伴随着密集的余震活动，据中国地震台网中心测定，截至2017年8月10日，共记录到九寨沟地震余震1700余次，其中最大的余震震级为4.8级.主震和密集发生的余震给灾区造成了严重的人员伤亡和财产损失.

地震发生后，国内外各个地震研究机构积极响应，一些学者利用远震地震资料快速给出了破裂过程反演的初步结果(张勇等，2017；张旭和许力生，2017；王卫民等，2017；罗翊菁和赵里，2017).总体特征上，多数远震反演结果较为一致地显示了此次地震的双侧破裂特征.另一方面，一些结果在破裂范围、最大滑动量和滑动区域的深度上，存在一些细节上的差异.造成这些差异的主要原因之一是震级较小.尽管九寨沟地震的面波震级较大，达到M_S7.0，但矩震级相对较小，只有M_W6.5，接近灾害性地震的震级下限，其震源尺度约20 km.通过位于数千公里外且信噪比不是很高的远震资料分辨M_W6.5地震约20 km范围内的破裂时空细节，具有较大的不确定性.

九寨沟地震发生后，中国地震局工程力学研究所“国家强震动台网中心”快速汇集了这次地震的强震记录，并根据强震记录的峰值加速度和峰值速度，确定和发布了地震的烈度分布(http://www.iem.net.cn/detail.html?id=880).作为地面运动最直接的反映，仪器烈度提供了震感强弱的第一手资料.然而，除峰值加速度和峰值速度等幅度信息外，强震波形记录中还包含了丰富的震源破裂信息.由于强震数据一般不易限幅，即使在近断层处也能记录到完整的波形信息.如果能够通过近场强震数据反演得到破裂过程，将有助于提高破裂模型的分辨能力，更好地保证其准确性和可靠性.

本文以此次九寨沟地震为例，采用基于近场数据自动反演的迭代反褶积和叠加方法(IDS)反演此次地震的强震数据，确定破裂过程.我们将讨论不同机构发布的定位结果对反演的影响，并简单展望强震数据反演在震后快速响应工作中的应用前景.

中国地震局工程力学研究所“国家强震动台网中心”在震后公布了66个台站的三分向强震动加速度记录.在反演前，我们将数据重采样为2 sps，并将加速度记录积分成速度记录.为尽可能避开强震数据中发生基线漂移和浅层结构的放大效应所对应的频带，并保证子断层点源近似条件成立，对数据进行了0.02~0.1 Hz的带通滤波.通过初次反演，对比观测数据与合成数据，发现66个台站中，有25个台站在到时和波形上可以得到较好的拟合，因此选择这25个台站的三分向记录，用于破裂模型的最终确定(图 1).

图 1

Fig. 1

图 1 2017年8月8日九寨沟地震震中附近区域构造背景 黄色和绿色五角星分别为CENC和USGS定位得到的震中位置；椭圆形区域为震后考察得到的烈度分布；红色线条为活动断层(TZF：塔藏断裂，MJF：岷江断裂，HYF：虎牙断裂)；青色小圆圈为九寨沟地震余震位置分布(房立华等，2017)1).左下角插图中的三角形为九寨沟地震的强震台站分布，其中粉色三角形为本研究选用的台站. Fig. 1 Regional tectonics of epicentral area of the 8 August 2017 Jiuzhaigou earthquake Yellow and green stars are the epicenters determined by CENC and USGS, respectively. Ellipses denote the earthquake intensities. Red lines are the active faults (TZF:Tazang fault, MJF:Minjiang fault, HYF:Huya fault). Cyan circles present the aftershocks (Fang et al., 2017)1). Triangles in inset on the bottom-left are the strong motion stations, and the pink triangles are those chosen for inversion.

1) 房立华, 吴建平等. 2017.四川九寨沟M_S7.0级地震主震及其余震序列精定位研究.

在近场强震数据反演中，由于断层到台站的距离较近，需要考虑不同子断层到台站的路径的差异，因此每一个子断层到台站处的格林函数都需要进行计算.为此，我们事先计算并构建了区域格林函数数据库.格林函数基于震中附近Shen等(2016)给出的地壳及上地幔顶部速度模型，采用Wang等(2017)的方法计算得到.

我们采用Zhang等(2014a)发展的IDS方法来反演近场强震数据.此方法综合了台网反演方法和反投影方法的优势，在数据条件和断层模型齐备的情况下，无需给定或尝试子断层破裂持续时间、最大破裂速度、时间和空间光滑权重等参数，因此能够最大限度地减少重复尝试和人为干预，以较高效率一次反演计算确定复杂破裂模型.目前，IDS方法的高效性、稳定性、以及对不同强震台网资料的适用性已在2008年汶川地震、2011年东日本地震和2015年南纳帕地震的应用中得到验证(Zhang et al., 2014a, 2015)，在针对土耳其马尔马拉(Marmara)海的潜在大地震研究中，这一方法也展现出了较强的应用价值(Diao et al., 2016).

本文采用美国地质调查局(United States Geological Survey，USGS)的体波震源机制解中沿北西方向延伸的断层节面(走向/倾角/滑动角=151°/85°/-4°)，作为实际发震断层面.在反演中，为保证所有可能的破裂行为不被遗漏，断层面的尺度需要给的足够大.鉴于九寨沟地区地壳厚度约为50 km(Bassin et al., 2000)，我们假定潜在破裂最大可抵达约一半地壳厚度，由此设断层面长为42 km，宽为24 km，分别沿走向和倾向方向离散成21×12=252个子断层，每个子断层的尺度为2 km×2 km.

作为破裂过程反演的基础震源参数之一，地震定位结果对强震数据反演的影响主要有两个方面(Zhang et al., 2014b).一是发震时刻决定了数据截取的时间；二是震源位置决定了断层的位置.不同的数据截取和断层位置有可能导致不同的破裂过程结果.其中，当震源位置的误差与破裂尺度相当时，不同的震源位置有可能导致完全不同的优势破裂方向(Zhang et al., 2014b).本文中，我们分别采用CENC和USGS的定位结果进行反演，并比较基于这两个定位结果得到的破裂模型之间的异同.

基于CENC定位结果，九寨沟地震的发震时间为北京时间2017年8月8日21:19:46.7，震源深度20 km；USGS定位结果则显示，发震时间为北京时间2017年8月8日21:19:49.5，震源深度9 km(表 1).不同定位结果中发震时间和震源深度之间存在互斥性(trade-off)，有可能给破裂过程反演造成一定干扰.比如对于远震地震资料反演，由于只能得到破裂和滑动相对于震源的位置分布，不同震源深度势必影响对滑动分布绝对位置的确定.

表 1 Table 1

表 1 不同机构给出的九寨沟地震震源时空参数 Table 1 The spatial and temporal parameters of the hypocenter released by different institutes 时间(BJT) 纬度 /(°N) 经度 /(°E) 深度 /km 机构 2017-08-08 21:19:46.7 33.200 103.820 20.0 CENC 2017-08-08 21:19:49.5 33.198 103.861 9.0 USGS 2017-08-08 21:19:53.8 33.210 103.900 14.9 gCMT 注：CENC和USGS给出的是发震时刻，gCMT给出的是矩心时间.

表 1 不同机构给出的九寨沟地震震源时空参数 Table 1 The spatial and temporal parameters of the hypocenter released by different institutes

反演中，我们将观测波形的起始时间固定在发震时刻，并把震源置于断层面上与震源深度相近的位置.在采用CENC定位结果的反演中，观测波形起始时间定为21:19:46.7，震源位于走向方向第11个、倾向方向第10个子断层处(震源深度约为20 km)；而在采用USGS定位结果的反演中，观测波形起始时间定为21:19:49.5，震源位于走向方向第11个、倾向方向第5个子断层处(震源深度约为9 km).

采用两个定位结果得到的破裂模型见图 2—3(以下分别简称为CENC模型和USGS模型).IDS方法经过22次和20次迭代后，地震矩分别收敛于6.89×10¹⁸ Nm和6.87×10¹⁸ Nm，矩震级皆为M_W6.5，在矩震级规模上与其他研究机构所得结果一致(张勇等，2017；张旭和许力生, 2017; 张旭等，2017；王卫民等，2017；罗翊菁和赵里，2017).

图 2

Fig. 2

图 2 分别采用CENC和USGS定位结果反演得到的破裂模型：(a) CENC结果；(b) USGS结果 从上到下依次为震源时间函数、断层面上静态滑动量分布、以及子断层震源时间函数.其中子断层震源时间函数图像中，横轴宽度为12 s，纵轴高度为0.4 m·s-1. Fig. 2 Rupture models obtained with the hypocenter location results of (a) CENC and (b) USGS From top to bottom are the source time function, the static slip distribution, and the sub-fault source time functions. The time widow of each sub-fault source time function is 12 s, and the peak slip rate is 0.4 m·s-1.

图 3

Fig. 3

图 3 断层面上破裂的时空分布 (a) CENC模型中滑动速率随时间的变化；(b) CENC模型中累积滑动量随时间的变化；(c) USGS模型中滑动速率随时间的变化；(d) USGS模型中累积滑动量随时间的变化. Fig. 3 Spatial and temporal variation of earthquake fault ruptures (a) Temporal variations of slip-rate in CENC model; (b) Temporal variations of cumulated slip in CENC model; (c) Temporal variations of slip-rate in USGS model; (b) Temporal variations of cumulated slip in USGS model.

图 2显示了反演得到的震源时间函数、断层面上滑动量分布和子断层震源时间函数.从整个地震的震源时间函数结果看，地震的时间过程比较简单，没有显示出多次子事件的破裂特征.相对于USGS模型，CENC模型的震源时间函数向后延迟了约2.5 s，接近但稍区别于发震时刻的差别(2.8 s)，这一区别可能是由于震中位置的差异(约3.8 km)所导致.整体上看，地震矩的释放过程在破裂开始后逐渐上升，在4~6 s达到峰值，与gCMT解中的矩心时间一致(表 1).所有地震矩在破裂开始后8~10 s内基本释放完毕.尽管采用了截然不同的震源深度(20 km和9 km)，CENC模型和USGS模型都显示了非常相似的静态滑动分布特征(图 2)，可见近场强震资料反演确定的滑动分布并不显著依赖于震源位置.滑动整体上呈椭圆型分布，主要集中在震中西北15 km至震中东南10 km的范围内，为不对称的双侧破裂模式.这一西北略占优势的双侧破裂模式与余震分布特征一致(图 1)，也与InSAR资料初步反演结果接近(单新建等，2017).深度方向上，滑动分布位于0~10 km处，其中在接近地表处(1 km)的滑动量约为0.3~0.4 m，这一地表破裂特征有待进一步研究确认.但另一方面，由于震中附近区域山高林密，人烟稀少，且可能存在一定厚度的土层，浅处基岩(深度1 km)上发生的0.3~0.4 m的破裂也不容易被野外考察发现.最大滑动量约为1.1 m，位于震源上方约5 km处.根据滑动分布结果计算得到的地震矩释放的矩心深度约为6 km.如果不考虑2.5 s左右的时间偏移，两个破裂模型的子断层震源时间函数特征也非常相近，滑动量大的子断层的峰值速率(震源时间函数幅度)大，持续时间长，形状简单，可能对应于滑动凹凸体的大规模破裂；而滑动量小的子断层的峰值速率小，持续时间短，形状复杂，可能对应于破裂遭受阻碍而发生的缓慢破裂.

图 3显示了滑动速率和滑动量的时空分布.从CENC模型的滑动时空分布图像可以看出，破裂的扩展和传播由3个阶段组成.第一阶段发生在地震发生后2 s，破裂主要沿断层面向上倾方向(地表方向)传播，对应的破裂速度约5.2 km·s^-1，明显偏大.需要说明的是，这一过程与给定的震源深度有一定联系，具有较强的不确定性.第二阶段发生在破裂开始后2~6 s，破裂在深度5 km为圆心，向四周扩展，平均破裂速度大约为2.5 km·s^-1.第三阶段发生在6 s之后，破裂开始朝西北方向扩展，但地震矩规模略有减弱.相比之下，USGS模型的破裂时空演化只包含2个阶段，其中第一阶段发生在破裂开始后0~4 s，第二阶段发生在4 s以后，分别与CENC模型的第二和第三阶段相对应.可见，CENC模型中第一阶段破裂由深部沿断层上倾方向向浅部传播的过程是由于给定20 km震源深度造成的，一旦震源深度变浅，这一破裂特征随即消失，因此具有较大的不确定性.

观测波形与合成波形的比较见图 4，用CENC定位结果和USGS定位结果得到合成波形几乎完全相同，反演残差皆约为0.5，属于尚可接受的范围(Zhang et al., 2014a).造成波形残差难以进一步减小的原因是波形中存在较多的误差干扰因素.由于绝大多数强震台站布置在较为松软的土层上，导致地震过程中台站下方发生复杂而明显的倾斜，导致较强的基线漂移，形成低频误差；另一方面，土层对高频地震波形具有明显的放大作用，且由于地下结构的三维非均匀性和各向异性很难在格林函数计算中完全考虑，观测数据和格林函数在较高频段也存在一定的误差.正是由于观测数据中的低频和高频误差，以及格林函数中的高频误差，限制了波形拟合的进一步提高.然而从反演结果来看，尽管我们给定了一个很大的断层面，但得到的破裂分布非常集中、且破裂传播具有很强的规律性，表明IDS方法具有很强的抗误差干扰的能力.

图 4

Fig. 4

图 4 九寨沟地震强震数据观测波形(黑色)与合成波形(红色和绿色)的比较 红色波形为CENC模型结果，绿色波形为USGS模型结果. Fig. 4 Comparisons of the observed (black) and synthetic (red and green) waves Red and green waves are for CENC and USGS models, respectively.

九寨沟地震之后，许多研究团队确定并发布了这次地震的破裂过程模型.其中大多为采用震后美国地震学研究联合会(Incorporated Research Institutions for Seismology, IRIS)数据中心提供的远震数据反演得到的远震模型.图 5比较了不同研究组确定破裂模型与本文强震模型的滑动分布.远震模型中，最显著的滑动区域大多集中在震源附近，且基本以震源为中心分布在断层面上.从走向上看，这一滑动分布特征表明破裂以双侧破裂为主，不存在明显的优势方向，这与本文强震数据反演结果基本一致.从深度方向看，远震模型的滑动分布依赖于给定的震源深度，由于多数震源深度给定在10 km左右，因此滑动分布也大多集中在此深度.相比之下，强震资料反演不依赖给定的震源深度，其滑动在深度方向上的分布非常稳定，集中在0~10 km范围内，据此我们认为强震模型结果，尤其是其显示的破裂深度分布，具有更高的可信度.

图 5

Fig. 5

图 5 不同研究组得到的断层面上滑动量分布 (a), (b)为本文强震数据反演结果；(c), (d), (e), (f)为远震数据反演结果；(g)为远震数据和InSAR数据联合反演结果. Fig. 5 Static slip distributions determined by different groups (a) and (b) show the strong-motion models determined in this study; (c), (d), (e) and (f) are teleseismic rupture models; (g) shows the joint model estimated from teleseismic and InSAR data.

为进一步检验本文得到的破裂模型的可靠性，我们分别基于CENC和USGS模型，分别正演了如图 6a所示的17个台站处的垂直向远震P波波形.其中远震格林函数基于AK135全球大陆速度结构模型，采用Wang等(2017)的方法计算得到.在正演中，对观测数据和格林函数进行了与强震数据反演相同的带通滤波(0.02~0.1 Hz).图 6b为观测波形与正演合成波形的比较，观测波形和合成波形拟合得相当好，且CENC模型和USGS模型正演得到的合成波形几乎完全相同.由此可见，本文得到的破裂模型尽管通过反演强震数据得到，但也可以很好地解释远震数据，从另一个侧面证明了模型可能具有较高的可靠性.

图 6

Fig. 6

图 6 九寨沟地震远震数据观测波形(黑色)和合成波形(红色和绿色)的比较 (a)蓝色三角形为远震台站分布，绿色五角星为CENC定位得到的震中位置；(b)红色为CENC模型正演结果，绿色为USGS模型正演结果. Fig. 6 Comparisons of the observed (black) and synthetic (red and green) waves (a) Blue triangles are the teleseismic stations, green star is the epicenter determined by CENC; (b) Red and green waves are synthesized based on the CENC and USGS rupture models, respectively.

一般地，近场强震资料和远震资料反演在波形数据拟合上存在差别.远震资料主要拟合波形，而不考虑拟合到时，因此往往会将数据在时间上进行适当偏移以实现波形拟合的进一步提升；强震资料反演一般既拟合波形，也拟合到时.由于拟合了到时，强震资料反演对破裂的绝对位置具有一定的定位能力(Zhang et al., 2014b).相比之下，远震资料反演不具备对破裂绝对位置的分辨能力，只能确定破裂相对于震源(初始破裂点)的相对位置.因此，相对于远震资料反演，强震资料反演对震源位置的依赖性较弱.本文开展的强震数据反演工作也证实了这一点，尽管采用了不同的地震定位结果进行反演，但得到的破裂模型在绝大部分特征上都非常一致.特别地，基于不同深度反演得到的滑动分布都集中在0~10 km范围内.这实际上体现了强震数据反演在震后应急工作中的一大优势.震后较短时间内，快速定位结果得到的发震时间和震源位置可能具有较大的不确定性.而强震数据反演可以在很大程度上无视这种不确定性，得到可靠的破裂分布和滑动位置，为极震区位置估计提供可靠的参数.

还需要提及的是，作为破裂过程反演的基础震源参数之一，震源机制结果对破裂过程的影响也同样重要.在这方面，Zhang等(2015)曾经做过测试，发现断层面参数的不确定性在±15°范围内时，不会对破裂模型的主要特征产生显著影响.因此作为面向震后快速响应的工作之一，强震数据反演具有较强的可行性.

除对破裂绝对位置的分辨能力之外，强震数据相对于远震宽频带数据在破裂模型确定的快速性上还具有其他一些优势.首先，宽频带资料由于动态范围的限制，在近距离上容易限幅，只能在较远距离上才能获取到可用的资料用于震源反演.由于距离较远，地震波传播需要一定时间.即使采用P波资料进行反演，在远震距离上，至少也需要6 min(对应于30°震中距，浅源地震的情况)的传播时间.强震资料由于分布在较近距离上，地震波传播时间相比之下要短得多，这有助于更早开展震源快速反演，为震后响应争取更多时间.

强震数据的一大劣势是台站效应比较复杂.如前文所述，强震数据的基线漂移和高频放大都容易引入误差，在一定程度上降低了数据质量.这也是为何有66个强震台站记录，但只有25个台站被选择用于反演的原因.此外，我们在反演中还采用了0.02~0.1 Hz的带通滤波，这一滤波频带实际上是在提取有效震源信息和克服噪声干扰之间寻求平衡的结果.反演结果表明，这一频带范围可以基本保证有效震源信息的提取，并得到可接受的波形拟合.

由于IDS方法无需尝试和给定一系列的反演参数，在很大程度上改进了传统工作中多次试错(trial-and-error)的反演模式，大大提高了反演的自动化程度，在完成数据筛选并确定好断层面参数的情况下，仅需一次计算即可确定破裂模型.本文针对九寨沟地震的反演应用中，在普通的微机上(CPU: interl i5-6400)22次迭代共耗时20 s.如果进一步使用高性能计算机进行并行计算，那么IDS方法的耗时可以进一步缩短至忽略不计.

本文采用强震数据在九寨沟地震中的应用研究表明，强震数据反演在很大程度上不依赖于地震定位结果，可以高效稳定地确定破裂过程，尤其是滑动区域的绝对位置.快速稳定的破裂模型可以为强地面运动模拟和烈度估计提供可靠的输入参数，对震后快速开展的应急响应工作具有重要价值.随着强震观测技术的日益完善和强震数据质量的不断提高，未来在强震监测重点区域，反演强震数据确定破裂模型的工作有可能成为震后最主要的应急内容之一.

采用的强震数据由中国地震局工程力学研究所“国家强震动台网中心”提供，论文的一些图件由GMT(Generic Mapping Tools)绘制.作者对此一并表示感谢.