基于微动方法研究五大连池火山区尾山火山锥浅层剪切波速度结构

为了研究五大连池火山区尾山火山锥浅层三维波速结构特征，在尾山火山锥附近区域布设了无线地震检波器密集台阵，记录连续地震背景噪声波形数据。基于微动方法（拓展空间自相关方法）提取了台站间2~5Hz频率范围的Rayleigh面波相速度频散曲线。利用面波层析成像方法反演获得2~5Hz Rayleigh面波二维相速度图像，基于每一个网格节点的频散曲线，进一步反演获得了尾山火山锥附近区域地表至700m深度的三维剪切波速度结构。成像结果显示：在0~150m较浅深度，靠近尾山火山锥区域显示为相对高速异常，远离火山锥区域则显示为相对低速异常。而至150~700m较深深度，波速异常特征与浅部相反，靠近尾山火山锥的区域显示为相对低速异常，远离火山锥的区域显示为相对高速异常。在远离尾山火山锥区域，浅层的相对低速异常可能与松散沉积层有关，深部的高速异常则反映了结晶变质岩的影响，在靠近尾山火山锥区域，浅部的相对高速异常应该反映了出露地表的玄武岩，而深部的相对低速异常则可能反映了火山通道周围广泛发育的破碎裂隙结构及其火山喷发后孔隙流体填充的影响。

数十米至数百米深度的地壳浅层剪切波速度结构可为地震工程场地评估、强地面运动模拟、活动断裂、地热资源勘查等提供重要的参考。对于具有较厚沉积层地区而言，可靠的浅层波速结构可为地震定位、地壳上地幔深部结构成像等研究提供可靠的浅层结构模型（Xiaetal，1999；何正勤等，2007；徐佩芬等，2009，2013；Xuetal，2012；Lietal，2014，2016；Nietal，2014；李志伟等，2015），因此，地壳浅层剪切波速度结构对于大型工程建设和基础研究均具有重要意义。目前地壳浅层剪切波结构主要的研究包括主动源和被动源方法。主动源方法主要为采用人工震源（如爆炸震源，气枪，重锤，震源车）激发地震面波信号，提取面波频散信息并反演剪切波速度结构。主动源方法具有信噪比高、探测精度高等优点（Xiaetal，1999；Wangetal，2016）然而，主动源探测方法的成本较高，同时基于炸药震源的方法对环境有一定破坏性，在城市等人口密集地区开展探测工作受到一定限制。被动源方法不依赖于人工震源，而是利用包括天然地震事件、地震背景噪声等在内的被动震源产生的信号探测浅层结构。例如，利用地震台站记录到的近震事件P波波形，基于初至P波垂向与径向振幅比，以及来自沉积基底的Ps、Sp转换震相的时差和振幅比，可以获得台站下方数十米至数百米深度的浅层剪切波速度结构（Lietal，2014；Nietal，2014）。

地球上时刻存在着由海洋、大气、地震等各种震源激发的地震波，携带了噪声源及地球内部结构的丰富信息。随着地震学理论和方法的发展，基于地震背景噪声数据的研究在地壳波速结构成像、地球内部波速结构随时间变化、地震水平位置测定等方面取得了丰富的成果（Shapiroetal，2005；Stehlyetal，2006；Zhanetal，2011）通过提取地震台阵记录的地震背景噪声数据中的面波信号，可以获得面波频散并反演地壳及上地幔的剪切波速度结构。而在地壳浅层结构研究中，也发展了包括空间自相关（Spatial Autocorrelation，SPAC）、拓展空间自相关（ESPAC）、频率波数法（F-K）、H/V谱比法等基于地震背景噪声数据的微动方法。在地壳浅部数十米至数百米深度的波速结构探测中，基于台阵的微动探测方法已在地热资源、重大工程场地等研究领域获得了广泛应用，取得了良好的应用效果（Tokszand Lacoss，1968；Ling，1994；Okada，2003；Xuetal，2012）。基于地震背景噪声的地壳浅层结构成像方法不受地震事件时空分布的限制，施工成本较低，对环境无破坏，因此已经取得了较快发展，成为研究地壳浅层波速结构的重要地震学方法。

Aki（1957）最早提出了SPAC方法的概念，并给出了平稳随机波动理论假设下一维和二维的空间自相关理论和应用实例。Aki（1965）,Asten（1976）,Ferrazzini等（1991）,Chouet等（1998）,Satoh等（2001）又对微动探测的原理和技术方法进行了研究，推动了微动方法的发展。何正勤等（2007）利用微动方法研究了0~3km深度的剪切波速度结构，为确定地热钻井位置提供了依据，并得到了后期钻孔资料的验证。徐佩芬等（2009，2013）基于微动方法开展了煤矿陷落柱和地铁沿线地下“孤石”等浅层波速结构相关的探测研究。

