CN104820242A - Trace gather amplitude frequency division compensation method for prestack inversion - Google Patents

Abstract

The invention brings forward a trace gather amplitude frequency division compensation method for prestack inversion. The method comprises the following steps: carrying out preprocessing including denoising, residual dynamic calibration and incident angle trace gather transferring on original prestack trace gather data; then performing frequency division processing based on wavelet transformation on the prestack trace gather data after the preprocessing; then extracting incident angles of middle angles and corresponding amplitude values from the prestack trace gather data with different frequencies by taking a forward modeling AVO characteristic trend as a constraint; performing AVO parabola fitting on the extracted incident angle and the extracted amplitude value of the trace gathers with the different frequencies, and calculating each coefficient of a parabola equation; and finally, according to a parabola equation obtained through the fitting, re-determining reflection coefficients of the different incident angles in the trace gathers with the different frequencies, and according to the re-determined reflection coefficients, performing frequency division amplitude compensation on channels with weak near-channel far-end energy, the change of amplitude along with the incident angles meets the rule of a parabola, the requirement of a prestack inversion method on input data is met, and the reliability and stability of inversion of prestack elastic parameters are improved.

Description

A kind of road collection amplitude frequency division compensation method towards prestack inversion

Technical field

The present invention relates to reservoir exploration technical field, especially a kind of road collection amplitude frequency division compensation method towards prestack inversion.

Background technology

In recent years, increasing along with lithological reservoir exploration demand, seismic reservoir prediction technology supporting with it, particularly prestack attribute analysis and prestack elastic parameter inversion technology obtain significant progress.These technology are used widely in lithology prediction and fluid detection, and its precision of prediction, except affecting by inversion method itself, is also subject to the impact of the input quality of data.Wherein, prestack road collection is as the input data of prestack attribute analysis and prestack elastic parameter inversion, and the impact of its quality on inversion result precision is larger.

At present, those skilled in the art generally recognize prestack road collection denoising, wavelet coherence, ranges of incidence angles, residue moves school and the pre-service such as distortion is very large on the impact of prestack inversion result a long way.Such as, but these pre-service mainly concentrate on the aspects such as improvement data signal to noise ratio (S/N ratio) and resolution, lack fidelity, road collection dynamic information---the research of amplitude, frequency, phase compensation, especially rare especially to wide-angle research.But prestack elastic parameter inversion will obtain stable density item, need incident angle wide as far as possible.In addition, owing to can affect by attenuation by absorption when seismic event is propagated in underground medium, therefore the high-frequency energy of seismic event has loss, and amplitude, frequency and phase place can distort.Although in the prior art, also can compensate various amplitude in process in the early stage of prestack road collection, these compensation methodes are just limited to geometrical attenuation longitudinally and absorption compensation, and do not consider collection compensation transversely.Xiong Dingyu etc. propose the compensation method meeting parabola rule based on road collection with incident angle change in " keeping the noise reduction method of seismologic record prestack AVO attribute " literary composition, but it is when calculating parabolic coefficient, have employed zero incident angle data to calculate, due to reasons such as nearly far away offset distance degree of covering when shortcut goes ground roll, skew are uneven, closely-energy is more weak a long way, all low-angle raw data are not suitable for carrying out asking for of parameter, therefore have certain limitation.In sum, although the amplitude of existing technical method to road collection makes moderate progress, but be all conventional single entire compensation method, do not consider that different frequency has different attenuation by absorption features and AVO feature, or the demand of pre-stack seismic inversion cannot be met, reliable data foundation can not be provided for prestack elastic parameter inversion.

Summary of the invention

For the problems referred to above, the present invention proposes a kind of road collection amplitude frequency division compensation method towards prestack inversion.

This compensation method, comprises the following steps:

Pre-treatment step, comprises denoising, the dynamic school of residue to original prestack road set information and proceeds to the pre-service of penetrating angle gathers;

Scaling down processing step, carries out scaling down processing to pretreated prestack road set information;

Parameter extracting step, extracts to meet and just drills the incident angle of AVO characteristic trend and the amplitude of correspondence from the prestack road set information of different frequency;

The Fitting Calculation step, carries out AVO Parabolic Fit to the incident angle of different frequency road collection extracted and amplitude, calculates each coefficient of parabolic equation;

Energy compensating step, the parabolic equation obtained according to matching redefines reflection coefficient when different incidence angles is concentrated in different frequency road, and the road weak to shortcut, distal energy according to the reflection coefficient redefined carries out frequency division amplitude compensation.

According to embodiments of the invention, in above-mentioned pre-treatment step, denoising is carried out in the beam focusing based on in-phase stacking.

According to embodiments of the invention, the denoising based on the beam focusing of in-phase stacking comprises the following steps:

For the offset distance that is determined, first define the contiguous range of shot point and geophone station, being extracted by the road collection within the scope of this is a beam;

In beam, choose sign back layer, when choosing near sign back layer, window carries out superposition and forms standard track, carries out correlation analysis with beam Nei Gedao;

If the time that the relevant time difference is corresponding more than 1/4 wavelet lengths, reject;

The each road meeting in-phase stacking principle is superposed.

According to embodiments of the invention, it is the method adopting axis of reflection automatic tracing that above-mentioned residue moves school process, with near, middle road, related coefficient is asked for a long way to each prestack given in timing window, obtain the time shift amount relative to zero-offset, then progressively sliding window obtains the time shift amount of whole road collection, and Jiang Daoji evens up.

According to embodiments of the invention, above-mentioned scaling down processing step, carries out scaling down processing based on wavelet transformation to pretreated prestack road set information.

According to embodiments of the invention, in above-mentioned parameter extraction step, preferably extract to meet from the prestack road set information of different frequency and just drill the incident angle of angle and the amplitude of correspondence in AVO characteristic trend.

According to embodiments of the invention, the incident angle of above-mentioned middle angle is preferably 10 ~ 25 degree.

According to embodiments of the invention, in above-mentioned the Fitting Calculation step, carry out AVO Parabolic Fit, by coefficients R, W, V of following formula determination parabolic equation according to the incident angle extracted and amplitude;

$\underset{i=1}{\overset{n}{\mathrm{\Σ}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i$

$\underset{i=1}{\overset{n}{\mathrm{\Σ}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)x_i=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i{x}_i$

$\underset{i=1}{\overset{n}{\mathrm{\Σ}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)x_i^2=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i{x}_i^2$

In formula, x _i=sin ²θ _i, θ _ii-th incident angle; y _iintermediate variable corresponding to i-th incident angle.

According to embodiments of the invention, in above-mentioned the Fitting Calculation step, least square method can be adopted to calculate each coefficient of parabolic equation.

Compared with prior art, one or more embodiment of the present invention can have the following advantages by tool:

The present invention carries out the pretreatment work such as necessary denoising, the dynamic school of residue to seismic traces collection, obtains signal to noise ratio (S/N ratio) high, the basic data that lineups are evened up; For different frequency seismic event, there is different attenuation by absorption features, introduce the thought that frequency division compensates, overcome the shortcoming that frequency factor is not considered in single compensation method; The AVO matching degree of well shake is improve in addition based on the AVO signature amplitude compensation method just drilling trend constraint, expand road collection incident angle range, improve the stability of prestack elastic parameter especially density item, and then improve reliability and the accuracy of prestack elastic parameter inversion.

Other features and advantages of the present invention will be set forth in the following description, and partly become apparent from instructions, or understand by implementing the present invention.Object of the present invention and other advantages realize by structure specifically noted in instructions, claims and accompanying drawing and obtain.

Accompanying drawing explanation

Accompanying drawing is used to provide a further understanding of the present invention, and forms a part for instructions, with embodiments of the invention jointly for explaining the present invention, is not construed as limiting the invention.In the accompanying drawings:

Fig. 1 is that neighborhood calculates schematic diagram;

Fig. 2 is the workflow diagram of road collection amplitude frequency division compensation method in the embodiment of the present invention.

Embodiment

In order to further illustrate object of the present invention, technical scheme and technique effect, first explain principle of work of the present invention below.

When seismic event is propagated in underground medium, owing to being subject to attenuation by absorption impact, therefore the high-frequency energy of seismic event incurs loss, and amplitude, frequency and phase place distort.Although existing Q compensation technique can compensate the amplitude caused because of seismic wave attenuation and phase distortion to a certain extent, but due to reasons such as nearly far away offset distance degree of covering when shortcut goes ground roll, skew are uneven, cause near-energy is more weak a long way, in addition due to can bring in field acquisition and indoor processing procedure the interference of a large amount of background noise and road collection unequal, therefore still cannot the demand of prestack inversion.First the present invention carries out pre-service to original road collection, adopts the beam focusing technology based on in-phase stacking, and the various random disturbance introduced in compacting gatherer process and all kinds of linear disturbance different from effective reflection ray parameter, as ground roll, multiple reflection etc.; Secondly consider that the seismic event of different frequency has different attenuation by absorption features, use for reference the thought of frequency division, scaling down processing is carried out to pretreated road collection, then on the seismic channel set of different frequency section, with the AVO trend just drilled for constraint, parabolical rule is met according to amplitude variation with Offset, AVO Parabolic Fit is carried out in preferred least square method pointwise, again reflection coefficient is asked for, shortcut in different frequency section, distal energy are compensated, finally the road collection of the different frequency section after compensation is merged.Wherein, AVO is the english abbreviation of amplitude variation with Offset.

Introduce correlation theory involved in the present invention below.

1) beam is theoretical

Seismic trace in one, earth's surface place shot point neighborhood and a geophone station neighborhood, no matter how subsurface velocities field changes, and the seismic event that these roads are corresponding all can be propagated near a beam or a central ray, still like this after reflection.According to reflection seismology, if select a Fresnel zone on reflecting interface, the reflection that secondary focal point then in Fresnel zone produces or diffracted wave can be produced the effect of in-phase stacking by same wave detector reception, similarly, received by one group of adjacent wave detector from the reflection of this Fresnel zone or diffracted wave, also can produce the effect of in-phase stacking, Here it is generates the basis of ray beam trace collection.

Suppose that R is the half of Fresnel zone, z is the degree of depth of target reflection horizon to ground, and r is the half of contiguous range, and h is half offset distance.The average velocity shown between inspection surface and target reflecting surface is v hypothetically simultaneously.For seismic event dominant frequency f _m, suppose that seismic reflection within Fresnel zone can in-phase stacking, as shown in Figure 1, fresnel zone reflected waves path difference should be less than that is:

$|\sqrt{z^2+{(h+R+r)}^2}-\sqrt{z^2+{(h-R-r)}^2}|\≤\frac{\mathrm{\λ}}{4}$

$r=R-\frac{\mathrm{\λ}}{2}\sqrt{z^2+h^2-\frac{{\mathrm{\λ}}^2}{4}/{16h}^2-\frac{{\mathrm{\λ}}^2}{4}}$

$r=R-\frac{v}{2f_m}\sqrt{z^2+h^2-\frac{v^2}{4f_m^2}/{16h}^2-\frac{v^2}{4f_m^2}}$

In formula, R is determined by resolution, roughly can calculate contiguous range thus.

2) wavelet transformation frequency division is theoretical

Wavelet transformation has multi-resolution characteristics, localization analysis can be carried out to signal in time domain and frequency field simultaneously, automatically window shape can be adjusted along with the change of frequency content, to meet the demand (namely in low frequency part, there is higher frequency resolution and lower temporal resolution, at HFS, there is lower frequency resolution and higher temporal resolution) of different frequency signals.This is just meeting low frequency signal change slowly and high-frequency signal changes feature rapidly.Continuous wavelet transform is equal to unlike signal employing time domain meticulous gradually or spatial domain sampling step length, can focus on any minor detail of signal.

The theory of wavelet transformation of signal f (t) is, to wavelet basis function { ψ _{σ, τ}(t) | σ, τ ∈ R} carries out translation (τ) and flexible (σ), obtains the wavelet ψ of f (t) _{σ, τ}(t):

${\mathrm{\ψ}}_{\mathrm{\σ},\mathrm{\τ}}\left(t\right)=\frac{1}{\sqrt{\mathrm{\σ}}}\mathrm{\ψ}\left(\frac{t-\mathrm{\τ}}{\mathrm{\σ}}\right)$

Again by ψ _{σ, τ}t () and f (t) do inner product, obtain the wavelet transformation F of f (t) _w(σ, τ):

$F_W(\mathrm{\&sigma;},\mathrm{\&tau;})=<f\left(t\right),{\mathrm{\&psi;}}_{\mathrm{\&sigma;},\mathrm{\&tau;}}>={\&Integral;}_{-\&infin;}^{+\&infin;}f\left(t\right)\frac{1}{\sqrt{\mathrm{\&sigma;}}}\stackrel{\&OverBar;}{\mathrm{\&psi;}}\left(\frac{t-\mathrm{\&tau;}}{\mathrm{\&sigma;}}\right)\mathrm{dt}$

In formula: τ is time shift method; σ is scale factor.Along with the increase of σ, the centre frequency of time frequency window reduces, and time width increases, but area remains unchanged.The synergy of σ and τ can make time frequency window cover whole phase space, be equivalent to go to analyze f (t) frequency content near moment τ, the window function that the f (t) different to frequency content uses width different by a series of continually varying Short Time Fourier Transform.That is, high frequency window time narrow, low frequency window time wide.Wavelet transformation is that a kind of window size is fixed, but shape, time window and frequency window can according to the self-adjusting localization analytical approach of the concrete form of signal.

3) amplitude compensation based on AVO feature is theoretical

Suppose that the incident angle of incident compressional angle and angle of transmission is no more than critical angle and is also not equal to 90 degree, then Zuo Pulizi equation can be reduced to following formula:

R(θ)＝P+Gsin ²θ+C(tan ²θ-sin ²θ) ①

Wherein $P=\frac{1}{2}(\frac{{\mathrm{\&Delta;v}}_p}{v_p}+\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})$

$G=\frac{1}{2}\frac{{\mathrm{\&Delta;v}}_p}{v_p}-2\frac{v_s^2}{v_p^2}(2\frac{{\mathrm{\&Delta;v}}_s}{v_s}+\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})$

$C=\frac{1}{2}\frac{{\mathrm{\&Delta;v}}_p}{v_p}$

Formula 1. in, V _pfor velocity of longitudinal wave; V _sfor shear wave velocity; ρ is rock density; θ is incident angle, and R (θ) is reflection coefficient; be respectively velocity of longitudinal wave, shear wave velocity and rate of change of the density.

Above formula contains little incident angle item, middle incident angle item and large incident angle item three parts, and Shuey carries out rearrangement to trigonometric function item wherein can obtain following formula:

$R\left(\mathrm{\&theta;}\right){\cos }^2\mathrm{\&theta;}\&ap;\frac{1}{2}(\frac{{\mathrm{\&Delta;v}}_p}{v_p}+\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})-(\frac{1}{2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}+4\frac{v_s^2}{v_p^2}\frac{{\mathrm{\&Delta;v}}_s}{v_s}+2\frac{v_s^2}{v_p^2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}){\sin }^2\mathrm{\&theta;}+(4\frac{v_s^2}{v_p^2}\frac{{\mathrm{\&Delta;v}}_s}{v_s}+2\frac{v_s^2}{v_p^2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}){\sin }^4\mathrm{\&theta;}$ ②

2. formula illustrates the relation that reflection R (θ) (i.e. amplitude) changes with incidence angle θ (i.e. geophone offset) well, wherein:

represent reflection coefficient item when zero degree is incident;

reflection coefficient item when angle is incident in expression;

represent reflection coefficient item during large angle incidence.

Suppose to there is following relation:

y＝R(θ)cos ²θ

$R=\frac{1}{2}(\frac{{\mathrm{\&Delta;v}}_p}{v_p}+\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})$

$W=-(\frac{1}{2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}+4\frac{v_s^2}{v_p^2}\frac{{\mathrm{\&Delta;v}}_s}{v_s}+2\frac{v_s^2}{v_p^2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})$

$V=(4\frac{v_s^2}{v_p^2}\frac{{\mathrm{\&Delta;v}}_s}{v_s}+2\frac{v_s^2}{v_p^2}\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}})$

x＝sin ²θ

Then 2. formula can be expressed as following formula:

y≈R+Wx+Vx ²③

3. formula is the common parabolic equation of a class, and Shuey approximate equation is become more directly perceived by this parabolic equation, and can solve with least square principle.

If make error be Δ y, then 3. formula can be expressed as

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2-y_i)}^2=\mathrm{\&Delta;y}$ ④

Formula 4. in, x _i=sin ²θ _i, θ _ii-th incident angle, y _ibe intermediate variable corresponding to i-th incident angle, n is the incident angle number of participating in calculating.Make Δ y=0, then respectively differential carried out to R, W, V, can following formula be obtained:

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i$ ⑤

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)x_i=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{x}_i$ ⑥

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+{\mathrm{Wx}}_i+{\mathrm{Vx}}_i^2)x_i^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{x}_i^2$ ⑦

In reflection R (θ) and incidence angle θ is known, namely can in the hope of x the Sine-squared of incidence angle θ (also), simultaneous 5. 6. 7. three formulas composition system of equations, try to achieve the value of coefficients R, W, V.

In actual treatment, first the process such as denoising, the dynamic school of residue is carried out to seismic channel set, then the value of R, W, V of asking for and x value are substituted into 3. formula try to achieve new y value, more again ask for new reflection R (θ) by y value.In physical record, due to factors such as nearly far away offset distance degree of covering when shortcut goes ground roll, skew are uneven, prestack road collection is caused to be fusiform reflection, closely, energy is weak a long way, show near-energy exception a long way, therefore when asking for R, W, V value, the road collection data can choosing more stable middle angle incidence pointedly carry out the Fitting Calculation, also namely choose the amplitude meeting more stable middle angle incident angle and correspondence that AVO is just drilling trend and carry out the Fitting Calculation.

Generally speaking, the concrete scope of low-angle, middle angle and wide-angle needs to determine according to data situation.The such as ranges of incidence angles of road collection is 0 ~ 30 degree, and so low-angle is exactly 0 ~ 10 degree, and middle angle is exactly 10 ~ 20 degree, and wide-angle is exactly 20 ~ 30 degree.When data is better, incident angle also may reach 40 more and spend.Therefore, in the present invention, the road collection data of so-called more stable middle angle incidence refer to that in data, Amplitudeversusangle situation can meet parabolical rule, and its error is less than certain threshold value.In the present embodiment, the amplitude choosing 10 ~ 25 degree of incident angles and correspondence carries out the Fitting Calculation.

The workflow of road provided by the invention collection vibration amplitude compensation method is explained below in conjunction with an embodiment.As shown in Figure 2, the method mainly comprises the following steps.

S100, pre-treatment step.

No matter be field acquisition or indoor process, all can be attended by a large amount of line noises and random noise.In prior art, be no matter the correlation analysis denoising of time domain, the various filtering and noise reductions of frequency field, although can suppress noise to a certain extent, but have the Noise Elimination poor effect of similar spectrum signature to useful signal.Therefore, in an embodiment of the present invention, the random noise and the line noise different from effective reflection ray parameter that adopt the beam focusing technology based on in-phase stacking to suppress to introduce in acquisition and processing process, as ground roll, multiple reflection etc.Beam is corresponding one section of linear event on seismic channel set, with concentrate in common offset road carry out the most directly perceived.For the offset distance that is determined, first the contiguous range of shot point and geophone station is defined, being extracted by road collection within the scope of this is a beam, then in beam, sign back layer is chosen, when choosing near sign back layer, window carries out superposition and forms standard track, carries out correlation analysis with beam Nei Gedao.If the relevant time difference, more than 1/4 wavelet lengths, rejects.Finally carry out each road meeting in-phase stacking principle to superpose the object reaching the various linear and random noise of compacting.

