A pre-stack seismic direct inversion method for P-wave and S-wave velocity ratios

Abstract

The present invention discloses a pre-stack seismic direct inversion method for P-wave and S-wave velocity ratios, which is applied to the field of geophysical inversion. Currently, P-wave and S-wave velocity ratios are mostly calculated based on the P-wave and S-wave velocities obtained according to the AVO approximate formula, but the inversion accuracy is inevitably reduced. In response to the above problems, the present invention proposes a pre-stack seismic method for directly inverting the P-wave and S-wave velocity ratios, the method comprising: constructing a Bayesian nonlinear inversion algorithm based on the PSR-EI equation by deriving the elastic impedance (EI) equation of the P-wave and S-wave velocity ratio (PSR), and obtaining a robust P-wave and S-wave velocity ratio estimation by repeatedly iterating the objective function of the algorithm. The technology provided by the present invention can improve the inversion accuracy of the P-wave and S-wave velocity ratios, thereby improving the processing and interpretation efficiency and accuracy of actual work.

Description

Technical Field

[0001] The present invention belongs to the field of geophysical seismic data processing, and in particular relates to a technology for extracting longitudinal and transverse wave velocity ratios. Background Art

[0002] The P-wave velocity ratio (PVR) is a crucial parameter for describing formations. It eliminates the influence of density and can be used to indicate lithology. For example, lithology classification can be achieved using the P-wave velocity ratio and natural gamma ray crossplots. Furthermore, the P-wave velocity ratio is sensitive to changes in oil and gas reservoirs. In practical applications, the P-wave velocity ratio and compensated neutron crossplots are often used to locate gas zones. Therefore, the P-wave velocity ratio is a crucial parameter for intra-stratum reservoir location and lithology classification, and research on techniques for extracting the P-wave velocity ratio has naturally become a hot topic.

[0003] Currently, AVO approximate equation inversion is the most commonly used method for extracting the P-wave velocity ratio in prestack seismic technology. The Aki-Richards equation is a representative example. The Aki-Richards equation can be used to invert P-wave and S-wave velocities, and then calculate the P-wave velocity ratio. However, this distributed inversion method can lead to cumulative errors in the calculation results, inevitably reducing the inversion accuracy. With the rise of precision oil and gas exploration, more accurate predictions of exploration "sweet spots" are required. Conventional AVO attribute inversion is no longer sufficient to achieve the required accuracy, and people are eager to seek new prestack seismic techniques to improve inversion results.

[0004] While studying elastic impedance (EI), Cambois discovered that EI inversion provides information describing intra-layer properties, which is more accurate than AVO inversion of inter-layer properties. Consequently, EI inversion technology has gained considerable attention in recent decades. Connolly, by introducing the concept of acoustic impedance, first proposed a specific expression for the elastic impedance equation. However, this elastic impedance cannot be directly and quantitatively compared with acoustic impedance. To address this issue, Whitcombe performed a normalization preprocessing on the equation. When using the elastic impedance equation to invert P- and S-wave impedance, Lu discovered that the P- and S-wave velocity ratio in the equation significantly affects the inversion accuracy. By definition, the velocity of adjacent dielectric layers should be used for calculation. However, in practice, to simplify the nonlinear EI inversion problem into a linear one, the P- and S-wave velocity ratio is often assumed to be constant, ultimately compromising the reliability of the inversion results. To address this issue, Su Yun used the Zoeppritz equation as a correction function and, by optimizing the P- and S-wave velocity ratio, proposed an EI inversion method based on the P- and S-wave velocity ratio scanning method. This method requires an accurate velocity field, which is difficult to obtain. Subsequently, Du Qizhen analyzed the relationship between P-waves, S-waves, density, and the P-wave velocity ratio in the elastic impedance equation and proposed a joint P-wave and S-wave inversion method that iterates the P-wave and S-wave velocity ratio. This method first inverts P-waves and density, then inverts S-waves and the P-wave and S-wave velocity ratio. This method improves the stability of EI inversion to a certain extent, but increases the time complexity of the inversion. Thanks to the tireless exploration and research of scholars at home and abroad, EI inversion technology has gradually matured. However, currently, when EI inversion is applied to extract the P-wave and S-wave velocity ratio, either EI inversion is used to obtain P-wave and S-wave velocities, and then the ratio is calculated, resulting in large errors in the inverted P-wave and S-wave velocity ratio; or the method uses a two-step inversion method using P-wave and S-wave joint inversion, which improves the inversion accuracy but is time-consuming and inefficient. Summary of the Invention

[0005] In order to solve the above technical problems, the present invention proposes a pre-stack seismic direct inversion method for P-wave velocity ratio and S-wave velocity ratio, and establishes a Bayesian nonlinear inversion algorithm based on the derived PSR-EI equation.

[0006] The technical solution adopted by the present invention is: a pre-stack seismic direct inversion method of longitudinal and transverse wave velocity ratios, comprising:

[0007] S1. By deriving the differential relationship between the P-wave-velocity to S-wave-velocity Ratio (PSR) and the product of the P-wave-velocity and S-wave-velocity ratios, an approximate equation for the reflection coefficient of the P-wave-velocity ratio is established. Furthermore, the PSR-EI equation is derived by combining it with the elastic impedance (EI).

[0008] S2. By establishing the logarithmic relationship between elastic impedance and seismic records, a Bayesian nonlinear inversion algorithm based on the PSR-EI equation is constructed, and the initial form of the objective function is obtained by maximizing the posterior probability function;

[0009] S3. The objective function of step S2 is iteratively solved by the least square method to obtain a robust longitudinal and transverse wave velocity ratio parameter.

[0010] Furthermore, step S1 specifically includes the following steps:

[0011] S11. The expression of the Aki-Richards equation is as follows:

[0012]

[0013] Among them, V P 、V S represents the longitudinal wave velocity and the shear wave velocity respectively; ρ is the density; θ is the incident angle; k is the background parameter,

[0014] ΔV P , ΔV S , Δρ represents the property difference between adjacent dielectric layers;

[0015] S12. According to the differential theorem, the differential forms of the longitudinal and transverse wave velocity products and the longitudinal and transverse wave velocity ratios are derived:

[0016] △(V P V S )=V P △V S +V S △V P ,

[0017] S13. Solve the two equations in step S12 simultaneously to obtain the longitudinal wave differential and the transverse wave differential. The expressions are as follows:

[0018]

[0019] S14. Substitute the expression of step S13 into the expression of step S11 to obtain the reflection coefficient of the longitudinal and transverse wave velocity ratio, which is expressed as follows:

[0020]

[0021] in: represents the reflection coefficient of the ratio of longitudinal and transverse wave velocities; represents the longitudinal and transverse wave velocity product reflection coefficient;

[0022] △ρ / ρ=2(ρ1-ρ2) / (ρ1+ρ2), which represents the density reflection coefficient; ρ1 and ρ2 are the longitudinal and transverse wave velocity ratio, longitudinal and transverse wave velocity product and density of adjacent dielectric layers respectively;

[0023] S15. According to the difference quotient of elastic impedance:

[0024]

[0025] Where: EI is elastic impedance;

[0026] S16. Combine the equation of step S15 with the reflection coefficient equation of step S14 regarding the ratio of longitudinal and transverse wave velocities to obtain the basic expression:

[0027]

[0028] S17. Remainder in the equation of step S16 Can be split into The PSR-EI equation is obtained:

[0029]

[0030] Where: a(θ)=8ksin 2 θ, b(θ) = 1 + tan 2 θ-8ksin 2 θ, c(θ) = 1-4ksin 2 θ.

[0031] Furthermore, step S2 specifically includes the following steps:

[0032] S21. Introducing vertical travel time parameters and seismic gather parameters into the PSR-EI equation, we obtain the expression:

[0033]

[0034] Where: t represents the vertical travel time parameter, x represents the seismic gather parameter;

[0035] S22. Synthesis of seismic data, the expression is as follows:

[0036] syn(θ,t n ,x m )=G(θ)R(θ,t n ,x m )

[0037] Where: syn is seismic data; m is the number of seismic gather sequences; n is the number of travel time sequences; G is the seismic wavelet matrix;

[0038] S23. Elastic impedance and reflection coefficient have a recursive form on the logarithm, and the expression is:

[0039]

[0040] S24. Combining the expression in step S23 with step S22, derive the direct relationship between elastic impedance and seismic data:

[0041] syn(θ,t n ,x m )=G(θ)·W·ln(EI(θ,t,x m ))+noise

[0042] Where: ln(EI(θ,t,x m )) is the logarithmic form of the PSR-EI equation, and noise is the noise matrix;

[0043]

[0044] S25. In the Bayesian theory framework, the posterior probability expression is:

[0045]

[0046] Where: P(EI|syn) is the posterior function; P(EI) is the elastic impedance prior function; P(syn|EI) is the likelihood function.

[0047] Furthermore, step S25 specifically includes the following steps:

[0048] A1. Let the prior function in step S25 be the Cauchy distribution function, expressed as:

[0049]

[0050] Where: EI represents the standard deviation of elastic impedance; N is the number of stratigraphic sequences;

[0051] A2. The likelihood function in step S25 obeys Gaussian distribution and is expressed as:

[0052]

[0053] Where: n is the noise mean square error;

[0054] A3. Integrate steps A1 and A2 into the posterior probability function of step S25 to obtain the expression:

[0055]

[0056] Where: M is a constant coefficient.

[0057] Furthermore, step S3 specifically includes the following steps:

[0058] S31, convert the posterior probability expression in step S25 into the objective function form:

[0059]

[0060] S32, using least squares iterative optimization EI for the objective function of step S31 to obtain a PSR-EI data volume;

[0061] S33. Use the PSR-EI equation to construct a three-parameter AX=B solution, where A is the coefficient matrix related to the PSR-EI equation, X is the parameter vector of the P-wave velocity ratio to be solved, and B is the PSR-EI data volume obtained by inversion of the pre-stack seismic trace gather. The P-wave velocity ratio is obtained by solving the matrix.

[0062] The beneficial effects of the present invention are as follows: the present invention proposes a pre-stack seismic direct inversion method for P-wave and S-wave velocity ratios, constructs a Bayesian nonlinear inversion algorithm by deriving the PSR-EI equation, and then iterates the objective function using the weighted least squares method to obtain a PSR-EI data volume, which can be solved to obtain a robust P-wave and S-wave velocity ratio; the present invention proposes a pre-stack direct inversion method based on the PSR-EI equation, which has better performance than the indirect inversion of the Aki-Richards equation, and is well applied in actual pre-stack seismic data, thereby improving the processing and interpretation efficiency and accuracy of actual work. BRIEF DESCRIPTION OF THE DRAWINGS

[0063] Figure 1 A schematic diagram of reflection and transmission of an elastic interface provided by an embodiment of the present invention;

[0064] Figure 2 Schematic diagram of Bayesian inversion PSR-EI provided by an embodiment of the present invention;

[0065] in, Figure 2 (a) is the prior distribution, Figure 2 (b) is the likelihood function distribution, Figure 2 (c) is the posterior distribution;

[0066] Figure 3 It is a flow chart of the solution of the present invention;

[0067] Figure 4 The effect diagram of the two-dimensional model for inverting the longitudinal and transverse wave velocity ratio provided by the embodiment of the present invention:

[0068] in, Figure 4 (a) is the longitudinal and transverse wave velocity ratio profile obtained by inversion of the Aki-Richards equation. Figure 4 (b) is the longitudinal and transverse wave velocity ratio profile obtained by inversion using the present invention;

[0069] Figure 5The original pre-stack seismic slice 1 provided in the embodiment of the present invention;

[0070] Figure 6 The second original pre-stack seismic slice provided by the embodiment of the present invention;

[0071] Figure 7 Convert the original pre-stack seismic slice 1 into angle gather slice 1;

[0072] Figure 8 The original pre-stack seismic slice 2 is converted into the angle gather slice 2;

[0073] Figure 9 This is an effect diagram of the inverted P-wave and S-wave velocity ratios of an angle gather slice provided by an embodiment of the present invention;

[0074] Figure 10 This is a diagram showing the effect of the P-wave and S-wave velocity ratios obtained by inverting the second slice of the angle gather provided by an embodiment of the present invention; DETAILED DESCRIPTION

[0075] To facilitate those skilled in the art to understand the content of the present invention, the following technical aspects are now explained:

[0076] Elastic impedance equation

[0077] When a P wave enters the stratum, a series of reflections and transmissions will occur. The model is as follows: Figure 1 As shown, the P-wave incident angle, P-wave transmission angle, SV-wave reflection angle, and SV-wave transmission angle all satisfy Snell's law.

[0078] According to Snell's law, the relationship can be expressed as:

[0079]

[0080] When we know the velocity parameters and incident angles of the adjacent strata, we can calculate other angles. Based on the conditions of continuity of stress and displacement at the interface, we can deduce:

[0081]

[0082] This equation is the famous Zoeppritz equation. Because the Zoeppritz equation is inconvenient to apply in practice, in recent decades, many types of Zoeppritz approximations have been proposed to make its physical meaning more intuitive. In the 1980s, Aki and Richards studied the previous equations and obtained the following formula:

[0083]

[0084] Where k is the background parameter, ΔV P , ΔVS , Δρ are the longitudinal wave velocity difference, shear wave velocity difference and density difference of adjacent dielectric layers respectively.

[0085] In 1999, Connolly first proposed the concept of elastic impedance by introducing the concept of acoustic impedance. First, elastic impedance is calculated according to the method of acoustic impedance as follows:

[0086]

[0087] Combining equations (3) and (4), we get

[0088]

[0089] The right side of the equation can be expressed by the approximate relationship get

[0090] △ln(EI)=(1-4ksin 2 θ)△lnρ+sec 2 θ△lnV P -8ksin 2 θ△lnV S (6) After integrating both sides of equation (6) and exponentially converting them, we get

[0091]

[0092] Equation (7) is the famous Connolly elastic impedance equation. When the incident angle is equal to 0°, this equation is the expression of acoustic impedance.

[0093] Bayesian nonlinear inversion

[0094] Bayesian nonlinear inversion is an inversion method in the mathematical statistics school. The algorithm optimizes the objective function iteratively until the posterior probability is maximized. According to Bayesian theory, the posterior function can be expressed as:

[0095]

[0096] Where m is the model parameter, d is the observed data, P(m) is the prior distribution, P(d|m) is the likelihood function, and the larger the likelihood function value, the closer it is to the true value, and vice versa. P(m|d) is the required posterior distribution function. If only the shape of the function is considered, the constant term P(d) can be ignored, and Equation (8) can be simplified as:

[0097] P(m|d)∝P(d|m)P(m)(9)

[0098] like Figure 4 This is the verification effect of the Marmousi2 local model, where Figure 4(a) is the longitudinal and transverse wave velocity ratio profile obtained by indirect inversion using the Aki-Richards equation. Figure 4 (b) is the profile of the P-wave velocity ratio obtained by direct inversion of the PSR-EI equation of the present invention. By comparison, it is found that the direct inversion of the P-wave velocity ratio by the present invention has a better effect. Figure 5-6 These are two sets of pre-stack seismic data from a real work area in China. Figure 7-8 It is the result of converting pre-stack seismic data into angle gather data. Figure 9-10 The two sets of longitudinal and transverse wave velocity ratio results obtained by inverting the technology of the present invention from the angle gathers are used to verify the practicality of the technology.

[0099] Those skilled in the art will appreciate that the embodiments described herein are intended to aid the reader in understanding the principles of the present invention, and it should be understood that the scope of the present invention is not limited to such specific descriptions and embodiments. Various modifications and variations are readily apparent to those skilled in the art. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention are intended to be included within the scope of the claims.

Claims

1. A pre-stack seismic direct inversion method for longitudinal and transverse wave velocity ratios, characterized in that: include: S1. By deriving the differential relationship between the P-wave-velocity to S-wave-velocity Ratio (PSR) and the product of the P-wave-velocity and S-wave-velocity ratios, an approximate equation for the reflection coefficient of the P-wave-velocity ratio is established. Furthermore, the PSR-EI equation is derived by combining it with the elastic impedance (EI). S2. By establishing the logarithmic relationship between elastic impedance and seismic records, a Bayesian nonlinear inversion algorithm based on the PSR-EI equation is constructed, and the initial form of the objective function is obtained by maximizing the posterior probability function; S3. The objective function of step S2 is iteratively solved by the least square method to obtain a robust longitudinal and transverse wave velocity ratio parameter.

2. The method for pre-stack seismic direct inversion of longitudinal and transverse wave velocity ratios according to claim 1, characterized in that: Step S1 specifically includes the following sub-steps: S11. The expression of the Aki-Richards equation is as follows: Among them, V P 、V S represents the longitudinal wave velocity and the shear wave velocity respectively; ρ is the density; θ is the incident angle; k is the background parameter, ΔV P , ΔV S , Δρ represents the property difference between adjacent dielectric layers; S12. According to the differential theorem, the differential forms of the longitudinal and transverse wave velocity products and the longitudinal and transverse wave velocity ratios are derived: S13. Solve the two equations in step S12 simultaneously to obtain the longitudinal wave differential and the transverse wave differential. The expressions are as follows: S14. Substitute the expression of step S13 into the expression of step S11 to obtain the reflection coefficient of the longitudinal and transverse wave velocity ratio, which is expressed as follows: in: represents the reflection coefficient of the ratio of longitudinal and transverse wave velocities; represents the longitudinal and transverse wave velocity product reflection coefficient; △ρ / ρ=2(ρ1-ρ2) / (ρ1+ρ2), which represents the density reflection coefficient; ρ1 and ρ2 are the longitudinal and transverse wave velocity ratio, longitudinal and transverse wave velocity product and density of adjacent dielectric layers respectively; S15. According to the difference quotient of elastic impedance: Where: EI is elastic impedance; S16. Combine the equation of step S15 with the reflection coefficient equation of step S14 regarding the ratio of longitudinal and transverse wave velocities to obtain the basic expression: S17. Remainder in the equation of step S16 Can be split into The PSR-EI equation is obtained: Among them: a(θ)=8ksin 2 θ, b(θ)=1+tan 2 θ-8ksin 2 θ, c(θ)=1-4ksin 2 I.

3. The method for pre-stack seismic direct inversion of longitudinal and transverse wave velocity ratios according to claim 2, characterized in that: Step S2 includes the following sub-steps: S21. Introducing vertical travel time parameters and seismic gather parameters into the PSR-EI equation, we obtain the expression: Where: t represents the vertical travel time parameter, x represents the seismic gather parameter; S22. Synthesis of seismic data, the expression is as follows: syn(θ,t n ,x m )=G(θ)R(θ,t n ,x m ) Where: syn is seismic data; m is the number of seismic gather sequences; n is the number of travel time sequences; G is the seismic wavelet matrix; S23. Elastic impedance and reflection coefficient have a recursive form on the logarithm, and the expression is: S24. Combining the expression in step S23 with step S22, derive the direct relationship between elastic impedance and seismic data: syn(θ,t n ,x m )=G(θ)·W·ln(EI(θ,t,x m ))+noise Where: ln(EI(θ,t,x m )) is the logarithmic form of the PSR-EI equation, and noise is the noise matrix; S25. In the Bayesian theory framework, the posterior probability expression is: Where: P(EI|syn) is the posterior function; P(EI) is the elastic impedance prior function; P(syn|EI) is the likelihood function.

4. The method for pre-stack seismic direct inversion of longitudinal and transverse wave velocity ratios according to claim 3, characterized in that: Step S25 includes the following steps: A1. Let the prior function in step S25 be the Cauchy distribution function, expressed as: Where: EI represents the standard deviation of elastic impedance; N is the number of stratigraphic sequences; A2. The likelihood function in step S25 obeys Gaussian distribution and is expressed as: Where: n is the noise mean square error; A3. Integrate steps A1 and A2 into the posterior probability function of step S25 to obtain the expression: Where: M is a constant coefficient.

5. The method for pre-stack seismic direct inversion of longitudinal and transverse wave velocity ratios according to claim 4, characterized in that: Step S3 solution includes the following sub-steps: S31, convert the posterior probability expression in step S25 into the objective function form: S32, using least squares iterative optimization EI for the objective function of step S31 to obtain a PSR-EI data volume; S33. Use the PSR-EI equation to construct a three-parameter AX=B solution, where A is the coefficient matrix related to the PSR-EI equation, X is the parameter vector of the P-wave velocity ratio to be solved, and B is the PSR-EI data volume obtained by inversion of the pre-stack seismic trace gather. The P-wave velocity ratio is obtained by solving the matrix.

Images