除了空间自相关方法外，还有一些基于地震背景噪声的地壳浅层结构探测方法，Capon（1970）利用最大似然F-K法成功地从地震背景噪声中提取了面波信号。Nakamura（1989）发展出新的台阵数据处理方法，即H/V谱比法研究场地的卓越频率及放大因子等参数，用于估计沉积层厚度。Claerbout（1968）提出了地震背景噪声互相关理论的雏形，证明了在水平层状且底层为均匀无限半空间模型中，自由地表记录的底部透射地震波自相关波形与其自激自收波形是等价的。Xu等（2013）在三峡库区巴东县黄土坡滑坡体隧道内用47个Texans 4.5Hz检波器组成的测线记录了28个小时微动数据，利用波形互相关方法得到了3~8Hz面波信号，并采用多道面波分析方法（MASW）提取相速度频散曲线，最后反演获得滑坡体200m深度内的剪切波速度结构。徐义贤和罗银河（2015）详细总结了噪声地震学在浅层勘探中的应用现状，并对噪声地震学的发展前景进行了展望。上述研究进一步表明微动方法在地壳浅层结构探测中的良好应用效果。

五大连池火山是我国最重要的活火山之一（其它重要的活火山包括腾冲火山、长白山火山及镜泊湖火山），其最近一次喷发发生在约300年前，喷发的火山口位于尾山西南侧的老黑山和火烧山。尾山位于五大连池火山群的北部，靠近整个五大连池火山区的中间区域。准确可靠的地壳浅表层剪切波结构对于火山地震准确定位和火山下方地壳深部三维结构成像具有重要意义，能够提高深部结构成像的可靠性。我们在五大连池火山区尾山火山锥附近区域布设了由29个地震检波器组成的密集台阵，对记录的连续地震背景噪声波形数据，使用微动方法提取Rayleigh面波相速度频散曲线，并采用面波层析成像方法反演了研究区地壳浅层剪切波速度结构。在此基础上，我们对五大连池火山区尾山火山锥附近区域的波速异常结构及其相关地质意义进行了初步讨论。

基于空间平稳随机波动理论和基阶面波是微动信号中的主要能量两个假设下，Aki（1957）对空间自相关理论进行了研究．在规则台阵的观测系统下，设圆心处的垂向分量波场为u（x,y,t）,而圆周上任一台站的垂向分量波场为u（x+ξ,y+η,t）,则空间自相关函数φ（ξ,η,t）可表示为：

上划线表示时间平均．假设波场以单一速度c在圆内传播，则空间自相关函数的方位平均为

上式为极坐标表达形式，ψ为圆形观测台阵中半径为r处台站与中心点台站的方位角，ξ=rcos（ψ）,η=rsin（ψ）.空间自相关函数方位平均与波场功率谱密度关系可表示为

Φ（ω）为波场的功率谱密度，ω为角频率，J0为第一类零阶贝塞尔函数，对上述信号应用窄带滤波，功率谱密度可表示为

P（ω0）为在频率为ω0时的功率谱密度，δ（ω-ω0）为狄拉克函数，此时方位平均的空间自相关函数可表示为

自相关系数定义为

由于功率谱密度与台站位置无关，因此可将方位平均后的自相关系数写为

Chávez-García等（2005）研究表明：长时间的地震时间序列可以替代台站的空间方位平均．即对于两个台站来说，通过延长观测记录时间，可以达到噪声源趋向于空间均匀分布的效果，因此长时间的观测与台站的空间方位平均可以达到相似的效果．对于台站间长时间的微动记录，自相关系数按照如下公式进行计算：

ρij（f）是自相关系数，是频率f的函数；Ui（f）与Uj（f）表示台站i与台站j的频率域微动数据，＊表示共轭，ω0为角频率，r为台站i与台站j间的距离，c（ω0）为Rayleigh面波相速度.

Asten（2006）提出当台站在圆周上充分多或为各向同性波场时，抑或台站间记录了足够长时间的微动信号时，自相关系数谱的虚部为零，其实部为第一类零阶贝塞尔函数．在实际计算中，虚部一般不为零，表明噪声源空间分布并非绝对均匀，或者地震噪声波场不是稳定不变的平面波，这会在自相关系数谱计算中引入一定误差．一般来说，将自相关系数谱的实部近似为第一类零阶贝塞尔函数，即

