CN103278849B - The method and system of higher-order spectra are carried out based on seismic data and well-log information - Google Patents

Abstract

The present invention provides a kind of method and system carrying out higher-order spectra based on seismic data and well-log information, and described method includes: gather the seismic data of current formation, well-log information；Minimum phase wavelet is estimated according to described seismic data；Carry out 01 codings according to described minimum phase wavelet, obtain multiple coded data；Well lie data are extracted according to described seismic data；Reflection coefficient is determined according to described well-log information；Well control evaluation model is determined according to described well lie data, reflection coefficient and described convolution model；Optimal wavelet is estimated according to described well control evaluation model and described multiple coded data；Carry out reservoir prediction according to described optimal wavelet, obtain the reservoir prediction result of current formation.Use genetic manipulation and combine well control as interpretational criteria to process sub-wave phase problem, efficiently reducing amount of calculation, accelerate the convergence rate of root varitron wave phase Power estimation, improve higher-order spectra precision.

Description

The method and system of higher-order spectra are carried out based on seismic data and well-log information

Technical field

The present invention is about geophysical prospecting for oil technical field, especially with regard to reservoir prediction based on higher-order spectra Technology, is concretely a kind of method and system carrying out higher-order spectra based on seismic data and well-log information.

Background technology

In geophysical prospecting for oil technical field, amplitude with angle change (Amplitude versus angle, AVA) inverting be utilize reflection seismic wave-amplitude with the Changing Pattern inverting subsurface lithologic parameter of angle of incidence, and then estimation reservoir rock Stone and the character of fluid also carry out the technology of petroleum-gas prediction.Prestack seismic data is frequently utilized that the amplitude change with offset distance (Amplitude Versus Offset, AVO) carries out inverting to stratum three parameter (velocity of longitudinal wave, shear wave velocity and density), from And obtain for identifying the basic lithological information such as the velocity of longitudinal wave of hydrocarbon-bearing formation, shear wave velocity and density.

Inversion technique of the prior art (all two kinds described above), when carrying out reservoir prediction, is required to be applied to convolution Model, and the seismic wavelet in convolution model is whether accurate, will directly affect the accuracy of follow-up reservoir prediction.In reality In seismic wavelet estimation, due to the impact of the factors such as noise, the phase spectrum of wavelet is difficult to accurately estimate.Based on inaccurate phason Ripple carries out deconvolution, and the resolution of its result will be somewhat limited.

In prior art, the estimation about seismic wavelet mainly has following several ways: (1), Ziolkowski et al. carry The method that geological data is carried out transform gone out；(2), the auto-correlation to geological data that Lindsey proposes carries out transform Method；(3), the method for Ursin proposition is to solve a series of extension Yule-Walker equation group, and passes through only to mixed The maximum phase component closing phase place wavelet carries out transform rooting to estimate optimal mixed-phase inverse filter；(4), utilization piece Act method, Li Lihong and Ma Zhaogui proposes the root to minimum phase wavelet transform and is scanned estimating optimal wavelet；(5), normal Phase place rotates wavelet estimation；(6), root conversion higher-order spectra mode.

The technological deficiency that in prior art, above-mentioned several higher-order spectra methods are primarily present is: (1), Ziolkowski side The higher-order spectra precision of method is the highest, and is vulnerable to effect of noise；(2), Lindsey method utilizes auto-correlation information to estimate son Ripple, it is difficult to control the phase spectrum of wavelet；(3), the precision and stability of Ursin method is easily subject to time window length and operator length Impact；(4), enumerative technique amount of calculation is huge, it is difficult to be applied to real data；(5), Chang Xiangwei rotation wavelet estimation can lead to The mode crossing phasescan realizes, and its amount of calculation is the most little, but application conditions is harsher；(6), root conversion higher-order spectra Although mode compares the phase condition meeting actual wavelet, but amount of calculation is huge, it is difficult to be applied in real data process.

Therefore, work area is explored for actual seismic, how accurately to estimate seismic wavelet according to different geological datas, and then Inverting multiple subsurface lithologic parameter, carrying out petroleum-gas prediction is this area technical barrier urgently to be resolved hurrily.

Summary of the invention

For the problems referred to above overcoming prior art to exist, the invention provides a kind of based on seismic data and well-log information Carrying out the method and system of higher-order spectra, by wavelet root conversion process being carried out specific coding and decoding, and introducing heredity plan Slightly, use genetic manipulation and combine well control as interpretational criteria to process sub-wave phase problem, efficiently reducing amount of calculation, add The fast convergence rate of root varitron wave phase Power estimation, improves higher-order spectra precision.

It is an object of the invention to provide a kind of method carrying out higher-order spectra based on seismic data and well-log information, Including: gather the seismic data of current formation, well-log information；Minimum phase wavelet is estimated according to described seismic data；According to Described minimum phase wavelet carries out 0-1 coding, obtains multiple coded data；Well lie number is extracted according to described seismic data According to；Reflection coefficient is determined according to described well-log information；According to described well lie data, reflection coefficient and described convolution Model determines well control evaluation model；Optimal wavelet is estimated according to described well control evaluation model and described multiple coded data； Carry out reservoir prediction according to described optimal wavelet, obtain the reservoir prediction result of current formation.

It is an object of the invention to provide and a kind of based on what seismic data and well-log information carried out higher-order spectra be System, including: harvester, for gathering the seismic data of current formation, well-log information；Minimum phase wavelet determines device, uses In estimating minimum phase wavelet according to described seismic data；Code device, for carrying out according to described minimum phase wavelet 0-1 encodes, and obtains multiple coded data；Well lie data extraction device, for extracting well lie according to described seismic data Data；Reflection coefficient determines device, for determining reflection coefficient according to described well-log information；Well control evaluation model determines dress Put, for determining well control evaluation model according to described well lie data, reflection coefficient and described convolution model；Optimal son Ripple estimation unit, for estimating optimal wavelet according to described well control evaluation model and described multiple coded data；Reservoir is pre- Survey device, for carrying out reservoir prediction according to described optimal wavelet, obtain the reservoir prediction result of current formation.

The beneficial effects of the present invention is, by wavelet root conversion process is carried out specific coding and decoding, and introduce something lost Pass strategy, use genetic manipulation and combine well control as interpretational criteria to process sub-wave phase problem, efficiently reducing calculating Amount, accelerates the convergence rate of root varitron wave phase Power estimation, improves higher-order spectra precision, solve of the prior art When root converter technique asks for mixed phase wavelet, sub-wave phase is had the technical problem of a certain degree of dependence.

For the above and other objects, features and advantages of the present invention can be become apparent, preferred embodiment cited below particularly, And coordinate institute's accompanying drawings, it is described in detail below.

Accompanying drawing explanation

In order to be illustrated more clearly that the embodiment of the present invention or technical scheme of the prior art, below will be to embodiment or existing In having technology to describe, the required accompanying drawing used is briefly described, it should be apparent that, the accompanying drawing in describing below is only this Some embodiments of invention, for those of ordinary skill in the art, on the premise of not paying creative work, also may be used To obtain other accompanying drawing according to these accompanying drawings.

In Fig. 1, (a), for without wavelet composite traces figure initial in the case of making an uproar, (b) is for without protocol figure, (c) in the case of making an uproar For the wavelet composite traces figure estimated without the present invention in the case of making an uproar；

Fig. 2 is the comparison diagram estimating wavelet without wavelet initial in the case of making an uproar, true wavelet and the present invention；

In Fig. 3, (a) be noisy in the case of initial wavelet composite traces figure, (b) be noisy in the case of protocol figure, (c) Wavelet composite traces figure is estimated for the present invention in the case of noisy；

Fig. 4 be noisy in the case of initial wavelet, true wavelet and the present invention estimate the comparison diagram of wavelet；

In Fig. 5, (a) is Chang Xiangwei wavelet composite traces figure, and (b) is well lie real data figure, and (c) is that the present invention estimates Wavelet composite traces figure；

Fig. 6 (a) is the wavelet figure that the present invention estimates；

Fig. 6 (b) is that Chang Xiangwei rotates the wavelet figure estimated；

Fig. 7 (a) is original seismic data figure；

Fig. 7 (b) for estimating the deterministic deconvolution result figure of wavelet based on the present invention；

Fig. 7 (c) is Chang Xiangwei wavelet deconvolution result figure；

A kind of method carrying out higher-order spectra based on seismic data and well-log information that Fig. 8 provides for the embodiment of the present invention Flow chart；

Fig. 9 is the particular flow sheet of step S103 in Fig. 8；

Figure 10 is the particular flow sheet of step S107 in Fig. 8；

Figure 11 is the particular flow sheet of step S303 in Figure 10；

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 12 provides for the embodiment of the present invention Structured flowchart；

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 13 provides for the embodiment of the present invention The concrete structure block diagram of middle code device；

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 14 provides for the embodiment of the present invention In the concrete structure block diagram of optimal higher-order spectra device；

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 15 provides for the embodiment of the present invention The concrete structure block diagram of middle decoder module；

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 16 provides for the embodiment of the present invention The concrete structure figure of middle evaluation module；

Figure 17 is the schematic diagram that minimum phase wavelet carries out 0-1 coding and decoding.

Detailed description of the invention

Below in conjunction with 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 Describe, it is clear that described embodiment is only a part of embodiment of the present invention rather than whole embodiments wholely.Based on Embodiment in the present invention, it is every other that those of ordinary skill in the art are obtained under not making creative work premise Embodiment, broadly falls into the scope of protection of the invention.

The present invention proposes a kind of method and system carrying out higher-order spectra based on seismic data and well-log information, passes through Wavelet root conversion process is carried out specific coding and decoding, and introduces genetic algorithm, use Genetic Strategies and combine well control conduct Interpretational criteria processes sub-wave phase problem, efficiently reduces amount of calculation, accelerates the receipts of root varitron wave phase Power estimation Hold back speed, improve higher-order spectra precision, reduce the purpose of prediction of oil-gas reserve risk.

A kind of method carrying out higher-order spectra based on seismic data and well-log information that Fig. 8 provides for the embodiment of the present invention Flow chart, as shown in Figure 8, the method specifically includes:

S101: gather the seismic data of current formation, well-log information.

In oil-gas exploration, seismic data, well-log information remain abundant information, how to use such information for obtaining Take Oil/gas Well parameter and for oil and gas detection and identify that providing data supporting is one of of paramount importance work of seismic prospecting.

S102: estimate minimum phase wavelet according to described seismic data.This step can be realized by prior art, all Such as spike deconvolution.

S103: carry out 0-1 coding according to described minimum phase wavelet, obtains multiple coded data.

Fig. 9 is the particular flow sheet of step S103 in Fig. 8, and as shown in Figure 9, this step specifically includes:

S201: described minimum phase wavelet is carried out transform, obtains the multinomial about variable z.

It is determined by temper ripple method of estimation or a minimum phase wavelet table that statistically higher-order spectra method estimates It is shown as w (t)=(w₀,w₁,…,w_n), n is the sampling number of seismic wavelet, wherein w₀w_n≠0.This wavelet sequence is carried out transform Obtain one only about the multinomial of variable z, be shown below:

Z(w(t))=w_n(z-α₁)(z-α₂)…(z-α_i)…(z-α_n) (1)

S202: solve the root that described multinomial is corresponding to determine described multinomial.Make (1) formula equal to 0, use multinomial It is (α that decomposition can obtain this root of polynomial₁,α₂,…,α_i,…,α_n)。

S203: the root that described multinomial is corresponding is carried out 0-1 coding, obtains multiple coded data.In order to be further described The cataloged procedure of the present invention, the most first sketches the coding principle of the present invention.

Assume that seismic data can be by Robinson(Lu Binxun) convolution model describes, then has

s(t)=w(t)*r(t)=Z^-1[Z(w(t))]*r(t) (2)

=w_n(-α₁,1)*(-α₂,1)*…*(-α_i,1)*…*(-α_n,1)*(r₀,r₁,…,r_m)

Wherein, t express time, s (t) represents earthquake record, w (t)=(w₀,w₁,…,w_n) expression seismic wavelet, r (t)= (r₀,r₁,…,r_m) representing reflection coefficient, m is the sampling number of reflection coefficient, Z^-1Represent Z anti transform, (α₁,α₂,…,α_i,…, α_n) represent seismic wavelet transform root, (-α_i, 1) represent w (t) 2 wavelets.Inventionbroadly, most of earthquake deconvolution Or Optimum Impedance Inversion Method is exactly based on what formula (2) realized.

2 wavelet (-α to w (t)_i, 1) and carry out Fourier transformation, and make α_i=x_i+jy_i, have

F[(-α_i,1)]=-α_i+e^-jwt=-(x_i+jy_i)+cos(wt)-jsin(wt) (3)

=[cos(wt)-x_i]-j[sin(wt)+y_i]

Wherein,Representing imaginary unit, w represents angular frequency.Formula (3) is taken amplitude spectrum and phase spectrum respectively, Have

$A\left[({-\mathrm{\α}}_i,1)\right]=\sqrt{{[\cos \left(\mathrm{wt}\right)-x_i]}^2+{[\sin \left(\mathrm{wt}\right)+y_i]}^2}$ (4)

$={[1+x_i^2+y_i^2-{2x}_i\cos \left(\mathrm{wt}\right)+{2y}_i\sin \left(\mathrm{wt}\right)]}^{1/2}$

$P\left[(-{\mathrm{\α}}_i,1)\right]=-\arctan \left[\frac{\sin \left(\mathrm{wt}\right)+y_i}{\cos \left(\mathrm{wt}\right)-x_i}\right]---\left(5\right)$

Due to w (t)=w_n(-α₁,1)*(-α₂,1)*…*(-α_i,1)*…*(-α_n, 1), so the amplitude spectrum of wavelet and phase place Spectrum can be expressed as

$A\left[w\left(t\right)\right]=\underset{i=1}{\overset{n}{\mathrm{\Π}}}A\left[({-\mathrm{\α}}_i,1)\right]=\underset{i=1}{\overset{n}{\mathrm{\Π}}}[1+x_i^2+y_i^2-{2x}_i\cos \left(\mathrm{wt}\right)+{2y}_i\sin \left(\mathrm{wt}\right){]}^{1/2}---\left(6\right)$

$P\left[w\left(t\right)\right]=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}P\left[({-\mathrm{\α}}_i,1)\right]=-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\arctan \left[\frac{\sin \left(\mathrm{wt}\right)+y_i}{\cos \left(\mathrm{wt}\right)-x_i}\right]---\left(7\right)$

2 wavelet (-α as w (t)_i, 1) and quiltAfter replacement, amplitude spectrum does not change, this is because

$A\left[\right|{\stackrel{\‾}{\mathrm{\α}}}_i\left|(-\frac{1}{{\stackrel{\‾}{\mathrm{\α}}}_i},1)\right]=\sqrt{x_i^2+y_i^2}\sqrt{{[\cos \left(\mathrm{wt}\right)-\frac{x_i}{x_i^2+y_i^2}]}^2+{[\sin \left(\mathrm{wt}\right)+\frac{y_i}{x_i^2+y_i^2}]}^2}---\left(8\right)$

$=[1+x_i^2+y_i^2-{2x}_i\cos \left(\mathrm{wt}\right)+{2y}_i\sin \left(\mathrm{wt}\right){]}^{1/2}=A[({-\mathrm{\α}}_i,1)]$

Wherein, x_iy_i≠0.2 waveletsPhase spectrum can be expressed as

$P\left[\right|{\stackrel{\‾}{\mathrm{\α}}}_i\left|(-\frac{1}{{\stackrel{\‾}{\mathrm{\α}}}_i},1)\right]=-\mathrm{arctyan}\left[\frac{\sin \left(\mathrm{wt}\right)+\frac{y_i}{x_i^2+y_i^2}}{\cos \left(\mathrm{wt}\right)-\frac{x_i}{x_i^2+y_i^2}}\right]---\left(9\right)$

It is apparent that work asI.e.Time,(-α_i, 1) phase spectrum just the same.This Time, (-α_i, 1) andSubstantially represent same binary wavelet.WhenTime,Therefore, by only convertingBinary wavelet (-α_i, 1) beJust can scan all equal length, the out of phase music score ripple of same-amplitude spectrum.Note working as α_iIt is multiple Number, i.e. y_iWhen ≠ 0, its Conjugate complex roots needs to convert simultaneously.

It follows that by moving part or allRoot α_i(i=1,2 ..., n) be, can scan Series of identical length, same-amplitude spectrum, the seismic wavelet of out of phase spectrum.Substantially, root α_iIt is moved intoIt is one to one Process.Root is from unit circleIn moving to unit circle at a certain defined location, or in unit circleMove at the outer a certain defined location of unit circle.

Binary the wavelet (-α of conversion w (t)_i, 1) beAfter, the amplitude spectrum of wavelet does not change, and is Phase spectrum there occurs change, therefore the scanning of sub-wave phase be circle outside root in circle some fixed position or circle in root to The process that circle some fixed position outer is moved.Based on this feature, the present invention proposes a kind of 0-1 coded method to complete this Plant and move root mode.

Figure 17 is the schematic diagram that minimum phase wavelet carries out 0-1 coding and decoding.As shown in Figure 17, chromosome is by necessarily The 0 of quantity or 1 is constituted, and the front end of chromosome represents real number gene, only maps real root；Tail end represents plural number gene, only maps multiple Root.The root that this gene position is mapped is not moved in 0 expression, and 1 represents the root that this gene position mobile is mapped.The length of chromosome only with Initial wavelet w (t) estimated is relevant, and after w (t) determines, the length of chromosome is fixed the most therewith.Additionally, plural number gene string is also Deterministically control another plural gene string (present invention is referred to as controlled plural number gene string).At controlled plural number In gene string, 0 ' represents the Conjugate complex roots not moving with moving the root that corresponding plural gene position is mapped respectively with 1 '.So Coded system us just can be made to directly utilize the selection of genetic algorithm, intersect and make a variation and update the phase spectrum of wavelet.

According to above-mentioned cataloged procedure, to minimum phase wavelet w (t)=(w₀,w₁,…,w_n) carry out 0-1 coding after, can be formed 2ⁿPlant coded data.

As shown in Figure 8, what the present invention provided carry out the method for higher-order spectra based on seismic data and well-log information also includes:

S104: extract well lie data according to described seismic data.

S105: determine reflection coefficient according to described well-log information.When first this step extracts sound wave from well-log information Difference curve and density curve, then calculate reflection coefficient according to interval transit time curve and density curve.

S106: determine well control evaluation model according to described well lie data, reflection coefficient and described convolution model. Evaluate the quality between individuality, need the size by fitness value to weigh.The present invention, using well as controlling standard, uses and closes The correlation coefficient of one-tenth earthquake record w (t) * r (t) and well lie earthquake record is as well control evaluation model, i.e.

O=max[xcorr(w(t)*r(t),s(t))] (10)

Wherein, max [] expression takes maximum, and s (t) is well lie data, and w (t) is the seismic wavelet that decoding obtains, r (t) For the reflection coefficient determined according to well-log information, xcorr represents between composite traces w (t) * r (t) and well lie record s (t) Correlation coefficient, O is object function.

S107: estimate optimal wavelet according to described well control evaluation model and described multiple coded data.Figure 10 is step The particular flow sheet of rapid S107, as shown in Figure 10, this step specifically includes:

S301: choose the coded data of predetermined number from described multiple coded data.Operation due to genetic algorithm To liking a colony, the coded data individual (chromosome) therefore randomly generating predetermined number forms population.The rule of population Mould is exactly number individual in population, is typically taken as even number.Each individuality represents an initial wavelet, and population represents one initially Wavelet collection.Concentrating in initial wavelet, each wavelet has identical length, identical amplitude spectrum and different phase spectrums.? In specific embodiment, predetermined number can preset, the even number being such as set to by predetermined number in 20-100 scope.

S302: be decoded the coded data of described predetermined number, obtains the wavelet sequence of predetermined number.This step has Body includes:

S401: described coded data is carried out 0-1 decoding, obtains the root of correspondence.It is described further below the present invention's Decoding process.According to described binary wavelet transform, become described wavelet multinomial the 0-1 coding mapping of coded data genetic algorithm The root of formula.

S402: ask for the multinomial about variable z according to described root.I.e. according to root (α₁,α₂,…,α_i,…,α_n) obtain About the multinomial of variable z, the multinomial solved such as formula (1).

S403: the multinomial of described variable z is carried out Z anti transform, obtains the wavelet sequence that described coded data is corresponding.Right Multinomial about variable z carries out Z anti transform, and the wavelet sequence obtained is w (t)=(w₀,w₁,…,w_n), n is seismic wavelet Sampling number, wherein w₀w_n≠0。

The i.e. decoding of basic genetic algorithmic is by { 0,1} binary string is transformed into decimal number, and then linear transformation becomes again One real number.And the decoding of the present invention is by { 0,1} binary string is transformed into n root (including compound radical and real root), then passes through Seek multinomial and take coefficient and obtain a sub-wave train.Assume that certain individuality after being decoded isThe individuality then decoded Have been assigned clear and definite physical significance.

As shown in Figure 10, S107 step also includes:

S303: whether exist optimal in the wavelet sequence according to the predetermined number described in described well control evaluation model judgement Wavelet sequence.Figure 11 is the particular flow sheet of step S303 in Figure 10, and as shown in Figure 11, step S303 is especially by including:

S401: obtain well control threshold coefficient set in advance.In a particular embodiment, well control threshold value can be preset Coefficient, such as 80%.

S402: described wavelet sequence carries out convolution with corresponding reflection coefficient, forms synthetic seismogram；

S403: determine the phase relation of described synthetic seismogram and well lie record according to described well control evaluation model Number；

Whether S404: the correlation coefficient described in judgement is more than or equal to described well control threshold coefficient.

As shown in Figure 10, S107 step also includes:

S304: when being judged as YES, the optimal wavelet sequence described in output, namely the maximum wavelet sequence of correlation coefficient i.e. For optimal wavelet sequence.

S305: when being judged as NO, the phase place of the wavelet described in renewal.Return according to the sub-wave phase after updating and perform step Rapid S301, chooses the coded data of new a collection of predetermined number, re-executes step S302 to S303, to choose optimal wavelet.

Namely utilize herein binary strings genetic algorithm to update the phase place of wavelet, and by making by synthetic seismogram and well Seismic channel reaches most preferably to mate to obtain optimal higher-order spectra.

When updating the phase place of wavelet, genetic manipulation can be used to update sub-wave phase, genetic manipulation generally comprises selection, intersects And variation, its main purpose is to transform individuality.In the present invention, their effect is the phase place updating wavelet.

Wherein " select " to be to select excellent wavelet, under making them have an opportunity to be to concentrate from current wavelet as parent Generation breeding descendants.Ideal adaptation degree is the highest, and its selected chance is the most.Selection is the key of genetic algorithm, and it embodies The thought of survival of the fittest in nature.

" intersect " is the part-structure that two fathers are individual to be replaced restructuring and generates new individual operation.This operation Purpose is able to produce new individuality the next generation, just as the marriage process of human society, is operated by intersection, loses Propagation algorithm search equal length, same-amplitude spectrum, the ability of out of phase music score ripple will be greatly improved." intersecting " is heredity Algorithm obtains optimal wavelet and the most important means of reflection coefficient.

" make a variation " is that the gene position to individual carries out inversion operation.By mutation operation, genetic algorithm can be made to keep planting The multiformity of group, to prevent premature convergence.It is said that in general, mutation probability all obtains less, the most even take 0.

Select, intersect and make a variation and be respectively arranged with the most different methods.For the particular problem of the present invention, in detailed description of the invention In, roulette wheel selection, binary system uniform crossover method and binary system alternative method can be chosen respectively.

As shown in Figure 8, a kind of of embodiment of the present invention offer carries out higher-order spectra based on seismic data and well-log information Method also includes:

S108: carry out reservoir prediction according to described optimal wavelet, obtains the reservoir prediction result of current formation.

When carrying out reservoir prediction, be required to be applied to convolution model, and the seismic wavelet in convolution model accurately with No, the accuracy of follow-up reservoir prediction will be directly affected.In actual seismic higher-order spectra, due to the shadow of the factors such as noise Ringing, the phase spectrum of wavelet is difficult to accurately estimate.Carrying out deconvolution based on inaccurate phase place wavelet, the resolution of its result will It is somewhat limited.Amplitude is such as utilized to carry out with angle change (Amplitude versus angle, AVA) inverting Petroleum-gas prediction, such as utilizes amplitude with the change (Amplitude Versus Offset, AVO) of offset distance to stratum three parameter (velocity of longitudinal wave, shear wave velocity and density) carries out the identification of hydrocarbon-bearing formation.Therefore, the method that the present invention proposes is directly to estimating The wavelet of meter carries out transform rooting, then by the phase place of binary system genetic manipulation scanning wavelet, and then evolves rapidly Good wavelet, improves higher-order spectra precision, and then improves the precision of reservoir prediction.

A kind of system carrying out higher-order spectra based on seismic data and well-log information that Figure 12 provides for the embodiment of the present invention Structured flowchart, as shown in Figure 12, shown system specifically includes:

Harvester 100, for gathering the seismic data of current formation, well-log information.In oil-gas exploration, earthquake provides Remaining abundant information in material, well-log information, how using such information for obtaining Oil/gas Well parameter is oil and gas detection and knowledge You can well imagine for data supporting is one of of paramount importance work of seismic prospecting.

Minimum phase wavelet determines device 200, for estimating minimum phase wavelet according to described seismic data.This device Can be realized by prior art, such as spike deconvolution.

Code device 300, for carrying out 0-1 coding according to described minimum phase wavelet, obtains multiple coded data.Figure 13 is the concrete structure block diagram of code device 300, and as shown in Figure 13, code device specifically includes:

Transform unit 301, for described minimum phase wavelet is carried out transform, obtains the multinomial about variable z.

It is determined by temper ripple method of estimation or a minimum phase wavelet table that statistically higher-order spectra method estimates It is shown as w (t)=(w₀,w₁,…,w_n), n is the sampling number of seismic wavelet, wherein w₀w_n≠0.This wavelet sequence is carried out transform Obtain a multinomial only about variable z, as shown in Equation (1).

Z(w(t))=w_n(z-α₁)(z-α₂)…(z-α_i)…(z-α_n) (1)

Polynomial solving unit 302, for solving the root that described multinomial is corresponding to determine described multinomial.Order is described (1) formula equal to 0, using Factoring Polynomials can obtain this root of polynomial is (α₁,α₂,…,α_i,…,α_n)。

Coding unit 303, for the root that described multinomial is corresponding is carried out 0-1 coding, obtains multiple coded data.Conversion Binary the wavelet (-α of w (t)_i, 1) beAfter, the amplitude spectrum of wavelet does not change, and simply phase spectrum there occurs Change, the root outside therefore the scanning of sub-wave phase is circle in circle the root in some fixed position or circle outside circle, some is solid The process of movement is put in location.Based on this feature, the present invention proposes a kind of 0-1 coded method to complete this shifting root mode.

Figure 17 is the schematic diagram that minimum phase wavelet carries out 0-1 coding and decoding.As shown in Figure 17, chromosome is by necessarily The 0 of quantity or 1 is constituted, and the front end of chromosome represents real number gene, only maps real root；Tail end represents plural number gene, only maps multiple Root.The root that this gene position is mapped is not moved in 0 expression, and 1 represents the root that this gene position mobile is mapped.The length of chromosome only with Initial wavelet w (t) estimated is relevant, and after w (t) determines, the length of chromosome is fixed the most therewith.Additionally, plural number gene string is also Deterministically control another plural gene string (present invention is referred to as controlled plural number gene string).At controlled plural number In gene string, 0 ' represents the Conjugate complex roots not moving with moving the root that corresponding plural gene position is mapped respectively with 1 '.So Coded system us just can be made to directly utilize the selection of genetic algorithm, intersect and make a variation and update the phase spectrum of wavelet.

According to above-mentioned cataloged procedure, to minimum phase wavelet w (t)=(w₀,w₁,…,w_n) carry out 0-1 coding after, just obtain 2ⁿPlant coded data.

As shown in Figure 12, a kind of of present invention offer carries out the system of higher-order spectra also based on seismic data and well-log information Including:

Well lie data extraction device 400, for extracting well lie data according to described seismic data；

Reflection coefficient determines device 500, for determining reflection coefficient according to described well-log information.This device is first from survey Well data is extracted interval transit time curve and density curve, then calculates instead according to interval transit time curve and density curve Penetrate coefficient.

Well control evaluation model determines device 600, for according to described well lie data, reflection coefficient and described pleat Product module type determines well control evaluation model.Evaluate the quality between individuality, need the size by fitness value to weigh.The present invention Using well as controlling standard, use synthetic seismogram w (t) * r (t) correlation coefficient with well lie earthquake record as well control Evaluation model, i.e.

O=max[xcorr(w(t)*r(t),s(t))]

Wherein, max [] expression takes maximum, and s (t) is well lie data, and w (t) is decoding seismic wavelet out, r (t) The reflection coefficient calculated for log, xcorr represents between composite traces w (t) * r (t) and well lie record s (t) Correlation coefficient, O is object function.

Optimal higher-order spectra device 700, for estimating according to described well control evaluation model and described multiple coded data Counting optimal wavelet, Figure 14 is the concrete structure block diagram of optimal higher-order spectra device, and as shown in Figure 14, optimal higher-order spectra device has Body includes:

Coded data acquisition module 701, for choosing the coded data of predetermined number from described multiple coded data. Owing to the operation of genetic algorithm is to liking a colony, the coded data individual (chromosome) therefore randomly generating predetermined number is come Composition population.The scale of population is exactly number individual in population, is typically taken as even number.Each individuality represents an initial son Ripple, population represents an initial wavelet collection.Concentrating in initial wavelet, each wavelet has identical length, identical amplitude Spectrum and different phase spectrums.In a particular embodiment, predetermined number can preset, and is such as set to by predetermined number Even number in 20-100 scope.

Decoder module 702, for being decoded the coded data of described predetermined number, obtains the wavelet sequence of predetermined number Row.The concrete structure block diagram of decoder module as shown in figure 15, specifically includes:

Decoding unit 7021, for described coded data is carried out 0-1 decoding, obtains the root of correspondence.Further below Introduce the decoding process of the present invention.According to described binary wavelet transform, coded data is encoded anti-according to the 0-1 of genetic algorithm It is mapped to described wavelet root of polynomial.

Multinomial asks for unit 7022, for asking for the multinomial about variable z according to described root.I.e. according to root (α₁, α₂,…,α_i,…,α_n) obtain the multinomial about variable z, the multinomial solved is as shown in Equation (1).

Z anti transform unit 7023, for the multinomial of described variable z is carried out Z anti transform, obtains described coded data pair The wavelet sequence answered.Multinomial about variable z is carried out Z anti transform, and the wavelet sequence obtained is w (t)=(w₀,w₁,…, w_n), n is the sampling number of seismic wavelet, wherein w₀w_n≠0。

The i.e. decoding of basic genetic algorithmic is by { 0,1} binary string is transformed into decimal number, and then linear transformation becomes again One real number.And the decoding of the present invention is by { 0,1} binary string is transformed into n root (including compound radical and real root), then passes through Seek multinomial and take coefficient and obtain a sub-wave train.Assume that certain individuality after being decoded is, then the individuality quilt that decodes Impart clear and definite physical significance.

As shown in Figure 14, optimal higher-order spectra device also includes:

Evaluation module 703, in the wavelet sequence of the predetermined number described in judging according to described well control evaluation model Whether there is optimal wavelet sequence, Figure 16 is the concrete structure figure of evaluation module, and as shown in Figure 16, evaluation module specifically includes:

Threshold value acquiring unit 7031, is used for obtaining well control threshold coefficient set in advance；

Convolution unit 7032, for described wavelet sequence carries out convolution with corresponding reflection coefficient, is formed synthetically Shake record；

Correlation coefficient determines unit 7033, for determining described synthetic seismogram according to described well control evaluation model Correlation coefficient with well lie record；

Judging unit 7034, for judging that whether described correlation coefficient is more than or equal to described well control threshold coefficient.

As shown in Figure 14, optimal higher-order spectra device also includes:

Optimal wavelet determines module 704, for the optimal wavelet when described evaluation module is judged as YES, described in output Sequence；

Sub-wave phase more new module 705, is used for when described evaluation module is judged as NO, the phase of the wavelet described in renewal Position.The coded data acquisition module described in performing is returned according to the sub-wave phase after updating.

Genetic manipulation can be used to update sub-wave phase, and genetic manipulation generally comprises selection, intersects and make a variation, its main purpose It is that individuality is transformed.In the present invention, their effect is the phase place updating wavelet.

Wherein " select " to be to select excellent wavelet, under making them have an opportunity to be to concentrate from current wavelet as parent Generation breeding descendants.Ideal adaptation degree is the highest, and its selected chance is the most.Selection is the key of genetic algorithm, and it embodies The thought of survival of the fittest in nature.

" intersect " is the part-structure that two fathers are individual to be replaced restructuring and generates new individual operation.This operation Purpose is able to produce new individuality the next generation, just as the marriage process of human society, is operated by intersection, loses Propagation algorithm search equal length, same-amplitude spectrum, the ability of out of phase music score ripple will be greatly improved." intersecting " is heredity Algorithm obtains optimal wavelet and the most important means of reflection coefficient.

" make a variation " is that the gene position to individual carries out inversion operation.By mutation operation, genetic algorithm can be made to keep planting The multiformity of group, to prevent premature convergence.It is said that in general, mutation probability all obtains less, the most even take 0.

Select, intersect and make a variation and be respectively arranged with the most different methods.For the particular problem of the present invention, in specific embodiment party In formula, roulette wheel selection, binary system uniform crossover method and binary system alternative method can be chosen respectively.Namely utilize binary system herein Genetic algorithm updates the phase place of wavelet, and by making synthetic seismogram reach most preferably to mate with seismic trace near well to obtain Good higher-order spectra.

As shown in Figure 12, a kind of of embodiment of the present invention offer carries out higher-order spectra based on seismic data and well-log information System also includes:

Reservoir prediction device 800, for carrying out reservoir prediction according to described optimal wavelet, obtains the reservoir of current formation Predict the outcome.When carrying out reservoir prediction, be required to be applied to convolution model, and the seismic wavelet in convolution model accurately with No, the accuracy of follow-up reservoir prediction will be directly affected.In actual seismic higher-order spectra, due to the shadow of the factors such as noise Ringing, the phase spectrum of wavelet is difficult to accurately estimate.Carrying out deconvolution based on inaccurate phase place wavelet, the resolution of its result will It is somewhat limited.Amplitude is such as utilized to carry out with angle change (Amplitude versus angle, AVA) inverting Petroleum-gas prediction, such as amplitude with the change (Amplitude Versus Offset, AVO) of offset distance to stratum three parameter (compressional wave Speed, shear wave velocity and density) carry out the identification of hydrocarbon-bearing formation.Therefore, the method that the present invention proposes is directly to estimating Wavelet carries out transform rooting, then by the phase place of binary strings genetic algorithm scanning wavelet, and then evolves optimal son rapidly Ripple, improves higher-order spectra precision, and then improves the precision of reservoir prediction.

Below in conjunction with specific embodiment, present invention offer is discussed in detail carries out son based on seismic data and well-log information The method of ripple estimation and the feasibility of system.

Fig. 1 (a) is for without wavelet composite traces figure initial in the case of making an uproar, and Fig. 1 (b) is for without protocol figure, Fig. 1 in the case of making an uproar C () is the wavelet composite traces figure estimated without the present invention in the case of making an uproar.Wherein, the phase of initial wavelet composite traces and protocol Closing coefficient is 46%, and the present invention estimates that the correlation coefficient of wavelet composite traces and protocol is 100%.Can be sent out by above-mentioned accompanying drawing Existing, in the case of without making an uproar, the well shake matching degree of initial wavelet is the lowest, after the sub-wave phase Power estimation of the inventive method, The well shake coupling of whole mixed phase wavelet has reached 100%.

Fig. 2 is the comparison diagram estimating wavelet without wavelet initial in the case of making an uproar, true wavelet and the present invention.As shown in Figure 2, exist In the case of geological data not Noise, the inventive method can accurately ask for mixed-phase seismic wavelet.

Fig. 3 (a) be noisy in the case of initial wavelet composite traces figure, Fig. 3 (b) be noisy in the case of protocol figure, Fig. 3 (c) be noisy in the case of the present invention estimate wavelet composite traces figure；Institute's plus noise be 50% random noise (meaning of 50% is to make an uproar The ratio of acoustic energy and signal energy is 1:2).Wherein, initial wavelet composite traces is 36% with the correlation coefficient of protocol, The present invention estimates that the correlation coefficient of wavelet composite traces and protocol is 82%.From above-mentioned accompanying drawing, disturb feelings at Noise Under condition, the well shake matching degree of the mixed phase wavelet of estimation has obtained bigger raising.

Fig. 4 be noisy in the case of initial wavelet, true wavelet and the present invention estimate the comparison diagram of wavelet, as shown in Figure 4, deposit When noise, estimate that wavelet is coincide preferably with true wavelet.

Fig. 5 (a) is Chang Xiangwei wavelet composite traces schematic diagram；Fig. 5 (b) is well lie real data schematic diagram；Fig. 5 (c) The wavelet composite traces schematic diagram estimated for the present invention.Wherein, the phase of Chang Xiangwei wavelet composite traces and well lie real data Closing coefficient is 61.8%, and the present invention estimates that the correlation coefficient of wavelet composite traces and well lie real data is 67.1%.For reality Border data, the well shake matching degree of the inventive method improves 5.3% than traditional method.

Fig. 6 (a) is the wavelet schematic diagram that the present invention estimates；Fig. 6 (b) is that Chang Xiangwei rotates the wavelet schematic diagram estimated.Fig. 7 A () is original seismic data schematic diagram；Fig. 7 (b) for estimating the deterministic deconvolution result schematic diagram of wavelet based on the present invention；Figure 7 (c) is Chang Xiangwei wavelet deconvolution result schematic diagram.From above-mentioned accompanying drawing, relatively conventional method, the anti-pleat of the inventive method The resolution of long-pending result is effectively improved, and has preferable dependency with well.

In sum, owing to mostly actual seismic wavelet is mixed-phase, for Chang Xiangwei rotates wavelet estimation, The present invention is the most advanced, and precision is higher.Additionally, phase estimation is forbidden also to be that in actual production, well shakes that stated accuracy is low Individual key reason.Rotate the well shake matching result obtained compared to Chang Xiangwei, this invention promotes precision and reaches more than 5%.

When existing converter technique asks for mixed phase wavelet, sub-wave phase is had a certain degree of dependence.Such as korean pine Relaxing and estimate the condition of sub-wave phase, its amount of calculation is huge, is difficult at present put in actual production.And the present invention is by wavelet Root conversion process carries out specific coding and decoding, and introduces genetic algorithm, uses genetic manipulation and combines well control as evaluating standard Then process sub-wave phase problem, efficiently reduce amount of calculation, the Fast Convergent that mixed phase wavelet is estimated can be realized.

The useful achievement of the present invention is:

1) present invention utilizes root in the feature of the inside and outside movement of unit circle, carries out specific coding to moving root mode, and then uses The selection of genetic algorithm, intersect and make a variation and heuristically update sub-wave phase, after repeatedly evolving, scan optimal wavelet.

2) present invention is by carrying out specific coding and decoding to wavelet root conversion process, and introduces genetic algorithm, uses and loses Pass operation and combine well control as interpretational criteria to process sub-wave phase problem, efficiently reducing amount of calculation, accelerating root and become Change the Fast Convergent of sub-wave phase Power estimation, improve higher-order spectra precision.

One of ordinary skill in the art will appreciate that all or part of flow process realizing in above-described embodiment method, Ke Yitong Crossing computer program and complete to instruct relevant hardware, described program can be stored in general computer read/write memory medium In, this program is upon execution, it may include such as the flow process of the embodiment of above-mentioned each method.Wherein, described storage medium can be magnetic Dish, CD, read-only store-memory body (Read-Only Memory, ROM) or random store-memory body (Random Access Memory, RAM) etc..

Those skilled in the art are it will also be appreciated that various functions that the embodiment of the present invention is listed are by hardware or soft Part realizes depending on specifically applying the design requirement with whole system.Those skilled in the art can be specific for every kind Application, it is possible to use the function described in the realization of various methods, but this realization is understood not to protect beyond the embodiment of the present invention The scope protected.

The present invention applies specific embodiment principle and the embodiment of the present invention are set forth, above example Explanation be only intended to help to understand method and the core concept thereof of the present invention；Simultaneously for one of ordinary skill in the art, According to the thought of the present invention, the most all will change, in sum, in this specification Hold and should not be construed as limitation of the present invention.

Claims (10)

1. the method carrying out higher-order spectra based on seismic data and well-log information, is characterized in that, described method includes:

Gather the seismic data of current formation, well-log information；

Minimum phase wavelet is estimated according to described seismic data；

Described minimum phase wavelet is carried out transform, obtains the multinomial about variable z；

Solve the root that described multinomial is corresponding to determine described multinomial；

Move according to described minimum phase wavelet corresponding multinomial location of root in transform territory and carry out 0-1 coding, obtain many Planting coded data, wherein, root is not moved in 0 expression, and 1 represents mobile root；

Well lie data are extracted according to described seismic data；

Reflection coefficient is determined according to described well-log information；

Well control evaluation model is determined according to described well lie data, reflection coefficient and convolution model；

Optimal wavelet is estimated according to described well control evaluation model and described multiple coded data；

Carry out reservoir prediction according to described optimal wavelet, obtain the reservoir prediction result of current formation.

Method the most according to claim 1, is characterized in that, according to described well lie data, reflection coefficient and described The well control evaluation model that determines of convolution model be:

O=max [xcorr (w (t) * r (t), s (t))]

Wherein, max [] expression takes maximum, and s (t) is well lie data, and w (t) is the seismic wavelet that decoding obtains, and r (t) is root The reflection coefficient determined according to described well-log information, xcorr is between synthetic seismogram w (t) * r (t) and well lie record s (t) Correlation coefficient, O is object function.

Method the most according to claim 2, is characterized in that, according to described well control evaluation model and described multiple coding The optimal wavelet of data estimation specifically includes:

The coded data of predetermined number is chosen from described multiple coded data；

The coded data of described predetermined number is decoded, obtains the wavelet sequence of predetermined number；

Whether the wavelet sequence according to the predetermined number described in described well control evaluation model judgement exists optimal wavelet sequence；

Optimal wavelet sequence when being judged as YES, described in output；

When being judged as NO, the phase place of the wavelet described in renewal；

The coded data performing to choose predetermined number from described multiple coded data is returned according to the sub-wave phase after updating.

Method the most according to claim 3, is characterized in that, is decoded the coded data of described predetermined number, obtains The wavelet sequence of predetermined number specifically includes:

Described coded data is carried out 0-1 decoding, obtains the root of correspondence；

The multinomial about variable z is asked for according to described root；

The multinomial of described variable z is carried out Z anti transform, obtains the wavelet sequence that described coded data is corresponding.

Method the most according to claim 3, is characterized in that, according to the predetermined number described in described well control evaluation model judgement Whether purpose wavelet sequence exists and meets the wavelet sequence of condition and specifically include:

Obtain well control threshold coefficient set in advance；

Described wavelet sequence is carried out with corresponding reflection coefficient convolution, forms synthetic seismogram；

The correlation coefficient of described synthetic seismogram and well lie record is determined according to described well control evaluation model；

Whether the correlation coefficient described in judgement is more than or equal to described well control threshold coefficient.

6. carrying out a system for higher-order spectra based on seismic data and well-log information, it is characterized in that, described system includes:

Harvester, for gathering the seismic data of current formation, well-log information；

Minimum phase wavelet determines device, for estimating minimum phase wavelet according to described seismic data；

Code device, specifically includes:

Transform unit, for described minimum phase wavelet is carried out transform, obtains the multinomial about variable z；

Polynomial solving unit, for solving the root that described multinomial is corresponding to determine described multinomial；

Coding unit, is carried out for moving according to described minimum phase wavelet corresponding multinomial location of root in transform territory 0-1 encodes, and obtains multiple coded data, and wherein, root is not moved in 0 expression, and 1 represents mobile root；

Well lie data extraction device, for extracting well lie data according to described seismic data；

Reflection coefficient determines device, for determining reflection coefficient according to described well-log information；

Well control evaluation model determines device, for determining well according to described well lie data, reflection coefficient and convolution model Control evaluation model；

Optimal higher-order spectra device, for estimating optimal son according to described well control evaluation model and described multiple coded data Ripple；

Reservoir prediction device, for carrying out reservoir prediction according to described optimal wavelet, obtains the reservoir prediction knot of current formation Really.

System the most according to claim 6, is characterized in that, described well control evaluation model determines that the well control that device determines is commented Valency model is:

O=max [xcorr (w (t) * r (t), s (t))]

Wherein, max [] expression takes maximum, and s (t) is well lie data, and w (t) is the seismic wavelet that decoding obtains, and r (t) is root The reflection coefficient determined according to well-log information, xcorr is the phase relation between composite traces w (t) * r (t) and well lie record s (t) Number, O is object function.

System the most according to claim 7, is characterized in that, described optimal higher-order spectra device specifically includes:

Coded data acquisition module, for choosing the coded data of predetermined number from described multiple coded data；

Decoder module, for being decoded the coded data of described predetermined number, obtains the wavelet sequence of predetermined number；

Whether evaluation module, exist in the wavelet sequence of the predetermined number described in judging according to described well control evaluation model Optimal wavelet sequence；

Optimal wavelet output module, for the optimal wavelet sequence when described evaluation module is judged as YES, described in output；

Sub-wave phase more new module, is used for when described evaluation module is judged as NO, the phase place of the wavelet described in renewal, according to Sub-wave phase after renewal returns the coded data acquisition module described in performing.

System the most according to claim 8, is characterized in that, described decoder module specifically includes:

Decoding unit, for described coded data is carried out 0-1 decoding, obtains the root of correspondence；

Multinomial asks for unit, for asking for the multinomial about variable z according to described root；

Z anti transform unit, for the multinomial of described variable z is carried out Z anti transform, obtains the wavelet that described coded data is corresponding Sequence.

System the most according to claim 8, is characterized in that, described evaluation module specifically includes:

Threshold value acquiring unit, is used for obtaining well control threshold coefficient set in advance；

Convolution unit, for described wavelet sequence carries out convolution with corresponding reflection coefficient, forms synthetic seismogram；

Correlation coefficient determines unit, for determining described synthetic seismogram and well lie according to described well control evaluation model The correlation coefficient of record；

Judging unit, for judging that whether described correlation coefficient is more than or equal to described well control threshold coefficient.