面波频散谱多模式高分辨率成像的多道信号比较法

易佳， 刘伊克 ， 胡昊， 张郁山， 杨宗奇

1 中国地震灾害防御中心， 北京 100029 2 中国科学院地质与地球物理研究所， 北京 100029 3 University of Houston, Houston Texas, 77204-5007, USA 4 中国石油集团东方地球物理勘探责任有限公司， 河北 涿州 072750

0 引言

在浅层地震勘探领域，由于面波在地震波场能量中占主导部分，而面波携带有丰富的地下介质的横波速度信息，因此横波速度被认为是一种十分重要的地震学参数(Aki and Richards, 1980).在地震工程领域，横波速度是确定场地土类别的重要参数和主要依据，从而用于场地地震反应分析计算，最终服务于场地的地震安全性评价(尤红兵等, 2009; 张郁山, 2009, 2010; Wang et al., 2019).因此，探明浅层精细的横波速度结构对于以上两个领域工作的开展意义重大.

自20世纪50年代面波的频散性质被发现后，在浅层速度结构的分析中，研究者们便利用面波的频散特性反演估计浅层横波速度结构.随着瞬态面波法的兴起和应用，面波频谱分析法(Stokoe and Nazarian, 1983; Stokoe et al., 1994)、多道面波分析法(Park et al., 1998; Xia et al., 1999, 2012)等通过提取面波的频散曲线，采用最小二乘局部优化算法或全局优化方法，反演频散曲线估计一系列的横波速度随深度变化的剖面(Beaty et al., 2002; Luo et al., 2007).通过对一维剖面的插值处理，可重建拟二维的横波速度结构(Boiero and Socco, 2010; Bergamo et al., 2012; Mi et al., 2017).随着波形反演(Pérez Solano et al., 2014; Liu et al., 2018)和波动方程走时反演(Yi et al., 2019)的飞速发展，基于波动方程的频散曲线反演(Li et al., 2017; Liu et al., 2019)和二维频散谱反演方法(Masoni et al., 2014; Zhang and Alkhalifah, 2019)随即发展起来，该类方法基于波形反演和波动方程走时反演的思想，通过建立频散曲线或二维的频散谱相关的目标函数，反传频散曲线或频散谱的残差直接反演获得二维的横波速度结构.

基于面波的频散特性反演获取横波速度结构的方法，均涉及一个关键步骤：对面波的频散能量进行成像，获得面波的频散谱.因此，高分辨率的面波频散谱成像对于准确提取频散曲线十分关键.此外，面波的高阶模式相比于基阶模式对于地层参数、软弱夹层等更加敏感，且高阶模式具有更深的勘探深度，因此，考虑精确的多模式的面波频散谱成像能有效提高反演的精度以及探测深度(Luo et al., 2007; Zhang and Alkhalifah, 2019).在主动源面波勘探中，常用的面波频散谱成像方法有：τ-p变换法(McMechan and Yedlin, 1981)，f-k变换法(Yilmaz, 1987)，相移法(Park et al., 1998)，频率分解倾斜叠加法(Xia et al., 2007)，高分辨率线性拉东变换法(Luo et al., 2008)，以及矢量波数变换法(杨振涛等, 2019).τ-p变换法将数据按不同的截距和斜率的直线进行叠加，对高阶模式成像质量较好，但对于基阶模式的成像分辨率较差，且易出现端点效应和假频问题 (邵广周和李庆春, 2010).f-k变换法对数据进行二维傅里叶变换，原理简单，但要求时间和空间采样等间隔，对采集方式要求较高，不能有空道和坏道，否则将大大影响频散谱的成像精度.相移法将各道数据沿时间方向做傅里叶变换，再沿空间方向进行积分，不要求检波器等间隔，但由于每个频率只能得到一个相位值，单个频率的频谱为基阶和高阶的耦合，因此相移法对于高阶模式不能很好地成像.频率分解倾斜叠加法首先通过频率扫描，得到伪炮集记录，再通过倾斜叠加的方式获得面波的频散谱，该方法在叠加前对频率进行了分解，因此，相比于τ-p变换法有更高的成像精度，但该方法在频率扫描阶段存在许多无效计算及频率偏差的问题(沈超, 2017).高分辨率线性拉东变换法将频散谱成像设计成一个反演过程，通过迭代可以逐渐提高成像的分辨率(Luo et al., 2008).矢量波数变换法通过对观测波场的频谱进行矢量波数变换，在频率波数域对变换后的波场进行扫描获得面波的频散谱.该方法对基阶和高阶模式均有较高的成像分辨率，但需要已知震源子波的信息，因此在实际资料震源未知的情况下，还需要对震源做近似处理(杨振涛等, 2019).