常规SPAC方法通常将观测系统布设成同心圆，以减小噪声源方位分布不均匀对面波频散曲线可靠性的影响．为了能够使空间自相关方法在较大区域开展研究，Ling和Okada（1993）发展出了基于不规则台阵的拓展空间自相关方法（ExtendSPAC,即ESPAC），ESPAC方法布置观测系统比较灵活，可以布设成一字型、十字型或呈网状的任意形状，并通过两两台站间的数据处理获得台站间的面波频散信息．随着地震背景噪声方法理论的发展，以及对于台阵布设无严格限制的优点，基于不规则台阵的ESPAC方法在地壳浅层结构研究得到广泛应用.Ohori等（2002）利用ESPAC和FK方法处理了T型台阵数据，并对获得的Rayleigh相速度频散曲线进行对比分析，发现基于不规则台阵ESPAC方法得到的Rayleigh面波相速度频散曲线更为可靠。Chávez-García等（2005,2006）和Ekstrm等（2009）研究表明在噪声来源方向随机分布的情况下，由两个台站获取的空间自相关函数为贝塞尔函数，与同心圆台站得到的结果相同，即台站对间长时间的地震时间序列可以替代Aki（1957）提出的规则台阵空间方位平均；而在新西兰开展的多次不规则台阵观测也表明反演得到的浅层剪切波速度结构与已有研究成果一致．不仅如此，Chávez-García和Rodríguez（2007）基于不规则台阵背景噪声数据的研究表明：频率域互相关（SPAC）与时间域波形互相关方法得到的Rayleigh相速度频散曲线是一致的.Tsai和Moschetti（2010）从理论上也证明了基于频率域互相关的SPAC方法与时间域波形互相关方法的一致性.Nakata等（2011）和何正勤等（2013）利用ESPAC方法从线性台阵记录的地震背景噪声中恢复了面波信号，测量了面波相速度频散曲线，发现与主动源多道面波分析方法的结果相比具有较好一致性．上述研究表明：基于不规则台阵的微动方法可以有效探测地壳浅层剪切波速度结构。

1、 微动观测系统

五大连池火山群由14个火山锥组成，实心正三角为各火山口位置，尾山位于该火山群最北侧，整体呈北西－南东走向，其东侧及南侧为五大连池堰塞湖（图1a）.研究区域广布农田和树林，高程在320~520m之间，西侧和北侧海拔较高，东南侧较低，整体地势较为平坦．尾山火山锥高程520m左右．我们将29个主频为2Hz的地震检波器组成密集观测台阵，布设在五大连池火山区尾山火山锥附近（图1b）.由于尾山火山口附近森林密布、地形复杂，地震仪器主要布设在尾山火山口的西侧．卢川等（2015）利用相同类型的地震检波器开展了沉积盆地的微动探测并取得了良好效果，表明此类地震检波器适用于微动探测研究．台站位置及时间服务则依据仪器内置GPS获得，采样率设为200Hz.仪器采用内置锂电池供电，每个仪器都记录了至少15h微动信号．由于白天人文干扰较大，因此挑选出2014年7月7号晚上6点至7月8号早上9点共15h数据进行进一步处理.

2、 数据处理

在计算台站对空间自相关系数前，我们分别对单个台站的原始数据进行了预处理，包括去均值、去倾斜、去尖脉冲等操作，并对数据进行频谱分析（图2）。数据振幅谱在1~5Hz频段内能量较强，带通滤波窗口选为0.1~10Hz，将每个台站垂向分量微动数据按60s分段，共计900段．首先对每段数据进行傅里叶变换，按公式（9）用一个台站的傅里叶谱与另一个台站的谱的共轭计算自相关系数，并按公式（10）取实部后进行叠加．然后采用滑动绝对平均方法对叠加后的自相关系数波形进行平滑处理，得到每个台站对的自相关系数曲线，并剔除了信噪比较低、未能呈现第一类零阶贝塞尔曲线形态的自相关系数曲线．从不同台站距的自相关系数波形可以看出，自相关系数曲线与第一类零阶贝塞尔曲线类似，且随着台站间距的增大，第一个零点对应的频率呈下降趋势（图2c）。

通过分析自相关系数波形的零点和极点，并与第一类零阶贝塞尔函数零点和极点进行拟合，得到了Rayleigh面波相速度频散曲线（图3）.即已知第一类零阶贝塞尔函数y=J0（x）的零点及极点对应的x值，由相关系数曲线

的零点及极点可得其对应频率f,由公式即可计算得到每个台站对的Rayleigh相速度频散曲线c（f），最终经过质量控制和挑选，共得到60条较高质量的Rayleigh波相速度频散曲线（图3）。

