2. 吉林大学地球探测科学与技术学院, 长春 130026
2. College of Geo-Exploration Science and Technology, Jilin University, Changchun 130026, China
被动源地震勘探,无需人工震源激发,利用背景噪声或者天然地震,节约大量的生产成本,迅速成为国际研究热点.在声波被动源数据重构的研究基础上,Vasconcelos和Snieder(2008)又将地震波干涉技术应用于弹性介质,给出弹性波干涉技术的基本理论,并将其应用于钻孔地震波成像.Gaiser和Vasconcelos(2010)对弹性波地震干涉技术进行深入的研究,并用于处理模拟的OBC(Ocean Bottom Cable)数据,将虚拟源重构到海底,同时获得纵波、横波和转换波信息.Snieder和Larose(2013)提出从实际噪声记录中提取弹性波场.与声波被动源数据相比,弹性波被动源数据不仅含有纵波(P波)信息,还含有横波(S波)信息,在传播过程中还会产生大量的转换波,包含更加丰富的地下地质构造信息.与主动源记录类似,弹性波被动源数据合成的虚拟炮记录中,同样包含大量的多次波(此处所指多次波,为表面相关多次波,下同),影响对虚拟炮记录的后续处理以及最终的地质构造解释和判断.直接对弹性波被动源数据合成的虚拟炮记录进行多次波的预测与匹配减去十分困难.
Van Groenestijn和Verschuur(2009a, 2009b)提出稀疏反演一次波估计方法(Estimating Primaries by Sparse Inversion, EPSI), 并将其应用于海上主动源地震记录.与表面相关多次波去除方法(Surface-Related Multiple Elimination, SRME)类似,均是基于表面相关多次波模型,完全数据驱动,无需知道地下任何先验信息.EPSI方法是在一次波反射系数稀疏的假设条件下,对数据驱动的波场信息进行反演,直接对一次波进行估计,避免从原始数据中预测和匹配减去多次波的过程,取而代之的是一个大规模的反演过程.此后,稀疏反演一次波估计方法又被推广到双检数据(Baardman et al., 2010)和OBC数据(Van Groenestijn and Ross, 2011).Lin和Herrmann(2013)提出L1范数约束的强劲稀疏反演一次波估计(Robust Estimation of Primaries by Sparse Inversion, REPSI),将求取稀疏解作为明确的目标,与传统算法相比,取得了更加准确的一次波响应.同年,冯飞等(2013)在L1范数约束稀疏反演过程中,结合3D Curvelet变换,使数据变得更加稀疏,提高了解的准确性.
Van Groenestijn和Verschuur(2010)根据被动源数据与主动源数据的差别,依照主动源EPSI理论,给出被动源数据稀疏反演一次波估计方法.在数据满足稀疏假设的条件下,直接从被动源数据中获得不含表面相关多次波的一次波响应,避免对虚拟炮记录进行表面相关多次波匹配减去的过程.程浩等(2015a)将REPSI方法应用于脉冲型被动源数据,在稀疏反演过程中,结合Curvelet-wavelet变换,进一步提高解的准确性.同年,提出噪声型被动源数据一次波估计方法(程浩等,2015b),并对比两种被动源数据类型一次波估计的特点(Cheng et al., 2015).
本文首先回顾声波被动源数据一次波估计的基本理论,在此基础上,推导弹性波被动源一次波估计方法的理论过程,并结合L1范数约束最优化求解方法;其次,通过合成的2D数据进行弹性波被动源数据一次波估计方法的验证,同时,利用地震波干涉技术合成弹性波被动源数据的虚拟炮记录;最后,将所得结果与弹性波主动源记录进行对比,验证本文提出方法的有效性与准确性,并给出相应的结论.
1 声波被动源数据EPSI在传统地震勘探中,地震记录通常包含一次波(注意,此处所谓的一次波是指没有经过地表反射,包含层间多次波)和表面相关多次波.以细节简化算子的形式(Berkhout, 1982),一个粗体符号代表一个频率域的叠前数据矩阵,每一列代表一个共炮点道集,每一行代表一个共检波点道集.通过这种数据矩阵的形式,采集到的上行可控源地震数据表示为
(1) |
式中,X0表示一次波反射系数序列;S表示震源子波特性矩阵函数;R表示自由表面反射算子.将一次波脉冲响应X0和震源子波特性矩阵函数S在频率域相乘,则等于一次波P0=X0S.将一次波脉冲响应X0和自由表面反射算子R与总的上行波场P三者在频率域相乘,则等于表面相关多次波M=X0RP.
声波被动源与主动源数据最大的区别,在于接收到的被动源数据中,没有类似于常规数据的一次波响应,而是一个来自于被动源的直达波响应.该直达波响应引起一系列与之相关的多次波响应,表示为下式:
(2) |
式中,Ppas表示在采集中观测到的声波被动源上行波场数据;方程(2)右侧第一项与方程(1)不同,Pdir表示地下声波被动源到检波器的直达波响应;方程(2)右侧第二项与方程(1)形式相同,均是利用一次波反射系序列与表面反射算子和声波被动源总体上行波场三者相乘,来表示由被动源直达波引起的一系列与之相关的多次波响应.
在方程(2)中包含两个未知数,Pdir和X0.在时间域,对X0的稀疏约束,不能同时求解两个未知数.因此,需要一个假设条件对其进行约束.将方程(2)进行数据重排,可得
(3) |
此处,假设Pdir具有最小能量,则得到一个趋于最小的目标函数J:
(4) |
式中,i表示迭代次数;
(5) |
式中,
由于L0范数约束在迭代过程中存在不收敛和不稳定问题,在足够稀疏的条件下,利用L1范数约束代替L0范数约束进行最优化问题的求解,得下式
(6) |
式中,τ为1范数值.求解方程(6)的最优化问题,可以转化为
(7) |
式中,σ表示实测被动源数据与估计数据的差值.程浩和王德利等利用凸优化问题求解方法SPGL1(Spectral Projected Gradient for L1 minimization, 谱梯度投影法)(Hennenfent et al., 2008)求解声波被动源数据稀疏反演一次波估计问题,得到了更加准确的一次波响应.
2 弹性波被动源数据EPSI与声波相比,弹性波被动源数据中,不仅包含P波和S波,在传播过程中还会产生大量转换波,如图 1所示.以数据矩阵的形式,表示为下式:
(8) |
式中,Pe-pas表示在采集过程中观测到的弹性波被动源总上行波场数据,PP-dire-pas、PS-dire-pas、PP-S-dire-pas、PS-P-dire-pas分别表示地下弹性震源激发到检波器的P波直达波响应,S波直达波响应,转换S波直达波响应,转换P波直达波响应;PP-dir-mule-pas和PS-dir-mule-pas分别表示由P波直达波响应和S波直达波响应引起的一系列与之相关的多次波响应,PP-dir-S-mule-pas和PS-dir-P-mule-pas分别表示由P波直达波响应和S波直达波响应产生的转换S波和P波引起的一系列与之相关的多次波响应,PP-S-dir-mule-pas和PS-P-dir-mule-pas分别表示转换S波直达波响应和转换P波直达波响应引起的一系列与之相关的多次波响应,PP-S-dir-P-mule-pas和PS-P-dir-S-mule-pas分别表示转换S波直达波响应和转换P波直达波响产生的转换P波和S波引起的一系列与之相关的多次波响应.
这些由直达波引起的多次波可以表示为下式:
(9) |
(10) |
(11) |
(12) |
(13) |
(14) |
(15) |
(16) |
式中,X为一次波反射系数;R为自由表面算子,表示所有的波近似完全反射,因此公式(9)到(16)中的自由表面算子一致.注意,XPS表示一次波的双程走时中,包含了一半P波和一半S波所组成的一次波反射系数,同理于XSP,如图 2所示.
由图 2可知,式(9)和(12)是由P波直达波响应引起的多次波和转换波多次波,所以在右侧都褶积上的是PP-dir-mule-pas,同理于式(10)和(11)、(13)和(16)、(14)和(15).由于波在传播过程中会产生能量损失,将多次转换的波略掉,并将(9)到(16)式代入到(8)式,可得:
(17) |
(17) 式的相同项合并,可得:
(18) |
由于XP、XS-P、XP-S-P都属于P波一次波反射系数,此处,用
(19) |
则,式(18)简化为
(20) |
(21) |
令
(22) |
因此,得到与声波类似的表达式,进行弹性波被动源数据一次波估计.假设
(23) |
同理,转化为L1范数约束的最优化问题求解,得下式
(24) |
由上述可知,估计出来的一次波反射系数不仅包含P波和S波,还包含转换波的一次波反射系数序列.
3 数据测试一由于弹性波传播的复杂性,首先利用一个二维的水平层状模型进行验证,如图 3所示.模型水平长度为3000 m,垂直深度为1000 m,仅包含一个水平界面,深度200 m,设置第一层的纵波速度为2000 m·s-1,横波速度为1200 m·s-1,密度为1000 kg·m-3,第二层的纵波速度为2500 m·s-1,横波速度为1500 m·s-1,密度为1200 kg·m-3.检波器设置在模型上表面,共300道,均匀分布,道间距为10 m,位置固定.随机震源分布在模型的第二层中,震源位置随机分布,进行激发,模拟弹性波被动源数据.
图 4a为水平层状模型弹性波被动源数据的正演模拟结果,与声波被动源相比,弹性波被动源数据含有大量的横波信息.本文同时模拟了对应的主动源记录,同样为300道,道间距10 m,300炮,炮间距10 m.图 4b为水平层状模型在自由表面条件下的第150炮主动源记录,可以发现,仅含一个界面的水平层状模型的主动源记录非常复杂,其中包括一次波PP、SS,多次波PPPP、SSSS、PPPPPP以及在传播过程中产生的转换波PS、PPPS、PPSS、PSSS、PPPPPS,命名规则如表 1所示.图 4c为水平层状模型在吸收表面条件下的第150炮主动源记录,三个同相轴依次为纵波反射波PP、纵波在分界面处产生的转换横波PS、横波反射波SS.
对比图 4(b, c)可知,图 4b中PPPP为纵波的一阶多次波;PPPS为一阶多次波在分界面处产生的转换横波;PPSS为一次波在自由表面产生的转换横波, 相当于一个纵波双程走时加上一个横波双程走时产生的响应;约0.6 s处两个反射同相轴重合,通过走时分析,PPPPPP为纵波的二阶多次波,PSSS为约0.3 s处的转换横波又经过地下一个双程走时被接收到的横波, 相当于一半的纵波双程走时加上1.5倍的横波双程走时所得的横波响应;约0.68 s处同样为两个反射同相轴的重合,SSSS为横波的一阶多次波,PPPPPS为纵波二阶多次波产生的转换横波.由于深部的同相轴能量较弱,不再做具体分析.综上分析,由于在自由表面处产生的多次波的影响以及多次波在后续传播中产生的转换波,使弹性波记录变得复杂多变,难以进行准确的分析和判断.
通过常规的互相关地震波干涉技术处理弹性波被动源数据,进行被动源数据的波场重构,合成虚拟炮记录.如图 5a所示,为虚拟炮记录的第150炮,与给出的主动源记录对应.与图 4b对比可知,图 5a中虚拟炮记录中的一次波PP、SS,多次波PPPP、SSSS、PPPPPP,以及转换波PS、PPPS、PPSS、PSSS、PPPPPS都得到准确的重构,与主动源记录的同相轴一一对应.和主动源记录相同,弹性波被动源数据合成的虚拟炮记录中仍然存在着多次波的影响以及由多次波产生的转换波.
利用本文提出的弹性波被动源数据一次波估计方法,处理弹性波被动源数据,直接从被动源数据中提取一次波信息.图 5b给出估计的弹性波一次波记录的第150炮,与公式(22)中
由此可见,通过本文所提出的弹性波被动源数据一次波估计方法,水平层状模型的纵波和横波所对应的多次波都能得到有效地衰减;纵波的多次波在地下层界面处产生的转换横波同样被衰减;保留了纵波在自由表面处产生的转换横波以及纵波在层界面处产生的转换波.此外,在反演过程中,产生的干扰噪声,影响了估计一次波的精度.
4 数据测试二在某些实际情况下,会含有由于断点形成的绕射波.利用断层模型讨论被动源一次波估计方法对于绕射波在地表多次反射的消除问题.如图 6所示,断层模型的水平长度和垂直深度为3000 m和1000 m,断层断距为100 m.和测试数据一类似,检波器设置在模型上表面,共300道,均匀分布,道间距为10 m,位置固定.随机震源分布在模型的第二层中,震源位置随机分布,进行激发,用来模拟弹性波被动源数据.
图 7a为利用断层模型模拟的弹性波被动源记录,随机震源点位于模型水平方向中间位置,显示出断层存在的迹象.图 7b和7c分别为利用断层模型模拟的自由表面和吸收表面条件下主动源弹性波记录.图 7b中对不同的波进行标记,表明各种波之间的对应关系.PPPP为纵波一阶多次波,PPPS为纵波一阶多次波在界面处形成的转变波,PPSS为纵波一次波在自由界面处形成的转换波,PPPPPP为纵波二阶多次波,PSSS为转换横波的一阶多次波.由于断点的存在,主动源记录中同样含有绕射波以及与其对应的绕射多次波.
图 8a为利用互相关方法重构的断层模型虚拟炮记录.与自由表面条件下的弹性波主动源记录图 7b对比,纵、横波的一次波、多次波以及转换波同相轴都得到准确地重构.对图 8a中形成的多次波及转换波进行标记,给出相互之间的对应关系,与自由表面弹性波主动源记录一致.图 8b为利用本文提出的方法估计出来的弹性波被动源数据一次波记录.对比图 8(a,b)可见,PPPP、PPPS和PPPPPP得到有效地衰减,转换波PPSS和PSSS作为一次波得到保留.
为了体现弹性波被动源一次波估计方法对绕射多次的压制能力,分别抽取自由表面条件下弹性波主动源记录、互相关记录和估计一次波记录的零偏移距剖面,如图 9所示.对比可知,圈定部位包含的断层右侧纵波二阶多次波对应的绕射多次以及纵波一阶多次波在界面处形成转变波对应的绕射多次都得到了充分地压制,表明本文所提方法同样能够压制绕射多次.
为了测试本文提出方法的广泛适用性,又利用同时含水平层和断层的模型进行验证,如图 10所示.模型水平长度和垂直深度同样为3000 m和1000 m,包含了一个水平层界面,深度200 m,和一个断层构造,断距为100 m.纵波速度、横波速度和密度从上到下依次增大.检波器设置在模型上表面,共300道,均匀分布,道间距为10 m,位置固定.随机震源分布在模型的底层,震源位置随机分布,进行激发,用来模拟弹性波被动源数据.
图 11a为含断层模型弹性波被动源数据的正演模拟结果,与上述水平层状模型相比,被动源数据变得更加复杂.图 11(b,c)分别为对应的自由表面和吸收表面条件下的弹性波主动源单炮记录.对比可知,图 11b比图 11c增加了纵波和横波的表面相关多次波,以及由多次波引起的转换波.图中蓝色字母表示与水平界面有关的地震波,红色字母表示与断层左侧有关的地震波,绿色字母表示与断层右侧有关的地震波.图 11b蓝PPPP为纵波的一阶多次波,它淹没了纵波在断层左侧界面产生的转换横波红PS;蓝PPPS为纵波一阶多次波在水平界面处产生的转换横波,它淹没了断层左侧界面产生的横波反射同相轴红SS;红蓝PPPP为断层左侧界面对应的纵波在第一层中产生的一阶多次波;蓝PSSS为纵波在第一层中产生的转换横波的一阶多次波;红蓝PPSS为断层左侧界面对应的纵波在第一层中产生的转换横波;红PPPP为断层左侧界面对应的纵波在断层左侧界面形成的一阶多次波;红PPPS为断层左侧界面对应的纵波在断层左侧界面形成的一阶多次波引起的转换横波;绿蓝PPPP为断层右侧界面对应的纵波在第一层中产生的一阶多次波;绿蓝PPSS为断层右侧界面纵波在自由界面引起的转换横波.由于其他的转换波和深部同相轴能量相对较弱,此处不再具体分析.
通过上述分析可知,由于自由表面的存在,使弹性波主动源变得十分复杂,产生大量的多次波,以及伴生的转换波,对有效信息的判断产生巨大的影响.
利用地震波干涉技术互相关算法,处理所得的弹性波被动源数据,进行波场重构,合成虚拟炮记录.图 12a给出了合成的虚拟炮记录的第150炮.与自由表面条件下的弹性波主动源记录图 11b进行对比,纵、横波的一次波、多次波以及转换波同相轴都得到有效地重构.蓝PP、红PP、绿PP为重构的纵波一次波,蓝PS、红PS、绿PS为与之对应的转换波,蓝SS、红SS、绿SS为重构的横波一次波,蓝PPPP、红蓝PPPP、绿蓝PPPP、红PPPP为重构的表面相关多次波,蓝PPPS、蓝PSSS、红蓝PPSS、绿蓝PPSS、红PPPS为重构的转换波.可见多次波和由多次波产生的转换波都清晰可见,严重影响对弹性波虚拟炮记录的后续处理和对地下地质构造信息的有效判断.
利用本文提出的一次波估计方法,处理含水平层和断层模型的被动源数据,提取弹性波被动源数据一次波信息.图 12b为估计出来的一次波响应的第150炮,与图 11b和图 12a对比可知,图 11b中纵波一阶多次波蓝PPPP被衰减掉,由断层左侧界面纵波产生的转换横波红PS凸显出来;由纵波一阶多次波在水平界面处产生的转换横波蓝PPPS被衰减掉,断层左侧界面对应的横波红SS得以凸显;断层左侧界面对应的纵波在第一层中产生的一阶多次波红蓝PPPP被衰减掉;纵波在第一层中形成的横波的一阶多次波蓝PSSS,作为一次波被保留;断层左侧界面对应的纵波在第一层中产生的一阶多次波红蓝PPSS被保留;断层左侧界面对应的纵波在断层左侧界面形成的一阶多次波红PPPP被衰减掉,与其对应的转换横波红PPPS同样被衰减掉;断层右侧界面对应的纵波在第一层中产生的一阶多次波绿蓝PPPP被衰减掉;断层右侧界面纵波在第一层中产生的转换横波绿蓝PPSS被保留.与公式(22)中
通过本文所提出的弹性波被动源数据一次波估计方法,含水平层和断层模型的纵波和横波所对应的多次波都得到有效地衰减;纵波的多次波在层界面处产生的转换横波同样被衰减掉;纵波在层界面处产生的转换横波的一阶多次波被作为新的一次波得到保留.同时,在一次波估计的过程中,会产生噪声干扰,影响一次波估计的精度.
6 结论弹性波被动源数据一次波估计方法,直接从弹性波被动源数据中进行一次波信息的提取,避免了对弹性波虚拟炮记录进行多次波衰减的一系列处理.通过水平层状模型和含断层模型,以及同时含水平层和断层的模型的验证,得到以下结论:
(1) 估计出的弹性波被动源数据一次波记录中,纵波和横波的一次波信息得到有效地重构,纵波和横波所对应的表面相关多次波都得到充分地压制;
(2) 在弹性波被动源数据一次波估计的过程中,由纵波的多次波在地下层界面处产生的转换横波,同样可以被衰减掉;
(3) 由纵波在自由界面处,产生的转换波,经过地下反射回自由界面,被作为新的一次波,得到保留;
(4) 由纵波在地下界面处产生的转换横波,又经过地表一次反射,所形成的一阶多次波,同样作为新的一次波得到保留.
综上所述,本文提出的弹性波被动源数据一次波估计方法,有效地衰减了与纵波和横波相关的多次波,以及与之相关的部分转换波.在估计出的一次波记录中,不但包含与地下地质构造相关的纵、横波一次波响应,而且包含一些新的转换波信息.如何对新的转换波信息进行有效的利用,更好地解读地下地质构造信息,将是后续需要解决的问题.本文理论数据的测试,将为实际弹性波被动源数据的处理提供良好的理论基础.
Baardman R H, Verschuur D J, Van Borselen R G, et al. 2010. Estimation of primaries by sparse inversion using dual-sensor data.//80th Annual International Meeting. SEG, 3468-3472.
|
Berkhout A J. 1982. Seismic Migration:Imaging of Acoustic Energy by Wave Field Extrapolation-Part A:Theoretical Aspects. New York: Elsevier.
|
Cheng H, Wang D L, Feng F, et al. 2015. Estimating primaries from passive seismic data. Exploration Geophysics, 46(2): 184-191. DOI:10.1071/EG14079 |
Cheng H, Wang D L, Feng F, et al. 2015. Estimating primaries by sparse inversion of passive-source seismic data with L1-norm constraint. Chinese J. Geophys. (in Chinese), 58(2): 674-684. DOI:10.6038/cjg20150228 |
Cheng H, Wang D L, Feng F, et al. 2015. Estimating primaries from passive seismic data activated by noise sources. Acta Petrolei Sinica (in Chinese), 36(2): 197-202, 223. |
Feng F, Wang D L, Zhu H, et al. 2013. Estimating primaries by sparse inversion of the 3D curvelet transform and the L1-norm constraint. Applied Geophysics (in Chinese), 10(2): 201-209. DOI:10.1007/s11770-013-0378-0 |
Gaiser J E, Vasconcelos I. 2010. Elastic interferometry for ocean bottom cable data:Theory and Examples. Geophysics Prospecting, 58(3): 347-360. DOI:10.1111/(ISSN)1365-2478 |
Hennenfent G, Van Den Berg E, Frendlander M P, et al. 2008. New insights into one-norm solvers from the Pareto curve. Geophysics, 73(4): A23-A26. DOI:10.1190/1.2944169 |
Lin T T Y, Herrmann F J. 2013. Robust estimation of primaries by sparse inversion via one-norm minimization. Geophysics, 78(3): R113-R150. DOI:10.1190/geo2012-0411.1 |
Snieder R, Larose E. 2013. Extracting Earth's elastic wave response from noise measurements. Ann. Rev. Earth Plant. Sci., 41(1): 183-206. DOI:10.1146/annurev-earth-050212-123936 |
Van Groenestijn G J A, Verschuur D J. 2009a. Estimating primaries by sparse inversion and application to near-offset data reconstruction. Geophysics, 74(3): A23-A28. DOI:10.1190/1.3111115 |
Van Groenestijn G J A, Verschuur D J. 2009b. Estimation of primaries and near-offset reconstruction by sparse inversion:marine data application. Geophysics, 74(6): R119-R128. DOI:10.1190/1.3213532 |
Van Groenestijn G J A, Verschuur D J. 2010. Estimation of primaries by sparse inversion from passive seismic data. Geophysics, 75(4): SA61-SA69. DOI:10.1190/1.3460431 |
Van Groenestijn G J A, Ross W. 2011. Primary estimation on OBC data by sparse inversion.//81th Annual International Meeting. SEG, 3531-3535.
|
Vasconcelos I, Snieder R. 2008. Interferometry by deconvolution, Part 2-Theory for elastic waves and applications to drill-bit seismic imaging. Geophysics, 73(3): S129-S141. DOI:10.1190/1.2904985 |
程浩, 王德利, 冯飞, 等. 2015a. L1范数约束被动源数据稀疏反演一次波估计. 地球物理学报, 58(2): 674-684. DOI:10.6038/cjg20150228 |
程浩, 王德利, 冯飞, 等. 2015b. 噪声被动源数据一次波估计. 石油学报, 36(2): 197-202, 223. |
冯飞, 王德利, 朱恒, 等. 2013. 三维曲波变换L1范数约束稀疏反演一次波估计方法研究. 应该地球物理, 10(2): 201-209. |