CN104316966B - A kind of Fluid Identification Method and system - Google Patents

Abstract

This application provides a kind of Fluid Identification Method, the method includes: will be converted into the data of prestack CRP road collection through protecting the data of the prestack CMP road collection that dampening information processes；The superposition of data of at least 3 angles is extracted according to described prestack CRP road collection；The viscoelastic fluid factor, modulus of shearing and density is calculated according to described superposition of data；Utilize the described viscoelastic fluid factor, described modulus of shearing and described density that reservoir fluid is identified.Present invention also provides a kind of fluid identification system.Can effectively improve the precision of fluid identification.

Description

A kind of Fluid Identification Method and system

Technical field

The application relates to field of geophysical exploration, particularly to a kind of Fluid Identification Method and system.

Background technology

Along with the intensification of geophysical exploration degree, underground condition becomes day by day complicated, carries out with seismic data for main body Fluid identification of reservoir can increase exploration success ratio to a certain extent, improves production efficiency and reduces cost of exploration & development, It it is one of the study hotspot of present stage petroleum exploration domain.

Existing Fluid Identification Method is that the dielectric resilient parameter by subsurface formations carries out reservoir prediction and fluid is known Other, but the ball medium of reality is viscoelastic media, existing Fluid Identification Method carries out the precision of fluid identification to be also had Wait to improve.

Therefore, the most effectively carrying out high-precision fluid identification of reservoir is that those skilled in the art are presently required solution Technical problem.

Summary of the invention

Technical problems to be solved in this application are to provide a kind of Fluid Identification Method and system, solve in prior art The problem that the precision of fluid identification has much room for improvement.

Its concrete scheme is as follows:

A kind of Fluid Identification Method, the method includes:

The data of prestack CRP road collection will be converted into through protecting the data of the prestack CMP road collection that dampening information processes；

The superposition of data of at least 3 angles is extracted according to described prestack CRP road collection；

The viscoelastic fluid factor, modulus of shearing and density is calculated according to described superposition of data；

Utilize the described viscoelastic fluid factor, described modulus of shearing and described density that reservoir fluid is identified.

Above-mentioned method, it is preferred that the described calculating viscoelastic fluid factor, modulus of shearing and density include:

Set up viscoplasticity impedance model；

Carry out, according to described viscoplasticity impedance inversion approach, the viscoelastic fluid factor inverting that decays, obtain described visco elastic fluids The body factor, described modulus of shearing and described density.

Above-mentioned method, it is preferred that described viscoplasticity impedance model of setting up includes:

Set up the longitudinal wave reflection coefficient equation of viscoelastic medium；

Utilize weak viscoplasticity and similar medium it is assumed that obtain described longitudinal wave reflection coefficient approximate expression；

Described longitudinal wave reflection coefficient approximate expression is utilized to set up the calculating viscoelastic fluid factor, modulus of shearing and volume close The equation of degree.

Above-mentioned method, it is preferred that the described longitudinal wave reflection coefficient equation setting up viscoelastic medium specifically includes:

Determine viscoelastic fluid factor expression；

Determine viscoelastic media longitudinal wave reflection coefficient expressions；

Described viscoelastic fluid factor expression is brought in described viscoelastic media longitudinal wave reflection coefficient expressions, Arrive with the longitudinal wave reflection coefficient equation of described viscoelastic fluid factor representation.

Above-mentioned method, it is preferred that described according to viscoplasticity impedance inversion approach carry out decay the viscoelastic fluid factor anti- Drill, obtain the described viscoelastic fluid factor, described modulus of shearing and described density and include:

Described viscoplasticity impedance model is carried out linear transformation；

The viscoplasticity impedance data that at least 3 angles described in acquisition are corresponding；

According to described viscoplasticity impedance data and log, it is calculated the described viscoelastic fluid factor, described shearing Modulus and described bulk density.

Above-mentioned method, it is preferred that viscoplasticity impedance data corresponding at least 3 angles of described acquisition includes:

Extract angle wavelet；

Carry out viscoplasticity Impedance Inversion, obtain at least 3 viscoplasticity impedance datas corresponding to angle.

Above-mentioned method, it is preferred that the extraction conditions of described extraction angle wavelet is:

E=∑ (s_i-d_i)²

Wherein, E is error energy, and s=r*w is synthetic seismogram, and r is reflectance factor, and w is the wavelet extracted, and d is ground Shake record, selects optimal wavelet, makes error energy minimum simultaneously.

Above-mentioned method, it is preferred that described in carry out viscoplasticity Impedance Inversion, obtain at least 3 viscoplasticity corresponding to angle Impedance data particularly as follows:

With the viscoplasticity impedance data body of each angle of Sparse Pulse Inversion inversion method, its objective optimization function is:

F=L_p(r(θ))+λL_q(s-d)+α^-1(ΔZ_trend)

Wherein, r (θ) is angle reflection coefficient sequence, and Δ Z is and the difference sequence of viscoplasticity impedance trend, and d is seismic channel sequence Row, s is synthetic seismogram sequence, and λ is residual error weight factor, and α is trend weight factor, p and q is the L mould factor.

Above-mentioned method, it is preferred that the described data through protecting the prestack CMP road collection that dampening information processes specially comprise The data of the prestack CMP road collection affected by viscoplasticity.

Above-mentioned method, it is preferred that the described superposition of data according to prestack CRP road collection at least 3 angles of extraction is specially Low-angle, middle angle and the superposition of data of wide-angle is extracted according to described prestack CRP road collection.

A kind of fluid identification system, including:

Conversion unit: for the data of the prestack CMP road collection through protecting dampening information process are converted into prestack CRP road collection Data；

First extraction unit: for extracting the superposition of data of at least 3 angles according to described prestack CRP road collection；

First computing unit: for calculating the viscoelastic fluid factor, modulus of shearing and density according to described superposition of data；

Recognition unit: be used for utilizing the described viscoelastic fluid factor, described modulus of shearing and described density to reservoir fluid It is identified.

Above-mentioned system, it is preferred that described first computing unit includes:

First sets up unit: be used for setting up viscoplasticity impedance model；

First inverting unit: anti-for carrying out the decay viscoelastic fluid factor according to described viscoplasticity impedance inversion approach Drill, obtain the described viscoelastic fluid factor, described modulus of shearing and described density.

Above-mentioned system, it is preferred that described first sets up unit includes:

Second sets up unit, for setting up the longitudinal wave reflection coefficient equation of viscoelastic medium；

Approximating unit: be used for utilizing weak viscoplasticity and similar medium it is assumed that obtain described longitudinal wave reflection coefficient approximate expression Formula；

3rd sets up unit: be used for utilizing described longitudinal wave reflection coefficient approximate expression set up calculate viscoelastic fluid because of Son, modulus of shearing and the equation of density.

Above-mentioned system, it is preferred that described second sets up unit includes:

First determines unit: be used for determining viscoelastic fluid factor expression；

Second determines unit: be used for determining viscoelastic media longitudinal wave reflection coefficient expressions；

Bring unit into: for described viscoelastic fluid factor expression being brought into described viscoelastic media longitudinal wave reflection system In number expression formula, obtain the longitudinal wave reflection coefficient equation with described viscoelastic fluid factor representation.

Above-mentioned system, it is preferred that described first inverting unit includes:

Converter unit: for described viscoplasticity impedance model is carried out linear transformation；