本文采用面波层析成像方法对2~5Hz的相速度频散数据进行了层析成像反演（Yanovskayaand Ditmar,1990）.面波射线路径分布显示频散数据对研究区有较好覆盖（图4）.反演采用0.1km×0.1km网格，并经过多次成像测试确定合适的正则化参数和阻尼系数，从而保证成像结果相对平滑且能够较好拟合观测数据（图5）.为了更清晰地分析层析成像结果随空间和周期的变化，以尾山火山口为原点O,分别沿西向OA,西南向OB,南向OC以及尾山南侧东西走向DD′（见图1b）获得4组相速度频散曲线变化剖面图（图6）.在研究区西侧、西南侧及南侧靠近尾山火山口的区域，在2.8~4.0Hz频段内相速度相对较高。

面波相速度频散数据不同周期的敏感核测试表明：2.0~3.0Hz频散数据对0~700m结构均有一定分辨能力，对100~400m结构较为敏感，其中2Hz频散对100~500m较浅深度具有较好分辨能力（图3）.为了获得研究区三维剪切波速度结构，将每个网格的相速度频散数据采用线性方法进行反演（Herrmannand Ammon,2004）,得到1-D随深度变化的剪切波速度结构，并通过插值平滑便可获得研究区0~700m的三维剪切波速度结构（图7和8）.100m深度的剪切波速度图像显示：剪切波波速变化范围为600~950m·s-1,在靠近尾山火山锥和研究区西北部为高速异常，研究区东南部则为低速异常；至200m以下深度，靠近尾山火山锥区域和研究区西北部则显示为低速异常，研究区东南部为高速异常．为了更好反映研究区波速结构变化特征，图8和图9分别给出了绝对和相对剪切波速度变化剖面．图8给出了由尾山火山口向西1.7km（OA）,向东南2.0km（OB）,向南1km（OC）及尾山南侧东西方向长度为2km的测线DD′等四组剖面，清晰地展示了0~700m深度绝对剪切波速度变化特征.OA,OB,OC三组由火山口为原点的测线海拔由480m逐渐降为360m,尾山近处100m以上浅层剪切波速度较远处高，并呈现出地表低速覆盖层越靠近尾山越薄的特征.DD′测线由东向西海拔由400m降至360m,整体较为平坦，同样地，100m深度以上的浅层，靠近尾山火山锥处的剪切波速度要较远处高，并呈现出低速覆盖层越靠近尾山越薄的特征．图9同样反映了OA,OB,OC,DD′四组剖面在150m深度以上的浅层地表尾山火山锥近处速度较远处偏高，而在150~700m深度范围内，尾山火山锥近处速度较远处偏低。

鉴于面波频散数据覆盖范围的限制，本文只展示和讨论了研究区主要的异常结构特征，而更小尺度异常特征可能具有较大不确定性不进行讨论．野外地质调查发现在尾山火山锥附近区域地表覆盖有松散沉积土壤层，且靠近火山锥附近广泛出露有多孔隙玄武岩．因此，推测尾山火山锥附近的剪切波速度异常特征可能由以下因素造成：尾山为整体走向为北西－南东，尾山火山口近处及研究区西北部地势较高，有黑色多孔隙玄武岩出露，东南部地势较低，地表分布为松散沉积层，推测靠近尾山火山锥处及西北部深部为火山喷发的玄武岩，在波速结构成像结果中显示为高速异常；而在远离尾山火山锥处，其浅部的低速异常反映了地表松散沉积层的影响；并且越靠近尾山火山锥，松散沉积层的厚度越薄（图7—9）.至200m以下深度，剪切波速度分布特征与火山喷发时火山口附近的岩浆分布密切相关：在靠近尾山火山锥区域，火山作用造成火山通道周围较为破碎多裂隙发育，加之火山岩为多孔隙玄武岩富含孔隙流体，使之显示为低速异常；而在远离尾山火山锥区域，其异常特征可能反映了具有较高剪切波速度的结晶变质岩的影响。

本文利用由29台检波器组成的无线地震检波器密集台阵，在五大连池火山区尾山火山锥附近开展了连续地震背景噪声观测，并基于微动和面波层析成像方法获得了研究区地表至700m深度的地壳浅层三维剪切波速度结构．五大连池火山区尾山火山锥附近的三维剪切波速度成像结果显示：在0~150m的浅部靠近尾山火山锥的区域显示为相对高速异常，远离火山锥的区域显示为相对低速异常．至150~700m深部，波速异常特征则明显与0~150m深度的波速特征相反，即靠近尾山火山锥的区域显示为相对低速异常，远离火山锥的区域显示为相对高速异常．在远离尾山火山锥区域，浅层的相对低速异常可能与松散沉积层有关，深部的高速异常则反映了结晶变质岩的影响；在靠近尾山火山锥区域，浅部的相对高速异常应该反映了出露地表的玄武岩，而深部的相对低速异常则可能反映了火山通道周围广泛发育的破碎裂隙结构及其火山喷发后孔隙流体填充的影响。