A method and device for prestack density stable inversion

Abstract

The present invention discloses a method and apparatus for stable prestack density inversion. The method comprises: preprocessing seismic data to convert the offset data volume of the seismic data into an angle gather data volume; calculating the angle gather data volume according to a first elastic impedance formula to obtain a first elastic impedance; introducing a low-frequency model into the first elastic impedance to obtain a second elastic impedance; and extracting elastic parameters from the second elastic impedance based on weighting coefficients. Compared with existing technologies, the present invention can produce stable, high-precision density inversion results.

Description

Technical Field

[0001] The present invention relates to the technical field of seismic data interpretation, and in particular to a pre-stack density stable inversion method and device. Background Art

[0002] Since Connolly first proposed the concept of elastic impedance (EI) in 1999, the theory has developed rapidly and become an important tool in seismic inversion. The elastic impedance equation is derived based on the reflection coefficient approximation and is a generalization of the post-stack wave impedance, although it is only a calculation quantity and has no clear physical meaning. The equation contains information on P-wave velocity, density, and incident angle, and can extract a variety of elastic parameters to reflect reservoir lithology and fluid information. VerWest et al. (2000) derived the elastic impedance equation in the ray domain and extended the elastic impedance theory to the ray domain. With the development of multi-wave and multi-component seismic exploration technology, the concept of elastic impedance has been extended to converted shear waves, and Duffaut et al. (2000) derived the converted shear wave elastic impedance (SEI) equation.

[0003] However, Connolly's EI equation is affected by the angle of incidence, resulting in significant dimensional variation in elastic impedance. Whitecombe (2002) therefore proposed the concept of extended elastic impedance (EEI) to standardize the equation. Martins (2002) generalized the elastic impedance equation to weakly anisotropic media and applied it to real-world data. Lu et al. (2004) applied the elastic impedance inversion method to extract target parameters from real-world data.

[0004] Domestic scholars have also conducted extensive research on elastic impedance inversion. Ni Yi (2003) proposed the normalized elastic impedance method, and Ma Jinfeng et al. (2003) proposed the concepts of generalized elastic impedance (GEI) and Zoeppritz elastic impedance (ZEI), which were successfully applied to hydrocarbon detection. Gan Lideng et al. (2003) conducted a comparative analysis of the accuracy of reflection coefficients. Wang Baoli (2007) proposed an elastic wave impedance inversion method based on the Gray approximation equation, which can directly obtain the Lame constant. Yin Xingyao et al. (2010) proposed an elastic wave impedance equation expressed in terms of fluid factor, Lame constant and density, and the corresponding inversion method based on the approximate equation of porous fluid-saturated elastic media. Zong Zhaoyun (2011) used the sparse pulse inversion method in the Bayesian framework to realize pre-stack elastic wave impedance inversion based on the Gray approximation. Liu Xiaojing et al. (2016) applied basis pursuit and reflection coefficient odd-even decomposition theory to prestack inversion, realized the inversion of Gassmann fluid terms and shear modulus, and analyzed the sensitivity of different elastic parameters to reservoir fluids. Zhang Fenglin (2017) introduced the Bayesian principle to construct regularization terms, and improved the inversion stability through sparse constraints and low-frequency soft constraints. Wang (2019) et al. used non-convex regularization to extract stable elastic impedance from seismic data, and then extracted longitudinal and shear wave velocities and densities from the elastic impedance under the Bayesian inversion framework. Hao Yaju (2020) projected the seismic signal onto the odd and even sub-wave library to obtain the relative reflection coefficient of each channel, and then performed channel integration to obtain the relative elastic wave impedance, thereby improving the thin layer resolution and inversion accuracy.

