CN106597532B - Pre-stack seismic data frequency band expanding method combining well data and horizon data - Google Patents
Pre-stack seismic data frequency band expanding method combining well data and horizon data Download PDFInfo
- Publication number
- CN106597532B CN106597532B CN201610998901.3A CN201610998901A CN106597532B CN 106597532 B CN106597532 B CN 106597532B CN 201610998901 A CN201610998901 A CN 201610998901A CN 106597532 B CN106597532 B CN 106597532B
- Authority
- CN
- China
- Prior art keywords
- data
- frequency
- seismic
- frequency band
- well
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 52
- 238000012545 processing Methods 0.000 claims abstract description 27
- 238000010183 spectrum analysis Methods 0.000 claims abstract description 9
- 230000002146 bilateral effect Effects 0.000 claims abstract description 8
- 238000001228 spectrum Methods 0.000 claims description 44
- 230000004927 fusion Effects 0.000 claims description 13
- 238000007500 overflow downdraw method Methods 0.000 claims description 6
- 238000011161 development Methods 0.000 claims description 5
- 230000007704 transition Effects 0.000 claims description 5
- 230000002087 whitening effect Effects 0.000 claims description 5
- 238000001914 filtration Methods 0.000 description 12
- 230000008569 process Effects 0.000 description 6
- 238000003672 processing method Methods 0.000 description 5
- 238000013459 approach Methods 0.000 description 4
- 239000002131 composite material Substances 0.000 description 3
- 230000009466 transformation Effects 0.000 description 3
- 241000587240 Cynanchum Species 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 238000007405 data analysis Methods 0.000 description 2
- 238000000354 decomposition reaction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 230000003595 spectral effect Effects 0.000 description 2
- 241000692870 Inachis io Species 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 230000000149 penetrating effect Effects 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 238000011160 research Methods 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
-
- 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/32—Transforming one recording into another or one representation into another
-
- 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/61—Analysis by combining or comparing a seismic data set with other data
- G01V2210/616—Data from specific type of measurement
- G01V2210/6169—Data from specific type of measurement using well-logging
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
The invention provides a pre-stack seismic data frequency band expanding method combining well data and horizon data, which comprises the following steps: performing signal-to-noise ratio spectral analysis and spectral analysis sub-band scanning on the seismic data to determine a dominant frequency band and a frequency band to be compensated; correcting the logging data and matching well seisms to construct a well-side expected seismic record; constructing a sparse reflection coefficient profile by using the reflection coefficient at the well point, the well side channel seismic record and the horizon interpretation data, and further obtaining the expected seismic record at the non-well point; solving corresponding bilateral inverse wavelets by making convolution of the original seismic record and the inverse wavelets approximate to the expected seismic record; using the frequency broadening operators obtained from the superposition records on the CRP gathers, and processing the data of each CRP gather by using the frequency broadening operators to obtain high-frequency and low-frequency broadening records of the CRP gathers; and carrying out data reconstruction on the low-frequency data obtained by frequency extension, the original dominant frequency band data and the high-frequency data obtained by frequency extension to obtain broadband seismic data.
Description
Technical Field
The invention relates to the technical field of seismic data processing, in particular to a prestack seismic data frequency band expanding method for improving resolution by combining well data and horizon data.
Background
At present, the common processing methods for improving the resolution of seismic data are inverse Q filtering, generalized S transformation, spectrum simulation deconvolution and the like. Inverse Q filtering is the most common and effective method for improving seismic data resolution. At present, anti-Q filtering methods at home and abroad mainly comprise three types: inverse Q filtering methods using series expansion for approximate high frequency compensation, inverse Q filtering methods based on wave field continuation, and other inverse Q filtering methods. The inverse Q filtering requires the quality factor Q value of the formation, and its accuracy directly affects the accuracy of the inverse Q filtering. However, Q is difficult to be accurate and computational efficiency is low.
Deconvolution is a conventional method for improving seismic data resolution, and has become an important component of seismic data processing since the beginning of Robinson's deconvolution, which worked by combining Wiener's theory. Robinson and Treitel argue and Peacock and Treitel and extend this theory. This algorithm relies on a simple seismic trace and two reasonable assumptions-seismic reflection coefficients and seismic wavelets, respectively. The seismic traces are treated as convolution of the unknown wavelets and the unknown reflection coefficients. This problem is not solved without some a priori information. The white noise of the core wavelet minimum phase reflection coefficient of the Robinson deconvolution method solves this problem. First, the wavelet minimum phase means that once the amplitude spectrum of the wavelet is determined, the phase spectrum is determined. The statistical white reflection coefficient then means that his amplitude spectrum is close to constant, so the autocorrelation is a Dirac delta function. The white noise assumption allows us to estimate the amplitude spectrum of a wavelet from the amplitude spectrum of a seismic trace. The minimum phase assumption allows us to determine his phase spectrum from the amplitude spectrum of the sub-wave. When the wavelet is known that the true deconvolution process is only a numerical calculation problem. Although a solution can be obtained by applying Robinson deconvolution, the constant multiple of this solution still satisfies the condition.
In recent years deconvolution methods have been further optimized, including sparse pulse deconvolution, Gabor deconvolution, and the like. Compared with the traditional wavelet estimation method based on the conventional convolution model, the attenuation wavelet estimation based on the Gabor deconvolution theoretical model is more in line with the actual propagation process of wavelets, and the adopted seismic data does not need amplitude recovery. However, in general, deconvolution makes some assumptions about reflection coefficients and wavelets, which is not consistent with the reality, deconvolution improves resolution and reduces signal-to-noise ratio of data, and specifically, the conditions of the development process from Robinson deconvolution, minimum entropy deconvolution, finite entropy deconvolution, and spectral simulation deconvolution are all changed, but there is a common point that the process of improving resolution has no evidence of other data, cannot be combined with well data, and is affected more seriously by high-frequency signal-to-noise ratio of seismic data. The recently developed well control deconvolution technology is usually to approach well side channel seismic channels to synthetic seismic records, and the similarity coefficient of the current well seismic calibration is only about 0.7, which means that the seismic data is greatly improved while frequency expansion is carried out.
Sunpai et al applied a combined process of broadband processing and frequency broadening processing, combined with the thin sheath tuning theory, to emphasize thin layer reflections by enhancing the tuning energy. Ringning and the like use a deconvolution processing method, obtain a longitudinal loss rule of seismic records by enabling well-side seismic channels to approach to the broadband seismic data obtained by well logging after correction, obtain a stratum transverse filtering effect by utilizing the comparison of the well-side seismic channels and any transverse seismic channels, and connect the well-side seismic channels and the transverse seismic channels in series to achieve the purpose of improving the seismic data resolution. The frequency band of seismic data in well data which can be recovered is extracted by applying a matching tracking algorithm to the Cynanchum Liuhuanense and the Zhang-traditional Chang, and then the frequency bands of the seismic data are widened by establishing a mapping relation between the amplitude spectra of the Cynanchum Liuhuanense and the amplitude spectra of the Zhang-traditional Chang through Gaussian fitting.
Besides, there are various frequency expansion processing methods such as multi-wavelet decomposition and reconstruction method HFE, wavelet analysis frequency expansion reconstruction method, spectrum decomposition method, and tuning energy enhancement method and point spectrum method based on time-frequency analysis. However, these methods have problems of low resolution and complicated processing method to some extent.
Disclosure of Invention
The invention aims to provide a pre-stack seismic data frequency band expanding method combining well data and horizon data, which has high resolution and simple processing method, aiming at the defects in the prior art.
The object of the invention can be achieved by the following technical measures:
a pre-stack seismic data frequency band expanding method combining well data and horizon data comprises the following steps:
step 1: carrying out signal-to-noise ratio spectrum and frequency spectrum analysis on the seismic data, and dividing a frequency band range;
step 2: constructing a high-resolution expected seismic channel;
and step 3: solving high-frequency seismic data;
and 4, step 4: solving low-frequency seismic data;
and 5: fusing the high-frequency data and the low-frequency data, and reconstructing to obtain broadband seismic data;
step 6: and solving a frequency extension operator, and performing frequency extension processing under the signal-to-noise ratio spectrum constraint to obtain the broadband pre-stack seismic data.
The scheme is further refined as follows:
in step 1, a signal-to-noise ratio spectrum of seismic data is obtained by utilizing the cross correlation between seismic channels and the autocorrelation of each seismic channel, a frequency spectrum of the seismic data is obtained by utilizing Fourier transform, a proper threshold is selected by combining the characteristics of the two, and a seismic data frequency band is divided into a low-frequency band, a dominant frequency band, a transition frequency band and a frequency band to be promoted;
in the step 2, selecting broadband wavelets, and utilizing convolution theory or finite difference forward modeling to match logging data with well-side seismic channels after pre-stack time offset to construct high-resolution expected seismic channels;
in step 3, firstly, a sparse pulse method is utilized at a well point, the development condition of a thin layer in a target layer and the characteristic waveform characteristics of the reflection coefficient are combined, and a corresponding relation rule of the thin layer waveform and the reflection coefficient of the target layer is established, namely a time-varying threshold is set, the extreme value point of the waveform is larger than the threshold, and the value is used as the reflection coefficient; then, the obtained reflection coefficient synthetic seismic record is calculated according to the new record and the seismic channel at the point, and the reflection coefficient is corrected by keeping the root mean square amplitude in each sliding time window consistent; thirdly, the seismic record synthesized by the corrected reflection coefficients must ensure that the record in the dominant frequency band corresponds to the wave crest and the wave trough of the seismic channel beside the well, and the high-frequency record of the seismic record corresponds to the wave crest and the wave trough of the high-frequency seismic record made by the logging data;
in step 4, the low-frequency seismic data are obtained by interpolation of wave impedance along the horizon obtained by a logging curve;
the fusion reconstruction method in step 5 includes: a time domain fusion, a frequency domain fusion or a fusion method based on wavelet transform;
in step 6, a bilateral frequency broadening operator is calculated, so that the zero offset distance of the leveled common reflection point gather and the convolution least square of the leveled common reflection point gather are close to the zero offset distance seismic channel obtained by frequency broadening, and the pre-whitening percentage is 1% -5%.
Further:
in step 1, a suitable threshold is selected by using a signal-to-noise ratio spectrum greater than 1.
The logging data in step 2 comprises logging wave impedance and logging density curves.
In step 4, the frequency band of the low frequency seismic data is in the range of 0-10 Hz.
The frequency band expanding method of the pre-stack seismic data combining well data and horizon data determines dominant frequency bands and frequency bands to be compensated by performing signal-to-noise ratio spectral analysis and frequency spectrum analysis and sub-band scanning on seismic data; correcting the logging data and matching well seisms to construct a well-side expected seismic record; constructing a sparse reflection coefficient profile by using the reflection coefficient at the well point, the well side channel seismic record and the horizon interpretation data, and further obtaining the expected seismic record at the non-well point; solving corresponding bilateral inverse wavelets by making convolution of the original seismic record and the inverse wavelets approximate to the expected seismic record; using a frequency broadening operator obtained from the superposition records on a common reflection point gather (CRP gather for short), and processing data of each CRP gather by using the frequency broadening operator to obtain high-frequency and low-frequency broadening records of the CRP gather; and carrying out data reconstruction on the low-frequency data obtained by frequency extension, the original dominant frequency band data and the high-frequency data obtained by frequency extension to obtain broadband seismic data. The method has relatively simple processing flow, the processing result is more practical, and the resolution of data processing is greatly improved.
Drawings
FIG. 1 is a flow chart of the present invention;
FIG. 2 is a data volume before processing in a work area;
FIG. 3 is a signal-to-noise ratio spectrum of actual data in a work area;
FIG. 4 is a division of a work area signal band;
FIG. 5 is a graph of a well side channel and a log synthetic seismic record comparison for a sparse impulsive method reconstruction seismic channel in a work area;
FIG. 6 is the results obtained from sparse pulsing;
FIG. 7 is a seismic section after frequency broadening processing;
FIG. 8 is a composite prestack AVA trace gather;
FIG. 9 is a composite pre-stack AVA high frequency (left image) and low frequency record (right image) gather;
FIG. 10 is a frequency fused post-stack AVA trace gather;
FIG. 11 is a comparison of CRP gathers from and to a certain processing area (where the left panel is data before processing and the right panel is data after processing);
fig. 12 shows the signal-to-noise ratio spectrum change after the frequency extension processing.
Detailed Description
In order to make the aforementioned and other objects, features and advantages of the present invention comprehensible, preferred embodiments accompanied with figures are described in detail below.
In embodiment 1, a method for expanding frequency bands of pre-stack seismic data by combining well data and horizon data includes:
And 2, selecting broadband wavelets, and matching logging data (logging wave impedance and logging density curves) with the well-side seismic channels after pre-stack time offset by utilizing convolution theory or finite difference forward modeling to construct the high-resolution expected seismic channels.
And 3, obtaining high-frequency seismic data (the frequency band range is related to the conditions of a work area and the horizon depth) by a sparse reflection coefficient profile construction method, wherein the sparse reflection coefficient profile construction method is combined with logging data (logging wave impedance and logging density curves), geological horizon data and post-stack seismic waveform data to construct a sparse reflection coefficient model at the target horizon.
And 4, when low frequency (the frequency range is 0-10Hz) is absent, interpolating the low frequency seismic data along the horizon by the wave impedance obtained by the logging curve.
And 5, a data fusion method comprehensively considers the signal-to-noise ratio spectrum of the seismic data, researches the characteristics of the target depth of the horizon, the high-low frequency energy relation and the like, and performs fusion reconstruction on the high-frequency and low-frequency components to obtain the broadband seismic data under the condition of ensuring that the dominant frequency band is not changed.
And 6, solving a frequency broadening operator, applying a bilateral deconvolution method to enable zero offset seismic channels of the CRP gather (common reflection point gather) with pre-stack seismic flattening and operator convolution to approximate to obtain broadband seismic data, applying the frequency broadening operator to process data in the same CRP gather (common reflection point gather), and fusing the generated high-frequency data with the original data frequency, wherein the method can maintain a relative amplitude relation while expanding the frequency band.
step 1: carrying out signal-to-noise ratio spectrum and frequency spectrum analysis on the seismic data, and dividing a frequency band range;
step 2: constructing a high-resolution expected seismic channel;
and step 3: solving high-frequency seismic data;
and 4, step 4: solving low-frequency seismic data;
and 5: fusing the high-frequency data and the low-frequency data, and reconstructing to obtain broadband seismic data;
step 6: and solving a frequency extension operator, and performing frequency extension processing under the signal-to-noise ratio spectrum constraint to obtain the broadband pre-stack seismic data.
Further:
in step 1, a signal-to-noise ratio spectrum of seismic data is obtained by using cross correlation between seismic traces and self correlation of each seismic trace, a frequency spectrum of the seismic data is obtained by using fourier transform, a proper threshold (for example, the signal-to-noise ratio spectrum is greater than 1) is selected by combining characteristics of the two, and a seismic data frequency band is divided into a low frequency band, a dominant frequency band, a transition frequency band and a frequency band to be promoted (as shown in fig. 4).
In step 2, selecting broadband wavelets, and utilizing convolution theory or finite difference forward modeling to match logging data (logging wave impedance and logging density curve) with the well-side seismic channels after pre-stack time migration to construct high-resolution expected seismic channels.
In step 3, firstly, a sparse pulse method is utilized at a well point, the development condition of a thin layer in a target layer and the characteristic waveform characteristics of the reflection coefficient are combined, a corresponding relation rule of the thin layer waveform and the reflection coefficient of the target layer is established, namely a time-varying threshold is set, and a value which is greater than the threshold in an extreme point of the waveform is taken as the reflection coefficient; then, the obtained reflection coefficient synthetic seismic record is calculated according to the new record and the seismic channel at the point, and the reflection coefficient is corrected by keeping the root mean square amplitude in each sliding time window consistent; the modified reflection coefficient synthesized seismic record must then ensure that the record at the dominant band corresponds to the well-side seismic trace peak-trough, while its high frequency record corresponds to the high frequency seismic record peak-trough made from the log data (as shown in FIG. 5).
In step 4, low frequency seismic data (frequency band range of about 0-10Hz) are interpolated along the horizon from the wave impedance values obtained from the log.
In step 5, the invention adopts three different fusion reconstruction methods which are respectively as follows: time domain fusion, frequency domain fusion reconstruction or fusion method based on wavelet transformation. Wherein: the time domain fusion method is basically characterized by high operation efficiency, but cannot control frequency components of different frequency bands. The frequency domain fusion reconstruction is characterized in that the seismic data phase spectrum can be kept unchanged. The fusion method based on wavelet transformation is characterized in that the wavelet obtained by reconstruction is good in smoothness. FIG. 8 is a diagram of a synthetic prestack AVA gather; FIG. 9 is a composite pre-stack AVA high frequency record (left image) and low frequency record (right image) gather; FIG. 10 is a frequency fused post-stack AVA trace gather.
In step 6, a bilateral frequency broadening operator is calculated, so that the zero offset distance of the leveled common reflection point gather and the convolution least square of the leveled common reflection point gather are close to the zero offset distance seismic channel obtained by frequency broadening, and the pre-whitening percentage is 1% -5%. Compared with the conventional deconvolution, the method has the advantages that the operator obtained in the step (1) has the value after the 0 moment and the value before the 0 moment, and the approximation effect on the target is better. (2) In the same CRP gather, each gather gathers the same operator. (3) And (3) signal-to-noise ratio spectrum constraint, wherein the signal-to-noise ratio of the seismic data in the frequency region to be expanded is ensured to be basically unchanged before and after frequency expansion. Fig. 12 shows the signal-to-noise ratio spectrum after the frequency extension processing, and it can be seen that the signal-to-noise ratio after the processing is basically kept unchanged.
in step 1, work area data analysis is carried out, wherein the work area data analysis comprises a work area seismic data frequency band range, signal-to-noise ratio spectrum analysis, a signal purity spectrum and a frequency division section, and the work area data frequency band is divided into a low frequency band, a transition frequency band, a dominant frequency band, a transition frequency band and a high frequency area to be promoted.
And step 2, analyzing the acoustic wave and density logging curves to remove singular values so as to eliminate the influence of correcting the borehole diameter. The method mainly comprises the steps of mainly adopting a longitudinal wave velocity curve and a density curve in the logging data and a time-depth relation curve derived from other commercial software, carrying out median filtering, Walsh conversion and layer merging processing (the speed selection threshold of general layering parameters is 50-200 m/s, the layer thickness parameters are 10-200 m), correcting according to the matched time-depth relation to obtain a logging reflection coefficient, establishing an expected seismic channel at a well point, selecting a broadband wavelet, carrying out finite difference forward synthetic recording or convolution synthetic seismic recording, and carrying out time-varying amplitude balance processing on the recording.
In step 3, a proper threshold is selected at the well point according to the corresponding relation of the corrected well thin layer characteristics, reflection coefficient characteristics and waveforms, the value at the extreme point of the seismic waveform is used as the reflection coefficient, the value smaller than the threshold is set to zero, the synthesized seismic record corresponds to the logging synthesized seismic record waveform, and the rule is applied along the interpretation horizon to obtain a high-frequency reflection coefficient section.
And 4, filtering the reflection coefficient obtained by the logging data (acoustic logging curve and density curve) to obtain low-frequency information of the well point, and interpolating the low-frequency information along the interpretation horizon to obtain low-frequency data.
And 5, performing Fourier transform on the low-frequency seismic data, the dominant frequency band seismic data and the high-frequency seismic data, keeping the phase spectrum of the dominant frequency band seismic data unchanged, fusing the amplitude spectrum with the amplitude spectrum of the low-frequency and high-frequency seismic data, and performing inverse transform reconstruction to obtain fused seismic data.
And 6, calculating a bilateral frequency broadening operator to enable the zero offset distance and the convolution least square of the leveled prestack CRP gather to approach to the zero offset distance seismic trace obtained by frequency broadening, wherein the dereferencing range of the pre-whitening parameter is 1% -5%, applying the operator to other seismic traces in the CRP gather, and performing frequency fusion to obtain the prestack seismic data subjected to frequency broadening.
(1) the method comprises the steps of obtaining a signal-to-noise ratio spectrum, a signal purity spectrum and a frequency division scanning section of a target work area, researching signal characteristics of the target work area, dividing an intermediate frequency band with high signal-to-noise ratio and large amplitude into a seismic data dominant frequency band, selecting a frequency band with weak effective signals appearing on the frequency division section but obviously lower than the dominant frequency band on the signal-to-noise ratio spectrum as a high frequency region to be compensated, and taking a low frequency region missing on a frequency spectrum as a low frequency region to be compensated (fig. 2 is original data, fig. 3 is signal-to-noise ratio spectrum distribution of the data and fig. 4 is a divided frequency band range).
(2) And analyzing the acoustic wave and density logging curves to remove singular values so as to eliminate the influence of correcting the borehole diameter. In the logging data, mainly adopting a sound wave time difference curve and a density curve in the data and a time-depth relation curve derived from other commercial software, carrying out median filtering and layer position merging processing (the speed of general layering parameters is selected to be 50-200 m/s, and the parameters after the layer are 10-200 m), and correcting according to the well-matched time-depth relation to obtain a logging reflection coefficient.
(3) An improved sparse pulse method is utilized at a well point, the thin layer development condition in the layer and the characteristic waveform characteristics of the reflection coefficient are combined, a corresponding relation rule of the thin layer waveform and the reflection coefficient in the layer is established, namely a time-varying threshold is set, and the value which is larger than the threshold in the extreme point of the waveform is used as the reflection coefficient (as shown in fig. 6, the established sparse pulse record is obtained).
(4) And selecting continuous horizons penetrating through a plurality of wells, and interpolating the low frequency of the corrected logging wave impedance in space according to the trend of layers and the relationship between the layers (as shown in figure 5).
(5) The method comprehensively considers the characteristics of a signal-to-noise ratio spectrum, a target depth, a high-low frequency energy relation and the like, fuses high-frequency and low-frequency components to reconstruct the seismic data under the condition of ensuring that the dominant frequency band is not changed, and mainly applies reconstruction methods such as time domain, frequency domain and wavelet transform-based reconstruction methods.
(6) And calculating a bilateral frequency broadening operator to enable the zero offset distance and the convolution least square of the zero offset distance of the flattened pre-stack CRP gather (common reflection point gather) to approach the zero offset distance seismic trace obtained by frequency broadening, wherein the pre-whitening percentage is 1% -5%. And (4) performing convolution (other fields are called convolution) on the obtained operator and the actual seismic trace prestack CRP data to obtain a high-frequency or low-frequency band seismic data volume (shown as data after frequency extension in FIG. 7). (for example, fig. 8 is an original common reflection point gather of the theory, the left side of fig. 9 is high-frequency data obtained by frequency extension of fig. 8, the right side of fig. 9 is low-frequency data obtained by frequency extension of fig. 8, fig. 10 is data after fusion, the left side of fig. 11 is an actual CRP gather (common reflection point gather) before processing, the right side of fig. 11 is a CRP gather after processing, it can be seen that a frequency band is extended, fig. 12 is signal-to-noise ratio spectral distribution before and after processing, it can be seen that the signal-to-noise ratio is not reduced after processing.
Claims (4)
1. A pre-stack seismic data frequency band expanding method combining well data and horizon data is characterized by comprising the following steps:
step 1: carrying out signal-to-noise ratio spectrum and frequency spectrum analysis on the seismic data, and dividing a frequency band range;
step 2: constructing a high-resolution expected seismic channel;
and step 3: solving high-frequency seismic data;
and 4, step 4: solving low-frequency seismic data;
and 5: fusing the high-frequency data and the low-frequency data, and reconstructing to obtain broadband seismic data;
step 6: obtaining a frequency extension operator, and performing frequency extension processing under the signal-to-noise ratio spectrum constraint to obtain broadband pre-stack seismic data;
in step 1, a signal-to-noise ratio spectrum of seismic data is obtained by utilizing the cross correlation between seismic channels and the autocorrelation of each seismic channel, a frequency spectrum of the seismic data is obtained by utilizing Fourier transform, a proper threshold is selected by combining the characteristics of the two, and a seismic data frequency band is divided into a low-frequency band, a dominant frequency band, a transition frequency band and a frequency band to be promoted;
in the step 2, selecting broadband wavelets, and utilizing convolution theory or finite difference forward modeling to match logging data with well-side seismic channels after pre-stack time offset to construct high-resolution expected seismic channels;
in step 3, firstly, a sparse pulse method is utilized at a well point, the development condition of a thin layer in a target layer and the characteristic waveform characteristics of the reflection coefficient are combined, a corresponding relation rule of the thin layer waveform and the reflection coefficient of the target layer is established, namely a time-varying threshold is set, and a value which is greater than the threshold in an extreme point of the waveform is taken as the reflection coefficient; secondly, the obtained reflection coefficient synthetic seismic record is calculated according to the new record and the seismic channel at the point, and the reflection coefficient is corrected by keeping the root mean square amplitude in each sliding time window consistent; thirdly, the seismic record synthesized by the corrected reflection coefficients must ensure that the record in the dominant frequency band corresponds to the wave crest and the wave trough of the seismic channel beside the well, and the high-frequency record of the seismic record corresponds to the wave crest and the wave trough of the high-frequency seismic record made by the logging data;
in step 4, the low-frequency seismic data are obtained by interpolation of wave impedance along the horizon obtained by a logging curve;
the fusion reconstruction method in step 5 includes: a time domain fusion, a frequency domain fusion or a fusion method based on wavelet transform;
in step 6, a bilateral frequency broadening operator is calculated, so that the zero offset distance of the leveled common reflection point gather and the convolution least square of the leveled common reflection point gather are close to the zero offset distance seismic trace obtained by frequency broadening, and the pre-whitening percentage is 1% -5%.
2. The method of claim 1, wherein in step 1, the selection of the appropriate threshold is performed by using a signal-to-noise ratio spectrum greater than 1.
3. The method of claim 1 or 2, wherein the logging data in step 2 comprises log wave impedance and log density curves.
4. The method of claim 1 or 2, wherein in step 4, the frequency band of the low frequency seismic data is in the range of 0-10 Hz.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610998901.3A CN106597532B (en) | 2016-11-14 | 2016-11-14 | Pre-stack seismic data frequency band expanding method combining well data and horizon data |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610998901.3A CN106597532B (en) | 2016-11-14 | 2016-11-14 | Pre-stack seismic data frequency band expanding method combining well data and horizon data |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106597532A CN106597532A (en) | 2017-04-26 |
| CN106597532B true CN106597532B (en) | 2020-06-30 |
Family
ID=58590163
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201610998901.3A Active CN106597532B (en) | 2016-11-14 | 2016-11-14 | Pre-stack seismic data frequency band expanding method combining well data and horizon data |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN106597532B (en) |
Families Citing this family (27)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107489418B (en) * | 2017-08-14 | 2020-06-09 | 中国石油天然气股份有限公司 | Density logging forward modeling method and device |
| CN109655883A (en) * | 2017-10-10 | 2019-04-19 | 中国石油化工股份有限公司 | A kind of earthquake dividing method and system for target |
| CN108363739B (en) * | 2018-01-22 | 2021-07-30 | 中国海洋大学 | Seismic data high-low frequency expansion method based on sparse acquisition |
| CN110967733B (en) * | 2018-09-28 | 2022-03-08 | 中国石油化工股份有限公司 | Improved seismic source reverse time positioning method and computer readable storage medium |
| CN109358365A (en) * | 2018-10-16 | 2019-02-19 | 中国海洋石油集团有限公司 | A kind of thin layer recognition methods based on well shake information matches |
| CN111077578B (en) * | 2018-10-22 | 2022-05-10 | 中国石油天然气股份有限公司 | Rock stratum distribution prediction method and device |
| CN111610562B (en) * | 2019-02-22 | 2023-10-27 | 中国石油天然气集团有限公司 | Method and device for separating strong reflection signals in seismic data |
| CN111665542B (en) * | 2019-03-06 | 2023-04-25 | 中国石油天然气集团有限公司 | Seismic data frequency extension method and system |
| CN111830559B (en) * | 2019-04-19 | 2023-01-10 | 中国石油天然气股份有限公司 | Method and device for frequency extension processing of seismic data |
| CN112034510B (en) * | 2019-06-04 | 2023-04-07 | 中国石油天然气股份有限公司 | Low-frequency compensation method and low-frequency compensation device |
| CN112764094B (en) * | 2019-10-21 | 2023-10-31 | 中国石油天然气股份有限公司 | Inversion method and device for seismic time-frequency reflection coefficient |
| CN112711067B (en) * | 2019-10-24 | 2023-10-31 | 中国石油天然气股份有限公司 | Thin reservoir prediction method and device |
| CN112925019A (en) * | 2019-12-06 | 2021-06-08 | 中国石油天然气股份有限公司 | Method and device for identifying pore type dolomite |
| CN112946750B (en) * | 2019-12-11 | 2024-01-30 | 中国石油天然气股份有限公司 | Well shock calibration method and system |
| CN111060960B (en) * | 2019-12-27 | 2022-03-18 | 恒泰艾普（北京）能源科技研究院有限公司 | FWI modeling method based on synthetic gun records |
| CN111242198B (en) * | 2020-01-06 | 2023-05-09 | 中国石油化工股份有限公司 | Gather tuning optimization dimension increasing method based on sparse representation |
| CN113391353B (en) * | 2020-03-11 | 2024-06-25 | 中国石油天然气集团有限公司 | Seismic data processing method and device |
| CN113625361B (en) * | 2020-05-09 | 2024-05-07 | 中国石油化工股份有限公司 | Well constraint-based seismic data quantitative frequency expansion method and device |
| CN113721294B (en) * | 2020-05-26 | 2023-06-23 | 中国石油化工股份有限公司 | Complex domain least square constraint spectrum bluing frequency-expanding method |
| CN113960671B (en) * | 2020-07-20 | 2024-04-16 | 中国石油化工股份有限公司 | Method, apparatus, computer device and storage medium for compressing frequency-variable wavelets |
| CN111765960B (en) * | 2020-07-23 | 2022-02-25 | 国网山西省电力公司太原供电公司 | Method for extracting seismic signals of OPGW (optical fiber composite overhead ground wire) optical cable based on distributed optical fiber sensing |
| CN112305623B (en) * | 2020-08-20 | 2022-06-17 | 中国地质科学院地球物理地球化学勘查研究所 | Potential field feature acquisition method and device based on spectrum fusion |
| CN112213775B (en) * | 2020-09-16 | 2023-01-24 | 中国石油天然气股份有限公司 | Fidelity frequency-boosting method for high-coverage-frequency pre-stack seismic data |
| CN112764099B (en) * | 2020-12-23 | 2023-03-31 | 中国海洋石油集团有限公司 | Seismic data frequency extension method based on seismic geometry information |
| CN112882099B (en) * | 2021-02-25 | 2024-01-12 | 中海石油深海开发有限公司 | Earthquake frequency band widening method and device, medium and electronic equipment |
| CN113341463B (en) * | 2021-06-10 | 2023-05-26 | 中国石油大学(北京) | Non-stationary blind deconvolution method for pre-stack seismic data and related components |
| CN115639603B (en) * | 2022-10-25 | 2024-02-23 | 电子科技大学 | Seismic frequency expansion method, device and storage medium based on superposition and interference removal of sampling points |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE3307575A1 (en) * | 1983-03-03 | 1984-09-06 | Erich Dipl.-Geophys. 8372 Lindberg Lippmann | Electrodynamic seismometer |
| CN103675899A (en) * | 2012-09-04 | 2014-03-26 | 中国石油天然气集团公司 | Method for expanding post-stack seismic data frequency band based on wavelet compression |
| CN105259580A (en) * | 2015-10-30 | 2016-01-20 | 中国石油大学(北京) | Low-frequency extension method for vibroseis signal |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9170343B2 (en) * | 2012-04-30 | 2015-10-27 | Conocophillips Company | Quasi-impulsive displacement source |
-
2016
- 2016-11-14 CN CN201610998901.3A patent/CN106597532B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE3307575A1 (en) * | 1983-03-03 | 1984-09-06 | Erich Dipl.-Geophys. 8372 Lindberg Lippmann | Electrodynamic seismometer |
| CN103675899A (en) * | 2012-09-04 | 2014-03-26 | 中国石油天然气集团公司 | Method for expanding post-stack seismic data frequency band based on wavelet compression |
| CN105259580A (en) * | 2015-10-30 | 2016-01-20 | 中国石油大学(北京) | Low-frequency extension method for vibroseis signal |
Non-Patent Citations (1)
| Title |
|---|
| "宽频带地震数据融合重构方法研究";李健 等;《中国地球科学联合学术年会2016》;20161015;第1517-1520页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN106597532A (en) | 2017-04-26 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106597532B (en) | Pre-stack seismic data frequency band expanding method combining well data and horizon data | |
| CN109061764B (en) | Frequency-division fusion wave impedance inversion method | |
| CN105353407B (en) | Post-stack seismic wave impedance inversion method | |
| US20150168573A1 (en) | Geologic quality factor inversion method | |
| CN104849756B (en) | A kind of seismic data resolution that improves strengthens the method for effective weak signal energy | |
| US11435493B2 (en) | Enhanced waveform analysis for target modes of borehole waves | |
| CN108845357B (en) | Method for estimating formation equivalent quality factor based on synchronous extrusion wavelet transform | |
| CN107678064B (en) | Real-time extraction method for sound wave time difference | |
| CN110596758B (en) | Seismic signal low-frequency energy compensation method | |
| CN111522062B (en) | Underburden amplitude compensation method based on volcanic shielding quantitative analysis | |
| CN109946740B (en) | Seismic resolution enhancement method based on wide flat spectrum seismic wavelet shaping | |
| CN106556865A (en) | A kind of tandem type seismic signal optimizes time-frequency conversion method | |
| CN112255690B (en) | Self-adaptive surrounding rock strong reflection separation method based on seismic phase decomposition | |
| CN113391353B (en) | Seismic data processing method and device | |
| CN114371505A (en) | Multi-wavelet inversion method and system based on seismic frequency division technology | |
| CN110673211B (en) | Quality factor modeling method based on logging and seismic data | |
| CN111427080B (en) | Method for extracting space-variant generalized wavelets of seismic data | |
| CN114137606A (en) | Stable spectrum simulation deconvolution method | |
| CN113933902B (en) | High-frequency expansion method based on stereoscopic space wavelet | |
| CN113589366B (en) | Broadband fusion modeling method based on full waveform inversion | |
| CN111665536A (en) | Well depth design method based on micro logging wavelet quantification constraint | |
| CN113419275B (en) | High-resolution seismic processing method based on sparse dictionary learning | |
| CN117368978A (en) | Thin interbed target processing and identifying method based on Q compensation compressed sensing | |
| CN117075195B (en) | Coal seam suppression method based on prestack gather frequency division | |
| CN113359187B (en) | Wavelet sidelobe elimination method for seismic data |
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 |