Seismic inversion method under joint constraints of physical model and prior information

Abstract

A seismic inversion method under joint constraints of a physical model and prior information, comprising: extracting seismic wavelets on the basis of seismic data, and determining a wavelet amplitude scaling factor; collecting statistics about priori information of impedance parameters; using seismic structure interpretation data and well logging data to establish an initial impedance parameter model; obtaining a linear approximate equation on the basis of an interface weak elasticity difference hypothesis, and using the linear approximate equation to forward simulate a seismic gather and calculate an inversion residual error; using a generalized linear inversion idea to rewrite a target function into a function related to the impedance parameters, using an iterative reweighted least square algorithm to obtain the impedance parameters, and updating the impedance parameters; and repeating the steps above until the inversion residual error reaches the requirement or reaches the maximum number of iterations, and outputting a final processing result. The inversion method can achieve stable and high-precision prediction of impedance parameters.

Description

A seismic inversion method constrained by a physical model and prior information Technical Field

[0001] This invention relates to the field of seismic exploration technology, and in particular to a seismic inversion method that uses physical models and prior information as constraints. Background Technology

[0002] Seismic exploration is a method of inferring subsurface geological structures by artificially generating seismic waves and recording their propagation characteristics underground. The propagation speed, reflection, and transmission characteristics of seismic waves in different media carry rich geological information, and seismic impedance is one of the key parameters describing this information. Seismic impedance is defined as the product of rock density and its P-wave velocity; it directly reflects the physical properties of the rock, such as porosity, saturation, and rock type. These properties are crucial for identifying oil and gas reservoirs and other economically valuable minerals. With the continued growth of global energy demand, the search for and development of new oil and gas resources has become increasingly urgent. The depletion of traditional oil and gas fields and the scarcity of new resources require exploration technologies to be more efficient and accurate. Seismic impedance inversion can provide high-resolution subsurface images, helping geologists accurately identify the boundaries, fluid properties, reservoir quality, and hydrocarbon potential of oil and gas reservoirs, thereby improving exploration success rates and reducing drilling risks and costs. The core objective of seismic impedance inversion is to convert seismic data into parameters that directly reflect the physical properties of rocks; this step is indispensable for reservoir prediction and reservoir characterization. As a comprehensive parameter, wave impedance can bridge the gap between the physical propagation characteristics of seismic waves and the actual geological properties of rocks. Through wave impedance inversion, geologists can indirectly estimate reservoir parameters such as porosity and permeability, which is crucial for assessing reservoir productivity and formulating development strategies. With in-depth research, various wave impedance inversion methods have emerged, including but not limited to trace integral inversion, generalized linear inversion, iterative inversion, and nonlinear inversion. Each method has its advantages and limitations. For example, trace integral inversion is simple and direct but has limited accuracy, while generalized linear inversion, although highly accurate, is susceptible to high-frequency noise. Model-based inversion techniques break through the traditional limitations of seismic resolution, theoretically achieving the same resolution as well logging data. However, the high and low frequency components provided by well logging information and the model in the inversion results lead to multiple solutions, and are constrained by the number of wells and the distribution of the well network. To reduce the multiple solutions, some prior information is usually introduced during the inversion process. However, the existing prior information (such as the initial model) constraint method directly calculates errors with the inversion results, leading to reduced accuracy.

[0003] In summary, current research on seismic impedance inversion methods suffers from the following problems: 1. Trace integral seismic impedance inversion is affected by seismic wavelet band limiting, resulting in low resolution and significant cumulative error due to the influence of initial impedance values; 2. The method of first inverting reflection coefficients and then recursively calculating impedance is greatly affected by seismic data noise, resulting in significant cumulative error; 3. Existing model-based seismic impedance inversion methods introduce prior constraints such as initial models to reduce the ambiguity of inversion results, but directly using these constraints reduces the resolution of the inversion results; 4. Deterministic seismic impedance inversion methods can only provide a unique inversion result, making it impossible to evaluate the reliability of the result and increasing the risk of subsequent interpretation. Summary of the Invention