在天然地震探测领域，双台法(Knopoff et al., 1966; Yao et al., 2006; Foster et al., 2014)又被称为LSC方法(Zheng and Hu, 2017)，被广泛用来对面波的频散能量进行成像，提取频散曲线并用于大尺度的面波层析成像.该方法通过对两个台站记录的来自同一方位的面波进行互相关信号比较，即可获得面波的频散谱，但互相关的成像方式对于低频信号极不敏感，会造成频散谱在低频端较低的分辨率，以及整个频段分辨率的不均匀，导致无法准确拾取频散曲线.NLSC方法(Zheng and Hu, 2017)利用指数函数克服了低频端分辨率较低的问题，同时引入了一个可调参数，极大地提高了频散谱的成像分辨率.主动源面波勘探观测系统与天然地震观测系统不同，通常沿着直线布设多个检波器接收地震信号，因此，在信号比较法的基础上，考虑多道的MSC方法(Hu et al., 2019)被提出用于多道地震数据的面波频散谱成像，对于缺道以及噪声问题，MSC方法也具有较强的稳健性，与此同时，MSC方法也被用于地震数据的面波分离当中，取得了良好的分离效果(Hu and Zheng, 2019).

在研究中我们进一步发现，LSC和NLSC方法利用两道地震记录，并不能将面波的频散特性完整地考虑在内，仅对基阶模式成像准确，高阶模式的成像存在较大的误差.MSC方法考虑了完整的面波频散信息，能获得准确且分辨率较高的面波频散谱多模式成像，但该方法需要计算任意两道的频散谱并进行叠加，对于多道记录，存在冗余的计算耗时较长.因此，我们对MSC方法进行了相应的改进，以求高效获得频散谱的多模式高分辨率精确成像，从而有助于多模式频散曲线的精确提取，以及利用多模式频散谱进行高精度的浅层地震勘探.本文系统性地对LSC、NLSC方法的理论、关键技术以及存在的问题进行了阐述分析，并给出了改进的MSC方法的原理.通过对含低速夹层的水平层状模型及一套实际地震资料的Rayleigh波频散谱的成像，深入对比研究了几种基于信号比较理论方法在面波频散谱成像中的应用效果，从而为不同的面波频散谱分析研究提供合适的解决方案.

1 基本原理

1.1 传统线性信号比较法 (LSC)

在天然地震面波频散性质的分析当中，LSC方法利用互相关的方式比较两个台站接收到的地震信号，获得它们之间的相位差，这个相位差与两个接收点之间的距离及面波的相速度有关，通过这种方式即可获得面波的频散谱.其方法原理如图1所示.震源激发的地震波，传播经过同一个方位的两个台站Rec 1和Rec 2被接收到，两者之间的距离为x.Rec 1处的地震记录表示为R1(t)，Rec 2处的地震记录表示为R2(t).如果能测得Rec 2处特定频率的面波到达的延迟时τ，那么就可以反推获得面波的相速度v.

图1 传统线性信号比较法原理示意图Fig.1 The conventional LSC method

然而，通常两个台站之间的距离是已知的，因此，只需给出一定的相速度范围，移动地震记录R2(t)，求取与R1(t)的互相关，即可获得面波的频散谱(Zheng and Hu, 2017)，记为ELSC(ω,v)：

(1)

1.2 非线性信号比较法 (NLSC)

为克服传统LSC方法分辨率的限制，并保证整个频带范围内均匀分布的分辨率.Zheng和Hu (2017)提出了NLSC方法：

