CN107132576A - The seismic data Time-Frequency Analysis Method of wavelet transformation is extruded based on second order - Google Patents
The seismic data Time-Frequency Analysis Method of wavelet transformation is extruded based on second order Download PDFInfo
- Publication number
- CN107132576A CN107132576A CN201710543582.1A CN201710543582A CN107132576A CN 107132576 A CN107132576 A CN 107132576A CN 201710543582 A CN201710543582 A CN 201710543582A CN 107132576 A CN107132576 A CN 107132576A
- Authority
- CN
- China
- Prior art keywords
- mrow
- mover
- msub
- mfrac
- omega
- 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
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. analysis, for interpretation, for correction
- G01V1/30—Analysis
- G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/624—Reservoir parameters
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
The present invention discloses a kind of seismic data Time-Frequency Analysis Method that wavelet transformation is extruded based on second order, first on the basis of continuous wavelet transform, introduce second order sync extruding operator, the conversion enters rearrangement by the coefficient to wavelet transformation in Time-Scale Domain, so as to suppress energy dissipation, this method can obtain a kind of analysis result with higher time frequency resolution.Compared to synchronous extruding conversion before, change method to be more suitable for this to seismic signal analyzing with relatively strong warbled signal, in the time frequency analysis for being further applied to actual seismic data, the border of the hydrocarbon reservoir in underground can more accurately be defined, finer portrays the geological informations such as river course and tomography, so as to be conducive to being explained further and well location determination for seismic data.
Description
【Technical field】
The invention belongs to seismic exploration technique field, it is related to a kind of seismic data time-frequency that wavelet transformation is extruded based on second order
Analysis method.
【Background technology】
Due to the decay and filter action during seimic wave propagation by underground medium, seismic signal is typical non-flat
Steady signal.Traditional Fourier conversion, which does not possess, portrays the part of signal frequency ability, and time frequency analysis can be for description
Seismic signal changes with time rule, the local feature of seismic signal is portrayed, so as to reflect the geology corresponding to these features
Structure and reservoir information, therefore time-frequency conversion method is to analyze the ideal tools of non-stationary signal.Construction or selection are suitable
Time frequency analyzing tool is one of key issue of seismic data processing and attributes extraction.With the hair of modern age signal processing theory
Exhibition, emerges substantial amounts of time frequency analyzing tool and Time-Frequency Analysis Method, wherein being much also applied in seismic signal analysis.
Nineteen forty-six, Gabor proposes Short Time Fourier Transform (Short-Time Fourier Transform, STFT),
The advantage of Short Time Fourier Transform is that algorithm realizes easy, explicit physical meaning, the overall time-frequency trend of reflection.Adding window Fourier becomes
The selection for changing middle window function is the compromise of temporal resolution and frequency resolution, once window function is selected, the time point of time-frequency figure
Resolution and frequency resolution are also determined therewith.As typical non-stationary signal, the frequency of seismic signal constantly becomes over time
Change, it is necessary to carry out multiscale analysis.In order to solve the deficiency of adding window Fourier transformation presence, wavelet transformation arises at the historic moment.20 generation
Record the later stage eighties, the perfect wavelet transformation theory such as Morlet.Continuous wavelet transform inherits the excellent of Short Time Fourier Transform
Gesture, realizes the requirement for becoming window processing, with good time frequency localization characteristic.Stockwell proposed S changes equal to 1996
Change, S-transformation is a kind of linear, Time-Frequency Analysis Method of multiresolution, lossless reciprocal, combines Short Time Fourier Transform and small echo
The advantage of conversion and the deficiency for avoiding them.But, the window function of S-transformation is to be changed with fixed trend with frequency, no
Window function can be changed according to application demand.By being replaced to the window function in standard S-transformation according to certain rule, just
The S-transformation expression formula of various generalized forms can be obtained.Generalized S-transform with different window functions portrays energy to unlike signal
Power is different, and good application has been obtained when solving including various practical problems including seism processing.
But, the Linear Time-Frequency Analysis method of the above by uncertainty principle due to being restricted, while separation is in time domain
Distance too small high fdrequency component and on frequency domain it is still highly difficult at a distance of too near low frequency signal.In order to overcome uncertainty principle
Existing Time-Frequency Analysis Method is developed and promoted by limitation, many scholars, it is proposed that many new Time-Frequency Analysis Methods.Its
Middle application more widely be shuffle algorithm, this method by the way that each time-frequency coefficients are reset to its corresponding energy barycenter,
So as to obtain the time-frequency result with higher resolution, but such algorithm amount of calculation than larger, and reset later time-frequency
As a result original signal can not be reconstructed.Daubechies et al. proposed synchronous extruding conversion in 2011
(Synchrosqueezing transform), changes method on the basis of wavelet transformation, introduces the weight based on instantaneous frequency
Algorithm is arranged, the wavelet coefficient of diffusion is expressed on its corresponding actual frequency in frequency direction, signal is not only substantially increased
Time frequency resolution, and can to signal carry out Accurate Reconstruction.But the algorithm is there is also certain limitation, only in analysis
Relatively good result can be obtained during the multicomponent data processing being superimposed by some narrow band signals.
【The content of the invention】
In view of the shortcomings of the prior art, it is a kind of high-precision based on second order extruding wavelet transformation present invention aims at providing
Seismic data Time-Frequency Analysis Method, analysis result is affected by noise smaller, and is easily achieved, and computational efficiency is high.
To reach above-mentioned purpose, the present invention, which is adopted the following technical scheme that, to be achieved：
Based on second order extrude wavelet transformation seismic data Time-Frequency Analysis Method, second order sync extruding wavelet transformation include with
Lower step：
1) continuous wavelet transform is carried out to signal, obtains wavelet transformation spectrum W (a, b), W (a, b) time domain and frequency domain presentation
Formula is respectively
Wherein, ψ (t) is mother wavelet function, and a and b represent yardstick and shift factor respectively, and X (ω) and Ψ (ω) are respectively x
(t) with ψ (t) Fourier transformation, * represents to take conjugation；Wavelet ψ (t) square integrables, no DC component, and meet with
Lower admissibility condition：
2) the corresponding reference frequency of each wavelet coefficient is calculated
Wherein
3) by when m- scale domain wavelet coefficient reset to time-frequency domain, second order sync extrudes the form of wavelet transformation
For
δ (g) represents uni-impulse function；By formula (4), when m- scale domain is all and instantaneous frequencyIt is right
The wavelet coefficient answered carries out longitudinal summation along dimension, and then the result of summation is placed on corresponding to time-frequency domain w
Position；
The seismic data Time-Frequency Analysis Method for extruding wavelet transformation based on above second order specifically includes following steps：
Step 1：A certain seismic profile in 3D data volume chooses typical earthquake road, and carries out two to the geological data
Rank extrudes wavelet transformation, and finds out the frequency content corresponding to abnormal area；
Step 2：Second order sync extruding wavelet transformation is carried out by road to whole seismic profile, extracted relative for abnormal area
The frequency slice answered；
Step 3:Second order extruding wavelet transformation is carried out to whole three-dimensional data, corresponding frequency data body is obtained, according to solution
The layer position information of personnel's offer is provided, horizon slice is made, is explained for geological personnel.
Further, the step 3) by step 2) obtain after the corresponding instantaneous frequency of each wavelet coefficient, pass through energy
Reset, by when m- scale domain wavelet coefficient be expressed in time-frequency domain；
Second order sync extruding wavelet transformation has following several forms：
Distribution form
Approximate form
In formula (5), h (g) represents function unlimited smooth and with local support, andε is expressed as
Prevent the unstable set threshold parameter of numerical computations, by formula (4) or (5), when m- scale domain it is all and
Instantaneous frequencyCorresponding wavelet coefficient carries out longitudinal summation along dimension, when then the result of summation is placed on
Position corresponding to m- frequency domain w.
Further, the step 3) second order sync extruding wavelet transformation discrete form be
The invention has the advantages that：
The present invention extrudes the seismic signal time-frequency analysis method of wavelet transformation based on second order sync, becomes first in continuous wavelet
On the basis of changing, introduce second order sync extruding operator, the conversion by the coefficient to wavelet transformation when m- scale domain weighed
Row, so as to suppress energy dissipation, this method can obtain a kind of analysis result with higher time frequency resolution.Compared to it
Preceding synchronous extruding conversion, change method be more suitable for it is this to seismic signal analyzed with relatively strong warbled signal,
In the time frequency analysis for being further applied to actual seismic data, the side of the hydrocarbon reservoir in underground can be more accurately defined
Boundary, finer portrays the geological informations such as river course and tomography, so as to be conducive to being explained further and well location determination for seismic data.
The present invention extrudes the seismic signal time-frequency analysis method of wavelet transformation based on second order sync, and composite signal shows, phase
Than in traditional Time-Frequency Analysis Method, such as wavelet transformation, synchronization extrude wavelet transformation, and it is right that second order sync extruding wavelet transformation passes through
The rearrangement of wavelet coefficient, obtains a kind of time frequency analysis result with higher time frequency resolution, and can accurately portray letter
Number frequency change with time relation.Cut into slices by extracting single-frequency to two-dimension earthquake section, second order sync extruding wavelet transformation
The result position and border that depict reservoir that can become apparent from.On the horizon slice of 3-d seismic data set, second order
The feature for portraying subsurface fault and river course that the result of synchronous extruding wavelet transformation can become apparent from.
【Brief description of the drawings】
Fig. 1 is several time frequency analysis results of test signal；
(a) test signal；(b) continuous wavelet transform is converted；(c) synchronous extruding wavelet transformation；(d) second order sync extruding is small
Wave conversion
Fig. 2 is certain oil gas field two dimensional cross-section；
Fig. 3 is the 35Hz frequency slices of Fig. 2 two dimensional cross-sections；
(a) wavelet transform result；(b) synchronous extruding transformation results；(c) second order sync extruding wavelet transform result；
Fig. 4 is certain oil field three-dimensional data volume 35Hz horizon slice schematic diagram；
(a) synchronous extruding wavelet transformation 35Hz horizon slices；(b) second order sync extruding wavelet transformation 35Hz horizon slices.
【Embodiment】
With reference to the accompanying drawing in the embodiment of the present invention, the technical scheme in the embodiment of the present invention is carried out clear, complete
Description, embodiment described herein is only a part of embodiment of the present invention, and not all embodiment.Based on this hair
Embodiment in bright, the every other embodiment that those of ordinary skill in the art are obtained under the premise of no creative work,
Belong to the scope of the present invention.
The seismic data Time-Frequency Analysis Method for extruding wavelet transformation based on second order sync proposed by the invention is specifically included
Following steps：
Second order sync extruding wavelet transformation comprises the following steps：
Step 1：Continuous wavelet transform is carried out to signal, wavelet transformation spectrum W (a, b), W (a, b) time domain and frequency domain is obtained
Expression formula is respectively
Wherein, ψ (t) is mother wavelet function, and a and b represent yardstick and shift factor respectively.X (ω) and Ψ (ω) is respectively x
(t) with ψ (t) Fourier transformation, * represents to take conjugation.
Step 2：Calculate the corresponding reference frequency of each wavelet coefficient
Wherein
We are with linear chirp signal h (t)=Acos (2 π f0t2) exemplified by, can be worked as | W (a, b) | when ≠ 0,Therefore, when signal Analysis is that linear chirp signals are, by the operation of formula (3), can to become from small echo
The instantaneous frequency of signal is demodulated in the coefficient changed.
Step 3：By when m- scale domain wavelet coefficient reset to time-frequency domain；
Obtained after the corresponding instantaneous frequency of each wavelet coefficient, reset by energy by step 2, can by when m- chi
The wavelet coefficient in degree domain is expressed in time-frequency domain.
Second order sync extruding wavelet transformation has following several forms：
Distribution form
Approximate form
In formula (4), δ (g) represents uni-impulse function；In formula (5), h (g) represents unlimited smooth and carries office
The function of portion's support, andε, which is represented, prevents the unstable set threshold parameter of numerical computations.It is logical
Cross formula (4) or (5), when m- scale domain is all and instantaneous frequencyCorresponding wavelet coefficient enters along dimension
The summation of row longitudinal direction, is then placed on the position corresponding to time-frequency domain w by the result of summation.
When carrying out numerical computations, it is necessary to discrete to yardstick a progress, discrete way is aj=a02j/nv, a0=1/T, T are
The time span of signal, nv typically takes 32 or 64. yardstick interval aj-aj-1=(Δ a)j；Then division w is carried out to frequencyj=Δ w
× j, wherein Δ w=1/ Δs t/N, Δ t, N are respectively sampling interval and the sampling number of signal.Second order sync extrudes wavelet transformation
Discrete form be
The material base of the present invention is 3-d seismic data set, using Dou Shizhu roads processing method.
The present invention is concretely comprised the following steps based on the seismic data Time-Frequency Analysis Method that second order sync extrudes wavelet transformation：
Step 1:Some typical earthquake roads are chosen in 3-d seismic data set, spectrum analysis is carried out first, the ground is found
Shake the corresponding dominant frequency of data volume.
Step 2:The a certain seismic profile of 3D data volume is chosen, second order sync extruding wavelet transformation is carried out by road, then
Make the frequency slice corresponding to dominant frequency.
Step 3：Second order sync extruding wavelet transformation is carried out by road to whole 3D data volume, the frequency near dominant frequency is extracted
Section, obtains frequency data body；
Step 4：The layer position information explained according to geological personnel, makes horizon slice.
Effect analysis
First, numerical experiment
First, contrast second order sync extrudes wavelet transformation with conventional wavelet transformation and synchronous extruding wavelet transformation to synthesis
The result of signal time frequency analysis.Composite signal such as Fig. 1 (a), this signal is made up of three below component：
s1(t)=sin (2 π (330t+16cos (3 π t)))
s2(t)=asin (160 π t)
Fig. 1 (b)-(d) corresponds to the knot that wavelet transformation, synchronous extruding wavelet transformation and second order sync extrude wavelet transformation respectively
Really.Although the result of these three conversion can be due to be limited by uncertainty principle the two components, small echo
The obvious energy dissipation that transformation results occur.Comparison diagram 1 (c) and (d) by energy as can be seen that reset, and synchronous extruding is small
The result of wave conversion and second order sync extruding wavelet transformation substantially increases time frequency resolution.It is worth noting that, in Fig. 1 (c)
In, in the Instantaneous frequency variations of signal faster local (black arrow), the time-frequency of second order sync extruding wavelet transformation gather (d)
Collection property is better than the result of synchronous extruding wavelet transformation, so as to more accurately portray the time-varying characteristics of signal.
2nd, actual seismic data
Below, second order sync is extruded wavelet transformation application actual seismic data time frequency analysis by us, and compared for small echo
Conversion and second order sync extrude the result of wavelet transformation.Scheme the two dimensional cross-section that (2) are certain offshore oil and gas field, the section contains 1000
Road, the sampled point of per pass is 751, and the sampling interval is 2ms.Time-frequency conversion is carried out to per pass geological data, 35Hz is therefrom extracted
Frequency slice, such as figure (3) shown in., can be more because second order sync extruding wavelet transformation has higher time frequency resolution
Accurately position the position of reservoir.The position of reservoir is indicated in figure (3) with black box.Finally, to three dimensions in certain work area
Time frequency analysis is carried out according to body.Synchronous extruding wavelet transformation and second order sync extruding wavelet transformation is respectively adopted to carry out the data volume
Analysis, obtains 35Hz frequency data body first, the layer position data provided according to the personnel of explanation, makes horizon slice.Such as
Scheme shown in (4), it can be seen that the feature of underground palaeostream becomes comparison on the synchronously horizon slice of extruding wavelet transformation and obscured,
Due to higher time frequency resolution, on the horizon slice that second order sync extrudes wavelet transformation, river course (black arrow)
Feature obtains very clearly portraying.
Above example is merely to illustrate technical scheme rather than its limitations, although with reference to above-described embodiment pair
The present invention is described in detail, and those of ordinary skill in the art can still be carried out to specific embodiments of the present invention
Modification or equivalent substitution, and these any modifications or equivalent substitution without departing from spirit and scope of the invention, it exists
Within the claims of the present invention.
Claims (3)
1. the seismic data Time-Frequency Analysis Method of wavelet transformation is extruded based on second order, it is characterised in that：
Second order sync extruding wavelet transformation comprises the following steps：
1) continuous wavelet transform is carried out to signal, obtains wavelet transformation spectrum W (a, b), W (a, b) time domain and frequency-domain expression point
It is not
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<msqrt>
<mi>a</mi>
</msqrt>
</mfrac>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mi>&infin;</mi>
</msubsup>
<mi>x</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>&psi;</mi>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mi>t</mi>
<mo>-</mo>
<mi>b</mi>
</mrow>
<mi>a</mi>
</mfrac>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<mfrac>
<msqrt>
<mi>a</mi>
</msqrt>
<mrow>
<mn>2</mn>
<mi>&pi;</mi>
</mrow>
</mfrac>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mi>&infin;</mi>
</msubsup>
<mi>X</mi>
<mrow>
<mo>(</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<msup>
<mi>&Psi;</mi>
<mo>*</mo>
</msup>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<msup>
<mi>e</mi>
<mrow>
<mi>i</mi>
<mi>b</mi>
<mi>&omega;</mi>
</mrow>
</msup>
<mi>d</mi>
<mi>&omega;</mi>
<mo>,</mo>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
Wherein, ψ (t) is mother wavelet function, and a and b represent yardstick and shift factor respectively, and X (ω) and Ψ (ω) are respectively x (t)
With ψ (t) Fourier transformation, * represents to take conjugation；Wavelet ψ (t) square integrables, no DC component, and meet following hold
Perhaps property condition：
<mrow>
<msub>
<mi>C</mi>
<mi>&psi;</mi>
</msub>
<mo>=</mo>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mi>&infin;</mi>
</msubsup>
<mfrac>
<msup>
<mrow>
<mo>|</mo>
<mi>&psi;</mi>
<mrow>
<mo>(</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
<mn>2</mn>
</msup>
<mrow>
<mo>|</mo>
<mi>&omega;</mi>
<mo>|</mo>
</mrow>
</mfrac>
<mi>d</mi>
<mi>&omega;</mi>
<mo><</mo>
<mi>&infin;</mi>
<mo>.</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mo>)</mo>
</mrow>
</mrow>
2) the corresponding reference frequency of each wavelet coefficient is calculated
<mrow>
<mi>&omega;</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mover>
<mi>&omega;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mover>
<mi>q</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>(</mo>
<mi>b</mi>
<mo>-</mo>
<mover>
<mi>&tau;</mi>
<mo>^</mo>
</mover>
<mo>(</mo>
<mrow>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
<mtd>
<mrow>
<msub>
<mo>&part;</mo>
<mi>b</mi>
</msub>
<mover>
<mi>&tau;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>&NotEqual;</mo>
<mn>0</mn>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mover>
<mi>&omega;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
<mtd>
<mrow>
<msub>
<mo>&part;</mo>
<mi>b</mi>
</msub>
<mover>
<mi>&tau;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mn>0</mn>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
Wherein
<mrow>
<mover>
<mi>&omega;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<msub>
<mo>&part;</mo>
<mi>b</mi>
</msub>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mi>i</mi>
<mn>2</mn>
<mi>&pi;</mi>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
</mrow>
<mrow>
<mover>
<mi>&tau;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mi>&infin;</mi>
</msubsup>
<mi>t</mi>
<mi>x</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mfrac>
<mn>1</mn>
<msqrt>
<mi>a</mi>
</msqrt>
</mfrac>
<mi>&psi;</mi>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mi>t</mi>
<mo>-</mo>
<mi>b</mi>
</mrow>
<mi>a</mi>
</mfrac>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
</mrow>
<mrow>
<mover>
<mi>q</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>Re</mi>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<msub>
<mo>&part;</mo>
<mi>b</mi>
</msub>
<mover>
<mi>&omega;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<msub>
<mo>&part;</mo>
<mi>b</mi>
</msub>
<mover>
<mi>&tau;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
3) by when m- scale domain wavelet coefficient reset to time-frequency domain, the form of second order sync extruding wavelet transformation is
<mrow>
<mi>T</mi>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mrow>
<mo>{</mo>
<mi>a</mi>
<mo>:</mo>
<mi>a</mi>
<mo>></mo>
<mn>0</mn>
<mo>,</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>w</mi>
<mo>,</mo>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>&NotEqual;</mo>
<mn>0</mn>
<mo>}</mo>
</mrow>
</msub>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mi>&delta;</mi>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>-</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mo>(</mo>
<mrow>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<msup>
<mi>a</mi>
<mrow>
<mo>-</mo>
<mn>3</mn>
<mo>/</mo>
<mn>2</mn>
</mrow>
</msup>
<mi>d</mi>
<mi>a</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
</mrow>
δ (g) represents uni-impulse function；By formula (4), when m- scale domain is all and instantaneous frequencyIt is corresponding small
Wave system number carries out longitudinal summation along dimension, and the result of summation then is placed on into the position corresponding to time-frequency domain w
Put；
The seismic data Time-Frequency Analysis Method for extruding wavelet transformation based on above second order specifically includes following steps：
Step 1：A certain seismic profile in 3D data volume chooses typical earthquake road, and it is crowded to carry out second order to the geological data
Wavelet transformation is pressed, and finds out the frequency content corresponding to abnormal area；
Step 2：Second order sync extruding wavelet transformation is carried out by road to whole seismic profile, extracted corresponding for abnormal area
Frequency slice；
Step 3:Second order extruding wavelet transformation is carried out to whole three-dimensional data, corresponding frequency data body is obtained, according to explanation people
The layer position information that member provides, makes horizon slice, is explained for geological personnel.
2. the seismic data Time-Frequency Analysis Method according to claim 1 that wavelet transformation is extruded based on second order, its feature is existed
In：The step 3) by step 2) obtain after the corresponding instantaneous frequency of each wavelet coefficient, reset by energy, by when it is m-
The wavelet coefficient of scale domain is expressed in time-frequency domain；
Second order sync extruding wavelet transformation has following several forms：
Distribution form
<mrow>
<mi>T</mi>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mrow>
<mo>{</mo>
<mi>a</mi>
<mo>:</mo>
<mi>a</mi>
<mo>></mo>
<mn>0</mn>
<mo>,</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>w</mi>
<mo>,</mo>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>&NotEqual;</mo>
<mn>0</mn>
<mo>}</mo>
</mrow>
</msub>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mi>&delta;</mi>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>-</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mo>(</mo>
<mrow>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<msup>
<mi>a</mi>
<mrow>
<mo>-</mo>
<mn>3</mn>
<mo>/</mo>
<mn>2</mn>
</mrow>
</msup>
<mi>d</mi>
<mi>a</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
</mrow>
Approximate form
<mrow>
<msubsup>
<mi>T</mi>
<mi>&epsiv;</mi>
<mi>&lambda;</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mrow>
<mo>{</mo>
<mi>a</mi>
<mo>:</mo>
<mi>a</mi>
<mo>></mo>
<mn>0</mn>
<mo>,</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>w</mi>
<mo>,</mo>
<mo>|</mo>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mo>></mo>
<mi>&epsiv;</mi>
<mo>}</mo>
</mrow>
</msub>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mfrac>
<mn>1</mn>
<mi>&lambda;</mi>
</mfrac>
<mi>h</mi>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mi>w</mi>
<mo>-</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mrow>
<mo>(</mo>
<mi>a</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
</mrow>
<mi>&lambda;</mi>
</mfrac>
<mo>)</mo>
</mrow>
<msup>
<mi>a</mi>
<mrow>
<mo>-</mo>
<mn>3</mn>
<mo>/</mo>
<mn>2</mn>
</mrow>
</msup>
<mi>d</mi>
<mi>a</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>5</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula (5), h (g) represents function unlimited smooth and with local support, andε represents anti-
Only there is the unstable set threshold parameter of numerical computations, by formula (4) or (5), when m- scale domain it is all and instantaneous
FrequencyCorresponding wavelet coefficient carries out longitudinal summation along dimension, m- when then the result of summation is placed on
Position corresponding to frequency domain w.
3. the seismic data Time-Frequency Analysis Method according to claim 1 that wavelet transformation is extruded based on second order, its feature is existed
In：The step 3) second order sync extruding wavelet transformation discrete form be
<mrow>
<mi>T</mi>
<mrow>
<mo>(</mo>
<mi>w</mi>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munder>
<mo>&Sigma;</mo>
<mrow>
<msub>
<mi>a</mi>
<mi>j</mi>
</msub>
<mo>:</mo>
<mo>|</mo>
<mover>
<mi>&omega;</mi>
<mo>~</mo>
</mover>
<mrow>
<mo>(</mo>
<msub>
<mi>a</mi>
<mi>j</mi>
</msub>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>w</mi>
<mo>|</mo>
<mo>&le;</mo>
<mi>&Delta;</mi>
<mi>w</mi>
<mo>/</mo>
<mn>2</mn>
</mrow>
</munder>
<mi>W</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>a</mi>
<mi>j</mi>
</msub>
<mo>,</mo>
<mi>b</mi>
<mo>)</mo>
</mrow>
<msub>
<mrow>
<mo>(</mo>
<mi>&Delta;</mi>
<mi>a</mi>
<mo>)</mo>
</mrow>
<mi>j</mi>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>6</mn>
<mo>)</mo>
</mrow>
<mo>.</mo>
</mrow>
2
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710543582.1A CN107132576B (en) | 2017-07-05 | 2017-07-05 | The seismic data Time-Frequency Analysis Method of wavelet transformation is squeezed based on second order sync |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710543582.1A CN107132576B (en) | 2017-07-05 | 2017-07-05 | The seismic data Time-Frequency Analysis Method of wavelet transformation is squeezed based on second order sync |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107132576A true CN107132576A (en) | 2017-09-05 |
| CN107132576B CN107132576B (en) | 2018-10-30 |
Family
ID=59737347
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710543582.1A Active CN107132576B (en) | 2017-07-05 | 2017-07-05 | The seismic data Time-Frequency Analysis Method of wavelet transformation is squeezed based on second order sync |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107132576B (en) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109884694A (en) * | 2019-02-19 | 2019-06-14 | 西安交通大学 | A kind of high-speed rail focus seismic signal time-frequency analysis method based on extruding adding window Fourier transformation |
| CN110941013A (en) * | 2018-09-21 | 2020-03-31 | 中国石油化工股份有限公司 | Time frequency domain energy focusing method and reservoir prediction method |
| CN111273347A (en) * | 2018-12-05 | 2020-06-12 | 中国石油天然气集团有限公司 | Seismic signal decomposition method and device |
| CN111366978A (en) * | 2020-04-29 | 2020-07-03 | 中国海洋石油集团有限公司 | Earthquake time-frequency analysis method and system based on multi-extrusion wavelet transform |
| CN111427091A (en) * | 2020-05-06 | 2020-07-17 | 芯元(浙江)科技有限公司 | Seismic exploration signal random noise suppression method by squeezing short-time Fourier transform |
| CN112394402A (en) * | 2019-08-19 | 2021-02-23 | 中国石油化工股份有限公司 | Method and system for detecting microseism signals based on synchronous extrusion wavelet transform |
| CN113435259A (en) * | 2021-06-07 | 2021-09-24 | 吉林大学 | Tensor decomposition-based satellite magnetic field data fusion seismic anomaly extraction method |
| CN114252915A (en) * | 2021-11-03 | 2022-03-29 | 成都理工大学 | Oil and gas reservoir identification method based on second-order horizontal multiple synchronous extrusion transformation |
| CN114563824A (en) * | 2022-02-25 | 2022-05-31 | 成都理工大学 | Identification method for second-order multiple synchronous extrusion polynomial chirp transform thin reservoir |
| CN112904412B (en) * | 2019-12-03 | 2024-02-23 | 北京矿冶科技集团有限公司 | Mine microseismic signal P-wave first arrival time extraction method and system |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160178772A1 (en) * | 2013-05-27 | 2016-06-23 | Statoil Petroleum As | High Resolution Estimation of Attenuation from Vertical Seismic Profiles |
| CN106291700A (en) * | 2016-09-28 | 2017-01-04 | 西安交通大学 | Based on the earthquake weighted average instantaneous frequency distilling method synchronizing extruding conversion |
| CN106556865A (en) * | 2016-11-25 | 2017-04-05 | 成都理工大学 | A kind of tandem type seismic signal optimizes time-frequency conversion method |
| CN104880730B (en) * | 2015-03-27 | 2017-04-26 | 西安交通大学 | Seismic data time-frequency analysis and attenuation estimation method based on Synchrosqueezing transform |
-
2017
- 2017-07-05 CN CN201710543582.1A patent/CN107132576B/en active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160178772A1 (en) * | 2013-05-27 | 2016-06-23 | Statoil Petroleum As | High Resolution Estimation of Attenuation from Vertical Seismic Profiles |
| CN104880730B (en) * | 2015-03-27 | 2017-04-26 | 西安交通大学 | Seismic data time-frequency analysis and attenuation estimation method based on Synchrosqueezing transform |
| CN106291700A (en) * | 2016-09-28 | 2017-01-04 | 西安交通大学 | Based on the earthquake weighted average instantaneous frequency distilling method synchronizing extruding conversion |
| CN106556865A (en) * | 2016-11-25 | 2017-04-05 | 成都理工大学 | A kind of tandem type seismic signal optimizes time-frequency conversion method |
Non-Patent Citations (2)
| Title |
|---|
| T OBERLIN: "The second-order wavelet synchrosqueezing transform", 《2017 IEEE INTERNATIONAL CONFERENCE ON ACOUSTIC，SPEECH AND SIGNAL PROCESSING（ICASSP）》 * |
| 游龙庭: "基于解析信号重构的同步挤压小波变换的时频谱影响因素分析", 《CT理论与应用研究》 * |
Cited By (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110941013A (en) * | 2018-09-21 | 2020-03-31 | 中国石油化工股份有限公司 | Time frequency domain energy focusing method and reservoir prediction method |
| CN111273347A (en) * | 2018-12-05 | 2020-06-12 | 中国石油天然气集团有限公司 | Seismic signal decomposition method and device |
| CN109884694A (en) * | 2019-02-19 | 2019-06-14 | 西安交通大学 | A kind of high-speed rail focus seismic signal time-frequency analysis method based on extruding adding window Fourier transformation |
| CN109884694B (en) * | 2019-02-19 | 2020-06-19 | 西安交通大学 | High-speed rail seismic source seismic signal time-frequency analysis method based on extrusion windowing Fourier transform |
| CN112394402A (en) * | 2019-08-19 | 2021-02-23 | 中国石油化工股份有限公司 | Method and system for detecting microseism signals based on synchronous extrusion wavelet transform |
| CN112904412B (en) * | 2019-12-03 | 2024-02-23 | 北京矿冶科技集团有限公司 | Mine microseismic signal P-wave first arrival time extraction method and system |
| CN111366978A (en) * | 2020-04-29 | 2020-07-03 | 中国海洋石油集团有限公司 | Earthquake time-frequency analysis method and system based on multi-extrusion wavelet transform |
| CN111427091A (en) * | 2020-05-06 | 2020-07-17 | 芯元(浙江)科技有限公司 | Seismic exploration signal random noise suppression method by squeezing short-time Fourier transform |
| CN113435259A (en) * | 2021-06-07 | 2021-09-24 | 吉林大学 | Tensor decomposition-based satellite magnetic field data fusion seismic anomaly extraction method |
| CN113435259B (en) * | 2021-06-07 | 2022-06-03 | 吉林大学 | Tensor decomposition-based satellite magnetic field data fusion earthquake anomaly extraction method |
| CN114252915A (en) * | 2021-11-03 | 2022-03-29 | 成都理工大学 | Oil and gas reservoir identification method based on second-order horizontal multiple synchronous extrusion transformation |
| CN114563824A (en) * | 2022-02-25 | 2022-05-31 | 成都理工大学 | Identification method for second-order multiple synchronous extrusion polynomial chirp transform thin reservoir |
| CN114563824B (en) * | 2022-02-25 | 2024-01-30 | 成都理工大学 | Second-order multiple synchronous extrusion polynomial chirp let transformation thin reservoir identification method |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107132576B (en) | 2018-10-30 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107132576A (en) | The seismic data Time-Frequency Analysis Method of wavelet transformation is extruded based on second order | |
| CN102879822B (en) | A kind of seismic multi-attribute fusion method based on contourlet transformation | |
| Taborda et al. | Ground‐motion simulation and validation of the 2008 Chino Hills, California, earthquake | |
| CN107817527A (en) | Seismic exploration in desert stochastic noise suppression method based on the sparse compressed sensing of block | |
| CN111505716B (en) | Seismic time-frequency analysis method for extracting generalized Chirplet transform based on time synchronization | |
| CN101545984A (en) | Seismic coherence algorithm based on wavelet transformation | |
| CN101545983A (en) | Multiattribute frequency division imaging method based on wavelet transformation | |
| CN102221708B (en) | Fractional-Fourier-transform-based random noise suppression method | |
| CN107632320A (en) | Seismic data Time-Frequency Analysis Method based on synchronous extraction S-transformation | |
| Pang et al. | Automatic passive data selection in time domain for imaging near-surface surface waves | |
| CN101923176B (en) | Method for oil and gas detection by utilizing seismic data instantaneous frequency attribute | |
| CN107272063B (en) | Anisotropism depicting method based on high-resolution time frequency analysis and consistency metric | |
| CN107356967A (en) | A kind of sparse optimization method suppressed seismic data and shield interference by force | |
| CN102692647B (en) | Stratum oil-gas possibility prediction method with high time resolution | |
| CN106707334B (en) | A method of improving seismic data resolution | |
| CN107255831A (en) | A kind of extracting method of prestack frequency dispersion attribute | |
| CN102323615B (en) | Method for reservoir prediction and fluid identification with earthquake data and device | |
| CN104880730B (en) | Seismic data time-frequency analysis and attenuation estimation method based on Synchrosqueezing transform | |
| CN105353408A (en) | Wigner higher-order spectrum seismic signal spectral decomposition method based on matching pursuit | |
| CN102305940B (en) | Method for extracting fluid factor | |
| CN105954801A (en) | Frequency-conversion-attribute-based reservoir permeability evaluation method | |
| CN104730576A (en) | Curvelet transform-based denoising method of seismic signals | |
| CN104007466B (en) | The reservoir that a kind of no restriction from borehole data prestack inversion based on P-wave amplitude realizes and fluid prediction method | |
| CN103645504A (en) | Weak earthquake signal processing method based on generalized instantaneous phase and P norm negative norm | |
| CN105005073A (en) | Time-varying wavelet extraction method based on local similarity and evaluation feedback |
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 |