CN104459770A - 一种高维地震数据规则化方法 - Google Patents
一种高维地震数据规则化方法 Download PDFInfo
- Publication number
- CN104459770A CN104459770A CN201310439558.5A CN201310439558A CN104459770A CN 104459770 A CN104459770 A CN 104459770A CN 201310439558 A CN201310439558 A CN 201310439558A CN 104459770 A CN104459770 A CN 104459770A
- Authority
- CN
- China
- Prior art keywords
- data
- frequency
- dimensional
- prior
- alias
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明提供了一种高维地震数据规则化方法,属于油气勘探开发中地震数据预处理领域。本方法包括:(1)输入原始三维炮道集;(2)利用高维傅里叶变换将时间空间域的所述原始三维炮道集变换到频率空间域中,然后进行分选得到五维炮道集;(3)对一个频率片进行最优化求解,得到无假频频谱;(4)判断是否达到迭代门槛值,如果是,则转入步骤(5),如果否,则返回步骤(3);(5)利用高维傅里叶反变换将所述无假频频谱变换到时间空间域,即得到规则炮道集数据;(6)输出所述规则炮道集数据。
Description
技术领域
本发明属于油气勘探开发中地震数据预处理领域,具体涉及一种高维地震数据规则化方法。
背景技术
在地震数据采集中,空间采样的不规则以及空间采样间隔过大(采集成本限制)导致地震数据的空间假频经常发生。由于偏移算子不能消除由空间采样引起的假频,以致这些假频能量以假象的形式出现在剖面中,这些假象会误导地震剖面的解释。因此为了实现真振幅成像,必须在偏移前通过地震数据插值等方法除去这些假频能量。在过去的十多年里,地球物理学家发展了很多解决地震数据非规则性和假频问题的方法,大体上可分为四类。描述如下:
第一类,基于褶积理论。这种方法的主要思想是,用地震数据的低频分量构造一个滤波器,并利用这个滤波器实现对地震数据的插值。属于这种类型的方法包括FX域预测滤波插值(Spitz,1991),TX域预测误差滤波插值(Claerbout,1992),FK域预测滤波插值(Gulunay,2003)。它们通常要求数据域的同相轴是线性的,只有满足线性假设条件,所构造的滤波器才有意义。这些方法的另一个要求是地震数据的空间采样方式必须是规则的,但是在数据采集过程中这个要求是很难实现的。
第二类,基于偏移算子和反偏移算子。这种方法先应用偏移算子将地震数据投影到成像域中,通过在成像域里完成地震数据的插值工作,然后将成像域的数据反偏移回数据域中。由Stolt(2002)和Trad(2003)介绍的方法都属于这种类型。这种方法的偏移算子种类繁多,包括正常时差校正(NMO),倾角时差校正(DMO),零偏移距偏移(MZO),方位角时差校正(AMO),叠前时间偏移(PSTM),叠前深度偏移(PSDM)等。这类方法虽然可以很好的保持地震数据中同相轴的相位信息,但是振幅信息往往是失真的。
第三类,基于FK频谱估计。这种方法可以另外分为两小类:一种基于反演理论,另一种基于相关理论。基于反演理论的主流方法可参考A.J.W.Duijndam(1999)和Sacchi(2004,2007)的文章。这类方法的主要思想是,在反演理论的框架下,应用合适的正演算子和约束条件,估计地震数据的FK频谱。这类反演问题的模型空间和数据空间通常是FX域和FK域的地震数据,而正算子是傅立叶变换。这些方法的约束条件差异很大,不同的约束条件会产生不同的结果,因此约束条件是该方法的关键因素。总所周知,一个反问题经常会对应着多个解,不同的限制条件意味着在不同的解空间寻求结果,因此对于不同的约束条件而言解是完全不同的。另一类基于相关理论的方法包括抗泄露傅立叶变换(ALFT)(Sheng Xu,2005),凸集投影(POCS)(Ray Abma,2006),匹配追踪(MP)(Al i Ozbek,)以及其他一些与基于地震数据相关的方法。在相关过程中,这些方法将相关能量最大的基作为地震数据中能量最大的分量,然后选择这个最大能量作为FK谱估计的一个分量。实际上如果将相关理论看作解反演问题的一种方法,那么每次选择的基就是梯度,而相关能量就是这个反演问题的步长。如果数据中的假频不严重,那么这类方法能得到很好的结果。但是如果假频严重,最大相关能量的基可能就不再是地震数据有效信号的分量而是假频,则此方法得到的插值结果就很差。
最后一种,基于矩阵分解。如果将地震数据看作一个高维矩阵,那么在奇异值分解(SVD)中,数据插值可以通过插值奇异向量得到。由Sacchi(2011)提出的方法就属于这一类。如果数据域中的同相轴是线性的,那这种方法将得到更好的结果,因为如果同相轴是线性的,那么大奇异值的所对应的奇异向量的数量将大大减少。但是其缺陷也是显而易见的,因为地震数据的随机性很强,因此在奇异值分解时就会产生奇异值从而导致矩阵误解。
上述四类方法中,最后两类是当今地震数据插值中最为流行的,因为它们都具备高维地震数据插值的能力,并且插值效果也比较理想。
发明内容
本发明的目的在于解决上述现有技术中存在的难题,提供一种高维地震数据规则化方法,解决地震数据采集中,空间采样的不规则以及空间采样间隔过大(采集成本限制)导致地震数据的空间假频的问题,这些假频能量以假象的形式出现在剖面中,这些假象会误导地震剖面的解释。本方法是在反演过程中以高维傅里叶变换为正演算子,数据线性拉东变换谱为约束的高维数据规则化方法。
本发明是通过以下技术方案实现的:
一种高维地震数据规则化方法,所述方法包括:
(1)输入原始三维炮道集;
(2)利用高维傅里叶变换将时间空间域的所述原始三维炮道集变换到频率空间域中,然后进行分选得到五维炮道集;
(3)对一个频率片(该频率片是由时间域通过傅里叶变换到频率域,有多少个时间采样点就有多少个频率片。)进行最优化求解,得到无假频频谱;
(4)判断是否达到迭代门槛值,如果是,则转入步骤(5),如果否,则返回步骤(3);
(5)利用高维傅里叶反变换将所述无假频频谱变换到时间空间域,即得到规则炮道集数据;
(6)输出所述规则炮道集数据。
所述步骤(3)是采用共轭梯度算法求解下面的公式实现的:
这个公式就是每个频率片要求解的,即公式是频率域的基础上求解。
其中,是无假频的规则数据频谱,即无假频频谱,mprior是先验模型,CM为数据的Radon谱,CD为一般约束条件(最小二乘反演的约束条件,一般给0.01),dobs为观测数据(观测数据就是野外得到的炮集数据经过傅里叶变到频率域的数据。),A表示采样算子矩阵(采样算子就是根据野外观测***定义的一个矩阵,采集到的地方为1,没有采集到的地方为0.),AH为A的共轭转置矩阵,mprior=0。
所述数据的Radon谱是这样得到的:
利用步骤(1)中的炮道集数据中同相轴的局部线性,进行线性拉东变换,得到拉东谱,具体实现公式如下:
R(f,p)=∫d(f,x)e-2πfpxdx
其中,f,p,x分别表示频率、射线参数和空间采样位置,d(f,x)是原始数据的频率空间表示,R(f,p)为所述数据的Radon谱。
所述步骤(4)中的迭代门槛值设置成10-6(即10的负六次方)。
与现有技术相比,本发明的有益效果是:
1)本发明方法使用了高精度,高效率的Intel编译库中高维傅里叶正反变换,其优势在于进行同等效果的傅里叶时的速度更加的快捷稳健,提高计算效率;
2)本发明方法利用原始数据同相轴线性拉东谱作为反演约束条件,其优势在于利用同相轴的拉东谱约束算法更加稳健,反演的结果更加的接近真解。
附图说明
图1本发明实施例中的规则采样的五个同相轴理论合成数据;
图2本发明实施例中图1的傅里叶变换规则数据频谱;
图3本发明实施例中图1随机抽取三分之一道后非规则采样数据;
图4本发明实施例中图3的傅里叶变换频谱;
图5本发明实施例中的通过本方法得到的频谱反变换得到的规则数据;
图6本发明实施例中的通过本方法得到的规则数据频谱;
图7本发明实施例中的合成的规则炮道集;
图8本发明实施例中的图7随机抽取三分之一道得到的非规则炮道集;
图9本发明实施例中的通过本方法得到的规则炮道集;
图10本发明实施例中的规则后的炮道集与原始规则炮道集的残差;
图11本发明实施例中随机抽取三维盐丘模型三分之一炮道集以及道集内随机抽取三分之一道得到的非规则炮道集;
图12本发明实施例中通过本方法得到的规则炮道集;
图13本发明实施例中图11局部放大非规则炮道集;
图14本发明实施例中图12局部放大规则炮道集;
图15本发明实施例中某野外实际非规则炮道集;
图16本发明实施例中通过本方法差值充零得到的炮道集;
图17本发明实施例中通过本方法得到的规则炮道集;
图18本发明方法的步骤框图。
具体实施方式
下面结合附图对本发明作进一步详细描述:
本发明通过对具有一定信噪比(即数据中没有强能量的野值和线性面波能量)的炮道集数据进行高维傅里叶变换,然后利用炮道集数据同相轴的局部线性拉东变换频谱作为约束,通过反演将FK谱(频率波数谱)的假频分量减弱得到规则数据的频谱,最后再利用高维傅里叶反变换得到规则数据。
如图18所示,所述方法具体包括以下步骤:
(1)准备具有高分辨率、高信噪比的原始地震炮集数据
(2)利用Intel编译库中高维傅里叶变换将时间空间域的炮集数据变换到频率空间域中,这个过程中将地震道头信息的炮点坐标记录下来形成道头信息索引;野外采集的地震数据是时间和空间的3维数据,而在数据规则化中要把3维数据转换到5维,即频率、x、y、hx和hy,hx和hy分别是x和y方向的偏移距,这样的数据称为高维数据。5维数据是通过高维地震数据插值产生的:输入的是原始3D炮道集,进入程序后进行傅里叶变换到频率空间域,之后再进行分选到5D进行插值,即x、y,hx,hy,t,这样得到的效果更加的接近实际采集的情况;
(3)利用(1)中原始炮道集数据中同相轴的局部线性,进行线性拉东变换,得到拉东谱,具体实现公式如下所示:
R(f,p)=∫d(f,x)e-2πfpxdx
其中,f,p,x分别表示频率、射线参数和空间采样位置,d(f,x)是原始数据的频率空间表示,R(f,p)为数据的Radon谱;
(4)反演过程中将拉东谱作为约束条件,对每个频率片进行反演得到规则数据的无假频频谱,具体如下:
其中是这个问题的最佳模型即无假频的规则数据频谱,mprior是先验模型,CM为数据的Radon谱,CD为一般约束条件,dobs为观测数据,A表示采样算子矩阵,AH为A的共轭转置矩阵。一般而言mprior=0,解最优化问题用共轭梯度算法求解。
(5)利用Intel编译库中高维傅里叶反变换将(4)中反演得到的无假频频谱变换到时间空间域得到规则采样的炮道集数据,即规则炮道集数据。
图1到图6是2维理论数据验证方法理论;图7到图10是3维理论数据进行算法测试;图11到图17是利用3维盐丘模型进行算法测试的结果。
本发明属于第三种类型中基于反演理论的方法。该方法的正过程即规则且无假频地震数据的傅立叶谱m(f,k)是本发明想要的模型。存在空间假频的且非规则的野外数据d(f,x)是已知的。则正演过程可以表示为d(f,x)=A(x,k)m(f,k)向量d(f,x),m(f,k)及矩阵A(x,k)的元素可以相应地作如下表达(下面3个公式表示每个矩阵的具体元素的含义,即原始观测数据,反演模型空间矩阵元素和采样矩阵算子的元素):
di(f,x)=d(f,xi)
mi(f,k)=m(f,ki)
amn=exp(i2πxm·kn)
其中d(f,xi)是野外数据,xi是其空间坐标。m(f,ki)是要估计模型的分量,其坐标是ki。xm·kn是xm和kn的内积。f是地震数据的频率。这个正过程与Dui jndam等人(1999)给出的定义相同,反过程利用了地震数据的拉东谱能量分布构造一个椭圆范数做为约束条件,所构造的椭圆范数能提高FK谱的相对能量分布。如果同相轴是线性的,那么拉东变换可以在FK域中有效地聚集这些同相轴的能量,相应地,FK谱中的假频分量将相对地被减弱。因此,在这样的约束条件下插值过程会变得更加稳定,反演问题的解也会变得更加稀疏且更合逻辑。
如果数据采样均匀,每个采样点在数据空间的元素中的地位是等价的,并且是相互独立的。这样数据的协方差矩阵主对角线上的值相等。其余非对角元素为0,即是一个单位阵I。这样起能量分布就是一个球面,实际数据采样是不等间隔的,这样数据的协方差矩阵主对角线值不再相等,也无法用单位阵表示。协方差矩阵不仅要描述各个采样点处能量的强弱关系,还要描述不同采样点之间的相互关系,协方差矩阵的值主要由观测***的形式决定。首先,协方差矩阵必须能反映各个采样点处如下的能量关系:在采样比较密集的方向上,对应着比较强的能量,因此在协方差矩阵对角线上与之相应的元素的值应该比较大;在采样比较稀疏的方向上,对应着比较弱的能量,因此在协方差矩阵对角线上与之相应的元素的值应该比较小。协方差矩阵主对角线上的元素是采样点的自相关,表示各个采样点处能量的大小;而非主对角线上的元素是不同采样点之间互相关,表示不同采样点之间相关性的大小。在正常情况下,均假设不同点处采集到的地震数据是不相关的,即各个采样点之间是互相独立的。因此,协方差矩阵非主对角线上的元素的值均默认为0。如果每个检波器所接收到的地震数据是相互独立的,那么协方差矩阵非对角线上的元素的值均为0,协方差矩阵变成了一个对角矩阵。即数据空间的单位圆对应的不再是一个球面,而是一个椭球面。
上述技术方案只是本发明的一种实施方式,对于本领域内的技术人员而言,在本发明公开了应用方法和原理的基础上,很容易做出各种类型的改进或变形,而不仅限于本发明上述具体实施方式所描述的方法,因此前面描述的方式只是优选的,而并不具有限制性的意义。
Claims (4)
1.一种高维地震数据规则化方法,其特征在于:所述方法包括:
(1)输入原始三维炮道集;
(2)利用高维傅里叶变换将时间空间域的所述原始三维炮道集变换到频率空间域中,然后进行分选得到五维炮道集;
(3)对一个频率片进行最优化求解,得到无假频频谱;
(4)判断是否达到迭代门槛值,如果是,则转入步骤(5),如果否,则返回步骤(3);
(5)利用高维傅里叶反变换将所述无假频频谱变换到时间空间域,即得到规则炮道集数据;
(6)输出所述规则炮道集数据。
2.根据权利要求1所述的高维地震数据规则化方法,其特征在于:所述步骤(3)是采用共轭梯度算法求解下面的公式实现的:
其中,是无假频的规则数据频谱,即无假频频谱,mprior是先验模型,CM为数据的Radon谱,CD为一般约束条件,dobs为观测数据,A表示采样算子矩阵,AH为A的共轭转置矩阵,mprior=0。
3.根据权利要求2所述的高维地震数据规则化方法,其特征在于:所述数据的Radon谱是这样得到的:
利用步骤(1)中的炮道集数据中同相轴的局部线性,进行线性拉东变换,得到拉东谱,具体实现公式如下:
R(f,p)=∫d(f,x)e-2πfpxdx
其中,f,p,x分别表示频率、射线参数和空间采样位置,d(f,x)是原始数据的频率空间表示,R(f,p)为所述数据的Radon谱。
4.根据权利要求1所述的高维地震数据规则化方法,其特征在于:所述步骤(4)中的迭代门槛值设置成10-6。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310439558.5A CN104459770B (zh) | 2013-09-24 | 2013-09-24 | 一种高维地震数据规则化方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310439558.5A CN104459770B (zh) | 2013-09-24 | 2013-09-24 | 一种高维地震数据规则化方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN104459770A true CN104459770A (zh) | 2015-03-25 |
CN104459770B CN104459770B (zh) | 2017-06-16 |
Family
ID=52906115
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201310439558.5A Active CN104459770B (zh) | 2013-09-24 | 2013-09-24 | 一种高维地震数据规则化方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN104459770B (zh) |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105549078A (zh) * | 2015-12-31 | 2016-05-04 | 中国石油天然气股份有限公司 | 不规则地震数据的五维插值处理方法及装置 |
CN106646612A (zh) * | 2016-12-20 | 2017-05-10 | 中国地质大学(北京) | 基于矩阵降秩的地震数据重建方法 |
CN107703539A (zh) * | 2017-09-18 | 2018-02-16 | 中国石油天然气股份有限公司 | 抗假频的地震数据插值方法和装置 |
CN108169795A (zh) * | 2018-02-11 | 2018-06-15 | 中国石油化工股份有限公司 | 基于随机采样的数据规则化方法 |
CN108345034A (zh) * | 2018-02-06 | 2018-07-31 | 北京中科海讯数字科技股份有限公司 | 一种地震数据规则化方法 |
CN108490486A (zh) * | 2018-02-01 | 2018-09-04 | 北京奥能恒业能源技术有限公司 | 一种地震资料反演的方法及装置、设备 |
CN110244353A (zh) * | 2019-06-25 | 2019-09-17 | 北京中科海讯数字科技股份有限公司 | 一种基于稀疏范数优化算法的地震数据规则化方法 |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN1332225C (zh) * | 2004-12-29 | 2007-08-15 | 中国石油天然气集团公司 | 一种采用深度域滤波消除线性与非线性干扰波的方法 |
US7751277B2 (en) * | 2008-03-17 | 2010-07-06 | Pgs Geophysical As | Method for interpolating seismic data by anti-alias, anti-leakage Fourier transform |
CN102692650B (zh) * | 2011-03-23 | 2014-10-29 | 中国石油天然气集团公司 | 一种具有假频压制功能的井筒波分离方法 |
-
2013
- 2013-09-24 CN CN201310439558.5A patent/CN104459770B/zh active Active
Non-Patent Citations (3)
Title |
---|
孟小红,等: "基于非均匀快速傅里叶变换的最小二乘反演地震数据重建", 《地球物理学报》 * |
张华,等: "基于jitter采样和曲波变换的三维地震数据重建", 《地球物理学报》 * |
王立明,等: "一种有效的地震道插值方法", 《地球物理学进展》 * |
Cited By (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105549078A (zh) * | 2015-12-31 | 2016-05-04 | 中国石油天然气股份有限公司 | 不规则地震数据的五维插值处理方法及装置 |
CN106646612A (zh) * | 2016-12-20 | 2017-05-10 | 中国地质大学(北京) | 基于矩阵降秩的地震数据重建方法 |
CN107703539A (zh) * | 2017-09-18 | 2018-02-16 | 中国石油天然气股份有限公司 | 抗假频的地震数据插值方法和装置 |
CN107703539B (zh) * | 2017-09-18 | 2019-05-07 | 中国石油天然气股份有限公司 | 抗假频的地震数据插值方法和装置 |
CN108490486A (zh) * | 2018-02-01 | 2018-09-04 | 北京奥能恒业能源技术有限公司 | 一种地震资料反演的方法及装置、设备 |
CN108345034A (zh) * | 2018-02-06 | 2018-07-31 | 北京中科海讯数字科技股份有限公司 | 一种地震数据规则化方法 |
CN108345034B (zh) * | 2018-02-06 | 2021-08-03 | 北京中科海讯数字科技股份有限公司 | 一种地震数据规则化方法 |
CN108169795A (zh) * | 2018-02-11 | 2018-06-15 | 中国石油化工股份有限公司 | 基于随机采样的数据规则化方法 |
CN110244353A (zh) * | 2019-06-25 | 2019-09-17 | 北京中科海讯数字科技股份有限公司 | 一种基于稀疏范数优化算法的地震数据规则化方法 |
CN110244353B (zh) * | 2019-06-25 | 2024-01-30 | 北京中科海讯数字科技股份有限公司 | 一种基于稀疏范数优化算法的地震数据规则化方法 |
Also Published As
Publication number | Publication date |
---|---|
CN104459770B (zh) | 2017-06-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN104459770A (zh) | 一种高维地震数据规则化方法 | |
Zhou et al. | Reverse time migration: A prospect of seismic imaging methodology | |
Fang et al. | Seismic data interpolation based on U-net with texture loss | |
Kumar et al. | Source separation for simultaneous towed-streamer marine acquisition—A compressed sensing approach | |
CN108037526B (zh) | 基于全波波场vsp/rvsp地震资料的逆时偏移方法 | |
Trad | Five-dimensional interpolation: Recovering from acquisition constraints | |
CA2793504C (en) | Windowed statistical analysis for anomaly detection in geophysical datasets | |
Cheng et al. | Fast dual-domain reduced-rank algorithm for 3D deblending via randomized QR decomposition | |
CN102636811B (zh) | 一种海上二维地震资料中多次波的消除方法 | |
CN106932824B (zh) | 陆地地震勘探资料的降维自适应层间多次波压制方法 | |
CN106133556A (zh) | 确定波场的分量 | |
CN105277985A (zh) | 一种基于图像处理的ovt域地震数据规则化方法 | |
CN103149592A (zh) | 一种变偏移距vsp波场分离方法 | |
CN109738952A (zh) | 基于全波形反演驱动的被动源直接偏移成像方法 | |
Kim et al. | Hypocentral parameter inversion for regions with poorly known velocity structures | |
CN105629299A (zh) | 角度域叠前深度偏移的走时、角度表获取方法及成像方法 | |
Yin et al. | Improving horizontal resolution of high-frequency surface-wave methods using travel-time tomography | |
CN102736108B (zh) | 基于样条拟合的真三维地震数据噪声压制方法 | |
Gan et al. | EWR‐net: Earthquake waveform regularization network for irregular station data based on deep generative model and RESNet | |
CN110850469A (zh) | 一种基于克希霍夫积分解的地震槽波深度偏移的成像方法 | |
Igel et al. | SANS: Publicly available daily multi‐scale seismic ambient noise source maps | |
Qu et al. | 3-D Least-Squares Reverse Time Migration in Curvilinear-τ Domain | |
Oropeza | The singular spectrum analysis method and its application to seismic data denoising and reconstruction | |
CN106054242A (zh) | 三维各向异性衰减介质波场模拟方法 | |
Montel et al. | Kinematics of common-image gathers—Part 1: Theory |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant |