CN113534250A - 一种基于快速匹配追踪的多尺度地震反演方法 - Google Patents
一种基于快速匹配追踪的多尺度地震反演方法 Download PDFInfo
- Publication number
- CN113534250A CN113534250A CN202010308485.6A CN202010308485A CN113534250A CN 113534250 A CN113534250 A CN 113534250A CN 202010308485 A CN202010308485 A CN 202010308485A CN 113534250 A CN113534250 A CN 113534250A
- Authority
- CN
- China
- Prior art keywords
- scale
- seismic
- inversion
- frequency
- longitudinal wave
- 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 63
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 25
- 230000004927 fusion Effects 0.000 claims abstract description 8
- 230000000903 blocking effect Effects 0.000 claims abstract description 6
- 239000011159 matrix material Substances 0.000 claims description 28
- 238000012937 correction Methods 0.000 claims description 9
- 238000004458 analytical method Methods 0.000 claims description 6
- 230000010354 integration Effects 0.000 claims description 6
- 238000001228 spectrum Methods 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 claims description 4
- 238000013507 mapping Methods 0.000 claims description 4
- 230000008569 process Effects 0.000 claims description 4
- 239000013598 vector Substances 0.000 claims description 4
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims description 3
- ONUFESLQCSAYKA-UHFFFAOYSA-N iprodione Chemical compound O=C1N(C(=O)NC(C)C)CC(=O)N1C1=CC(Cl)=CC(Cl)=C1 ONUFESLQCSAYKA-UHFFFAOYSA-N 0.000 claims description 2
- 230000001105 regulatory effect Effects 0.000 claims description 2
- 230000017105 transposition Effects 0.000 claims description 2
- 230000015572 biosynthetic process Effects 0.000 claims 1
- 238000005070 sampling Methods 0.000 description 9
- 238000012360 testing method Methods 0.000 description 5
- 238000012512 characterization method Methods 0.000 description 3
- 238000009499 grossing Methods 0.000 description 3
- 239000004215 Carbon black (E152) Substances 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 229930195733 hydrocarbon Natural products 0.000 description 2
- 150000002430 hydrocarbons Chemical class 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000007493 shaping process Methods 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 230000008033 biological extinction Effects 0.000 description 1
- 230000001276 controlling effect Effects 0.000 description 1
- 238000012217 deletion Methods 0.000 description 1
- 230000037430 deletion Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000036039 immunity Effects 0.000 description 1
- 230000009191 jumping Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 238000000638 solvent extraction Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Images
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
- G01V1/30—Analysis
- G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
-
- 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
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/40—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/622—Velocity, density or impedance
- G01V2210/6226—Impedance
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开了一种基于快速匹配追踪的多尺度地震反演方法。包括：进行地震信号多尺度分解的分析，对不同尺度地震信号如何划分进行分析，并同时对不同尺度地震信号进行分析。快速匹配追踪块化地层纵波阻抗反演方法，对不同信噪比情况下的迭代阈值与迭代次数设置进行分析用以获得最佳反演结果。多尺度快速匹配追踪纵波阻抗反演方法，在低频模型约束基础上进行大尺度快速匹配追踪反演以获得富含低频信息的纵波阻抗稀疏反演结果，在此基础上，利用大尺度纵波阻抗反演结果作为低频模型进行小尺度纵波阻抗稀疏反演，最后对不同尺度反演结果进行多尺度融合以获得多尺度纵波阻抗稀疏反演结果。大幅提高了反演结果的高低频信息，提升了反演方法的分辨率。
Description
技术领域
本发明涉及地球物理地震反演技术领域，具体是一种基于快速匹配追踪的多尺度地震反演方法，通过对地震数据进行多尺度分解后利用多尺度地震数据与匹配一种对地震和测井数据波形极值特征点分离与波形重构方法
背景技术
地震反演问题的核心在于何如利用已获得的观测数据来反推地球内部的地球物理模型特征。但在实际的地震信号传播过程中，可将地球内部看作一带通滤波器，即人们在地表观测所得的地震信号呈现带限特征，因此地球物理反演方法呈现出普遍的病态性，这中病态直接表现为观测数据远大于未知模型数据，而观测数据为非线性独立，使得传统地球物理反演方法呈现出普遍的多解性。此外，对于勘探目标的逐渐减小，勘探对象逐渐由构造油气藏变换为隐藏油气藏，分辨率的问题严重制约着薄储层的刻画。
反演多解性及薄储层的刻画一直以来困扰着业界学者，为了精确的描述薄储层边界，近年来基追踪、匹配追踪等地震反演方法成为业界研究重点对象。1993年Mallat等人首次提出全局搜索匹配追踪算法，通过对时频原子字典进行全域搜索以达到重构地震信号的目的，但耗时长、计算效率低等问题严重制约了其发展。在此基础上，张翻昌(2010)等人研究了利用地震三瞬属性对最佳匹配原子进行限制从而缩小搜索范围以提高匹配追踪计算效率。除此之外，匹配追踪算法在地震勘探领域中具有广泛的应用，其具体为利用匹配追踪算法进行谱分解(WANG Y，2010；张繁昌等，2012)以及去强反射(李海山等，2014)，以及利用瞬时谱进行砂体尖灭点识别(赵天姿等，2008)。
CN103698808B共开了《一种对地震和测井数据波形极值特征点分离与波形重构方法》。其技术方案包括：
1)野外地震、测井数据，得到经过叠后或叠前处理地震波形和预处理后的测井曲线，截取目标时窗和深窗分析范围，并对地震或测井采样数据进行多点平滑、整形处理；所述的平滑、整形处理包括波形上的棘刺、高频跳跃点、拐点和拐折点形成待分析的初始级波形信号。
2)确定一维信号的中央基准线；所述的中央基准线对地震道记录数据，直接以零基线作为中央基准线或铅垂线，零基线两侧区域为地震道反射波的正、负极性，分别对应波峰与波谷反射能量的左、右两个区带；沿初始波形曲线深度方向、逐点自上向下先区分出左、右极值点集并采用多项式平滑拟合出左、右初始极值包络线；然后按照测井原始采样间隔连续计算每个采样位置上穿越左、右包络线两极值点的中点值，连接所有垂向中点值、采取多项式平滑拟合获得中央基准线；中央基准线为基准划分出测井两个物理值左、右分布区域，完成对测井信号物理属性左、右采样极值点分区。所述的区分是：设函数y＝f(x)在x0的一个邻域内有定义，若对于该邻域内异于x0的x恒有：(1)f(x0)>f(x)，则f(x0)为函数f(x)的极大值，x0为f(x)的极大值点；(2)f(x0)<f(x)，则f(x0)为函数f(x)的极小值，x0为f(x)的极小值点；极大值点、极小值点为极值点。
3)分别求取初始级波形上位于中央基准线左侧的极大值点集(绝对值)和其右侧的极大值点集(绝对值)，通过相邻左侧极大值点和右侧极大值点的线性插值连续求取出与中央基准线相交的地震波形回零点集、测井曲线物理平衡点集，再对垂向上相邻左、右极值点集和回零点集或物理平衡点集进行波形重构获得一级波形；将波形重构后的一级波形进行三次样条函数或多项式插值方法拟合,并使得一级波形曲线采样点数与初始波形曲线采样点数相等，获得光滑的、合理的1D一级波形序列及对应的波形剖面；所述的波形重构按照采样点的深度值或反射时间由大到小、从下往上逐个将垂向相邻的特征点连接；所述的连接按照地震或测井信号特点分以下两种情形：①当相邻两极值特征点为异号交替出现时，两点直接连线并与中央基准线相交，线性插值后可获得多个新的回零点或物理属性变化平衡点，重构波形是单波；②当局部相邻两极值特征点为同号时，两点连线时，需要在它们之间构建一个左或右极小值点，并用两同号极值特征点值的算术和的三分之一进行赋值，重构的波形是复合波。
4)重复步骤3)，把一级波形曲线当作初始级波形再进行极值点集筛选、插值与二级新波形的重构；多次重复，获得多级一维新波形；所述的重复次数依据原始信号的频率、采样率和被分析目标的大小而定，通常情况下进行3次以上。对地震单道逐道地进行上述相同的变换方式，多道间采用横向多点平滑技术，就可获得2D分频波形剖面及3D分频波形数据体。
但是该反演方法由于频带带限特征导致地震信号(也称地震数据)高频下降以及低频缺失问题，此外，反演多解性也依然存在。
发明内容
本发明的目的为缓解反演方法多解性同时增加反演高低频信息，而提出一种基于快速匹配追踪的多尺度地震反演方法。本方法充分考虑地震信号中蕴含的“频带”特征，将地震信号进行多尺度分解，并采取逐级迭代的思想进行反演。
本发明的技术方案是：
一种基于快速匹配追踪的多尺度地震反演方法，包括：
(1)进行地震信号多尺度分解，包括对不同尺度地震信号划分，对划分后的不同尺度地震信号进行分析；
(2)建立快速匹配追踪块化地层纵波阻抗反演方法，包括对不同信噪比情况下的迭代阈值与迭代次数设置；
(3)建立多尺度快速匹配追踪纵波阻抗反演方法，包括在低频模型约束基础上进行大尺度快速匹配追踪反演以获得富含低频信息的纵波阻抗稀疏反演结果，在此基础上，利用大尺度纵波阻抗反演结果作为低频模型进行小尺度纵波阻抗稀疏反演，最后对不同尺度反演结果进行多尺度融合。
上述方案进一步包括：
所述步骤(1)地震信号多尺度分解是：根据地下介质岩性不同导致地震数据表现出不同的频率成分，利用多尺度地震数据对地震层序进行有效刻画，对于多尺度分解后的地震信号进行分析可知，大尺度地震信号为地震信号低频信息，小尺度地震信号为地震信号高频信息，且小尺度地震信号镶嵌于大尺度地震信号中。
所述步骤(2)快速匹配追踪块化地层纵波阻抗反演方法设置包括：对于信噪比较高的地震数据，减小迭代硬阈值并同时提高迭代次数；对于信噪比较低的地震数据，减小迭代次数并同时增大迭代阈值。其中步骤(2)信噪比分别按照5:1、2:1和1:1进行划分。
所述步骤(3)最后对不同尺度反演结果进行多尺度融合是对不同尺度反演结果进行加权平均，以获得多尺度纵波阻抗稀疏反演结果。
所述多尺度地震反演方法是依据地震信号的时频谱将地震信号划分为不同频带对应的地震信号：
利用最小二乘将目标泛函表示为公式(6)，
其中，s为实际地震信号，sfi为多尺度地震信号，公式1具体物理含义为实际地震信号表示为所有尺度地震信号之和；公式6中，J(R)f为在不同尺度条件下的频率域目标泛函，用L2范数表示，其中α1为频率域反演调控参数，Yf为在不同尺度下，由频率域地震信号虚部和实部组成的列向量，G1f为在不同尺度下由频率域子波实部与虚部组成的核矩阵，R为反射系数序列；
所述快速匹配追踪的算法是通过在抄完被子波字典中寻找最佳匹配原子，在残差到达允许误差范围内将非平稳信号表示为：
匹配追踪公式：
匹配追踪字典为：
DN＝[g1,g2,…,gm] (8)
获取字典后进行振幅修正：
上述公式(7)-(9)中：s为实际地震信号，Rk代表第k次迭代所得到的残差，Rn为残差，g1-gm为不同时频子波，m为第m个子波，其具体由实际地震数据决定，T为矩阵转置符号，σ为gγk为第k次迭代的最佳匹配原子，DN为匹配追踪字典，aN为修正振幅，P为残差信号；
为实现稀疏反演，加入约束相，快速匹配追踪稀疏反演方法最终目标泛函表达为：
引入低频模型约束公式为
且：
其中，ψ为低频约束的阻抗模型参数，K为迭代次数，IPi与IP0分别为实际与相对阻抗，通过对反射系数进行道积分获得纵波阻抗信息：
其中，公式10为最终目标反演，Y为由频率域地震信号虚部和实部组成的列向量，G1为频率域子波实部与虚部组成的核矩阵，R为反射系数序列，ψ为低频约束的阻抗模型参数，Rit匹配追踪迭代次数，ε为收敛硬阈值，该部分物理含义当迭代次数小于设定的迭代次数即Rit＜K，或残差小于收敛阈值则迭代停止；公式12为纵波阻抗计算公式，其中Ip(t)为纵波阻抗，Ip(t0)为低频模型的第一个值R(τ)为反射系数，也即反演结果，利用反射系数反演结果进行道积分得到最终纵波阻抗反演结果。
所述公式(1)至公式(6)的演化过程为：
根据Robinson提出的自激自收褶积模型
D＝G*m+N (2)
和Margrave提出的联合域褶积模型
其中，D为地震信号，G为映射核矩阵，即为子波矩阵，m为反射系数序列，N为噪声；公式(3)为联合域褶积算子，其中s(f)为频率域地震信号，sf为对应不同尺度地震信号，ωf(f)为不同尺度的频率域地震子波，r(t)为反射系数序列，N(f)为随机噪声频率域表达式，
其中ωf为尺度子波，fi为部分频率，则联合域褶积模型公式改写为：Sf＝GfR+N，将该式改写为虚部与实部之和为：
对于公式(4)，ω为频率域子波矩阵(也称核矩阵)，ωf为在不同尺度下的频率域子波矩阵，ωf(2πf)为在不同尺度下的频率域子波；公式5中，Re[Sf]表示频率域地震信号的实部，Im[Sf]表示频率域地震信号的虚部，Re[Gf]表示频率域核矩阵的实部，Im[Gf]表示频率域核矩阵的虚部，R为反射稀疏序列，Re[N]表示频率域噪声的实部，Im[N]表示频率域噪声的虚部；
若G1f＝[Re[Gf] Im[Gf]]T，Yf＝[Re[Sf] Im[Sf]]T，Ν1＝[Re[N] Im[N]]T，则公式(5)改写为Yf＝G1fR+N1，此时利用最小二乘将目标泛函表示为：
通过对地震信号多尺度分解，大尺度地震信号对应大而厚的地层，即为低频信息，小尺度地震信号对应小而薄的地层，即为高频信息。
所述快速匹配追踪算法中，公式(11)为反演所加入的低频模型约束，也就是公式(10)中的||ψ-MR||改写为
其中H为时频原子字典，gi为最佳匹配原子，αi为振幅。
本发明是利用不同尺度(频带)地震数据进行纵波阻抗反演，通过不断更新总部阻抗低频模型进行反演，首先，利用大尺度地震信号进行纵波阻抗反演得到大尺度纵波阻抗反演结果，在此基础上，利用大尺度总部阻抗反演结果作为低频模型进行中尺度总部阻抗反演，依次类推利用中尺度总部阻抗反演结果作为低频模型进行小尺度总部阻抗反演，这种方法除了可以对曾为进行精确刻画，其另一优点在于对于初始模型精度依赖程度大大降低，即使初始模型精度不足抑或十分平滑，利用本专利方法依然可以进行有效的总部阻抗反演，并得到较好的反演结果，此外，本发明也缓解了反演多解性问题。
附图说明
图1地震信号多尺度分解：(a)实际地震信号，(b)大尺度地震信号，(c)中尺度地震信号，(d)小尺度地震信号；其中黑色箭头表示实际地震信号，黑色空心箭头表示尺度地震信号。
图2纵波阻抗反演结果：(a)大尺度纵波阻抗反演结果，(b)中尺度纵波阻抗反演结果，(c)小尺度纵波阻抗反演结果，(d)多尺度纵波阻抗融合结果。
图3抗噪性测试：(a)无噪声情况下多尺度纵波阻抗反演结果，(b)信噪比为5:1情况下多尺度纵波阻抗反演结果，(c)信噪比为2:1情况下多尺度纵波阻抗反演结果，(d)信噪比为1情况下多尺度纵波阻抗反演结果
图4实际资料测试：(a)中国东部SL勘探工区内的某条地震剖面，(b)常规反演方法纵波阻抗反演结果，(c)本发明反射系数反演结果，(d)本发明纵波阻抗反演结果。
具体实施方式
下面结合附图和实施例对本发明做进一步说明。
实施例1
一种基于快速匹配追踪的多尺度地震反演方法，首先(1)进行了地震信号多尺度分解的分析，对不同尺度地震信号如何划分进行分析，并同时对不同尺度地震信号进行分析(其中地震尺度是由地震信号进行频带分析并进行相应频带分解得到，不同地震数据划分结果不同)。(2)快速匹配追踪块化地层纵波阻抗反演方法，对不同信噪比(信噪比分别为5:1，2:1以及1:1)的情况下的迭代阈值与迭代次数设置进行分析用以获得最佳反演结果，其结论为：在信噪比较高的情况下，可以提高匹配次数以及减小迭代阈值以追求精确的反演结果，当信噪比较低时，应适当减小匹配次数增大阈值来缓解噪声影响，有效的突出了薄层边界。(3)多尺度快速匹配追踪纵波阻抗反演方法，在低频模型约束基础上进行大尺度快速匹配追踪反演以获得富含低频信息的纵波阻抗稀疏反演结果，在此基础上，利用大尺度纵波阻抗反演结果作为低频模型进行小尺度纵波阻抗稀疏反演，最后对不同尺度反演结果进行多尺度融合(对不同尺度反演结果进行加权平均)以获得多尺度纵波阻抗稀疏反演结果。大幅提高了反演结果的高低频信息，同时有效的提升了反演方法的分辨率。
地震信号多尺度分解，由于地下介质岩性不同会导致地震数据表现出不同的频率成分，因此利用多尺度地震数据可以对地震层序进行有效刻画，对于多尺度分解后的地震信号进行分析可知，不同尺度地震数据对应不同的频率成分，其具体为大尺度地震信号为地震信号低频信息，小尺度地震信号为地震信号高频信息，对于不同尺度地震信号，其频率变化不大，较为稳定且小尺度镶嵌于大尺度中。
所述快速匹配追踪块化地层纵波阻抗反演方法的参数设置，对于信噪比较高的地震数据，应减小迭代硬阈值并同时提高迭代次数。对于信噪比较低的地震数据，应减小迭代次数并同时增大迭代阈值，以避免对噪声进行匹配，其具体数值应根据地震数据进行调整。
(1)多尺度地震反演方法
本实施例首先介绍了时频域反演问题，然后介绍了快速匹配追踪算法，最终由公式(10)给出最终目标泛函。
依据地震信号的时频谱将地震信号划分为不同频带对应的地震信号：
其中，sfi为多尺度地震信号，若G1f＝[Re[Gf] Im[Gf]]T，Yf＝[Re[Sf] Im[Sf]]T，Ν1＝[Re[N] Im[N]]T则式5可改写为Yf＝G1fR+N1，此时利用最小二乘可将目标泛函表示为公式6，1954年经典自激自收褶积模型被Robinson提出：
D＝G*m+N (2)
其中D为地震信号，G为映射核矩阵，即为子波矩阵，m为反射系数序列，N为噪声。
1998年Margrave等人提出了联合域褶积模型：
其中，ωf(f)为频率域尺度子波矩阵，r(t)为时间域反射系数矩阵，N(f)为随机噪声频率域表达式。且：
其中ωf为尺度子波，fi为部分频率，则褶积公式可改写为：Sf＝GfR+N，将该式改写为虚部与实部之和为：
若G1f＝[Re[Gf] Im[Gf]]T，Yf＝[Re[Sf] Im[Sf]]T，Ν1＝[Re[N] Im[N]]T则式5可改写为Yf＝G1fR+N1，此时利用最小二乘可将目标泛函表示为：
对于地震信号多尺度分解，如图1所示为地震信号的多尺度分解。图1(a)为实际地震信号。图1(b)为实际地震信号及其对应的大尺度地震信号。图1(c)为实际地震信号及其对应的中尺度地震信号。图1(d)为实际地震信号及其对应的小尺度地震信号。由图可知，大尺度地震信号对应大而厚的地层，即为低频信息。小尺度地震信号对应小而薄的地层，即为高频信息。
(2)快速匹配追踪算法
匹配追踪算法通过在抄完被子波字典中寻找最佳匹配原子，在残差到达允许误差范围内可将非平稳信号表示为：
式7为匹配追踪公式，其中s为地震信号。gγk为第k次迭代的最佳匹配原子，Rn为残差。公式8为匹配追踪字典的表示方式，公式9为进行的振幅修正，公式(11)为反演所加入的低频模型约束，也就是公式(10)中的||ψ-MR||，公式(2)还可以改写为
其中H为时频原子字典gi为最佳匹配原子，αi为振幅。为增加匹配追踪计算效率，本发明发展了快速匹配追踪算法，通过同时对多个原子进行匹配搜索来减少计算时间。其中字典为：
DN＝[g1,g2,…,gm] (8)
获取字典后还应进行振幅修正：
其中aN为修正振幅，P为残差信号。为实现稀疏反演，需要加入约束相，快速匹配追踪稀疏反演方法最终目标泛函可表达为：
为了缓解超定问题给反演方法带来的多解性，本文引入低频模型约束来缓解该问题，其中低频模型约束公式为：
其中ψ为低频约束的阻抗模型参数，K为迭代次数，且：
其中IPi与IP0为实际与相对阻抗，通过对反射系数进行道积分可获得纵波阻抗信息：
实施例2
结合附图，一种基于快速匹配追踪的多尺度地震反演方法，首先，图1(a)所示为原始地震数据，图1(b)-图1(b)分别为大尺度地震信号，中尺度地震信号以及小尺度地震信号，其中黑色箭头表示原始地震信号，白色箭头为多尺度地震信号。由图1可知，大尺度地震信号对应地震信号的低频信息，小尺度地震信号对应地震信号的高频信息，本发明在此基础上进行逐级迭代反演，可以有效缓解由于地震信号高频下降以及低频缺失所带来的反演多解性问题。
图2所示为纵波阻抗反演结果，图2(a)为大尺度纵波阻抗反演结果，由图可知其代表了纵波阻抗低频信息。图2(b)为中尺度纵波阻抗反演结果，其在大尺度纵波阻抗反演结果的基础上进行中尺度反演，从图中可以看出，中尺度纵波阻抗反演结果在大尺度纵波阻抗反演的结果上增加了部分高频信息。图2(c)为小尺度纵波阻抗反演结果，该反演结果较好的反应了纵波阻抗的高频信息。图2(d)为多尺度纵波阻抗融合结果，由图可知反演结果与原始纵波阻抗具有较高的相似性，说明了该方法的有效性。值得说明的是，本发明是利用了迭代反演的算法，其具体为利用大尺度纵波阻抗反演结果作为中尺度纵波阻抗反演的初始模型，在此基础上利用中尺度地震信号的反演结果作为小尺度纵波阻抗反演的初始模型。
图3为抗噪性测试，其中图3(a)-(d)分别为在无噪声、信噪比为5:1、信噪比为2:1、信噪比为1的情况下进行的多尺度纵波阻抗反演结果，由图可知，该方法具有较强的抗噪能力，进一步说明了该方法的实用性。
图4为实际数据测试结果，其中(a)为中国东部SL勘探工区内的某条地震剖面，(b)为常规反演方法纵波阻抗反演结果，(c)本发明反射系数反演结果，(d)本发明纵波阻抗反演结果。由实际数据测试结果可知，本发明方法相较于常规反演方法分辨率有显著的提升，与此同时可以在连续变化的地层进行精确有效的纵波阻抗刻画，这也进一步说明了该方法的可行性及有效性。
根据发明内容主要流程如下：
(1)地震信号多尺度分解
依据地震信号的时频谱可将地震信号划分为不同频带对应的地震信号：
其中，sf为多尺度地震信号。1954年经典自激自收褶积模型被Robinson提出：
D＝G*m+N (2)
其中D为地震信号，G为映射核矩阵，即为子波矩阵，m为反射系数序列，N为噪声。
1998年Margrave等人提出了联合域褶积模型：
其中，ωf(f)为频率域尺度子波矩阵，r(t)为时间域反射系数矩阵，N(f)为随机噪声频率域表达式。且：
其中ωf为尺度子波，fi为部分频率，则褶积公式可改写为：Sf＝GfR+N，将该式改写为虚部与实部之和为：
若G1f＝[Re[Gf] Im[Gf]]T，Yf＝[Re[Sf] Im[Sf]]T，Ν1＝[Re[N] Im[N]]T则式5可改写为Yf＝G1fR+N1，此时利用最小二乘可将目标泛函表示为：
对于地震信号多尺度分解，如图1所示为地震信号的多尺度分解。图1(a)为实际地震信号。图1(b)为实际地震信号及其对应的大尺度地震信号。图1(c)为实际地震信号及其对应的中尺度地震信号。图1(d)为实际地震信号及其对应的小尺度地震信号。由图可知，大尺度地震信号对应大而厚的地层，即为低频信息。小尺度地震信号对应小而薄的地层，即为高频信息。
(2)快速匹配追踪算法
匹配追踪算法通过在抄完被子波字典中寻找最佳匹配原子，在残差到达允许误差范围内可将非平稳信号表示为：
其中s为非平稳信号。gγk为第k次迭代的最佳匹配原子，Rn为残差。公式(2)还可以改写为
其中H为时频原子字典gi为最佳匹配原子，αi为振幅。为增加匹配追踪计算效率，本发明发展了快速匹配追踪算法，通过同时对多个原子进行匹配搜索来减少计算时间。其中字典为：
DN＝[g1,g2,…,gm] (8)
获取字典后还应进行振幅修正：
aN＝[(DN)T(DN)+σ2E]-1(DN)TP (9)
其中aN为修正振幅，P为残差信号。为实现稀疏反演，需要加入约束相，快速匹配追踪稀疏反演方法目标泛函可表达为：
其中ψ为低频约束的阻抗模型参数，K为迭代次数，且：
其中IPi与IP0为实际与相对阻抗，通过对反射系数进行道积分可获得纵波阻抗信息：
本发明充分考虑地震信号中蕴含的“频带”特征，将地震信号进行多尺度分解，并采取逐级迭代的思想进行反演，首先将地震信号进行时频分析获得有效频带范围后进行多尺度分解，其中大尺度信号对应低频信息即厚地层，小尺度信号对应高频信息即薄层，大尺度地震反演可以获得波阻抗的低频信息，在此基础上，利用大尺度纵波阻抗反演结果作为初始模型进行小尺度纵波阻抗反演，使得纵波阻抗反演结果兼具高低频信息。除此之外，为提升薄层边界识别能力，本发明同时考虑通过控制匹配次数以及迭代硬阈值来实现快速匹配追踪块化地层反演。
Claims (8)
1.一种基于快速匹配追踪的多尺度地震反演方法，其特征在于包括：
(1)进行地震信号多尺度分解，包括对不同尺度地震信号划分，对划分后的不同尺度地震信号进行分析；
(2)建立快速匹配追踪块化地层纵波阻抗反演方法，包括对不同信噪比情况下的迭代阈值与迭代次数设置；
(3)建立多尺度快速匹配追踪纵波阻抗反演方法，包括在低频模型约束基础上进行大尺度快速匹配追踪反演以获得富含低频信息的纵波阻抗稀疏反演结果，在此基础上，利用大尺度纵波阻抗反演结果作为低频模型进行小尺度纵波阻抗稀疏反演，最后对不同尺度反演结果进行多尺度融合。
2.根据权利要求1所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述步骤(1)地震信号多尺度分解是：根据地下介质岩性不同导致地震数据表现出不同的频率成分，利用多尺度地震数据对地震层序进行有效刻画，对于多尺度分解后的地震信号进行分析可知，大尺度地震信号为地震信号低频信息，小尺度地震信号为地震信号高频信息，且小尺度地震信号镶嵌于大尺度地震信号中。
3.根据权利要求1所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述步骤(2)快速匹配追踪块化地层纵波阻抗反演方法设置包括：对于信噪比较高的地震数据，减小迭代硬阈值并同时提高迭代次数；对于信噪比较低的地震数据，减小迭代次数并同时增大迭代阈值。
4.根据权利要求3所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述步骤(2)信噪比分别按照5:1、2:1和1:1进行划分。
5.根据权利要求4所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述步骤(3)最后对不同尺度反演结果进行多尺度融合是对不同尺度反演结果进行加权平均，以获得多尺度纵波阻抗稀疏反演结果。
6.根据权利要求1-5任一所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于，
所述多尺度地震反演方法是依据地震信号的时频谱将地震信号划分为不同频带对应的地震信号：
利用最小二乘将目标泛函表示为公式(6)，
其中，s为实际地震信号，
为多尺度地震信号，公式1具体物理含义为实际地震信号表示为所有尺度地震信号之和；公式6中，J(R)f为在不同尺度条件下的频率域目标泛函，用L2范数表示，其中α1为频率域反演调控参数，Yf为在不同尺度下，由频率域地震信号虚部和实部组成的列向量，G1f为在不同尺度下由频率域子波实部与虚部组成的核矩阵，R为反射系数序列；
所述快速匹配追踪的算法是通过在抄完被子波字典中寻找最佳匹配原子，在残差到达允许误差范围内将非平稳信号表示为：
匹配追踪公式：
匹配追踪字典为：
DN＝[g1,g2,…,gm] (8)
获取字典后进行振幅修正：
aN＝[(DN)T(DN)+σ2E]-1(DN)TP (9)
上述公式(7)-(9)中：s为实际地震信号，Rk代表第k次迭代所得到的残差，Rn为残差，g1-gm为不同时频子波，m为第m个子波，其具体由实际地震数据决定，T为矩阵转置符号，σ为gγk为第k次迭代的最佳匹配原子，DN为匹配追踪字典，aN为修正振幅，P为残差信号；
为实现稀疏反演，加入约束相，快速匹配追踪稀疏反演方法最终目标泛函表达为：
引入低频模型约束公式为
且：
其中，ψ为低频约束的阻抗模型参数，K为迭代次数，IPi与IP0分别为实际与相对阻抗，通过对反射系数进行道积分获得纵波阻抗信息：
其中，公式10为最终目标反演，Y为由频率域地震信号虚部和实部组成的列向量，G1为频率域子波实部与虚部组成的核矩阵，R为反射系数序列，ψ为低频约束的阻抗模型参数，Rit匹配追踪迭代次数，ε为收敛硬阈值，该部分物理含义当迭代次数小于设定的迭代次数即Rit＜K，或残差小于收敛阈值则迭代停止；公式12为纵波阻抗计算公式，其中Ip(t)为纵波阻抗，Ip(t0)为低频模型的第一个值R(τ)为反射系数，也即反演结果，利用反射系数反演结果进行道积分得到最终纵波阻抗反演结果。
7.根据权利要求6所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述公式(1)至公式(6)的演化过程为：
根据Robinson提出的自激自收褶积模型
D＝G*m+N (2)
和Margrave提出的联合域褶积模型
其中，D为地震信号，G为映射核矩阵，即为子波矩阵，m为反射系数序列，N为噪声；公式(3)为联合域褶积算子，其中s(f)为频率域地震信号，sf为对应不同尺度地震信号，ωf(f)为不同尺度的频率域地震子波，r(t)为反射系数序列，N(f)为随机噪声频率域表达式，
其中ωf为尺度子波，fi为部分频率，则联合域褶积模型公式改写为：Sf＝GfR+N，将该式改写为虚部与实部之和为：
对于公式(4)，ω为频率域子波矩阵(也称核矩阵)，ωf为在不同尺度下的频率域子波矩阵，ωf(2πf)为在不同尺度下的频率域子波；公式5中，Re[Sf]表示频率域地震信号的实部，Im[Sf]表示频率域地震信号的虚部，Re[Gf]表示频率域核矩阵的实部，Im[Gf]表示频率域核矩阵的虚部，R为反射稀疏序列，Re[N]表示频率域噪声的实部，Im[N]表示频率域噪声的虚部；
若G1f＝[Re[Gf] Im[Gf]]T，Yf＝[Re[Sf] Im[Sf]]T，Ν1＝[Re[N] Im[N]]T，则公式(5)改写为Yf＝G1fR+N1，此时利用最小二乘将目标泛函表示为：
通过对地震信号多尺度分解，大尺度地震信号对应大而厚的地层，即为低频信息，小尺度地震信号对应小而薄的地层，即为高频信息。
8.根据权利要求7所述的基于快速匹配追踪的多尺度地震反演方法，其特征在于所述快速匹配追踪算法中，公式(11)为反演所加入的低频模型约束，也就是公式(10)中的||ψ-MR||改写为
其中H为时频原子字典，gi为最佳匹配原子，αi为振幅。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010308485.6A CN113534250A (zh) | 2020-04-18 | 2020-04-18 | 一种基于快速匹配追踪的多尺度地震反演方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010308485.6A CN113534250A (zh) | 2020-04-18 | 2020-04-18 | 一种基于快速匹配追踪的多尺度地震反演方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN113534250A true CN113534250A (zh) | 2021-10-22 |
Family
ID=78123505
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010308485.6A Pending CN113534250A (zh) | 2020-04-18 | 2020-04-18 | 一种基于快速匹配追踪的多尺度地震反演方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113534250A (zh) |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105388527A (zh) * | 2015-11-30 | 2016-03-09 | 中国石油大学(北京) | 一种基于复数域匹配追踪算法的油气检测方法 |
| CN106291677A (zh) * | 2015-05-22 | 2017-01-04 | 中国石油化工股份有限公司 | 一种基于匹配追踪方法的叠后声波阻抗反演方法 |
| CN106842298A (zh) * | 2015-12-04 | 2017-06-13 | 中国石油化工股份有限公司 | 一种基于匹配追踪的不整合强反射自适应分离方法 |
| CN109143328A (zh) * | 2017-06-19 | 2019-01-04 | 中国石油化工股份有限公司 | 一种叠后地震反演方法 |
| US20190113641A1 (en) * | 2017-10-12 | 2019-04-18 | Southern University Of Science And Technology | Seismic travel time tomographic inversion method based on two point ray tracing |
| CN110516358A (zh) * | 2019-08-28 | 2019-11-29 | 电子科技大学 | 一种基于稀疏表示的地震各向异性参数反演方法 |
| CN110687597A (zh) * | 2019-10-22 | 2020-01-14 | 电子科技大学 | 一种基于联合字典的波阻抗反演方法 |
-
2020
- 2020-04-18 CN CN202010308485.6A patent/CN113534250A/zh active Pending
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106291677A (zh) * | 2015-05-22 | 2017-01-04 | 中国石油化工股份有限公司 | 一种基于匹配追踪方法的叠后声波阻抗反演方法 |
| CN105388527A (zh) * | 2015-11-30 | 2016-03-09 | 中国石油大学(北京) | 一种基于复数域匹配追踪算法的油气检测方法 |
| CN106842298A (zh) * | 2015-12-04 | 2017-06-13 | 中国石油化工股份有限公司 | 一种基于匹配追踪的不整合强反射自适应分离方法 |
| CN109143328A (zh) * | 2017-06-19 | 2019-01-04 | 中国石油化工股份有限公司 | 一种叠后地震反演方法 |
| US20190113641A1 (en) * | 2017-10-12 | 2019-04-18 | Southern University Of Science And Technology | Seismic travel time tomographic inversion method based on two point ray tracing |
| CN110516358A (zh) * | 2019-08-28 | 2019-11-29 | 电子科技大学 | 一种基于稀疏表示的地震各向异性参数反演方法 |
| CN110687597A (zh) * | 2019-10-22 | 2020-01-14 | 电子科技大学 | 一种基于联合字典的波阻抗反演方法 |
Non-Patent Citations (7)
| Title |
|---|
| SONG PEI ET AL: "《Fast matching pursuit based multi-scale seismic inversion》", SEG TECHNICAL PROGRAM EXPANDED ABSTRACTS 2019, pages 1881 - 1821 * |
| 刘晓晶等: "《基于正交匹配追踪算法的叠前地震反演方法》", 石油地球物理勘探, vol. 50, no. 5, 31 October 2015 (2015-10-31), pages 925 - 935 * |
| 印兴耀等: "《多尺度快速匹配追踪多域联合地震反演方法》", 地球物理学报, vol. 63, no. 9, pages 3431 - 3441 * |
| 李坤: "《非平稳地震频率域稀疏反演方法研究》", 中国优秀硕士学位论文全文数据库 基础科学辑, no. 2018, 15 June 2018 (2018-06-15), pages 21 - 22 * |
| 李坤等: "《基于匹配追踪谱分解的时频域FAVO流体识别方法》", 石油学报, vol. 37, no. 6, pages 777 - 786 * |
| 李坤等: "《基于快速匹配追踪的混合域地震稀疏反演方法》", 中国石油大学学报( 自然科学版), vol. 42, no. 1, 28 February 2018 (2018-02-28), pages 50 - 59 * |
| 李坤等: "基于快速匹配追踪的混合域地震稀疏反演方法", 中国石油大学学报(自然科学版), vol. 42, no. 1, pages 50 - 59 * |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104849756B (zh) | 一种提高地震数据分辨率增强有效弱信号能量的方法 | |
| CN104237945B (zh) | 一种地震资料自适应高分辨处理方法 | |
| CN109541685B (zh) | 一种河道砂体识别方法 | |
| CN110895348B (zh) | 一种地震弹性阻抗低频信息提取方法、系统及存储介质 | |
| CN111142156B (zh) | 基于相控的地震强屏蔽时频信息提取及剥离方法 | |
| CN104730576A (zh) | 基于Curvelet变换的地震信号去噪方法 | |
| CN113077386A (zh) | 基于字典学习和稀疏表征的地震资料高分辨率处理方法 | |
| CN108226997A (zh) | 一种基于叠前地震数据的地震相划分方法 | |
| CN116520419A (zh) | 一种热流体裂缝通道识别方法 | |
| CN109655883A (zh) | 一种针对目标的地震分频方法及系统 | |
| CN111090117B (zh) | 一种相控正演约束下的有效储层预测方法及系统 | |
| CN113534250A (zh) | 一种基于快速匹配追踪的多尺度地震反演方法 | |
| CN110673211B (zh) | 一种基于测井与地震数据的品质因子建模方法 | |
| CN110865409B (zh) | 一种基于波阻抗低秩正则化的地震波阻抗反演方法 | |
| CN106291675A (zh) | 一种基于基追踪技术的地震数据重构方法 | |
| CN105891882A (zh) | 一种基于裂缝时频表征的匹配追踪分频方法 | |
| CN113589365B (zh) | 基于时频域信息的储层尖灭线描述方法 | |
| Liu et al. | An improved Gaussian frequency domain sparse inversion method based on compressed sensing | |
| CN114859404A (zh) | 超采样地震波形匹配方法及装置 | |
| CN115561814A (zh) | 断陷盆地缓坡带中浅层地层-岩性油藏的勘探方法 | |
| CN111239816B (zh) | 一种基于匹配追踪的超低信噪比高精度速度谱生成方法 | |
| CN113419275B (zh) | 一种基于稀疏字典学习的高分辨率地震处理方法 | |
| US11960043B2 (en) | Data driven method to invert for the formation anisotropic constants using borehole sonic data | |
| CN113009579B (zh) | 地震数据反演方法及装置 | |
| CN114428283B (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 |