《叠前地震数据特征波场成像与层析反演.docx》由会员分享,可在线阅读,更多相关《叠前地震数据特征波场成像与层析反演.docx(7页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、叠前地震数据特征波场成像与层析反演1局部平面波的成像与成像地震资料包含丰富有效的信号,但时间间隔内的信号表达方法无法直接描述地震资料的特 征。因此,可以通过特定的数学变换来提取地震信号的主要特征,进行地震数据的稀疏表 示,这有助于地震数据处理、波场分析和反演图像。事实上,地震数据可以在转换区域中 稀疏表示。常用的转换区域包括foyer、小波场、curvelet rando gaor等。一般来说, 转换域中的地震数据可以是稀疏的,数据在转换区域中的稀疏性可以是稀疏的。利用数据 在转换区域中的稀疏性,可以压缩和重建地震数据,并具有插值的能力。xu等人(2005) 反演不规则地震数据的频谱,并进行数
2、据重建。将震源和卡盘上的地震数据稀疏度(2008) 和白兰淑(2014)用于地震数据重建。然而,Curvelet变换无法直观描述地震波的传播特征,因而没有明确的地球物理含义.本 文中,我们把叠前地震数据抽象为不同出射角度(即出射慢度矢量)以及不同到达时的局 部平面波的线性叠加,而接收到的地震记录可以看作局部平面波波前在空间离散的检波器 处的采样.在立体层析理论(Billette和LambarG, 1998)中,地震记录被简化为炮点、 检波点坐标,炮点、检波点处的射线参数以及射线的双程走时.上述参数构成了立体层析 中的数据空间,但忽略了地震子波的波形信息.若用Gabor框架作为基函数,则可以地震
3、数据精确重构,cek (2003)用相干态变换 (coherentstate transform)方法构造一组完备的Gabor框架函数,实现了地震数据的高 斯波包分解,由于Gabor框架函数的冗余度高(Sondergaard, 2007),该方法无法实现数 据的稀疏表达.Geng (2011)忽略了高斯波包中一些不重要的参数,压缩了重构观测波场 的高斯波包泛函空间,从而减少框架函数的数乱李辉等(2014)提出了地震数据的特征 高斯波包(CGP)分解方法,在Radon域选取出地震数据的时空和方向特征,大幅减少了 高斯波包框架的数量,实现了 Gabor域的稀疏化表达.基于线性Radon变换(LRT
4、)的局部平面波分解的方法是主要的地震信号的线性表达方法, 以其明确的物理含义(可以描述局部平面波的传播方向)而被广泛应用.非线性的Radon 变换,如抛物 Radon 变换(Sacchi 和 Ulrych, 1995),双曲 Radon 变换(Liu 和 Sacchi, 2002),多项式Radon变换(牛滨华等,2001)等,则主要用于多次波压制、去噪、数据 重构等方面,不同于Fourier变换和正交小波变换,LRT的基函数是非完备正交的,因而导 致LRT结果难以实现稀疏表达.王雄文和王华忠(2014)将常规的LRT转化为压缩感知问 题,在压缩感知的理论框架下实现了高分辨率的LRT.在Rad
5、on域中,可以测量局部平面波 在炮点处的入射射线参数P木文中,首先考虑地震数据在局部平面波域中的稀疏表达,并引入地震数据的特征分解方 法,随后给出基于特征波数据的偏移成像与反演成像方法.2局部平面波的成像在天然地震学中,通常用检波器排列的某种组合来分离相干信号和噪声.最基本的方法就 是波束合成方法(Beam Forming, Rost 和 Thomas, 2002;Kr(iger 等,1993),它利用局部 平面波波前到达不同检波器的时差,对检波器接收到的地震记录进行时移叠加,使局部平 面波相干叠加而其他信号相互抵消,从而提取某种特定的波型(震相).波束合成的思想也可应用于勘探地震学中,如高斯
6、束偏移方法(Hill, 1990;Hale, 1992),在炮集中对不同检波器中心进行加窗局部倾斜叠加(即波束合成),利用炮集的 丁-P谱进行偏移.在射线束偏移(Sun等,2001)和控制束偏移(Gao, 2006)中,通过选 择t-P谱的能量,实现有选择性的波场传播与成像.波束合成算法需要估计局部平面波的方向信息.本文中,我们考虑Radon谱约束下的局部 平面波反演方法(胡江涛等,2013):式中,X(3)为局部空间窗内的地震信号,X (a) =x上述Radon谱约束下的局部平面波反演方法可以提高Radon谱的精度并压制泄漏的噪声. 更进一步地,为了得到更稀疏的局部平面波解(王雄文和王华忠,
7、2014),可以把式(1) 改写为压缩感知问题,在L0范数意义下求解如下反问题:利用匹配追踪方法求解(2)式,可以得到局部空间窗内的稀疏平面波解.通过传统的线性Radon变换或者高分辨率的平面波反演方法,可以得到局部平面波及其慢 度矢量(或射线参数P).此时,我们可以考虑对炮集的波束合成.在共炮点集中,对中心检波点X其中,d (x我们把在中心检波点处进行一次波束合成之后的地震数据记为d若将地震数据分选为共检波点道集,此时可以再次在中心炮点X其中,(X同公式(3),经过中心炮点处的再次波束合成后,实现了入射方向为p将(3)式中的波束合成代入(4)式中,可得同时在中心炮点和检波点处的波束合成算法:
8、在三维时,上述公式中的坐标与射线参数均为矢量,表示为x从公式(5)可以看出,特征波场合成与局部平面波的入射与出射射线参数有关,其本质 是利用局部平面波的入射和出射慢度计算道间时差(At=pAx),并将波束合成孔径内 的地震道按照道间时差进行时移叠加.在时移叠加的过程中,只有给定慢度矢量(局部平 面波入射与出射射线参数)的局部平面波波前实现了同相叠加,而其他的信号相互抵消, 从而实现不同震相(波型)的分离,即特征波场合成.事实上,上述方法也可以推广到任意线性排列的观测系统中(如VSP或井间观测系统等), 只要满足局部的炮点/检波点排列具有线性特征.相应的,局部平面波的慢度矢量(射线参 数)定义为
9、地震波走时对地震排列(炮点或检波点坐标)的方向导数.此外,在频率域特征波场合成公式(5)对应的相移矩阵为A3特征波成像优势地震数据偏移成像的过程可以描述为地表波场的反传播加上成像条件.以运动学射线追踪 为传播算子的Kirchhoff叠前深度偏移(PSDM)技术在横向变速剧烈情况下存在多路径和 振幅焦散等问题,基于射线束的成像方法(刘少勇等,2012)在一定程度上解决了这些问 题.高斯束偏移技术(Hill, 1990;Hale, 1992)已经成为工业界一个重要的叠前深度偏移 成像工具.Sun和Schuster (2001)提出菲涅尔带控制的波路径偏移方法,将积分法偏移 在成像空间中沿等时面的投
10、影转化为沿胖射线路径的投影,在计算效率上有一定优势.Liu 和Palacharla (2011)对高斯束偏移进行简化,只保留其运动学特征,推导出简单的解析 表达式来控制射线束的宽度和振幅,理论上说该实现方式效率要优于高斯束偏移,成像精 度基本有保障.快速射线束(Fast-beam)偏移(Gao等,2006)便是该思想的一种实现,特 征波偏移成像利用射线束波传播算子来反传局部特征波场,既保留了射线类成像的灵活性, 又保证了成像精度.特征波场是地震波场在高维空间中的稀疏表达,且局部平面波的入射与出射慢度矢量已知. 利用这个性质,射线束类的偏移方法可以实现从“画弧”向“搬家”的转变,进一步提高 成像
11、效率,与上节的特征波场合成方法结合,匹配合适的局部射线束传播算子,可以实现 特征波叠前深度偏移成像(characteristic wave imaging,简记为CWI).更进一步的,“搬家”或一对一映射的偏移方法与基于角度道集的层析速度反演可以自然的融合成一体, 共享同样的数据空间和成像空间的数据体.成像与层析速度反演一体化也是特征波成像的 一大优势.传统基于炮道集的射线束(包括高斯束)成像条件可以描述为其中:x表示成像点坐标,W其中,W传统成像条件和特征波成像条件的几何解释可由图1表示.图中红色反射层为地层真实位 置,传统射线束成像方法是将局部平面波投影到整个椭圆(三维情况为椭球)上.由于
12、特 征波场在地表的入射与出射慢度已知,因而特征波成像只对特征波场进行“能量搬家”. 红色小椭圆的位置代表一对特征波射线束在成像域的菲涅尔带范围.从图1可以看出,基 于特征波场表达,实现了成像方式从“画弧”到“搬家”的转变.传统积分类成像的画弧成像过程必然会带来画弧噪声,而搬家的方式可最大限度减少此类噪声同时能兼顾菲涅尔 带的成像分辨率,其成像效率也会得到很大的提高.4特征波场走时反演传统的波动方程走时反演方法(Lu。和Schuster, 1991)利用互相关方法估计观测数据与 模拟数据的时差,需要对两者进行加窗处理,从而提取初至波或者折射波.该方法需要拾 取特定的波型并在时间域加窗,繁重的拾取
13、工作限制了该方法的应用.利用特征波场合成,单道地震记录被分离为一系列单独的震相.此时,可直接利用独立的 震相进行反演,无需对地震记录加窗处理.更重要的是,特征数据使得多尺度、逐层反演 成为可能,如先利用折射波反演浅层速度,再利用由浅到深的反射波反演中深层速度.由于特征波场只含有单独的波型(震相),可用互相关函数(Luo和Schuster, 1991)测 量观测数据与模拟数据的时差,无需对地震数据加窗.本文中,我们将特征波场合成算法 应用初至波场分解,并应用于透射波走时反演中.至于特征波场分解在反射波走时反演中 的应用,可以参考冯波(2015),本文不再赘述.若观测数据d与模拟数据U经过特征波场
14、合成后,得到的初至震相分别可以表示为d观测数据与模拟数据的初至走时残差At定义为基于走时残差L2范数的目标函数可以表示为经过推导,目标泛函梯度(冯波,2015)为其中,g对于初至走时反演,d梯度类的迭代算法为:其中,a5单独震相的成像反演根据上文的分析,我们可以利用特征波场合成实现不同波型(震相)的分离,而分离后的 单独震相可以用于后续的成像及反演.在数值实验中,我们首先分析不同信噪比的数据对射线参数估计以及特征波场合成的影响, 然后,用一个简单模型(凹陷模型)展示特征波场成像的效果及其优势,最后,用特征波 场合成方法提取初至波场,并用于波动方程初至走时反演中.5. 1特征波场成算法我们用Ma
15、rmousi模型正演观测数据.速度模型为:X方向网格数目nx二737, Z方向网格数目 nz=250;网格间隔均为12. 5 m.观测系统设置为:第一炮的坐标为625m,炮间隔25叫总共 301炮.检波器为固定排列,第一道的坐标为625m,检波器间隔12. 5叫总共601道.用10 Hz主频的Ricker子波作为震源.首先,选取中心地震道的炮点和检波点坐标为(1850, 3000) .图2a中,从左到右分别是为了测量局部平面波的入射射线参数,需要抽取共检波点道集.图2b中,从左到右分别是利用图2中拾取的入射和出射射线参数,根据公式(5),可进行特征波场合成.图3a中, 第14道为合成的特征波场
16、.其中,第1道代表对直达波的特征波场合成,它将直达波从 地震记录中提取.第24道代表对反射波的特征波场合成,它将反射波从地震记录中提取. 第5道为重构的中心道(第14道叠加),与原始中心地震道(第6道)相比,只保留 了我们提取的震相,提高了特征波的信噪比,同时,考虑到特征波场(第14道)的时间 局部特征,可以(在变换域)进一步压缩存储.接着,为了测试噪声对特征波场合成的影响,我们对原始数据加入了高斯白噪声 (S/N=100) .图4中展示了相应的道集及其T-p谱(所有参数含义同图2).可以看出,高 斯随机噪声(S/N=100)降低了 T-p谱的分辨率,因而可能会对T-p谱的自动拾取有一 定的影
17、响,但对人工拾取的影响较小,图3b中,展示了对图4中的数据进行特征波场合成 的结果(所有参数同图2a).显然,特征波场合成之后,仍然能够得到高信噪比的单独震 相.这是由于特征波场合成的过程中,同相的信号相互增强,而其他信号相互干涉.若进一 步增强噪声强度(S/N=20),此时,随机噪声完全压制了有效信号,即使在t-p谱也无 法拾取射线参数.若我们仍然用真实的射线参数(图2中拾取的结果),特征波场合成结 果如图4c所示,其中,直达波的波形仍然具有较为保真的形态,而后续的几个反射波则受 到了一定程度的噪声干扰,略有震荡.上述试验说明,本文中的特征波场合成算法十分稳健,随机噪声对其影响较小.当噪声能
18、 量逐渐增强时,工-p谱的分辨率会逐渐降低,直至无法拾取准确的射线参数,进而导致特 征波场合成算法失效.5.2特征波场成像优势本文用洼陷模型测试特征波场成像方法的有效性,其速度模型如图5a所示.模拟数据共 201炮,炮间距20叫 每一炮301道接收,道间距10m.本节分别对该模型数据进行传统的 Kirchhoff PSDM和特征波成像来进行效率效果的对比.图5b是Kirchhoff PSDM的成像剖 面,由于其成像过程是通过等时面叠加收敛绕射波,成像剖面上留下了画弧噪声.图5c是 特征波场成像得到的成像剖面,由于其通过对特征域数据进行“搬家”,其画弧范围控制 在局部范围内,其画弧噪声明显弱于传
19、统Kirchhoff PSDM.由于特征波的成像过程中可以 得到射线方向,张角道集可以方便的输出.传统Kirchhoff PSDM一般只能输出地表偏移距 道集,在对后续的速度分析或层析成像的适应性方面上来说,张角道集质量要优于地表偏 移距道集.图5 (d, e)分别为洼陷模型3500 m处的Kirchhoff PSDM地表偏移距道集和特 征波成像的张角道集,对比两图可见,张角道集的成像质量要高于地表偏移距道集.特征波成像由于只对特征数据进行反投影(特征数据是对原地震数据的一个压缩处理), 其成像效率要远远高于传统的Kirchhoff PSDM或者传统的射线束成像.本文对比了 Kirchhoff
20、 PSDM和特征波成像在洼陷模型上的成像效率,具体数据见表1,数据得到了有 效的压缩,同时成像效率也有近一个数量级的提高.由于特征波场包含了数据方向信息, 特征波场成像方法可以方便的推广到目标体成像,图5f展示了仅仅对最后一个反射层的 面向目标成像结果.对特定目标成像将会进一步提高效率,缩短处理周期.5. 3初至走时反演速度模型在本节中,我们将特征波场合成算法分离初至波型,并用于波动方程初至走时层析中.速 度模型(BP2004速度模型的右端部分,横向约26k叫 深度约12km)如图6a所示.观测系 统为:炮点深度12. 5%炮点间隔为200% 共107炮,检波器深度12. 5% 间隔400叫
21、每炮 54道.观测系统总共5778道,初始速度模型采用一维模型,其中,水体速度以及水底深度 已知.总共进行18次反演迭代,反演得到的速度模型见图6b.反演得到了整体背景速度的 变化趋势,其中模型右侧的高速隆起反演得较好,但没有反演出浅层小尺度的低速异常. 图7中为初至波走时对比曲线.其中,红色、绿色和蓝色曲线分别是真实模型、初始模型 和第18次反演的速度模型中的初至走时,可以看出,真实的初至走时与反演得到的速度模 型中的初至走时基本吻合,实现了初至走时反演的目标.6特征波场辅助偏移成像本文提出了一种地震数据的稀疏表达方法.不同于传统的变换域中的数据稀疏表达理论, 该方法利用局部平面波的传播方向
22、(慢度矢量),在中心炮检点处同时进行波束合成,从 而将地震数据投影到局部平面波域(高维空间,二维观测对应5D数据;三维观测对应9D 数据)中.由于波束合成后的地震数据描述了局部平面波的方向特征,因此我们称之为特 征波场,相应的,我们称这种波场合成(分离)算法为特征波场合成方法.特征波场合成的思想是利用局部平面波到达不同检波器的时差,对地震记录进行时移叠加, 使得给定的入射/出射慢度的局部平面波相干加强而其他信号相互抵消,从而分离地震记 录中的不同的波型(震相).因此,该方法对含有噪声的地震数据有较好的适应能力.文中 的特征波场合成正问题需要估计射线参数,当地震数据受噪声干扰严重时,难以在常规
23、t-P谱中自动估计局部平面波的射线参数(慢度矢量).此时,可以考虑基于反演理论的 特征波场合成方法.在压缩感知理论框架下同时反演局部平面波及其传播方向,从而提高 特征波合成的自动化程度并保持方法的稳健性.然而,特征波场合成正问题要求局部炮/检排列在同一水平面上,因此,对于地表高程剧烈 变化的起伏地表,该方法不再成立此外,若炮点(检波器)采样过于稀疏,也会影响射 线参数估计的精度.此时,可以仅考虑单次波束合成.通过特征波场合成,可以将地震数据分解为单独的震相.此时,传统的Kirchhoff及高斯 束偏移方法可以实现从“沿等时面的画弧”到向“向反射点归位搬家投影”的转变.相应 的,常速介质中的偏移脉冲响应由“椭圆/球”退化为点响应.同时, 由于在偏移过程中不 需要对炮/检射线参数扫描,因而叠前偏移效率将会大幅提高,若炮点的慢度矢量(入射射 线参数)难以估计,特征波场的偏移成像可以自动退化为高斯束类的成像方法.更进一步地,特征波偏移成像这种一对一映射的偏移方法非常适合与基于角度道集的速度 层析反演方法融合成一体.它们具有同样的数据空间和成像空间.这种特征波场驱动的高效 偏移成像和速度层析反演融为一体的技术体系会改变目前的地震成像流程,与大规模计算 机集群和现代图像处理方法结合,在很大程度上可以实现半自动或自动的地震波成像.