CN109375265A - A kind of ideal earthquake spectral factorization method based on changeable phases Ricker wavelet match tracing - Google Patents
A kind of ideal earthquake spectral factorization method based on changeable phases Ricker wavelet match tracing Download PDFInfo
- Publication number
- CN109375265A CN109375265A CN201810962986.9A CN201810962986A CN109375265A CN 109375265 A CN109375265 A CN 109375265A CN 201810962986 A CN201810962986 A CN 201810962986A CN 109375265 A CN109375265 A CN 109375265A
- Authority
- CN
- China
- Prior art keywords
- wavelet
- phase
- spectrum
- signal
- frequency
- 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 33
- 230000003595 spectral effect Effects 0.000 title abstract 2
- 238000001228 spectrum Methods 0.000 claims abstract description 95
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 44
- 238000004458 analytical method Methods 0.000 claims description 23
- 238000004364 calculation method Methods 0.000 claims description 9
- 238000010606 normalization Methods 0.000 claims description 2
- 230000009466 transformation Effects 0.000 abstract description 3
- 230000002123 temporal effect Effects 0.000 abstract 1
- 239000006185 dispersion Substances 0.000 description 3
- 238000001514 detection method Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 description 1
- 208000035126 Facies Diseases 0.000 description 1
- 241000446313 Lamella Species 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000001427 coherent effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 238000011065 in-situ storage Methods 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/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/32—Transforming one recording into another or one representation into another
- G01V1/325—Transforming 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/40—Transforming data representation
- G01V2210/48—Other transforms
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 discloses a kind of ideal earthquake spectral factorization methods based on changeable phases Ricker wavelet match tracing, first building changeable phases Ricker wavelet library；Then matching wavelet and corresponding peak factor are sought according to match tracing method, and the subsignal obtained as decomposition that peak factor is multiplied with wavelet atom；Then Fourier transformation is done to the corresponding subsignal in single reflection interface, calculates the amplitude spectrum and rotatable phase spectrum of the interface reflection seismic waves, Fourier transform spectrum was placed at corresponding delay time, ideal time-frequency spectrum decomposition result is finally obtained.The Time-frequency Decomposition result includes amplitude spectrum and rotatable phase composes two parts, and rotatable phase spectrum is capable of providing more abundant geology and reservoir information.Amplitude spectrum decomposition result accurately reflects delay time and its energy distribution varying with frequency that seismic wavelet reaches interface, it is not limited by the uncertainty principle limit, it is a kind of ideal time-frequency spectrum for seismic signal with high temporal resolution and frequency resolution.
Description
Technical Field
The invention belongs to the field of oil and gas geophysical exploration, and particularly relates to an ideal earth spectrum decomposition method based on phase-varying Rake wavelet matching pursuit.
Background
The seismic spectrum decomposition technology converts time domain seismic signals into a time-frequency domain for analysis, can acquire richer stratum information from seismic data, and is widely applied to the aspects of determining stratum thickness, sedimentary facies interpretation, oil gas detection and the like.
Common seismic spectrum decomposition methods mainly comprise linear time-frequency analysis, bilinear time-frequency analysis and self-adaptive time-frequency analysis methods. A linear time-frequency analysis method, namely a kernel function decomposition method, realizes the time-frequency locality of signals by using a window function, and mainly comprises short-time Fourier transform (STFT), Continuous Wavelet Transform (CWT), S Transform (ST) and Generalized S Transform (GST). The linear time-frequency analysis method is limited by a Heisenberg uncertainty principle, the time resolution and the frequency resolution of the linear time-frequency analysis method are mutually restricted and contained, and the time resolution and the frequency resolution cannot be simultaneously reduced to any degree. The distribution of signal energy in a time-frequency domain is obtained by the bilinear time-frequency analysis method, and mainly comprises Cohen time-frequency distribution and affine time-frequency distribution. The common Wigner-Ville distribution (WVD) is a Cohen-like time frequency distribution, and compared with a linear time frequency analysis method, the WVD takes a self signal as a window, so that the self-adaptive matching is realized, the limitation of an inaccurate principle is avoided, and the time frequency resolution is higher. Cross terms exist in the WVD decomposition, serious interference is generated on a frequency spectrum decomposition result, and the cross interference terms can be reduced by applying a smooth function to carry out two-dimensional low-pass filtering on the WVD. The higher the smoothing degree of the time-frequency domain, the smaller the influence of the interference term, and the lower the resolution of the time-frequency spectrum, and a trade-off between the joint time-frequency resolution and the influence of the coherent term is needed.
A seismic spectrum decomposition method based on Matching Pursuit (MP) is a self-adaptive time-frequency analysis method, firstly an over-complete wavelet base is constructed, then sparse decomposition of signals is realized based on matching pursuit, and finally time-frequency distribution of each decomposed sub-component signal is obtained by utilizing WVD, and superposition and summation are carried out. The seismic spectrum decomposition result based on matching pursuit has high time-frequency resolution, and avoids the influence of cross terms, thereby being widely applied to the field of seismic data interpretation and inversion. Wavelet atoms commonly used in matching pursuit algorithms include Gabor wavelets, Morlet wavelets, and rake wavelets. The matching tracking algorithm based on the Rake wavelets usually selects the Rake wavelets with zero phase to establish a wavelet base, and calculates the amplitude spectrum and the phase spectrum of the matched wavelets as a whole, wherein the phase spectrum obtained by the method is the linear phase change caused by time delay. However, the zero-phase rake wavelet base has poor redundancy compared to the actual seismic wavelets, and phase-varying rake wavelets need to be constructed. Meanwhile, in the actual seismic recording, not only linear phase change caused by time delay exists, but also phase rotation caused by lamella, velocity dispersion, energy attenuation and the like exists, and conventional matching tracking (no matter zero-phase or non-zero-phase Rake wavelets, real domain or complex domain matching tracking) does not separately consider the rotating phase spectrum.
Disclosure of Invention
In order to overcome the defects of the prior art, the invention provides an ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit, which not only has high time-frequency resolution, but also can obtain a rotating phase spectrum, and the rotating phase spectrum contains richer geological information and is expected to play an important role in reservoir prediction and oil gas detection.
The technical scheme adopted by the invention for solving the technical problems is as follows: an ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit is constructed, and the method comprises the following steps:
(1) inputting a path of seismic data as an initial residual signal;
(2) constructing a complex analysis signal of a current residual signal, then calculating an instantaneous envelope of the complex analysis signal, and determining the delay time of a matched wavelet atom according to the instantaneous envelope;
(3) constructing a zero-phase Rake wavelet according to the delay time of the matched wavelet atoms, scanning the dominant frequency and the phase of the zero-phase Rake wavelet, and constructing a phase-variable Rake wavelet base;
(4) extracting a wavelet atom which is most matched with a current residual signal from a phase-changing Rake wavelet base, solving an amplitude coefficient corresponding to the most matched wavelet atom, multiplying the amplitude coefficient and the most matched wavelet atom to obtain a sub-signal obtained by decomposition, then subtracting the sub-signal obtained by decomposition from the current residual signal, and taking the obtained residual as a new residual signal;
(5) repeating the steps (2) to (4) until the iteration number reaches a preset value or the residual energy is less than a certain threshold value;
(6) after iteration is finished, obtaining a plurality of phase-variable Rake matching wavelets and corresponding amplitude coefficients, and respectively and correspondingly multiplying the matching wavelets and the amplitude coefficients to obtain single-interface reflection seismic signals serving as a sub-signal set;
(7) and performing Fourier transform on each sub-signal in the sub-signal set, and corresponding the obtained Fourier transform spectrum to the delay time corresponding to the interface to finally obtain an ideal time-frequency spectrum decomposition result.
Further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (2) specifically comprises:
constructing a complex analysis signal of the residual signal s (t):
S(t)＝s(t)+isH(t)，
wherein s isH(t) represents the Hilbert transform of the real signal s (t), the instantaneous amplitude envelope of the complex analysis signal being:
and taking the time corresponding to the maximum value of the instantaneous amplitude envelope as the delay time of the matched wavelet atom.
Further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (3) specifically comprises:
(31) constructing a zero-phase Rake wavelet according to the delay time of the matched wavelet atoms:
where τ is the delay time of the matched wavelet atom, fdDominant frequency of the Rake wavelet;
(32) by scanning the time-frequency parameter gamma ═ fdPhi, constructing a phase-variable Rake wavelet atom library, wherein the scanning range of the wavelet main frequency is 0 to Nyquist frequency, the scanning range of the phase is 0-180 degrees, and the constructed phase-variable Rake wavelet library is as follows:
wγ(t)＝w(t-τ,fd)cos(φ)-wH(t-τ,fd)sin(φ)，
where φ is the phase of the Rake wavelet, wH(t) is the Hilbert transform of the zero-phase Rake wavelet w (t).
Further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (4) specifically comprises:
for the first random iteration number k, all wavelet atoms in the phase-varying Rake wavelet base and the current residual signal sres k-1(t) performing correlation calculation to obtain a correlation coefficient corr:
wherein,representing the current residual signal and the wavelet atom wγ(t) the inner product of the (t),representing wavelet atoms wγ(t) a signal normalization operator;
then, the maximum value of the correlation coefficient and the wavelet atom w corresponding to the maximum value are obtainedγk(t)，
The wavelet atom wγk(t) is the wavelet atom which is most matched with the current residual signal, and then the wavelet atom w is obtained through calculationγk(t) amplitude coefficient C corresponding tok：
The sub-signal decomposed by the kth iteration is Ckwγk(t) subtracting the k-th decomposed sub-signal from the current residual to obtain a new residual signal sres k(t)：
sres k(t)＝sres k-1(t)-Ckwγk(t)。
Further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (6) specifically comprises:
will sk(t)＝Ckwγk(t) as the sub-signals generated in the kth iteration, all sub-signals are formed into a sub-signal set.
Further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (7) specifically comprises:
evaluating a subsignal sk(t) Fourier transform spectrum Sk(f)：
Wherein f represents the frequency of the fourier transform;
then Fourier transform spectrum Sk(f) Corresponding to the delay time corresponding to the interface to obtain the ideal time spectrum IDFTs(t,f)：
Obtaining the amplitude spectrum and the phase spectrum of the reflected seismic waves of the interface,
further, in the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit, the ideal time frequency spectrum IDFTsThe amplitude spectrum a (t, f) and the phase spectrum Φ (t, f) of (t, f) are:
further, in the seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit of the invention, the step (7) further comprises obtaining an ideal time frequency spectrum IDFTsRotational phase spectrum Θ (t, f) of (t, f):
wherein phik(f) Is a propagation time of τkBest match the phase of the rake wavelet.
Compared with the prior art, the invention has the beneficial effects that:
the amplitude spectrum is formed by Fourier transform of a single-interface fixed propagation time sub-signal, the energy of a reflection wavelet of a single interface is reflected in time, and matching pursuit is a sparse decomposition method, so that the amplitude spectrum has ideal time resolution; in terms of frequency, the Fourier transform is directly carried out on the sub-signals with the same length as the seismic record length, the frequency resolution cannot be reduced due to the influence of the window length, and therefore the ideal frequency resolution is achieved.
The phase spectrum contains two types of phase information, the first is linear phase change due to time delay, and the second is phase rotation due to thin layer, velocity dispersion and energy attenuation. The phase spectrum obtained by the conventional matching tracking based on the zero-phase Rake wavelet only reflects the linear phase change caused by time delay, and the phase spectrum obtained by the matching tracking based on the non-zero-phase Rake wavelet reflects the comprehensive effect of two influences of the linear phase change and the phase rotation. The linear phase change caused by the time delay has no direct relation with the lithology and the reservoir, so for the conventional time frequency analysis, only amplitude spectrum information is generally utilized, and phase spectrum information is less utilized. In fact, in seismic exploration, the rotating phase spectrum can play an important role in formation identification and reservoir prediction because the rotating phase contains information of thin layers, velocity dispersion and energy attenuation.
Drawings
The invention will be further described with reference to the accompanying drawings and examples, in which:
FIG. 1 is a flow chart of an ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit according to the present invention;
FIG. 2 is a seismic trace of the model of example 1;
FIG. 3 is a sub-signal set decomposed by the model seismic trace based on the matching pursuit of the phase-varying Rake wavelets of example 1;
FIG. 4(a) is an amplitude spectrum of an ideal seismic spectrum decomposition result obtained by using the ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit according to the invention for a model seismic trace of example 1;
FIG. 4(b) is a graphical representation of an ideal amplitude spectrum;
FIG. 4(c) is an ideal rotating phase spectrum;
FIGS. 5(a) and 5(b) are an amplitude spectrum and a phase spectrum obtained by the S-transform time-frequency analysis method for the model seismic traces of example 1;
FIG. 6 shows the result of the matching pursuit-WVD time-frequency decomposition of the model seismic traces of example 1.
Detailed Description
For a more clear understanding of the technical features, objects and effects of the present invention, embodiments of the present invention will now be described in detail with reference to the accompanying drawings.
Referring to fig. 1, fig. 1 is a flowchart of an embodiment of an ideal seismic spectrum decomposition method based on phase-varying raeket wavelet matching pursuit of the present invention, and the method includes the following steps:
and S1, inputting single-channel seismic data as an initial residual signal.
S2, constructing a complex analysis signal of the current residual signal S (t):
S(t)＝s(t)+isH(t)
wherein s isH(t) represents the Hilbert transform of the real signal s (t), then the instantaneous amplitude A (t) of the complex analysis signal is:
and taking the time corresponding to the maximum value of the instantaneous amplitude envelope as the delay time tau of the matched wavelet atom.
S3, constructing a zero-phase Rake wavelet according to the delay time tau of the matched wavelet atom calculated in the step S2
w(t)＝{1-2[πfd(t-τ)]2}exp{-[πfd(t-τ)]2}
Wherein, the delay time tau of the matching wavelet atom is the propagation time of the Rake wavelet in the underground medium, fdTo match the dominant frequencies of the rake wavelet. To obtain a Rake wavelet with any constant phase, an analytic signal expression of a Rake wavelet with a zero phase is obtained
W(t)＝w(t)+iwH(t)
Wherein wH(t) is the Hilbert transform of a zero-phase Rake wavelet. After rotating its phase by phi, the new non-zero phase complex field Rake wavelet is
W'(t)＝W(t)exp(iφ)
Namely, it is
W'(t)＝[w(t)+iwH(t)][cos(φ)+isin(φ)]
Then the time-frequency atoms of the constant-phase Rake wavelet can be represented as
wγ(t)＝w(t-τ,fd)cos(φ)-wH(t-τ,fd)sin(φ)
If the sampling frequency of the signal is fsThe maximum frequency of the atoms of the Rake wavelet is f according to the Nyquist's laws2, in practice, the primary frequency of seismic waves is usually within 100 Hz; the phase interval of the Rake wavelets is 0-180 degrees. Scanning time-frequency parameter gamma ═ tau, fdPhi, then obtaining the Rake wavelet base. The smaller the scanning interval, the more wavelet atoms are generated, the more sparse the decomposition result is, but the calculation amount required by time-frequency decomposition is increased.
In the invention, the phase of the Rake wavelet is assumed to be phi, and the time-frequency parameter gamma is changed to { f through scanningdPhi scans the dominant frequency and phase of the zero-phase Rake wavelet to construct a phase-variable Rake wavelet atom library, wherein the scanning range of the dominant frequency of the wavelet is 0 to Nyquist frequency, the scanning range of the phase is 0 DEG to 180 DEG, and the constructed phase-variable Rake wavelet library wγ(t) is:
wγ(t)＝w(t-τ,fd)cos(φ)-wH(t-τ,fd)sin(φ)
wherein wH(t) is the Hilbert transform of the zero-phase Rake wavelet w (t). Generally, if the sampling frequency of the signal is fsThe frequency sweep interval may be the minimum frequency interval in the fourier transform, i.e., △ f-fsand/N, wherein N is the data length. The scan interval of the phase angle may take pi/36 radians. The smaller the scanning interval, the more wavelet atoms are generated, the more sparse the decomposition result is obtained, but the more time is required for calculation, and a compromise needs to be made between the calculation accuracy and the calculation time.
S4, combining all wavelets in the Rake wavelet base with the current residual signal Sres k-1(t) performing a correlation calculation, the correlation thereofCoefficient of
Calculating maximum value of correlation coefficient and corresponding wavelet atom
Where D is a library of matching Rake wavelet atoms, the wavelet atom wγk(t) is the wavelet atom that matches best with the current residual signal, and its corresponding amplitude coefficient is:
the signal decomposed by the kth iteration isSubtracting the subsignal decomposed at the k time from the current residual to obtain a new residual signal:
and S5, repeating the steps S2 to S4 for the new residual signal until the energy of the residual signal is less than a certain threshold value or the iteration number exceeds the set maximum iteration number, and ending the iteration.
S6, multiplying the decomposed matching Rake wavelets by the corresponding amplitude coefficients to obtain sub signals:
the original seismic signal may be decomposed into sub-signals and decomposed residuals using matching pursuits:
where K is the total number of iterations performed, RK(t) is the matching pursuit decomposition residue, which can be ignored, and the in-situ seismic signal is approximately written as:
all sub-signals can be regarded as reflection signals of a single interface, and the actual seismic record can be regarded as superposition of a plurality of reflection interface reflection waves.
S7, for all the decomposed sub-signals, the amplitude spectrum and the phase spectrum are obtained by fourier transform. For any sub-signal, its propagation time τkIs calculated according to the extreme value of the instantaneous amplitude of the analytic signal, so that the propagation time of all the sub-signals in the underground medium is known. The conventional time-frequency decomposition method based on matching pursuit usually performs WVD on the decomposed sub-signals to obtain a high-resolution time-frequency spectrum. In fact, according to the stacking principle, the seismic record can be regarded as stacking of a plurality of single interface reflected waves, so that the decomposed wavelets are subjected to fourier transform:
and Fourier transform the result Sk(f) The ideal time frequency spectrum IDFT can be obtained by placing the IDFT on the corresponding propagation time of the interfaces(t,f)：
Ideally, the spectrum is a complex matrix, whose elements are written in exponential form,
the amplitude spectrum a (t, f), the phase spectrum Φ (t, f), and the rotational phase spectrum Θ (t, f) are:
in one embodiment of the invention, a model seismic trace (FIG. 2) is match tracked (FIG. 3) and its ideal amplitude spectrum and rotational phase spectrum are calculated (see FIGS. 4(a) - (c)). Comparing the spectrum decomposition result (fig. 4) based on the present invention with the spectrum decomposition results of other time-frequency analysis methods (fig. 5(a), 5(b), 6), it can be seen that: the time-frequency resolution of S transformation is relatively poor, and the phase spectrum of the S transformation comprises a linear phase and a rotating phase, so that the practicability is poor; the amplitude spectrum in the WVD decomposition result based on matching pursuit represents the energy of the signal but not the amplitude, and no phase spectrum exists; the spectrum decomposition result of the invention has higher time and spectrum resolution, can obtain the information of the rotation phase and has higher application value.
While the present invention has been described with reference to the embodiments shown in the drawings, the present invention is not limited to the embodiments, which are illustrative and not restrictive, and it will be apparent to those skilled in the art that various changes and modifications can be made therein without departing from the spirit and scope of the invention as defined in the appended claims.
Claims (8)
1. An ideal seismic spectrum decomposition method based on phase-varying Rake wavelet matching pursuit is characterized by comprising the following steps:
(1) inputting a path of seismic data as an initial residual signal;
(2) constructing a complex analysis signal of a current residual signal, then calculating an instantaneous envelope of the complex analysis signal, and determining the delay time of a matched wavelet atom according to the instantaneous envelope;
(3) constructing a zero-phase Rake wavelet according to the delay time of the matched wavelet atoms, scanning the dominant frequency and the phase of the zero-phase Rake wavelet, and constructing a phase-variable Rake wavelet base;
(4) extracting a wavelet atom which is most matched with a current residual signal from a phase-changing Rake wavelet base, solving an amplitude coefficient corresponding to the most matched wavelet atom, multiplying the amplitude coefficient and the most matched wavelet atom to obtain a sub-signal obtained by decomposition, then subtracting the sub-signal obtained by decomposition from the current residual signal, and taking the obtained residual as a new residual signal;
(5) repeating the steps (2) to (4) until the iteration number reaches a preset value or the residual energy is less than a certain threshold value;
(6) after iteration is finished, obtaining a plurality of phase-variable Rake matching wavelets and corresponding amplitude coefficients, and respectively and correspondingly multiplying the matching wavelets and the amplitude coefficients to obtain single-interface reflection seismic signals serving as a sub-signal set;
(7) and performing Fourier transform on each sub-signal in the sub-signal set, and corresponding the obtained Fourier transform spectrum to the delay time corresponding to the interface to finally obtain an ideal time-frequency spectrum decomposition result.
2. The method for decomposing an ideal seismic spectrum based on the phase-varying Ricker wavelet matching pursuit as claimed in claim 1, wherein the step (2) comprises:
constructing a complex analysis signal of the residual signal s (t):
S(t)＝s(t)+isH(t)，
wherein s isH(t) represents the Hilbert transform of the real signal s (t), the instantaneous amplitude envelope of the complex analysis signal being:
and taking the time corresponding to the maximum value of the instantaneous amplitude envelope as the delay time of the matched wavelet atom.
3. The method for decomposing an ideal seismic spectrum based on the phase-varying Ricker wavelet matching pursuit as claimed in claim 1, wherein the step (3) comprises:
(31) constructing a zero-phase Rake wavelet according to the delay time of the matched wavelet atoms:
w(t)＝{1-2[πfd(t-τ)]2}exp{-[πfd(t-τ)]2}，
where τ is the delay time of the matched wavelet atom, fdDominant frequency of the Rake wavelet;
(32) by scanning the time-frequency parameter gamma ═ fdPhi, constructing a phase-variable Rake wavelet atom library, wherein the scanning range of the wavelet main frequency is 0 to Nyquist frequency, the scanning range of the phase is 0-180 degrees, and the constructed phase-variable Rake wavelet library is as follows:
wγ(t)＝w(t-τ,fd)cos(φ)-wH(t-τ,fd)sin(φ)，
where φ is the phase of the Rake wavelet, wH(t) is the Hilbert transform of the zero-phase Rake wavelet w (t).
4. The method for decomposing an ideal seismic spectrum based on the phase-varying Ricker wavelet matching pursuit as claimed in claim 3, wherein the step (4) comprises:
for the first random iteration number k, all wavelet atoms in the phase-varying Rake wavelet base and the current residual signal sres k-1(t) performing correlation calculation to obtain a correlation coefficient corr:
wherein,representing the current residual signal and the wavelet atom wγ(t) the inner product of the (t),representing wavelet atoms wγ(t) a signal normalization operator;
then, the obtained value is obtainedMaximum value of correlation coefficient and corresponding wavelet atom
The wavelet atomI.e. the wavelet atom which is most matched with the current residual signal, and then the wavelet atom is obtained by calculationCorresponding amplitude coefficient Ck：
The sub-signal decomposed by the k iteration isSubtracting the sub-signal decomposed at the k time from the current residual to obtain a new residual signal sres k(t)：
5. The method for decomposing an ideal seismic spectrum based on phase-varying Ricker wavelet matching pursuit according to claim 4, wherein the step (6) comprises:
will be provided withAs the sub-signals generated in the kth iteration, all the sub-signals constitute a set of sub-signals.
6. The method for decomposing an ideal seismic spectrum based on phase-varying Ricker wavelet matching pursuit according to claim 5, wherein the step (7) comprises:
evaluating a subsignal sk(t) Fourier transform spectrum Sk(f)：
Wherein f represents the frequency of the fourier transform;
then Fourier transform spectrum Sk(f) Corresponding to the delay time corresponding to the interface to obtain the ideal time frequency spectrum IDFTs(t,f)：
And obtaining the amplitude spectrum and the phase spectrum of the reflected seismic waves of the interface.
7. The method of claim 5, wherein the IDFT is an ideal time spectrum IDFTsThe amplitude spectrum a (t, f) and the phase spectrum Φ (t, f) of (t, f) are:
8. the method of claim 5, wherein the step (7) further comprises obtaining an IDFT (inverse discrete Fourier transform) of the ideal time-frequency spectrumsRotational phase spectrum Θ (t, f) of (t, f):
wherein phik(f) Is a propagation time of τkBest match the phase of the rake wavelet.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810962986.9A CN109375265B (en) | 2018-08-22 | 2018-08-22 | Ideal seismic spectrum decomposing method based on phase-varying Rake wavelet matching pursuit |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810962986.9A CN109375265B (en) | 2018-08-22 | 2018-08-22 | Ideal seismic spectrum decomposing method based on phase-varying Rake wavelet matching pursuit |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109375265A true CN109375265A (en) | 2019-02-22 |
| CN109375265B CN109375265B (en) | 2020-01-17 |
Family
ID=65404522
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810962986.9A Expired - Fee Related CN109375265B (en) | 2018-08-22 | 2018-08-22 | Ideal seismic spectrum decomposing method based on phase-varying Rake wavelet matching pursuit |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109375265B (en) |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110133717A (en) * | 2019-04-15 | 2019-08-16 | 长江大学 | Determine the method and apparatus of regional earthquake wave phase |
| CN111142156A (en) * | 2020-01-06 | 2020-05-12 | 中国石油化工股份有限公司 | Phase control-based seismic strong shielding time-frequency information extraction and stripping method |
| CN112213782A (en) * | 2020-09-29 | 2021-01-12 | 中国石油大学(北京) | Processing method and device for sub-phase seismic data and server |
| CN112666609A (en) * | 2020-11-30 | 2021-04-16 | 中国石油大学(华东) | Seismic data low-frequency compensation method and device based on sparse envelope |
| CN113703049A (en) * | 2021-08-27 | 2021-11-26 | 中铁二十局集团安哥拉国际有限责任公司 | Wavelet estimation method based on n-time Fourier spectrum |
| CN114384580A (en) * | 2021-12-31 | 2022-04-22 | 同济大学 | Vibroseis-based ideal wavelet customizing method |
| CN116626760A (en) * | 2023-06-02 | 2023-08-22 | 四川中成煤田物探工程院有限公司 | Stratum discontinuity detection method and device based on self-adaptive high-order maximum entropy WVD |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105022090A (en) * | 2015-07-14 | 2015-11-04 | 北京博达瑞恒科技有限公司 | Wavelet decomposition-based earthquake spectrum decomposition method |
| CN105353408A (en) * | 2015-11-20 | 2016-02-24 | 电子科技大学 | Wigner higher-order spectrum seismic signal spectral decomposition method based on matching pursuit |
| CN105388527A (en) * | 2015-11-30 | 2016-03-09 | 中国石油大学(北京) | Hydrocarbon detection method based on complex domain matching pursuit algorithm |
| WO2016063125A1 (en) * | 2014-10-23 | 2016-04-28 | Cgg Services Sa | Imaging the near subsurface with surface consistent deconvolution operators |
| CN106855638A (en) * | 2016-12-19 | 2017-06-16 | 中国石油天然气股份有限公司 | Matching pursuit seismic spectrum decomposition method and device |
| CN107255831A (en) * | 2017-06-14 | 2017-10-17 | 中国石油化工股份有限公司 | A kind of extracting method of prestack frequency dispersion attribute |
-
2018
- 2018-08-22 CN CN201810962986.9A patent/CN109375265B/en not_active Expired - Fee Related
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2016063125A1 (en) * | 2014-10-23 | 2016-04-28 | Cgg Services Sa | Imaging the near subsurface with surface consistent deconvolution operators |
| CN105022090A (en) * | 2015-07-14 | 2015-11-04 | 北京博达瑞恒科技有限公司 | Wavelet decomposition-based earthquake spectrum decomposition method |
| CN105353408A (en) * | 2015-11-20 | 2016-02-24 | 电子科技大学 | Wigner higher-order spectrum seismic signal spectral decomposition method based on matching pursuit |
| CN105388527A (en) * | 2015-11-30 | 2016-03-09 | 中国石油大学(北京) | Hydrocarbon detection method based on complex domain matching pursuit algorithm |
| CN106855638A (en) * | 2016-12-19 | 2017-06-16 | 中国石油天然气股份有限公司 | Matching pursuit seismic spectrum decomposition method and device |
| CN107255831A (en) * | 2017-06-14 | 2017-10-17 | 中国石油化工股份有限公司 | A kind of extracting method of prestack frequency dispersion attribute |
Non-Patent Citations (1)
| Title |
|---|
| 李坤 等: "基于匹配追踪谱分解的时频域FAVO流体识别方法", 《石油学报》 * |
Cited By (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110133717A (en) * | 2019-04-15 | 2019-08-16 | 长江大学 | Determine the method and apparatus of regional earthquake wave phase |
| CN111142156A (en) * | 2020-01-06 | 2020-05-12 | 中国石油化工股份有限公司 | Phase control-based seismic strong shielding time-frequency information extraction and stripping method |
| CN111142156B (en) * | 2020-01-06 | 2021-07-20 | 中国石油化工股份有限公司 | Phase control-based seismic strong shielding time-frequency information extraction and stripping method |
| CN112213782A (en) * | 2020-09-29 | 2021-01-12 | 中国石油大学(北京) | Processing method and device for sub-phase seismic data and server |
| CN112213782B (en) * | 2020-09-29 | 2022-03-04 | 中国石油大学(北京) | Processing method and device for sub-phase seismic data and server |
| CN112666609A (en) * | 2020-11-30 | 2021-04-16 | 中国石油大学(华东) | Seismic data low-frequency compensation method and device based on sparse envelope |
| CN113703049A (en) * | 2021-08-27 | 2021-11-26 | 中铁二十局集团安哥拉国际有限责任公司 | Wavelet estimation method based on n-time Fourier spectrum |
| CN113703049B (en) * | 2021-08-27 | 2024-04-23 | 中铁二十局集团安哥拉国际有限责任公司 | Wavelet estimation method based on n times of Fourier spectrum |
| CN114384580A (en) * | 2021-12-31 | 2022-04-22 | 同济大学 | Vibroseis-based ideal wavelet customizing method |
| CN114384580B (en) * | 2021-12-31 | 2023-05-02 | 同济大学 | Ideal wavelet customizing method based on controllable focus |
| CN116626760A (en) * | 2023-06-02 | 2023-08-22 | 四川中成煤田物探工程院有限公司 | Stratum discontinuity detection method and device based on self-adaptive high-order maximum entropy WVD |
| CN116626760B (en) * | 2023-06-02 | 2024-05-10 | 四川省自然资源投资集团物探勘查院有限公司 | Stratum discontinuity detection method and device based on self-adaptive high-order maximum entropy WVD |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109375265B (en) | 2020-01-17 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109375265B (en) | Ideal seismic spectrum decomposing method based on phase-varying Rake wavelet matching pursuit | |
| Puryear et al. | Layer-thickness determination and stratigraphic interpretation using spectral inversion: Theory and application | |
| US11333783B2 (en) | Integrated method for estimation of seismic wavelets and synthesis of seismic records in depth domain | |
| Fomel | Local seismic attributes | |
| Wang | Multichannel matching pursuit for seismic trace decomposition | |
| US8265875B2 (en) | Interpolation of periodic data | |
| Liu et al. | Time-frequency analysis of seismic data using local attributes | |
| Ji et al. | Source description of the 1999 Hector Mine, California, earthquake, part I: Wavelet domain inversion theory and resolution analysis | |
| Velis | Stochastic sparse-spike deconvolution | |
| AU2023202187B2 (en) | Use nuos technology to acquire optimized 2d data | |
| Zhang et al. | Seismic sparse-layer reflectivity inversion using basis pursuit decomposition | |
| US10169652B2 (en) | Spatial expansion seismic data processing method and apparatus | |
| EP2689271B1 (en) | Simultaneous wavelet extraction and deconvolution in the time domain | |
| Pérez et al. | High-resolution prestack seismic inversion using a hybrid FISTA least-squares strategy | |
| Hale | A method for estimating apparent displacement vectors from time-lapse seismic images | |
| Li et al. | Plane-wave least-squares reverse time migration with a preconditioned stochastic conjugate gradient method | |
| Duchkov et al. | Discrete almost-symmetric wave packets and multiscale geometrical representation of (seismic) waves | |
| van der Baan et al. | Nonstationary phase estimation using regularized local kurtosis maximization | |
| Sacchi et al. | Minimum weighted norm wavefield reconstruction for AVA imaging | |
| EP3273274A1 (en) | Device and method for estimating pre-stack wavelet model from seismic gathers | |
| Matos et al. | Detecting stratigraphic discontinuities using time-frequency seismic phase residues | |
| Wason et al. | Low-cost time-lapse seismic with distributed compressive sensing—Part 2: Impact on repeatability | |
| Queißer et al. | Localizing CO 2 at Sleipner—Seismic images versus P-wave velocities from waveform inversion | |
| Zhang et al. | Depth-domain seismic reflectivity inversion with compressed sensing technique | |
| Cao et al. | Joint deblending and data reconstruction with focal transformation |
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 | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20200117 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |