WO2013152468A1 - 一种地层品质因子反演方法 - Google Patents

一种地层品质因子反演方法 Download PDF

Info

Publication number
WO2013152468A1
WO2013152468A1 PCT/CN2012/001686 CN2012001686W WO2013152468A1 WO 2013152468 A1 WO2013152468 A1 WO 2013152468A1 CN 2012001686 W CN2012001686 W CN 2012001686W WO 2013152468 A1 WO2013152468 A1 WO 2013152468A1
Authority
WO
WIPO (PCT)
Prior art keywords
frequency
amplitude spectrum
vertical seismic
seismic section
wave
Prior art date
Application number
PCT/CN2012/001686
Other languages
English (en)
French (fr)
Inventor
张固澜
王熙明
张庆红
李彦鹏
彭继新
赵予凤
容娇君
李可恩
金其虎
郭晓玲
Original Assignee
中国石油天然气集团公司
中国石油集团东方地球物理勘探有限责任公司
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by 中国石油天然气集团公司, 中国石油集团东方地球物理勘探有限责任公司 filed Critical 中国石油天然气集团公司
Priority to RU2014145635/28A priority Critical patent/RU2579164C1/ru
Priority to EP12874053.7A priority patent/EP2837953A4/en
Priority to US14/394,100 priority patent/US20150168573A1/en
Publication of WO2013152468A1 publication Critical patent/WO2013152468A1/zh

Links

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
    • 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/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • 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
    • G01V1/42—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging using generators in one well and receivers elsewhere or vice versa
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2200/00—Details of seismic or acoustic prospecting or detecting in general
    • G01V2200/10—Miscellaneous details
    • G01V2200/14—Quality control
    • 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
    • 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/63—Seismic attributes, e.g. amplitude, polarity, instant phase
    • 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/67—Wave propagation modeling
    • G01V2210/677—Spectral; Pseudo-spectral

Definitions

  • the invention relates to a seismic exploration data processing technology, and is a formation quality factor inversion method which utilizes the vertical wave amplitude spectrum property of vertical seismic section (VSP) data and has good stability.
  • VSP vertical seismic section
  • the absorption attenuation of the stratum is mainly manifested by the attenuation of the amplitude of the seismic wave during the propagation process, the distortion of the phase, and the lowering of the frequency, and the high-frequency part decays faster than the low-frequency part, and the attenuation is faster in the shallow layer than in the deep layer, which seriously reduces the resolution of the seismic data. rate.
  • the inversion method of formation quality factor Q mainly uses the logarithmic spectrum ratio method, the center frequency offset method and the peak frequency offset method, the combination of scanning analysis technology and time-frequency analysis, and the multi-window spectrum analysis method for the amplitude spectrum of seismic wavelets. Wait.
  • the center frequency offset method and the peak frequency offset method assume seismic wave vibration
  • the spectrum can be represented by a Gaussian spectrum; the time-frequency analysis assumes that the seismic wavelet is zero-phase.
  • Mathneey and Nowack proposed the instantaneous frequency matching method, which uses an iterative process to modify the causal decay operator, so that the operator acts on the weighted instantaneous frequency at the envelope peak after the reference pulse and the weighted instantaneous frequency at the peak of the target pulse envelope.
  • the closest thus the quality factor of the media, they used this method to estimate the attenuation of the crustal diffraction seismic data;
  • Dasios et al. used the instantaneous frequency matching method to estimate the attenuation of the full-wave array sonic logging.
  • This method overcomes Some disadvantages of the log-ratio method, such as the need to select a variable frequency band range, etc., but this method requires the Hilbert transform method to calculate the instantaneous frequency, and also uses a complex iterative process to match the instantaneous frequency. It is well known that the Hilbert transform pair Noise sensitive, so the use of instantaneous frequency matching in noise-containing seismic signals is limited. Barnes assumes that the source wavelet is an ideal bandpass wavelet, and gives a relationship between instantaneous frequency and Q value and transmission time, but the actual source wavelet and the ideal bandpass wavelet are quite different.
  • VSP vertical seismic profile
  • the downhole detector receives the vertical seismic section data, and the detector close to the source receives the monitoring wavelet signal corresponding to each vertical seismic section record;
  • step 5 In the frequency-wavenumber (FK) spectrum obtained in step 4), multiply the frequency-wavenumber (FK) spectrum corresponding to the up-wave by zero; then perform inverse Fourier transform in the wavenumber direction to obtain the amplitude spectrum; Performing an inverse Fourier transform on the obtained amplitude spectrum in the frequency direction to obtain a wave field two in the time domain;
  • a time window is opened from the first sample point in each downlink wave, and the signal in the time window is Fourier transformed to obtain an amplitude spectrum of each frequency; And the amplitude spectrum corresponding to each frequency is divided by the square of its corresponding frequency value to obtain an amplitude spectrum in the exponential form.
  • step 6) Repeat step 6) to obtain the amplitude spectrum of the exponential form of each frequency in all the down-going waves.
  • step 7 After taking the natural logarithm of the amplitude spectrum 2 in step 7), the quadratic function fitting related to the frequency is performed by the least square method, and the primary term coefficient and the quadratic term coefficient corresponding to the downlink wavelet are obtained;
  • step 8) to obtain the primary and quadratic coefficients corresponding to the descending wavelets of all the tracks in the vertical seismic section record;
  • step 10 Picking up the first to second of the monitoring wavelet recorded in step 1), and recording the corresponding monitoring wavelet signal for each vertical seismic section, starting from the beginning to the second, opening a time window, within the time window The signal is Fourier transformed to obtain the amplitude spectrum of each frequency of the monitored wavelet; and the amplitude spectrum corresponding to each frequency is divided by the square of its corresponding frequency value to obtain the amplitude spectrum of all the track wavelet indices. four; 11) After taking the natural logarithm of the amplitude spectrum in step 10), the frequency-dependent quadratic function fitting is performed by the least squares method, and the primary coefficient and the quadratic coefficient of the monitored wavelet spectrum are obtained;
  • step 11) Repeat step 11) to obtain the primary and secondary coefficients corresponding to the monitored wavelets of all the tracks in the vertical seismic section record;
  • step 12) obtaining an average value of the quadratic coefficient of all the tracks in the vertical seismic section record and the corresponding monitoring wavelet quadratic coefficient;
  • step 14 After taking the natural spectrum of the amplitude spectrum in step 7), subtract the product of the mean value of the quadratic coefficient of the seismic trace obtained in step 13) and the square of the frequency, and then use the square method to perform frequency-dependent Quadratic function fitting, obtaining primary term coefficients and quadratic coefficients;
  • step 15 using the first-order time of each track in the vertical seismic section record divided by the first-order coefficient of the track obtained in step 14) to obtain an equivalent Q (stratigraphic quality factor) value of one;
  • step 15 Repeat step 15) to obtain the equivalent formation quality factor value of all vertical seismic section tracks, and statistically smooth the equivalent Q (stratigraphic quality factor) values of all the tracks to obtain the equivalent Q (stratigraphic quality factor). Value two;
  • the layer Q ratio quality factor
  • step 18 Repeat step 18) until the corresponding layer Q (stratigraphic factor) values for each vertical seismic section record are inverted.
  • FIG. 1 Schematic diagram of the downlink wave
  • Figure 2 is a schematic diagram of the downlink wave intercepted
  • Figure 3 shows the amplitude spectrum of the downlink wave intercepted
  • the invention provides a method for inverting the formation quality factor of the downlink wave using the vertical seismic profile (VSP) data and having good stability, and the specific implementation steps are as follows:
  • the downhole detector receives the vertical seismic section data, and the detector close to the source receives the monitoring wavelet signal corresponding to each vertical seismic section record;
  • step 5 In the frequency-wavenumber (FK) spectrum obtained in step 4), multiply the frequency-wavenumber (FK) spectrum corresponding to the up-wave by zero; then perform inverse Fourier transform in the wavenumber direction to obtain the amplitude spectrum; Performing an inverse Fourier transform on the obtained amplitude spectrum in the frequency direction to obtain a wave field two in the time domain;
  • step 6) Repeat step 6) to obtain an amplitude spectrum of the exponential form of each frequency in all the down-going waves.
  • step 7 After taking the natural logarithm of the amplitude spectrum in step 7), the frequency-dependent quadratic function fitting is performed by the least square method to obtain the primary term coefficient and the quadratic coefficient corresponding to the downlink wavelet of the channel;
  • step 8) to obtain the primary and quadratic coefficients corresponding to the descending wavelets of all the tracks in the vertical seismic section record;
  • step 1) Picking up the first to second of the monitoring wavelet recorded in step 1), and recording the corresponding monitoring wavelet signal for each vertical seismic section, starting from the beginning to the second, opening a time window, within the time window
  • the signal is Fourier transformed to obtain the amplitude spectrum of each frequency of the monitored wavelet; and the amplitude spectrum corresponding to each frequency is divided by the square of its corresponding frequency value to obtain the amplitude spectrum of all the track wavelet indices.
  • step 11 After taking the natural logarithm of the amplitude spectrum in step 10), the frequency-dependent quadratic function fitting is performed by the least squares method, and the primary coefficient and the quadratic coefficient of the monitoring wavelet spectrum of the channel are obtained;
  • step 11) Repeat step 11) to obtain the primary and secondary coefficients corresponding to the monitored wavelets of all the tracks in the vertical seismic section record;
  • Step 12) Obtain the average of the quadratic coefficients of all the tracks in the vertical seismic section record and the corresponding monitoring wavelet quadratic coefficients;
  • step 14 After taking the natural spectrum of the amplitude spectrum in step 7), subtract the product of the mean value of the quadratic coefficient of the seismic trace obtained in step 13) and the square of the frequency, and then use the square method to perform frequency-dependent Quadratic function fitting, obtaining primary term coefficients and quadratic coefficients;
  • step 15 using the first-order time of each track in the vertical seismic section record divided by the primary term coefficient of the track obtained in step 14) to obtain an equivalent Q (stratigraphic quality factor) value of one;
  • step 15 Repeat step 15) to obtain the equivalent formation quality factor of all vertical seismic section tracks
  • the value is one, and the equivalent Q (stratigraphic quality factor) value of all the tracks is statistically smoothed to obtain the equivalent Q (stratigraphic quality factor) value of two;
  • the layer Q ratio quality factor
  • step 18 Repeat step 18) until the corresponding layer Q (stratigraphic factor) values for each vertical seismic section record are inverted. As shown in Fig. 4, as the depth increases, the layer Q also tends to increase.
  • the invention has strong anti-random interference capability, can eliminate the difference of the excitation wave, and not only the algorithm is simple but also greatly saves the workload, and the layer Q value of the inversion has good stability and high precision.

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因子反演的方法。该方法首先利用F-K（频率-波数）方法对VSP原始数据进行波场分离，得到下行波；然后选取下行子波和监控子波进行傅立叶变换得到振幅谱，对振幅谱进行多项式拟合从而得到等效Q，利用等效Q和层Q之间的公式，反演得到层Q。该方法具有很强抗随即干扰能力，能消除激发子波的差异，不仅算法简单还能大大节约工作量，反演的层Q值具有很好的稳定性、且精度很高。

Description

一种地层品质因子反演方法
技术领域
本发明涉及地震勘探数据处理技术， 是一种利用垂直地震剖面 （VSP) 数据的下行波振幅谱属性， 且具有很好的稳定性的地层品质因子反演方 法。
背景技术
随着地震勘探精度要求的提高， 需要高分辨率的地震资料进行油气储 层的详细描述， 而地层的吸收衰减作用是影响地震资料分辨率的一个主要 因素。 地层的吸收衰减主要表现为地震波在传播过程振幅发生衰减、 相位 发生畸变、 频率变低， 且高频部分比低频部分衰减更快， 在浅层比深层衰 减更快， 严重降低了地震资料的分辨率。 准确地估算地层的品质因子 （Q 因子） 值， 然后对叠前或叠后地震记录进行有效的反 Q滤波补偿， 可以使 地震剖面的浅、 中、 深层反射波波形基本一致， 中、 深层的高频部分得到 加强， 频谱被拓宽， 使其恢复原来的地震波形态， 消除子波时变影响， 从 而满足反褶积和子波估计要求子波时不变的假设前提， 并能有效地提高地 震剖面质量， 从而更有利于地震资料的处理和解释。
零偏 VSP (垂直地震剖面） 资料采集中， 炮点离井口很近， 不同深度 接收到的下行直达波具有相同的传播路径， 可直接利用不同深度地震记录 中下行直达波反演地层品质因子 (Q因子)， 并进行反 Q滤波， 从而提高 VSP 资料的分辨率及驱动地面地震处理提高分辨率； 因此如何利用零偏 VSP资 料进行精确的 Q提取， 具有很重要的实际应用价值。
地层品质因子 Q 的反演方法主要对地震子波的振幅谱利用对数谱比 法、中心频率偏移法和峰值频率偏移法、扫描分析技术与时频分析结合法、 多窗谱分析法等。 其中， 中心频率偏移法和峰值频率偏移法假设地震波振 幅谱可由高斯谱来表示； 时频分析法假设地震子波是零相位的。
Mathneey和 Nowack提出了瞬时频率匹配法， 即采用一个迭代过程修 改因果衰减算子， 使该算子作用于基准脉冲后的包络峰值处的加权瞬时频 率和目标脉冲包络峰值处的加权瞬时频率最接近， 由此反演出介质的品质 因子，他们用此方法估计了地壳绕射地震资料的衰减; Dasios等用瞬时频率 匹配法估计了全波列声波测井记录的衰减.这种方法克服了对数谱比法的 一些缺点， 比如不需要选择可变的频带范围等.但该方法需要利 Hilbert变 换法计算瞬时频率， 而且还要用复杂的迭代过程来匹配瞬时频率.众所周 知， Hilbert变换对噪声敏感， 因此瞬时频率匹配法在含噪地震信号中使用 受到限制。 Barnes假设震源子波为理想的带通子波， 给出了一个瞬时频率 和 Q值以及传输时间的关系， 但实际震源子波和理想带通子波差异较大。
以上所有方法对实际资料几乎没有应用，也没有公开如何利用 VSP资 料的下行波， 反演的层 Q值与地层的层速度值也几乎没有对应性， 无法评 价反演的 Q值的合理性。 另外， 以上所有方法都没有考虑采集中激发环境 造成的激发子波差异， 几乎不具有适用性和推广性， 势必将影响品质因子 Q的稳定性。
发明内容
本发明目的是提供一种利用垂直地震剖面 （VSP ) 数据的下行波的振 幅谱属性， 且具有很好的稳定性的地层品质因子反演方法。
本发明具体步骤包括：
1 ) 地面震源激发， 井下检波器接收得到垂直地震剖面数据， 靠近震 源的检波器接收得到每道垂直地震剖面记录对应的监控子波信号；
2 ) 拾取每道垂直地震剖面记录的初至一， 以及每道垂直地震剖面记 录对应的监控子波信号的初至二； 3 ) 对每道垂直地震剖面记录各个样点的时间都减去该道的初至时间 一， 从而将下行波拉平得到第一波场；
4)对第一波场先进行时间方向的傅里叶变换， 从而变换到频率域， 得 到整个垂直地震剖面记录的振幅谱； 然后对振幅谱在道号方向进行傅里叶 变换， 从而变换到波数域， 得到其频率 -波数 (F-K)谱；
5)在步骤 4)得到的频率 -波数 (F-K)谱中， 将上行波对应的频率 -波数 (F-K) 谱乘以零； 然后进行波数方向的反傅里叶变换， 得到其振幅谱； 再 对得到的振幅谱进行频率方向的反傅里叶变换， 得到时间域的波场二；
6 ) 在波场二中， 在每道的下行波中从第一个样点开始往后开一时窗， 对该时窗内的信号进行傅里叶变换， 得到每个频率的振幅谱一； 并将每个 频率对应的振幅谱都除以其对应的频率值的平方， 得到指数形式的振幅谱
7 )重复步骤 6) ， 得到所有道下行波中每个频率的指数形式的振幅谱
8 )对步骤 7) 中的振幅谱二取自然对数后， 利用最小二乘法进行与频 率有关的二次函数拟合， 得到该道下行子波对应的一次项系数和二次项系 数；
9)重复步骤 8)得到垂直地震剖面记录中所有道的下行子波对应的一 次项和二次项系数；
10) 拾取步骤 1 ) 中记录的监控子波的初至二， 并对每道垂直地震剖 面记录对应的监控子波信号， 从初至二开始往后开一个时窗， 对该时窗内 的信号进行傅里叶变换， 得到监控子波每个频率的振幅谱三； 并将每个频 率对应的振幅谱都除以其对应的频率值的平方， 得到所有道监控子波指数 形式的振幅谱四； 11 )对步骤 10 ) 中的振幅谱四取自然对数后， 利用最小二乘法进行与 频率有关的二次函数拟合， 得到该道监控子波频谱的一次项系数和二次项 系数；
12)重复步骤 11 )得到垂直地震剖面记录中所有道对应的监控子波对 应的一次项系数和二次项系数；
13 )求取步骤 12 )得到垂直地震剖面记录中所有道的二次项系数和其 对应的监控子波二次项系数的平均值；
14) 将步骤 7) 中的振幅谱二取自然对数后， 减去步骤 13) 中得到的 该地震道二次项系数的平均值与频率平方的乘积后， 利用二乘法进行与频 率有关的二次函数拟合， 得到一次项系数和二次项系数；
15)利用垂直地震剖面记录中每道的初至时间一除以步骤 14)中得到 的该道的一次项系数， 得到等效 Q (地层品质因子） 值一；
16)重复步骤 15 )得到所有垂直地震剖面记录道的等效地层品质因子 值一,并对所有道的等效 Q (地层品质因子）值一进行统计平滑， 得到等效 Q (地层品质因子） 值二；
17 ) 用垂直地震剖面记录每道的初至时间一除以该道记录对应的等效 Q (地层品质因子） 值二， 得到该道的吸收系数；
18 ) 利用相邻垂直地震剖面记录道的初至一之间的差值除以相邻道的 吸收系数之间的差值， 得到该道垂直地震剖面记录对应的层 Q (地层品质 因子） 值；
19)重复步骤 18) ， 直到每道垂直地震剖面记录对应的层 Q (地层品 质因子） 值都反演完。
附图说明
图 1 下行波示意图； 图 2 截取的下行波示意图；
图 3 截取的下行波的振幅谱；
图 4 本发明反演的层 Q值。
具体实施方式
以下结合附图详细说明本发明。
本发明一种利用垂直地震剖面 （VSP ) 数据的下行波的振幅谱属性， 且具有很好的稳定性的地层品质因子反演方法， 具体实现步骤如下：
1 ) 地面震源激发， 井下检波器接收得到垂直地震剖面数据， 靠近震 源的检波器接收得到每道垂直地震剖面记录对应的监控子波信号；
2 ) 拾取每道垂直地震剖面记录的初至一， 以及每道垂直地震剖面记 录对应的监控子波信号的初至二；
3 ) 对每道垂直地震剖面记录各个样点的时间都减去该道的初至时间 一， 从而将下行波拉平得到第一波场；
4)对第一波场先进行时间方向的傅里叶变换， 从而变换到频率域， 得 到整个垂直地震剖面记录的振幅谱； 然后对振幅谱在道号方向进行傅里叶 变换， 从而变换到波数域， 得到其频率 -波数 (F-K)谱；
5)在步骤 4)得到的频率 -波数 (F-K)谱中， 将上行波对应的频率 -波数 (F-K) 谱乘以零； 然后进行波数方向的反傅里叶变换， 得到其振幅谱； 再 对得到的振幅谱进行频率方向的反傅里叶变换， 得到时间域的波场二；
6 ) 在波场二中， 如图 1 所示， 在每道的下行波中从第一个样点开始 往后开一时窗， 如图 2所示， 对该时窗内的信号进行傅里叶变换， 得到每 个频率的振幅谱一， 如图 3所示； 并将每个频率对应的振幅谱都除以其对 应的频率值的平方， 得到指数形式的振幅谱二；
7 )重复步骤 6 ) ， 得到所有道下行波中每个频率的指数形式的振幅谱 8)对步骤 7 ) 中的振幅谱二取自然对数后， 利用最小二乘法进行与频 率有关的二次函数拟合， 得到该道下行子波对应的一次项系数和二次项系 数；
9)重复步骤 8)得到垂直地震剖面记录中所有道的下行子波对应的一 次项和二次项系数；
10) 拾取步骤 1 ) 中记录的监控子波的初至二， 并对每道垂直地震剖 面记录对应的监控子波信号， 从初至二开始往后开一个时窗， 对该时窗内 的信号进行傅里叶变换， 得到监控子波每个频率的振幅谱三； 并将每个频 率对应的振幅谱都除以其对应的频率值的平方， 得到所有道监控子波指数 形式的振幅谱四；
11 )对步骤 10) 中的振幅谱四取自然对数后， 利用最小二乘法进行与 频率有关的二次函数拟合， 得到该道监控子波频谱的一次项系数和二次项 系数；
12)重复步骤 11 )得到垂直地震剖面记录中所有道对应的监控子波对 应的一次项系数和二次项系数；
13)求取步骤 12)得到垂直地震剖面记录中所有道的二次项系数和其 对应的监控子波二次项系数的平均值；
14) 将步骤 7 ) 中的振幅谱二取自然对数后， 减去步骤 13 ) 中得到的 该地震道二次项系数的平均值与频率平方的乘积后， 利用二乘法进行与频 率有关的二次函数拟合， 得到一次项系数和二次项系数；
15 )利用垂直地震剖面记录中每道的初至时间一除以步骤 14)中得到 的该道的一次项系数， 得到等效 Q (地层品质因子） 值一；
16 )重复步骤 15 )得到所有垂直地震剖面记录道的等效地层品质因子 值一，并对所有道的等效 Q (地层品质因子）值一进行统计平滑， 得到等效 Q (地层品质因子） 值二；
17 ) 用垂直地震剖面记录每道的初至时间一除以该道记录对应的等效 Q (地层品质因子） 值二， 得到该道的吸收系数；
18 ) 利用相邻垂直地震剖面记录道的初至一之间的差值除以相邻道的 吸收系数之间的差值， 得到该道垂直地震剖面记录对应的层 Q (地层品质 因子） 值；
19)重复步骤 18 ) ， 直到每道垂直地震剖面记录对应的层 Q (地层品 质因子） 值都反演完。 如图 4， 随着深度的增大， 层 Q也呈增大的趋势。 工业实用性
本发明具有很强抗随机干扰能力， 能消除激发子波的差异， 不仅算法 简单还能大大节约工作量，反演的层 Q值具有很好的稳定性、且精度很高。

Claims

权利要求书
1、 一种地层品质因子反演方法， 特点是具体步骤包括：
1 ) 地面震源激发， 井下检波器接收得到垂直地震剖面数据， 靠近震 源的检波器接收得到每道垂直地震剖面记录对应的监控子波信号；
2 ) 拾取每道垂直地震剖面记录的初至一， 以及每道垂直地震剖面记 录对应的监控子波信号的初至二；
3 ) 对每道垂直地震剖面记录各个样点的时间都减去该道的初至时间 一， 从而将下行波拉平得到第一波场；
4)对第一波场先进行时间方向的傅里叶变换， 从而变换到频率域， 得 到整个垂直地震剖面记录的振幅谱； 然后对振幅谱在道号方向进行傅里叶 变换， 从而变换到波数域， 得到其频率-波数 (F-K)谱；
5)在步骤 4)得到的频率 -波数 (F-K)谱中， 将上行波对应的频率 -波数 (F-K) 谱乘以零； 然后进行波数方向的反傅里叶变换， 得到其振幅谱； 再 对得到的振幅谱进行频率方向的反傅里叶变换， 得到时间域的波场二；
6) 在波场二中， 在每道的下行波中从第一个样点开始往后开一时窗， 对该时窗内的信号进行傅里叶变换， 得到每个频率的振幅谱一； 并将每个 频率对应的振幅谱都除以其对应的频率值的平方， 得到指数形式的振幅谱
7 )重复步骤 6) ， 得到所有道下行波中每个频率的指数形式的振幅谱
8)对步骤 7) 中的振幅谱二取自然对数后， 利用最小二乘法进行与频 率有关的二次函数拟合， 得到该道下行子波对应的一次项系数和二次项系 数；
9)重复步骤 8)得到垂直地震剖面记录中所有道的下行子波对应的一 次项和二次项系数；
10) 拾取步骤 1 ) 中记录的监控子波的初至二， 并对每道垂直地震剖 面记录对应的监控子波信号， 从初至二开始往后开一个时窗， 对该时窗内 的信号进行傅里叶变换， 得到监控子波每个频率的振幅谱三； 并将每个频 率对应的振幅谱都除以其对应的频率值的平方， 得到所有道监控子波指数 形式的振幅谱四；
11 )对步骤 10) 中的振幅谱四取自然对数后， 利用最小二乘法进行与 频率有关的二次函数拟合， 得到该道监控子波频谱的一次项系数和二次项 系数；
12)重复步骤 11 )得到垂直地震剖面记录中所有道对应的监控子波对 应的一次项系数和二次项系数；
13)求取步骤 12 )得到垂直地震剖面记录中所有道的二次项系数和其 对应的监控子波二次项系数的平均值；
14) 将步骤 7 ) 中的振幅谱二取自然对数后， 减去步骤 13) 中得到的 该地震道二次项系数的平均值与频率平方的乘积后， 利用二乘法进行与频 率有关的二次函数拟合， 得到一次项系数和二次项系数；
15)利用垂直地震剖面记录中每道的初至时间一除以步骤 14)中得到 的该道的一次项系数， 得到等效地层品质因子值一；
16)重复步骤 15 )得到所有垂直地震剖面记录道的等效地层品质因子 值一，并对所有道的等效地层品质因子值一进行统计平滑， 得到等效地层 品质因子值二；
17 ) 用垂直地震剖面记录每道的初至时间一除以该道记录对应的等效 地层品质因子值二， 得到该道的吸收系数；
18 ) 利用相邻垂直地震剖面记录道的初至一之间的差值除以相邻道的 吸收系数之间的差值， 得到该道垂直地震剖面记录对应的层地层品质因子 值；
19)重复步骤 18 ) ， 直到每道垂直地震剖面记录对应的层地层品质因 子值都反演完。
PCT/CN2012/001686 2012-04-13 2012-12-11 一种地层品质因子反演方法 WO2013152468A1 (zh)

Priority Applications (3)

Application Number Priority Date Filing Date Title
RU2014145635/28A RU2579164C1 (ru) 2012-04-13 2012-12-11 Способ обращения для определения добротности геологической среды
EP12874053.7A EP2837953A4 (en) 2012-04-13 2012-12-11 METHOD FOR INVERTING A GEOLOGICAL QUALITY FACTOR
US14/394,100 US20150168573A1 (en) 2012-04-13 2012-12-11 Geologic quality factor inversion method

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201210109416.8 2012-04-13
CN201210109416.8A CN103376464B (zh) 2012-04-13 2012-04-13 一种地层品质因子反演方法

Publications (1)

Publication Number Publication Date
WO2013152468A1 true WO2013152468A1 (zh) 2013-10-17

Family

ID=49326985

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2012/001686 WO2013152468A1 (zh) 2012-04-13 2012-12-11 一种地层品质因子反演方法

Country Status (5)

Country Link
US (1) US20150168573A1 (zh)
EP (1) EP2837953A4 (zh)
CN (1) CN103376464B (zh)
RU (1) RU2579164C1 (zh)
WO (1) WO2013152468A1 (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2016008105A1 (zh) * 2014-07-15 2016-01-21 杨顺伟 一种基于柯西分布的叠后波阻抗反演方法
CN107179544A (zh) * 2017-07-11 2017-09-19 中国石油集团川庆钻探工程有限公司地球物理勘探公司 转换波宽频延拓及提高转换波数据分辨率的方法
CN109471162A (zh) * 2018-10-08 2019-03-15 中国石油天然气集团有限公司 层间多次波处理方法、系统、电子设备及可读介质
CN111965695A (zh) * 2020-08-17 2020-11-20 山西潞安环保能源开发股份有限公司五阳煤矿 一种基于反射槽波的小断层落差探测方法
CN112162314A (zh) * 2020-09-25 2021-01-01 武汉市工程科学技术研究院 一种人工地震信号剖面的二维插值方法
CN112305587A (zh) * 2019-08-02 2021-02-02 中国石油化工股份有限公司 恢复地震数据分辨率的方法、存储介质及计算机设备
CN112578436A (zh) * 2019-09-27 2021-03-30 中国石油化工股份有限公司 一种子波提取方法及系统
RU2751088C2 (ru) * 2016-09-09 2021-07-08 Чайна Петролеум Энд Кемикл Корпорейшн Способ и система обработки сейсмических данных
CN116840916A (zh) * 2023-07-04 2023-10-03 成都理工大学 一种地震速度信号和加速度信号联合子波提取方法

Families Citing this family (35)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103698812B (zh) * 2013-12-26 2016-06-29 中国石油天然气集团公司 利用叠前地震道集计算地层品质因数的方法及装置
CN105095559A (zh) * 2014-05-09 2015-11-25 中央大学 实施全息希尔伯特频谱分析的方法与系统
CN105388523A (zh) * 2014-09-04 2016-03-09 中国石油化工股份有限公司 一种高精度的品质因子提取方法
CN104502965B (zh) * 2014-12-22 2017-04-05 中国石油天然气集团公司 一种振幅补偿因子的反演方法
CN105182416A (zh) * 2015-09-06 2015-12-23 中国石油天然气股份有限公司 基于分频数据的地震反演方法及其装置
CN106814393B (zh) * 2015-11-27 2019-07-12 中国石油化工股份有限公司 一种地层品质因子q的估算方法
CN107300718B (zh) * 2016-04-14 2019-11-12 中国石油天然气股份有限公司 一种品质因子三维衰减模型的建立方法
CN107340538B (zh) * 2016-05-03 2019-02-01 中国石油化工股份有限公司 基于混频处理的储层预测方法和装置
CN107544087B (zh) * 2016-06-23 2019-02-15 中国石油天然气股份有限公司 一种测量近地表地层品质因子的方法及装置
US10481287B2 (en) 2016-08-05 2019-11-19 Saudi Arabian Oil Company Surface consistent statics solution and amplification correction
US20200003923A1 (en) * 2016-12-01 2020-01-02 Arkady Yurievich Segal Method for determining physical characteristics of a homogeneous medium and its boundaries
CN107219553B (zh) * 2017-06-06 2019-11-08 中国石油化工股份有限公司 基于gr分频反演的暗河充填预测方法
CN107356964B (zh) * 2017-07-05 2018-10-30 西安交通大学 S变换域基于变分原理的q值估计与补偿方法
CN109425903A (zh) * 2017-08-21 2019-03-05 中国石油天然气股份有限公司 一种近地表地层品质因子的获取方法
CN108469633B (zh) * 2018-02-07 2019-10-11 中国石油天然气集团有限公司 一种地层品质因子的计算方法及装置
CN108845357B (zh) * 2018-06-13 2020-12-22 成都信息工程大学 一种基于同步挤压小波变换估计地层等效品质因子的方法
CN108919354B (zh) * 2018-09-27 2019-09-27 中国科学院地质与地球物理研究所 近地表q偏移方法及装置
CN109765615A (zh) * 2019-01-10 2019-05-17 中国石油天然气股份有限公司 一种地层品质因子反演方法及装置
CN110456415A (zh) * 2019-07-17 2019-11-15 中国石油大港油田勘探开发研究院 一种基于局部峰值频率的陆相薄储层解释方法及系统
CN110988986B (zh) * 2019-12-25 2021-01-01 成都理工大学 改善深层碳酸盐岩储层刻画精度的地震资料低频增强方法
CN113138419B (zh) * 2020-01-20 2022-05-10 中国石油天然气集团有限公司 提取下行子波和衰减参数的方法、装置
CN113341457A (zh) * 2020-02-18 2021-09-03 中国石油天然气集团有限公司 一种时频域等效q场的获取方法及装置
US11391855B2 (en) 2020-03-13 2022-07-19 Saudi Arabian Oil Company Developing a three-dimensional quality factor model of a subterranean formation based on vertical seismic profiles
CN111596350B (zh) * 2020-04-20 2023-03-24 江苏省地震局 一种地震台网波形数据质量监控方法和装置
US11703607B2 (en) 2020-06-15 2023-07-18 Saudi Arabian Oil Company Determining a seismic quality factor for subsurface formations from a seismic source to a first VSP downhole receiver
CN112099083B (zh) * 2020-08-26 2023-10-13 中化地质矿山总局地质研究院 一种基于双谱谱比对数的品质因子估计方法及系统
CN114152983B (zh) * 2020-09-08 2024-05-07 中国石油化工股份有限公司 地震子波提取方法、系统、存储介质以及电子设备
CN112099088B (zh) * 2020-09-16 2022-04-12 中油奥博(成都)科技有限公司 一种基于高密度光纤地震数据的油气指示及表征方法
US11859472B2 (en) 2021-03-22 2024-01-02 Saudi Arabian Oil Company Apparatus and method for milling openings in an uncemented blank pipe
US11573346B2 (en) 2021-04-15 2023-02-07 Saudi Arabian Oil Company Determining a seismic quality factor for subsurface formations for marine vertical seismic profiles
CN113589381B (zh) * 2021-08-09 2023-06-27 成都理工大学 一种基于压缩感知的相位与反射系数同时反演方法
CN113777650B (zh) * 2021-08-12 2022-10-25 西安交通大学 一种基于混合范数和小波变换的稀疏时频谱分解方法、装置、设备及存储介质
US11788377B2 (en) 2021-11-08 2023-10-17 Saudi Arabian Oil Company Downhole inflow control
CN114861563B (zh) * 2022-04-27 2022-12-13 中国石油大学(华东) 物理嵌入深度学习地层压力预测方法、装置、介质及设备
CN115902528B (zh) * 2023-02-21 2023-05-26 华东交通大学 一种直流牵引网振荡与短路故障辨识方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030086335A1 (en) * 2001-11-07 2003-05-08 Charles Naville Method for absolute preserved amplitude processing of seismic well data
WO2005019868A1 (en) * 2003-08-23 2005-03-03 Westerngeco Seismic Holdings Ltd. Multiple attenuation method
US20060265132A1 (en) * 2005-05-13 2006-11-23 Chevron U.S.A. Inc. Method for estimation of interval seismic quality factor
CN101046516A (zh) * 2007-04-22 2007-10-03 罗仁泽 利用垂直地震剖面和微测井进行地震信号补偿方法
CN101071175A (zh) * 2006-05-11 2007-11-14 中国石油集团东方地球物理勘探有限责任公司 零井源距垂直地震剖面纵横波数据深度域走廊叠加剖面处理方法
CN101630017A (zh) * 2008-07-16 2010-01-20 中国石油天然气集团公司 二维垂直地震剖面不同类型地震波场分离方法
CN102023311A (zh) * 2010-08-10 2011-04-20 中国石油大学(华东) 地层的品质因子谱及其求取方法
CN102269822A (zh) * 2010-06-02 2011-12-07 中国石油天然气集团公司 一种混合的地层吸收补偿方法

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5067112A (en) * 1991-01-04 1991-11-19 Mobil Oil Corporation Method for removing coherent noise from seismic data through f-x filtering
RU2153182C1 (ru) * 1998-12-30 2000-07-20 Московский государственный университет леса Способ оценки нефтегазового месторождения
GB0018480D0 (en) * 2000-07-27 2000-09-13 Geco Prakla Uk Ltd A method of processing surface seismic data
US20040122596A1 (en) * 2002-12-19 2004-06-24 Core Laboratories, Inc. Method for high frequency restoration of seismic data
US6931324B2 (en) * 2003-10-16 2005-08-16 Rdspi, L.P. Method for determining formation quality factor from seismic data
CN102884447B (zh) * 2010-05-05 2015-08-19 埃克森美孚上游研究公司 Q层析成像方法
US9291733B2 (en) * 2011-01-31 2016-03-22 Cggveritas Services Sa Device and method for determining S-wave attenuation in near-surface condition

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030086335A1 (en) * 2001-11-07 2003-05-08 Charles Naville Method for absolute preserved amplitude processing of seismic well data
WO2005019868A1 (en) * 2003-08-23 2005-03-03 Westerngeco Seismic Holdings Ltd. Multiple attenuation method
US20060265132A1 (en) * 2005-05-13 2006-11-23 Chevron U.S.A. Inc. Method for estimation of interval seismic quality factor
CN101071175A (zh) * 2006-05-11 2007-11-14 中国石油集团东方地球物理勘探有限责任公司 零井源距垂直地震剖面纵横波数据深度域走廊叠加剖面处理方法
CN101046516A (zh) * 2007-04-22 2007-10-03 罗仁泽 利用垂直地震剖面和微测井进行地震信号补偿方法
CN101630017A (zh) * 2008-07-16 2010-01-20 中国石油天然气集团公司 二维垂直地震剖面不同类型地震波场分离方法
CN102269822A (zh) * 2010-06-02 2011-12-07 中国石油天然气集团公司 一种混合的地层吸收补偿方法
CN102023311A (zh) * 2010-08-10 2011-04-20 中国石油大学(华东) 地层的品质因子谱及其求取方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
AN, HUI ET AL.: "distance varying inverse q filtering in vsp", COMPUTING TECHNIQUES FOR GIOPHYSICAL AND GEOCHENICAL EXPLORATION, vol. 20, no. 1, February 1998 (1998-02-01), pages 19 - 24, XP008175238 *
See also references of EP2837953A4 *
ZHU, GUANGMING ET AL.: "quality factor inversion using section modeling method", GEOPHYSICAL PROSPECTING FOR PETROLEUM, vol. 31, no. 3, September 1992 (1992-09-01), pages 55 - 65, XP008175236 *

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2016008105A1 (zh) * 2014-07-15 2016-01-21 杨顺伟 一种基于柯西分布的叠后波阻抗反演方法
RU2751088C2 (ru) * 2016-09-09 2021-07-08 Чайна Петролеум Энд Кемикл Корпорейшн Способ и система обработки сейсмических данных
CN107179544A (zh) * 2017-07-11 2017-09-19 中国石油集团川庆钻探工程有限公司地球物理勘探公司 转换波宽频延拓及提高转换波数据分辨率的方法
CN107179544B (zh) * 2017-07-11 2019-03-26 中国石油集团东方地球物理勘探有限责任公司 转换波宽频延拓及提高转换波数据分辨率的方法
CN109471162A (zh) * 2018-10-08 2019-03-15 中国石油天然气集团有限公司 层间多次波处理方法、系统、电子设备及可读介质
CN112305587A (zh) * 2019-08-02 2021-02-02 中国石油化工股份有限公司 恢复地震数据分辨率的方法、存储介质及计算机设备
CN112578436A (zh) * 2019-09-27 2021-03-30 中国石油化工股份有限公司 一种子波提取方法及系统
CN111965695A (zh) * 2020-08-17 2020-11-20 山西潞安环保能源开发股份有限公司五阳煤矿 一种基于反射槽波的小断层落差探测方法
CN112162314A (zh) * 2020-09-25 2021-01-01 武汉市工程科学技术研究院 一种人工地震信号剖面的二维插值方法
CN112162314B (zh) * 2020-09-25 2024-01-02 武汉市工程科学技术研究院 一种人工地震信号剖面的二维插值方法
CN116840916A (zh) * 2023-07-04 2023-10-03 成都理工大学 一种地震速度信号和加速度信号联合子波提取方法
CN116840916B (zh) * 2023-07-04 2024-03-26 成都理工大学 一种地震速度信号和加速度信号联合子波提取方法

Also Published As

Publication number Publication date
CN103376464B (zh) 2016-04-06
US20150168573A1 (en) 2015-06-18
RU2579164C1 (ru) 2016-04-10
EP2837953A4 (en) 2016-04-06
EP2837953A1 (en) 2015-02-18
CN103376464A (zh) 2013-10-30

Similar Documents

Publication Publication Date Title
WO2013152468A1 (zh) 一种地层品质因子反演方法
Edgar et al. How reliable is statistical wavelet estimation?
CN109669212B (zh) 地震数据处理方法、地层品质因子估算方法与装置
WO2004059342A1 (en) A method for high frequency restoration of seimic data
US7616524B1 (en) Wavelet based intercept attribute for seismic exploration
Wang et al. Improving the resolution of seismic traces based on the secondary time–frequency spectrum
Gao et al. Estimation of quality factor Q from the instantaneous frequency at the envelope peak of a seismic signal
CN110967749A (zh) 一种vsp地震资料频变q值估计与反q滤波方法
CN105092343B (zh) 去除薄层调谐效应的方法和识别预测薄储层及气层的方法
Neep et al. Measurement of seismic attenuation from high-resolution crosshole data
CN100412569C (zh) 利用地震微测井进行地震信号高频补偿方法
Cheng et al. Q estimation by a match-filter method
Lideng et al. Breakthrough and significance of technology on internal multiple recognition and suppression: a case study of Sinian Dengying Formation in Central Sichuan Basin, SW China
Zhang et al. Interval Q inversion based on zero-offset VSP data and applications
Herman et al. Analysis and removal of multiply scattered tube waves
CN110568491B (zh) 一种品质因子q的估算方法
CN111505707B (zh) 一种从垂直地震剖面数据中提取频散曲线的方法
Wu et al. Iterative deblending based on the modified singular spectrum analysis
Yang et al. Seismic attenuation estimation from instantaneous frequency
Zhang et al. Effect of Scholte wave on rotation of multi-component OBC seismic data in shallow water environment of the Arabian Gulf
Pramanik et al. Estimation of Q from borehole data and its application to enhance surface seismic resolution: A case study
Zhao et al. Absorption-Constrained Wavelet Power Spectrum Inversion for Enhancing Resolution of Nonstationary Seismic Data
Zhang et al. Improving OBC data quality of the geophone components in shallow-water Persian Gulf through advanced time-frequency analysis
Cao et al. Joint inversion method for interval quality factor based on amplitude and phase information
Qi Signal enhancement and impedance inversion for highly cyclically stratified sedimentation

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 12874053

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

WWE Wipo information: entry into national phase

Ref document number: 2012874053

Country of ref document: EP

ENP Entry into the national phase

Ref document number: 2014145635

Country of ref document: RU

Kind code of ref document: A

WWE Wipo information: entry into national phase

Ref document number: 14394100

Country of ref document: US