Acquiring unit: be used for the viscoplasticity impedance data that at least 3 angles described in obtaining are corresponding；

Second computing unit: for according to described viscoplasticity impedance data and log, be calculated described viscoplasticity Fluid factor, described modulus of shearing and described bulk density.

Above-mentioned system, it is preferred that described acquiring unit includes:

Second extraction unit: be used for extracting angle wavelet；

Second inverting unit: be used for carrying out viscoplasticity Impedance Inversion, obtains at least 3 viscoplasticity impedance numbers corresponding to angle According to

In a kind of Fluid Identification Method that the application provides, by the number through protecting the prestack CMP road collection that dampening information processes According to the data being converted into prestack CRP road collection；The superposition of data of at least 3 angles is extracted according to described prestack CRP road collection；Calculate viscous The elastic fluid factor, modulus of shearing and density；Utilize the described viscoelastic fluid factor, described modulus of shearing and described density to storage Layer fluid is identified, and can effectively improve the precision of fluid identification.

Accompanying drawing explanation

For the technical scheme being illustrated more clearly that in the embodiment of the present application, in embodiment being described below required for make Accompanying drawing be briefly described, it should be apparent that, the accompanying drawing in describing below is only some embodiments of the application, for From the point of view of those of ordinary skill in the art, on the premise of not paying creative work, it is also possible to obtain it according to these accompanying drawings His accompanying drawing.

Fig. 1 is the flow chart of a kind of Fluid Identification Method embodiment 1 of the application；

Fig. 2 is the flow chart of a kind of Fluid Identification Method embodiment 2 of the application；

Fig. 3 is the flow chart of a kind of Fluid Identification Method embodiment 3 of the application；

Fig. 4 is the flow chart of a kind of Fluid Identification Method embodiment 4 of the application；

Fig. 5 is the flow chart of a kind of Fluid Identification Method embodiment 5 of the application；

Fig. 6 is the flow chart of a kind of Fluid Identification Method embodiment 6 of the application；

Fig. 7 is the flow chart that a kind of Fluid Identification Method of the application is specifically applied；

Fig. 8 is the flow chart of the another concrete application of a kind of Fluid Identification Method of the application；

Fig. 9 is the structural representation of a kind of fluid identification system embodiment 1 of the application；

Figure 10 is the structural representation of a kind of fluid identification system embodiment 2 of the application；

Figure 11 is the structural representation of a kind of fluid identification system embodiment 3 of the application；

Figure 12 is the structural representation of a kind of fluid identification system embodiment 4 of the application；

Figure 13 is the structural representation of a kind of fluid identification system embodiment 5 of the application；

Figure 14 is the structural representation of a kind of fluid identification system embodiment 6 of the application.

Detailed description of the invention

The core of the present invention is to provide the method and system of a kind of fluid identification, it is possible to the effective essence improving fluid identification Degree.

Below in conjunction with the accompanying drawing in the embodiment of the present application, the technical scheme in the embodiment of the present application is carried out clear, complete Describe, it is clear that described embodiment is only some embodiments of the present application rather than whole embodiments wholely.Based on Embodiment in the application, 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 the application protection.

With reference to Fig. 1, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 1, may comprise steps of:

Step S101: the number of prestack CRP road collection will be converted into through protecting the data of the prestack CMP road collection that dampening information processes According to.

In the application, the described data through protecting the prestack CMP road collection that dampening information processes specially contain by viscoelastic Property impact the data of prestack CMP road collection.

Owing to CRP road collection being overlapped or uses other monitoring means, the whether accurate of migration velocity, institute can be checked To need the data of prestack CMP road collection to be converted into the data of prestack CRP road collection.

Concrete conversion process is:

1) analyze the offset distance distribution of seismic channel in whole district CMP road collection, determine maximum, smallest offset away from and neighboring track skew Away from increment.

2) utilize final stacking velocity field, the road in CMP road collection is carried out the dynamic correction of difference, all roads is moved school to above The offset distance position determined.

3) in CMP road collection all have same offset away from road be overlapped, and all scarce offset distances are supplemented null value Road, is that the road in the most individual CMP has the offset distance homogenizing sequence.

4) using through interpolation, smoothing processing stacking velocity field as initial offset velocity field, to input CMP road collection altogether Offset distance territory bias internal, obtains initial CRP road collection.

5) utilize identical velocity field, CRP road collection is carried out reaction correction, reaction correction Hou CRP road is re-started Velocity analysis, obtains new migration velocity field.

6) utilize new migration velocity field, input CMP road collection is re-started pre-stack time migration, obtains Xin CRP road Collection.

Step S102: extract the superposition of data of at least 3 angles according to described prestack CRP road collection.

Before fetching portion angular stack data, needing to divide angular-trace gather, general division methods is:

Firstly, it is necessary to determine maximum incident angle degree that main purpose interval can completely receive and to be disturbed wave action little Minimum incident angle.The maximum incident angle degree that main purpose interval to be determined can completely receive, we should adopt with reference to earthquake The maximum offset of collection and the degree of depth of target zone, the seismic interval velocity simultaneously process obtained, CheckShot or compressional wave time difference are bent Line is demarcated, and carries out offset distance angular transition, determines the maximum angle that interval of interest can completely receive, and exceedes this angle When spending, the degree of covering of target zone can be not.The determination of minimum angles is mainly come by the shortcut Analysis signal-to-noise ratio (SNR) of road collection data Determine.As in the processing of ocean data owing to repeatedly POP is all over growing, have, near migration range, many subwaves residual that comparison is many, fit When by the impact adjusting the many subwaves of reduction of minimum angles, bigger minimum angles can reduce many subwaves or the shadow of interference ripple Ring.

Between maximum angle and minimum angles, how many angular stack data volumes it are divided into it addition, be also predefined.From reason Say in opinion that superposition of data body angular divisions fine degree is directly connected to the holding of seismic AVO feature, but in view of actual seismic The factors such as workload such as the signal to noise ratio of prestack road collection, the internal memory of computer, disk space and man-machine interaction wavelet extraction demarcation, General recommendations is superimposed as 3-5 angular stack data volume.Such as, the maximum angle general stackable one-tenth of observation system less than 40 degree Three angular stack data volumes, maximum angle four or more angular stack data volume of stackable one-tenth more than 45 degree.

In this application, the described superposition of data according to prestack CRP road collection at least 3 angles of extraction is specially according to described Prestack CRP road collection extracts low-angle, middle angle is folded and the superposition of data of wide-angle.

Step S103: calculate the viscoelastic fluid factor, modulus of shearing and density according to described superposition of data.

Step S104: utilize the described viscoelastic fluid factor, described modulus of shearing and described density that reservoir fluid is carried out Identify.

Think that attenuation when seismic wave is propagated in real medium is owing to viscosity pore fluid flow causes, viscous Elasticity can characterize this attenuation apparently, and the viscoelastic fluid factor higher can carry out fluid identification, improves stream The precision of body identification.

The actual seismic wave received not only contains the impact of elastic property, also contains the impact of inelastic nature, This stiff impact is attributed to pore-fluid and there is the Viscoelastic effect caused, according to more stiff reservoir predictions Carry out fluid identification with fluid identification technology, improve the precision of identification.

With reference to Fig. 2, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 2, described calculating visco elastic fluids The body factor, modulus of shearing and density include:

Step S201: set up viscoplasticity impedance model.

Step S202: carry out, according to described viscoplasticity impedance inversion approach, the viscoelastic fluid factor inverting that decays, obtain institute State the viscoelastic fluid factor, described modulus of shearing and described density.

With reference to Fig. 3, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 3, described viscoplasticity of setting up hinders Anti-model includes:

Step S301: set up the longitudinal wave reflection coefficient equation of viscoelastic medium.

Step S302: utilize weak viscoplasticity and similar medium it is assumed that obtain described longitudinal wave reflection coefficient approximate expression.

Step S303: utilize described longitudinal wave reflection coefficient approximate expression to set up and calculate the viscoelastic fluid factor, shearing mould Amount and the equation of density.

With reference to Fig. 4, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 4, described viscoplasticity of setting up is situated between The longitudinal wave reflection coefficient equation of matter specifically includes:

Step S401: determine viscoelastic fluid factor expression.

Step S402: determine viscoelastic media longitudinal wave reflection coefficient expressions.

Step S403: described viscoelastic fluid factor expression is brought into described viscoelastic media longitudinal wave reflection coefficient table Reach in formula, obtain the longitudinal wave reflection coefficient equation with described viscoelastic fluid factor representation.

With reference to Fig. 5, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 5, described according to viscoplasticity resistance Anti-inversion method carries out the viscoelastic fluid factor inverting that decays, and obtains the described viscoelastic fluid factor, described modulus of shearing and institute State density to include:

Step S501: described viscoplasticity impedance model is carried out linear transformation.

Step S502: the viscoplasticity impedance data that at least 3 angles described in acquisition are corresponding.

Step S503: according to described viscoplasticity impedance data and log, be calculated the described viscoelastic fluid factor, Described modulus of shearing and described density.

With reference to Fig. 6, it is shown that the flow chart of the application a kind of Fluid Identification Method embodiment 6, at least 3 angles of described acquisition The viscoplasticity impedance data that degree is corresponding includes:

Step S601: extract angle wavelet.

The extraction conditions of described extraction angle wavelet is:

E=∑ (s_i-d_i)²

Wherein, E is error energy, and s=r*w is synthetic seismogram, and r is reflectance factor, and w is the wavelet extracted, and d is ground Shake record, selects optimal wavelet, makes error energy minimum simultaneously.

Step S602: carry out viscoplasticity Impedance Inversion, obtains at least 3 viscoplasticity impedance datas corresponding to angle.

Described carry out viscoplasticity Impedance Inversion, obtain at least 3 viscoplasticity impedance datas corresponding to angle particularly as follows:

With the viscoplasticity impedance data body of each angle of Sparse Pulse Inversion inversion method, its objective optimization function is:

F=L_p(r(θ))+λL_q(s-d)+α^-1(ΔZ_trend)

Wherein, r (θ) is angle reflection coefficient sequence, and Δ Z is and the difference sequence of viscoplasticity impedance trend, and d is seismic channel sequence Row, s is synthetic seismogram sequence, and λ is residual error weight factor, and α is trend weight factor, p and q is the L mould factor.

With reference to Fig. 7, it is shown that the flow chart that a kind of Fluid Identification Method of the application is specifically applied, describe concrete calculating and glue The process of the elastic fluid factor, modulus of shearing and density.

First the method carries out subangle overlap-add procedure, rule of thumb formula to seismic data, utilizes log to set up The dampening fluid factor and modulus of shearing model, utilize angle part superposition road collection to extract angle wavelet, finally retrain at model Under the conditions of utilize non-resilient impedance inversion approach to carry out dampening fluid factor inverting.Concrete process is divided into following main step Rapid:

1) carry out earthquake data before superposition, on the basis of corresponding relative amplitude preserved processing and denoising, carrying out subangle superposition Obtain the subangle superposition seismic data of at least three angle.

2) calculate quality factor, viscoelastic fluid factor parameter f and modulus of shearing μ from existing well-log information, and utilize The pseudo-well curve calculated sets up corresponding model with layer position explanation results.

3) the angle wavelet of response is extracted according to different partial angle stack data volumes and well-log information.

4) under Sparse Pulse Inversion framework, add model constraint and carry out the viscoelastic fluid factor inverting that decays.

On the basis of each embodiment above-mentioned, the concrete computational algorithm of the viscoelastic fluid factor, modulus of shearing and density As follows:

1) according to earthquake Rock physical analysis and calculate each medium parameter (p-wave impedance I_p, velocity of longitudinal wave V_p, Poisson ratioσ, Bulk modulus K, Lame parameter λ, elastic pore-fluid parameter f, viscoelastic fluid factor f_ane), the fluid comparing this group parameter is quick Perception, it is known that, in numerous attributes, the sensitiveness of viscoelastic fluid factor pair pore-fluid is the strongest, can be used for reservoir fluid Identify, improve accuracy of identification.

2) elastic fluid fluid factor is expressed as:OrSimilarly, by the stream of viscoelastic media The body factor is referred to as the viscoelastic fluid factor, and its form is:

$\begin{array}{c}{f}_{\mathrm{anela}}=\mathrm{\ρ}{\mathrm{\α}}^2-\mathrm{c\ρ}{\mathrm{\β}}^2\\ =\mathrm{\ρ}{V}_P^2(1+\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P})-\mathrm{c\ρ}{V}_S^2(1+\frac{2}{\mathrm{\π}{Q}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S})\\ =\mathrm{\ρ}{V}_P^2-\mathrm{c\ρ}{V}_S^2+\mathrm{\ρ}{V}_P^2(\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P})-\mathrm{c\ρ}{V}_S^2(\frac{2}{\mathrm{\π}{Q}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S})\\ =f_{\mathrm{ela}}+\mathrm{\Δ}{f}_Q\end{array}---\left(1\right)$

Wherein, corresponding with elastic part fluid factor isThe non-resilient impact on fluid factor Be equivalent to the small sample perturbations applied in elastic background, therefore use Δ f_QRepresent the non-resilient disturbance to fluid factor, such as following formula:

$\begin{array}{c}\mathrm{\Δ}{f}_Q=\mathrm{\ρ}{V}_P^2(\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P})-\mathrm{c\ρ}{V}_S^2(\frac{2}{\mathrm{\π}{Q}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S})\\ =\mathrm{\ρ}{V}_P^2\left(\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)\right)-\mathrm{c\ρ}{V}_S^2\left(\frac{2}{\mathrm{\π}{Q}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)\right)-(\mathrm{\ρ}{V}_P^2\frac{i}{Q_P}-\mathrm{c\ρ}{V}_S^2\frac{i}{Q_S})\\ =\mathrm{\ρ}\frac{2}{\mathrm{\π}}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)[\frac{V_P^2}{Q_P}-c\frac{V_S^2}{Q_S}]-\mathrm{i\ρ}[\frac{V_P^2}{Q_P}-c\frac{V_S^2}{Q_S}]\end{array}---\left(2\right)$

Therefore, the viscoelastic fluid factor is ultimately expressed as:

f_anela=f_ela+Δf_Q (3)

Wherein,

$f_{\mathrm{ela}}=\mathrm{\ρ}{V}_P^2-\mathrm{c\ρ}{V}_S^2$

$\begin{array}{c}\mathrm{\Δ}{f}_Q=\mathrm{\ρ}{V}_P^2(\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P})-\mathrm{c\ρ}{V}_S^2(\frac{2}{\mathrm{\π}{Q}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S})\end{array}$

Viscoelastic media elastic impedance and the viscoelastic fluid factor provide for extracting viscoelastic media parameter from geological data New approach.

3) the reflectance factor equation of viscoelastic media is:

$R_{\mathrm{PP}}^{\mathrm{ane}}\left(\mathrm{\θ}\right)=a\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{V}_P}{V_P}+b\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{V}_S}{V_S}+c\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}+d\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{Q}_P}{Q_P}+e\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{Q}_S}{Q_S}---\left(4\right)$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{2}(1+{\tan }^2\mathrm{\θ})$

$b\left(\mathrm{\θ}\right)=-4{\sin }^2\mathrm{\θ}\frac{V_S^2}{V_P^2}$

$c\left(\mathrm{\θ}\right)=\frac{1}{2}(1-4{\sin }^2\mathrm{\θ}\frac{V_S^2}{V_P^2})$

$d\left(\mathrm{\θ}\right)=-\frac{1}{2}(1+{\tan }^2\mathrm{\θ})\frac{\mathrm{\Γ}}{Q_P+\mathrm{\Γ}}$

$e\left(\mathrm{\θ}\right)=4{\sin }^2\mathrm{\θ}\frac{V_S^2}{V_P^2}\frac{\mathrm{\Γ}}{Q_S+\mathrm{\Γ}}$

$\mathrm{\Γ}=\frac{1}{\mathrm{\π}}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{2}$

In above formula, first three items is in reference frequency ω about P-wave And S_rThe phase velocity at place and the reflection system of density composition Number, latter two is the reflectance factor change caused due to inelastic nature.When only considering the decay that pore-fluid causes, due to Shear wave is the most fluid affected, therefore, containing Q_SItem can omit, formula (4) becomes:

$R_{\mathrm{PP}}^{\mathrm{ane}}\left(\mathrm{\θ}\right)=a\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{V}_P}{V_P}+b\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{V}_S}{V_S}+c\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}+d\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{Q}_P}{Q_P}---\left(5\right)$

4) the reflectance factor equation that with fluid factor represent corresponding with elastic fluid:

$R_{\mathrm{PP}}^{\mathrm{ela}}\left(\mathrm{\θ}\right)=a\left(\mathrm{\θ}\right)\frac{\mathrm{\Δf}}{f}+b\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}+c\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}---\left(6\right)$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{4}(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}){\sec }^2\mathrm{\θ}$

$b\left(\mathrm{\θ}\right)=\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}$

$c\left(\mathrm{\θ}\right)=\frac{1}{2}(1-4{\sin }^2\mathrm{\θ}\frac{V_S^2}{V_P^2})$

At this for building by the reflectance factor equation of viscoelastic fluid factor representation, according to a upper joint reflection coefficient formula Deriving, the form of observation type (5), containing elastic parameter i.e. reference velocity V in formula_PAnd V_S, utilize the expression of the elastic fluid factor Formula, by V_PIt is expressed as the elastic fluid factor and V_S, the function of ρ, substitute into formula (5), be use elastic parameter build elastic fluid because of Son, obtains the viscoelastic media longitudinal wave reflection coefficient expressions of following form:

$R_{\mathrm{PP}}^{\mathrm{ane}}\left(\mathrm{\θ}\right)=a\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{f}_{\mathrm{ela}}}{f_{\mathrm{ela}}}+b\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}+c\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}+d\left(\mathrm{\θ}\right)\frac{\mathrm{\Δ}{Q}_P}{Q_P}---\left(7\right)$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{4}(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}){\sec }^2\mathrm{\θ}$

$b\left(\mathrm{\θ}\right)=\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}$

$c\left(\mathrm{\θ}\right)=\frac{1}{2}(1-4{\sin }^2\mathrm{\θ}\frac{V_S^2}{V_P^2})$

$d\left(\mathrm{\θ}\right)=-\frac{1}{2}(1+{\tan }^2\mathrm{\θ})\frac{\mathrm{\Γ}}{Q_P+\mathrm{\Γ}}$

$\mathrm{\Γ}=\frac{1}{\mathrm{\π}}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{2}$

By f_ela=f_anela-Δf_Q, substitute into (7) formula, obtain by viscoelastic fluid factor f_anelaWith Δ f_QAnd other parameter lists The viscoplasticity reflectance factor equation shown.

Then further according to V_SWith the relation of μ, through certain mathematical derivation, obtain indulging with viscoelastic fluid factor representation Wave reflection coefficient equation, at this, we directly give this form:

$R_{\mathrm{pp}}\left(\mathrm{\θ}\right)\≈\left[(\frac{1}{4}-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^{\′2}}){\sec }^2\mathrm{\θ}\right]\frac{\mathrm{\Δ}{f}_{\mathrm{ane}}}{{\stackrel{\‾}{f}}_{\mathrm{ane}}}+\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^{\′2}}[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4}{\sec }^2\mathrm{\θ}-2{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ\μ}}{\stackrel{\‾}{\mathrm{\μ}}}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\stackrel{\‾}{\mathrm{\ρ}}}---\left(8\right)$

Wherein, $f_{\mathrm{ane}}=\mathrm{\ρ}{\mathrm{\α}}^2-{\mathrm{\γ}}_{\mathrm{dry}}^2\mathrm{\μ}=\mathrm{\ρ}{\mathrm{\α}}^2-{\mathrm{\γ}}_{\mathrm{dry}}^2\mathrm{\ρ}{\mathrm{\β}}^2.$

Complex velocity phase velocity is represented, substitutes into above formula,

$\frac{1}{\mathrm{\α}}=\frac{1}{V_P}(1-\frac{1}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)+\frac{i}{2Q_P})---\left(9\right)$

$\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^{\′2}}=\frac{{\mathrm{\α}}^2}{{\mathrm{\β}}^2}\≈\frac{V_P^2}{V_S^2}[1+\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P}]---\left(10\right)$

$\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^{\′2}}=\frac{{\mathrm{\β}}^2}{{\mathrm{\α}}^2}\≈\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^2}[1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)+\frac{i}{Q_P}]---\left(11\right)$

Assume that dry rock P-S wave velocity ratio is not by influence of fading.The P-S wave velocity ratio of saturated rock is substituted into reflection system Number expression formula,

$\begin{array}{c}{R}_{\mathrm{pp}}\left(\mathrm{\θ}\right)\≈\frac{1}{4}[1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}(1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)+\frac{i}{Q_P})]{\sec }^2\mathrm{\θ}\frac{\mathrm{\Δ}{f}_{\mathrm{ane}}}{f_{\mathrm{ane}}}\\ +\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^2}[1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)+\frac{i}{Q_P}](\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4}{\sec }^2\mathrm{\θ}-2{\sin }^2\mathrm{\θ})+\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}(\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4})\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}\end{array}---\left(12\right)$

It can be seen that this derivation result is equivalent to carry out elasticity of substitution fluid factor by the inelastic fluid factor, be only by In introducing decay, difference on parameter characterization.The reason obtaining this result is: we only assume that compressional wave has decay, Shear wave is not decayed, so so that modulus of shearing and density item in reflectance factor are unaffected, and only fluid factor occurs Change.