[0005] The elastic impedance equation enables the inversion of far-offset angular stacking data using acoustic impedance inversion techniques, but its dimension varies with the angle of incidence, complicating the comprehensive analysis of AI and EI. In the absence of noise, the elastic impedance equation based on Connolly's method can accurately invert elastic parameters. However, as noise increases, the errors in the resulting P- and S-wave velocities also increase. Mallick et al. have shown that in the presence of 2% random noise, this elastic impedance equation is no longer sufficient to invert reasonable elastic parameters. Consequently, many researchers have derived elastic impedance equations for various elastic parameter data volumes. In 2008, Wang Baoli proposed and standardized an elastic impedance equation based on the Fatti approximation, enabling direct inversion of the impedance data volume, avoiding cumulative computational errors and resulting in more accurate inversion results. Long-term practice has demonstrated that direct inversion of elastic parameters is difficult and subject to significant instability, while inverting elastic impedance is a relatively stable process. Consequently, parameter inversion methods based on elastic impedance have long been used in seismic inversion.

[0006] Density, a key indicator of reservoir properties, is closely related to parameters such as porosity, fluid type, saturation, and mineral content, and can intuitively reflect the changing patterns of reservoirs and fluids. Improving the stability of density inversion is an important research direction in the field of prestack inversion. In AVO (Amplitude Variation with Offset inversion) inversion, the commonly used reflection coefficient approximation equation generally assumes that the seismic wavelet does not change with offset. However, in actual processing, factors such as dynamic correction stretching can cause wavelet variations, complicating AVO attribute extraction. In contrast, elastic impedance inversion uses different wavelet inversions for stacked gathers at different angles, taking into account the variation of the wavelet with offset, thereby providing more reliable elastic parameters. However, the existing elastic impedance equation deviates from the exact value when calculating the reflection coefficient, and this deviation increases with increasing angle, affecting the accuracy of the inversion results.

[0007] In addition, the traditional pre-stack parameter inversion method based on elastic impedance faces two major problems: first, the Aki-Richards reflection coefficient approximation formula, after linearization into the expression of elastic impedance, is insufficiently accurate, which is inconsistent with the development trend of high-precision inversion; second, in the process of using elastic impedance to extract elastic parameters, the strong nonlinearity between elastic parameters and impedance increases the difficulty of parameter estimation, making it difficult to achieve the goal of stable inversion. Summary of the Invention

[0008] The technical problem to be solved by the present invention is to address at least one defect of the related technology mentioned in the above background technology: in the existing technology, it is difficult to obtain stable and high-precision density inversion results, the conventional reflection coefficient approximation equation has a large deviation from the exact value, and the deviation gradually increases with the increase of the incident angle, and a pre-stack density stable inversion method and device are provided.

[0009] The technical solution adopted by the present invention to solve the technical problem is to construct a pre-stack density stable inversion method, which includes the following steps:

[0010] Preprocessing the seismic data, converting the offset data volume of the seismic data into the angle gather data volume;

[0011] Calculate the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance;

[0012] Introducing the low-frequency model into the first elastic impedance to obtain the second elastic impedance;

[0013] An elastic parameter is extracted from the second elastic impedance according to a weighting coefficient.

[0014] In some embodiments, preprocessing the seismic data to convert the offset data volume of the seismic data into an angle gather data volume, followed by: obtaining well logging data through well-seismic calibration, and extracting angle wavelets from the angle gather data volume of the seismic data and the well logging data;

[0015] Calculating the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance includes: using the first elastic impedance formula to calculate the angle gather data volume of seismic data and logging data corresponding to the angle wavelet to obtain the first elastic impedance.

[0016] In some embodiments, the angle gather data volume includes angle gathers at three different angles.

[0017] In some embodiments, the low-frequency model is introduced into the first elastic impedance to obtain the second elastic impedance, and the method further includes:

[0018] The low-frequency components in the logging data are weighted and calculated to establish a low-frequency model.

[0019] In some embodiments, extracting elastic parameters from the second elastic impedance according to the weighting coefficients further includes: solving a coefficient matrix using a coefficient estimation method with spatial variation characteristics to obtain the weighting coefficients.

[0020] In some embodiments, a coefficient estimation method with spatial variation characteristics is used to solve a coefficient matrix to obtain weighted coefficients, including: performing statistical analysis on the second elastic impedance, establishing a coefficient matrix for the second elastic impedance with bimodal or multimodal characteristics, and solving the weighted coefficients using singular value decomposition.

[0021] In some embodiments, the first elastic impedance and the second elastic impedance include longitudinal wave velocity, shear wave velocity, and density, respectively.

[0022] In some embodiments, the first elastic impedance formula is an elastic impedance formula based on power series approximation.

[0023] In some embodiments, preprocessing the seismic data includes fitting missing data based on existing data.

[0024] The present invention also constructs a pre-stack density stable inversion device, which includes:

[0025] one or more processors;

[0026] A storage device is used to store one or more programs, which, when executed by the one or more processors, enable the one or more processors to implement the pre-stack density stable inversion method as described in any of the above embodiments.

[0027] By implementing the present invention, the following beneficial effects are achieved:

[0028] The present invention uses a pre-stack density stability inversion method to preprocess seismic data, convert the offset data volume into an angle gather data volume, calculate the angle gather data volume according to the first elastic impedance formula to obtain a first elastic impedance, introduce a low-frequency model into the first elastic impedance to obtain a second elastic impedance, and extract elastic parameters from the second elastic impedance based on weighting coefficients. Compared with existing technologies, the present invention can achieve stable and high-precision density inversion results. BRIEF DESCRIPTION OF THE DRAWINGS

[0029] The present invention will be further described below with reference to the accompanying drawings and embodiments, in which:

[0030] Figure 1 A flow chart of an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0031] Figure 2 A comparison diagram of the reflection coefficients of the first type AVO model in one embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0032] Figure 3 A comparison diagram of the reflection coefficient errors of the first type AVO model in one embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0033] Figure 4 A comparison diagram of reflection coefficients of the second type AVO model in an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0034] Figure 5 A comparison diagram of reflection coefficient errors of the second type AVO model in an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0035] Figure 6 A comparison diagram of reflection coefficients of the third type AVO model in an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0036] Figure 7 A comparison diagram of reflection coefficient errors of the third type AVO model in an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0037] Figure 8 Graphs showing elastic impedance inversion results at three angles when the signal-to-noise ratio is 0 in one embodiment of the pre-stack density stability inversion method of the present invention are shown;

[0038] Figure 9 Graphs showing elastic impedance inversion results at three angles when the signal-to-noise ratio is 10 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0039] Figure 10Graphs showing elastic impedance inversion results at three angles when the signal-to-noise ratio is 5 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0040] Figure 11 FIG1 shows the inversion result of elastic parameters of the present invention when the signal-to-noise ratio is 0 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0041] Figure 12 The figure shows the inversion result of the elastic parameters of the Connolly elastic impedance equation when the signal-to-noise ratio is 0 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0042] Figure 13 FIG1 shows the inversion result of elastic parameters of the present invention when the signal-to-noise ratio is 10 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0043] Figure 14 The figure shows the inversion result of the elastic parameters of the Connolly elastic impedance equation when the signal-to-noise ratio is 10 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0044] Figure 15 FIG1 shows the inversion result of elastic parameters of the present invention when the signal-to-noise ratio is 5 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0045] Figure 16 The figure shows the inversion result of the elastic parameters of the Connolly elastic impedance equation when the signal-to-noise ratio is 5 in one embodiment of the pre-stack density stable inversion method of the present invention;

[0046] Figure 17 The figure shows the inversion result of the elastic impedance of the 10° well bypass channel in one embodiment of the pre-stack density stability inversion method of the present invention;

[0047] Figure 18 The figure shows the inversion result of the elastic impedance of the 20° well bypass channel in one embodiment of the pre-stack density stability inversion method of the present invention;

[0048] Figure 19 The figure shows the inversion result of the elastic impedance of the 30° well bypass channel in one embodiment of the pre-stack density stability inversion method of the present invention;

[0049] Figure 20 The figure shows the inversion results of the actual well-side channel P-wave velocity, S-wave velocity and density in one embodiment of the pre-stack density stable inversion method of the present invention;

[0050] Figure 21 A 10° angle stack gather diagram of an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0051] Figure 22A 20° angle stack gather diagram of an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0052] Figure 23 A 30° angle stack gather diagram of an embodiment of the pre-stack density stable inversion method of the present invention is shown;

[0053] Figure 24 FIG1 shows the inversion result of the longitudinal wave velocity in one embodiment of the pre-stack density stable inversion method of the present invention;

[0054] Figure 25 FIG1 shows the inversion result of shear wave velocity in one embodiment of the prestack density stable inversion method of the present invention;

[0055] Figure 26 The figure shows the density inversion result of one embodiment of the pre-stack density stable inversion method of the present invention. DETAILED DESCRIPTION

[0056] In order to have a clearer understanding of the technical features, purposes and effects of the present invention, specific embodiments of the present invention are now described in detail with reference to the accompanying drawings.

[0057] It should be noted that the flowcharts shown in the accompanying drawings are for illustrative purposes only and do not necessarily include all content and operations / steps, nor must they be executed in the order described. For example, some operations / steps may be decomposed, while others may be combined or partially combined. Therefore, the actual execution order may vary depending on the actual situation.

[0058] The block diagrams shown in the accompanying drawings are merely functional entities and do not necessarily correspond to physically separate entities. That is, these functional entities may be implemented in software, in one or more hardware modules or integrated circuits, or in different networks and / or processor devices and / or microcontroller devices.

[0059] Some embodiments of the present invention disclose a method for prestack density stability inversion, such as Figure 1 As shown, the method includes the following steps:

[0060] Preprocessing the seismic data, converting the offset data volume of the seismic data into the angle gather data volume;

[0061] Calculate the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance;

[0062] Introducing the low-frequency model into the first elastic impedance to obtain the second elastic impedance;

[0063] An elastic parameter is extracted from the second elastic impedance according to a weighting coefficient.

[0064] In some embodiments, preprocessing the seismic data includes fitting missing data based on existing data.

[0065] In some embodiments, the offset data volume refers to the data set involved in the process of "deflecting" the recorded data back to its original real spatial position in seismic data processing, and the angle gather data volume is a data set that classifies and organizes seismic data according to specific incident angles or reflection angles.

[0066] In some embodiments, the first elastic impedance is the well bypass elastic impedance (EI) or the relative elastic impedance, and the second elastic impedance is the absolute elastic impedance.

[0067] In some embodiments, the first elastic impedance formula is: an elastic impedance formula based on power series approximation, and the power series approximation is specifically as follows: first, Taylor expansion is used to express the reflection coefficient formula as a power series form with respect to the ray parameter p, and then approximation is performed by retaining only the quadratic term of the ray parameter p. The four-term approximation formula obtained by this method has shown higher accuracy than the traditional linear approximation formula in the test of four types of AVO models, and in terms of the sensitivity matrix, the condition number with respect to the three elastic parameters has also been significantly improved. The following is the derivation process of the elastic impedance formula based on power series approximation, and the specific power series approximation reflection coefficient expression of the three elastic parameters is shown as follows:

[0068]

[0069] Where R(θ) is the formation reflection coefficient, α, β, and ρ are the P-wave velocity, S-wave velocity, and density, respectively, and Δα, Δβ, and Δρ are the differences in P-wave velocity, S-wave velocity, and density between the lower and upper layers on both sides of the interface, respectively. are the average values ​​of the longitudinal wave velocity, shear wave velocity and density on both sides of the interface, θ is the average value of the incident angle and the reflection angle, and γ is the ratio of the average values ​​of the longitudinal and shear wave velocities of the upper and lower layers.

