中国地震  2016, Vol. 32 Issue (2): 426-437
用于龙门山断裂南段深部介质探测的精密可控震源
崔仁胜1,2, 孙安辉1,2, 王洪体1,2, 刘宁1,2, 陈阳1,2, 林湛1,2     
1. 中国地震局地震预测重点实验室, 北京市复兴路63号 100036;
2. 中国地震局地震预测研究所, 北京市复兴路63号 100036
摘要:汶川地震后在龙门山断裂带南段建立了精密可控震源发射台,结合布设的28套流动地震仪的信号记录,尝试监测地下介质的动态应力变化。为了提取记录中的主动源发射的信号,本文采用了匹配滤波方法来实现弱信噪比条件下的信号检测并进行了波形变换。在提取出每个周期的主动源信号后,利用信号叠加,分别获得了1、5、10天时间的波形,有效地提高了记录波形的信噪比。然后通过互相关检测精确得出地震波的走时变化,并分析了计算误差。最后,针对实验期间发生的绵竹MS5.6地震进行分析,获得了地震前后的走时剖面和信号变化,结果表明,距离较远的龙门山断裂北段地震未能对南段实验区域的地下介质波速产生影响。
关键词龙门山断裂带    精密可控震源    走时变化    信号提取    主动监测    
The controlled accurate seismic source for exploration of deep medium below south part of the Longmenshan fault zone
Cui Rensheng1,2, Sun Anhui1,2, Wang Hongti1,2, Liu Ning1,2, Chen Yang1,2, Lin Zhan1,2     
1. Key Laboratory of Earthquake Prediction, CEA, Beijing 100036, China;
2. Institute of Earthquake Science, CEA, Beijing 100036, China
Abstract: In order to monitor temporal stress variation of underground medium, we built a signal transmitters based on the Controlled Accurate Seismic Source(CASS)near south part of Longmenshan Fault zone after the Wenchuan MS8.0 earthquake in 2008. 28 seismometers were set up to record the elastic waves emitted from the CASS. The source signal was detected in the records with low signal to noise ratio of the seismic stations, and transformed to seismic wavelet by matched filter technique. Then the SNR of the seismic wavelet was improved by stacking technique, and seismic wavelet of every day, every 5 days and every 10 days were obtained. In addition, the variation of seismic travel time was extracted by cross-correlation and the error range was assessed. Using CASS system is feasible for monitor variation of seismic travel time from the result. Finally, we compared and analyzed the variation of seismic travel time before and after the Mianzhu MS5.6 earthquake in 2009. The result showed that the earthquake with long distance had no effort on the seismic velocity of underground medium below our experimental region.
Key words: Longmenshan Fault zone     CASS     Variation of travel time     Signal extraction     Active monitoring    
0 引言

2008年5月12日汶川8.0级地震是一场巨大灾难,给震区人民带来了巨大的生命财产损失,因此,开展地震发生机理的研究具有重要的意义。众所周知,地震是在构造应力作用下由地下介质发生破裂而引起的。地下介质的应力状态与地震的孕育和发生密切相关,研究地球深部应力随时间的变化是理解地震发生机理的有效途径。然而直接测量地下介质的应力变化比较困难,但穿过地下介质的弹性波的波速变化可以有效反映地下介质因构造应力变化所导致的改变。本宗充等(2005)通过重复爆破方法证明了大震孕育区内的应力变化可以通过观测地震波速的变化得到监测;刘志坤等(2010)利用背景噪声互相关技术研究发现,汶川地震震源区地震波速度最大降幅达0.4%;其他学者通过重复地震、人工震源等手段也监测到地下介质的波速变化(Poupient et al,1984Kunitomo et al,2004Wang et al,2006Ryoya et al,2007Niu et al,2008Wang et al,2008)。因此,测量地震波的走时可以实现对地下介质波速变化的测量,进而通过对地下介质应力变化的动态监测来达到探索地震孕育环境的目的。

除天然地震和背景噪声等观测方法以外,人工震源也是探测地下介质结构的有力手段,由其可以获得不同尺度的地球深部图像(陈颙等,20052007)。本文采用的人工震源是精密可控震源,属机电可控的旋转式人工震源。国际上与之相类似的有俄罗斯的可控震源和日本精密控制常时运行震源(ACROSS)。它们各有特点,其中,俄罗斯的可控震源是垂直作用力的线性震源,以作用力大、传播距离远见长;日本的ACROSS为水平旋转震源,其显著特点是精密控制、长时间运行。利用这2种人工震源已分别在俄罗斯和日本开展了大量的实验(Alekseev et al,2005Ryoya et al,2007)。中国的机械旋转震源研究始于2002年,由最初的水平旋转震源,到输出水平线性作用力的第2代可控震源,再到输出垂直线性作用力的第3代可控震源,逐步改善了对地耦合水平;后结合俄罗斯和日本可控震源的优点,在大吨位、垂直输出力的同时,还实现了发射信号的高精度控制,并在此基础上研制出了第4代大吨位精密可控震源,即由北京港震机电技术有限公司生产的40t精密可控震源(廖成旺等, 2003, 王洪体等,2009)。

用于深部结构探测的40t精密可控震源是一种大功率机电可控的旋转式人工震源,具有绿色、环保、发射频率可控制、信号高度重复和可连续长时间运行等优点,是地下介质结构探测的有力工具。汶川大地震发生后,中国地震局地震预测研究所于2008年9月~2009年10月间组织实施了龙门山断裂南段的精密可控震源主动观测实验,尝试开展对该区域固定地震波射线路径走时变化的监视与研究。

1 实验概况

对GPS观测资料和汶川地震库仑应力触发(唐文清等,2004安其美等,2004)等的分析显示,龙门山断裂带南段应力应变水平上升。安其美等(2004)对龙门山断裂带两侧的近地表(100~500m深的钻孔)水压致裂地应力测量结果表明,南段的应力积累已经达到临界值。故我们认为,在汶川震区的南端存在的历史地震空区是未来可能发生地震的危险区域,为此,选定了位于汶川震区南端的雅安、宝兴一带(图 1)开展实验研究。野外工作共分为2期进行,其中第1期为2008年9月~2009年1月,第2期为2009年5~10月,本文处理分析的是第2期实验数据。本次实验研究区域位于龙门山断裂南端,精密可控震源信号发射台位于四川省芦山县太平镇,共布设流动测线2条。如图 1所示,A测线由雅安市雨城区望鱼镇开始向西北方向展布长度约120km,架设了23套流动数字地震仪;B测线由5个台组成,分布在可控震源发射台两侧。实验使用40t大功率精密可控震源,用于野外观测的地震计由10台FSS-3短周期地震计和18台BBVS-1宽频带地震计组成,采集器为EDAS-24 IP型数据采集器。实验中利用精密可控震源信号发射台每天连续6hr发射线性扫频信号,通过流动地震计记录观测数据。2009年5~10月开展了为期6个月的实际观测,获得了大量有效的可控震源观测数据。

图 1 精密可控震源发射台及测点分布 五角星为精密可控震源信号发射台;三角形为观测点
2 资料处理

精密可控震源发射的是一种连续信号(调制正弦波,频率为4~10Hz),观测台站记录数据的信噪比低于0dB,与脉冲震源发出的脉冲信号有所不同,不能由其直接分析和提取地下介质信息。而目前解释介质特性的技术,都是基于震源信号为脉冲信号的波形,因此,需要对观测数据进行波形变换和信号叠加(王洪体等,2009)。另外,由于精密可控震源发射连续的线性调频信号,故可以使用数字信号处理技术来提取有效信号。

实验期间,精密可控震源发射信号的频率为4~10Hz,信号激发自每天22点到次日凌晨4点共计6hr,每小时1个周期,共6个周期。每个周期内有5次升频-降频过程,具体变频模式如图 2所示。通过记录近场可控震源信号并计算其相关性发现,精密可控震源发射信号相关系数高达99.995%以上,发射信号的高度重复保证了测量相同路径地震波走时变化的可靠性。

图 2 发射信号模式
2.1 波形变换

波形变换的目的是把湮没在噪声中的可控震源信号提取出来,并转换成可以进一步提取震相信息的地震子波(王洪体等,2009崔仁胜等,2011)。王洪体等(2009)杨微等(2010)崔仁胜等(2011)刘希康等(2013)张正帅等(2015)分别用匹配滤波、反褶积、加权匹配滤波、Wigner-Hough变换和时变窄带滤波等提取精密可控震源信号。本文采用了匹配滤波方法来实现弱信噪比条件下的信号检测并进行波形变换。其原理如下:

s(t)为可控震源发射信号;x(t)为台站记录信号;n(t)为台站噪声;r(t)为信号传播路径的传递函数,则

$ x\left(t \right) = r\left(t \right) * s\left(t \right) + n\left(t \right) $ (1)

其中,*代表褶积。根据匹配滤波理论存在这样一个滤波器,它使得台站记录信号x(t)通过滤波后输出信号y(t)达到最大信噪比,即

$ y\left( t \right) = \int\limits_{ - \infty }^{ + \infty } {x\left( {t + \tau } \right)s\left( \tau \right){\rm{d}}\tau } $ (2)

这等价于计算x(t)和s(t)的互相关。根据匹配滤波理论,通过式(2) 计算出的滤波输出的最大信噪比可达到

$ {\rho _{\max }} = \frac{E}{{{N_0}/2}} $ (3)

其中,N0为噪声的平均功率谱值;E为震源发射信号的能量。精密可控震源发射信号的能量与信号激发时间正相关,因此能通过增加发射信号的持续时间来增加匹配滤波输出的信噪比。

由式(2) 可见,台站记录的匹配滤波输出就是传递函数与发射信号自相关波形的褶积。单个传递函数波形形态接近于发射信号的自相关波形,匹配滤波的输出看起来就是发射信号的自相关波形组成的序列,序列中的每一个子波就是震源信号的自相关波形,它对应着1个ρmax,出现的时间就是对应传递函数的时间延迟,即不同震相的观测走时。

图 3为实验中2009年7月21日0时布设在信号发射台附近的近场观测记录的1个周期(1hr)的源信号(图 3(a)),以及距发射台18.3km的lmsb05台的记录波形(图 3(b))。对源信号和记录波形4~10Hz滤波预处理后再进行匹配滤波得到转换波形(图 3(c)),即加窗传递函数,其峰值对应的是Pg等多个震相。

图 3 源信号(a)、lmsb05台观测记录波形(b)及转换波形(c)的对比
2.2 波形叠加

由于在部分观测时段可能受到强烈的短时环境噪声的影响,部分台站记录信号的单个周期的转换波形信噪比较低,因此本文通过波形叠加提高信噪比,如下所示

$ y\left(t \right) = \frac{1}{n}\sum\limits_{i = 1}^n {{x_i}\left(t \right)} \;\;\;i = 1, 2, \cdots, 6 $ (4)

其中,xi为每个周期的转换波形。如式(4) 所示,叠加就是把每一个测点每次计算得到的波形严格按照同时间点累加起来,这主要是利用了信号的相关性和噪声的不相关性,从理论上讲,信噪比的改善与叠加次数的平方根正相关。在实验中,精密可控震源每天运行6hr(周期),分别对6个周期的观测记录进行波形转换,得到1天的6条转换波形。通过叠加获得日记录后,分别叠加5天日记录和10天日记录来获得以5、10天为等时间间隔的波形时间序列。

图 4为第2期实验中2009年5月17日~8月11日lmsb05台站的日记录波形序列(图 4(a))、5天记录波形序列(图 4(b))和10天记录波形序列(图 4(c))。由图 4可见,波形序列各震相波组清晰明确,到时稳定,叠加次数越多,震相信噪比越高,越有利于同相位震相追踪及走时提取。

图 4 lmsb05测点的1天(a)、5天(b)及10天(c)记录波形序列
3 走时变化提取和分析 3.1 走时变化提取方法

对于精密可控震源信号发射台和同一个观测台站,假设信号传播的路径固定,可控震源信号激发的时刻为ts,某一震相到达接收点的时刻为te,则沿着固定地震波的走时τ

$ \tau = {t_{\rm{e}}} - {t_{\rm{s}}} $ (5)

图 4可见,震相走时τ可以在转换波形中直接拾取读出,而不同时间2次激发的信号传递到同一台站时同一震相的走时分别为τ1τ2,则地震波走时变化为τ1τ2之差。精密可控震源信号的2次激发时间已知,观测波形经处理后可提取震相到达时间,对2条波形进行互相关检测则可精确得出地震波走时的变化(罗桂纯等,2008王伟涛,2009)。

3.2 走时变化计算及误差分析

由于精密可控震源发射信号的高重复性可以减少可控震源对走时变化的影响,且互相关方法能够有效排除干扰和不确定因素,故可对地震波走时变化进行精确测量。对某一接收点的某个人震相,取2个时间点同一震相波组进行信号重采样,计算2个震相波组之间互相关,通过检测互相关波形中最大值出现的位置,则可确定2个波组的走时差。

计算中首先对提取的波组数据线性内插,得到200000sps的波形数据,再对插值后的波形进行零相位低通滤波,得到重采样的200000sps的波数据。分别选取时间窗口1s的重采样后的参考波形和当前波形进行互相关检测,得到相对参考震相走时的变化,参考波形由所有波形叠加平均得到。按照上述的计算过程,走时变换计算分辨率为5μs,通过限定允许的最大走时变化为100ms所得结果如图 5所示。如前所述,在得到3组等间隔波形时间序列的基础上,通过互相关检测参考波形和当前波形的震相窗口,计算得到2009年5月17日~8月11日间lmsb05测点的1、5、10天尺度的走时变化曲线,其中1天尺度的走时变化范围为-5~2ms,5天尺度和10天尺度的走时变化范围为-2~2ms。由图 5可见,同一组震相在不同尺度下的走时变化趋势一致,但小尺度能将变化刻画得更细微些,时间分辨率更高。

图 5 lmsb05台(震中距18.30km)走时变化 图(a)、(b)、(c)上图均为参考波形,箭头处为检测窗起始位置,窗长1s;下图均为提取的走时变化曲线

互相关检测的误差下限可由下式所示的Cramer-Rao Lower Bound法则(Silver et al,2007)计算

$ {\sigma _\tau } \ge \sqrt {\frac{3}{{2f_0^3{\pi ^2}T\left({{B^2} + 12B} \right)}}\left[ {\frac{1}{{{\rho ^2}}}{{\left({1 + \frac{1}{{{\rm{SN}}{{\rm{R}}^2}}}} \right)}^2} - 1} \right]} $ (6)

其中,f0为震源的中心频率;B为频带宽度比;ρ为波形的相关系数;SNR为信噪比;T为所选相关窗口的长度。实验中,f0≈7Hz,ρ≈1,SNR≈20,T=1s,B≈2.5,代入式(6) 计算可得误差约为0.2ms。

3.3 龙门山断裂北段地震对地下介质波速的影响

据中国地震台网中心测定,2009年6月30日02时03分四川省绵竹市(31.4°N,104.1°E)发生MS5.6地震。该地震震中距精密可控震源信号发射台约170km(图 6),本文尝试计算B测线的走时变化,以分析该地震对龙门山断裂南段地下介质波速的影响。

图 6 绵竹地震位置示意图