Viscoplasticity longitudinal wave reflection coefficient has two problems, adds the difficulty of actual application.First, containing with reference to frequency Rate, second, containing imaginary part item.For first problem, we choose reference frequency is 1Hz, and frequency is the dominant frequency of seismic wavelet 30Hz, the fluid parameter so characterized is corresponding with the reference frequency of 1Hz.For the imaginary part item in reflectance factor, as it is assumed that Weak inelastic nature, therefore, imaginary part item is that we are ignored the most in a small amount.

Carry out two above it is assumed that obtain reflectance factor and be:

$\begin{array}{c}{R}_{\mathrm{pp}}\left(\mathrm{\θ}\right)\≈\frac{1}{4}[1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}(1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(30\right))]{\sec }^2\mathrm{\θ}\frac{\mathrm{\Δ}{f}_{\mathrm{ane}}}{f_{\mathrm{ane}}}\\ +\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^2}[1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(30\right)](\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4}{\sec }^2\mathrm{\θ}-2{\sin }^2\mathrm{\θ})+\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}(\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4})\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}\end{array}---\left(13\right)$

Formula (13) is the reflectance factor equation comprising the viscoelastic fluid factor that the present invention derives.

For studying the inelastic nature impact on reflectance factor, we carry out analysis below, and the parameter in use table 1 compares Viscoelastic media reflectance factor and the difference of elastic fluid reflectance factor.It can be seen that inelastic nature causes reflectance factor inclined From elastic fluid reflectance factor.

Table 1 elastic fluid parameter

According to reflectance factor approximate expression, with reference to the derivation of elastic impedance, can obtain compressional wave quality factor is Q Viscoplasticity impedance corresponding to viscoelastic media, be expressed as with QEI:

$\mathrm{QEI}\left(\mathrm{\θ}\right)=f_{\mathrm{ane}}^{a\left(\mathrm{\θ}\right)}{\mathrm{\μ}}^{b\left(\mathrm{\θ}\right)}{\mathrm{\ρ}}^{c\left(\mathrm{\θ}\right)}---\left(14\right)$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{2}[1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}(1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(30\right))]{\sec }^2\mathrm{\θ}$

$b\left(\mathrm{\θ}\right)=\frac{1}{{\mathrm{\γ}}_{\mathrm{sat}}^2}[1-\frac{2}{\mathrm{\π}{Q}_P}\log \left(30\right)](\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{2}{\sec }^2\mathrm{\θ}-4{\sin }^2\mathrm{\θ})$

$c\left(\mathrm{\θ}\right)=1-\frac{{\sec }^2\mathrm{\θ}}{2}$

This expression formula will be used for calculating the viscoelastic fluid factor.

5) in order to obtain the viscoelastic fluid factor, modulus of shearing and density, it is necessary to obtain at least three from inverting different The viscoplasticity elastic impedance body of angle, it may be assumed that QEI (θ₁)、QEI(θ₂)、QEI(θ₃), utilize these three viscoplasticity elastic impedance body just The viscoelastic fluid factor, modulus of shearing and density can be extracted.Non-linear due to equation (14), if directly calculating, certainly will bring Many troubles, therefore carry out linear transformation by model.Are taken the logarithm in equation both sides, can obtain such as drag:

Ln (QEI)=a (θ) ln (f_ane)+b(θ)ln(μ)+c(θ)ln(ρ) (15)

In order to obtain ln (f_ane), ln (μ) and ln (ρ), need three different angles viscoplasticity impedance body.By three angles Angle value substitutes into theoretical model expression formula (15) respectively can obtain following equation group:

$\left[\begin{array}{ccc}a\left({\mathrm{\θ}}_1\right)& b\left({\mathrm{\θ}}_1\right)& c\left({\mathrm{\θ}}_1\right)\\ a\left({\mathrm{\θ}}_2\right)& b\left({\mathrm{\θ}}_2\right)& c\left({\mathrm{\θ}}_3\right)\\ a\left({\mathrm{\θ}}_3\right)& b\left({\mathrm{\θ}}_3\right)& c\left({\mathrm{\θ}}_3\right)\end{array}\right]\left[\begin{array}{c}\ln \left(f_{\mathrm{ane}}\right)\\ \ln \left(\mathrm{\μ}\right)\\ \ln \left(\mathrm{\ρ}\right)\end{array}\right]=\left[\begin{array}{c}\ln \left(\mathrm{QEI}\left({\mathrm{\θ}}_1\right)\right)\\ \ln \left(\mathrm{QEI}\left({\mathrm{\θ}}_2\right)\right)\\ \ln \left(\mathrm{QEI}\left({\mathrm{\θ}}_3\right)\right)\end{array}\right]---\left(16\right)$

Equation (16) solve shape as Ax=B is solved, due to angle it is known that have only to know three separate Viscoplasticity impedance data body, just can be in the hope of with conventional method for solving.

In actual applications, if equation group (16) is directly solved, obtained ln (f_ane), ln (μ) and ln (ρ) Value there will be the phenomenon of instability in some sample point, and may run counter to actual physics and geological Significance.To this end, also Above formula need to be made conversion further.

Knowable to formula (16): when angle is identical, the same parameters (ln (f of each sampled point_ane), ln (μ) and ln (ρ)) institute right The coefficient value answered is identical, and they do not change over time.Therefore, it can be write as (16) formula:

Ln (QEI (t, θ))=a (θ) ln (f_ane(t))+b(θ)ln(μ(t))+c(θ)ln(ρ(t)) (17)

For with along with different sampled points i.e. different time t, have following formula to set up:

$\left[\begin{array}{ccc}\ln \left(f_{\mathrm{ane}}\left(t_1\right)\right)& \ln \left(\mathrm{\μ}\left(t_1\right)\right)& \ln \left(\mathrm{\ρ}\left(t_1\right)\right)\\ \ln \left(f_{\mathrm{ane}}\left(t_2\right)\right)& \ln \left(\mathrm{\μ}\left(t_2\right)\right)& \ln \left(\mathrm{\ρ}\left(t_2\right)\right)\\ M& M& M\\ \ln \left(f_{\mathrm{ane}}\left(t_n\right)\right)& \ln \left(\mathrm{\μ}\left(t_n\right)\right)& \ln \left(\mathrm{\ρ}\left(t_n\right)\right)\end{array}\right]\left[\begin{array}{c}a\left(\mathrm{\θ}\right)\\ b\left(\mathrm{\θ}\right)\\ c\left(\mathrm{\θ}\right)\end{array}\right]=\left[\begin{array}{c}\ln \left(\mathrm{QEI}(t_1,\mathrm{\θ})\right)\\ \ln \left(\mathrm{QEI}(t_2,\mathrm{\θ})\right)\\ M\\ \ln \left(\mathrm{QEI}(t_n,\mathrm{\θ})\right)\end{array}\right]---\left(18\right)$

Owing to the viscoplasticity impedance body used when extracting the parameter such as the viscoelastic fluid factor is from inverting income value, because of This, optional negating drills obtained well lie viscoplasticity impedance data and log (P-and S-wave velocity curve and density song Line) set up relation of plane, the relation between this rock parameter and viscoplasticity impedance is the closest, a (θ), the b obtained (θ), c (θ) is the most representative.

If three angles that viscoplasticity impedance body is corresponding are respectively θ₁、θ₂And θ₃, the most just can get three equations:

ln(QEI(t,θ₁))=a (θ₁)ln(f_ane(t))+b(θ₁)ln(μ(t))+c(θ₁)ln(ρ(t)) (19)

ln(QEI(t,θ₂))=a (θ₂)ln(f_ane(t))+b(θ₂)ln(μ(t))+c(θ₂)ln(ρ(t)) (20)

ln(QEI(t,θ₃))=a (θ₃)ln(f_ane(t))+b(θ₃)ln(μ(t))+c(θ₃)ln(ρ(t)) (21)

These three equation is substituted into respectively formula (18) formula and can get nine constant (a (θ₁)、b(θ₁)、c(θ₁)；a(θ₂)、b (θ₂)、c(θ₂)；a(θ₃)、b(θ₃)、c(θ₃)), (17) formula that they substituted into respectively can obtain following equations group:

$\left\{\begin{array}{c}\ln \left(\mathrm{QEI}(t,{\mathrm{\θ}}_1)\right)=a\left({\mathrm{\θ}}_1\right)\ln \left(f_{\mathrm{ane}}\left(t\right)\right)+b\left({\mathrm{\θ}}_1\right)\ln \left(\mathrm{\μ}\left(t\right)\right)+c\left({\mathrm{\θ}}_1\right)\ln \left(\mathrm{\ρ}\left(t\right)\right)\\ \ln \left(\mathrm{QEI}(t,{\mathrm{\θ}}_2)\right)=a\left({\mathrm{\θ}}_2\right)\ln \left(f_{\mathrm{ane}}\left(t\right)\right)+b\left({\mathrm{\θ}}_2\right)\ln \left(\mathrm{\μ}\left(t\right)\right)+c\left({\mathrm{\θ}}_2\right)\ln \left(\mathrm{\ρ}\left(t\right)\right)\\ \ln \left(\mathrm{QEI}(t,{\mathrm{\θ}}_3)\right)=a\left({\mathrm{\θ}}_3\right)\ln \left(f_{\mathrm{ane}}\left(t\right)\right)+b\left({\mathrm{\θ}}_3\right)\ln \left(\mathrm{\μ}\left(t\right)\right)+c\left({\mathrm{\θ}}_3\right)\ln \left(\mathrm{\ρ}\left(t\right)\right)\end{array}\right.---\left(22\right)$

Each angle viscoplasticity impedance body and each constant calculated by log of inverting gained are substituted in the lump Equation group (22), thus obtain the viscoelastic fluid factor, modulus of shearing and the density at each road arbitrary sampled point t.

With reference to Fig. 8, it is shown that the flow chart of the another concrete application of a kind of Fluid Identification Method of the application, by above-mentioned viscoplasticity The computational methods of fluid factor are applied to fluid identification of reservoir.

The method that Fluid Identification Method embodiment 1 a kind of with above-mentioned the application is provided is corresponding, sees Fig. 9, the application Additionally providing a kind of fluid identification system embodiment 1, in the present embodiment, this system includes:

Conversion unit 701: for the data of the prestack CMP road collection through protecting dampening information process are converted into prestack CRP The data of road collection.

First extraction unit 702: for extracting the superposition of data of at least 3 angles according to described prestack CRP road collection.

First computing unit 703: be used for calculating the viscoelastic fluid factor, modulus of shearing and density.

Recognition unit 704: be used for utilizing the described viscoelastic fluid factor, described modulus of shearing and described density to reservoir stream Body is identified.

On the basis of Fig. 7, the application provide a kind of fluid identification system embodiment 2 as shown in Figure 10, described first Computing unit 703 includes:

First sets up unit 801: be used for setting up viscoplasticity impedance model.

First inverting unit 802: for carrying out, according to described viscoplasticity impedance inversion approach, the viscoelastic fluid factor that decays Inverting, obtains the described viscoelastic fluid factor, described modulus of shearing and described density.

As shown in figure 11, described first sets up unit 801 wraps a kind of fluid identification system embodiment 3 that the application provides Include:

Second sets up unit 901, for setting up the longitudinal wave reflection coefficient equation of viscoelastic medium.

Approximating unit 902: be used for utilizing weak viscoplasticity and similar medium it is assumed that obtain described longitudinal wave reflection coefficient approximation table Reach formula.

3rd sets up unit 903: is used for utilizing described longitudinal wave reflection coefficient approximate expression to set up and calculates viscoelastic fluid The equation of the factor, modulus of shearing and density.

As shown in figure 12, described second sets up unit 901 wraps a kind of fluid identification system embodiment 4 that the application provides Include:

First determines unit 1001: be used for determining viscoelastic fluid factor expression.

Second determines unit 1002: be used for determining viscoelastic media longitudinal wave reflection coefficient expressions.

Bring unit 1003 into: anti-for described viscoelastic fluid factor expression is brought into described viscoelastic media compressional wave Penetrate in coefficient expressions, obtain the longitudinal wave reflection coefficient equation with described viscoelastic fluid factor representation.

As shown in figure 13, described first inverting unit 802 wraps a kind of fluid identification system embodiment 5 that the application provides Include:

Converter unit 1101: for described viscoplasticity impedance model is carried out linear transformation.

Acquiring unit 1102: be used for the viscoplasticity impedance data that at least 3 angles described in obtaining are corresponding.

Second computing unit 1103: for according to described viscoplasticity impedance data and log, be calculated described viscous The elastic fluid factor, described modulus of shearing and described density.

As shown in figure 14, described acquiring unit 1102 includes a kind of fluid identification system embodiment 6 that the application provides:

Second extraction unit 1201: be used for extracting angle wavelet.

Second inverting unit 1202: be used for carrying out viscoplasticity Impedance Inversion, obtains at least 3 viscoplasticity resistances corresponding to angle Anti-data

It should be noted that each embodiment in this specification all uses the mode gone forward one by one to describe, each embodiment weight Point explanation is all the difference with other embodiments, and between each embodiment, identical similar part sees mutually. For device class embodiment, due to itself and embodiment of the method basic simlarity, so describe is fairly simple, relevant part ginseng See that the part of embodiment of the method illustrates.

Finally, in addition it is also necessary to explanation, in this article, the relational terms of such as first and second or the like be used merely to by One entity or operation separate with another entity or operating space, and not necessarily require or imply these entities or operation Between exist any this reality relation or order.And, term " includes ", " comprising " or its any other variant meaning Containing comprising of nonexcludability, so that include that the process of a series of key element, method, article or equipment not only include that A little key elements, but also include other key elements being not expressly set out, or also include for this process, method, article or The key element that equipment is intrinsic.In the case of there is no more restriction, statement " including ... " key element limited, do not arrange Except there is also other identical element in including the process of described key element, method, article or equipment.

Above a kind of Fluid Identification Method provided herein and system are described in detail, used herein Principle and the embodiment of the application are set forth by specific case, and the explanation of above example is only intended to help to understand this The method of application and core concept thereof；Simultaneously for one of ordinary skill in the art, according to the thought of the application, specifically All will change on embodiment and range of application, in sum, this specification content should not be construed as the application's Limit.

Claims (13)

1. a Fluid Identification Method, it is characterised in that the method includes:

The data of prestack CRP road collection will be converted into through protecting the data of the prestack CMP road collection that dampening information processes；

The superposition of data of at least 3 angles is extracted according to described prestack CRP road collection；

The viscoelastic fluid factor, modulus of shearing and density is calculated according to described superposition of data；

Utilize the described viscoelastic fluid factor, described modulus of shearing and described density that reservoir fluid is identified；

Wherein, the described viscoelastic fluid factor is by formula

$\begin{array}{c}{f}_{anela}={\mathrm{\ρ\α}}^2-{\mathrm{c\ρ\β}}^2\\ ={\mathrm{\ρV}}_P^2\left(1+\frac{2}{{\mathrm{\πQ}}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P}\right)-{\mathrm{c\ρV}}_S^2\left(1+\frac{2}{{\mathrm{\πQ}}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S}\right)\\ ={\mathrm{\ρV}}_P^2-{\mathrm{c\ρV}}_S^2+{\mathrm{\ρV}}_P^2\left(\frac{2}{{\mathrm{\πQ}}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P}\right)-{\mathrm{c\ρV}}_S^2\left(\frac{2}{{\mathrm{\πQ}}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S}\right)\\ =f_{ela}+{\mathrm{\Δf}}_Q\end{array}$

Determine；

In formula: f_anelaFor the described viscoelastic fluid factor；f_elaFor the fluid factor corresponding with elastic part,Δf_QFor viscoplasticity to f_elaDisturbance；V_P、V_SWith ρ be respectively medium velocity of longitudinal wave, shear wave velocity and Density；α is multiple velocity of longitudinal wave；β is multiple shear wave velocity；ω is earthquake dominant frequency；ω_rFor reference frequency；Q_PFor compressional wave quality factor； Q_SFor shear wave quality factor；C be dry rock P-S wave velocity ratio square；I is imaginary unit.

Fluid Identification Method the most according to claim 1, it is characterised in that the described calculating viscoelastic fluid factor, shearing Modulus and density include:

Set up viscoplasticity impedance model；

According to described viscoplasticity impedance inversion approach carry out decay viscoelastic fluid factor inverting, obtain described viscoelastic fluid because of Modulus of shearing sub, described and described density；

Wherein, described viscoplasticity impedance model of setting up includes: set up the longitudinal wave reflection coefficient equation of viscoelastic medium；Utilize weak viscous Elasticity and similar medium are it is assumed that obtain described longitudinal wave reflection coefficient approximate expression；Utilize described longitudinal wave reflection coefficient approximation table Reach formula and set up the calculating viscoelastic fluid factor, modulus of shearing and the equation of bulk density；

Wherein, described longitudinal wave reflection coefficient approximate expression is

$\begin{array}{c}{R}_{pp}\left(\mathrm{\θ}\right)\≈\frac{1}{4}\[1-\frac{{\mathrm{\γ}}_{dry}^2}{{\mathrm{\γ}}_{sat}^2}\left(1-\frac{2}{{\mathrm{\πQ}}_P}\log \left(30\right)\right)\]{\sec }^2\mathrm{\θ}\frac{{\mathrm{\Δf}}_{ane}}{f_{ane}}\\ +\frac{1}{{\mathrm{\γ}}_{sat}^2}\[1-\frac{2}{{\mathrm{\πQ}}_P}\log \left(30\right)\]\left(\frac{{\mathrm{\γ}}_{dry}^2}{4}{\sec }^2\mathrm{\θ}-2{\sin }^2\mathrm{\θ}\right)\frac{\mathrm{\Δ}\mathrm{\μ}}{\mathrm{\μ}}+\left(\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}\right)\frac{\mathrm{\Δ}\mathrm{\ρ}}{\mathrm{\ρ}}\end{array}$

In formula: γ_dryRepresent dry rock P-S wave velocity ratio；γ_satRepresent fluid saturated rocks P-S wave velocity ratio；f_ane, μ and ρ The mean value of the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides respectively；Δf_ane, Δ μ and Δ ρ table respectively Show the difference of the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides；Q_PFor compressional wave quality factor；θ is incident Angle；

Wherein, described viscoplasticity impedance model is

$QEI\left(\mathrm{\θ}\right)=f_{ane}^{a\left(\mathrm{\θ}\right)}{\mathrm{\μ}}^{b\left(\mathrm{\θ}\right)}{\mathrm{\ρ}}^{c\left(\mathrm{\θ}\right)}$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{2}\[1-\frac{{\mathrm{\γ}}_{dry}^2}{{\mathrm{\γ}}_{sat}^2}(1-\frac{2}{{\mathrm{\πQ}}_P}log(30\left)\right)\]{\sec }^2\mathrm{\θ}$

$b\left(\mathrm{\θ}\right)=\frac{1}{{\mathrm{\γ}}_{sat}^2}\[1-\frac{2}{{\mathrm{\πQ}}_P}log\left(30\right)\](\frac{{\mathrm{\γ}}_{dry}^2}{2}{\sec }^2\mathrm{\θ}-4{\sin }^2\mathrm{\θ})$

$c\left(\mathrm{\θ}\right)=1-\frac{{\sec }^2\mathrm{\θ}}{2}$

In formula: f_ane, the mean value of μ and ρ respectively the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides；γ_dry Represent dry rock P-S wave velocity ratio；γ_satRepresent fluid saturated rocks P-S wave velocity ratio；θ is incident angle；Q_PFor compressional wave Quality factor.

Fluid Identification Method the most according to claim 2, it is characterised in that the described longitudinal wave reflection setting up viscoelastic medium Coefficient equation specifically includes:

Determine viscoelastic fluid factor expression；

Determine viscoelastic medium longitudinal wave reflection coefficient expressions；

Described viscoelastic fluid factor expression is brought in described viscoelastic medium longitudinal wave reflection coefficient expressions, is used The longitudinal wave reflection coefficient equation of described viscoelastic fluid factor representation.

Fluid Identification Method the most according to claim 2, it is characterised in that described enter according to viscoplasticity impedance inversion approach Row decay viscoelastic fluid factor inverting, obtains the described viscoelastic fluid factor, described modulus of shearing and described density and includes:

Described viscoplasticity impedance model is carried out linear transformation；

The viscoplasticity impedance data that at least 3 angles described in acquisition are corresponding；

According to described viscoplasticity impedance data and log, it is calculated the described viscoelastic fluid factor, described modulus of shearing With described density.

Fluid Identification Method the most according to claim 4, it is characterised in that corresponding the gluing of at least 3 angles of described acquisition Elastic impedance data include:

Extract angle wavelet；

Carry out viscoplasticity Impedance Inversion, obtain at least 3 viscoplasticity impedance datas corresponding to angle.

Fluid Identification Method the most according to claim 5, it is characterised in that the extraction conditions of described extraction angle wavelet For:

E=∑ (s_i-d_i)²

Wherein, E is error energy, and s=r*w is synthetic seismogram, and r is reflectance factor, and w is the wavelet extracted, and d is earthquake note Record, selects optimal wavelet, makes error energy minimum simultaneously.

Fluid Identification Method the most according to claim 6, it is characterised in that described in carry out viscoplasticity Impedance Inversion, obtain Viscoplasticity impedance data that at least 3 angles are corresponding particularly as follows:

With the viscoplasticity impedance data body of each angle of Sparse Pulse Inversion inversion method, its objective optimization function is:

F=L_p(r(θ))+λL_q(s-d)+α^-1(ΔZ_trend)

Wherein, r (θ) is angle reflection coefficient sequence, and Δ Z is and the difference sequence of viscoplasticity impedance trend, and d is seismic channel sequence, s For synthetic seismogram sequence, λ is residual error weight factor, and α is trend weight factor, p and q is the L mould factor.

Fluid Identification Method the most according to claim 1, it is characterised in that described through protecting the prestack that dampening information processes The data of CMP road collection specially contain the prestack CMP road collection data affected by viscoplasticity.

9. a fluid identification system, it is characterised in that including:

Conversion unit: for the number by being converted into prestack CRP road collection through protecting the data of the prestack CMP road collection that dampening information processes According to；

First extraction unit: for extracting the superposition of data of at least 3 angles according to described prestack CRP road collection；

First computing unit: for calculating the viscoelastic fluid factor, modulus of shearing and density according to described superposition of data；

Recognition unit: be used for utilizing the described viscoelastic fluid factor, described modulus of shearing and described density that reservoir fluid is carried out Identify；

Wherein, the described viscoelastic fluid factor is by formula

$\begin{array}{c}{f}_{anela}={\mathrm{\ρ\α}}^2-{\mathrm{c\ρ\β}}^2\\ ={\mathrm{\ρV}}_P^2\left(1+\frac{2}{{\mathrm{\πQ}}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P}\right)-{\mathrm{c\ρV}}_S^2\left(1+\frac{2}{{\mathrm{\πQ}}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S}\right)\\ ={\mathrm{\ρV}}_P^2-{\mathrm{c\ρV}}_S^2+{\mathrm{\ρV}}_P^2\left(\frac{2}{{\mathrm{\πQ}}_P}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_P}\right)-{\mathrm{c\ρV}}_S^2\left(\frac{2}{{\mathrm{\πQ}}_S}\log \left(\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right)-\frac{i}{Q_S}\right)\\ =f_{ela}+{\mathrm{\Δf}}_Q\end{array}$

Determine；

In formula: f_anelaFor the described viscoelastic fluid factor；f_elaFor the fluid factor corresponding with elastic part,Δf_QFor viscoplasticity to f_elaDisturbance；V_P、V_SWith ρ be respectively medium velocity of longitudinal wave, shear wave velocity and Density；α is multiple velocity of longitudinal wave；β is multiple shear wave velocity；ω is earthquake dominant frequency；ω_rFor reference frequency；Q_PFor compressional wave quality factor； Q_SFor shear wave quality factor；C be dry rock P-S wave velocity ratio square；I is imaginary unit.

Fluid identification system the most according to claim 9, it is characterised in that described first computing unit includes:

First sets up unit: be used for setting up viscoplasticity impedance model；

First inverting unit: for carrying out, according to described viscoplasticity impedance inversion approach, the viscoelastic fluid factor inverting that decays, To the described viscoelastic fluid factor, described modulus of shearing and described density；

Wherein, described first sets up unit includes: second sets up unit, for setting up the longitudinal wave reflection coefficient side of viscoelastic medium Journey；Approximating unit, is used for utilizing weak viscoplasticity and similar medium it is assumed that obtain described longitudinal wave reflection coefficient approximate expression；The Three set up unit, are used for utilizing described longitudinal wave reflection coefficient approximate expression to set up and calculate the viscoelastic fluid factor, modulus of shearing Equation with density；

Wherein, described longitudinal wave reflection coefficient approximate expression is

$\begin{array}{c}{R}_{pp}\left(\mathrm{\θ}\right)\≈\frac{1}{4}\[1-\frac{{\mathrm{\γ}}_{dry}^2}{{\mathrm{\γ}}_{sat}^2}\left(1-\frac{2}{{\mathrm{\πQ}}_P}\log \left(30\right)\right)\]{\sec }^2\mathrm{\θ}\frac{{\mathrm{\Δf}}_{ane}}{f_{ane}}\\ +\frac{1}{{\mathrm{\γ}}_{sat}^2}\[1-\frac{2}{{\mathrm{\πQ}}_P}\log \left(30\right)\]\left(\frac{{\mathrm{\γ}}_{dry}^2}{4}{\sec }^2\mathrm{\θ}-2{\sin }^2\mathrm{\θ}\right)\frac{\mathrm{\Δ}\mathrm{\μ}}{\mathrm{\μ}}+\left(\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}\right)\frac{\mathrm{\Δ}\mathrm{\ρ}}{\mathrm{\ρ}}\end{array}$

In formula: γ_dryRepresent dry rock P-S wave velocity ratio；γ_satRepresent fluid saturated rocks P-S wave velocity ratio；f_ane, μ and ρ The mean value of the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides respectively；Δf_ane, Δ μ and Δ ρ table respectively Show the difference of the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides；Q_PFor compressional wave quality factor；θ is incident Angle；

Wherein, described viscoplasticity impedance model is

$QEI\left(\mathrm{\θ}\right)=f_{ane}^{a\left(\mathrm{\θ}\right)}{\mathrm{\μ}}^{b\left(\mathrm{\θ}\right)}{\mathrm{\ρ}}^{c\left(\mathrm{\θ}\right)}$

Wherein,

$a\left(\mathrm{\θ}\right)=\frac{1}{2}\[1-\frac{{\mathrm{\γ}}_{dry}^2}{{\mathrm{\γ}}_{sat}^2}(1-\frac{2}{{\mathrm{\πQ}}_P}log(30\left)\right)\]{\sec }^2\mathrm{\θ}$

$b\left(\mathrm{\θ}\right)=\frac{1}{{\mathrm{\γ}}_{sat}^2}\[1-\frac{2}{{\mathrm{\πQ}}_P}log\left(30\right)\](\frac{{\mathrm{\γ}}_{dry}^2}{2}{\sec }^2\mathrm{\θ}-4{\sin }^2\mathrm{\θ})$

$c\left(\mathrm{\θ}\right)=1-\frac{{\sec }^2\mathrm{\θ}}{2}$

In formula: f_ane, the mean value of μ and ρ respectively the viscoelastic fluid factor, modulus of shearing and the density of interface media of both sides；γ_dry Represent dry rock P-S wave velocity ratio；γ_satRepresent fluid saturated rocks P-S wave velocity ratio；θ is incident angle；Q_PFor compressional wave Quality factor.

11. fluid identification systems according to claim 10, it is characterised in that described second sets up unit includes:

First determines unit: be used for determining viscoelastic fluid factor expression；

Second determines unit: be used for determining viscoelastic medium longitudinal wave reflection coefficient expressions；

Bring unit into: for described viscoelastic fluid factor expression being brought into described viscoelastic medium longitudinal wave reflection coefficient table Reach in formula, obtain the longitudinal wave reflection coefficient equation with described viscoelastic fluid factor representation.

12. fluid identification systems according to claim 10, it is characterised in that described first inverting unit includes:

Converter unit: for described viscoplasticity impedance model is carried out linear transformation；

Acquiring unit: be used for the viscoplasticity impedance data that at least 3 angles described in obtaining are corresponding；

Second computing unit: for according to described viscoplasticity impedance data and log, be calculated described viscoelastic fluid The factor, described modulus of shearing and described density.

13. fluid identification systems according to claim 12, it is characterised in that described acquiring unit includes:

Second extraction unit: be used for extracting angle wavelet；

Second inverting unit: be used for carrying out viscoplasticity Impedance Inversion, obtains at least 3 viscoplasticity impedance datas corresponding to angle.