[0070] For small to moderate changes in wave impedance, the reflection coefficient expressed as the logarithm of the wave impedance is accurate:

[0071]

[0072] Where EI is the elastic impedance, and ΔEI represents the change in elastic impedance. Based on the relationship in equation (2), we can rewrite equation (1) to obtain:

[0073]

[0074] Performing Taylor expansion on the fourth term on the right side of the equation and ignoring higher-order terms, we can obtain:

[0075]

[0076] Integrate both sides of the equation, set the constant to 0, and let is a constant, we can get:

[0077]

[0078] Since the value of formula (5) varies greatly with the incident angle, it is not conducive to comparing the elastic wave impedance values ​​at different angles. Therefore, formula (5) is further standardized, that is,

[0079]

[0080] Where α, β, and ρ are the longitudinal wave velocity, shear wave velocity, and density, respectively. are the average values ​​of the corresponding parameters in the logging curve, is the ratio of the average transverse and longitudinal wave velocities of the upper and lower layers, and formula (6) is the first elastic impedance formula in the present invention. By using the reflection coefficient approximation equation of the power series approximation to make the elastic impedance body, the problem of insufficient accuracy of the traditional elastic impedance equation can be well solved. The synthetic data test shows that the reflection coefficient calculated by the method of the present invention is compared with the exact value, and it can be found that this elastic impedance equation can effectively improve the accuracy of the reflection coefficient. Therefore, by using this formula to extract elastic parameters from the elastic impedance body, the elastic parameters obtained by inversion can be made more accurate and stable, and better describe the reservoir information. As shown in Table 1 and Figure 2-Figure 7 As shown, the present invention calculates the reflection coefficient (El-Quadratic) and compares it with the precise Zoeppritz reflection coefficient (Zoe) and the reflection coefficient calculated by the classic Connolly elastic impedance equation (Connally), demonstrating the significant advantages of the present invention in terms of calculation accuracy. In particular, in the calculation of the AVO model, the method of the present invention performs particularly well in the first and second types of AVO models. Its calculation results are not only closer to the exact value, but also have a smaller deviation from the exact value when the incident angle increases, which further demonstrates the effectiveness of the present invention in improving calculation accuracy. In the third type of AVO model, the accuracy of the reflection coefficient of the present invention is not much different from that of the Connolly elastic impedance equation.

[0081] Table 1 Formation model parameters

[0082]

[0083] In some embodiments, preprocessing the seismic data to convert the offset data volume of the seismic data into an angle gather data volume, followed by: obtaining well logging data through well-seismic calibration, and extracting angle wavelets from the angle gather data volume of the seismic data and the well logging data;

[0084] Calculating the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance includes: using the first elastic impedance formula to calculate the angle gather data volume of seismic data and logging data corresponding to the angle wavelet to obtain the first elastic impedance.

[0085] In some embodiments, the angle gather data body includes angle gathers at three different angles, specifically small, medium and large angle gathers, for example: 10° angle stacking gather, 20° angle stacking gather and 30° angle stacking gather.

[0086] In some embodiments, the low-frequency model is introduced into the first elastic impedance to obtain the second elastic impedance, and the method further includes:

[0087] The low-frequency components in the logging data are weighted to establish a low-frequency model. Specifically, the logging data is parameter-fitted, and the low-frequency components between wells are weighted according to distance. The low-frequency model for the entire inversion block is then established using geological interpretation data. Parameter fitting involves fitting missing parameters based on known parameters. For example, if shear wave velocity is missing, a shear wave velocity curve can be obtained by calculating and fitting the compressional wave velocity and density.

[0088] In some embodiments, extracting elastic parameters from the second elastic impedance according to the weighting coefficients further includes: solving a coefficient matrix using a coefficient estimation method with spatial variation characteristics to obtain the weighting coefficients.

[0089] In some embodiments, a coefficient estimation method with spatial variation characteristics is used to solve the coefficient matrix to obtain the weighted coefficient, including: performing statistical analysis on the second elastic impedance, establishing a coefficient matrix for the second elastic impedance with bimodal or multimodal characteristics, and solving the weighted coefficient using singular value decomposition. The singular value decomposition is specifically as follows: let the matrix D n×m The singular value decomposition of is:

[0090] D T =VΛU T (7)

[0091] Where r = rank(D). When r is much smaller than the number of rows and columns of the matrix D, the first r singular values ​​can be used to approximate the matrix, which is called partial singular value decomposition. That is, Λ is an r×r order diagonal matrix, with all elements except the main diagonal being 0, and each element σ on the main diagonal being i are all called eigenvalues, and σ1≥σ2≥L≥σ i ≥0;D T To represent the transpose of matrix D, U is an n×r matrix, U T To represent the transpose of matrix U, V T is an r×m-order matrix, V TTo represent the transpose of matrix V, when r is closer to m, the result of multiplication on the right side is closer to the left side.

[0092] D T D=(VΛU T )(U T ΛV)=VΛ(U T U)ΛV T =VΛ 2 V T (8)

[0093] (D T D) -1 =(VΛ 2 V T ) -1 =VΛ -2 V T (9)

[0094] Among them, Λ 2 represents the square of each singular value in Λ, Λ -2 Represents Λ 2 The inverse matrix of , that is, the inverse of the square of each singular value.

[0095] Then the least square generalized inverse can be transformed into the following form:

[0096] (D T D) -1 D T =(VΛ 2 V T ) -1 (VΛU T )=VΛ -1 U T (10)

[0097] So we can get the least square solution b of the inversion as:

[0098] b=(D T D) -1 D T m=VΛ -1 U T m(11)

[0099] Substituting the coefficient matrix D established by the second elastic impedance into formula (11) can obtain the weighted coefficient.

[0100] In some embodiments, the first elastic impedance and the second elastic impedance include longitudinal wave velocity, shear wave velocity, and density, respectively.

[0101] Some embodiments of the present invention further disclose a pre-stack density stabilization inversion device, comprising:

[0102] one or more processors;

[0103] The storage device is used to store one or more programs. When the one or more programs are executed by one or more processors, the one or more processors implement the pre-stack density stable inversion method as described in any of the above embodiments, which will not be repeated here.

[0104] The following are the effects of the embodiments of the present invention:

[0105] like Figures 8-10 As shown in the figure, the present invention inverts synthetic seismic data with signal-to-noise ratios (SNRs) of 0, 10, and 5, yielding elastic impedance values ​​at three angles. The blue solid line represents the elastic impedance calculated from the well bypass curve, the orange solid line represents the elastic impedance inverted from the synthetic seismic data, and the black curve represents the initial value. The results show that elastic impedance values ​​at small, medium, and large angles can be accurately inverted.

[0106] like Figures 11-16 As shown in the figure, the reflection coefficients calculated by the method of the present invention are compared with the accurate Zoeppritz reflection coefficient and the reflection coefficients calculated by the classical Connolly elastic impedance equation. In the synthetic data test, we observed that the blue curve represents the true value of each elastic parameter, the orange curve represents the inverted value of each elastic parameter, and the black curve represents the initial value of each elastic parameter. Under different signal-to-noise ratio (SNR) conditions, including no noise (SNR=0), SNR of 10 and SNR of 5, the inversion results of the elastic impedance equation of the present invention and the conventional Connolly elastic impedance equation are compared. The results show that under noise-free conditions, the elastic parameter inversion results of both elastic impedance equations are very accurate. When the SNR is 10, the inversion results of shear wave velocity and density show a certain degree of deviation. However, when the SNR is reduced to 5, the elastic impedance equation proposed by the present invention shows better results in the inversion of shear wave velocity and density than the conventional elastic impedance equation, which proves its advantage in stability.

