地球物理学报  2018, Vol. 61 Issue (11): 4436-4447   PDF    
基于背景噪声成像方法的新疆呼图壁储气库地区近地表速度结构研究
王娟娟1, 姚华建1,2, 王伟涛3, 王宝善3, 李成1, 魏斌4, 冯磊4     
1. 中国科学技术大学地球和空间科学学院, 合肥 230026;
2. 蒙城地球物理国家野外科学观测研究站, 安徽蒙城;
3. 中国地震局地球物理研究所, 北京 100081;
4. 新疆维吾尔自治区地震局, 乌鲁木齐 830011
摘要:近年来,背景噪声成像方法在恢复高频面波信号及获取近地表速度结构方面得到了广泛的应用,本文将该方法应用于准噶尔盆地南缘的呼图壁背斜地区的呼图壁储气库.采用储气库及其周边区域22个台站记录的连续背景噪声数据的垂直分量,通过噪声互相关方法获得了台站对之间的瑞利面波经验格林函数,并进一步提取了0.5~1.5 s的基阶瑞利面波群速度频散曲线.首先根据区域平均频散曲线得到了该地区地下数百米的平均一维横波速度结构,然后利用基于面波射线路径追踪的面波频散直接成像方法得到该地区深度为500 m以上的三维横波速度结构.反演结果显示该地区沉积层较厚,整体横波速度值较小(0.4~0.9 km·s-1).储气库在地表投影区域的横波速度值较小,这可能是由于抽注水、气引起的沉积岩石裂隙所导致.储气库东北和东南方向均有明显的相对高速区,推测是区域地下水位和地形起伏综合作用的结果.本研究获得的近地表三维速度结构为呼图壁储气库地区的上覆地层物性研究、区域微震精定位、场地效应的评估和去除浅层影响的深部介质成像等研究提供了重要基础.
关键词: 背景噪声      群速度      面波层析成像      近地表横波速度结构      储气库     
Study of the near-surface velocity structure of the Hutubi gas storage area in Xinjiang from ambient noise tomography
WANG JuanJuan1, YAO HuaJian1,2, WANG WeiTao3, WANG BaoShan3, LI Cheng1, WEI Bin4, FENG Lei4     
1. School of Earth and Space Sciences, University of Science and Technology of China, Hefei 230026, China;
2. Mengcheng National Geophysical Observatory, University of Science and Technology of China, Mengcheng Anhui, China;
3. Institute of Geophysics, China Earthquake Administration, Beijing 100081, China;
4. Earthquake Administration of Xinjiang Uygur Autonomous Region, Urumqi 830011, China
Abstract: Recently, ambient noise tomography has been widely used in retrieving high frequency surface waves and getting near-surface velocity models. Here we applied this method to the Hutubi gas storage area, which is located at the Hutubi anticline in the southern margin of Junggar basin in Xinjiang. We calculated the cross-correlation functions of the vertical component of continuous ambient noise data of 22 stations in the Hutubi gas storage and its neighboring regions, and measured the fundamental-mode Rayleigh wave group velocity dispersion curves in the period band 0.5~1.5 s. Firstly, using the average dispersion curve in this area, we obtained the one-dimensional average shear-velocity model, and then we used a direct three-dimensional surface-wave tomographic method based on period-dependent raytracing to invert for the 3-D near-surface shear-velocity structure to the maximum depth about 500 meters underground. The results show that the sedimentary layers in this area are thick and have low shear-wave velocity (about 0.4~0.9 km·s-1) in general. The low velocity structure in the gas storage region is probably caused by the fractured sedimentary rocks due to gas and water extraction and injection. In the southeast and northeast of the gas storage, there are two zones with relatively high velocities, which are probably due to the integrated effects of the local groundwater level and topographic relief. The near-surface velocity model can provide an important basis for understanding the overlying formation characteristics, precise location of local microearthquakes, assessing near-surface site effect, and removing the shallow structural effect in imaging the deep structures.
Keywords: Ambient noise    Group velocity    Surface-wave tomography    Near-surface shear-wave velocity structure    Gas storage    
0 引言

呼图壁储气库位于新疆维吾尔族自治区呼图壁县城东4.5 km处,在该处新疆油田公司发现了准噶尔盆地南缘第一个气田.一方面,这个气田油气已经枯竭, 另一方面北方地区冬夏两季用气量差别较大,且每一天的不同时段用气量也差别较大,所以中石油将其改建为储气库,用于季节调峰和战略储备.储气库的建设世界上早有先例,加拿大在1915年建成了世界上第一个枯竭油气田储气库,随后美国和前苏联等都建造了储气库.中国于1999年在大港油田大张坨建造了第一个地下储气库.储气库具有储藏量大、调峰能力强、安全性高和建造运行成本低等优势,因此在众多国家得到了大力发展.目前世界上的储气库主要由枯竭油气田,含水层,盐穴和废弃矿井等改造而成(丁国生和谢萍,2006曹锡秋,2013).

呼图壁储气库(如图 1)位于准噶尔盆地南缘(北天山山前坳陷),构造形态是近东西向的长轴背斜,沿长轴方向分别被两个近平行的逆掩断层(呼图壁断裂和呼图壁北断裂)切割(曲国胜等,2008).呼图壁背斜是发育在吐谷鲁背斜东北方向呼图壁以东到昌吉之间戈壁平原上的一个轻微隆起,该背斜下部7~8 km处侏罗纪煤系地层中滑脱面开始向上扩展,形成一个断坡,该背斜应属于逆断裂扩展或断裂弯曲型背斜(陈立春,2011).呼图壁背斜属天山山前坳陷第三排构造带东端(邓起东等,1999).北天山坳陷是因为海西运动期天山地槽褶皱回返形成的,是一个持续性的沉积坳陷.印支运动使该坳陷区的三叠系和侏罗系不整合接触.燕山运动期,天山强烈上升,山前坳陷沉降.喜马拉雅期,印度板块向欧亚板块俯冲,造成盆地基底区域性向南倾斜,山前坳陷沉降.中新世晚期山前一带遭强挤压,形成由山前到盆地腹部方向的3排背斜和两排向斜带.呼图壁背斜主要形成于喜马拉雅运动期,可能是第二排逆断层-背斜带向东扩展的产物(康竹林,1997).北天山山前坳陷沉积岩最大厚度达15000 m,速度较低(宋元林和胡新平,2001).

图 1 (a) 研究区域地质背景,黑色方框表示研究区域;(b)研究区域的台站分布及构造简图,灰色阴影是呼图壁储气库在地表的投影,在灰色阴影上的红色细线表示断裂,1、2和3分别对应呼图壁断裂、呼图壁北断裂和呼001井北断裂,三角形表示使用的22个台站,蓝色多边形表示反演的核心区域,粗绿线连接XJ01-XJ06台站对,粗红线连接XJ01-XJ19台站对 Fig. 1 (a) The geological background of the Hutubi gas storage. The area in the black box is the study region; (b) The distribution of seismic stations (triangles) and tectonic sketch map of the study area. The projection of the Hutubi gas storage on the surface is shown as the shaded area. Three faults 1, 2, and 3 (the red thin lines in the shaded area) are the Hutubi fault, Hutubi North fault, and Hutubi No.001 North Well fault, respectively. The main region of the inversion is given by the blue polygon. The two-station pairs XJ01-XJ06 and XJ01-XJ19 are shown as the bold green line and the bold red line, respectively

浅层沉积盆地的三维速度结构会对地震波的传播造成很大的影响,当地震波从地球内部传播到近地表区域时,经历了从高速区到低速区的变化,会引起比较强烈的振幅放大、共振、振动持续时间长等场地效应,可能会造成一些不可估量的地震灾害,所以有必要对场地效应进行评估,而评估的关键参数就是近地表三维速度结构(Flores-Estrella et al., 2006).同时,当研究深部结构时,浅层结构的干扰比较大,通常需要去除.所以研究该地区的近地表速度结构对于储气库区域的地震场地效应评估及深部结构和构造研究具有重要的意义.由于注采气过程会引发区域小地震的活动,可靠的浅层三维结构模型也为提升该地区微震定位的精度提供模型基础.

本文将使用背景噪声成像来研究呼图壁储气库地区的三维近地表速度结构.背景噪声成像是近年来发展迅速的一种成像方法,它的基本思想是通过计算台站对之间的噪声互相关函数来近似得到台站对间的经验格林函数.20世纪50年代Aki(1957)提出SPAC方法,通过环形台阵记录的噪声信号的互相关计算提取面波格林函数,这通常被认为是这个方法的一个开端.Weaver和Lobkis在声学研究领域通过实验发现两个传感器间的格林函数可以通过计算两点之间的热起伏噪声的波形互相关得到(Lobkis and Weaver, 2001; Weaver, 2005). Campillo和Paul建立了两点观测到的地震尾波互相关函数和格林函数之间的联系,推动了这个方法在地震学领域的发展(Campillo and Paul, 2003).鉴于地震尾波依赖于地震的发生强度和地下散射体分布等一些不足之处,Shapiro和Campillo等合作研究提出可以从连续背景噪声的互相关函数中提取基阶瑞利面波来反演地壳速度结构(Shapiro and Campillo, 2004Shapiro et al., 2005).Sabra等人对南加州地区超过150个台站的连续30天的背景噪声记录的处理,并提出了时间域的格林函数(TDGF)可以通过两个台站之间长时间平均的背景噪声互相关函数对时间求导得到(NCF)(Sabra et al., 2005).

随着数据处理流程的发展和各种数据处理技术的完善,背景噪声成像方法在大区域(如美国、欧洲和中国大陆等)和小区域(如美国南加州、新西兰、青藏高原东南缘和华北等)的应用,均得到了高信噪比的经验格林函数及较可靠的面波层析成像结果(例如Sabra et al., 2005; Yao et al., 2006, 2008, 2010; Yang et al., 2007; Lin et al., 2007; Zheng et al., 2008; Bensen et al., 2008; 房立华和吴建平, 2009).一般情况下,高频面波对浅层结构比较敏感,低频面波对深部结构比较敏感.先前的背景噪声成像研究面波频带主要集中在5~40 s周期段,用于反演地壳和上地幔顶部的结构.近年来应用背景噪声提取的高频面波成分来获取近地表速度结构的研究取得了显著的进展(例如Huang et al., 2010; Lin et al., 2013; Fang et al., 2015; Li et al., 2016).

本研究利用呼图壁储气库周边22个台站记录的连续波形数据的垂直分量,从背景噪声的互相关函数中提取了0.5~1.5 s周期段的基阶瑞利面波的群速度频散曲线.根据平均频散曲线得到了该地区地下数百米的一维横波速度结构,并将其作为反演的初始模型,然后利用所有路径的群速度频散曲线及基于射线路径追踪的面波频散直接反演方法(Fang et al., 2015),得到了研究区域更为精细的近地表三维横波速度结构.

1 数据

本研究挑选了呼图壁储气库附近22个流动台站的数据.其中15个台站,如图 1b中的红色三角形所示,最小台间距为0.6 km,在2015年7月到2016年6月,均采用Trillium地震计, 频带为40 s~50 Hz, 数据记录连续性比较好;从2016年6月至2016年12月,这15个台站采用CMG40T地震计,频带为2 s~100 Hz, 数据记录连续性较好.另外7个台站仪器类型为CMG-6TD一体化地震计,频带为30 s~100 Hz, 如图 1b中的黑色三角形所示,所使用的数据记录时间从2013年8月到2016年7月.由于部分仪器安放场地温度的原因,数据记录连续性比较差(王芳,2017).因为本次研究测量的是群速度,即台站对之间的距离除以不同周期的波包的到时.由于仪器本身的相位几乎不影响波包的到时,因此在后续数据处理中,没有去除仪器响应.

