CN101446645B - Method for determining fluid by utilizing seismic fluid impedance - Google Patents

Abstract

The invention discloses a method for determining a fluid by utilizing seismic fluid impedance. The method comprises the following steps: exciting and recording a seismic wave to obtain an angle gather; eliminating signal frequency distortion caused by a normal moveout correction stretching effect; extracting longitudinal wave reflectivity attributes and weighted transverse wave reflectivity attributes from the compensated angle gather; computing improved fluid factors; obtaining logging data to obtain a fluid impedance curve; performing iterative inversion to obtain a section of the fluid impedance; qualitatively analyzing fluid anomalies in pores by the fluid factors, and describing spatial distribution of the fluid by the fluid impedance. The method helps accurately extract the fluid factors and directly identify components of the fluid in the pores without an assumption of a velocity ratio of the longitudinal wave to the transverse wave as well as the relation between density and the longitudinal wave velocity, which reduces ambiguity of gas layer identification, realizes direct conversion from the fluid factors to the fluid impedance and improves the precision of gas reservoir identification.

Description

A kind ofly utilize the earthquake fluid impedance to carry out the method that fluid is determined

Technical field

The invention belongs to the hydro carbons detection technique of geophysical prospecting for oil, is a kind ofly to utilize the earthquake fluid impedance to carry out the method that fluid is determined.

Background technology

Geological data just plays an increasingly important role in petroleum exploration and development, and it is from being used for seeking the simple function of structure type hydrocarbon-bearing pool to lithology identification and the multi-functional transformation of fluid detection in the past.Utilizing seismic data to carry out the most frequently used technology of fluid detection is to utilize the entrained amplitude information of seismic signal that penetrates underground medium.Amplitude information typically refers to amplitude for the pre-stack seismic data and detects gas-bearing formation with variation (AVO) feature of offset distance; And for the poststack seismic data, gas-bearing formation shows as on section: bright spot, dim spot, flat spot and phase reversal etc." bright spot " technology starts from the seventies in last century, first milestone of the direct hydrocarbon detection of being known as (Hilterman, 2001).Behind the gassiness of oil gas field reservoir speed and density descend very fast, on stacked seismic data, have very big reflection amplitude promptly " bright spot " occur.But be not all to show as " bright spot " after all reservoir gas-bearing, and some non-gas-bearing formations also can form " bright spot " phenomenon, utilize " bright spot " technology for detection gas-bearing formation that very big multi-solution is arranged therefore merely.

The AVO technology that occur the eighties in last century is second milestone (Hilterman, 2001) on the direct hydrocarbon detection history.On prestack CDP road collection, the AVO of gas sand response shows as amplitude and increases (III class), reversal of poles (II class) and reduce (I class) three category features (Ostrander, 1984 with offset distance; Rutherford etc., 1989).P, G ripple attribute according to the theoretical inverting of AVO can be indicated gas reservoir, have bigger negative value as III class gas sand P, G value, show very strong on the occasion of feature on the P*G section; And II class gas sand approaches zero because of its P, G value, makes that P, G value and combination section thereof are very little with respect to background value difference, and therefore this method lost efficacy for the identification of II class sandstone.

" the fluid factor " is a kind of special Fluid Identification Method based on the AVO attribute, Smith and Gidlow (1987) combination Aki-Richards equation (Aki﹠amp; Richards, 1980) and mud stone line (Catagna etc., 1985) discern in the geological data unusual as hydrocarbon reservoir, they disclose can generate P ripple and S wave reflection rate road by the AVO approximate equation, obtain " the fluid factor " road by means of mud stone line combination P ripple and S wave reflection rate road, the existence of gas reservoir can be indicated in this fluid factor road.Gas sand in the clastic sequence produces higher fluid factor amplitude, and other reflection has lower amplitude.But the fluid factor formula of Smith and Gidlow (1987) definition need provide accurate transverse wave reflection rate and suppose longitudinal and transverse wave velocity ratio in advance, and this is difficult to satisfy in practical operation.

By the attribute that the theoretical inverting of AVO obtains, no matter be that P ripple, G involve its combination, the fluid factor, although they under specific situation, can both discern gas-bearing reservoir effectively, but these method majorities be qualitatively, descriptive.Amplitude waveform attributes as other hydro carbons indications, subsurface interface thicker for reservoir, that do not have to interfere reflects, they can indicate the vertical and cross direction profiles of gas reservoir effectively, but discern gas reservoir with the amplitude of wave form size and can bring error when gas-bearing formation is thin, the gas reservoir border can be difficult to accurately determine owing to thin-bed effect.

The elastic parameter inversion method that Goodway etc. (1997) propose has been avoided the influence of thin layer tuning effect convection cell recognition methods.The elastic parameter method is referred to as λ-μ-ρ method again, this method is at first concentrated P ripple and the S wave reflection rate extracted according to the AVO approximate formula of Fatti (1994) from pre-stack seismic P radio frequency channel, then these two reflectivity are carried out inverting, obtain P ripple and S wave impedance, product by P ripple and S wave impedance utilization elastic parameter formula conversion generation Lame's constant and density is λ ρ and μ ρ again, utilizes λ ρ and parameters such as μ ρ and Poisson ratio to discern gas reservoir.Russell etc. (2003) utilization Biot-Gassmann theory is extracted the fluid factor of their definition of institute and is indicated gas reservoir from the P ripple of seismic inversion and S Acoustic Impedance Data.Although from P ripple and S wave impedance, extract the qualitative question that elastic parameter has avoided amplitude waveform Attribute Recognition fluid to exist, but because seismic data is the band limit, need provide the P ripple simultaneously during by P ripple and S wave reflection rate road inverting P ripple and S wave impedance and S ripple well-log information sets up initial model and the initial model that provides is provided inversion result, the error of inverting certainly will influence the calculating of behindness parameter.And λ ρ and Poisson ratio still have multi-solution when test fluid, when particularly II class AVO water bearing sand and I class gas sand occur together, are difficult to differentiate.In addition, during by Fatti formulas Extraction S wave reflection rate, need suppose longitudinal and transverse wave velocity ratio in advance, this is difficult to satisfy in practical operation equally.

Summary of the invention

The object of the invention is to provide a kind of and extracts and the closely-related fluid factor of pore fluid and fluid impedance property parameters from earthquake and well-logging, determine that oil-gas space distributes, improve the earthquake fluid impedance that utilizes of oil gas detecting accuracy and carry out the method that fluid is determined.

The present invention is achieved through the following technical solutions, and concrete steps comprise:

1) excite and write down seismic event, earthquake data is routinely carried out high-fidelity and is handled, and is formed for common midpoint (CMP) the road collection of amplitude after with the normal-moveout correction of offset distance mutation analysis;

The described seismic data high-fidelity of step 1) is handled and is comprised geometry diffusion compensation, earth's surface-consistent amplitude restoration, earth's surface-consistent wavelet shaping deconvolution, the denoising of the prestacks guarantor width of cloth, refractive wave static correction, the residual static correction of reflection wave earth's surface-consistent, prestack road collection regularization processing.

2) utilize following formula that the common midpoint in the step 1) (CMP) road collection is transformed in the angle domain, obtain angular-trace gather,

$\sin \mathrm{\θ}=\frac{{\mathrm{xv}}_p}{v_{\mathrm{rms}}^2{(t_0^2+\frac{x^2}{v_{\mathrm{rms}}^2})}^{1/2}}---\left(1\right)$

In the formula, x is a geophone offset, V _pBe the interval velocity on stratum, t ₀Be zero-offset whilst on tour, V _RmsRoot-mean-square velocity for the stratum;

3) distort to eliminate the signal frequency that is caused owing to the NMO stretching effect in order to the collection working frequency compensation of following degree road, method diagonal angle;

The described distortion of step 3) is meant that seismic signal is along with incident angle increases and the artificial relaxation phenomenon of frequency reduction.

The described frequency compensation of step 3) is meant carries out the frequency division amplitude compensation to the earthquake reflective data on the different angles territory.

The amplitude variations of the different offset distances of following formula correction makes statistical average spectral amplitude coupling zero-offset road a long way,

${\mathrm{\α}}_x=\frac{\∂t}{{\∂t}_0}\≈\frac{t_0}{t}(1-\frac{2(t-t_0)}{V_{\mathrm{NMO}}}\frac{\∂V_{\mathrm{NMO}}}{{\∂t}_0})---\left(2\right)$

α in the formula _xBe compressibility coefficient, t ₀Be the zero-offset round trip travel-time, t is the round trip travel-time, V _NMOBe the normal moveout correction root-mean-square velocity;

Following formula is eliminated offset distance far away and because frequency reducing that nmo stretching caused,

$w_{\mathrm{\θ}}\left(t\right)=w\left(\frac{t}{\cos \mathrm{\θ}}\right)\stackrel{\mathrm{FT}}{\→}\tilde{w}(\mathrm{\ω}\·\cos \mathrm{\θ})---\left(3\right)$

W in the formula _θBe the back wavelet that stretches, w is the zero-offset wavelet, and θ is a reflection angle, and t is the travel-time, and ω is an angular frequency.

4) utilize the angle road of following Fatti formula after the compensation concentrate to extract longitudinal wave reflection rate attribute

With add

Power transverse wave reflection rate attribute

$R\left(\mathrm{\θ}\right)=\frac{1}{2}\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}(1+{\tan }^2\mathrm{\θ})-4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}{\sin }^2{\mathrm{\θ}}^{\′}---\left(4\right)$

The leaching process of two attributes in the formula is found the solution following system of equations:

$\left[\begin{array}{c}\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}\\ 4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}\end{array}\right]={\left[\begin{array}{cc}\frac{1+{\tan }^2{\mathrm{\θ}}_1}{2}& -{\sin }^2{\mathrm{\θ}}_1\\ \frac{1+{\tan }^2{\mathrm{\θ}}_2}{2}& -{\sin }^2{\mathrm{\θ}}_2\\ \mathrm{MM}& \\ \frac{1+{\tan }^2{\mathrm{\θ}}_n}{2}& -{\sin }^2{\mathrm{\θ}}_n\end{array}\right]}^{-1}\left[\begin{array}{c}R\left({\mathrm{\θ}}_1\right)\\ R\left({\mathrm{\θ}}_2\right)\\ M\\ R\left({\mathrm{\θ}}_n\right)\end{array}\right]---\left(5\right)$

θ in the formula _i, i=1,2 ..., n is the incident angle of road, the angle collection geological data that obtains in the step (3), R (θ _i) be road, the angle collection P wave reflection geological data that obtains in the step (3), V _p, V _sAnd I _p, I _sBe respectively the average speed value and the average impedance value of the P ripple and the S ripple of reflection aspect levels, Δ I _pAnd Δ I _sIt is respectively the variable quantity of P ripple, S ripple and the density of reflection aspect levels;

5) the longitudinal wave reflection rate attribute of above inverting and weighting transverse wave reflection rate attribute are utilized the fluid factor after the following formula computed improved:

$R_f=\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}-c4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}---\left(6\right)$

R in the formula _fBe the fluid reflectivity, be also referred to as the fluid factor, c is that constant is 2-3;

6) adopt conventional method to obtain well-log information, handle the back and from longitudinal and transverse ripple and densimetric curve, obtain fluid impedance curve;

Described step 6 is that the well-log information that longitudinal and transverse ripple and density are arranged is directly changed.When not having the SWAL data, can utilize rock physics dielectric model and Gassmann equation to construct plan shear wave data and be converted to fluid impedance curve again.

7) fluid impedance curve of extracting with the fluid factor and the step 6 of step 5 extraction is carried out iterative inversion, obtains the fluid impedance section;

Described step 7 iterative inversion algorithm is:

(1) utilizes fluid impedance curve Fluid Computation impedance rate of change R _f(t), again it and seismic wavelet w (t) convolution, obtain synthetic fluid factor record d (t)=R _f(t) * w (t);

(2) demarcate writing down with the well lie of the fluid factor of from seismic data, extracting by the synthetic fluid factor of well,

Above-mentioned correlation calibration is that the well lie of the synthetic fluid factor road that obtains of well with the fluid factor of extraction is presented at, more known geology reference lamina is labeled on the fluid factor road of synthetic and actual extracting,

Above-mentioned correlation calibration is that the geologic horizon in the well logging is mapped on the earthquake fluid factor road;

(3) adopt distance weighted usually interpolation method to set up initial fluid impedance model by the well logging fluid impedance;

(4) from fluid factor section, obtain final fluid impedance section, on the fluid impedance section, carry out the horizontal tracking that fluid distributes;

8) utilize the unusual variation of fluid in the fluid factor attribute qualitative analysis hole, utilize the space distribution of fluid impedance quantitative description fluid

The described fluid qualitative analysis of step 8) is on fluid factor section, and the absolute value of fluid factor attribute amplitude is relatively large to be gas-bearing reservoir, and amplitude goes to zero and then is non-gas-bearing reservoir.

The described fluid quantitative of step 8) is described and is: on the fluid impedance section, and the fluid impedance value minimum of gas-bearing reservoir, but not the fluid impedance value maximum of gas-bearing reservoir.Determine that according to the size of fluid impedance value reservoir is an oily or moisture, laterally tracks out the cross direction profiles and the thickness distribution of oil gas then.

The present invention compares with traditional method, does not need to make hypothesis to the p-and s-wave velocity ratio and to the relation between density and the velocity of longitudinal wave, has guaranteed the accuracy of the fluid factor attribute of extraction;

The new attribute of fluid impedance of the present invention is very responsive to the fluid in the hole, can Direct Recognition pore fluid composition, and reduced the explanation multi-solution of gas-bearing formation identification;

Fluid impedance inversion technique of the present invention has realized by the fluid factor having avoided the error transfer effect directly to the conversion of fluid impedance, has improved the precision of gas reservoir identification.

Description of drawings

Description of drawings of the present invention is as follows:

Fig. 1 is the simplified flow chart of fluid impedance inversion technique of the present invention;

Fig. 2 is fluid impedance and p-wave impedance, shear wave impedance, Poisson ratio σ and the elastic constant λ ρ X plot of three kinds of AVO gas reservoir types;

Fig. 3 is fluid impedance, λ ρ, p-wave impedance and the shear wave impedance of embodiment of the invention well objective interval and the X plot of WELL LITHOLOGY curve (Gr);

Fig. 4 is the impedance of another embodiment of the present invention resident fluid and p-wave impedance and λ ρ X plot;

Fig. 5 a is the fluid factor section that the CDP road collection after the normal moveout correction extracts, and T81 reflection seismic layer position is equivalent to Cretaceous System gas-bearing formation top reflectogram;

Fig. 5 b is that (Fig. 5 is the fluid impedance section that obtains of inverting a) by the fluid factor;

Fig. 6 a is the zero-offset compressional wave sectional view that the CDP road collection after the normal moveout correction extracts;

Fig. 6 b is that (Fig. 6 is the p-wave impedance section that obtains of inverting a) by compressional wave record.

Embodiment

The present invention fully utilizes earthquake, well-log information data, and by extracting fluid impedance from well-log information, and the fluid factor of extracting from seismic data comes the inverting fluid impedance, carried out the lateral prediction of reservoir fluid by the fluid impedance section.It requires the pre-stack seismic data is carried out the high-fidelity processing that amplitude keeps.

The present invention is achieved through the following technical solutions, and concrete steps comprise:

1) excite and write down seismic event, earthquake data is routinely carried out high-fidelity and is handled, and is formed for common midpoint (CMP) the road collection of amplitude after with the normal-moveout correction of offset distance mutation analysis;

The described seismic data high-fidelity of step 1) is handled and is comprised geometry diffusion compensation, earth's surface-consistent amplitude restoration, earth's surface-consistent wavelet shaping deconvolution, the denoising of the prestacks guarantor width of cloth, refractive wave static correction, the residual static correction of reflection wave earth's surface-consistent, prestack road collection regularization processing.

2) utilize following formula that the common midpoint in the step 1) (CMP) road collection is transformed in the angle domain, obtain angular-trace gather,

$\sin \mathrm{\θ}=\frac{{\mathrm{xv}}_p}{v_{\mathrm{rms}}^2{(t_0^2+\frac{x^2}{v_{\mathrm{rms}}^2})}^{1/2}}---\left(1\right)$

In the formula, x is a geophone offset, V _pBe the interval velocity on stratum, t ₀Be zero-offset whilst on tour, V _RmsRoot-mean-square velocity for the stratum;

3) distort to eliminate the signal frequency that is caused owing to the NMO stretching effect in order to the collection working frequency compensation of following degree road, method diagonal angle;

The described distortion of step 3) is meant that seismic signal is along with incident angle increases and the artificial relaxation phenomenon of frequency reduction.

The described frequency compensation of step 3) is meant carries out the frequency division amplitude compensation to the earthquake reflective data on the different angles territory.

The amplitude variations of the different offset distances of following formula correction makes statistical average spectral amplitude coupling zero-offset road a long way,

${\mathrm{\α}}_x=\frac{\∂t}{{\∂t}_0}\≈\frac{t_0}{t}(1-\frac{2(t-t_0)}{V_{\mathrm{NMO}}}\frac{\∂V_{\mathrm{NMO}}}{{\∂t}_0})---\left(2\right)$

α in the formula _xBe compressibility coefficient, t ₀Be the zero-offset round trip travel-time, t is the round trip travel-time, V _NMOBe the normal moveout correction root-mean-square velocity;

Following formula is eliminated offset distance far away and because frequency reducing that nmo stretching caused,

$w_{\mathrm{\θ}}\left(t\right)=w\left(\frac{t}{\cos \mathrm{\θ}}\right)\stackrel{\mathrm{FT}}{\→}\tilde{w}(\mathrm{\ω}\·\cos \mathrm{\θ})---\left(3\right)$

W in the formula _θBe the back wavelet that stretches, w is the zero-offset wavelet, and θ is a reflection angle, and t is the travel-time, and ω is an angular frequency.

4) utilize the angle road of following Fatti formula after the compensation concentrate to extract longitudinal wave reflection rate attribute

With add

Power transverse wave reflection rate attribute

$R\left(\mathrm{\θ}\right)=\frac{1}{2}\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}(1+{\tan }^2\mathrm{\θ})-4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}{\sin }^2{\mathrm{\θ}}^{\′}---\left(4\right)$

The leaching process of two attributes in the formula is found the solution following system of equations:

$\left[\begin{array}{c}\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}\\ 4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}\end{array}\right]={\left[\begin{array}{cc}\frac{1+{\tan }^2{\mathrm{\θ}}_1}{2}& -{\sin }^2{\mathrm{\θ}}_1\\ \frac{1+{\tan }^2{\mathrm{\θ}}_2}{2}& -{\sin }^2{\mathrm{\θ}}_2\\ \mathrm{MM}& \\ \frac{1+{\tan }^2{\mathrm{\θ}}_n}{2}& -{\sin }^2{\mathrm{\θ}}_n\end{array}\right]}^{-1}\left[\begin{array}{c}R\left({\mathrm{\θ}}_1\right)\\ R\left({\mathrm{\θ}}_2\right)\\ M\\ R\left({\mathrm{\θ}}_n\right)\end{array}\right]---\left(5\right)$

θ in the formula _i, i=1,2 ..., n is the incident angle of road, the angle collection geological data that obtains in the step (3), R (θ _i) be road, the angle collection P wave reflection geological data that obtains in the step (3), V _p, V _sAnd I _p, I _sBe respectively the average speed value and the average impedance value of the P ripple and the S ripple of reflection aspect levels, Δ I _pAnd Δ I _sIt is respectively the variable quantity of P ripple, S ripple and the density of reflection aspect levels;

5) the longitudinal wave reflection rate attribute of above inverting and weighting transverse wave reflection rate attribute are utilized the fluid factor after the following formula computed improved:

$R_f=\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}-c4\frac{V_s^2}{V_p^2}\frac{{\mathrm{\ΔI}}_s}{{\stackrel{\‾}{I}}_s}---\left(6\right)$

R in the formula _fBe the fluid reflectivity, be also referred to as the fluid factor, c is that constant is 2-3;

Conventional fluid factor notion is proposed by Smith and Gidlow (1987), the velocity of longitudinal wave of water bearing sand, mud stone and shale that they provide according to (1985) such as Castagna and " mud stone line " equation: the V between the shear wave velocity _p=1360+1.16V _sDerive fluid factor formula, promptly $R_f=R_p-c\frac{V_s}{V_p}{R}_s,$ R in the formula _pAnd R _sBe P ripple and S wave reflection rate, constant c=1.16.Need suppose longitudinal and transverse wave velocity ratio in advance when utilizing this formula to extract the fluid factor, this is difficult to satisfy in practical operation." mud stone line " of the present invention equation: $I_p^2=a+{\mathrm{cI}}_s^2,$ Derive formula 6, corresponding in the constant c of formula 6 and this formula.

Not only need not suppose longitudinal and transverse wave velocity ratio in advance with the formula 6 Fluid Computation factors time, and fluid impedance inverting quantitative description fluid that can performing step 7 distributes, this is that conventional fluid factor analysis is not available.

6) adopt conventional method to obtain well-log information, handle the back and from longitudinal and transverse ripple and densimetric curve, obtain fluid impedance curve;

Described step 6 is that the well-log information that longitudinal and transverse ripple and density are arranged is directly changed.When not having the SWAL data, can utilize rock physics dielectric model and Gassmann equation to construct plan shear wave data and be converted to fluid impedance curve again.

7) fluid impedance curve of extracting with the fluid factor and the step 6 of step 5 extraction is carried out iterative inversion, obtains the fluid impedance section;

Described step 7 iterative inversion algorithm is:

(1) utilizes fluid impedance curve Fluid Computation impedance rate of change R _f(t), again it and seismic wavelet w (t) convolution, obtain synthetic fluid factor record d (t)=R _f(t) * w (t);

(2) demarcate writing down with the well lie of the fluid factor of from seismic data, extracting by the synthetic fluid factor of well,

Above-mentioned correlation calibration is that the well lie of the synthetic fluid factor road that obtains of well with the fluid factor of extraction is presented at, more known geology reference lamina is labeled on the fluid factor road of synthetic and actual extracting,

Above-mentioned correlation calibration is that the geologic horizon in the well logging is mapped on the earthquake fluid factor road;

(3) adopt distance weighted usually interpolation method to set up initial fluid impedance model by the well logging fluid impedance;

(4) from fluid factor section, obtain final fluid impedance section, on the fluid impedance section, carry out the horizontal tracking that fluid distributes;

8) utilize the unusual variation of fluid in the fluid factor attribute qualitative analysis hole, utilize the space distribution of fluid impedance quantitative description fluid

The described fluid qualitative analysis of step 8) is on fluid factor section, and the absolute value of fluid factor attribute amplitude is relatively large to be gas-bearing reservoir, and amplitude goes to zero and then is non-gas-bearing reservoir.

The described fluid quantitative of step 8) is described and is: on the fluid impedance section, and the fluid impedance value minimum of gas-bearing reservoir, but not the fluid impedance value maximum of gas-bearing reservoir.Determine that according to the size of fluid impedance value reservoir is an oily or moisture, laterally tracks out the cross direction profiles and the thickness distribution of oil gas then.

Fig. 2 is fluid impedance and p-wave impedance, shear wave impedance, Poisson ratio σ and the elastic constant λ ρ X plot of three kinds of AVO gas reservoir types, and the typical models parameter is taken from even each homogeny model that Blangy (1994) provides.Can see that for three types of (III class, II class and I class) AVO reservoirs, the fluid impedance of gas sand all shows as low value, less than 110; And the fluid impedance of water sand and mud stone all shows as high value, greater than 110.(Fig. 2 a) has only the III class to separate with mud stone and water sand, and II class and I class gas sand are positioned at mud stone and water sand interval, is difficult to identification for p-wave impedance (Ip).For shear wave impedance (Is) (Fig. 2 b), III class, II class are identical with the distributed area of mud stone and water sand with I class gas sand distributed area, are difficult to identification.For Poisson ratio σ (Fig. 2 c), III class, II class and I class gas sand can separate with mud stone well, but II class and I class gas sand and water sand twists together, and are difficult to identification.For λ ρ (Fig. 2 d), III class and II class can be separated with mud stone and water sand, and I class gas sand is positioned at mud stone and water sand interval, are difficult to identification.Therefore, conventional gas reservoir recognition methods such as wave impedance, Poisson ratio and λ ρ can not determine uniquely that fluid, testing result depend on AVO gas reservoir type when using separately, and each parameter carries out all having multi-solution when fluid is explained; And fluid impedance can identify fluid well, and for three class AVO response, fluid impedance all shows as low value, and the precision of gas reservoir identification obviously improves.

Fig. 3 is fluid impedance, λ ρ, p-wave impedance and the shear wave impedance of embodiment of the invention well objective interval and the X plot of WELL LITHOLOGY curve (Gr).No matter can see, be p-wave impedance (I _p), shear wave impedance (I _s) or λ ρ (Fig. 3 a-c), gas sand, dry sand rock twist together mutually, are difficult to make a distinction.And at fluid impedance (I _Fluid) with lithology Gr X plot (Fig. 3 d), find out that obviously the fluid impedance of gas sand is minimum, can make a distinction mutually between it and dry sand and the mud stone.Therefore, utilize fluid impedance can identify fluid in the petroclastic rock well.

Fig. 4 is another enforcement resident fluid impedance and p-wave impedance and λ ρ X plot of the present invention.Can see that gassiness and vaughanite are at p-wave impedance (I _p) and λ ρ distributed area on overlapped, carry out having multi-solution when gas reservoir detects with them.And at fluid impedance (I _Fluid) having monambiguity, the dolomitic fluid impedance of gassiness is minimum, and it and vaughanite can make a distinction mutually.Therefore, utilize fluid impedance can identify fluid in the carbonate reservoir well.

Fig. 5 a is the fluid factor section that the CDP road collection after the normal moveout correction extracts, and T81 reflection seismic layer position is equivalent to the reflection of Cretaceous System gas-bearing formation top.Bored 2 mouthfuls of wells on this section, wherein being positioned at the crow of structure lower position, to join 1 well be a high gas rate well, and what be positioned at the higher position of structure is dry-well according to drawing 2 wells, and crow ginseng 1 well is corresponding to the high value of the fluid factor as can be seen, and dry-well is corresponding to fluid factor low value.Fig. 5 b is that (Fig. 5 is the fluid impedance section that obtains of inverting a) by the fluid factor.Can see that the fluid impedance value of gas-bearing formation shows as low value (non-green stripes), to according to draw 2 wells gradually transition be high value, to according to drawing 2 gas well gas layers disappearance (bottle green background).As seen the fluid impedance section has reduced the multi-solution of sound impedance section explanation gas-bearing formation, has improved horizontal tracking precision.

Fig. 6 a is the zero-offset compressional wave section that the CDP road collection after the normal moveout correction extracts, the dry-well that is positioned at the higher position of structure is much better than according to gas well crow ginseng 1 well that draws near the reflection amplitude of 2 wells zone of interest than being positioned at the structure lower position, this strong amplitude is owing to lithological change causes, and nonfluid causes unusually.As seen, the method for conventional " bright spot " identification gas reservoir here is impracticable.Fig. 6 b is that (Fig. 6 is the p-wave impedance section that obtains of inverting a) by compressional wave record.Can see, near near the wave impedance distribution similar characteristics of (comply with and draw 2 wells) the sound impedance value of (crow ginseng 1 well) and the non-gas-bearing formation gas-bearing formation, thereby, be the cross direction profiles that is difficult to track out gas reservoir from wave impedance section.

Claims (10)

1. one kind is utilized the earthquake fluid impedance to carry out the method that fluid is determined, it is characterized in that realizing by following concrete steps:

1) excite and write down seismic event, earthquake data is routinely carried out high-fidelity and is handled, and is formed for the common midpoint CMP road collection of amplitude after with the normal-moveout correction of offset distance mutation analysis;

2) utilize following formula that the common midpoint CMP road collection in the step 1) is transformed in the angle domain, obtain angular-trace gather,

$\sin \mathrm{\θ}=\frac{{\mathrm{xv}}_p}{v_{\mathrm{rms}}^2{(t_0^2+\frac{x^2}{v_{\mathrm{rms}}^2})}^{1/2}}---\left(1\right)$

In the formula, x is a geophone offset, V _pBe the interval velocity on stratum, t ₀Be zero-offset travel-time, V _RmsRoot-mean-square velocity for the stratum;

3) compensation of angular-trace gather working frequency is distorted with the signal frequency of eliminating owing to NMO stretching was caused; The amplitude variations of the different offset distances of following formula correction makes statistical average spectral amplitude coupling zero-offset road a long way,

${\mathrm{\α}}_x=\frac{\∂t}{{\∂t}_0}\≈\frac{t_0}{t}(1-\frac{2(t-t_0)}{V_{\mathrm{NMO}}}\frac{{\∂V}_{\mathrm{NMO}}}{\∂t_0})---\left(2\right)$

α in the formula _xBe compressibility coefficient, t ₀Be the zero-offset travel-time, t is the travel-time, V _NMOBe the normal moveout correction root-mean-square velocity;

Following formula is eliminated offset distance far away and because frequency reducing that NMO stretching caused,

$w_{\mathrm{\θ}}\left(t\right)=w\left(\frac{t}{\cos \mathrm{\θ}}\right)\stackrel{\mathrm{FT}}{\→}w(\mathrm{\ω}\·\cos \mathrm{\θ})---\left(3\right)$

W in the formula _θBe the back wavelet that stretches, w is the zero-offset wavelet, and θ is a reflection angle, and t is the travel-time, and ω is an angular frequency;

4) utilize following Fatti formula from the angular-trace gather after the compensation, to extract longitudinal wave reflection rate attribute

With weighting transverse wave reflection rate attribute

(4) leaching process of two attributes in the formula is found the solution following system of equations:

$\left[\begin{array}{c}\frac{{\mathrm{\ΔI}}_p}{{\stackrel{\‾}{I}}_p}\\ 4\frac{V_s^2}{V_p^2}\frac{\mathrm{\Δ}{I}_s}{{\stackrel{\‾}{I}}_s}\end{array}\right]{\left[\begin{array}{cc}\frac{1+{\tan }^2{\mathrm{\θ}}_1}{2}& -{\sin }^2{\mathrm{\θ}}_1\\ \frac{1+{\tan }^2{\mathrm{\θ}}_2}{2}& {-\sin }^2{\mathrm{\θ}}_2\\ \·\·\·\·\·\·& \\ \frac{1+{\tan }^2{\mathrm{\θ}}_n}{2}& {-\sin }^2{\mathrm{\θ}}_n\end{array}\right]}^{-1}\left[\begin{array}{c}R\left({\mathrm{\θ}}_1\right)\\ R\left({\mathrm{\θ}}_2\right)\\ \·\·\·\\ R\left({\mathrm{\θ}}_n\right)\end{array}\right]---\left(5\right)$

θ in the formula _i, i=1,2 ..., n is the incident angle of the angular-trace gather geological data that obtains in the step 3), R (θ _i) be the angular-trace gather P wave reflection geological data that obtains in the step 3), V _p, V _sAnd

Be respectively the average speed value and the average impedance value of the P ripple and the S ripple of reflection aspect levels, Δ I _pAnd Δ I _sBe respectively the P ripple of reflection aspect levels, the variable quantity of S ripple;

5) above longitudinal wave reflection rate attribute and weighting transverse wave reflection rate attribute are utilized the fluid factor after the following formula computed improved:

$R_f=\frac{\mathrm{\Δ}{I}_p}{{\stackrel{\‾}{I}}_p}-c4\frac{V_s^2}{V_p^2}\frac{\mathrm{\Δ}{I}_s}{{\stackrel{\‾}{I}}_s}---\left(6\right)$

R in the formula _fBe the fluid reflectivity, be also referred to as the fluid factor, c is that constant is 2-3;

6) adopt conventional method to obtain well-log information, handle the back and from longitudinal and transverse ripple and densimetric curve, obtain fluid impedance curve;

7) fluid impedance curve of extracting with the fluid factor and the step 6) of step 5) extraction is carried out iterative inversion, obtains the fluid impedance section;

The iterative inversion algorithm is:

(1) utilizes fluid impedance curve Fluid Computation factor R _f(t), again it and seismic wavelet w (t) convolution, obtain synthetic fluid factor record d (t)=R _f(t) * w (t);

(2) synthetic fluid factor record is demarcated with the well lie of the fluid factor of extracting from seismic data,

(3) adopt distance weighted usually interpolation method to set up initial fluid impedance model by the well logging fluid impedance;

(4) from fluid factor section, obtain final fluid impedance section, on the fluid impedance section, carry out the horizontal tracking that fluid distributes;

8) utilize the unusual variation of fluid in the fluid factor attribute qualitative analysis hole, utilize the space distribution of fluid impedance quantitative description fluid.

2. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, it is characterized in that the described seismic data high-fidelity of step 1) is handled and comprised geometry diffusion compensation, earth's surface-consistent amplitude restoration, earth's surface-consistent wavelet shaping deconvolution, the denoising of the prestacks guarantor width of cloth, refractive wave static correction, the residual static correction of reflection wave earth's surface-consistent, prestack road collection regularization processing.

3. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, it is characterized in that, the described distortion of step 3) is meant that seismic signal increases along with reflection angle θ and artificial relaxation phenomenon that frequency reduces.

4. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, it is characterized in that the described frequency compensation of step 3) is meant carries out the frequency division amplitude compensation to the earthquake reflective data on the different angles territory.

5. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, it is characterized in that the described impedance curve of step 6) is directly converted by the well-log information of longitudinal and transverse ripple and density.

6. utilize the earthquake fluid impedance to carry out the method that fluid is determined according to claim 1 or 5, it is characterized in that the described impedance curve of step 6) is to be converted to by utilizing rock physics dielectric model and Gassmann equation to construct plan shear wave data again.

7. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, the described demarcation of step 7) is that the geologic horizon in the well logging is mapped on the earthquake fluid factor record.

8. according to claim 1ly utilize the earthquake fluid impedance to carry out the method that fluid is determined, the described demarcation of step 7) is that the well lie of the synthetic fluid factor record that obtains with the fluid factor of extracting is presented at, and more known geology reference lamina is labeled on the fluid factor record of synthetic and actual extracting.

9. describedly utilize the earthquake fluid impedance to carry out the method that fluid is determined according to claim 1 or 7, it is characterized in that, the described qualitative analysis of step 8) is on fluid factor section, and the absolute value of fluid factor attribute amplitude is relatively large to be gas-bearing reservoir, and amplitude goes to zero and then is non-gas-bearing reservoir.

10. describedly utilize the earthquake fluid impedance to carry out the method that fluid is determined according to claim 1 or 7, it is characterized in that, the described quantitative description of step 8) is: on the fluid impedance section, the fluid impedance value minimum of gas-bearing reservoir, but not the fluid impedance value maximum of gas-bearing reservoir, determine that according to the size of fluid impedance value reservoir is an oily or moisture, laterally tracks out the cross direction profiles and the thickness distribution of oil gas then.

Images