[0004] The purpose of this invention is to overcome the shortcomings of existing technologies and provide a seismic inversion method that combines physical models and prior information constraints. It adopts a model-based inversion strategy, introduces a physical model to ensure that the prediction results meet the observed seismic data, designs prior information constraint terms for the initial model to calculate the error between the low frequency of the inversion results and the prior information of the initial model, reduces the ambiguity of the inversion results and improves the resolution of the inversion results, and gives the uncertainty of the inversion results while predicting the inversion results based on a Bayesian framework, thus meeting the requirements of high-precision seismic exploration and fine reservoir characterization.

[0005] The objective of this invention is achieved through the following technical solution: a seismic inversion method jointly constrained by a physical model and prior information, comprising the following steps:

[0006] Step 110: Extract wavelets using actual seismic data, perform forward modeling of seismic gathers based on well logging data and post-stack seismic reflection coefficient equations, and determine the amplitude scaling factor by combining actual well-side seismic data;

[0007] Step 120: Extract the impedance parameters and their mean values ​​based on all well logging data in the work area, and statistically calculate their variance and vertical variation function matrix;

[0008] Step 130: Use seismic data layer interpretation data and well logging data to establish an initial impedance parameter model in the time domain;

[0009] Step 140: Starting from the accurate post-stack reflection coefficient equation, and based on the assumption of weak elasticity difference at the interface, a linear approximate equation is derived. Based on the initial impedance parameter model and simplified equation in the time domain, the post-stack seismic gather is simulated in forward modeling. The inversion residual is directly calculated from the forward modeling record and the actual record.

[0010] Step 150: Based on Bayes' principle, introduce sparse prior constraints and initial model prior information constraints to construct the inversion objective function in the sense of maximum a posteriori probability. Solve the inversion objective function using the generalized linear inversion idea to obtain the solution expression for the impedance parameter.

[0011] Step 160: Based on the solution formulas for the model parameters and the inversion residuals, calculate the model parameters using the iterative reweighted least squares algorithm;

[0012] Step 170: Repeat steps 140, 150 and 160 above, and control the maximum number of iterations by the inversion residual to obtain the optimal impedance parameter inversion result.

[0013] Specifically, in step 110, the reflection coefficient is calculated using the post-stack seismic reflection coefficient equation with well logging data as input.

[0014] In the formula, m represents the logging impedance parameter data, and r represents the calculated reflection coefficient;

[0015] The reflection coefficient is then convolved with the extracted seismic wavelet to obtain a seismic gather. This gather is compared with the actual well-side seismic gather, and the amplitude scaling factor is calculated and applied to the extracted seismic wavelet to achieve amplitude matching between the simulated and actual records.

[0016] In the formula, x represents the synthetic seismic record. Let w represent the convolution operator, and w represent the extracted wavelet.

[0017] Specifically, in step 120, the required impedance parameters are obtained through well logging data analysis, the autocorrelation coefficient of the impedance parameters is calculated, the variance matrix is ​​constructed, and the vertical variation function is calculated based on the impedance parameters to form a prior distribution function of the impedance parameters that conforms to the work area.

[0018] The formula for calculating the vertical variation function is as follows:

[0019] Where υ represents the variogram value calculated based on well data, h represents the hysteresis distance, N represents the length of the logging data, and z(x i ) represents the impedance parameter value at a certain position i.

[0020] Specifically, in step 130, seismic tectonic interpretation data is used to establish a geological model based on the sedimentary model. Well logging data is interpolated and extrapolated according to the tectonic model to obtain the initial impedance parameter model for each survey line. The impedance parameter model is established using spatial interpolation methods. First, the data of each layer is interpolated using the scattered interpolation method to complete the geological layer modeling. Then, the impedance parameters are interpolated laterally according to the geological layers to calculate the impedance parameter value at each underground point, thus completing the task of initial impedance parameter modeling.

[0021] Specifically, in step 140, based on the assumption of weak elasticity difference, a linear approximation equation is derived:

[0022] In the formula, m represents the impedance parameter data, r represents the calculated reflection coefficient, and ln represents taking the logarithm of the data;

[0023] Then, using the initial impedance parameter model in the time domain as input, the reflection coefficient vector is directly calculated using simplified equations. The seismic wavelet is convolved with the reflection coefficient to obtain the post-stack seismic gather, and the difference between this and the actual seismic gather is used to obtain the inversion residual.

[0024] In the formula, K represents the wavelet matrix constructed based on w, D represents the differential operator, and m represents the impedance parameter data.

[0025] Specifically, in step 150, the initial model is used to constrain the inversion results, and it is assumed that the seismic data noise follows a Gaussian distribution, the prior information constraint terms of the initial model follow a Gaussian distribution, and the model follows a sparse distribution. Then, the inversion likelihood function and the prior probability distribution satisfy the Gaussian distribution and the joint distribution of Gaussian and sparse, respectively. According to Bayes' principle, the posterior probability distribution function is obtained by combining the inversion likelihood function and the prior distribution function. The objective function of the inversion is determined according to the posterior probability. The iterative solution formula is obtained by differentiating the objective function with respect to the model parameters.

[0026] Let the impedance model parameters be z. T =(z1,z2,…,z) n ) T The observed earthquake data is x T =(x1,x2,…,x n ) T According to Bayesian theory, given post-stack seismic data, the problem of inverting subsurface impedance parameters can be reduced to solving a posterior probability function:

[0027] Where P(x)=∫p(x|z)p(z)dz is the normalization factor, which can be regarded as a constant, P(x|z) is the likelihood function, and P(z) is the prior probability distribution. It is assumed that the likelihood function follows a Gaussian distribution:

[0028] In the above formula, Γ is the forward modeling operator, and C X It is the covariance matrix of the noise, N x Given the length of the observed data, based on Bayesian theory and combined with the prior model, after canceling the normalization constant, the posterior probability distribution of the inverted post-stack seismic data can be expressed as:

[0029] In the formula, Γ is the forward modeling operator, Φ(z) represents the sparse constraint term, and z0 represents the initial impedance parameter. Represents the smoothing operator, C X It is the covariance matrix of the noise, C z The model covariance matrix can be represented by the variance statistically derived from well logging data. α sum of variation function matrix υ Φ(z) is obtained through the Kronecker product. The specific expression is as follows:

[0030] Φ(z)=||Dz||0 (9)

[0031] In the formula, D represents the difference operator. If three-point smoothing is used, The expression can be written as:

[0032] Depending on the specific degree of smoothness, it can be modified according to equation (10). The expression,

[0033] The optimal solution can be obtained by finding the maximum value of equation (8), which is equivalent to finding the solution corresponding to the minimum value of the objective function in the following equation.

[0034] Applying the above objective function to the impedance parameter z Taking the derivative and setting it to zero, since the 0-norm is not differentiable, we use the iterative reweighted least squares algorithm to solve for the model parameters, which yields the iterative formula for updating the model parameters:

[0035] In the formula, L = 0.5·W·D represents the forward modeling operator; y = ln(z) represents the logarithm of the impedance parameter; W represents the wavelet matrix; X represents the seismic record; and I is N. x ×N x The unit matrix, N x λ is the length of the observation data; diag represents the operator for constructing the diagonal matrix; λ1 and λ2 represent the weights of the prior information constraint term and the sparse constraint term of the initial model, respectively; μ is a constant to ensure that the matrix is ​​diagonally dominant. Specifically, in step 170, based on the Bayesian framework, while calculating the optimal impedance parameter inversion result, the uncertainty information of the inversion result is also calculated. The seismic inversion method based on the joint constraints of the physical model and prior information can predict the inversion result while giving the uncertainty of the inversion result, realizing a quantitative reliability evaluation of the inversion result.

[0036] ∑=C z -((LC z ) T / (LC z L T +μI))·LC z(13)

[0037] In the formula, C z L represents the model covariance matrix; μ represents the forward modeling operator; I is a constant that ensures the matrix is ​​diagonally dominant; and I is the N... x ×N x The unit matrix, N x It is the length of the observed data.