(2)

其中，ENL表示未归一化的频散谱，ε∈[0,+∞)是一个非负实数，用来控制频散谱的整体分辨率.Zheng和Hu(2017)定义了一个归一化算子Eπ：

Eπ(ω,v;ε)=

(3)

对ENL进行归一化处理，即可得到归一化的高分辨率频散谱，记为ENLSC：

(4)

归一化算子Eπ可认为是两个台站记录的地震信号的相位差为π时，由(2)式计算得到.NLSC方法巧妙地运用指数函数的性质，克服了传统LSC方法生成的频散谱在整个频带范围分布不均匀的问题，提高了低频端的分辨率，同时引入参数ε可以有效提高频散谱的整体成像分辨率.

1.3 多道信号比较法(MSC)

由于探测精度的需求，主动源面波勘探观测系统与天然地震观测系统不同，通常沿着直线布设多个检波器接收地震信号，从而可利用多道信息反演获得小尺度的更为精细的地下结构.在研究中我们发现，只利用两道地震记录的LSC和NLSC方法，并不能将面波的频散特性完整地考虑在内，只对基阶模式有较准确的成像，对于高阶模式的成像存在较大的误差.因此，Hu等(2019)基于勘探地震中多道采集的方式提出了MSC方法：

(5)

其中，i,j表示地震道，n表示总的接收道数，Ei,j表示利用信号比较法计算的任意两道的频散谱，E′MSC表示将所有两道频散谱叠加的频散谱，采用归一化方法，对频散谱进行归一化处理，即可得到MSC方法计算的面波频散谱：

EMSC(ω,v)=

(6)

图2 改进的多道信号比较法原理示意图Fig.2 The modified multichannel signal comparison method

(7)

同理，改进的多道非线性信号比较法(MNLSC)的可表示为

(8)

改进的MSC方法追踪炮集记录上的面波波组，考虑了完整的面波频散信息同时不增加冗余的信息，可以有效提高计算效率.

2 数值实验

2.1 传统LSC和NLSC方法的频散谱成像分析

为了研究几种信号比较法在浅层地震勘探中的成像能力，采用含低速夹层的水平层状模型(LVL)进数值实验，模型参数如表1所示.利用垂直震源激发地震波，子波是主频为20 Hz的Ricker子波，自由边界条件，采用SPECFEM2D(Komatitsch and Tromp, 1999)模拟含有Rayleigh波的地震波场.线性观测系统，100个检波点均匀分布在地表，道间距1 m，采样间隔0.1 ms，一共记录时长为1 s.获得的炮集地震记录的垂直分量如图3所示，可以看到地震记录中Rayleigh波的能量占了主导成分，且可观测到明显的面波的频散特性.

表1 含低速夹层的水平层状模型参数Table 1 Parameters of the LVL model

图3 含低速夹层的水平层状模型的炮集地震记录的垂直位移分量Fig.3 The vertical displacement component seismograms of the LVL model

抽取其中的两道地震记录(如图3中红色地震道所示)，利用LSC和NLSC方法进行面波频散谱的成像计算.此外，还利用Knopoff算法(Schwab and Knopoff, 1970)计算得到该模型的理论频散曲线叠加在频散谱上，用以判断数值方法计算的频散谱的准确性.所得的Rayleigh波频散谱能量图如图4所示，可以看到LSC方法(图4a)的频散谱在给出的频带范围内呈现出非均匀的分辨率，高频端的分辨率较高，低频端分辨率较低，在小于10 Hz区域成像的分辨率尤其较低.而NLSC方法(图4c，ε=0.01时)成像的分辨率得到了有效地提高，尤其在低频端得到了显著地提升.当调整分辨率的参数ε=0.001时(图4d)，成像分辨率达到极高.在频散谱的准确性上，LSC和NLSC方法均呈现出强而清晰的基阶模式Rayleigh波频散能量，并且与理论频散曲线的基阶模式有较好地吻合.然而，两种方法的成像结果，在高频端均出现了扭曲，与理论频散曲线吻合较差，且高阶模式虽有着清晰的成像，却与理论频散曲线相差甚远.因此，只利用两道地震记录的信号比较法，无法对面波的多模式进行准确地成像，因此无法利用其进行高精度的浅层结构反演.

图4 LSC和NLSC方法的Rayleigh波频散谱能量图(a) LSC方法； (b) NLSC方法,ε=0.01； (c) NLSC方法,ε=0.001，图中白色圆点为Knopoff算法计算获得的理论频散曲线.Fig.4 The Rayleigh wave dispersion spectrum maps of the LVL model using (a) LSC method an (b) NLSC method with ε=0.01 and (c) NLSC method with ε=0.001. The white dots in the figure are the theoretical dispersion curves calculated by the Knopoff′s algorithm

2.2 改进的MSC方法频散谱成像分析

分别采用24和100道的地震数据，应用改进的MLSC和MNLSC方法对该地震记录(图3)进行频散谱的成像计算，为了对比展示成像效果，利用相移法也获得了面波的频散谱，将三种方法与理论频散曲线叠加，结果如图5和图6所示.可以看到，由于低速夹层的存在，三种方法均呈现出了清晰的多模式面波能量.图5显示，当参与计算的道数较少时，频散谱的成像效果均较差，相移法在低频端与理论频散曲线吻合较好，但分辨率较低，而MLSC和MNLSC方法的成像分辨率以及高阶模式的成像精度均明显高于相移法.当参与频散谱成像的地震道数增加时(图6)，三种方法的成像质量都获得了显著地提升，表明面波频散谱的分辨率和准确性会随着道数的增加而提高.图6a显示，相移法的成像分辨率在整个频段上分布极不均匀，低频端(<10 Hz)的成像分辨率较低，使得低频端能量峰值与理论频散曲线出现较大偏差.与相移法相比，MLSC(图6b)和MNLSC方法(图6c，6d)的成像质量具有明显优势，整个频段的分辨率均有所提升，尤其在低频端有显著的提高，与理论频散曲线吻合较好.MNLSC方法相比于MLSC方法在成像分辨率上则有更近一步的提升.因此，改进的MLSC和MNLSC方法，利用多道地震信号，在提升计算效率的同时，不仅能有效地对面波的多阶模式进行精确地成像，而且能获得均匀且分辨率较高的频谱，这在面波频散谱的成像上具有巨大的潜力.

图5 24道数据的Rayleigh波频散谱能量图(a) 相移法； (b) MLSC方法； (c) MNLSC方法,ε=0.03； (d) MNLSC方法,ε=0.01.Fig.5 The Rayleigh wave dispersion spectrum maps with 24 traces using (a) phase shift method (b) MLSC method (c) MNLSC method with ε=0.03 and (d) MNLSC method with ε=0.01

图6 100道数据的Rayleigh波频散谱能量图(a) 相移法； (b) MLSC方法； (c) MNLSC方法,ε=0.03； (d) MNLSC方法,ε=0.01.Fig.6 The Rayleigh wave dispersion spectrum maps with 100 traces using (a) phase shift method (b) MLSC method (c) MNLSC method with ε=0.03 and (d) MNLSC method with ε=0.01

3 实际地震资料测试

为研究改进的MSC方法在实际地震资料中的适用性，将其应用在一套红海海岸跨越Qademah断层采集的含Rayleigh波的实际地震资料中(Hanafy et al., 2015; Li et al., 2017; Zhang and Alkhalifah, 2019).采集参数为：垂直震源，120个检波器，道间距5 m，采样间隔 0.05 ms, 一共记录时长1.2 s.图7a展示了该资料沿着时间方向进行能量补偿(Claerbout, 1985)之后的典型单炮记录，我们对其进行了包括剔除坏道，开窗切除保留面波成分等预处理工作，预处理之后的炮集记录如图7b所示.由于该炮集记录远偏移距的地震道未能记录到完整的面波，因此，我们仅利用预处理之后的前64道进行多道面波频散谱的计算.

