地球物理学报  2018, Vol. 61 Issue (4): 1225-1237   PDF    
利用走时和波形拟合迭代反演阿拉善地壳速度结构及2015年阿拉善左旗5.8级地震震源机制
宋超, 盖增喜     
北京大学地球与空间科学学院地球物理系, 北京 100871
摘要:据中国地震台网测定,北京时间2015年4月15日15时39分,在内蒙古自治区阿拉善左旗(39.8°N,106.3°E)发生MS5.8地震,震源深度为10 km.地震发生后多家机构对其开展了研究,本文使用喜马拉雅Ⅱ期布设在南北地震带北段的台站观测数据,通过走时反演和波形拟合反演的迭代,获得了该地区地壳一维速度结构,接着利用直达P波观测与理论走时差对震中位置重定位,然后反演地震的最佳双力偶解以及震源深度,最终得到了区域速度结构、地震的三维坐标、发震时刻以及震源机制解.结果显示,此次地震发生于世界时2015年4月15日7时39分26.718s,震中(39.7663°N,106.4304°E),震源矩心深度18 km,矩震级MW5.25,节面Ⅰ走向176°,倾角85°,滑动角-180°,节面Ⅱ走向86°,倾角90°,滑动角-5°.结合该区域断裂带构造运动分析,本文认为此次地震是左旋走滑破裂,略带正断分量,断层面是节面Ⅱ,走向为NEE(近E-W)向,发震构造为震中附近的E-W向隐伏断裂.
关键词: 阿拉善地震      相位加权叠加      走时反演      速度结构      波形拟合      震源机制     
An iterative travel time inversion and waveform modeling method to determine the crust structure and the focal mechanism of 2015 Alxa Left Banner MS5.8 Earthquake
SONG Chao, GE ZengXi     
Institute of Geophysics, School of Earth and Space Sciences, Peking University, Beijing 100871, China
Abstract: According to CENC, an MS5.8 earthquake occurred in Alxa Left Banner, Inner Mongolia Autonomous Region, China (39.8°N, 106.3°E) on 15:39, April 15, 2015(GMT+8), with a source depth of 10 km. After that, several researches were carried on this earthquake. Based on these previous results, we use observation data received by stations built on the northern part of Chinese North-South Seismic Belt during the Project Himalaya Ⅱ, to obtain a better regional 1-D crustal velocity structure successfully via an iterative travel time inversion and waveform modeling method. Then we relocate the epicenter with the travel time difference between Pg observation and synthesis, and re-inverse the best double-couple focal mechanism and depth. Finally, we determine the regional velocity model, 3-D location of hypocenter, origin time and focal mechanism. Our result shows that, the epicenter occurred at 07:39:26.718 on April 15, 2015 (UTC) is located at (39.7663°N, 106.4304°E) with the centroid depth of 18 km and of MW5.25. The nodal plane Ⅰ has a strike of 176°, dip angle of 85° and slip angle of -180°, while plane Ⅱ has a strike of 86°, dip angle of 90° and slip angle of -5°. Considering the dynamic structure of local fault zone, we believe this earthquake is caused by a nearly pure left-lateral strike-slip fault with a few normal component. The nodal plane Ⅱ is the fault plane, which strikes towards NEE (nearly E-W). The seismogenic structure is likely to be an E-W striking buried fault near the epicenter.
Key words: Alxa earthquake    Phase-weighted stack    Travel time inversion    Velocity structure    Waveform misfit    Focal mechanism    
0 引言

研究表明,青藏高原东北地区的活动块体间以复杂的变形带接触,甘-青块体、阿拉善块体以及鄂尔多斯块体之间由活动断裂带、活动褶皱带和活动盆地带构成边界(许忠淮等,2000邓起东等,2003).其中,阿拉善和鄂尔多斯块体内部结构完整,强度较大,块内地震活动较少(王萍和王增光,1997邓起东等,2003张培震等,2003),但两块体边界汇聚部位却是应力集中区,活动性较强,历史上发生过1739年银川8级大地震(汪素云等,2000),自1970年有现代仪器记录起,该区及周围共发生M4.0以上地震11次,包括1976年巴音木仁MS6.2地震(许忠淮等,2000曹刚,2001韩晓明等,2015魏建民等,2017).从科学研究角度看,较大地震的震源机制解可显示断层运动性质,直接反映近期构造变形特征(许忠淮等,2000),但由于之前在上述地区观测台站少、技术受限,震源机制解等信息缺失或可靠性差;并且,尽管微震活动频繁,但中强地震发生次数少,不利于深入研究.

2015年阿拉善左旗MS5.8主震则是一个契机,它是继1976年巴音木仁地震以来该地区发生的最大地震,并且,喜马拉雅Ⅱ期在该地区布设的大量台站对该次地震序列有充足的波形记录,这些数据对研究地震本身,以及阿拉善和鄂尔多斯边界地带地震活动、现代构造变形特征等都具有重大意义.

图 1所示,地震发生在贺兰山隆起与巴彦乌拉山山前断裂之间的吉兰泰断陷带北端,鄂尔多斯块体西北边缘,是阴山东西向构造带、阿拉善弧形构造带和南北构造带等多组活动构造带的交接复合部位(吴桂桔等,2015魏建民等,2017).

图 1 震中区域地形构造示意图 黄色五角星表示台网测定的震中位置,红色实线表示活动及隐伏断层,灰色实线为活动块体边界,右上角的大尺度地图用蓝色矩形表示研究区域的位置,灰色实线仍为活动块体边界.块体及断层数据来自Zhang等(2003)吴桂桔等(2015)魏建民等(2017). Fig. 1 Tectonic structure of the epicenter region The yellow star denotes the epicenter. The red solid lines denote active or hidden faults. The grey solid line denote the boundaries of active blocks. The large-scale map on the upper left corner uses a blue rectangle to denote the relative location of studied region. The block and fault data come from Zhang et al. (2003), Wu et al. (2015) and Wei et al. (2017).

图 2所示,多家机构(GFZ,USGS,GCMT,中国地震局地球物理研究所)和文献(韩晓明等,2015郝美仙等, 2015, 2016)给出了主震震源机制解.其中,三篇文章采用了“CAP”(Cut and Paste)方法(Zhao and Helmberger, 1994Zhu and Helmberger, 1996Zhu and Ben-Zion, 2013Zhu and Zhou, 2016),都得到了走滑型震源机制,且NE-SW向均为压应力,NW-SE方向为张应力,节面参数也接近,说明即使数据量不同、假定的区域速度模型不同,该方法测定中等地震的震源机制时仍较为稳定.然而,他们反演的震源深度范围在13~30 km不等,而且对于哪个节面是断层面说法不一,这对我们认识真正的震源机制以及发震构造仍然构成困难.本文初步分析认为,产生差异的一部分原因是数据的质量和分布,更大程度上则可能源于不准确的速度模型.为了弥补较差的模型对结果的影响,CAP在互相关时允许波形的各个部分移动不同的时间以实现最大的波形拟合度.但这种情况下,各段的偏移时间较大,收敛性差,波形拟合度依然不理想,即便震源机制解可信,过大的深度误差、波形拟合残差无法合理解释观测数据.因此,要想得到真实的震源机制,首先要获取更准确的区域速度结构.相比于前人,我们最大的优势在于数据量充足,台站方位角覆盖度好(图 3),具备先反演出速度结构的条件.结合速度模型和观测波形,再对震中重定位,最后反演此次地震的震源机制和深度,而对于一维速度结构和震中的联合反演,前人已经有过成功的尝试(田玥和陈晓非,2006).最佳的模型才能最大程度解释观测数据,并得到可靠的震源机制解,才能帮助我们准确认识地震的发震构造以及震区断层活动特征.

图 2 前人和本文给出的2015年阿拉善地震震源机制解 黄色五角星为ANSS目录的震中位置,震源机制球颜色表示深度,部分参数来源于韩晓明等(2015). Fig. 2 Focal mechanisms of 2015 Alxa earthquake given by previous works and this study The yellow star denotes the epicenter location from ANSS catalog. The color of beach balls denotes the depth of source. Some parameters are cited from Han et al. (2015).
图 3 本文中所用观测台站分布及原始数据 (a)黄色五角星为ANSS目录震中位置,灰色实线为块体边界,黑色三角为观测台站;(b)原始数据波形Z(上)和T分量(下),红蓝色箭头指向对应震相. Fig. 3 Distribution of observation stations and original data used in this study (a) The yellow star denotes the epicenter location from ANSS catalog. The grey solid line denote the boundaries of active blocks. The black triangles denote the stations; (b) Z components of original data (up) and T components (down) of original data. The red and blue arrows point to the corresponding phases.
1 数据和方法

本文数据来自喜马拉雅Ⅱ期布设在南北地震带北段的宽频带地震台站和固定台站.根据数据质量和震中距,从614个有效台站中挑选了震中距在400 km以内155个台站的Z分量和T分量(图 3),并经过常规的预处理.实际反演中采用的发震时刻(UTC时间2015年4月15日7时39分27.210秒)和震中位置(39.7528°N,106.4034°E)信息来源于ANSS(Advanced National Seismic System)地震目录.

为了获得研究区域的速度结构,本文提出了一种地震走时和波形迭代反演的方法,图 4展示了该迭代求解方法流程图.从图 3中不难看出,数据中存在两类震相,一种是Pg或Sg这样的强信号,振幅大初动明显,可人工拾取走时后拟合走时曲线;另一种是Pn或者Sn这样的弱信号,由于振幅较小,无法准确辨别到时,为了得到对应的慢度和走时,我们通过倾斜叠加以增强信号.有了震相走时曲线,可以通过最小二乘反演得到速度模型,基于该模型,再通过CAP方法进行波形反演,得到震源机制、深度以及波形拟合结果,把截取理论和观测波形做互相关得到的偏移时间用来判断模型优劣并校正走时:如偏移时间较大,或出现系统性的偏差,说明模型较差,则计算该模型最佳深度下的理论直达波到时,加上偏移时间,作为新的直达波到时以取代拾取到时,再用此到时重新进行模型反演、波形反演,构成迭代;一旦偏移时间不再显著收敛减小、波形拟合度无明显提高,就退出迭代,得到最优模型.继续利用此时的直达波观测和理论走时差,对震中位置进行重定位.根据定位结果重新处理数据(震中距等相应改变)后再进行波形反演,得到最终的震源机制解.

图 4 走时与波形反演迭代方法流程示意图 Fig. 4 Flow diagram of the iterative travel time inversion and waveform modeling method

本文在处理弱震相的方法是倾斜叠加.一方面,叠加可以增强相干的弱信号;另一方面,Pn、Sn的走时曲线是直线,则用不同的速度或慢度对数据进行倾斜叠加后,在τ-pτ-v域直线会成为一个极值点,可直接读取对应慢度和到时.但由于干扰震相和噪声的存在,直接叠加出现的假象会影响判断.为此,本文采用了相位加权叠加(phase-weighted stack, PWS)方法(Schimmel and Paulssen, 1997).它主要利用与振幅无关的瞬时相位相关性对用于叠加的波形加权,微弱的相干信号因相位相同而被增强,噪声则被更好地压制.

对于直达波到时,本文采用最小二乘法拟合得到双曲线型走时曲线,结合叠加得到的首波慢度和走时,得到五个震相(Pn、pPn、Pg、Sn、Sg)的走时曲线作为反演数据.反演的目的是得到包含层厚、P波和S波波速的一维水平层状均匀速度模型.简化后的模型有四层,分别是沉积层、上地壳、下地壳以及地幔,地幔及以下为均匀半空间,震源位于上地壳,并将它分为上下两部分,两部分只有厚度差别.基于这五个震相,定义如下目标函数:

(1)

(2)

其中,例如,(TPg)mod表示模型(mod)确定的Pg震相到时,(TPg)fit表示观测数据拟合(fit)确定的Pg震相到时,(TPn)PWS表示观测数据叠加(PWS)确定的Pg震相到时,(TPn|x=0)|PWS表示观测数据叠加(PWS)确定的Pg震相震中距(x)为0处的截距时间.式(1)、(2)实际上就是震相模型走时与观测走时的残差平方和,其中EP用于反演各层P波波速及层厚,而ES用于反演各层S波波速.加上深度,待解的参数为12个,叠加得到的地幔P、S波速可以减少两个参数.对余下的参数给定合理范围,在所有的参数空间内搜索目标函数最小的一组参数,即为所求.

基于反演得到的模型,采用频率—波数(FK)方法计算格林函数(Haskell,1964Hudson,1969Wang and Herrmann, 1980Zhu and Rivera, 2002),然后采用“CAP”方法(Zhao and Helmberger, 1994; Zhu and Helmberger, 1996; Zhu and Ben-Zion, 2013Zhu and Zhou, 2016)反演震源机制,并拟合波形.用不同的深度参与反演,选取解的误差最小的深度作为最佳深度.

2 结果分析

为减少噪声干扰,数据在叠加之前经过了低通0.5 Hz滤波;为使不同震相完全分开,避免相互干扰,选择250~400 km之间的数据用于叠加(参考震中距326.0027km);用Pn震相校正极性,提高叠加效果.如图 5所示,通过叠加得到,极值点处Pn速度是8.12 km·s-1,到时是50.54 s;pPn速度7.98 km·s-1,到时56.48 s(图 5a黑色点);Sn速度4.48 km·s-1,到时85.96 s(图 5b黑色点).除此之外,我们还得到了直达波的慢度和走时(最值点)的估计值:Pg速度6 km·s-1,到时58.73 s(图 5a黑色星号);Sg速度3.69 km·s-1,到时96.39 s(图 5b黑色星号).因为直达波走时曲线在远震中距趋于直线,且慢度接近于震源所在层慢度,故会在图上形成亮点,但由于其能量大,续至震相干扰强,所以仅能为震源所在层慢度提供一定参考.作为叠加结果的检验,将对应的走时曲线标记在数据上,可以发现叠加得到的的确是真实的震相(图 5c图 5d).

图 5 相位加权叠加结果及震相检验 (a)和(c)分别为Z分量叠加结果和对应震相在波形上的走时曲线;(b)和(d)分别为T分量叠加结果和对应震相在波形上的走时曲线. Fig. 5 Phase-weighted stacking results and phase check (a) and (c) Stacking result and corresponding phase travel time curve on the waveform of Z components, respectively; (b) and (d) Stacking result and corresponding phase travel time curve on the waveform of T components, respectively.

在迭代反演的过程中,首先,人工拾取的到时经过走时曲线反演得到初始模型(表 1),用该模型参与波形反演得到的结果显示,最佳震源矩深度27 km,震级MW5.34,对应的震源机制节面Ⅰ走向267°,倾角90°,滑动角-4°,节面Ⅱ走向357°,倾角86°,滑动角180°(如图 6a所示).作为判据,此时理论波形和观测波形振幅差异大,相关系数不高,截取直达波窗作互相关得到的偏移时间出现系统性偏差(-2~0 s),说明模型不准,则用以反演模型的拾取到时误差大.然后,我们计算初始模型27 km深度处震源到各台站的理论Pg和Sg走时,分别加上对应的互相关偏移时间,作为新的Pg和Sg走时,再进行走时曲线反演得到下一个速度模型(同见表 1),该模型的波形反演结果显示,最佳震源矩深度18 km,震级MW5.25,节面Ⅰ走向86°,倾角90°,滑动角-5°,节面Ⅱ走向176°,倾角85°,滑动角180°(如图 6b所示).继续检验理论和观测波形,此时它们的振幅相当,互相关系数较上一代普遍提高,偏移时间无系统偏差,收敛性也变好.重复上述迭代过程,计算该模型18 km深度处震源到各台站的理论Pg和Sg走时并加上偏移时间,但此时在最小二乘意义下已经拟合不出一条双曲线型走时曲线,说明迭代已经收敛,因此终止迭代,将上代模型作为迭代反演的最优模型.

表 1 初始和最优模型参数 Table 1 Parameters of the initial model and best model
图 6 各模型反演的震源机制解及其相对误差随深度的变化和最佳矩心深度 (a)初始模型结果;(b)最优模型并重定位的最终结果. Fig. 6 Focal mechanisms, relative misfit errors variation with depth and the best focal depth of each model (a) Results of initial model; (b) Final results of best model after relocation.

与此同时,最后一次迭代中波形反演获得的互相关偏移时间是最优模型下Pg理论走时和观测走时之差,使用200 km内的值对震中位置重新定位,重定位后的震中坐标为(39.7663°N,106.4304°E).

使用重定位后的震中,我们对数据重新处理,再做一次波形反演,偏移时间收敛性更好,得到最终的震源机制(图 6b),震源矩心深度18 km,矩震级MW5.25,属走滑型地震,节面Ⅰ的走向176°,倾角85°,滑动角-180°,节面Ⅱ的走向86°,倾角90°,滑动角-5°.根据残差结果,得到的新发震时刻为世界时2015年4月15日7时39分26.718 s.鉴于本文反演得到的模型已经最大程度地减少了观测和理论波形间的走时差,且拟合残差小,据此得到的震源机制也比较可靠,所以我们进一步直接对比绘出两种波形,而不按照互相关偏移时间移动,以体现我们在CAP拟合方法上的提高,拟合结果如图 7所示.从图上可以看到,我们的模型和结果能较好地拟合波形前半部分和走时反演中使用的震相到时,而尾波部分存在一定的系统偏差,且后半部分拟合不好.因为本文使用波形拟合互相关的走时差来校正走时,是两者的结合,需要权衡数据的高频和低频信息,存在一定的折衷(tradeoff):高频信息对走时约束好,分辨率高;低频信息更容易拟合,互相关的准确率高,相关系数高.为了得到可信的波形反演的结果,本文在进行CAP反演的时候对Pnl和面波部分分别采用0.05~0.2 Hz和0.05~0.1 Hz的滤波频带,该频带的有效性经过了吕坚等(2008)的检验,在此频带下我们进行了最大程度的走时校正.从最终的CAP的拟合结果可以看到,加上了偏移时间的拟合效果很好,但是图 7中也显示,低频带来的一定的系统偏差是无法避免的.另外,尾波的后半部分拟合不理想,原因则来自实际速度结构的横向不均匀性等,它导致的波形的复杂散射部分,难以用一维模型简单解释.

图 7 使用最优模型下并重定位后反演得到的最佳震源位置处产生的理论波形和观测波形的完整拟合 黑色为观测波形,红色为理论波形.(a) Z分量,滤波至0.05~0.20 Hz; (b) T分量,滤波至0.05~0.10 Hz. Fig. 7 The complete waveform fitting of observation and synthesis generated at the best centroid location using best model after relocation The black lines denote the observation waveforms. The red lines denote the synthetical ones. (a) Z component, bandpass filtered by 0.05~0.2 Hz; (b) T component, bandpass filtered by 0.05~0.10 Hz.
3 讨论

本文利用初始和最优模型(表 1)来说明模型不同对结果造成的影响,图 6给出了两者波形反演结果的对比.从参数上看,两模型差异明显,最优模型比初始模型厚度大,P、S波波速也更大.从结果上看,最佳深度处的震源机制解的节面走向、倾角、滑动角都十分接近,但是深度值不同、波形拟合度不同,波形各段偏移时间大小也不同.最优模型下解的误差随深度的变化梯度更大,对深度敏感;最优深度观测波形和理论波形相关系数高,且振幅相当,互相关偏移时间绝对值较小,且更加集中,基本没有系统偏差,从这个角度来讲,本文反演的模型是比较可靠的.值得注意的是,杨宜海等(2015)也进行过先反演模型,再反演震源机制的尝试.他们采用了Herrmann等(2011)的全波形反演方法,并测试得出速度模型对结果的影响有限,对深度的影响大多在1 km以内.但是他们使用的滤波频带为芦山主震滤波频带为0.01~0.05 Hz;M≥4为0.02~0.08 Hz;M3~4为0.02~0.1 Hz,比CAP方法的滤波频带更低,因而更容易进行波形拟合;将观测和理论波形按照互相关偏移时间移动后进行对比,尽管相关系数较高,但两者仍存在着明显的正的系统偏差,说明理论波形到得快,表明了模型仍然不够准确.同为全波形反演方法,本文使用的CAP方法区别在于,它将波形分成Pnl和面波段,并进行独立的反演.它的一大优势在于对速度模型和横向的地壳变化不敏感,原因在于CAP方法对P,S和面波分别采用不同的时间位移,以减少模型不准确性带来的误差,从而保证震源机制解的稳定性.但是当模型不准时也会出现上述问题,为了解决它,我们才设计了本文的方法,即走时反演和波形反演同时进行,目的就是为了让观测和数据之间的走时差也能够显著降低,这种降低意味着模型准确性的提高,因此除了给出CAP反演的移动后的对比来表明波形拟合的良好效果外,我们也画出移动前的结果(图 7)来直观地展示这种提高.但是由于波形拟合是在低频下进行,频率越低,波长越长,想要使用低频下的互相关来校正高频意义下的震相到时,就必然要做出权衡,因此,图 7已经是最大程度校正的结果,后半部分仍旧存在一定的系统偏差是难以避免的.

图 3的台站分布可以看到,台站多分布在震中南部,鄂尔多斯块体权重略大一些,所以我们参考杨明芝等(2007)给出的鄂尔多斯块体的速度结构,并选择Crust1.0(Laske et al., 2013)在台站集中区的一点(38°N,107°E)作为代表,加上本文的最佳模型,对比了它们对于实际数据的解释情况,其中杨明芝等(2007)的模型只有P波速度且无莫霍面以下的速度,如图 8所示,为了展示迭代的效果,初始模型也加入了对比.(a)首先给出各模型的参数,可以看到,本文以及杨的模型较为简化,最佳模型和Crust1.0在深度7 km和38 km左右的界面十分接近,且莫霍面以下的速度几乎一致.而杨的模型则和Crust1.0在4 km的节面接近,但速度差别大.四个模型在10~20 km深度内的P波、S波速度差别很小,本文和Crust1.0的S波速度更为接近.它们最大的差别在于莫霍面的深度,最佳模型最大,达到了54 km左右,而杨和Crust1.0都给出了相近的42 km左右.(b)是四个模型对于直达波Pg和Sg的拟合效果,黑色星号是经过校正的到时标记,可以看到本文最佳模型的拟合效果最好,Crust1.0对Sg的拟合也比较好.最后我们以实际数据的倾斜叠加结果为背景,检验各模型在参考震中距下Pn、pPn、Sn震相的慢度和截距是否符合观测,(c)和(d)分别表示P和S波结果,可以看到,尽管慢度接近,但是最佳模型比Crust1.0更能拟合到时,而这种到时截距的差很大程度上来源于两模型莫霍面深度的不同,因此我们尽管得到的莫霍面深度更大,但更能解释数据.徐果明等(2000)得到的瑞利面波反演结果表明鄂尔多斯块体西缘以西地区,地壳厚度约为50 km;并且,童蔚蔚等(2007)得到的接收函数结果表明青藏高原东北缘和鄂尔多斯的接触过渡带接收函数震相复杂, 鄂尔多斯西南缘地壳平均厚度约为50 km;也有深地震探测表明鄂尔多斯盆地的北侧和西侧, 莫霍面深度可达45~49 km(熊小松等,2011).结合本文研究区域台站较多的南侧有山脉分布,地形高(图 3),我们认为54 km左右的莫霍面深度是一个在误差允许范围内的合理值.图 2给出了本文得到的震源机制与前人结果的对比,可以看到,本文的结果节面参数比较接近USGS和GCMT,深度则介于两者之间.

图 8 本文得到的模型与其他模型的对比结果 (a)模型参数的对比,其中杨明芝等(2007)的模型只有P波速度,且无莫霍面以下速度,虚线标记莫霍面深度;(b)各模型下的Pg走时曲线对比,黑色星号是校正后的直达波到时;(c)和(d)各模型在参考震中距的莫霍面首波震相速度和截距对倾斜叠加结果解释的对比,(c)表示P波,(d)表示S波. Fig. 8 Comparison between models inverted by this study and other models (a) Comparison of model parameters. The model in Yang et al. (2007) only contains P velocity and no velocity information below the Moho interface. The dashed lines denote the Moho depth; (b) Comparison of the Pg phase travel time curves based on each model. The black asterisks denote the corrected direct wave arrival time; (c), (d) Comparison of the interpretation to the stacking results by the Moho head wave velocities and intercepts at reference epicentral distance based on each model. (c) P wave, (d) S wave.

本文之所以要进行走时反演和波形反演的迭代,一方面,震相标记的人工误差(滤波后波长变大)、震相间干扰和噪声等会使得用于反演的走时数据不佳,即使用完美的方法也会产生不合理的模型;另一方面,走时拟合时整个参数域内可能存在多个局部极值,不同的参数范围会得到不同的结果,如果只考虑最小值,可能对应的模型拟合残差很小,但是波速比不合理,波形拟合结果也不够好,选择合适的参数范围既重要,又带有经验性.因此,波形反演不仅为了得到震源机制,也是检验模型准确性的必要手段,只有既能拟合走时,又能拟合波形的模型才是我们的目标.如果理论和观测波形拟合互相关系数低,各部分偏移时间绝对值较大且分散,或者出现系统性的偏差,则表明模型与真实结构之间存在较大差距,需要修正.而用偏移时间作为判据,并据此校正直达波走时,简单有效,也刚好利用了CAP方法对模型不敏感的特性.研究中只对200 km的台站做校正也经过了反复验证,当震中距在200 km以内时,直达波是记录中第一个震相,用波形互相关校正到时的准确性高,而在200 km之后,首波震相开始超前,经滤波后,直达波震相会和它们产生干扰,降低互相关的可靠性.

震源深度是走时反演的一个参数,同时CAP反演也得到了震源深度,但这两个值并不一致.原因在于,深度震相对于震源深度的约束比较强,走时拟合反演时,数据的滤波频带是低通0.5Hz,在此频带下,我们可以清楚地辨别pPn震相,并且倾斜叠加的结果也支持了我们的判断.但是在进行CAP反演的时候,Pnl部分滤波频带是0.05~0.2 Hz,此时已经基本看不到明显的pPn.因此在互相关的过程中,振幅微弱的pPn对整个Pnl部分的拟合影响远小于直达波和续至震相,拟合的结果更多地决定于能量大的震相,故CAP得到的深度实际上对应的是矩心深度.而走时反演的结果只依赖到时,反映的是起始破裂深度.杨宜海等(2015)中的图14给出了不同模型对深度的影响,但是并未同时给出对应的反演拟合残差、偏移时间和误差随深度的变化,并且,Herrmann等(2011)的方法拟合的主要是面波,对深部结构不够敏感,实际测试也发现了滤波频率更低(0.05~0.1 Hz)的面波部分总是更容易拟合,所以在CAP中,我们增大了体波的反演权重以更好地约束矩心深度.

理论上来说,Pn和pPn震相的慢度是相同的,都对应莫霍面以下的地幔速度,但实际倾斜叠加得到的速度有微小差异(图 5).其一,从实际数据中可以看出pPn附近其他震相以及噪声的干扰更强,即使用了相位加权叠加,仍无法忽略噪声对叠加分辨率的影响,相比而言Pn震相更加干净、突出(图 5).其二,由于pPn比Pn盲区大,且到时晚,数据范围内pPn得到叠加的台站数远小于Pn,也影响了它的叠加效果.其三,从原理上来说,经过极性校正后,叠加考察的是震相对应的波峰叠加,而实际的震相形态并不是理想的脉冲,不同震相在不同震中距具有不同的波长,随后的滤波操作更会放大这种差异,所以即使理论上慢度相同,传播路径不同的震相也不太可能叠加得到完全相同的慢度.综合考虑后,本文使用的是Pn叠加的慢度,而pPn震相只使用它和Pn的走时差.另外值得一提的是,最优模型的震源所在层P、S波速也确实接近图 5的叠加结果(参数见表 1).

本文的结果给出,地震的两个节面分别为NEE(近E-W)向和NNW(近N-S)向,最大压应力轴P为NE向,最大张应力轴T为NW向,而该区域一直受到来自甘青块体NE向的挤压作用,说明结果符合当地应力环境(国家地震局,1988邓起东等,1999高熹微等,2015).有研究表明,青藏块体东北缘活动断裂带地质构造活动以NE向挤压-逆冲和走滑-旋转(左旋走滑)为主要特征,南北地震带北段为NE向走滑类型的应力状态(王椿镛等,2015武永彩等,2016).最新的GPS流动观测资料显示,相对于稳定的华南块体,吉兰泰地区水平运动为SWW和W向,其中2009—2011年贺兰山东麓断裂带西侧平行断裂带分量表现为很强的左旋走滑,垂直分量表现为拉张变形,2011—2013年走滑分量变小(马海萍,2014).此次地震的烈度Ⅵ、Ⅶ等震线分布呈椭圆形,长轴方向NEE(杨锐等,2016魏建民等,2017).图 9b援引了韩晓明等(2015)进行重定位后的余震分布结果,从图上可以看出余震在主震东西两侧分布更集中,优势方向为NEE,接近E-W方向(韩晓明等,2015魏建民等,2017).再结合当地的地形构造(图 1),本文认为,此次地震为左旋走滑破裂,略带正断分量,断层面走向为NEE(近E-W)向.虽然位于震中东侧的磴口—本井隐伏断裂是距离震中最近的构造,但是它的走向为NE,与震源机制解确定的断裂节面差距较大,且震区布格重力异常结果并不完全支持它是此次地震的发震构造(吴桂桔等,2015).本文分析认为,真正的发震构造更像是附近的E-W走向隐伏断裂,最终结果如图 9所示.该研究区发育了多组NE、NEE和E-W走向的断层,尽管一部分已被地质调查发现,但对于断层性质仍知之甚少.就我们的研究来说,此次事件对应的断层或许意味着所有的E-W向断层都具有共同的左旋走滑性质;也有一定的可能是,当地存在新生不久的隐伏断层,或者早就存在但未被发现.这些结果也许对当地未来的地震活动性和现代断层构造结构的研究有所启发.

图 9 本文确定的此次地震的震中、震源机制以及发震断层 (a)当地构造背景下确定的重定位后的震中以及震源机制,绿色五角星表示1976年巴音木仁地震,黄色五角星表示本次地震,红色虚线表示可能的发震断层,蓝色线框表示余震重定位的区域; (b)是(a)蓝框放大的区域,表示修改自韩晓明等(2015)的余震重定位结果.圆圈代表余震的位置,圆圈填充颜色的变化表示余震相对主震的流逝时间,大小对应震级的不同. Fig. 9 The final epicenter, focal mechanism and seismogenic fault (a) The determined epicenter and focal mechanism after relocation under the background of regional tectonics. The green star denotes the 1976 Bayinmuren Earthquake. The yellow star denotes the Alxa Earthquake. The red dashed line denotes possible seismogenic fault. The blue rectangle denotes the region of relocated aftershocks. (b) The amplified region in (a) and denotes the aftershocks relocation results modified from Han et al. (2015). The circles around the star denote the aftershocks, which are coded by color, whose variation represents the elapsed time. The size of circles corresponds to different magnitudes.

我们的结果还显示,此次地震与1976年巴音木仁MS6.2地震震中位置(39.90°N,106.40°E)相当接近,并且震源机制(节面Ⅰ走向3°,倾角87°,滑动角176°;节面Ⅱ走向93°,倾角86°,滑动角3°)相似(许忠淮等,2000王晓山等,2015),可能是相近断层应力积累下的再释放(图 9a).

4 结论

本文提出了一种地震走时和波形迭代反演的方法,使用喜马拉雅Ⅱ期台站数据,获得了阿拉善及鄂尔多斯局部区域的一维地壳速度结构.基于该速度结构,反演得到2015年4月15日阿拉善左旗5.8级地震发震时刻为世界时2015年4月15日7时39分26.718s,震中位置(39.7663°N,106.4304°E),震源深度18 km,矩震级MW5.25,最佳双力偶解的节面Ⅰ走向176°,倾角85°,滑动角-180°,节面Ⅱ走向86°,倾角90°,滑动角-5°.通过研究,本文得到了如下结论:

(1) 如果数据质量较高,数量充足,台站分布均匀,速度模型对波形反演获取震源深度和波形拟合存在明显的影响,但是对于震源机制解影响不大,即使是稍差的速度模型,仍可以得到稳定的震源机制.

(2) 通过走时和波形拟合的迭代反演,可以得到并优化速度模型.进一步地,过程中得到的理论与观测走时差还可以实现震中位置的重定位,最后再用波形反演能够得到更准确的深度及震源机制,最大程度拟合观测波形.

(3) 结合当地构造活动特征,本文认为,此次地震为左旋走滑破裂,略带正断分量,断层面是节面Ⅱ,走向为NEE近E-W向,真正的发震构造是附近的E-W向隐伏断裂.由于研究区发育了多组E-W走向的断层,这可能意味着它们都具有类似的大倾角左旋走滑的性质.

致谢

本文使用的fk3.2、gCAP3D1.0程序由美国圣路易斯大学的朱露培教授提供,本文数据来自喜马拉雅Ⅱ期布设在南北地震带北段的宽频带地震台站和固定台站.文中图 2绘制的震源机制解来源于GFZ、GCMT和USGS等机构的结果.文中大多数图件的绘制使用了GMT绘图软件,在此一并提出感谢.同时,作者也感谢两位匿名审稿人和编辑对本文提出的宝贵意见,你们的修改建议使得论文增色.