[0038] A computer device includes a memory and a processor, the memory storing computer-executable instructions and the processor executing the computer-executable instructions, which, when executed by the processor, implement the steps of the method described above.

[0039] A computer-readable storage medium storing computer-executable instructions that, when executed by a processor, implement the steps of the above-described method.

[0040] The present invention has the following advantages:

[0041] 1. This invention estimates impedance directly from seismic data, avoiding the cumulative errors caused by trace integrals and recursive inversion methods;

[0042] 2. This invention introduces a physical model to ensure that the prediction results meet the observed seismic data, and designs an initial model prior information constraint term to reduce the multiple solutions of the inversion results while improving the accuracy of the inversion results;

[0043] 3. This invention introduces sparse prior constraints and vertical constraints, which can effectively highlight the boundary features of the inversion results and ensure the vertical continuity of the inversion results;

[0044] 4. Based on the Bayesian framework, the inversion results are predicted and the uncertainty of the inversion results are given, which can realize the quantitative reliability evaluation of the inversion results. Attached Figure Description

[0045] Figure 1 is a flowchart of the seismic inversion method of the present invention, which is jointly constrained by the physical model and prior information.

[0046] Figure 2 shows the noiseless post-stack seismic data and initial impedance parameter model input in this invention;

[0047] Figure 3 shows the input noisy post-stack seismic data and initial impedance parameter model of this invention;

[0048] Figure 4 shows the impedance and its uncertainty information obtained by the seismic inversion method of the present invention under noise-free conditions, which uses physical models and prior information as joint constraints.

[0049] Figure 5 shows the impedance and its uncertainty information obtained by the seismic inversion method of the present invention, which uses physical model and prior information as joint constraints under noisy conditions. Detailed Implementation

[0050] The present invention will be further described below with reference to the accompanying drawings, but the scope of protection of the present invention is not limited to the following description. As shown in Figures 1 to 5, a seismic inversion method jointly constrained by a physical model and prior information includes the following steps:

[0051] Step 110: Assuming the seismic wavelet is known before inversion, it is necessary to extract the wavelet using statistical methods based on actual seismic gathers and well logging data. Actual seismic amplitudes are often relative values, and the amplitudes of the seismic data simulated using the post-stack seismic reflection coefficient equation differ somewhat from the actual amplitudes. Using well logging data as input, the reflection coefficient is calculated using the post-stack seismic reflection coefficient equation.

[0052] In the formula, m represents the logging impedance parameter data, and r represents the calculated reflection coefficient;

[0053] Then, the reflection coefficient is convolved with the extracted seismic wavelet to obtain the seismic gather, which is compared with the actual well-side seismic gather. The amplitude scaling factor is calculated and applied to the extracted seismic wavelet to achieve amplitude matching between the simulated record and the actual record.

[0054] In the formula, x represents the synthetic seismic record. Let w represent the convolution operator, and w represent the extracted wavelet;

[0055] Step 120: Extract impedance parameters and their mean values ​​from all well logging data within the work area, and statistically calculate their variance and vertical variation function matrix; obtain the required impedance parameters through well logging data analysis, calculate the autocorrelation coefficients of the impedance parameters, construct the variance matrix, and calculate the vertical variation function based on the impedance parameters to form a prior distribution function of impedance parameters conforming to the work area; the calculation formula for the vertical variation function is as follows:

[0056] Where υ represents the variogram value calculated based on well data, h represents the hysteresis distance, N represents the length of the logging data, and z(x i ) represents the impedance parameter value at a certain location i.

[0057] Step 130: Using seismic data stratigraphic interpretation data and well logging data, establish an initial impedance parameter model in the time domain; using seismic tectonic interpretation data, establish a geological model based on the sedimentary model, and interpolate and extrapolate the well logging data according to the tectonic model to obtain the initial impedance parameter model for each survey line; the impedance parameter model is established using spatial interpolation methods. First, the data of each stratigraphic level is interpolated using scattered interpolation to complete the geological stratigraphic model. Then, the impedance parameters are interpolated laterally according to the geological stratigraphic level to calculate the impedance parameter value at each underground point, thus completing the initial impedance parameter model.

[0058] Step 140: Starting from the accurate post-stack reflection coefficient equation, and based on the assumption of weak elastic difference at the interface, a linear approximate equation is derived:

[0059] In the formula, m represents the impedance parameter data, r represents the calculated reflection coefficient, and ln represents taking the logarithm of the data;

[0060] The post-stack seismic gather is simulated using a time-domain initial impedance parameter model and simplified equations. Taking the time-domain initial impedance parameter model as input, the reflection coefficient vector is directly calculated using simplified equations. The seismic wavelet is convolved with the reflection coefficients to obtain the post-stack seismic gather, and the difference between this and the actual seismic gather is used to obtain the inversion residual.

[0061] In the formula, K represents the wavelet matrix constructed based on w, D represents the differential operator, and m represents the impedance parameter data.

[0062] Step 150: Based on Bayes' principle, and by introducing sparse prior constraints and initial model prior information constraints, construct the inversion objective function in the sense of maximum a posteriori probability. Solve the inversion objective function using the generalized linear inversion method to obtain the solution expression for the impedance parameter. Constrain the inversion results using the initial model, and assume that the seismic data noise follows a Gaussian distribution, the initial model prior information constraints follow a Gaussian distribution, and the model follows a sparse distribution. Then, the inversion likelihood function and the prior probability distribution satisfy the Gaussian distribution and the joint distribution of Gaussian and sparse, respectively. According to Bayes' principle, synthesize the inversion likelihood function and the prior distribution function to obtain the posterior probability distribution function. Determine the inversion objective function based on the posterior probability, and differentiate the objective function with respect to the model parameters to obtain the iterative solution formula.

[0063] Let the impedance model parameters be z. T =(z1,z2,…,z) n ) T The observed earthquake data is x T =(x1,x2,…,x n) T According to Bayesian theory, given post-stack seismic data, the problem of inverting subsurface impedance parameters can be reduced to solving a posterior probability function:

[0064] Where P(x)=∫p(x|z)p(z)dz is the normalization factor, which can be regarded as a constant, P(x|z) is the likelihood function, and P(z) is the prior probability distribution. It is assumed that the likelihood function follows a Gaussian distribution:

[0065] In the above formula, Γ is the forward modeling operator, and C X It is the covariance matrix of the noise, N x Given the length of the observed data, based on Bayesian theory and combined with the prior model, after canceling the normalization constant, the posterior probability distribution of the inverted post-stack seismic data can be expressed as:

[0066] In the formula, Γ is the forward modeling operator, Φ(z) represents the sparse constraint term, and z0 represents the initial impedance parameter. Represents the smoothing operator, C X It is the covariance matrix of the noise, C z The model covariance matrix can be represented by the variance statistically derived from well logging data. α sum of variation function matrix υ Φ(z) is obtained through the Kronecker product. The specific expression is as follows:

[0067] Φ(z)=||Dz||0 (9)

[0068] In the formula, D represents the difference operator. If three-point smoothing is used, The expression can be written as:

[0069] Depending on the specific degree of smoothness, it can be modified according to equation (10). The expression,

[0070] The optimal solution can be obtained by finding the maximum value of equation (8), which is equivalent to finding the solution corresponding to the minimum value of the objective function in the following equation.

[0071] Applying the above objective function to the impedance parameter z Taking the derivative and setting it to zero, since the 0-norm is not differentiable, we use the iterative reweighted least squares algorithm to solve for the model parameters, which yields the iterative formula for updating the model parameters:

[0072] In the formula, L = 0.5·W·D represents the forward modeling operator; y = ln(z) represents the logarithm of the impedance parameter; W represents the wavelet matrix; X represents the seismic record; and I is N. x ×N x The unit matrix, N x λ is the length of the observed data; diag represents the operator for constructing the diagonal matrix; λ1 and λ2 represent the weights of the prior information constraint term and the sparse constraint term of the initial model, respectively; μ is a constant that ensures the matrix is ​​diagonally dominant.

[0073] Step 160: Based on the solution formulas for the model parameters and the inversion residuals, calculate the model parameters using the iterative reweighted least squares algorithm;

[0074] Step 170: Repeat steps 140, 150, and 160, and control the maximum number of iterations by using the inversion residuals to obtain the optimal impedance parameter inversion result; as shown in Figures 4 and 5, these are the impedance parameters obtained by the seismic inversion method based on the joint constraints of physical models and prior information in this embodiment of the invention. In the figures, the vertical axis represents time in seconds, and the horizontal axis represents impedance (unit: g / cm). 3 With a speed of ·km / s), this invention can predict impedance parameter information with high accuracy. As shown in Figure 5, the inversion result is obtained when random noise with a signal-to-noise ratio of 4 is added. The introduction of the prior model plays a key role in maintaining the stability of the inversion process and improving the accuracy of the inversion result. Based on the Bayesian framework, the uncertainty information of the inversion result can be calculated at the same time as the inversion result (Equation 13). The uncertainty is shown in the right of Figure 4 and the right background color of Figure 5. This invention can predict the inversion result and give the uncertainty of the inversion result at the same time, and can realize the quantitative reliability evaluation of the inversion result.

[0075] ∑=C z -((LC z ) T / (LC z L T +μI))·LC z (13)

[0076] In the formula, C z L represents the model covariance matrix; μ represents the forward modeling operator; I is a constant that ensures the matrix is ​​diagonally dominant; and I is the N... x ×N x The unit matrix, N x It is the length of the observed data.

[0077] A computer device includes a memory and a processor, the memory storing computer-executable instructions and the processor executing the computer-executable instructions, which, when executed by the processor, implement the steps of the method described above.

[0078] A computer-readable storage medium storing computer-executable instructions that, when executed by a processor, implement the steps of the above-described method.

[0079] The present invention can achieve stable and high-precision prediction of impedance parameters through the above method.

Claims

1. A seismic inversion method jointly constrained by a physical model and prior information, characterized in that: Includes the following steps: Step 110: Extract wavelets using actual seismic data, perform forward modeling of seismic gathers based on well logging data and post-stack seismic reflection coefficient equations, and determine the amplitude scaling factor by combining actual well-side seismic data; Step 120: Extract the impedance parameters and their mean values ​​based on all well logging data in the work area, and statistically calculate their variance and vertical variation function matrix; Step 130: Use seismic data layer interpretation data and well logging data to establish an initial impedance parameter model in the time domain; Step 140: Using the accurate post-stack reflection coefficient equation, based on the assumption of weak elasticity difference at the interface, a linear approximate equation is derived. Based on the initial impedance parameter model and simplified equation in the time domain, the post-stack seismic gather is simulated in forward modeling. The inversion residual is directly calculated from the forward modeling record and the actual record. Step 150: Based on Bayes' principle, introduce sparse prior constraints and initial model prior information constraints to construct the inversion objective function in the sense of maximum a posteriori probability. Solve the inversion objective function using the generalized linear inversion idea to obtain the solution expression for the impedance parameter. Step 160: Based on the solution formulas for the model parameters and the inversion residuals, calculate the model parameters using the iterative reweighted least squares algorithm; Step 170: Repeat steps 140, 150 and 160, and control the maximum number of iterations by the inversion residual, and output the inversion result of the optimal impedance parameters.

2. The seismic inversion method based on joint constraints of physical model and prior information according to claim 1, characterized in that: In step 170, based on the Bayesian framework, while calculating the optimal impedance parameter inversion result, the uncertainty information of the inversion result is also calculated to achieve a quantitative reliability evaluation of the inversion result: ∑=C z -((LC z ) T / (LC z L T +μI))·LC z (13) In the formula, C z L represents the model covariance matrix; μ represents the forward modeling operator; I is a constant that ensures the matrix is ​​diagonally dominant; and I is the N... x ×N x The unit matrix, N x It is the length of the observed data.

3. The seismic inversion method based on joint constraints of physical model and prior information according to claim 1, characterized in that: In step 110, the reflection coefficient is calculated using the post-stack seismic reflection coefficient equation with well logging data as input model: In the formula, m represents the logging impedance parameter data, and r represents the calculated reflection coefficient; The reflection coefficient is then convolved with the extracted seismic wavelet to obtain a seismic gather. This gather is compared with the actual well-side seismic gather, and the amplitude scaling factor is calculated and applied to the extracted seismic wavelet to achieve amplitude matching between the simulated and actual records. In the formula, x represents the synthetic seismic record. Let w represent the convolution operator, and w represent the extracted wavelet.

4. The seismic inversion method based on joint constraints of physical model and prior information according to claim 3, characterized in that: In step 120, the required impedance parameters are obtained through well logging data analysis, the autocorrelation coefficient of the impedance parameters is calculated, the variance matrix is ​​constructed, and the vertical variation function matrix is ​​calculated based on the impedance parameters to form the a priori distribution function of the impedance parameters that conforms to the work area. The formula for calculating the vertical variation function is as follows: Where υ represents the variogram value calculated based on well data, h represents the hysteresis distance, N represents the length of the logging data, and z(x i ) represents the impedance parameter value at a certain position i.

5. The seismic inversion method based on joint constraints of physical model and prior information according to claim 4, characterized in that: In step 130, seismic tectonic interpretation data are used to establish a geological model based on the sedimentary model. Well logging data are then interpolated and extrapolated according to the tectonic model to obtain the initial impedance parameter model for each logging line.

6. The seismic inversion method based on joint constraints of physical model and prior information according to claim 5, characterized in that: In step 130, the impedance parameter model is established using a spatial interpolation method. First, the data of each layer is interpolated using the scattered interpolation method to complete the geological layer modeling. Then, the impedance parameters are interpolated laterally according to the geological layers to calculate the impedance parameter values ​​at each underground point, thus completing the initial impedance parameter modeling.

7. The seismic inversion method based on joint constraints of physical model and prior information according to claim 6, characterized in that: In step 140, based on the assumption of weak elasticity difference, a linear approximation equation is derived: In the formula, m represents the impedance parameter data, r represents the calculated reflection coefficient, and ln represents taking the logarithm of the data; Then, using the initial impedance parameter model in the time domain as input, the reflection coefficient vector is directly calculated using simplified equations. The seismic wavelet is convolved with the reflection coefficient to obtain the post-stack seismic gather, and the difference between this and the actual seismic gather is used to obtain the inversion residual. In the formula, K represents the wavelet matrix constructed based on w, D represents the differential operator, and m represents the impedance parameter data.

8. The seismic inversion method based on joint constraints of physical model and prior information according to claim 7, characterized in that: In step 150, the initial model is used to constrain the inversion results. It is assumed that the seismic data noise follows a Gaussian distribution, the prior information constraint terms of the initial model follow a Gaussian distribution, and the model itself follows a sparse distribution. Therefore, the inversion likelihood function and the prior probability distribution satisfy the Gaussian distribution and the joint distribution of Gaussian and sparse distribution, respectively. According to Bayes' principle, the posterior probability distribution function is obtained by combining the inversion likelihood function and the prior distribution function. The objective function for inversion is determined based on the posterior probability distribution function. The iterative solution formula is obtained by differentiating the objective function with respect to the model parameters. In the formula, L = 0.5·W·D represents the forward modeling operator; D represents the difference operator; W represents the wavelet matrix; and X represents the seismic record. Let y = ln(z) represent the smoothing operator, and let I represent the logarithm of the impedance parameter. x ×N x The unit matrix, N x λ is the length of the observed data; diag represents the operator for constructing the diagonal matrix; λ1 and λ2 represent the weights of the prior information constraint term and the sparse constraint term of the initial model, respectively; μ is a constant that ensures the matrix is ​​diagonally dominant.

9. A computer device, characterized in that: The method includes a memory and a processor, the memory being used to store computer-executable instructions, and the processor being used to execute the computer-executable instructions, which, when executed by the processor, implement the steps of the method according to any one of claims 1 to 8.

10. A computer-readable storage medium storing computer-executable instructions, characterized in that: When the computer-executable instructions are executed by a processor, they implement the steps of the method according to any one of claims 1 to 8.

Images