图7 实际资料的典型单炮地震记录(a) 能量补偿后的单炮记录； (b) 预处理之后的单炮记录.Fig.7 The common shot gathers of the field data after (a) energy compensation, and (b) after preprocessing

分别采用相移法、改进的MLSC和MNLSC方法对该实际地震资料进行Rayleigh波的频散谱成像计算，结果如图8所示，可以看到MLSC(图8b)和MNLSC(图8c和图8d)频散谱成像的分辨率，在高频端和低频端，均优于相移法(图8a).此外，对比图中红色箭头处，MLSC和MNLSC相比于相移法对高阶模式有更清晰的成像，对比白色箭头处，MLSC和MNLSC则清晰地指示了模式的分离.相比于相移法，MLSC和MNLSC方法频散谱的成像细节，能更好地指示地下横波速度的改变.

图8 实际地震资料Rayleigh波频散谱能量图(a) 相移法； (b) MLSC方法； (c) MNLSC方法,ε=0.05； (d) MNLSC方法,ε=0.02.Fig.8 The Rayleigh wave dispersion spectrum maps of the field data using (a) phase shift method (b) MLSC method (c) MNLSC method with ε=0.05 and (d) MNLSC method with ε=0.02.

4 讨论

在数值试验和实际地震资料的测试中我们可以看到，改进的MSC方法不需要计算每两道地震信号的频散谱，提高了计算效率，但同时充分利用接收的地震信息，获得了准确的多模式高分辨率的频散谱成像.在抗噪性方面，从方法原理可以看出，改进的MSC方法的成像效果较为依赖参考道的数据质量，若参考道数据质量较差，则会对成像结果有较大的影响，因此，在实际资料的处理中，可在近偏移距中选取数据质量较高的地震道作为参考道进行成像计算.对于MNLSC方法中的参数ε可依据反演需求来设置：当采用直接基于面波频散谱的反演方法获得二维横波速度结构时(Zhang and Alkhalifah, 2019)，可将ε设置为一个稍大的数值，或者直接采用MLSC方法来获得频散能量清晰且分布均匀的频散谱；当采用传统提取频散曲线的方法，估计一维横波速度(Xia et al., 1999; Boiero and Socco, 2010) 或基于波动方程走时反演理论，利用频散曲线直接反演获得二维速度结构时(Li et al., 2017)，可将ε设置为一个较小的数值，以获得高分辨率的频散谱.因此，改进的MSC方法的另一个优势是，可依据反演需求获得合适的面波频散谱的分辨率.理论上，当ε趋于零时，MNLSC方法的频散谱分辨率将达到极致，但当ε取值较小时，低频端的成像质量也会有所损失，如图9所示.因此，在利用MNLSC方法计算时，合理地对参数ε取值，对于频散谱的成像质量也十分关键.

图9 LVL模型MNLSC方法 (ε=0.002)Rayleigh波频散谱能量图Fig.9 The Rayleigh wave dispersion spectrum maps of the LVL model using MNLSC method with ε=0.002

5 结论

天然地震探测中的两道信号比较法(LSC方法和NLSC方法)在面波频散谱的计算中存在一定的局限性：低频端较低的分辨率，整个频段非均匀的分辨率，高阶模式不准确的成像，这些问题导致无法准确拾取频散曲线，以及无法有效利用面波的高阶模式反演获得较为准确的横波速度结构.本文根据主动源面波勘探中多道采集的方式，提出的改进的多道信号比较法(MSC方法)通过追踪面波波组，考虑了完整的面波频散信息的同时不增加冗余的计算，可以获得高精度、高分辨率的面波多模式频散谱的成像，该方法原理简单，易于实现，并且理论上不需要观测系统等间隔排列，降低了实际数据采集时对场地的要求.因此，在浅地表地震勘探以及地震工程的横波速度结构探测中，多道信号比较法将具有广阔的应用前景.

致谢感谢两位评审专家提出的宝贵意见，感谢沙特阿卜杜拉国王科技大学Schuster教授及其课题组成员在采集地震数据过程中的辛苦奉献，感谢Zhengdong Zhang博士提供的实际地震数据以及相关讨论.