2 背景噪声数据处理

背景噪声数据处理流程主要参照Bensen等(2007)介绍的背景噪声数据提取面波频散曲线的方法, 简要概括为单台数据预处理,计算噪声互相关函数和提取面波频散曲线三个步骤.具体分为:

(1) 单台数据预处理:原始的垂直分量连续波形数据是mseed格式,采样率为100 Hz.首先进行数据格式转换、插值拼接、去均值、去趋势,再重采样到10 Hz;对数据做谱白化处理,拓宽频谱,防止单一频谱对噪声信号的干扰(Bensen et al., 2007房立华和吴建平, 2009).将数据滤波到0.5~2 s周期段,然后对数据进行时域归一化(Bensen et al., 2007)来减小地震信号、仪器噪声和台站附近不稳定噪声对噪声互相关函数质量的影响.

(2) 计算噪声互相关函数:将每个台站对、每天预处理后的0.5~2 s频带的数据按照记录时间对齐并做互相关计算,最大互相关时间延迟取为150 s,最后通过线性叠加得到每个台站对在该频带的所有天叠加的噪声互相关函数,互相关时间从100天到520天不等.我们挑选出了信噪比较高的噪声互相关函数,并按照台间距排列在一起,如图 2所示,可以看出正负半轴均有比较明显的信号,如图中黑色虚线标出的范围,正负分支并非完全对称.另外,从图中也可以估计该地区的瑞利面波群速度大小基本介于0.3~0.9 km·s-1之间,其中较慢的基阶面波比较靠近0.3 km·s-1的群速度线.

图 2 所有双台路径0.5~2 s周期的噪声互相关函数(信噪比大于5),其中对称的黑色虚线表示群速度分别为0.3 km·s-1和0.9 km·s-1的时距线 Fig. 2 All the interstation cross-correlation functions with the signal-to-noise ratio (SNR) greater than 5 in the 0.5~2 s period band. The black dashed lines give the group velocity range between 0.3 km·s-1 and 0.9 km·s-1 as indicated

(3) 提取面波频散曲线:计算的噪声互相关函数有正负分支,分别反映信号从作为虚拟源的A台站到接收器B台站和作为虚拟源的B台站到接收器A台站的经验格林函数.理想情况下,如果噪声源分布均匀,正负两支完全对称.而在实际处理中,则常常需要通过平均正负分支把双边信号转换成单边信号(对称分量信号)的方法来减轻因噪声源分布不均匀带来的影响(如图 2)(Yang et al., 2007; 房立华和吴建平, 2009).Lin等(2007)通过单独计算正负分支的信噪比谱发现对称分量的信噪比在各个周期都相对较高.本次数据处理采用的是对称分量信号.瑞利面波群速度频散曲线的测量是基于多时窗的时频滤波方法(Dziewonski et al., 1969; Levshin et al., 1992; Yao et al., 2011).我们仅挑选出信噪比大于5的噪声互相关函数进行频散分析.信噪比的定义如下:首先确定信号窗和噪声窗,如图 3a所示,黑色波形对应于信号窗(对应于群速度窗口大于0.3 km·s-1),红色波形对应于噪声窗(信号窗后30 s),信噪比定义为信号窗内波形的波包最大振幅值与噪声窗内波形振幅绝对值的平均值的比值,在实际分析中通过高斯滤波的方法分别计算每个周期的信噪比.在频散曲线提取时,要求台间距大于2倍的波长,以满足面波传播远场近似的假设,尽管有学者提出台间距大于1倍或1.5倍的波长即可(Yao et al., 2011; Luo et al., 2015).以XJ01与XJ06台站对间的频散曲线为例,从图 3a可以看出这个台站对的面波信号大概到时为20 s,已知台间距约为8 km,由此估算面波群速度为0.4 km·s-1.从图 3b群速度频散分析的时频能量分布图可以看出,这个台站对的面波能量比较集中,提取的频散曲线周期范围大概为0.5~1.4 s,提取的群速度频散数值(红点所示)与图 3a反映的结果基本一致.

图 3 XJ01-XJ06台站对的频散曲线测量示意图 (a)该台站对的噪声互相关函数,黑色波形表示信号时窗波形,红色波形表示噪声时窗的部分波形;(b)多重滤波分析后得到的时频能量分布图,黑线表示提取的群速度频散曲线,红点表示最终保留的对应周期信噪比大于5的频散. Fig. 3 Group velocity dispersion measurement of the station pair XJ01-XJ06 (a) The noise cross-correlation function of XJ01-XJ06. The black waveform shows for the signal window and red for the noise window; (b) Time-period energy distribution diagram after multiple filtering analysis. The black curve shows the traced group velocity dispersion curve with the red solid circles for the periods with the signal to noise ratio (SNR) above 5 that are finally saved.

0.5~2 s周期段的所有基阶瑞利面波的群速度频散曲线如图 4所示.结果显示该地区的速度比较低,群速度随周期的变化不大;但频散曲线值较分散,从0.4 km·s-1到1 km·s-1都有分布,在后续的一维平均横波速度结构反演中我们将其分为两个部分.图 4中绿色的频散曲线是图 1b中绿线连接起来的台站对XJ01-XJ06之间的频散曲线,它代表穿过储气库或者储气库附近的群速度较低的频散曲线.图 4中红色的频散曲线是图 1b中粗红线连接起来的台站对XJ01-XJ19之间的频散曲线,它代表远离储气库的群速度较高的频散曲线.当周期超过1.5 s时,路径数不足30条,所以接下来只考虑0.5~1.5 s周期段的频散情况,最终得到202条基阶瑞利面波的群速度频散曲线,用于反演研究区域的三维近地表横波速度结构.

图 4 0.5~1.7 s周期段的基阶瑞利面波的群速度频散曲线(黑色实线)及每个周期的频散曲线的数量(黑色虚线),其中红色和绿色频散曲线分别对应于图 1中红色和绿色双台路径之间的测量频散曲线 Fig. 4 Fundamental mode Rayleigh wave group velocity dispersion curves (black solid lines) in the 0.5~1.7 s period band and the number of dispersion measurements (black dashed line) at every period. The dispersion curves in red and green correspond to the two-station pairs connected with the red and blue ray paths in Fig. 1, respectively
3 区域平均一维横波速度结构的反演

为了得到较好的三维横波速度结构反演的初始模型,本研究将整个区域的平均频散曲线(图 5a)反演得到研究区域平均一维横波速度模型,这里采用线性最小二乘迭代反演方法,初始速度模型由CRUST1.0模型(Laske et al., 2012)插值得到,密度模型通过经验关系式(Brocher, 2005)得到.提取频散曲线时发现,当某些台站对之间的射线路径离储气库较远(不经过储气库)或者极少部分路径经过储气库时,得到的群速度明显高于射线路径主要穿过储气库地区的台站对的群速度.所以在反演一维横波速度结构时,首先用整个区域的平均频散(图 5a)反演得到了区域平均的一维横波速度结构(图 5d),然后对一条群速度比较低的频散曲线(图 4中绿色的频散曲线,图 5b,对应于图 1b中XJ01-XJ06绿色路径)和一条群速度比较高的频散曲线(图 4中红色的频散曲线,图 5c,对应于图 1中XJ01-XJ19红色路径)分别来反演,获得了两条不同路径下方平均的一维横波速度结构.反演得到的台站对XJ01-XJ06之间的平均一维横波速度结构(图 5e)较低,而台站对XJ01-XJ19之间的一维横波速度结构较高,但是对比发现整个区域的平均速度结构(图 5d)和XJ01-XJ06之间的速度结构(图 5e)相似.

图 5 (a) 所有路径的平均群速度频散曲线; (b) XJ01-XJ06台站对(图 1b中绿线所示)的频散曲线; (c) XJ01-XJ19台站对(图 1b中红线所示)的频散曲线; (d)、(e)和(f)分别为(a)、(b)和(c)中频散曲线反演得到的平均一维横波速度结构模型 Fig. 5 (a) The average dispersion curve of all two-station pairs; (b) The dispersion curve of the station pair XJ01-XJ06 (green path in Fig. 1b); (c) The dispersion curve of the station pair XJ01-XJ19 (red in Fig. 1b); (d, e, f) give the average one-dimensional shear velocity model obtained from the inversion of the dispersion curve in (a, b, c), respectively

为了解面波群速度在深度方向对横波速度的约束效果,需要计算群速度对横波速度的深度敏感核.图 6是基于该地区平均一维横波速度模型计算得到的基阶瑞利面波群速度对横波速度的深度敏感核,从结果可以看出在选取的面波频带内,提取的频散曲线可以对500 m以上的地壳浅层结构有较好的约束.

图 6 不同周期(0.5 s, 0.8 s, 1.1 s和1.4 s)的基阶瑞利面波群速度对横波速度的深度敏感核() Fig. 6 Depth sensitivity kernels () of the fundamental mode Rayleigh wave group velocities to shear velocities at different periods (0.5 s, 0.8 s, 1.1 s, and 1.4 s)
4 三维横波速度结构反演方法及结果

本研究使用由所有路径基阶瑞利面波频散走时反演近地表三维横波速度结构的面波直接成像方法(Fang et al., 2015; Li et al., 2016),它考虑了复杂介质中不同频率面波射线弯曲的影响, 且避免了反演二维群速度和相速度图.对于正演问题,首先把模型离散为K个规则的水平网格节点, 每个网格点下对应着一个由垂直方向离散格点控制的一维层状模型,用于计算该格点处对应的群速度频散曲线,然后通过双线性插值得到角频率为ω的二维水平面上的群速度分布图, 最后使用快速行进法(Rawlinson and Sambridge, 2004)计算所有源和台站之间的不同频率的面波二维射线路径和走时.

在高频近似的情况下,面波在源A和接收器B之间沿路径的传播时间可以离散表示为

(1)

Sp是沿着面波传播路径的ΔlAB段的群慢度(群速度的倒数),P是射线路径总段数.

对于反演问题,需要找到一个模型m,使得不同路径所有频率ω的观测到时和模型的正演到时之间的差值最小.在角频率为ω时,参考模型的射线路径i的传播时间扰动可以表示为

(2)

tiobs(ω)是观测的面波传播时间,ti(ω)是由初始参考模型或每一步迭代之后得到的模型正演计算走时,vik是沿着与第i个传播时间相关的射线路径的双线性插值系数,Uk(ω)和δUk(ω)分别是二维平面上第k个网格点的群速度及其扰动.其中群速度扰动可以用一维层状模型的横波速度、纵波速度和密度的相对扰动来表示.由于面波频散对横波速度β比较敏感,所以我们通过经验关系式(Brocher, 2005),把纵波速度α和密度ρ的扰动通过横波速度的扰动来表达(Fang et al., 2015).最后第i条路径的传播时间扰动可以写为

(3)

θk代表平面上第k个网格点的一维模型, J表示深度方向的网格点数,M=KJ.Rα(Z)和Rρ(Z)分别表示纵波速度和密度相对于横波速度的拟合多项式系数(Fang et al., 2015). 分别为群速度对纵波速度、密度和横波速度的一维深度敏感核.

方程(3)可以写成d=GM的形式,d是所有路径所有频率的传播时间残差矩阵,G是数据敏感核矩阵,M是待求解的三维横波速度扰动向量.通过反演可得到速度模型的更新量,然后在新的模型下正演计算射线路径和敏感核矩阵,并进一步迭代更新模型,直到反演收敛.这里为了反演的稳定性,我们对模型参数施加了平滑约束.

反演过程中模型的网格参数设置如下:南北方向的网格间距是0.01°,东西方向的网格间距是0.02°,水平方向共有网格27×25个.垂直方向的网格点是10个, 深度范围为0~530 m.

为了检测反演结果对初始模型依赖程度,本研究使用了两种初始一维模型:①由该地区平均频散曲线反演得到的平均一维横波速度结构模型;②对XJ01-XJ06台站对间的频散曲线反演得到的一维横波速度结构模型和XJ01-XJ19台站对间的频散曲线反演得到的一维横波速度结构模型求平均,得到一个一维横波速度结构模型.对以上两种一维横波速度结构模型分别进行三维延拓作为初始模型.采用初始模型①得到的三维横波速度结构模型如图 7所示.采用不同初始模型的纵剖面结果如图 8a图 8b所示, 分别对应①和②模型,结果基本一致,但是初始模型①残差较小,所以后续只以初始模型①的反演结果为例进行讨论.

图 7 三维横波速度结构反演结果 (a), (b), (c)和(d)分别是60 m, 180 m, 300 m和420 m深度的横波速度结构剖面图, 近椭圆黑色曲线表示储气库的边界,红色的线段表示断裂,黑色三角形表示台站,(a)中两条近垂直的紫线(pp′, qq′)表示图 8中的纵剖面位置,(c)中两条蓝色线段bo1-bo2, bo2-bo3和图形南侧边界之间的区域是工业园区的大概位置,两个黑色圆点分别表示HK19井和HUK20井的位置. Fig. 7 The three-dimensional Vs model after the inversion (a—d) The Vs structure at four depths: 60 m, 180 m, 300 m and 420 m. The Elliptical shape shows the boundary of Hutubi gas storage, the red lines show the three faults, the black triangles represent the seismic stations, the two nearly perpendicular purple lines (pp′, qq′) in (a) show the two vertical cross-sections in Fig. 8, the area between two blue lines (bo1-bo2, bo2-bo3) in (c) and the southern boundary of the graph shows the location of the industry park, and the two black dots represent the location of the HK19 and HUK20 wells.
图 8 沿图 7a中pp′, qq′中纵剖面的横波速度结果 图中的黑色+号表示储气库的边界,黑色三角形表示台站,黑点表示HK19井的位置,(a)和(b)分别为初始模型①和②的反演结果. Fig. 8 Vertical Vs cross-sections (pp′, qq′) with the location marked in Fig. 7a The black plus represents the boundary of the gas storage and the black triangle represents the seismic station, the black dot represents the HK19 well. (a) and (b) are the results from the initial models ① and ②, respectively.
5 结果的分辨率与讨论

图 9反映了三个深度的棋盘测试恢复情况,总体来说在反演区域恢复较好.当深度逐渐增加时,在储气库东南方向的恢复会逐渐变差.图 10展示反演区域(蓝色方框)有较好的路径覆盖,但是在储气库东南方向路径稍微稀疏一些,这也与棋盘测试结果吻合较好.反演过程是比较稳定的,总的走时残差的均方差随着迭代次数的增加而不断下降(图 11a),最终收敛.反演后走时残差的平均值由0.17 s降为0.091 s,基本呈高斯分布(图 11a).初始模型每个周期的走时残差平均值(图 11b中蓝色线)与反演后模型每个周期的走时残差平均值(图 11b中红色线)相差较小,表明反演的初始模型比较合理.

图 9 反演结果分辨率的棋盘测试 (a)设置的棋盘测试的输入模型; (b), (c)和(d)分别是60 m, 180 m和360 m的棋盘测试恢复结果. Fig. 9 Checkerboard resolution tests of the inversion results (a) The input checkerboard model; (b—d) The recovered checkerboard results at depths 60 m, 180 m and 360 m.
图 10 两个周期(0.6 s,1.2 s)的路径覆盖,红色三角形表示台站,蓝色多边形表示反演结果较可靠的区域(即射线路径覆盖较好的区域) Fig. 10 Ray-path coverage for the two selected periods: 0.6 s and 1.2 s. The red triangles represent the stations, and the blue polygon shows the boundary of region with more reliable inversion results (i.e., with better path coverage)
图 11 (a) 反演过程中总走时残差的均方差的变化及反演前后的总走时残差分布情况(蓝色及红色柱状图所示); (b)反演前后每个周期的走时残差的均值(蓝线和红线所示)和标准偏差(蓝色和红色误差棒所示) Fig. 11 (a) Variation of the standard deviation of all travel-time residuals with the iteration number and the histograms of all travel-time residuals before inversion (blue bars) and after 10 iterations (red bars); (b) The distribution of the mean travel-time residuals (blue and red lines) and the corresponding standard deviations (blue and red error bars) of different periods before inversion (blue) and after 10 iterations (red)

由于台站在气库区域较多,在外围区域(尤其是其南部)台站较少.图 10所示蓝色多边形区域外侧台站的加入,增加了反演核心区域的路径覆盖,有助于增加核心地区反演结果的可靠性.但是蓝色多边形(图 10)外侧区域路径覆盖较少,这可能导致反演结果对初始模型的依赖程度较大,可能会导致反演出现一些虚假异常,所以本研究仅讨论蓝色多边形区域内的速度结构异常.

曹锡秋(2013)分析了HUK20井(位置如图 7所示)的岩石分层情况,李延军等(2004)给出了呼2井(储气库西北段)的岩石分层情况,两个井的岩石分层基本一致,但是前者在本次研究深度范围内分层更精细.李一峰(2014)李杰(2016)等人给出了呼图壁储气库较为清晰的井位信息,再结合前人(陈华慧等,1994宋元林和胡新平,2001李学义等,2003)的一些研究结果,呼图壁储气库地区发育的地层从上到下依次是:第四系西域组(Q1x),平均厚度429 m;新近系独山子组(N2d),平均厚度1341 m;塔西河组(N1t),平均厚度440 m;沙湾组(N1s),平均厚度297 m;古近系安集海河组(E2-3a),平均厚度850 m;紫泥泉子组,厚度575 m(E1-2z).图 7中结果显示呼图壁储气库区域的横波速度很低,我们认为导致这个低速的主要原因是呼图壁储气库位于准噶尔盆地南缘,沉积层厚度较大,且较年轻,也有可能与储气库区域抽注水与抽放气导致该区沉积岩层裂隙发育有关.图 7显示了三条断裂,从南到北依次是呼图壁断裂,呼图壁北断裂,呼001井北断裂(结合图 1).呼图壁断裂将地层分为上下盘,在下盘发育呼图壁北断裂,储气库主要是被呼图壁断裂和呼图壁北断裂控制.呼图壁断裂是最主要的且断距最大的断裂带,地震反射资料揭示,呼图壁断裂断开的层系最多,从侏罗系(J)到新近系塔西河组(N1t),而其余两条断裂带均在侏罗系(J)和古近系安集海河组(E2-3a)发育(曹锡秋,2013康竹林,1997邓起东等,1999曲国胜等,2008).由于断裂发育的深度深于反演深度,所以速度结构异常受断裂带的分布影响可能很小,图 7的反演结果也显示断裂与速度结构异常之间并不存在空间上明显的相关性.

图 12展示了储气库东南段的一个测井HUK20 (位置如图 7所示)的地质资料(曹锡秋,2013),研究区域的地层主要是第四系西域组,它的主要岩性是灰色砂砾岩和砂质小砾岩,中间夹层为褐灰色泥岩.该地区的近地表区域岩性变化比较大,精细的地层划分是有必要的.图 8展示了两条分别沿着储气库走向和垂直储气库走向的纵向剖面的横波速度结果.本次反演的最大深度在500 m左右,这大约到了第四系西域组的底部.总体来说,浅层速度值较小,深部速度值较大,但是在不同的深度速度值有一些横向变化,可能是受泥岩夹层的影响.HK19气井(位置如图 7所示)的南北两侧分别是工业园区和大面积的农田.储气库区所在的呼图壁县二十里店镇,农田水利灌溉大多以地下水为主,地下水超采引起强烈的地面沉降(李杰等,2016).在卫星影像地图上可以清楚地发现位于蓝色线段bo1-bo2和bo2-bo3之间的区域是工业园区(图 7c),储气库的南侧海拔高于北侧海拔约50 m.杨宏伟和朱瑾(2012)认为呼图壁河流域下游和玛纳斯河流域中游2009年的平均地下水位均超过400 m,玛纳斯河流域下游多年平均地下水位在360 m左右.结合图 12分析,如图 8a中黑框区域所示,qq′剖面右侧(储气库的东南方向)和pp′剖面的右侧(储气库东北方向)速度值较高,推测可能是地下水位差异和地形起伏综合作用的结果, 也有可能受地下岩性差异的影响.

图 12 研究区域HUK20井的钻井地质资料(修改自曹锡秋,2013).Q1x:第四系西域组, N2d:新近系独山子组.黑点区域表示灰色砂砾岩,砂质小砾岩和泥质小砾岩.褐灰色,浅棕色和黑色区域分别表示褐灰色泥岩,浅棕色泥岩和砂泥岩 Fig. 12 Borehole data of the HUK20 well in study region (modified from Cao, 2013). Q1x: Xiyu group of the Quaternary, N2d: Dushanzi group of the Neogene. The patches with the black dots are gray sandy gravel, sandy breccia, and mud conglomerate. And the patches in taupe gray, light brown, and black represent brown mudstone, light brown mudstone, and sand-mudstone, respectively
6 结论

本文使用呼图壁储气库周边22个台站的连续波形数据,利用背景噪声成像方法获得了呼图壁储气库地区的浅层横波速度结构.近地表横波速度值整体较小,大致介于0.4~0.9 km·s-1之间,可能与该区域较厚、较年轻的松散沉积层有关.储气库区域横波速度值较低,有可能是抽注水和抽注气引起的岩石裂隙导致的.从纵剖面结果来看,浅层速度值较小,深部速度值较大,在不同的深度速度值有一些变化,这可能是由于泥岩夹层的影响.纵剖面240 m以下,储气库的东南方向和东北方向速度值较高,我们认为这可能主要是不同区域地下水和地形起伏综合作用的结果.近地表结构复杂,未来更加精细的成像工作需要与密集台站和岩石地化资料相结合,以更好的认识储气库近地表的结构与构造特征.

致谢  感谢两位评审专家提出的修改完善建议.
References
Aki K. 1957. Space and time spectra of stationary stochastic waves, with special reference to microtremors. Bulletin Earthquake Research Institute Tokyo University, 35: 415-456.
Bensen G D, Ritzwoller M H, Barmin M P, et al. 2007. Processing seismic ambient noise data to obtain reliable broad-band surface wave dispersion measurements. Geophysical Journal International, 169(3): 1239-1260. DOI:10.1111/j.1365-246X.2007.03374.x
Bensen G D, Ritzwoller M H, Shapiro N M. 2008. Broadband ambient noise surface wave tomography across the united states. Journal of Geophysical Research:Solid Earth, 113(B5): B05306. DOI:10.1029/2007jb005248
Brocher T M. 2005. Empirical relations between elastic wavespeeds and density in the earth's crust. Bulletin of the Seismological Society of America, 95(6): 2081-2092. DOI:10.1785/0120050077
Campillo M, Paul A. 2003. Long-range correlations in the diffuse seismic coda. Science, 299(5606): 547-549. DOI:10.1126/science.1078551
Cao X Q. 2013. A reaserch on reservoir geomechanic features of a gas storage in Xinjiang after natural depletion[Ph. D. thesis] (in Chinese). Beijing: China University of Geosciences (Beijing).
Chen H H, Lin X L, Guan K N, et al. 1994. Early Pleistocene deposits and its lower boundary (Q/N) in Tianshan MT. , Xinjiang region. Quaternary Sciences (in Chinese), 14(1): 38-47.
Chen L C. 2011. Late Quaternary behavior of the active tectonic system in the Urumqi transform region of the north Tianshan[Ph. D. thesis] (in Chinese). Beijing: Institute of Geology, China Earthquake Administration.
Deng Q D, Feng X Y, Zhang P Z, et al. 1999. Reverse fault and fold zone in the Urumqi range-front depression of the northern Tianshan and its genetic mechanism. Earth Science Frontiers (in Chinese), 6(4): 191-201.
Ding G S, Xie P. 2006. Current situation and prospect of Chinese underground natural gas storage. Natural Gas Industry (in Chinese), 26(6): 111-113.
Dziewonski A, Bloch S, Landisman M. 1969. A technique for the analysis of transient seismic signals. Bulletin of the Seismological Society of America, 59(1): 427-444.
Fang H J, Yao H J, Zhang H J, et al. 2015. Direct inversion of surface wave dispersion for three-dimensional shallow crustal structure based on ray tracing:Methodology and application. Geophysical Journal International, 201(3): 1251-1263. DOI:10.1093/gji/ggv080
Fang L H, Wu J P. 2009. Measurement of Rayleigh wave dispersion from ambient seismic noise and its application in North China. Acta Seismologica Sinica (in Chinese), 31(5): 544-554.
Flores-Estrella H, Yussim S, Lomnitz C. 2006. Seismic response of the mexico city basin:A review of twenty years of research. Natural Hazards, 40(2): 357-372. DOI:10.1007/s11069-006-0034-6
Huang Y C, Yao H J, Huang B S, et al. 2010. Phase velocity variation at periods of 0.5-3 seconds in the Taipei basin of Taiwan from correlation of ambient seismic noise. Bulletin of the Seismological Society of America, 100(5A): 2250-2263. DOI:10.1785/0120090319
Kang Z L. 1997. Petroleum prospect in the southern margin of Junggar Basin. Petroleum Explorationist (in Chinese), 2(4): 31-34.
Laske G, Masters G, Ma Z, et al. 2012. CRUST1.0:An updated global model of Earth's crust. EGU General Assembly, 14: 3743.
Levshin A, Ratnikova L, Berger J. 1992. Peculiarities of surface-wave propagation across central Eurasia. Bulletin of the Seismological Society of America, 82(6): 2464-2493. DOI:10.1785/0120130319
Li C, Yao H J, Fang H J, et al. 2016. 3D Near-surface shear-wave velocity structure from ambient-noise tomography and Borehole data in the Hefei Urban Area, China. Seismological Research Letters, 87(4): 882-892. DOI:10.1785/0220150257
Li J, Li R, Wang X Q, et al. 2016. Research on surface vertical deformation in the Hutubi underground gas storage reservoir. Earthquake Research in China (in Chinese), 32(2): 407-416.
Li X Y, Shao Y, Li T M. 2003. Three oil-reservoir combinations in south marginal of Jungar Basin, Northwest China. Petroleum Exploration & Development (in Chinese), 30(6): 32-34.
Li Y F, Li Y H, Gao Q, et al. 2014. New understanding of reservoir of Z2 sand layer of Ziniquanzi formation in Hutubi gas storage. Xinjiang Petroleum Geology (in Chinese), 35(2): 182-186.
Li Y J, Wang T D, Zhang Y Y. 2004. Natural gas genesis and formation of gas pools in the South margin of Junggar Basin. Acta Sedimentologica Sinica (in Chinese), 22(3): 529-534. DOI:10.3969/j.issn.1000-0550.2004.03.023
Lin F C, Ritzwoller M H, Townend J, et al. 2007. Ambient noise rayleigh wave tomography of new zealand. Geophysical Journal International, 170(2): 649-666. DOI:10.1111/j.1365-246X.2007.03414.x
Lin F C, Tsai V C, Schmandt B, et al. 2013. Extracting seismic core phases with array interferometry. Geophysical Research Letters, 40(6): 1049-1053. DOI:10.1002/grl.50237
Lobkis O I, Weaver R L. 2001. On the emergence of the green's function in the correlations of a diffuse field. The Journal of the Acoustical Society of America, 110(6): 3011-3017. DOI:10.1121/1.1417528
Luo Y H, Yang Y J, Xu Y X, et al. 2015. On the limitations of interstation distances in ambient noise tomography. Geophysical Journal International, 201(2): 652-661. DOI:10.1093/gji/ggv043
Qu G S, Ma Z J, Zhang N, et al. 2008. Fault structures in andaround Junggar Basin. Xinjiang Petroleum Geology (in Chinese), 29(3): 290-295.
Rawlinson N, Sambridge M. 2004. Wave front evolution in strongly heterogeneous layered media using the fast marching method. Geophysical Journal International, 156(3): 631-647. DOI:10.1111/j.1365-246X.2004.02153.x
Sabra K G, Gerstoft P, Roux P, et al. 2005. Extracting time-domain green's function estimates from ambient seismic noise. Geophysical Research Letters, 32(3): L03310. DOI:10.1029/2004GL021862
Shapiro N M, Campillo M. 2004. Emergence of broadband rayleigh waves from correlations of the ambient seismic noise. Geophysical Research Letters, 31(7): L07614. DOI:10.1029/2004GL019491
Shapiro N M, Campillo M, Stehly L, et al. 2005. High-resolution surface-wave tomography from ambient seismic noise. Science, 307(5715): 1615-1618. DOI:10.1126/science.1108339
Song Y L, Hu X P. 2001. Characteristics and comprehensive evaluation of Ziniquanzi reservoir. Journal of Xi'an Petroleum Institute (Natural Science Edition) (in Chinese), 16(4): 15-19.
Wang F. 2017. Shallow structure and temporal variation monitoring around Hutubi underground gas storage in Xinjiang based on seismic ambient noise[Ph. D. thesis] (in Chinese). Beijing: Institute of Geophysics, China Earthquake Administration.
Weaver R L. 2005. Information from seismic noise. Science, 307(5715): 1568-1569. DOI:10.1126/science.1109834
Yang H W, Zhu J. 2012. The groundwater dynamic and sustainable utilization countermeasuress of groundwater resources in southernplain of Junggar basin in recent 30 years. Xinjiang Geology (in Chinese), 30(3): 350-354.
Yang Y J, Ritzwoller M H, Levshin A L, et al. 2007. Ambient noise rayleigh wave tomography across Europe. Geophysical Journal International, 168(1): 259-274. DOI:10.1111/j.1365-246X.2006.03203.x
Yao H J, Van Der Hilst R D, De Hoop M V. 2006. Surface-wave array tomography in SE Tibet from ambient seismic noise and two-station analysis-Ⅰ. Phase velocity maps. Geophysical Journal International, 166(2): 732-744. DOI:10.1111/j.1365-246X.2006.03028.x
Yao H J, Beghein C, Van Der Hilst R D. 2008. Surface wave array tomography in SE Tibet from ambient seismic noise and two-station analysis-Ⅱ.Crustal and upper-mantle structure. Geophysical Journal International, 173(1): 205-219. DOI:10.1111/j.1365-246X.2007.03696.x
Yao H J, Van Der Hilst R D, Montagner J P. 2010. Heterogeneity and anisotropy of the lithosphere of SE Tibet from surface wave array tomography. Journal of Geophysical Research:Solid Earth, 115(B12): B12307. DOI:10.1029/2009JB007142
Yao H J, Gouédard P, Collins J A, et al. 2011. Structure of young east pacific rise lithosphere from ambient noise correlation analysis of fundamental- and higher-mode scholte-rayleigh waves. Comptes Rendus Geoscience, 343(8-9): 571-583. DOI:10.1016/j.crte.2011.04.004
Zheng S H, Sun X L, Song X D, et al. 2008. Surface wave tomography ofchina from ambient seismic noise correlation. Geochemistry, Geophysics, Geosystems, 9(5): Q05020. DOI:10.1029/2008GC001981
曹锡秋. 2013.新疆某地衰竭气藏地下储气库地应力特征研究[博士论文].北京: 中国地质大学(北京). http://cdmd.cnki.com.cn/Article/CDMD-11415-1013261800.htm
陈华慧, 林秀伦, 关康年, 等. 1994. 新疆天山地区早更新世沉积及其下限. 第四纪研究, 14(1): 38-47. DOI:10.3321/j.issn:1001-7410.1994.01.004
陈立春. 2011.北天山乌鲁木齐转换区构造系晚第四纪活动性[博士论文].北京: 中国地震局地质研究所. http://kns.cnki.net/KCMS/detail/detail.aspx?filename=GJZT201210012&dbname=CJFD&dbcode=CJFQ
邓起东, 冯先岳, 张培震, 等. 1999. 乌鲁木齐山前坳陷逆断裂褶皱带及其形成机制. 地学前缘, 6(4): 191-201. DOI:10.3321/j.issn:1005-2321.1999.04.001
丁国生, 谢萍. 2006. 中国地下储气库现状与发展展望. 天然气工业, 26(6): 111-113. DOI:10.3321/j.issn:1000-0976.2006.06.036
房立华, 吴建平. 2009. 背景噪声频散曲线测定及其在华北地区的应用. 地震学报, 31(5): 544-554. DOI:10.3321/j.issn:0253-3782.2009.05.007
康竹林. 1997. 准噶尔盆地南缘油气勘探前景. 勘探家, 2(4): 31-34.
李杰, 李瑞, 王晓强, 等. 2016. 呼图壁地下储气库部分区域地表垂直形变机理研究. 中国地震, 32(2): 407-416. DOI:10.3969/j.issn.1001-4683.2016.02.024
李学义, 邵雨, 李天明. 2003. 准噶尔盆地南缘三个油气成藏组合研究. 石油勘探与开发, 30(6): 32-34. DOI:10.3321/j.issn:1000-0747.2003.06.009
李延钧, 王廷栋, 张艳云, 等. 2004. 准噶尔盆地南缘天然气成因与成藏解剖. 沉积学报, 22(3): 529-534. DOI:10.3969/j.issn.1000-0550.2004.03.023
李一峰, 李永会, 高奇, 等. 2014. 呼图壁储气库紫泥泉子组紫二砂层组储集层新认识. 新疆石油地质, 35(2): 182-186.
曲国胜, 马宗晋, 张宁, 等. 2008. 准噶尔盆地及周缘断裂构造特征. 新疆石油地质, 29(3): 290-295.
宋元林, 胡新平. 2001. 呼图壁气田紫泥泉子储层特征及综合评价. 西安石油学院学报(自然科学版), 16(4): 15-19. DOI:10.3969/j.issn.1673-064X.2001.04.005
王芳. 2017.利用背景噪声研究新疆呼图壁储气库周边浅层介质结构及其变化[博士论文].北京: 中国地震局地球物理研究所. http://kns.cnki.net/KCMS/detail/detail.aspx?filename=GJZT201801017&dbname=CJFD&dbcode=CJFQ
杨宏伟, 朱瑾. 2012. 准噶尔盆地南缘平原区地下水动态及地下水资源可持续利用对策. 新疆地质, 30(3): 350-354. DOI:10.3969/j.issn.1000-8845.2012.03.022