CN107272062A - 一种数据驱动的地下介质q场估计方法 - Google Patents
一种数据驱动的地下介质q场估计方法 Download PDFInfo
- Publication number
- CN107272062A CN107272062A CN201710543599.7A CN201710543599A CN107272062A CN 107272062 A CN107272062 A CN 107272062A CN 201710543599 A CN201710543599 A CN 201710543599A CN 107272062 A CN107272062 A CN 107272062A
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- window
- msubsup
- molecule
- 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
- 238000000034 method Methods 0.000 title claims abstract description 57
- 238000001228 spectrum Methods 0.000 claims abstract description 34
- 238000013507 mapping Methods 0.000 claims abstract description 16
- 230000003044 adaptive effect Effects 0.000 claims abstract description 15
- 230000015572 biosynthetic process Effects 0.000 claims abstract description 15
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 6
- 238000001914 filtration Methods 0.000 claims description 14
- 238000004458 analytical method Methods 0.000 claims description 6
- 239000012634 fragment Substances 0.000 claims description 6
- 230000011218 segmentation Effects 0.000 claims description 5
- 238000006243 chemical reaction Methods 0.000 claims description 4
- 230000006835 compression Effects 0.000 claims description 4
- 238000007906 compression Methods 0.000 claims description 4
- 238000012545 processing Methods 0.000 claims description 4
- 239000000284 extract Substances 0.000 claims description 3
- 238000000605 extraction Methods 0.000 claims description 2
- 230000008569 process Effects 0.000 claims description 2
- 230000014509 gene expression Effects 0.000 claims 1
- 230000008859 change Effects 0.000 abstract description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000012216 screening Methods 0.000 description 3
- 238000013459 approach Methods 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 239000004215 Carbon black (E152) Substances 0.000 description 1
- 230000002159 abnormal effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 229930195733 hydrocarbon Natural products 0.000 description 1
- 150000002430 hydrocarbons Chemical class 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 230000011514 reflex Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000011664 signaling Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 230000006641 stabilisation Effects 0.000 description 1
- 238000011105 stabilization Methods 0.000 description 1
- 238000004454 trace mineral analysis 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
- 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
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/624—Reservoir parameters
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (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
一种数据驱动的地下介质Q场估计方法,采用地震包络局部峰值作为地层结构约束来构造自适应于地层结构的分子时窗,然后采用非线性压缩映射子波振幅谱估计方法从地震道分子分解时频谱中估计正比于时变子波振幅谱的分量,接着计算质心频率,并对质心频率做筛选,减少子波干涉和噪声的影响,最后用质心频率偏移法估计得到稳健的地下介质Q场。所述分子时窗的构造方法,可以保证每个时窗内至少有一个反射子波,在减少时窗数量、提高计算效率的同时,使得这些时窗在横向上与地层结构有关,有利于保持所估计Q场的横向连续性。本发明能够有效克服复杂地质结构变化和噪声干扰的影响,得到稳定的地下介质Q场。
Description
技术领域
本发明属于地球物理勘探领域,涉及一种地下介质物性参数估计方法,特别涉及一种数据驱动的地下介质Q场估计方法。
背景技术
地震波在地下介质中传播时会随着传播时间衰减,主要表现为振幅衰减、相位畸变、主频降低,且高频部分比低频部分衰减快,浅层比深层衰减快。引起衰减的因素包括几何扩散、散射、介质的非完全弹性等。其中由介质的非完全弹性引起的衰减常用品质因子Q来描述。如果知道地下介质Q值,可对地震记录进行有效的反Q滤波补偿,从而使地震记录的高频成分得到增强,地下浅、中、深层反射地震记录的纵向分辨率都得到提高;另外,Q值也与岩性、饱和度、孔隙率等有关,可用于储层识别和烃类检测。因此,估计地下介质Q值具有重要的实际意义。
迄今,国内外学者提出了多种直接估计Q值的方法。时间域方法有子波模拟法、解析信号法和脉冲上升时间法等,这些方法都需要振幅保真度高的数据,而实际记录常受到几何扩散、散射等因素的影响,致使这些时间域方法精度降低。频率域常用的方法有对数谱比值法(LSR,log-spectrum ratio)、谱匹配法、中心频率偏移法等,这些方法需要用时窗截取地震记录,再对截出的地震记录进行分析。这类方法时窗的选择对估计结果影响较大,如何最佳地选择窗的类型及长度是一个难题。利用脉冲瞬时属性的变化来估计Q值也是一条重要的途径。Mathneey和Nowack提出了瞬时频率匹配法,高静怀和杨森林等人提出利用小波域地震子波包络峰值处瞬时频率估计Q值的方法WEPIF(wavelet envelope peakinstantaneous frequency)。这类方法的优点是时窗易于选取或者不需要加时窗,估计得到的Q值分辨率较高,但由于只利用了脉冲包络峰值处单个或相邻几个数据值,估计得到的Q值受随机噪声的影响很大。Tonn对近十种直接估计Q值的方法进行了比较,结果表明每种方法都有一定的适用条件。
由于垂直地震剖面(VSP)资料可对穿过介质的波进行直接测量,这为上述直接估计介质Q值的方法以及利用波形反演估计介质Q值及其它参数的方法(如:Stewart提出的同时利用相邻四道记录的全波形在频率域反演速度、衰减因子、上行波和下行波等的方法;Amundsen和Mittet只利用初至下行波和一次反射波在频率域反演相速度和品质因子Q的方法)提供了途径。VSP是在地震测井基础上发展起来的,而测井数量有限,故由VSP资料估计得到的Q值不能反映整个地下介质的Q场分布。为此,学者们也尝试研究利用反射地震资料估计地下介质Q场的方法。目前常用的各种Q值估计(或反演)方法在用于反射地震资料时往往会碰到许多问题,如:复杂的地质结构使得反射波互相交叠,不能直接从地震记录获得时变子波的波形、振幅谱以及瞬时频率等。另外,噪声干扰影响瞬时属性的提取,进而影响Q值估计的准确性。也就是说,复杂的地质结构和噪声干扰常常使得估计出的地下介质的Q值不稳定。
发明内容
本发明的目的在于克服上述现有技术的缺点,提供一种数据驱动的地下介质Q场估计方法,该方法为一种数据驱动的基于自适应分子分解利用非线性压缩映射提取时变子波振幅谱并用质心频率偏移法估计地下介质Q场的方法,该方法采用地震包络局部峰值作为地层结构约束来构造自适应于地层结构的分子时窗,然后采用非线性压缩映射子波振幅谱估计方法从地震道分子分解时频谱中估计正比于时变子波振幅谱的分量,接着计算质心频率,并对质心频率做筛选,减少子波干涉和噪声的影响,最后用质心频率偏移法估计得到稳健的地下介质Q场。
为实现上述目的,本发明的目的是通过以下技术方案来解决的:
一种数据驱动的地下介质Q场估计方法,采用地震包络局部峰值作为地层结构约束来构造自适应于地层结构的分子时窗,然后采用非线性压缩映射子波振幅谱估计方法从地震道分子分解时频谱中估计正比于时变子波振幅谱的分量,接着计算质心频率,并对质心频率做筛选,减少子波干涉和噪声的影响,最后用质心频率偏移法估计得到稳健的地下介质Q场。
本发明进一步的改进在于,具体包括以下步骤:
1)提取地震道包络局部峰值:
2)生成满足单位分割的原子窗集:
选择基本原子窗函数G(t),令Gj(t)=G(t-jΔt)表示中心位于第j个采样点的原子窗,对原子窗族{Gj:1≤j≤N}按下式归一化
得一组单位分割原子窗集{gj:1≤j≤N},这里N为地震道采样点的个数;
3)构造自适应分子时窗:
选择每相邻两个包络局部峰值点的中点作为各分子窗的边界点,将边界点间的满足单位分割的原子窗叠加起来就形成了分子窗;设第k个分子窗对应的包络局部峰值点位于Pk,前一个包络局部峰值点位于Pk-1(P0=-P1),后一个包络局部峰值点位于Pk+1,则此分子窗ψk(t),对应的第一个原子窗的中心位于Mk-1+1=(Pk-1+Pk)/2+1(M0=0),对应的最后一个原子窗的中心位于Mk=(Pk+Pk+1)/2,分子窗由这中间的Mk-Mk-1个原子窗叠加得到,即分子窗ψk(t)由下式表示
令L为分子窗的个数,分子窗族{ψk(t):1≤k≤L}也构成单位分割;
4)对步骤3)得到的自适应分子时窗进行能量归一化:
令Ek表示第k个分子窗的能量,即
能量归一化以后的分子分析时窗为{ψk(t)/Ek:1≤k≤L};
对分子分析时窗进行平移和调制后,得到一组分子标架;此时,分子分解时频变换定义为
其中f为频率;
5)利用非线性压缩映射提取时变子波振幅谱:
设x0∈(a,b),取δ0>0,使得对于任意给定的δ>0,定义L2[a,b]空间中的子集
在[a,x0-δ0]上单增,且在[x0+δ0,b]上单减}
如果且q≤2,则有任取定义函数
此处0<q≤1;显然,0≤Fq(u;x)≤1;设cq>0,α,β>1,对于定义上的算子P为
P是到的非线性压缩映射算子;
对于第k个分子窗中的地震记录片段,其振幅谱为设定迭代初值为:
其中0<q≤1,
其中fM为振幅谱的截止频率;建立迭代uk=P(uk-1),其中算子参数由最小二乘得到;由于算子是压缩映射,得到不动点u*,则估计的子波振幅谱为记为Lk(f);
6)计算质心频率:
对于第k个分子窗中的地震记录片段,质心频率fc,k为
式中Fc为子波振幅谱的截止频率;对计算得到的质心频率处理后得到最终的质心频率
7)估算地下介质Q场:
用质心频率偏移法估算Q值,相应的估算公式如下
式中:分别为t1时刻频谱的质心频率和方差;是t2时刻频谱的质心频率;时差Δt=t2-t1;为t1和t2之间地层介质的品质因子;进一步推导,得到t时刻地下介质的等效Q值
式中fc,0和为初始时刻地震子波振幅谱的质心频率和方差,即
用步骤6)中计算得到的质心频率代替fc(t),得到各个窗中心点处的Q值
Tk表示第k个分子窗的中心点所对应的时间位置;对窗中心点处估得的Q值插值后即得到稳定的地下介质Q场。
本发明进一步的改进在于,步骤1)中提取地震道包络局部峰值具体过程为:设s*(t)为地震道s(t)的Hilbert变换,则
a(t)=[s(t)2+s*(t)2]1/2
为地震道s(t)的包络;由上式计算地震道的包络,并提取包络局部峰值点。
本发明进一步的改进在于,步骤6)中对计算得到的质心频率做多次光滑滤波,并计算滤波前后的绝对误差,去掉绝对误差较大的一些质心频率后,再次拟合滤波,得到最终的质心频率
与现有技术相比,本发明具有以下有益效果:
本发明利用地震包络峰值作为地层结构约束来构造自适应于地层结构的分子时窗,由于地震数据的包络峰值可以大致反映地层的层序结构,因此构造的时窗在横向上与地层结构有关,有利于保持估计的Q场的横向连续性,并且这种自适应分子时窗构造方法可以保证每个分子窗内至少有一个反射子波,在减少时窗数量、提高计算效率的同时,可以有效减少窗端点对子波的截断;采用非线性压缩映射子波振幅谱估计方法可以从较短数据中估计子波振幅谱,这为本发明从地震道分子分解时频谱中估计正比于时变子波振幅谱的分量,从而得到稳定的质心频率提供了保证;在计算Q值之前对质心频率做了滤波筛选,可以有效减少子波干涉和随机噪声的影响,得到稳定的地下介质Q场。
附图说明
图1是合成地震道时反射系数示意图;
图2是合成地震道的包络局部峰值示意图;其中,实线表示非平稳合成地震道,虚线表示地震道包络,星号表示包络局部峰值;
图3是构成单位分割的原子Gabor时窗示例图;
图4是地震道包络局部峰值约束下寻找分子Gabor时窗边界点示例图;
图5是将边界点间的原子Gabor时窗叠加形成分子Gabor时窗示例图;
图6是能量归一化后得到的分子Gabor时窗示例图;
图7是一段实际地震资料剖面图;
图8是从地震资料剖面中提取出的质心频率图;
图9是对质心频率做滤波筛选后得到的质心频率图;
图10是基于滤波筛选后的质心频率计算得到的地下介质Q场图。
具体实施方式
下面结合附图对本发明做进一步详细描述:
本发明针对复杂的地质结构和噪声干扰常常使得估计出的地下介质的Q值不稳定的问题,以地震道包络局部峰值作为地层结构,提出基于自适应分子分解利用非线性压缩映射提取时变子波振幅谱并用质心频率偏移法估计得到稳健的地下介质Q场的方法。
本发明的物质基础是叠后反射地震数据,即在地表设置震源,在地表安置检波器接收地震波,经过一系列处理,最终叠加或偏移成像得到的地震数据。本发明的具体步骤为:
1)提取地震道包络局部峰值;
根据复道分析原理,设s*(t)为地震道s(t)的Hilbert变换,则
a(t)=[s(t)2+s*(t)2]1/2
为地震道s(t)的包络。对比图1中合成地震道所用的反射系数和图2中合成地震道包络局部峰值点可见,地震道的包络局部峰值与反射界面有一定的对应关系。
2)生成满足单位分割的原子窗集;
恰当地选择基本原子窗函数G(t),令Gj(t)=G(t-jΔt)表示中心位于第j个采样点的原子窗,对原子窗族{Gj:1≤j≤N}按下式归一化
可得一组单位分割原子窗集{gj:1≤j≤N},这里N为地震道采样点的个数。如图3中的一组窗就是满足单位分割的原子窗集。
3)构造自适应分子时窗;
选择每相邻两个包络局部峰值点的中点作为各分子窗的边界点,将边界点间的满足单位分割的原子窗叠加起来就形成了分子窗。如图4所示,设第k个分子窗对应的包络局部峰值点位于Pk,前一个包络局部峰值点位于Pk-1(P0=-P1),后一个包络局部峰值点位于Pk+1,则此分子窗ψk(t),对应的第一个原子窗的中心位于Mk-1+1=(Pk-1+Pk)/2+1(M0=0),对应的最后一个原子窗的中心位于Mk=(Pk+Pk+1)/2,分子窗由这中间的Mk-Mk-1个原子窗叠加得到,即分子窗ψk(t)可由下式表示
令L为分子窗的个数,易知分子窗族{ψk(t):1≤k≤L}也构成单位分割,如图5所示。这样,可以保证每个分子窗内至少有一个反射子波,在减少时窗数量的同时,可以有效减少窗端点对子波的截断效应。并且由于包络峰值可以大致反映地层的层序结构,构成的分子窗在横向上与地层结构有关,有利于保持处理后资料的横向连续性。
4)对步骤3)得到的自适应分子时窗进行能量归一化;
分子分解时频变换是保能量的变换,由图5可以看出,在步骤3)得到的分子窗族中,不同的窗能量不同。如果用这些窗直接构造分子Gabor标架,则一个地震道的时频能量表示不仅与该地震道有关,还与分子窗有关。希望这一时频能量表示只与地震道有关,为此,需要对每个分子窗进行能量归一化。令Ek表示第k个分子窗的能量,即
能量归一化以后的分子分析时窗为{ψk(t)/Ek:1≤k≤L},如图6所示。
对分子分析时窗进行平移和调制后,可得到一组分子标架。此时,分子分解时频变换可定义为
其中f为频率。
5)利用非线性压缩映射提取时变子波振幅谱;
设x0∈(a,b),取δ0>0足够小,使得对于任意给定的δ>0,定义L2[a,b]空间中的子集
在[a,x0-δ0]上单增,且在[x0+δ0,b]上单减}
如果且q≤2,则有任取定义函数
此处0<q≤1。显然,0≤Fq(u;x)≤1。设cq>0,α,β>1,对于定义上的算子P为
显然,P[u;x]在区间[a,b]上是连续的,并且δ≤P[u;x]≤cq+δ,由此可得P[u;x]∈L2[a,b]。因此,P是到L2[a,b]的非线性算子。根据Banach空间算子不动点定理,可以证明:P是到的压缩算子,并且在中存在唯一不动点。
对于第k个分子窗中的地震记录片段,其振幅谱为设定迭代初值为:
其中0<q≤1,
其中fM为振幅谱的截止频率。建立迭代uk=P(uk-1),其中算子参数由最小二乘得到。由于算子是压缩映射,可以得到不动点u*,则估计的子波振幅谱为记为Lk(f)。
6)计算质心频率;
对于第k个分子窗中的地震记录片段,质心频率为
式中Fc为子波振幅谱的截止频率。理论上,地震子波在传播的过程中质心频率是随着传播时间衰减的,然而实际中(如图7所示),当地层较薄时,地震波之间的互相干涉会使得质心频率突然增大或者突然减小(如图8所示)。为此,本发明在计算Q值之前对计算得到的质心频率做多次光滑滤波,并计算滤波前后的绝对误差,去掉绝对误差较大的一些质心频率后,再次拟合滤波,得到最终的质心频率以减少地震波干涉的影响(如图9所示)。
7)估算地下介质Q场;
用质心频率偏移法估算Q值,相应的估算公式如下
式中:分别为t1时刻频谱的质心频率和方差;是t2时刻频谱的质心频率;时差Δt=t2-t1;为t1和t2之间地层介质的品质因子。进一步推导,可得t时刻地下介质的等效Q值
式中fc,0和为初始时刻地震子波振幅谱的质心频率和方差,即
用步骤6)中计算得到的质心频率代替fc(t),可得到各个窗中心点处的Q值
Tk表示第k个分子窗的中心点所在的时间位置。对窗中心点处估得的Q值插值后,就可以得到稳定的地下介质Q场(如图10所示)。
该方法用地震道包络局部峰值作为地层结构约束来构造自适应分子时窗,并用这些时窗生成分子标架对地震记录做分子分解;接着在分子分解时频域利用非线性压缩映射提取时变子波振幅谱,并计算质心频率;然后对质心频率做滤波筛选,以减少子波干涉和噪声影响;最后用质心频率偏移法计算每个时窗对应的等效衰减因子Q,并插值得到稳定的地下介质Q场。所述分子时窗的构造方法,可以保证每个时窗内至少有一个反射子波,在减少时窗数量、提高计算效率的同时,使得这些时窗在横向上与地层结构有关,有利于保持所估计Q场的横向连续性。本发明能够有效克服复杂地质结构变化和噪声干扰的影响,得到稳定的地下介质Q场。
Claims (4)
1.一种数据驱动的地下介质Q场估计方法,其特征在于,采用地震包络局部峰值作为地层结构约束来构造自适应于地层结构的分子时窗,然后采用非线性压缩映射子波振幅谱估计方法从地震道分子分解时频谱中估计正比于时变子波振幅谱的分量,接着计算质心频率,并对质心频率做筛选,减少子波干涉和噪声的影响,最后用质心频率偏移法估计得到稳健的地下介质Q场。
2.根据权利要求1所述的一种数据驱动的地下介质Q场估计方法,其特征在于,具体包括以下步骤:
1)提取地震道包络局部峰值:
2)生成满足单位分割的原子窗集:
选择基本原子窗函数G(t),令Gj(t)=G(t-jΔt)表示中心位于第j个采样点的原子窗,对原子窗族{Gj:1≤j≤N}按下式归一化
<mrow>
<msub>
<mi>g</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>G</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>/</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>N</mi>
</munderover>
<msub>
<mi>G</mi>
<mi>i</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
得一组单位分割原子窗集{gj:1≤j≤N},这里N为地震道采样点的个数;
3)构造自适应分子时窗:
选择每相邻两个包络局部峰值点的中点作为各分子窗的边界点,将边界点间的满足单位分割的原子窗叠加起来就形成了分子窗;设第k个分子窗对应的包络局部峰值点位于Pk,前一个包络局部峰值点位于Pk-1(P0=-P1),后一个包络局部峰值点位于Pk+1,则此分子窗ψk(t),对应的第一个原子窗的中心位于Mk-1+1=(Pk-1+Pk)/2+1(M0=0),对应的最后一个原子窗的中心位于Mk=(Pk+Pk+1)/2,分子窗由这中间的Mk-Mk-1个原子窗叠加得到,即分子窗ψk(t)由下式表示
<mrow>
<msub>
<mi>&psi;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<msub>
<mi>M</mi>
<mrow>
<mi>k</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>+</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>M</mi>
<mi>k</mi>
</msub>
</munderover>
<msub>
<mi>g</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
令L为分子窗的个数,分子窗族{ψk(t):1≤k≤L}也构成单位分割;
4)对步骤3)得到的自适应分子时窗进行能量归一化:
令Ek表示第k个分子窗的能量,即
<mrow>
<msub>
<mi>E</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>&psi;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mo>|</mo>
<mo>=</mo>
<msup>
<mrow>
<mo>&lsqb;</mo>
<msubsup>
<mo>&Integral;</mo>
<mi>&infin;</mi>
<mrow>
<mo>+</mo>
<mi>&infin;</mi>
</mrow>
</msubsup>
<mo>|</mo>
<msub>
<mi>&psi;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<msup>
<mo>|</mo>
<mn>2</mn>
</msup>
<mi>d</mi>
<mi>t</mi>
<mo>&rsqb;</mo>
</mrow>
<mrow>
<mn>1</mn>
<mo>/</mo>
<mn>2</mn>
</mrow>
</msup>
</mrow>
能量归一化以后的分子分析时窗为{ψk(t)/Ek:1≤k≤L};
对分子分析时窗进行平移和调制后,得到一组分子标架;此时,分子分解时频变换定义为
<mrow>
<msub>
<mover>
<mi>s</mi>
<mo>^</mo>
</mover>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<msub>
<mi>E</mi>
<mi>k</mi>
</msub>
</mfrac>
<munderover>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mi>&infin;</mi>
</munderover>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<msub>
<mi>&psi;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mn>2</mn>
<mi>i</mi>
<mi>&pi;</mi>
<mi>f</mi>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
其中f为频率;
5)利用非线性压缩映射提取时变子波振幅谱:
设x0∈(a,b),取δ0>0,使得对于任意给定的δ>0,定义L2[a,b]空间中的子集
如果且q≤2,则有任取定义函数
<mrow>
<msub>
<mi>F</mi>
<mi>q</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>u</mi>
<mo>;</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mi>a</mi>
<mi>x</mi>
</msubsup>
<msup>
<mi>u</mi>
<mi>q</mi>
</msup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mi>a</mi>
<mi>b</mi>
</msubsup>
<msup>
<mi>u</mi>
<mi>q</mi>
</msup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>,</mo>
<mo>&ForAll;</mo>
<mi>x</mi>
<mo>&Element;</mo>
<mo>&lsqb;</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>&rsqb;</mo>
</mrow>
此处0<q≤1;显然,0≤Fq(u;x)≤1;设cq>0,α,β>1,对于定义上的算子P为
<mrow>
<mi>P</mi>
<mo>&lsqb;</mo>
<mi>u</mi>
<mo>;</mo>
<mi>x</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<msub>
<mi>c</mi>
<mi>q</mi>
</msub>
<msubsup>
<mi>F</mi>
<mi>q</mi>
<mi>&alpha;</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>u</mi>
<mo>;</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<msup>
<mrow>
<mo>&lsqb;</mo>
<mn>1</mn>
<mo>-</mo>
<msub>
<mi>F</mi>
<mi>q</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>u</mi>
<mo>;</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mi>&beta;</mi>
</msup>
<mo>+</mo>
<mi>&delta;</mi>
<mo>,</mo>
<mo>&ForAll;</mo>
<mi>x</mi>
<mo>&Element;</mo>
<mo>&lsqb;</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>&rsqb;</mo>
</mrow>
P是到的非线性压缩映射算子;
对于第k个分子窗中的地震记录片段,其振幅谱为设定迭代初值为:其中0<q≤1,
其中fM为振幅谱的截止频率;建立迭代uk=P(uk-1),其中算子参数由最小二乘得到;由于算子是压缩映射,得到不动点u*,则估计的子波振幅谱为记为Lk(f);
6)计算质心频率:
对于第k个分子窗中的地震记录片段,质心频率fc,k为
<mrow>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<mi>f</mi>
<mo>|</mo>
<msub>
<mi>L</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<mo>|</mo>
<msub>
<mi>L</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
</mfrac>
</mrow>
式中Fc为子波振幅谱的截止频率;对计算得到的质心频率处理后得到最终的质心频率
7)估算地下介质Q场:
用质心频率偏移法估算Q值,相应的估算公式如下
<mrow>
<msub>
<mi>Q</mi>
<mrow>
<msub>
<mi>t</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<msub>
<mi>t</mi>
<mn>2</mn>
</msub>
</mrow>
</msub>
<mo>=</mo>
<mi>&pi;</mi>
<mi>&Delta;</mi>
<mi>t</mi>
<mfrac>
<msubsup>
<mi>&sigma;</mi>
<mn>1</mn>
<mn>2</mn>
</msubsup>
<mrow>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<msub>
<mi>t</mi>
<mn>1</mn>
</msub>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<msub>
<mi>t</mi>
<mn>2</mn>
</msub>
</mrow>
</msub>
</mrow>
</mfrac>
</mrow>
式中:分别为t1时刻频谱的质心频率和方差;是t2时刻频谱的质心频率;时差Δt=t2-t1;为t1和t2之间地层介质的品质因子;进一步推导,得到t时刻地下介质的等效Q值
<mrow>
<msub>
<mi>Q</mi>
<mi>e</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mi>&pi;&sigma;</mi>
<mn>0</mn>
<mn>2</mn>
</msubsup>
<mi>t</mi>
</mrow>
<mrow>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mn>0</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>f</mi>
<mi>c</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
</mrow>
式中fc,0和为初始时刻地震子波振幅谱的质心频率和方差,即
<mrow>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mn>0</mn>
</mrow>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<mi>f</mi>
<mo>|</mo>
<mover>
<mi>w</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<mo>|</mo>
<mover>
<mi>w</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
</mfrac>
</mrow>
2
<mrow>
<msubsup>
<mi>&sigma;</mi>
<mn>0</mn>
<mn>2</mn>
</msubsup>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<msup>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>-</mo>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mn>0</mn>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>|</mo>
<mover>
<mi>w</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<msub>
<mi>F</mi>
<mi>c</mi>
</msub>
</msubsup>
<mo>|</mo>
<mover>
<mi>w</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>d</mi>
<mi>f</mi>
</mrow>
</mfrac>
</mrow>
用步骤6)中计算得到的质心频率代替fc(t),得到各个窗中心点处的Q值
<mrow>
<msub>
<mi>Q</mi>
<mrow>
<mi>e</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mi>&pi;&sigma;</mi>
<mn>0</mn>
<mn>2</mn>
</msubsup>
<msub>
<mi>T</mi>
<mi>k</mi>
</msub>
</mrow>
<mrow>
<msub>
<mi>f</mi>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mn>0</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mover>
<mi>f</mi>
<mo>~</mo>
</mover>
<mrow>
<mi>c</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
</mrow>
</mfrac>
</mrow>
Tk表示第k个分子窗的中心点所对应的时间位置;对窗中心点处估得的Q值插值后即得到稳定的地下介质Q场。
3.根据权利要求1所述的一种数据驱动的地下介质Q场估计方法,其特征在于,步骤1)中提取地震道包络局部峰值具体过程为:设s*(t)为地震道s(t)的Hilbert变换,则
a(t)=[s(t)2+s*(t)2]1/2
为地震道s(t)的包络;由上式计算地震道的包络,并提取包络局部峰值点。
4.根据权利要求1所述的一种数据驱动的地下介质Q场估计方法,其特征在于,步骤6)中对计算得到的质心频率做多次光滑滤波,并计算滤波前后的绝对误差,去掉绝对误差较大的质心频率后,再次拟合滤波,得到最终的质心频率
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710543599.7A CN107272062B (zh) | 2017-07-05 | 2017-07-05 | 一种数据驱动的地下介质q场估计方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710543599.7A CN107272062B (zh) | 2017-07-05 | 2017-07-05 | 一种数据驱动的地下介质q场估计方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107272062A true CN107272062A (zh) | 2017-10-20 |
CN107272062B CN107272062B (zh) | 2018-12-07 |
Family
ID=60072351
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710543599.7A Active CN107272062B (zh) | 2017-07-05 | 2017-07-05 | 一种数据驱动的地下介质q场估计方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107272062B (zh) |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108333630A (zh) * | 2018-02-02 | 2018-07-27 | 林腾飞 | 一种利用数据自适应性时窗估算几何属性的方法 |
CN108646289A (zh) * | 2018-03-19 | 2018-10-12 | 中国海洋石油集团有限公司 | 一种估计地震品质因子的方法 |
CN110515127A (zh) * | 2019-09-26 | 2019-11-29 | 中国石油大学(北京) | 一种地震品质因子确定方法、装置、设备、介质 |
CN110579805A (zh) * | 2019-10-17 | 2019-12-17 | 西南石油大学 | 一种基于自适应增益限反q滤波的地震资料处理方法 |
CN110824564A (zh) * | 2018-08-08 | 2020-02-21 | 中国石油化工股份有限公司 | 用于近地表品质因子q值反演的衰减曲线层析剥离方法 |
CN111381107A (zh) * | 2020-06-01 | 2020-07-07 | 成都市易冲半导体有限公司 | 一种无线充电高精度q值检测方法和电路 |
CN112099088A (zh) * | 2020-09-16 | 2020-12-18 | 中油奥博(成都)科技有限公司 | 一种基于高密度光纤地震数据的油气指示及表征方法 |
CN112526611A (zh) * | 2019-09-18 | 2021-03-19 | 中国石油天然气集团有限公司 | 表层地震波品质因子的提取方法及装置 |
WO2023173399A1 (en) * | 2022-03-18 | 2023-09-21 | Saudi Arabian Oil Company | Holographic inversion for hydrocarbon indicator |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102288992A (zh) * | 2011-04-26 | 2011-12-21 | 中国海洋石油总公司 | 利用地震信号包络峰值瞬时频率估计介质品质因子的方法 |
CN103412324A (zh) * | 2013-07-17 | 2013-11-27 | 西安交通大学 | 一种估计介质品质因子的epifvo 方法 |
CN106324663A (zh) * | 2015-06-17 | 2017-01-11 | 中国石油化工股份有限公司 | 一种品质因子的获取方法 |
-
2017
- 2017-07-05 CN CN201710543599.7A patent/CN107272062B/zh active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102288992A (zh) * | 2011-04-26 | 2011-12-21 | 中国海洋石油总公司 | 利用地震信号包络峰值瞬时频率估计介质品质因子的方法 |
CN103412324A (zh) * | 2013-07-17 | 2013-11-27 | 西安交通大学 | 一种估计介质品质因子的epifvo 方法 |
CN106324663A (zh) * | 2015-06-17 | 2017-01-11 | 中国石油化工股份有限公司 | 一种品质因子的获取方法 |
Non-Patent Citations (3)
Title |
---|
余连勇: "时间域质心频移法估算品质因子Q", 《西南石油大学学报(自然科学版)》 * |
李霞: "《油气地球物理技术新进展-第84届SEG年会论文摘要》", 31 December 2015 * |
王清振: "基于变时窗Gabor变换的高分辨率处理技术研究与应用", 《中国海上油气》 * |
Cited By (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108333630A (zh) * | 2018-02-02 | 2018-07-27 | 林腾飞 | 一种利用数据自适应性时窗估算几何属性的方法 |
CN108646289A (zh) * | 2018-03-19 | 2018-10-12 | 中国海洋石油集团有限公司 | 一种估计地震品质因子的方法 |
CN108646289B (zh) * | 2018-03-19 | 2019-09-13 | 中国海洋石油集团有限公司 | 一种估计地震品质因子的方法 |
CN110824564A (zh) * | 2018-08-08 | 2020-02-21 | 中国石油化工股份有限公司 | 用于近地表品质因子q值反演的衰减曲线层析剥离方法 |
CN112526611A (zh) * | 2019-09-18 | 2021-03-19 | 中国石油天然气集团有限公司 | 表层地震波品质因子的提取方法及装置 |
CN110515127A (zh) * | 2019-09-26 | 2019-11-29 | 中国石油大学(北京) | 一种地震品质因子确定方法、装置、设备、介质 |
CN110579805A (zh) * | 2019-10-17 | 2019-12-17 | 西南石油大学 | 一种基于自适应增益限反q滤波的地震资料处理方法 |
CN110579805B (zh) * | 2019-10-17 | 2021-03-12 | 西南石油大学 | 一种基于自适应增益限反q滤波的地震资料处理方法 |
CN111381107A (zh) * | 2020-06-01 | 2020-07-07 | 成都市易冲半导体有限公司 | 一种无线充电高精度q值检测方法和电路 |
CN112099088A (zh) * | 2020-09-16 | 2020-12-18 | 中油奥博(成都)科技有限公司 | 一种基于高密度光纤地震数据的油气指示及表征方法 |
CN112099088B (zh) * | 2020-09-16 | 2022-04-12 | 中油奥博(成都)科技有限公司 | 一种基于高密度光纤地震数据的油气指示及表征方法 |
WO2023173399A1 (en) * | 2022-03-18 | 2023-09-21 | Saudi Arabian Oil Company | Holographic inversion for hydrocarbon indicator |
Also Published As
Publication number | Publication date |
---|---|
CN107272062B (zh) | 2018-12-07 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN107272062B (zh) | 一种数据驱动的地下介质q场估计方法 | |
US11333783B2 (en) | Integrated method for estimation of seismic wavelets and synthesis of seismic records in depth domain | |
US10920585B2 (en) | Determining sand-dune velocity variations | |
CN103376464B (zh) | 一种地层品质因子反演方法 | |
EP2567259B1 (en) | Q tomography method | |
US20230083651A1 (en) | Method and system for analyzing filling for karst reservoir based on spectrum decomposition and machine learning | |
CN105093294B (zh) | 基于可变模态分解的地震波衰减梯度估计方法 | |
CN104849756B (zh) | 一种提高地震数据分辨率增强有效弱信号能量的方法 | |
CN105388518A (zh) | 一种质心频率与频谱比联合的井中地震品质因子反演方法 | |
CN107132579A (zh) | 一种保地层结构的地震波衰减补偿方法 | |
US11181653B2 (en) | Reservoir characterization utilizing ReSampled seismic data | |
CN104237945B (zh) | 一种地震资料自适应高分辨处理方法 | |
CN103163554A (zh) | 利用零偏vsp资料估计速度和q值的自适应波形的反演方法 | |
Khosro Anjom et al. | Full-waveform matching of VP and VS models from surface waves | |
CN106019376B (zh) | 一种频率驱动空变q值模型构建的地震波补偿方法 | |
CN107179550B (zh) | 一种数据驱动的地震信号零相位反褶积方法 | |
CN110471113A (zh) | 基于非稳态地震资料的反演动校正方法、装置及存储介质 | |
CN103439740A (zh) | 基于偶极地震子波多重积分的相对阻抗预测的方法及装置 | |
CN106707334A (zh) | 一种提高地震资料分辨率的方法 | |
WO2005026776A1 (en) | Wide-offset-range pre-stack depth migration method for seismic exploration | |
CN102854530B (zh) | 基于对数时频域双曲平滑的动态反褶积方法 | |
CN108646290A (zh) | 一种基于模型定量补偿的薄层反演方法 | |
Moon et al. | Collocated cokriging and neural-network multi-attribute transform in the prediction of effective porosity: A comparative case study for the Second Wall Creek Sand of the Teapot Dome field, Wyoming, USA | |
CN112147700A (zh) | 速度异常区的低频模型构建方法及系统 | |
Gao et al. | Acquisition and processing pitfall with clipped traces in surface-wave analysis |
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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |