CN110516358A - 一种基于稀疏表示的地震各向异性参数反演方法 - Google Patents
一种基于稀疏表示的地震各向异性参数反演方法 Download PDFInfo
- Publication number
- CN110516358A CN110516358A CN201910800069.5A CN201910800069A CN110516358A CN 110516358 A CN110516358 A CN 110516358A CN 201910800069 A CN201910800069 A CN 201910800069A CN 110516358 A CN110516358 A CN 110516358A
- Authority
- CN
- China
- Prior art keywords
- dictionary
- parameter
- sparse
- inversion
- atom
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 107
- 238000012549 training Methods 0.000 claims abstract description 45
- 230000008569 process Effects 0.000 claims abstract description 23
- 230000014509 gene expression Effects 0.000 claims abstract description 18
- 239000011159 matrix material Substances 0.000 claims description 25
- 238000012545 processing Methods 0.000 claims description 12
- 230000015572 biosynthetic process Effects 0.000 claims description 11
- 238000004364 calculation method Methods 0.000 claims description 10
- 235000013399 edible fruits Nutrition 0.000 claims description 4
- FESBVLZDDCQLFY-UHFFFAOYSA-N sete Chemical compound [Te]=[Se] FESBVLZDDCQLFY-UHFFFAOYSA-N 0.000 claims description 2
- 230000000750 progressive effect Effects 0.000 abstract description 2
- 230000006870 function Effects 0.000 description 12
- 238000000354 decomposition reaction Methods 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 5
- 238000005070 sampling Methods 0.000 description 5
- 238000003786 synthesis reaction Methods 0.000 description 5
- 241001269238 Data Species 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 239000002131 composite material Substances 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000005457 optimization Methods 0.000 description 2
- 238000003908 quality control method Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 102000008297 Nuclear Matrix-Associated Proteins Human genes 0.000 description 1
- 108010035916 Nuclear Matrix-Associated Proteins Proteins 0.000 description 1
- 238000012952 Resampling Methods 0.000 description 1
- 239000008186 active pharmaceutical agent Substances 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000007935 neutral effect Effects 0.000 description 1
- 210000000299 nuclear matrix Anatomy 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000010008 shearing Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000010998 test method Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/282—Application of seismic models, synthetic seismograms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Remote Sensing (AREA)
- General Physics & Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Environmental & Geological Engineering (AREA)
- Geophysics (AREA)
- Acoustics & Sound (AREA)
- Theoretical Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Artificial Intelligence (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开一种基于稀疏表示的地震各向异性参数反演方法,应用于地震数据处理领域,各向异性参数反演中的稳定性问题,本发明根据实测数据获得训练样本,基于K‑SVD方法从实际测井数据中学习六个参数各自的在线稀疏字典作为初始字典;再在反演过程中根据在线字典学习方法逐步更新字典,并将反演的各参数在对应的稀疏字典上的稀疏表达作为约束条件,从而进行各参数的反演,使得反演结果更趋向于测井数据,进而保证纵向分辨率的提高,使得反演结果更准确。
Description
技术领域
本发明属于地震数据解释领域,特别涉及一种地震各向异性参数反演技术。
背景技术
地震波各向异性反演是利用地表观测的叠前方位角道集地震资料,以已知地质规律和钻井、测井资料为约束,对地下岩层裂缝及各向异性特征的空间结构和进行成像的过程。各向异性反演为油气勘探中的储层预测提供了重要的依据。
各向异性参数需要和弹性参数同步反演,这些参数(后称六参数)同时反演时,其反演框架同叠前反演一致,是一种基于模型的叠前反演方法,是预测地层裂缝和各向异性参数的中有手段。该方法从初始地质模型出发,利用合成记录与实际地震数据之间的残差计算出的模型修正量来不断迭代更新初始模型,使模型正演合成记录与实际地震记录达到最佳吻合,最终得到的模型数据便是反演结果。
各向异性参数需要和弹性参数同步反演,这些参数同时反演时,由于其反演参数多、参数对地震数据影响微弱等原因导致反演结果存在极大地不稳定性。通常解决反演稳定性最直接方法就是加入先验信息约束,但是传统的先验约束信息都需要假设模型具有特定的结构(稀疏、光滑等)或具有一个唯一确定的数学表达形式(高斯或者柯西先验分布)。比如,增强反演结果光滑性最流行的正则化方法是Tikhonov型正则化,该类方法以最小化模型参数的低阶空间导数的l2范数作为先验约束;以及,2012年,Anagaw等基于全变差(TV)方法设计出能生成块状结构的弹性参数的反演方法,而这种方法是以最小化模型参数的梯度的l1范数为目标。这些先验信息通常只是一种单纯的数学抽象,只能适用于特定的工区,对地下信息的刻画并不完备,如果将这种先验信息用于各向异性反演必然导致反演结果存在误差。
由于传统的正则化方法没有考虑到测井数据中内含的地层特征,针对这一缺陷,基于稀疏表征的反演方法因此产生。稀疏表示技术是将信号在预先获得的一组基向量中进行分解,使得信号可以尽可能地由少量几个基向量线性表示。如今,稀疏表征理论在图像超分辨率重建,信号去噪,等领域中都取得成功的应用。
稀疏表征技术包括字典学习和稀疏重构。1996年,Olshausen等在《Nature》上发表的著名的Sparsenet字典学习方法奠定了字典学习的理论基础。其主要思想是从自然图像库中抽取大量的小图像块作为训练集,采用交替优化方法不断的学习训练更新字典,以获得信号的稀疏表示。2006年,Aharon提出了K-SVD方法,可很好地用于字典学习中;该方法使用正交匹配追踪技术,逐个考虑原子信号对数据拟合的贡献并逐个更新原子信号,目前,该方法已经成为字典学习应用最广泛的方法。稀疏重构就是将信号在预先学好的字典中进行稀疏分解,使得信号可以尽可能地由少量原子线性表示。稀疏分解实质是求一个l0范数的问题,求解方法主要有匹配追踪,正交匹配追踪等。
字典学习最早用于主成分分析方法中,2006年,Elad成功的将字典学习应用于图像去噪。2012年,Tang等引入了学习型超完备字典,将字典学习与总变差(TV)最小化结合起来用于地震数据去噪,并证明非均匀原子的字典更适合于地震数据去噪。2015年,Zhu改进了以往的字典方法,引入基于参数的字典学习方法用于地震数据去噪,在合成地震记录中有很好的效果。2016年,Chen将字典模型的分析方法和学习方法结合起来提出双稀疏字典学习方法用于地震数据去噪。在2017年,Nazari等通过多任务策略将二维的字典学习模型扩展到三维并用于三维地震数据去噪。同样地,2018年,Liu等将字典学习用于地震数据的稀疏表示,进而作为一个正则项用于地震数据去噪。前述在地震勘探中的字典学习方法均是用于数据去噪。2018年,Li首次将字典学习方法用于地震的全波形反演,提出基于多类字典学习的正则化框架,并在数值模拟中有较好的成果。同样地,在2018年,She等首次将字典学习方法运用于叠前AVO弹性参数反演,并提出将稀疏表示作为正则化项这一新方法,不管是在数值模拟还是在实际数据上均有很好的效果,这是反演历史的一个突破。但是She的方法只能利用初始的井数据学习字典,无法在反演过程中对字典进行更新。
发明内容
为解决各向异性参数反演中的稳定性问题,本发明提出了一种基于稀疏表示在线学习字典的地震各向异性参数反演方法,降低反演不稳定性。
本发明采用的技术方案为:一种基于稀疏表示的地震各向异性参数反演方法,根据同一工区中地层在横向上具有连续性,将各向异性介质中地震训练样本在字典上的稀疏表达作为约束条件进行反演,得到地层各向异性参数的反演结果,同时在反演过程中字典原子能够根据反演结果进行同步更新,避免因原始井数据较少导致的样本不足等问题。
进一步地,包括以下步骤:
S1、根据实际测井数据训练,具体的:根据实际的测井数据提取训练集,对训练集进行稀疏表示,基于K-SVD方法从实际测井数据中学习各参数对应的稀疏字典作为反演计算中所需要的初始字典;所述训练集中包括若干训练样本;
S2、通过对实际测井数据进行插值处理,建立各向异性介质中地震数据反演的初始模型;
S3、根据步骤S2建立的模型建立各向异性参数的反演目标函数,将一个地震剖面反演过程中,已反演好的结果再作为训练样本,进而实现在线更新字典,最后得到地层各向异性参数的反演结果;所述反演目标函数在同一工区中,地层在横向上具有连续性的假设下,将各参数在对应的稀疏字典上的稀疏表达作为约束条件。
进一步地,所述各参数包括各向异性参数及弹性参数;具体的各向异性参数包括ε,δ,γ,弹性参数包括:纵波速度Vp,横波速度Vs,密度ρ。
进一步地,步骤S3所述更新字典,具体为:
A1、将六个参数的反演结果进行分块处理,得到各参数的N个长度为M的训练样本集;
A2、设置正则化参数以及最大迭代次数;
A3、以步骤S1中从实际测井数据中学习得到的各参数对应的稀疏字典作为各参数对应的初始字典,并根据各参数的初始字典确定对应的初始化的稀疏系数矩阵;
A4、假设各字典中的原子固定,采用OMP方法对参数的训练样本进行稀疏编码,得到参数对应的稀疏系数矩阵;
A5、对字典中的原子进行逐个迭代更新,从而得到各参数对应的稀疏字典以作为反演计算中所需要的字典。
更进一步地,步骤S1包括以下分步骤:
S11、将实际测井数据进行分块处理,得到各参数的N个长度为M的训练样本集;
S12、设置正则化参数以及最大迭代次数;
S13、初始化各参数对应的字典和稀疏系数矩阵;
S14、假设各字典中的原子固定,采用OMP方法对参数的训练样本进行稀疏编码,得到参数对应的稀疏系数矩阵;
S15、对字典中的原子进行逐个迭代更新,从而得到各参数对应的稀疏字典以作为反演计算中所需要的初始字典。
进一步地,步骤S15所述逐个原子更新具体为:
B1、将当前要更新的第k列原子记为dk,记Ik={i|αi(k)≠0,1≤i≤N},i表示第i个测井数据的小块,αi(k)为αi中的第k个元素,则Ik表示所有样本中使用到第k列原子dk的索引集;
B2、用表示去除字典中第k列原子后的字典,表示去除第k行系数以后的系数矩阵,然后计算与与该参数训练样本集之间的残差矩阵Ek表示除去第k个原子后对所有样本造成的误差;
B3、根据Ik中的索引值选取Ek中相应的列向量,构成新的误差矩阵并对进行奇异值分解
B4、取矩阵U中的第一列,作为更新后的原子,并更新相应的系数。
本发明的有益效果:本发明将图像处理领域中应用广泛的在线字典学习与稀疏表征方法应用于地震记录叠前宽方位角道集的各向异性参数反演,相比于传统反演方法本发明方法利用反演结果对字典在线更新,使得工区具备自适应性的优点;本发明方法采用在线字典稀疏重构作为反演的约束项以代替传统正则化约束项,使得反演结果更趋向于实际观测数据,并能够根据实际反演数据结果进行在线调整,进而保证纵向分辨率的提高,使得本发明的反演结果可以更好的为储层预测提供帮助。
附图说明
图1为本发明的方案流程图;
图2为本发明实施例一提供的根据测井数据训练的字典DP、Dε的前30个原子;
其中,图2(a)为字典DP的前30个原子,图2(b)为字典Dε的前30个原子;
图3为本发明实施例一提供的方位角为30度的合成角道集和具有带限噪声的道集;
其中,图3(a)为EAGE模型#179在方位角等于30度时的无噪角道集地震记录,图3(b)为含噪地震记录;
图4为本发明实施例一提供的传统Tikhnov正则化和本发明方法的反演结果对比图;
其中,图4(a)为传统Tikhnov正则化的弹性参数Vp反演结果,图4(b)为传统Tikhnov正则化的弹性参数Vs反演结果,图4(c)为传统Tikhnov正则化的弹性参数ρ反演结果,图4(d)为本发明方法的弹性参数Vp反演结果,图4(e)为本发明的弹性参数Vs反演结果,图4(f)为本发明的弹性参数ρ反演结果,图4(g)为传统Tikhnov正则化的各向异性参数ε反演结果,图4(h)为传统Tikhnov正则化的各向异性参数δ反演结果,图4(i)为传统Tikhnov正则化的各向异性参数γ反演结果,图4(j)为本发明方法各向异性参数ε反演结果,图4(k)为本发明方法各向异性参数δ反演结果,图4(l)为本发明方法各向异性参数γ反演结果;
图5为本发明实施例二提供的根据测井数据训练的字典DP、Dε的前30个原子;
其中,图5(a)为字典DP的前30个原子,图5(b)为字典Dε的前30个原子;
图6为本发明实施例二提供的各向异性参数反演的过井剖面;
其中,图6(a)为参数ε反演的过井剖面,图6(b)为参数δ反演的过井剖面,图6(c)为参数γ反演的过井剖面。
具体实施方式
为便于本领域技术人员理解本发明的技术内容,下面结合附图对本发明内容进一步阐释。
现有的反演方法采用正则化约束,正则化的假设需要专家对工区地层特征的先验信息,通常这些方法不具有普遍的适应性;并且没有考虑参数本身所蕴含的地层特征。本发明提出的基于测井数据驱动的各向异性参数和弹性参数(共六参数)同时反演方法,不需要更多的先验信息,只需要同一工区中有足够的测井数据用于训练字典,并且再反演计算中,运用在线字典学习得方法逐步更新字典,获得更多得地层特征。
如图1所示,本发明的方法包括以下步骤:
S1、根据实际测井数据训练,具体的:根据实际的测井数据提取训练集,对训练集进行稀疏表示,基于K-SVD方法从实际测井数据中学习各参数对应的稀疏字典作为反演计算中所需要的初始字典;所述训练集中包括若干训练样本;
S2、通过对实际测井数据进行插值处理,建立各向异性介质中地震数据反演的初始模型;
S3、根据步骤S2建立的模型建立各向异性参数的反演目标函数,将一个地震剖面反演过程中,已反演好的结果再作为训练样本,进而更新字典,再得到地层各向异性参数的反演结果;所述反演目标函数在同一工区中,地层在横向上具有连续性的假设下,将六个参数在各参数对应的稀疏字典上的稀疏表达作为约束条件。
本发明方法分为两个步骤:字典学习和稀疏反演;
字典学习阶段分为两个部分,一个是通过测井数据学习初始字典,即测井数据分块得到训练数据,再运用K-SVD等方法学得初始字典,该字典可表述测井数据中蕴含的地层特征。(She et al.,2018)对用测井数据做字典学习的过程进行了详细的描述。第二部分,在反演过程中,将已反演的道的结果去替换原来的样本集,根据初始字典再重新更新字典,利用该字典对待反演的弹性参数和各向异性参数的稀疏表达做为约束再加入到反演过程中。
字典学习和稀疏表示在学术界的称谓是稀疏字典学习。稀疏表示理论最早是在研究信号处理应用中发展起来的。其基础是多尺度分析理论,在此基础上拓展,形成了相应的理论框架。字典学习的最终任务是从许多信号中推知原子信号以及各个信号的原子信号权重,它的实质是对大量数据信号的一种降维表示,它尝试找出隐藏在信号数据中每个信号所共有的特征。稀疏表示则是可以从原子信号中重构出原始信号,即用尽可能少的原子信号来表示原始信号,这种表示还能带来一个附加的好处,即计算速度快。近年来稀疏字典学习在图像处理领域中得到广泛应用。本发明在反演过程中,继续根据已反演好的道的结果作为样本集,用之前训练好的字典作为初始字典,再用K-SVD方法得到更新的字典,用更新的字典来做反演计算。以下为本发明采用在线稀疏字典学习的具体实现过程:
在根据测井数据训练字典时,本发明需要将重采样后的测井数据分成小块(是一个一维的数据),以作为原始信号,每个信号的长度为M,M的取值范围一般在30到100之间,可根据实际情况确定;然后用这些小块合并作为训练样本集。设训练样本集为R为Y中训练样本的个数。字典学习过程的目标便是找到一个字典能以最稀疏的方式表达训练集Y,若记此时的稀疏系数集为则该过程可以被表示为如下能量最小化问题:
其中,每一个长度为L的列向量αi为样本yi在字典D下的稀疏表示,||·||0是l0-范数,||·||0用以计算某一个向量的含零元素的个数,K是表示稀疏程度的常量。在式(1)中,其代价函数:确立了对训练集的分解的一致性,而约束项则限制了表达结果的稀疏性。
下面给出基于K-SVD的测井数据的初始字典学习方法,在反演过程中的在线字典学习流程与该流程一致,所不同的是样本数据需要逐步增加反演迭代的结果。
(1)对已知的每条测井曲线的数据分别进行分块处理,得到不同参数的R个长度为M的训练样本Y=[y1,y2,…,ym];
(2)设置正则化参数λ,最大迭代次数T;
(3)初始化字典D和稀疏系数矩阵α,置当前迭代次数t=1;
(4)稀疏编码:假设原子D固定,用OMP方法对训练样本集Y进行稀疏编码,得到稀疏系数矩阵α。
(5)字典更新:对字典中的原子进行逐个更新。对于第k(k=1,2,…K)个原子,按照下述步骤更新:
①dk为当前要进行更新的原子,记Ik={i|αi(k)≠0,1≤i≤N},αi(k)为αi中的第k个元素,则Ik表示所有样本中使用到第k列原子dk的索引集。
②用表示去除字典中第k个原子后的字典,表示去除第k行系数以后的系数矩阵,然后计算新的字典与系数矩阵与原信号之间的残差矩阵Ek表示除去第k个原子后对所有样本造成的误差;
③根据Ik中的索引值选取Ek中相应的列向量,构成新的误差矩阵并对进行奇异值分解
④取矩阵U中的第一列,即最大特征值对应的特征向量,作为更新后的原子,并更新相应的系数。
(6)置t=t+1。检查t是否大于T。如果t≤T,则返回步骤(4);如果t>T,则停止迭代,输出字典D。
在反演迭代计算过程中,需要将已反演好的道的结果替换样本集,将测井数据训练的字典作为初始字典,再根据上述K-SVD过程更新字典,过程通步骤(1)至(6),并将更新后的字典来约束反演。
AVAZ反演框架的构建
各向异性假设下双层HTI介质中纵波反射系数近似表达式;
其中,i、φ分别是入射角和方位角;α、β分别为纵波和横波速度;Z为纵波波阻抗Z=ρα;G为横波剪切模量,G=ρβ2;运算符Δ[·]为上下界面的参数之差,为上下界面参数之均值。当φ=π/2时,该公式就可退化成HTI介质中的与方位角无关的各向同性反射系数计算公式(3),Rpiso称为各向同性项:
对于一列地质模型,各向同性介质中的反射系数计算公式(3),利用Aki-Richards近似方程替换为:
其中,c1=1+tan2i,c2=-8ν2tan2i,ν=β/α,c3=-0.5tan2i+2ν2sin2i
同样地,对于与方位角有关的项可写成式(5),Rpani称为各向异性项;
Rpani(φ)=c4ε+c5δ+c6γ (5)
其中以及
因此,HTI介质中,反射系数可写为:
Rp=Rpiso+Rpani (6)
单道地震记录d(t)可以通过子波w(t)与反射系数卷积而得:
式(7)中前一项描述各向同性介质下,反射系数可以写成一个线性方程组(Hampson et al.,2005),并考虑到P波阻抗的对数和S波阻抗、密度的对数之间存在着线性关系,建立数学模型如下:
dpiso=Gisomiso (8)
其中,Giso表示由子波和权系数组成的核矩阵,而miso表示模型参数。
其中,N为地质模型采样点数,表示P波阻抗的对数,k表示第k层,和是如下式中对参数进行线性拟合后的剩余偏差:
其中,和分别表示S波阻抗和密度的对数,a,b,c,e为常线性拟合系数,可以从测井曲线中学习得到。
对于公式(7)中后一项,同理可将其写为矩阵形式;
dpani=Ganimani (12)
因此,各向异性介质中地震记录为:
d=dpiso+dpani=Gm (13)
其中
G=[Giso,Gani] (14)
m=[miso,mani]T (15)
根据(12),已知叠前方位纵波地震数据d,和矩阵G,估算地下参数m,此过程称为反演,通常,该过程目标函数如下:
其中,λ是正则化参数;
现有技术中采用正则化方法来约束反演。比如,Tikhonov正则化和全变差正则化。Tikhonov正则化通过最小化模型m的梯度的二范使得解具有光滑性:
全变差正则化通过假设模型m具有稀疏的梯度,使得解呈块状特征:
式(17)(18)中T为一阶差分算子。
现有技术采用正则化约束,存在各向异性参数反演中的稳定性差的问题,本发明提出一种能根据数据本身构建自适应的先验模型的方法以实现更各向异性参数更精确的反演,降低反演不稳定性。
由于地壳运动的尺度很大,故同一工区中,地层通常在横向上具有连续性。上述通过测井学得的字典使得每条测井数据都能在该字典上稀疏表达。由于同一工区中,地层通常在横向上具有连续性,故对于每个地震道的反演结果通常也能在该字典上稀疏表达。根据这条思路,将模型在字典上的稀疏表达作为约束条件来约束反演。
将(16)式重写如下:
其中记为:
其中,Vp,Vs,ρ分别表示纵波速度,横波速度,以及密度。n是经过块处理后模型的块的数量,下标P,S,ρ分别表示纵波,横波和密度;ε,δ,γ为各向异性参数。以Vp为例,表示从上述稀疏字典学习过程学得的具有L个原子的字典,表示训练样本的第i个小块,它是长度为M的一个列向量,相应地,每一个列向量是通过字典DP对的稀疏表示结果。上述描述对Vs,ρ,ε,δ,γ也是一样的。此外,μ1,μ2,μ3μ4,μ5和μ6表示六个合适的权重参数。算子“PC”是“正变换”,只是对Vp,Vs,ρ,对于各向异性参数不需要做变换。Ri是一个6M×6N的矩阵,用于提取模型参数m的第i个小块。
分析目标函数,第一项表示反演数据和真实数据的误差,后面六项则为反演结果在字典上稀疏表达的整体误差。约束项是使得反演参数能在各自对应字典上进行稀疏表示。等式约束则为保证反演模型m能转化成弹性参数。
求解(19)式过于困难,第一,参数太多;第二,等式约束复杂,弹性参数和各向异性参数组成的模型m,使得各参数必须得同时考虑。因此,本发明引入中间变量x,将目标函数写成式(21),约束项不变,式中β是一个新的正则化参数。
采用坐标下降方法,将求解式(21),求解过程分为以下两步:
1)求解目标函数(22)
当k=0时,设置x0=m0(初始模型),该步骤是一个凸优化问题,求解方法众多,如梯度下降方法,拟牛顿方法等。
2)求解目标函数(23)
式(23)等价于对每个参数求形如式(23)的目标函数:
其中j=1,2,3,4,5,6,依次代表需要反演的六个参数。
综上,本发明的方法流程为:
1.提取地震子波并进行井震标定,即,将测井数据和地震数据进行匹配,保证测井数据和地震数据的具有相同的采样间隔以及采样时间;
2.选取合适的原子长度和原子个数,然后运用K-SVD方法从实际测井数据中学习出稀疏字典D以作为初始字典;
3.通过测井数据来插值,建立反演计算过程中所要用到的初始值;
4.选取合适的正则化参数λ,和稀疏度K;
5.在反演计算过程中,将已反演道得结果作为样本集,重新更新字典,使得字典能逐步学习到储层参数结构特征。
6.利用拟牛顿方法和OMP方法最小化目标函数,得到地层各向异性参数的反演结果。
原子长度和正则化参数的选取需要根据质控结果来确定,所谓质控,即质量控制,用不同的参数去反演,根据反演效果评估确定最优参数。
以下结合具体实例对本发明的内容做进一步阐述
实施例1
本发明用EAGE模型来测试,另外,本发明还对比了由1th-Tikh方法生成的反演结果。
给定一个雷克子波,再利用(3)式,便可以生成方位角0度、30度、60度、90度、入射角从0-36度的对应的合成角道集。对着该模型的每一道而言,本发明都添加了SNR=4:1的高斯白噪声,然后,本发明通过一个带通滤波器获得了最终的带限的合成观测数据。
通过采样和时深转换处理后,本发明采用时间域的大小为200×399的EAGE模型来测试方法的有效性。按照上述合成方式,本发明用主频为15Hz的雷克子波生成了不同方位的角道集。与其他方法不同,本发明的方法首先需要学习字典DP、DS、Dρ、Dε、Dδ、Dγ,因此本发明从399道数据中等间距地提取了20道井数据作为训练数据。此次训练的矩阵的大小是50×200,即,每个字典有200个原子,每个原子的采样点数为50。图2展示了字典DP、Dε的前30个原子,从该图中可以看出,每一个原子都继承了当前地质模型的块状特征。
接下来,本发明首先专注于那些参与训练的道。以道#80为例,图3展示了如图3(a)所示方位角为30度的合成角道集和如图3(b)所示具有带限噪声的道集。图4(a)、(b)(c)展示了传统Tikhnov正则化的各项同性参数(VP,Vs,ρ)的反演结果,图4(d)、(e)(f)展示了本发明方法的各项同性参数(VP,Vs,ρ)的反演结果,图4(g)、(h)、(i)展示了传统Tikhnov正则化的各向异性参数(ε,δ,γ)的反演结果,图4(j)、(k)、(l)展示了本发明方法的各向异性参数(ε,δ,γ)的反演结果;显然,显然本发明方法比传统方法有更精确的反演结果。
实施例2
本发明采用中国西部某工去的实际数据测试本发明的发明,根据测井数据分别训练弹性参数和各向异性参数共六个参数的稀疏字典,图5(a)和图5(b)分别表示用测井数据训练得纵波速度和ε的前40个原子,并且将发明的方法应用于该工区的过井剖面,反演结果如图6所示,从图6(a)、(b)、(c)中可以看出,该技术能较好的反演出地下各向异性的情况,说明了本发明的方法的有效性。
本领域的普通技术人员将会意识到,这里所述的实施例是为了帮助读者理解本发明的原理,应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说,本发明可以有各种更改和变化。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的权利要求范围之内。
Claims (6)
1.一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,根据同一工区中,地层在横向上具有连续性,将各向异性介质中各参数在字典上的稀疏表达作为约束条件进行反演,得到地层各向异性参数的反演结果。
2.根据权利要求1所述的一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,包括以下步骤:
S1、根据实际测井数据训练,具体的:根据实际的测井数据提取训练集,对训练集进行稀疏表示,基于K-SVD方法从实际测井数据中学习各参数对应的稀疏字典作为一个初始的字典供反演计算中使用;所述训练集中包括若干训练样本;
S2、通过对实际测井数据进行插值处理,建立各向异性介质中地震数据反演的初始模型;
S3、根据步骤S2建立的模型建立各向异性参数的反演目标函数,将一个地震剖面反演过程中,已反演好的结果再作为训练样本,进而更新字典,最后得到地层各参数的反演结果;所述反演目标函数在同一工区中,地层在横向上具有连续性的假设下,将各参数在对应的在线稀疏字典上的稀疏表达作为约束条件。
3.根据权利要求2所述的一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,所述各参数包括各向异性参数及弹性参数;具体的各向异性参数包括ε,δ,γ,弹性参数包括:纵波速度Vp,横波速度Vs,密度ρ。
4.根据权利要求3所述的一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,步骤S3所述更新字典,具体为:
A1、将六个参数的反演结果进行分块处理,得到各参数的N个长度为M的训练样本集;
A2、设置正则化参数以及最大迭代次数;
A3、以步骤S1中从实际测井数据中学习得到的各参数对应的稀疏字典作为各参数对应的初始字典,并根据各参数的初始字典确定对应的初始化的稀疏系数矩阵;
A4、假设各字典中的原子固定,采用OMP方法对参数的训练样本进行稀疏编码,得到参数对应的稀疏系数矩阵;
A5、对字典中的原子进行逐个迭代更新,从而得到各参数对应的稀疏字典以作为反演计算中所需要的初始字典。
5.根据权利要求3所述的一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,步骤S1包括以下分步骤:
S11、将六个参数的实际测井数据进行分块处理,得到各参数的N个长度为M的训练样本集;
S12、设置正则化参数以及最大迭代次数;
S13、初始化各参数对应的字典和稀疏系数矩阵;
S14、假设各字典中的原子固定,采用OMP方法对参数的训练样本进行稀疏编码,得到参数对应的稀疏系数矩阵;
S15、对字典中的原子进行逐个迭代更新,从而得到各参数对应的稀疏字典以作为反演计算中所需要的初始字典。
6.根据权利要求4或5任一权利要求所述的一种基于稀疏表示的地震各向异性参数反演方法,其特征在于,所述逐个原子更新具体为:
B1、将当前要更新的第k列原子记为dk,记Ik={i|αi(k)≠0,1≤i≤N},i表示第i个测井数据的小块,αi(k)为αi中的第k个元素,则Ik表示所有样本中使用到第k列原子dk的索引集;
B2、用表示去除字典中第k列原子后的字典,表示去除第k行系数以后的系数矩阵,然后计算与与该参数训练样本集之间的残差矩阵Ek表示除去第k个原子后对所有样本造成的误差;
B3、根据Ik中的索引值选取Ek中相应的列向量,构成新的误差矩阵并对进行奇异值分解
B4、取矩阵U中的第一列,作为更新后的原子,并更新相应的系数。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910800069.5A CN110516358A (zh) | 2019-08-28 | 2019-08-28 | 一种基于稀疏表示的地震各向异性参数反演方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910800069.5A CN110516358A (zh) | 2019-08-28 | 2019-08-28 | 一种基于稀疏表示的地震各向异性参数反演方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN110516358A true CN110516358A (zh) | 2019-11-29 |
Family
ID=68628383
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910800069.5A Pending CN110516358A (zh) | 2019-08-28 | 2019-08-28 | 一种基于稀疏表示的地震各向异性参数反演方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110516358A (zh) |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111368247A (zh) * | 2020-03-12 | 2020-07-03 | 电子科技大学 | 基于快速正交字典的稀疏表征正则化叠前avo反演方法 |
CN111965697A (zh) * | 2020-08-18 | 2020-11-20 | 电子科技大学 | 基于联合字典学习和高频预测的高分辨率地震反演方法 |
CN113534250A (zh) * | 2020-04-18 | 2021-10-22 | 中国石油化工股份有限公司 | 一种基于快速匹配追踪的多尺度地震反演方法 |
CN114063156A (zh) * | 2020-07-30 | 2022-02-18 | 中国石油天然气股份有限公司 | 基于稀疏表征的地震各向异性梯度反演方法及系统 |
CN115166822A (zh) * | 2022-08-02 | 2022-10-11 | 中国矿业大学(北京) | 储层弹性参数的预测方法、装置和电子设备 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105487110A (zh) * | 2014-09-16 | 2016-04-13 | 中国石油化工股份有限公司 | 一种基于透射方程的各向异性参数反演方法 |
WO2016065356A1 (en) * | 2014-10-24 | 2016-04-28 | Ion Geophysical Corporation | Methods for seismic inversion and related seismic data processing |
CN106597541A (zh) * | 2017-02-22 | 2017-04-26 | 中国石油大学(华东) | 基于Shearlet变换的地震数据重构方法 |
US20170145793A1 (en) * | 2015-08-20 | 2017-05-25 | FracGeo, LLC | Method For Modeling Stimulated Reservoir Properties Resulting From Hydraulic Fracturing In Naturally Fractured Reservoirs |
CN109143356A (zh) * | 2018-08-29 | 2019-01-04 | 电子科技大学 | 一种自适应混合范数字典学习地震波阻抗反演方法 |
-
2019
- 2019-08-28 CN CN201910800069.5A patent/CN110516358A/zh active Pending
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105487110A (zh) * | 2014-09-16 | 2016-04-13 | 中国石油化工股份有限公司 | 一种基于透射方程的各向异性参数反演方法 |
WO2016065356A1 (en) * | 2014-10-24 | 2016-04-28 | Ion Geophysical Corporation | Methods for seismic inversion and related seismic data processing |
US20170145793A1 (en) * | 2015-08-20 | 2017-05-25 | FracGeo, LLC | Method For Modeling Stimulated Reservoir Properties Resulting From Hydraulic Fracturing In Naturally Fractured Reservoirs |
CN106597541A (zh) * | 2017-02-22 | 2017-04-26 | 中国石油大学(华东) | 基于Shearlet变换的地震数据重构方法 |
CN109143356A (zh) * | 2018-08-29 | 2019-01-04 | 电子科技大学 | 一种自适应混合范数字典学习地震波阻抗反演方法 |
Non-Patent Citations (4)
Title |
---|
BIN SHE,ET AL: "Seismic Impedance Inversion using Dictionary Learning and Sparse Representation", 《2018 SEG INTERNATIONAL EXPOSITION AND 88TH ANNUAL MEETING》 * |
JULIEN MAIRAL,ET AL: "Online Learning for Matrix Factorization and Sparse Coding", 《JOURNAL OF MACHINE LEARNING RESEARCH》 * |
程三: "基于构造约束的地震信号多道反演方法研究", 《中国优秀硕士学位论文全文数据库 基础科学辑》 * |
陈怀震: "基于岩石物理的裂缝型储层叠前地震反演方法研究", 《中国博士学位论文全文数据库 基础科学辑》 * |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111368247A (zh) * | 2020-03-12 | 2020-07-03 | 电子科技大学 | 基于快速正交字典的稀疏表征正则化叠前avo反演方法 |
CN111368247B (zh) * | 2020-03-12 | 2021-11-30 | 电子科技大学 | 基于快速正交字典的稀疏表征正则化叠前avo反演方法 |
CN113534250A (zh) * | 2020-04-18 | 2021-10-22 | 中国石油化工股份有限公司 | 一种基于快速匹配追踪的多尺度地震反演方法 |
CN114063156A (zh) * | 2020-07-30 | 2022-02-18 | 中国石油天然气股份有限公司 | 基于稀疏表征的地震各向异性梯度反演方法及系统 |
CN111965697A (zh) * | 2020-08-18 | 2020-11-20 | 电子科技大学 | 基于联合字典学习和高频预测的高分辨率地震反演方法 |
CN111965697B (zh) * | 2020-08-18 | 2021-09-28 | 电子科技大学 | 基于联合字典学习和高频预测的高分辨率地震反演方法 |
CN115166822A (zh) * | 2022-08-02 | 2022-10-11 | 中国矿业大学(北京) | 储层弹性参数的预测方法、装置和电子设备 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
EP3894907B1 (en) | Machine learning-augmented geophysical inversion | |
CN110516358A (zh) | 一种基于稀疏表示的地震各向异性参数反演方法 | |
US11397272B2 (en) | Data augmentation for seismic interpretation systems and methods | |
US10996372B2 (en) | Geophysical inversion with convolutional neural networks | |
Yang et al. | Random noise attenuation based on residual convolutional neural network in seismic datasets | |
US11435498B2 (en) | Subsurface models with uncertainty quantification | |
CN106405651B (zh) | 一种基于测井匹配的全波形反演初始速度模型构建方法 | |
CN109143356A (zh) | 一种自适应混合范数字典学习地震波阻抗反演方法 | |
WO2020123099A2 (en) | Automated seismic interpretation-guided inversion | |
US11181653B2 (en) | Reservoir characterization utilizing ReSampled seismic data | |
Yang et al. | Denoising of distributed acoustic sensing data using supervised deep learning | |
CN113077386A (zh) | 基于字典学习和稀疏表征的地震资料高分辨率处理方法 | |
Min et al. | D 2 UNet: Dual decoder U-Net for seismic image super-resolution reconstruction | |
Wu et al. | Variable seismic waveforms representation: Weak-supervised learning based seismic horizon picking | |
Khosro Anjom et al. | Machine learning for seismic exploration: Where are we and how far are we from the holy grail? | |
Li et al. | Removing abnormal environmental noise in nodal land seismic data using deep learning | |
Bocheng et al. | Wavefield separation of VSP data based on deep neural network | |
CN114063156B (zh) | 基于稀疏表征的地震各向异性梯度反演方法及系统 | |
Athmer et al. | Integrating Seismic Interpretation, Classification and Geologic Process Modeling for Shale Reservoir Characterization | |
Chen et al. | Self-supervised Multistep Seismic Data Deblending | |
CN109471167A (zh) | 一种针对多震源缺失数据的波场重构反演方法 | |
Pérez et al. | Three-Dimensional Initial Salt Body Building for Full-Waveform Inversion Assisted by Deep Learning | |
Xu et al. | DAS-VSP coupled noise suppression based on U-Net network | |
CN114063156A (zh) | 基于稀疏表征的地震各向异性梯度反演方法及系统 | |
CN109212602A (zh) | 一种改进地震数据分辨率的反射系数反演方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
RJ01 | Rejection of invention patent application after publication | ||
RJ01 | Rejection of invention patent application after publication |
Application publication date: 20191129 |