我们对该地震前后5天的数据通过匹配滤波和叠加处理并进行4~10Hz的滤波获得B测线走时剖面(图 7)。图 8为各测点信号局部放大,除lms07台因信噪比较低引起的差异较大外,其他台站地震前后波形相位吻合较好,lmsb01~lms09台P波初至估算直达P波的视速度为5.06km/s,根据纵、横波波速比为1.73,利用P波波速推算S波波速为2.93km/s。lmsb02~lmsb04台估算直达P波的视速度为5.53km/s,推算S波波速为3.2km/s。为获得精确的走时延迟,将地震前、后的波形数据进行重采样(200000sps),互相关检测发现各观测点(除lms07台)S波走时(表 1)相对变化较小,其中,lms07台因信噪比较低而导致延迟较大。由此可见,绵竹MS5.6地震对龙门山断裂南段各观测点地震前、后走时延迟的影响不大。

图 7 绵竹地震前后B测线波形 蓝色为地震前; 红色为地震后

图 8 绵竹地震前后B测线各观测点信号局部放大 蓝色为地震前;红色为地震后

表 1 B测线各观测点绵竹地震前后走时延迟

研究表明,地震波的走时变化与大气压力之间有很好的相关性(Wang et al,20062008),此外还受温度、降雨等因素的影响。考虑到数据采集器的GPS服务时间误差一般为μs量级,比走时变化的ms量级小2~3个数量级,另外通过长时间的实验记录也可消除GPS授时异常产生的影响(杨微等,2009)。

汶川地震后杨微等(2009)在龙门山断裂带的东北缘布设了主动源动态监测系统,在测线附近观测到绵竹MS5.6地震前、后穿过断裂带的直达S波走时发生了时延为5~9ms的微弱变化。相比上述测线,我们的观测系统位于龙门山断裂南段,与绵竹地震震中平均距离约为170km,方位角大致与龙门山断裂带平行。本文通过分析实验数据未能观测到明显的地震波走时变化。故推测认为,一方面可能是由于该地震产生的震源破裂区未及龙门山断裂带南段,未能对我们观测系统所在区域的地下介质产生足够的影响;另一方面也可能实验所在区域受当地环境(如温度、降雨等)的影响而未能检测到信号变化。总之,还需要进一步优化观测方案和信号提取技术以综合提高观测精度,从而获得对于地下介质应力变化的更高分辨率的观测结果。

4 结语

利用基于精密可控震源的主动地震监测系统在龙门山断裂南段开展了为期1年的连续观测实验。本文利用匹配滤波方法在接收到的弱信号记录中提取出了有效信号并转换为地震子波,再通过信号叠加方法提高了信噪比,并按等时间间隔构成波形时间序列并提取走时变化信息,通过上述处理手段可使观测数据更好地用于地震学研究。

我们利用互相关方法对同一台站同一震相波形进行走时变化检测,提取的走时变化精度优于1ms,误差在0.2ms左右。本文还以绵竹MS5.6地震为例,分析获得了走时剖面和信号变化,结果表明,远距离的龙门山断裂北段地震未对南段实验区域的地下介质波速产生影响。

实验表明,借助于精密可控震源系统,并对观测数据进行处理和信噪比增强后,可提取精确的地震波走时并分析波速变化。这为实际监视地下介质物性的动态变化提供了新的有效技术手段,对了解地震孕育的物理过程具有重要的探索意义。