In addition, the process of dynamic school can also be remained to prestack road collection.In the prior art, the dynamic school of prestack residue relates generally to static correction method and speed adjusting method.The former supposes that the fluctuation of collection lineups is caused by the local velocity's disturbance on seismic raypaths, global speed field therefore can not be adopted to adjust, but the mode by road determination static correction drift value can be adopted to carry out correction error; The latter supposes that the fluctuation of collection lineups is caused by dynamic school velocity error, therefore adopts velocity analysis, recalculates velocity field high order component to correct error.These two kinds of methods operate all more loaded down with trivial details, and treatment cycle is longer.Therefore, in an embodiment of the present invention, on the basis of prestack road collection denoising, adopt the method for axis of reflection automatic tracing, ask for related coefficient with near, middle road a long way to each prestack given in timing window, obtain the time shift amount relative to zero-offset, then progressively sliding window obtains the time shift amount of whole road collection, road collection is evened up, realizes without the residual NMO correction under speed.The method does template with near, middle trace-stacking road, can reduce a long way nmo stretching on the impact in template road, and along time shaft slide, consider the sequence constraint relation that seismic event is forward and backward in time, and process forward and backward signal amplitude to remain unchanged, maintain AVO feature.

In addition, also need original seismic data to change into incident angle gathers during this, because this process belongs to the common practise of this area, therefore do not run business into particular one and state.

S200, scaling down processing step

Seismic data is the comprehensive response of all frequency informations in effective band, conventional compensation method is all compensate on this comprehensively response, and do not consider that the seismic event of different frequency has different attenuation by absorption features, in the embodiment of the present invention, scaling down processing based on wavelet transformation is carried out to the road collection entering denoising, residue moves school, the seismic event that the method is different to frequency uses different window functions, overcome in conventional method and consolidate timing window to the impact of result, road collection seismic channel set being resolved into different frequency composition carries out subsequent treatment.

S300, parameter extracting step.

As previously mentioned, in physical record, due to factors such as nearly far away offset distance degree of covering when shortcut goes ground roll, skew are uneven, cause prestack road collection be fusiform reflection, energy is weak closely, a long way, show near-energy is extremely a long way.Therefore, overall single compensation can be carried out to seismic data in Normal practice.But in an embodiment of the present invention, on the seismic channel set of different frequency, preferred intermediate stable reliable road collection data calculate respectively targetedly.In addition, the AVO forward simulation based on well data is also carried out before this, to consider the constraint just drilling AVO characteristic trend when choosing interim orbit collection.That is, just to drill AVO characteristic trend for constraint, from the prestack road set information of different frequency, the incident angle of middle angle and the amplitude of correspondence is chosen.

S400, the Fitting Calculation step.

Just as previously mentioned, in physical record, because when shortcut goes ground roll, skew, nearly far away offset distance degree of covering is uneven, prestack road collection be that fusiform reflects, and closely, a long way energy is weak.For revising fusiform road collection reflection-closely, energy is weak a long way problem, increase the fidelity of road collection, it is made to meet the requirement of prestack inversion, have based on stable parabola characteristic by road collection amplitude variation with Offset, AVO Parabolic Fit is carried out to the incident angle of different frequency road collection extracted and amplitude, such as, adopts least square method pointwise to carry out AVO Parabolic Fit.Ask for R, W, V value by matching, and 3. their substitution formulas are tried to achieve new y value, then recalculate reflection R (θ) according to y value.

S500, energy compensating step.

According to R, W, V value of the Fitting Calculation, 3. they are substituted into formula, recalculate the reflection R (θ) when different incidence angles is concentrated in different frequency road, the road weak to shortcut, distal energy according to the reflection coefficient recalculated carries out frequency division amplitude compensation, finally the road collection that frequency division compensates is merged, obtain the road collection of earthquake full frequency band.

Last it is noted that above each embodiment is only in order to illustrate technical scheme of the present invention, be not intended to limit; Although with reference to foregoing embodiments to invention has been detailed description, those of ordinary skill in the art is to be understood that: it still can be modified to the technical scheme described in foregoing embodiments, or carries out equivalent replacement to wherein some or all of technical characteristic; And these amendments or replacement, do not make the essence of appropriate technical solution depart from the protection domain of application claims.

Claims (9)

1., towards a road collection amplitude frequency division compensation method for prestack inversion, it is characterized in that, comprise

Pre-treatment step, comprises denoising, the dynamic school of residue to original prestack road set information and proceeds to the pre-service of penetrating angle gathers;

Scaling down processing step, carries out scaling down processing to pretreated prestack road set information;

Parameter extracting step, extracts to meet and just drills the incident angle of AVO characteristic trend and the amplitude of correspondence from the prestack road set information of different frequency;

The Fitting Calculation step, carries out AVO Parabolic Fit to the incident angle of different frequency road collection extracted and amplitude, calculates each coefficient of parabolic equation;

Energy compensating step, the parabolic equation obtained according to matching redefines reflection coefficient when different incidence angles is concentrated in different frequency road, carries out frequency division amplitude compensation according to the reflection coefficient redefined to shortcut and energy is weak a long way road.

2. road collection amplitude frequency division compensation method as claimed in claim 1, it is characterized in that, in described pre-treatment step, denoising is carried out in the beam focusing based on in-phase stacking.

3. road collection amplitude frequency division compensation method as claimed in claim 2, it is characterized in that, the denoising based on the beam focusing of in-phase stacking comprises the following steps:

For the offset distance that is determined, first define the contiguous range of shot point and geophone station, being extracted by the road collection within the scope of this is a beam;

In beam, choose sign back layer, when choosing near sign back layer, window carries out superposition and forms standard track, carries out correlation analysis with beam Nei Gedao;

Time window within the scope of, the road collection of the relevant time difference more than 1/4 wavelet lengths is rejected;

The each road meeting in-phase stacking principle is superposed.

4. the road collection amplitude frequency division compensation method as described in claims 1 to 3 any one, is characterized in that:

It is the method adopting axis of reflection automatic tracing that described residue moves school process, with near, middle road, related coefficient is asked for a long way to each prestack given in timing window, obtain the time shift amount relative to zero-offset, then progressively sliding window obtains the time shift amount of whole road collection, and Jiang Daoji evens up.

5. the road collection amplitude frequency division compensation method as described in claims 1 to 3 any one, is characterized in that:

Described scaling down processing step, carries out scaling down processing based on wavelet transformation to pretreated prestack road set information.

6. the road collection amplitude frequency division compensation method as described in claims 1 to 3 any one, is characterized in that, further,

In described parameter extracting step, extract to meet from the prestack road set information of different frequency and just drill the incident angle of angle and the amplitude of correspondence in AVO characteristic trend.

7. road collection amplitude frequency division compensation method as claimed in claim 6, is characterized in that:

The scope of the incident angle of described middle angle is 10 ~ 25 degree.

8. road collection amplitude frequency division compensation method as claimed in claim 1, is characterized in that, in described the Fitting Calculation step, carries out AVO Parabolic Fit, by coefficients R, W, V of following formula determination parabolic equation according to the incident angle extracted and amplitude;

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+W{x}_i+V{x}_i^2)=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i$

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+W{x}_i+V{x}_i^2)x_i=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{x}_i$

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(R+W{x}_i+V{x}_i^2)x_i^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{x}_i^2$

In formula, x _i=sin ²θ _i, θ _ii-th incident angle; y _iintermediate variable corresponding to i-th incident angle.

9. road collection amplitude frequency division compensation method as claimed in claim 8, is characterized in that:

In the Fitting Calculation step, least square method is adopted to calculate each coefficient of parabolic equation.