[0107] like Figure 17-Figure 19 As shown, to verify the effectiveness of the present invention in practical applications, the accuracy of the wellside channel application was first tested. Using the wellside channel logging data, the method of the present invention was applied to calculate the initial model and true value of the elastic impedance. Next, elastic impedance inversion was performed on three data sets at different angles. In the results, the blue solid line represents the true value of the elastic impedance, the red solid line represents the result obtained through inversion, and the black solid line represents the initial model of the inversion process. This intuitive comparison can be used to evaluate the accuracy and reliability of the present method in elastic impedance inversion.

[0108] like Figure 20As shown in the figure, the elastic impedance values ​​of the small, medium and large angle gathers obtained by inversion are used to successfully invert the elastic parameter values ​​of the wellside channel. In the result display, the black solid line represents the true value of each elastic parameter, the red solid line represents the parameter value obtained by inversion, and the blue dotted line represents the initial value of the inversion. Figure 20 As shown in the figure, the degree of match between the actual data of the well bypass and the inversion results is high, which verifies the feasibility of the method of the present invention. Among the various parameters, the inversion accuracy of the longitudinal wave velocity is the highest, while the shear wave velocity and density can basically invert the trend of the entire curve. Specifically, the correlation coefficients between the inversion results of the longitudinal wave velocity, shear wave velocity and density and the true values ​​are 0.7239, 0.1775 and 0.3209 respectively, which shows that the inversion effect of density is more ideal than that of the shear wave velocity. The coefficient matrix obtained based on the results of the well bypass can be further combined with the elastic impedance body obtained by inversion to obtain the inversion result of the entire profile.

[0109] like Figure 21-23 As shown, in order to verify the applicability of the algorithm in the inversion of actual seismic data, the present invention was applied to the inversion process of actual seismic data. First, based on the seismic data collected in the field, angular channel data containing 400 sampling points and 500 channels were extracted. These angular channel sets are divided into three angle ranges of small, medium and large for superposition, thereby obtaining the sub-waves of each partially stacked data body. Then, the well data model is established using logging data and layer interpretation, and the elastic impedance values ​​of different angle channel sets are obtained by conventional elastic impedance inversion methods. Then, combining the relationship between the matrix coefficients obtained from the wellside channel and the angle channel sets, the longitudinal wave velocity, shear wave velocity and density of each channel are inverted. During the inversion process, the parameters are continuously adjusted to make the inversion results closer to the true value. The final inversion effect is shown in the figure. Figure 21-23 As shown in Figure 2, the P-wave velocity profile shows good correspondence with the well logging data, while the density profile inverted by the method of the present invention not only has clear lateral continuity but also has a good correspondence with the well logging data. This demonstrates the high applicability and accuracy of the present invention in the inversion of actual seismic data.

[0110] like Figure 24-26 As shown, the inversion results of the entire profile show that the inversion results of P-wave velocity, S-wave velocity, and density have significant continuity. Among these parameters, the inversion effect of P-wave velocity is the most outstanding, and the match with the well data is the highest, which shows the efficiency and accuracy of the inversion method in P-wave velocity. At the same time, the density data obtained by inversion also shows good continuity and a high degree of match with the well data, which further confirms the effectiveness and high stability of this method in inverting density. These results show that the inversion method of the present invention can provide reliable and consistent underground medium parameter information in seismic data interpretation.

[0111] The implementation of the present invention has the following beneficial effects: by comparing the calculated reflection coefficient with the precise value, it is found that the elastic impedance equation of the present invention can significantly improve the accuracy of the reflection coefficient calculation. In the comparative analysis of the reflection coefficient, the calculation results of the Zoeppritz equation, the Connolly elastic impedance equation and the present invention in four types of AVO models are presented. In the first and second types of AVO models, the method of the present invention is closer to the precise value than the Connolly equation, especially when the incident angle increases, the deviation from the precise value is smaller, which highlights the superiority of the present invention in improving accuracy. In the third type of AVO model, the accuracy of the reflection coefficients obtained by different methods is not much different.

