A Kalman Filter Lateral Constraint Method and Apparatus Based on Accurate Zoeppritz Inversion
By introducing the Kalman filter algorithm and the exact Zoeppritz equation into AVO inversion, combined with adaptive incremental Kalman filtering, the problem of lack of lateral constraints in Bayesian inversion is solved, and the accuracy and stability of seismic inversion are improved.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-03-20
- Publication Date
- 2026-04-03
AI Technical Summary
In the existing technology, the AVO inversion method based on the Beyesian framework lacks lateral constraints, resulting in limited reliability of the inversion results. Furthermore, the state noise and measurement noise are unknown, affecting the inversion accuracy and stability.
A Kalman filter algorithm framework is introduced for lateral constraints, and combined with the exact Zoeppritz equation, an adaptive incremental Kalman filter algorithm is used to reduce the impact of noise. By combining adaptive incremental Kalman filtering with Bayesian nonlinear AVO inversion, longitudinal and lateral constraints on elastic parameters are achieved.
This improves the accuracy and stability of seismic inversion, reduces the impact of state noise and measurement noise, and achieves better inversion results.
Smart Images
Figure CN116125538B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of seismic inversion of geophysical oil and gas reservoirs, specifically to a Kalman filter lateral constraint method and apparatus based on accurate Zoeppritz inversion. Background Technology
[0002] While the AVO inversion method based on the Bayesian framework can achieve good inversion results for a single trace, its reliability is relatively limited due to the lack of consideration for lateral constraints. The Earth is laterally correlated, and adjacent seismic traces are interconnected. Based on this correlation, Kalman filtering can be used to laterally constrain the elastic parameters of the seismic traces. Combining Kalman filtering with Bayesian nonlinear AVO inversion constrains the elastic parameters of the seismic traces from two directions, achieving better inversion results.
[0003] Kalman filtering uses the mathematical model of the system to recursively estimate the motion state. Traditional Kalman filtering requires that the state noise and measurement noise are known, but in actual seismic AVO inversion projects, the state noise and measurement noise are mostly unknown. How to reduce the impact of state noise and measurement noise on the accuracy of AVO inversion results and improve the stability of AVO inversion has become an urgent problem to be solved. Summary of the Invention
[0004] The purpose of this invention is to address the lack of lateral constraints in current Bayesian AVO inversion methods by introducing a Kalman filter algorithm framework to implement lateral constraints on the AVO inversion results. Furthermore, it establishes measurement equations based on the accurate Zoeppritz equations and employs an optimized adaptive incremental Kalman filter algorithm to further reduce the impact of state noise and measurement noise on the accuracy and stability of the AVO inversion results.
[0005] To achieve the above-mentioned objectives, the present invention provides the following technical solution:
[0006] A Kalman filter lateral constraint method based on accurate Zoeppritz inversion includes:
[0007] Step 1: Establish state transition equations based on the lateral adjacency relationships of seismic traces;
[0008] Step 2: Construct a Bayesian inversion objective function based on the forward modeling operator formed by the exact Zoeppritz equation, and establish the measurement equation based on the Bayesian maximum a posteriori probability density function;
[0009] Step 3: Predict the elastic parameters of the current seismic trace from the elastic parameters of adjacent seismic traces based on the state transition equation;
[0010] Step 4: Correct the elastic parameters of the current seismic trace using the measured seismic data and the measurement equation.
[0011] This method integrates disciplines such as seismic inversion and mathematical statistics, organically combining mathematical calculation techniques and error statistical analysis techniques to effectively improve the accuracy of seismic inversion and overcome the shortcomings of Bayesian inversion technology in terms of lateral constraints.
[0012] In step 1, the state transition equation is:
[0013]
[0014] Where x k Represents the k-th seismic data, A k-1 W represents the transformation matrix of the state equation. k-1 Noise representing the state equation.
[0015] In step 2, the inversion objective function constructed using the forward operators based on the exact Zoeppritz equations is specifically as follows:
[0016] Min J(m) =
[0017] Where m is the elastic parameter vector to be inverted, d is the measured seismic data, and G is the nonlinear forward modeling operator constructed from the exact Zoeppritz equations. It is the noise covariance matrix. Let be the covariance matrix of the prior model, and μ be the mean vector.
[0018] The forward model corresponding to the nonlinear forward operator G formed by the exact Zoeppritz equation is:
[0019] = Simplified to: ,
[0020] Among them, the elastic parameter vector to be inverted V P1 V P2 Vs1 and Vs2 represent the longitudinal wave velocities of the upper and lower interfaces, respectively. , The density of the upper and lower interfaces. For observation data, The number of incident angles R is the number of parameters. pp Let W be the reflection coefficient and W be the wavelet matrix.
[0021] As a preferred embodiment of this application, when performing inversion using the exact Zoeppritz equation, intermediate variables are defined. ,in, , , , Reflection coefficient for intermediate variables ( The partial derivative of ) is:
[0022]
[0023]
[0024]
[0025]
[0026]
[0027] .
[0028] In step 2, the measurement equation is:
[0029]
[0030] Among them, the measured value ,
[0031] Transformation matrix of measurement equation ,
[0032] In the formula, I is the identity matrix; It is the prediction error of the state equation; is the forward operator; The expression is: R k For measurement error covariance; It is a priori model.
[0033] In step 3, predicting the elastic parameters of adjacent seismic traces includes calculating the prediction error covariance matrix:
[0034]
[0035] in Systematic errors representing the state transition process.
[0036] Step 4, correcting the elastic parameters of the adjacent seismic traces includes:
[0037] Calculate the Kalman gain:
[0038]
[0039] Calculate the measurement increment:
[0040] ;
[0041] Calculate the measurement residuals:
[0042] ;
[0043] Corrected forecast:
[0044] ;
[0045] Adaptive correction of error covariance in state equations:
[0046] ;
[0047] Update the error covariance matrix:
[0048]
[0049] In the formula The expression is:
[0050] .
[0051] This application also provides a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps of the above-described Kalman filter lateral constraint method based on accurate Zoeppritz inversion.
[0052] This application also provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the steps of the Kalman filter lateral constraint method based on accurate Zoeppritz inversion described above.
[0053] Compared with the prior art, the beneficial effects of the present invention are as follows:
[0054] In the scheme of this application, adaptive incremental Kalman filtering is used to constrain the elastic parameters of the seismic trace in the lateral direction based on the correlation of the elastic parameters of the subsurface medium. The adaptive incremental Kalman filtering is combined with Bayesian nonlinear AVO inversion to constrain the elastic parameters of the seismic trace from two directions, thereby achieving better inversion results.
[0055] To address the systematic error problem in the Kalman filter state equation, an adaptive method is used for correction. During the iterative process of the Kalman filter, the systematic error in the Kalman filter state equation is reduced, improving the lateral continuity of the inversion. Simultaneously, this adaptive method can effectively solve the problem of unknown state noise.
[0056] Incremental processing was performed on the measured values. The increment refers to the difference between adjacent seismic traces. Since the error values between adjacent seismic traces are similar, the error can be reduced by subtraction, thereby improving the inversion accuracy and solving the problem of unknown systematic errors caused by the stability of the measurement tools, the environment, and other factors. Attached Figure Description
[0057] Figure 1 Here is a flowchart of the Kalman filter lateral constraint method based on accurate Zoeppritz inversion in Example 1;
[0058] Figure 2 This is a flowchart of the adaptive incremental Kalman filter lateral constraint method in Example 1. Detailed Implementation
[0059] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments.
[0060] Therefore, the following detailed description of embodiments of the present invention is not intended to limit the scope of the claimed invention, but merely illustrates some embodiments of the invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.
[0061] It should be noted that, unless otherwise specified, the embodiments and features and technical solutions in the present invention can be combined with each other.
[0062] Example 1: See Figure 1 As shown in this embodiment, a Kalman filter lateral constraint method based on accurate Zoeppritz inversion is provided, including:
[0063] Step 1: Establish state transition equations based on the lateral adjacency relationships of seismic traces;
[0064] Step 2: Construct a Bayesian inversion objective function based on the forward modeling operator formed by the exact Zoeppritz equation, and establish the measurement equation based on the Bayesian maximum a posteriori probability density function;
[0065] Step 3: Predict the elastic parameters of the current seismic trace from the elastic parameters of adjacent seismic traces based on the state transition equation;
[0066] Step 4: Correct the elastic parameters of the current seismic trace using the measured seismic data and the measurement equation.
[0067] In step 1, the state transition equation is:
[0068]
[0069] Where x k Represents the k-th seismic data, A k-1 W represents the transformation matrix of the state equation. k-1 Noise representing the state equation.
[0070] Step 2 consists of two parts: Zoeppritz-based forward modeling and Bayesian-based AVO inversion. These two parts will be explained in detail below.
[0071] Forward modeling theory:
[0072] According to seismic wave field theory, when seismic waves are incident non-perpendicularly on the interface of an elastic medium, the waves not only generate converted transverse waves but also refract back into the incident medium and simultaneously transmit into the next layer of medium. This phenomenon can be described by the Zoeppritz equation:
[0073] , recorded as B; where the parameters of the upper and lower interfaces for non-perpendicular P-wave incidence are respectively and ; , The angles of reflection and transmission of the incident P-wave are denoted as . , To convert the reflection angle and transmission angle of the S-wave; , , These represent the reflection and transmission coefficients of the converted transverse and longitudinal waves, respectively.
[0074] The AVO inversion form with noise can be obtained using the seismic convolution model:
[0075]
[0076] Where α, β, and ρ correspond to the P-wave velocity, S-wave velocity, and density parameters, respectively; The observation data... In theoretical synthetic data ( , Expanding the Taylor series at (where is the wavelet matrix), and discarding the second-order higher-order terms, gives:
[0077] = Simplified to: ,
[0078] Among them, the elastic parameter vector to be inverted V P1 V P2Vs1 and Vs2 represent the longitudinal wave velocities of the upper and lower interfaces, respectively. , The density of the upper and lower interfaces. For observation data, The number of incident angles Let Rpp be the number of parameters, Rpp be the reflection coefficient, and W be the wavelet matrix. Expression:
[0079]
[0080] This is the forward operand G of the exact Zoeppritz equation.
[0081] As can be seen from the above equation, when performing inversion using the exact Zoeppritz equation, this application requires obtaining the reflection coefficient with respect to the elastic parameters. The partial derivative. To facilitate the solution, intermediate variables are defined. ,in, , , , Reflection coefficient for intermediate variables ( The partial derivative of ) can be obtained by the following formula:
[0082]
[0083]
[0084]
[0085]
[0086]
[0087]
[0088]
[0089]
[0090]
[0091] In the above equation, A, B, and R all correspond to matrices in the Zoeppritz equation, and a and b correspond to elements in A and B, respectively. Reflection coefficient For parameters The differential can be obtained by using the chain rule. and Based on the formula Ultimately, the partial derivatives of the reflection coefficient with respect to each inversion parameter can be obtained:
[0092]
[0093]
[0094]
[0095]
[0096]
[0097] .
[0098] Pre-stack AVO Bayesian inversion:
[0099] Due to the influence of factors such as the uncertainty of observational data, theoretical accuracy, and data completeness, the inversion process inevitably suffers from multiple solutions. Geophysical inversion aims to find the optimal solution under different objective functions by utilizing different inversion strategies and standards. Bayesian theory uses prior distributions and likelihood functions to represent the posterior probability density distribution of model parameters, combining known prior information and observational data to obtain the best-fit solution for the model parameters.
[0100]
[0101] Statistical inversion based on Bayesian theory treats the inversion parameters m as following a certain prior probability distribution. The random variable is obtained by maximizing the posterior probability density function. To obtain the corresponding solution. Among them, This is called the prior probability. The prior probabilities of model parameters are usually obtained through methods such as well logging or geostatistics. This is equivalent to introducing additional constraint information into the inversion process, acting as a regularization term and effectively reducing the ambiguity of the inversion. Meanwhile, P(d|m) is called the likelihood function, describing the conditional probability density of observed data d given model parameters m. It represents the degree of approximation between the model value calculated through forward modeling and the true value given model parameters m.
[0102] With the forward operand G of the exact Zoeppritz equation:
[0103]
[0104] We construct a nonlinear inversion objective function within a Bayesian theoretical framework. We assume that both the prior probability and the likelihood function follow a Gaussian distribution.
[0105]
[0106]
[0107] in, For the nonlinear forward modeling operator constructed from the new exact Zoeppritz equation, This is the noise covariance matrix. Assuming the noise data is uncorrelated with the observed data, the noise covariance matrix can be simplified to an identity matrix. = I. Where I is the identity matrix, This represents the noise mean square error. Let be the covariance matrix of the prior model, and μ be the mean vector.
[0108] Substituting the prior distribution expression and likelihood function expression into Bayes' theorem, the problem of finding the maximum a posteriori probability solution can be transformed into finding the minimum value of the objective function shown below:
[0109] Min J(m) = .
[0110] Adaptive incremental Kalman filtering:
[0111] Nonlinear AVO inversion based on a Bayesian framework provides vertical constraints on the inversion results for a single trace, but it does not consider the lateral continuity of the subsurface geological body. Adaptive incremental Kalman filtering is employed to impose lateral constraints on the elastic parameters of the seismic traces. By combining adaptive incremental Kalman filtering with nonlinear AVO inversion within a Bayesian framework, both vertical and lateral constraints are achieved, resulting in better inversion results.
[0112] Adaptive incremental Kalman filtering can effectively address the problem of inaccurate process noise estimation. The filtering process comprises two parts: state transition and measurement equations, such as... Figure 2 As shown, the specific steps include the following detailed steps:
[0113] Step 1: Obtain measurement values through Bayesian inversion Transformation matrix of measurement equation :
[0114] a. Measured value:
[0115]
[0116] b. Transformation matrix of the measurement equation:
[0117]
[0118] In the formula: I is the identity matrix; It is the prediction error of the state equation; is the forward operator; The expression is: Measurement error covariance R k Obtained through well logging data; It is a priori model.
[0119] Step 2: Set up the state equations:
[0120]
[0121] in Represents the k-th seismic data. The transformation matrix representing the state equation. Noise representing the state equation.
[0122] Step 3: Calculate the error covariance matrix :
[0123]
[0124] in This represents systematic error.
[0125] Step 4: Calculate the Kalman gain:
[0126] ;
[0127] Step 5: Calculate the measurement increment:
[0128]
[0129] in This represents the increment of the measured value.
[0130] Step Six: Calculate the measurement residuals:
[0131]
[0132] Step 7: Correct the predicted value:
[0133]
[0134] At this point, incremental processing has been performed on the measured values. The increment refers to the difference between adjacent seismic traces. Since the errors between adjacent seismic traces are similar, the error can be reduced and the inversion accuracy can be improved by subtracting the difference. This can effectively solve the problem of unknown systematic errors caused by the stability of the measuring tools, the environment, and other factors.
[0135] Step 8: Adaptive Correction of Error Covariance in State Equation:
[0136]
[0137] Step 9: Update the error covariance matrix:
[0138]
[0139] In the formula The expression is as follows:
[0140] .
[0141] By using an adaptive method to correct the systematic error problem of the Kalman filter state equation, the systematic error of the Kalman filter state equation is reduced during the continuous iteration of the Kalman filter, thereby improving the lateral continuity of the inversion. This adaptive method can effectively solve the problem of unknown state noise.
[0142] The above embodiments are only used to illustrate the present invention and are not intended to limit the technical solutions described herein. Although the present invention has been described in detail with reference to the above embodiments, the present invention is not limited to the specific embodiments described above. Therefore, any modifications or equivalent substitutions to the present invention, as well as all technical solutions and improvements that do not depart from the spirit and scope of the invention, are covered within the scope of the claims of the present invention.
Claims
1. A Kalman filter lateral constraint method based on accurate Zoeppritz inversion, characterized in that, include: Step 1: Establish state transition equations based on the lateral adjacency relationships of seismic traces; Step 2: Construct a Bayesian inversion objective function based on the forward modeling operator formed by the exact Zoeppritz equation, and establish the measurement equation based on the Bayesian maximum a posteriori probability density function; Step 3: Predict the elastic parameters of the current seismic trace from the elastic parameters of adjacent seismic traces based on the state transition equation; Step 4: Correct the elastic parameters of the current seismic trace using measured seismic data and the measurement equation; wherein: The forward model corresponding to the nonlinear forward operator G formed by the exact Zoeppritz equation is: = Simplified to: , Among them, the elastic parameter vector to be inverted V P1 、V P2 The longitudinal wave velocities are at the upper and lower interfaces. Vs 1 Vs 2 The transverse wave velocities are at the upper and lower interfaces. , The density of the upper and lower interfaces. For observation data, The number of incident angles For the number of parameters, R pp Let W be the reflection coefficient and W be the wavelet matrix.
2. The method as described in claim 1, characterized in that: The state transition equation is: ， in x k express k-th earthquake data, x k-1 This represents the (k-1)th seismic trace. A k-1 The transformation matrix representing the state equation. W k-1 This represents the noise in the state equation.
3. The method as described in claim 2, characterized in that: The inversion objective function constructed using the forward operators based on the exact Zoeppritz equations is as follows: Min J(m)= , in, m Let be the vector of elastic parameters to be inverted, d be the measured seismic data, and G be the nonlinear forward modeling operator constructed from the exact Zoeppritz equations. It is the noise covariance matrix. Let be the covariance matrix of the prior model, and μ be the mean vector.
4. The method as described in claim 1, characterized in that: When performing inversion using the exact Zoeppritz equation, define intermediate variables. ,in, , , , Reflection coefficient for intermediate variables The partial derivatives are: ， ， ， ， ， ， in, .
5. The method as described in claim 2, characterized in that: The measurement equation is: ； Among them, the measured value ; Transformation matrix of measurement equation , In the formula, I is the identity matrix; It is the prediction error of the state equation; is the forward operator; The expression is: , R k For measurement error covariance; It is a priori model; v k This is the error term.
6. The method as described in claim 5, characterized in that: Predicting the elastic parameters of adjacent seismic traces includes calculating the prediction error covariance matrix: ； in Systematic errors representing the state transition process.
7. The method as described in claim 6, characterized in that: The elastic parameters of the adjacent seismic traces are corrected as follows: Calculate the Kalman gain: ； Calculate the measurement increment: ； Calculate the measurement residuals: ； Corrected forecast: ； Adaptive correction of error covariance in state equations: ； Update the error covariance matrix: ； In the formula The expression is: 。 8. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by a processor, it implements the steps of the method according to any one of claims 1-7.
9. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the steps of the method according to any one of claims 1-7.