致谢: 中国科学院地质与地球物理研究所陈棋福研究员为龙门山南段精密可控震源监测实验提出了很多宝贵的指导意见。野外工作区域内条件十分艰苦,工作任务量大,且地形条件极差,高山大川延绵不断,在项目组其他成员(赵翠萍、李建飞,张妍、余伟、朱祖扬、周连庆、梁鸿森、齐诚、王力伟、周青云、苏金蓉、李乐、罗艳、王曙光、彭菲、谭毅培、张项、石玉涛、太龄雪、杨峰)的共同努力下以及四川省地震局、雅安测绘工程院和北京港震机电技术有限公司的协助下,顺利地完成了跨龙门山断裂带的野外主动源观测实验工作。两位匿名评审老师提出了中肯的意见,在此一并致谢。
参考文献
安其美, 丁立丰, 王海忠, 等. 2004, 龙门山断裂带的性质与活动性研究. 大地测量与地球动力学, 24(2): 115–119.
陈颙, 李宜晋. 2007, 地震波雷达研究展望:用人工震源探测大陆地壳结构. 中国科学技术大学学报, 37(8): 813–819.
陈颙, 朱日祥. 2005, 设立"地下明灯计划"的建议. 地球科学进展, 20(5): 485–489.
崔仁胜, 王洪体. 2011, 提取精密控制震源信号的自适应加权匹配滤波方法. 地震, 31(4): 133–139.
古本宗充, 平松良浩, 佐藤隆司. 2005, 地震波速变化和地壳内应力的测定. 世界地震译丛, 3(2): 28–35.
廖成旺, 庄灿涛, 梁鸿森. 2003, 精密可控常时震源系统(ACROSS)的初步实验. 中国地震, 19(1): 89–96.
刘希康, 崔仁胜, 王洪体, 等. 2013, 用Wigner-Hough变换检测精密可控震源信号. 地震, 33(3): 33–42.
刘志坤, 黄金莉. 2010, 利用背景噪声互相关研究汶川地震震源区地震波速度变化. 地球物理学报, 53(4): 853–863.
罗桂纯, 葛洪魁, 王宝善. 2008, 利用相关检测进行地震波速变化精确测量研究进展. 地球物理学进展, 23(1): 26–62.
唐文清, 刘宇平, 陈智梁, 等. 2004, 龙门山断裂构造带GPS研究. 大地测量与地球动力学, 24(3): 57–59.
王洪体, 庄灿涛, 薛兵, 等. 2009, 精密主动地震监测. 地球物理学报, 52(11): 1808–1815.
王伟涛, 2009, 基于人工震源的区域尺度介质波速探测研究, 博士学位论文, 82~84, 合肥: 中国科学技术大学.
杨微, 葛洪魁, 王宝善, 等. 2010, 由精密控制人工震源观测到的绵竹5.6级地震前后波速变化. 地球物理学报, 53(5): 1149–1157.
张正帅, 崔仁胜, 薛兵, 等. 2015, 基于时变窄带滤波技术提取可控震源扫频信号方法研究. 地震, 35(3): 44–56.
Alekseev A S, Chichinin I S, Korneev V A, et al. 2005, Powerful low-frequency vibrators for active seismology. Bull Seism Soc Am, 95(1): 1–17. DOI:10.1785/0120030261.
Poupient G, Ellsworth W L, Frechet J, et al. 1984, Monitoring velocity variations in the crust using earthquake doublets:an application to the Calaveras Fault. J Geophys Res, 89(B7): 5719–5731. DOI:10.1029/JB089iB07p05719.
Kunitomo T, Kumazawa M, 2004, Active monitoring of the earth's structure by the seismic ACROSS-Transmitting and receiving technologies of the seismic ACROSS, First International Workshop on Active Monitoring in the Solid Earth Geophysics, S4-04.
Niu F L, Silver P G, Daley T M, et al. 2008, Preseismic velocity changes observed from active source monitoring at the Parkfield SAFOD drill site. Nature, 454: 204–208. DOI:10.1038/nature07111.
Ryoya I, Koshun Y, Koji M, et al. 2007, Continuous monitoring of propagation velocity of seismic wave using ACROSS. Geophys Res Lett, 1(5): 681–689.
Silver P G, Daley T M, Niu F L, et al. 2007, Active source monitoring of cross-well seismic travel time for stress induced changes. Bull Seism Soc Am, 97: 281–293. DOI:10.1785/0120060120.
Wang B, Zhu P, Chen Y, et al. 2006, Continuous in-situ measurement of stress-induced travel time variation with coda interferometry. AGU Fall Meeting Abstracts.
Wang B, Zhu P, Chen Y, et al. 2008, Continuous subsurface velocity measurement with coda wave interferometry. J Geophys Res, 113(B12). DOI:10.1029/2007JB005023.