[0112] Furthermore, the effectiveness and applicability of the proposed method were verified by the elastic parameter inversion results of synthetic data. Compared with the results of conventional inversion methods, the proposed method showed higher accuracy, which further demonstrated its potential in practical applications.

[0113] Secondly, elastic impedance inversion generally offers greater stability and noise immunity than traditional AVO inversion. This means that the elastic parameters extracted from elastic impedance are more stable and accurate, which is crucial for accurately characterizing reservoir information. Therefore, this invention not only improves the calculation accuracy of reflection coefficients but also facilitates a more accurate description and understanding of reservoir properties by providing more stable elastic parameters.

[0114] It is understandable that the above embodiments only express some of the implementation methods of the present invention, and their descriptions are relatively specific and detailed, but they cannot be understood as limiting the scope of the patent of the present invention. It should be pointed out that for ordinary technicians in this field, without departing from the concept of the present invention, the above embodiments or technical features can be freely combined, and several deformations and improvements can be made, which all fall within the scope of protection of the present invention, that is, the embodiments described in "some embodiments" can be freely combined with any of the above and below embodiments. Therefore, all equivalent transformations and modifications made to the scope of the claims of the present invention should fall within the scope of the claims of the present invention.

Claims

1. A prestack density stability inversion method, characterized in that: The method comprises the following steps: Preprocessing the seismic data, converting the offset data volume of the seismic data into the angle gather data volume; Calculate the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance; Introducing the low-frequency model into the first elastic impedance to obtain the second elastic impedance; The coefficient estimation method with spatial variation characteristics is used to solve the coefficient matrix to obtain the weighted coefficients; extracting elastic parameters from the second elastic impedance according to a weighting coefficient; The method of using a coefficient estimation method with spatial variation characteristics to solve a coefficient matrix to obtain weighted coefficients includes: performing statistical analysis on the second elastic impedance, establishing a coefficient matrix for the second elastic impedance with multi-peak characteristics, and solving the weighted coefficients using singular value decomposition; The first elastic impedance formula is an elastic impedance formula based on power series approximation, and the specific expression is: ； In the formula EI is the first elastic impedance, are the longitudinal wave velocity, shear wave velocity and density, are the average values ​​of the corresponding parameters in the logging curve, is the ratio of the average transverse and longitudinal wave velocities of the upper and lower layers, is the average of the angles of incidence and reflection.

2. The prestack density stable inversion method according to claim 1, characterized in that: Preprocessing the seismic data to convert the offset data volume of the seismic data into an angle gather data volume, followed by: obtaining well logging data through well-seismic calibration, and extracting angle wavelets from the angle gather data volume of the seismic data and the well logging data; Calculating the angle gather data volume according to the first elastic impedance formula to obtain the first elastic impedance includes: using the first elastic impedance formula to calculate the angle gather data volume of seismic data and logging data corresponding to the angle wavelet to obtain the first elastic impedance.

3. The prestack density stable inversion method according to claim 1, characterized in that: The angle gather data volume includes three angle gathers at different angles.

4. The prestack density stable inversion method according to claim 2, characterized in that: The low-frequency model is introduced into the first elastic impedance to obtain the second elastic impedance, which also includes: The low-frequency components in the logging data are weighted and calculated to establish a low-frequency model.

5. The prestack density stable inversion method according to claim 1, characterized in that: The first elastic impedance and the second elastic impedance include longitudinal wave velocity, shear wave velocity and density, respectively.

6. The prestack density stable inversion method according to claim 1, characterized in that: Preprocess the seismic data, including fitting missing data based on existing data.

7. A pre-stack density stable inversion device, characterized in that: include: one or more processors; A storage device for storing one or more programs, which, when executed by the one or more processors, enables the one or more processors to implement the prestack density stable inversion method according to any one of claims 1 to 6.

Images