参考文献
Cao G. 2001. Earthquake Research in Inner Mongolia. Beijing: Seismological Press.
Deng Q D, Cheng S P, Min W, et al. 1999. Discussion on Cenozoic tectonics and dynamics of Ordos block. Journal of Geomechanics, 5(3): 13-21.
Deng Q D, Zhang P Z, Ran Y K, et al. 2003. Active tectonics and earthquake activities in China. Earth Science Frontiers, 10(S1): 66-73.
Gao X W, Wan Y G, Huang J C, et al. 2015. Tectonic stress field analysis and static coulomb stress changes of the Ms5.8 Inner Mongolia's Alxa left banner earthquake. North China Earthquake Sciences, 33(2): 48-54.
Han X M, Liu F, Zhang F, et al. 2015. Source mechanism of the 2015 Alxa Zuoqi Ms5.8 earthquake and its relocation. Acta Seismologica Sinica, 37(6): 1059-1063.
Hao M X, Han X M, Zhang F, et al. 2015. Digital waveform recording characteristics and focal mechanism solution of Alashanzuoqi 5.8 earthquake. Plateau Earthquake Research, 27(4): 6-10.
Hao M X, Yin Z J, Ding F H, et al. 2016. The analysis on focal mechanism solution of Alashan Zuoqi M5.8 earthquake, 2015. Seismological and Geomagnetic Observation and Research, 37(4): 1-5.
Haskell N A. 1964. Total energy and energy spectral density of elastic wave radiation from propagating faults. Bulletin of the Seismological Society of America, 54(6A): 1811-1841.
Herrmann R B, Malagnini L, Munafò I. 2011. Regional moment tensors of the 2009 L'Aquila earthquake sequence. Bulletin of the Seismological Society of America, 101(3): 975-993. DOI:10.1785/0120100184
Hudson J A. 1969. A quantitative evaluation of seismic signals at teleseismic distances-I radiation from point sources. Geophysical Journal International, 18(3): 233-249. DOI:10.1111/gji.1969.18.issue-3
Lü J, Zheng Y, Ni S D, et al. 2008. Focal mechanisms and seismogenic structures of the Ms5.7 Ms4.8 Jiujiang-Ruichang earthquakes of Nov. 26, 2005. Chinese Journal of Geophysics, 51(1): 158-164.
Laske G, Masters G, Ma Z T, et al. 2013. Update on CRUST1.0-a 1-degree global model of Earth's crust.//EGU General Assembly. Vienna, Austria:EGU, EGU2013-2658.
Ma H P. 2014. Study on dynamics characteristics of GPS crustal deformation in the northern section of North-South Seismic Belt (in Chinese). Lanzhou:Lanzhou Institute of Seismology, CEA.
National Seismological Bureau. 1988. Ordos Peripheral Active Fault System. Beijing: Seismological Press.
Schimmel M, Paulssen H. 1997. Noise reduction and detection of weak, coherent signals through phase-weighted stacks. Geophysical Journal International, 130(2): 497-505. DOI:10.1111/gji.1997.130.issue-2
Tian Y, Chen X F. 2006. Simultaneous inversion of hypocenters and velocity using the quasi-Newton method and trust region method. Chinese Journal of Geophysics, 49(3): 845-854.
Tong W W, Wang L S, Mi N, et al. 2007. Study of crust and upper mantle structure in Liupanshan region by receiver functions. Science in China (Series D:Earth Science), 37(S1): 193-198.
Wang C Y, Herrmann R B. 1980. A numerical study of P-, SV-, and SH-wave generation in a plane layered medium. Bulletin of the Seismological Society of America, 70(4): 1015-1036.
Wang C Y, Yang W C, Wu J P, et al. 2015. Study on the lithospheric structure and earthquakes in North-South Tectonic Belt. Chinese Journal of Geophysics, 58(11): 3867-3901. DOI:10.6038/cjg20151101
Wang P, Wang Z G. 1997. Division of the Alxa Block and its attribution. Earthquake, 17(1): 103-112.
Wang S Y, Gao A J, Xu Z H, et al. 2000. Relocation of earthquakes in northeastern region of Qinghai-Xizang Plateau and characteristics of earthquake activity. Acta Seismologica Sinica, 22(3): 241-248.
Wang X S, Lü J, Xie Z J, et al. 2015. Focal mechanisms and tectonic stress field in the field in the North-South Seismic Belt of China. Chinese Journal of Geophysics, 58(11): 4149-4162. DOI:10.6038/cjg20151122
Wei J M, Han X M, Zhang F, et al. 2017. The sequence character of Alashan Left Banner Ms5.8 earthquake on April 15, 2015. Seismological and Geomagnetic Observation and Research, 38(1): 32-37.
Wu G J, Shen C Y, Tan H B, et al. 2015. Tectonic implication of images of Bouguer gravity Anomaly and its normalized total gradient in 2015 Alashanzuoqi Ms5.8 Area. Journal of Geodesy and Geodynamics, 35(6): 936-940.
Wu Y C, Bao C Y, Tang H T. 2016. Strain field analysis of northeast margin of Qinghai-Tibet Block based on GPS velocity fields. Journal of Geodesy and Geodynamics, 36(11): 981-984.
Xiong X S, Gao R, Zhang X Z, et al. 2011. The Moho depth of north China and northeast China revealed by seismic detection. Acta Geoscientica Sinica, 32(1): 46-56.
Xu G M, Li G P, Zhou H S, et al. 2000. The 3-D structure of shear waves in the crust and mantle of east continental China inverted by Rayleigh wave data. Chinese Journal of Geophysics, 43(3): 366-376.
Xu Z H, Wang S Y, Gao A J. 2000. Present-day tectonic movement in the northeastern margin of the Qinghai-Xizang (Tibetan) Plateau as revealed by earthquake activity. Acta Seismologica Sinica, 22(5): 472-481.
Yang M Z, Ma H Q, Liao Y H. 2007. Earthquake Seismicity and Research in Ningxia. Beijing: Seismological Press.
Yang R, Xie Q C, Yang C, et al. 2016. Analysis of strong motion records of Alax Zuoqi M5.8 earthquake. Seismological and Geomagnetic Observation and Research, 37(5): 35-40.
Yang Y H, Liang C T, Su J R. 2015. Focal mechanism inversion based on regional model inverted from receiver function and its application to the Lushan earthquake sequence. Chinese Journal of Geophysics, 58(10): 3583-3600. DOI:10.6038/cjg20151013
Zhang P Z, Deng Q D, Zhang G M, et al. 2003. Active tectonic blocks and strong earthquakes in the continent of China. Science in China Series D:Earth Sciences, 46(S2): 13-24.
Zhang P Z, Wang M, Gan W J, et al. 2003. Slip rates along major active faults from GPS measurements and constraints on contemporary continental tectonics. Earth Science Frontiers, 10(S1): 81-92.
Zhao L S, Helmberger D V. 1994. Source estimation from broadband regional seismograms. Bulletin of the Seismological Society of America, 84(1): 91-104.
Zhu L P, Helmberger D V. 1996. Advancement in source estimation techniques using broadband regional seismograms. Bulletin of the Seismological Society of America, 86(5): 1634-1641.
Zhu L P, Rivera L A. 2002. A note on the dynamic and static displacements from a point source in multilayered media. Geophysical Journal International, 148(3): 619-627. DOI:10.1046/j.1365-246X.2002.01610.x
Zhu L P, Ben-Zion Y. 2013. Parametrization of general seismic potency and moment tensors for source inversion of seismic waveform data. Geophysical Journal International, 194(2): 839-843. DOI:10.1093/gji/ggt137
Zhu L P, Zhou X F. 2016. Seismic moment tensor inversion using 3D velocity model and its application to the 2013 Lushan earthquake sequence. Physics and Chemistry of the Earth, Parts A/B/C, 95: 10-18. DOI:10.1016/j.pce.2016.01.002
曹刚. 2001. 内蒙古地震研究. 北京: 地震出版社.
邓起东, 程绍平, 闵伟, 等. 1999. 鄂尔多斯块体新生代构造活动和动力学的讨论. 地质力学学报, 5(3): 13–21.
邓起东, 张培震, 冉勇康, 等. 2003. 中国活动构造与地震活动. 地学前缘, 10(S1): 66–73.
高熹微, 万永革, 黄骥超, 等. 2015. 内蒙古阿拉善左旗MS5.8地震的构造应力场和静态库伦应力变化分析. 华北地震科学, 33(2): 48–54.
国家地震局. 1988. 鄂尔多斯周缘活动断裂系. 北京: 地震出版社.
韩晓明, 刘芳, 张帆, 等. 2015. 2015年阿拉善左旗MS5.8地震的震源机制和重新定位. 地震学报, 37(6): 1059–1063. DOI:10.11939/jass.2015.06.015
郝美仙, 韩晓明, 张帆, 等. 2015. 阿拉善左旗5.8级地震数字化波形记录特征和震源机制解. 高原地震, 27(4): 6–10.
郝美仙, 尹战军, 丁风和, 等. 2016. 2015年阿拉善左旗5.8级地震震源机制解分析. 地震地磁观测与研究, 37(4): 1–5.
吕坚, 郑勇, 倪四道, 等. 2008. 2005年11月26日九江-瑞昌Ms5.7、Ms4.8地震的震源机制解与发震构造研究. 地球物理学报, 51(1): 158–164.
马海萍. 2014. 南北地震带北段GPS地壳形变动态特征研究[硕士论文]. 兰州: 中国地震局兰州地震研究所.
田玥, 陈晓非. 2006. 利用拟牛顿法和信赖域法联合反演震中分布与一维速度结构. 地球物理学报, 49(3): 845–854.
童蔚蔚, 王良书, 米宁, 等. 2007. 利用接收函数研究六盘山地区地壳上地幔结构特征. 中国科学D辑:地球科学, 37(S1): 193–198.
王椿镛, 杨文采, 吴建平, 等. 2015. 南北构造带岩石圈结构与地震的研究. 地球物理学报, 58(11): 3867–3901. DOI:10.6038/cjg20151101
王萍, 王增光. 1997. 阿拉善活动块体的划分及归宿. 地震, 17(1): 103–112.
汪素云, 高阿甲, 许忠淮, 等. 2000. 青藏高原东北地区地震重新定位及其活动特征. 地震学报, 22(3): 241–248.
王晓山, 吕坚, 谢祖军, 等. 2015. 南北地震带震源机制解与构造应力场特征. 地球物理学报, 58(11): 4149–4162. DOI:10.6038/cjg20151122
魏建民, 韩晓明, 张帆, 等. 2017. 2015年4月15日阿拉善左旗MS5.8地震序列特征. 地震地磁观测与研究, 38(1): 32–37.
吴桂桔, 申重阳, 谈红波, 等. 2015. 2015阿拉善左旗Ms5.8地震区布格重力异常及归一化总梯度的构造意义. 大地测量与地球动力学, 35(6): 936–940.
武永彩, 保长燕, 唐红涛. 2016. 基于GPS速度场的青藏块体东北缘应变场分析. 大地测量与地球动力学, 36(11): 981–984.
熊小松, 高锐, 张兴洲, 等. 2011. 深地震探测揭示的华北及东北地区莫霍面深度. 地球学报, 32(1): 46–56.
徐果明, 李光品, 周虎顺, 等. 2000. 用瑞利面波资料反演中国大陆东部地壳上地幔横波速度的三维构造. 地球物理学报, 43(3): 366–376.
许忠淮, 汪素云, 高阿甲. 2000. 地震活动反映的青藏高原东北地区现代构造运动特征. 地震学报, 22(5): 472–481.
杨明芝, 马禾青, 廖玉华. 2007. 宁夏地震活动与研究. 北京: 地震出版社.
杨锐, 解全才, 杨程, 等. 2016. 内蒙古阿拉善左旗M5.8地震强震记录分析. 地震地磁观测与研究, 37(5): 35–40.
杨宜海, 梁春涛, 苏金蓉. 2015. 用接收函数建立区域模型的震源机制反演及其在芦山地震序列研究中的应用. 地球物理学报, 58(10): 3583–3600. DOI:10.6038/cjg20151013
张培震, 王敏, 甘卫军, 等. 2003. GPS观测的活动断裂滑动速率及其对现今大陆动力作用的制约. 地学前缘, 10(S1): 81–92.