CN115857004A - High-resolution multi-wave joint prestack inversion method for VTI medium of shale reservoir - Google Patents

Abstract

The invention discloses a high-resolution multi-wave joint prestack inversion method for a VTI medium of a shale reservoir. In the past, VTI medium inversion focuses more on the precision of an inversion result rather than the resolution, a Gaussian constraint is mostly adopted to establish a target function, the resolution of the obtained result is insufficient, and the requirement for a high-resolution result in shale reservoir exploration cannot be met. The method introduces the basis pursuit decomposition into the VTI medium anisotropic inversion, and restrains the target function by adopting the sparse norm so as to improve the vertical resolution of the inversion result to a great extent. And considering that the anisotropy multi-parameter inversion has strong unsuitability, the multi-wave joint inversion combining PP data and PS converted wave data is adopted to improve the accuracy of the result. Compared with the conventional VTI inversion method, the method has higher vertical resolution and good noise immunity.

Description

High-resolution multi-wave joint prestack inversion method for VTI medium in shale reservoirs

Technical Field

The present application relates to the field of seismic forward and inversion of anisotropic media such as shale reservoirs, and in particular to a high-resolution multi-wave joint prestack inversion method for VTI media in shale reservoirs.

Background Art

Prestack seismic inversion can make full use of the AVO/AVA (amplitude variation with offset/incident angle) information in seismic data and is currently one of the most effective seismic techniques for reservoir characterization and fluid identification. It helps to depict reservoir rock boundaries, development ranges, and fluid identification. In the field of high-resolution seismic inversion, many methods have been proposed. A series of earlier proposed random inversion methods used logging and geostatistical information to achieve high vertical resolution inversion prediction, but their strong uncertainty has not been resolved. Subsequent scholars have introduced sparse constraints or edge protection regularization in seismic inversion to improve the resolution of the results ^[1] , in order to improve the problems of smoothing and edge blurring caused by conventional Gaussian constraints. Zhang and Castagna realized post-stack basis pursuit inversion that can obtain high-resolution sparse results by performing odd-even pair decomposition on the reflection coefficient and using the _L1 sparse norm to constrain ^[2] . Based on this, many scholars have carried out related research to expand basis pursuit inversion from post-stack to pre-stack inversion, and have achieved good results in practical applications. Although the above research on high-resolution seismic inversion is relatively mature, it is all carried out on isotropic media, and the high-resolution inversion work on anisotropic media is relatively weak.

At present, unconventional oil and gas reservoirs such as shale oil and gas have become the focus and difficulty of international exploration. This type of rock exhibits strong VTI anisotropy characteristics (transverse isotropy with vertical axis of symmetry, VTI), and the previous inversion technology based on the isotropy assumption is no longer applicable. Research on anisotropic inversion methods for VTI media is of great significance for improving the exploration accuracy of shale reservoirs. The forward modeling research on VTI media is relatively complete ^[3-6] , but there are not many studies on prestack inversion methods for VTI media, mainly using the reflection coefficient approximation as the forward operator ^[5,6] . Hou Dongjia et al. (2014) proposed to use the Bayesian method to establish the objective function and realized the simultaneous inversion of five parameters of VTI media based on the reflection coefficient approximation ^[7] . Zhang et al. (2019) proposed a step-by-step inversion method for VTI media, which reduced the five-parameter inversion process to three parameters to improve the stability of the inversion, and then obtained the final five-parameter result through indirect conversion. However, these studies are all aimed at the accuracy of anisotropic inversion, not the vertical resolution. The anisotropic high-resolution inversion method for VTI media needs to be studied.

[1]Guo Q, Zhang H, Cao H, et al., Hybrid seismic inversion based on multi-order anisotropic markov random field[J]. IEEE Transactions onGeoscience and Remote Sensing, 2019, vol.58, no.1, pp. 1-14.

[2]Zhang R,Sen M,and Srinivasan S,Prestack basis pursuit seismicinversion[J]Geophysics,2013,vol.78,no.1,pp.R1-R11.

[3]Hron F,Daley P F.Reflection and transmission coefficients for transversely isotropic media[J].Bulletin of the Seismological Society ofAmerica,1977,67(3):661-675.

[4]Thomsen L A.Weak anisotropic reflections,Offset-DependentReflectivity: Theory and Practice of AVO Analysis[M].Library.seg.org,1993:103–114

A.Variation of P-wave reflectivity with offset and azimuth inanisotropic media[J].Geophysics,1998,63(3):935-947.[5]

A.Variation of P-wave reflectivity with offset and azimuth inanisotropic media[J].Geophysics,1998,63(3):935-947.

A.P-wave reflection coefficients for transversely isotropicmodels with vertical and horizontal axis of symmetry[J].Geophysics,1997,62(3):713-722.[6]

AP-wave reflection coefficients for transversely isotropicmodels with vertical and horizontal axis of symmetry[J]. Geophysics, 1997, 62(3):713-722.

[7] Hou Dongjia, Liu Yang, Ren Zhiming, et al. Multi-wave prestack joint inversion in VTI media based on Bayesian theory [J]. Geophysical Prospecting for Petroleum, 2014, 53(3): 294-303.

[8] Zhang, F., Zhang, T., and Li, X. Seismic amplitude inversion for the transversely isotropic media with vertical axis of symmetry: Geophysical Prospecting, 2019, 67, 2368-2385.

Summary of the invention

In view of the defects of the prior art, the object of the present invention is to provide a high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media.

To achieve the above object, the present invention provides the following technical solution: a high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media, comprising the following steps:

S1. Based on the PP wave prestack angle gathers and PS wave prestack angle gathers of the target study area, the PP seismic wavelet and PS seismic wavelet for inversion are extracted respectively, and the PP wave wavelet matrix and PS wave wavelet matrix are established respectively;

S2, obtaining various elastic parameters at the wellhead position of the study area, respectively calculating the reflectivity of each elastic parameter, and constructing an initial inversion model using the reflectivity of each elastic parameter at the wellhead position;

S3, according to the PP wave sub-wave matrix, the PS wave sub-wave matrix, and the inversion initial model, using the VTI medium reflection coefficient approximation formula, respectively obtain the PP wave synthesis record and the PS wave synthesis record;

S4, performing basis pursuit decomposition on the inversion initial model to obtain the odd and even component reflection coefficient vectors;

S5. For the inversion initial model, the partial derivative of the reflectivity parameter is obtained using the VTI medium reflection coefficient approximation formula;

S6. For the inversion initial model, establish a sparse constrained inversion objective function, use the alternating direction multiplier method to solve the objective function, update the model parameters, obtain the inversion results of the odd and even component reflection coefficient pair vectors, and then use the odd and even component reflection coefficient pair vectors to calculate the reflectivity results of the elastic parameters;

S7, calculating and obtaining a total error value by using the PP wave prestack angle gather, the PP wave synthetic record, the PS wave prestack angle gather, and the PS wave synthetic record;

S8, determine whether the total error value is less than the preset error value, if so, output the reflectivity result of the elastic parameter, and then perform channel integration processing on the reflectivity result of the elastic parameter to obtain the final elastic parameter result, otherwise execute step S9; S9, return to execute steps S3 to S8 until the total error value is less than the preset error value, and obtain the reflectivity result of the elastic parameter corresponding to the total error value, and then perform channel integration processing on the reflectivity result of the elastic parameter to obtain the final elastic parameter result.

Furthermore, the aforementioned step S2 is specifically as follows:

S2.1. Obtain various elastic parameters at the wellhead position in the study area: P-wave velocity V _P , S-wave velocity V _S , density DEN , anisotropy parameter EPS , anisotropy parameter DEL ;

S2.2. Calculate the reflectivity of each elastic parameter according to the following formula:

l _VP =2(V _P2 -V _P1 )/(V _P2 +V _P1 ),

l _VS =2(V _S2 -V _S1 )/(V _S2 +V _S1 ),

l _DEN = 2(DEN ₂ -DEN ₁ )/(DEN ₂ +DEN ₁ ),

_EPS = EPS ₂ - EPS ₁ ,

l _DEL = DEL ₂ - DEL ₁ ,

Among them, l _VP , l _VS , l _DEN , l _EPS and l _DEL represent the reflectivity of longitudinal wave velocity _VP , the reflectivity of shear wave velocity _VS , the reflectivity of density DEN, the reflectivity of anisotropy parameter EPS and the reflectivity of anisotropy parameter DEL, respectively; _VP1 , _VS1 , DEN ₁ , EPS ₁ and DEL ₁ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the upper medium of the reflection interface, respectively; _VP2 , _VS2 , DEN ₂ , EPS ₂ and DEL ₂ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the lower medium of the reflection interface, respectively; S2.3. Extract the low-frequency information from each reflectivity result at the wellhead, use the preset layer interpretation data as the lateral constraint, and use the reflectivity of each elastic parameter to construct the initial inversion model.

Furthermore, in the aforementioned step S3, the PP wave synthesis record is obtained according to the following formula:

D _PP =W _PP ·R _PP =W _PP ·A _PP ·L,

Among them, D _PP is the synthetic record of PP wave, R _PP is the reflection coefficient of PP wave, W _PP is the wavelet matrix of PP wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and A _PP is the PP wave reflection coefficient parameter matrix;

The PS wave synthesis record is obtained according to the following formula:

D _PS =W _PS ·R _PS =W _PS ·A _PS ·L,

Among them, D _PS is the synthetic record of PS wave, R _PS is the reflection coefficient of PP wave, W _PS is the wavelet matrix of PS wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and A _PS is the PS wave reflection coefficient parameter matrix.

Furthermore, in the aforementioned step S3, the parameter vector L composed of reflectivity is calculated according to the following formula:

分别为第1层介质对应的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS的反射率、各向异性参数DEL的反射率；

分别为第n层介质对应的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS的反射率、各向异性参数DET的反射率，T为转置符号。Among them, l _VP , l _VS , l _DEN , l _EPS , and l _DEL are parameter vectors composed of the reflectivity of longitudinal wave velocity, the reflectivity of shear wave velocity, the reflectivity of density, the reflectivity of anisotropy parameter EPS, and the reflectivity of anisotropy parameter DEL, respectively;

They are respectively the reflectivity of the longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DEL corresponding to the first layer of the medium;

They are respectively the reflectivity of the longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DET corresponding to the nth layer of medium, and T is the transposition sign.

Furthermore, in the aforementioned step S3, the PP wave reflection coefficient parameter matrix A _PP is obtained according to the following formula:

A _PP = [A _PP1 (i,θ _k ),A _PP2 (i,θ _k ),A _PP3 (i,θ _k ),A _PP4 (i,θ _k ),A _PP5 (i,θ _k )],

Wherein, A _PP1 , A _PP2 , A _PP3 , A _PP4 , and A _PP5 represent the PP wave reflection coefficient parameters corresponding to the longitudinal wave velocity, the shear wave velocity, the density, the anisotropy parameter EPS, and the anisotropy parameter DEL, respectively; i represents the i-th reflection interface; and θ _k represents the k-th incident angle;

The PS wave reflection coefficient parameter matrix A _PS is obtained according to the following formula:

A _PS =[A _PS1 (i,θ _k ),A _PS2 (i,θ _k ),A _PS3 (i,θ _k ),A _PS4 (i,θ _k ),A _PS5 (i,θ _k )],

Among them, A _PS1 , A _PS2 , A _PS3 , A _PS4 , and A _PS5 represent the PS wave reflection coefficient parameters corresponding to the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS, and anisotropy parameter DEL, respectively; i represents the i-th reflection interface, and θ _k represents the k-th incident angle.

Furthermore, in the aforementioned step S4, according to the inversion initial model L = [l _VP ,l _VS ,l _DEN ,l _EPS ,l _DEL ] ^T , basis pursuit decomposition is performed on it to obtain the even and odd component reflection coefficient pair vectors, as shown in the following formula: L = H·F

F＝[f _VP,odd ,f _VP,even ,f _VS,odd ,f _VS,even ,f _DEN,odd ,f _DEN,even ,f _EPS,odd ,f _EPS,even ,f _DEL,odd ,f _{DEL ,even} ]T,

In the formula, H is a matrix composed of odd and even components, h _odd and h _even represent the odd and even components of the decomposition dictionary; F represents the reflection coefficient pair vector composed of odd and even components; f _VP,odd , f _VS,odd , f _DEN,odd , f _EPS,odd , and f _DEL,odd are the odd component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively; f _VP,even , f _VS,even , f _DEN,even , f _EPS,even , and f _DEL,even are the even component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively.

Furthermore, in the aforementioned step S5, for the inversion initial model, the partial derivative of the VTI medium reflection coefficient approximation formula to the reflectivity parameter is used, as shown in the following formula:

分别表示PP波反射系数和PS波反射系数对反射率参数向量L的偏导数；VTI介质反射系数近似式为：in,

They represent the partial derivatives of the PP wave reflection coefficient and the PS wave reflection coefficient with respect to the reflectivity parameter vector L respectively; the approximate formula of the VTI medium reflection coefficient is:

R _PP =A _PP ·L, R _PS =A _PS ·L,

Among them, R _PP and R _PS represent the VTI medium reflection coefficient of PP wave and PS wave respectively; A _PP and A _PS represent the PP wave reflection coefficient parameter matrix and PS wave reflection coefficient parameter matrix respectively, and L is the parameter vector composed of reflectivity, that is, the inversion initial model.

Furthermore, the aforementioned step S6 includes the following sub-steps:

S6.1. For the inversion initial model L, the L ₁ -L ₂ sparse norm is selected as the constraint to establish the inversion objective function, as follows:

分别表示PP波、PS波的叠前角度道集，D_PP、D_PS分别表示PP波、PS波的合成记录；W_PP、W_PS分别表示PP波、PS波对应的子波矩阵，A_PP、A_PS分别表示PP波、PS波反射系数参数矩阵；H为奇偶分量组成的矩阵；α₁、α₂分别为PP波、PS波数据的权重参数，κ为正则化参数，β为调节因子；D_PP＝W_PP·A_PP·H·F，D_PS＝W_PS·A_PS·H·F；Wherein, J(F) represents the inversion objective function; F represents the pair vector of odd and even component reflection coefficients, i.e., the inversion objective parameter; the symbol ||*|| ₁ represents the L ₁ norm, and the symbol ||*|| ₂ represents the L ₂ norm;

denote the prestack angle gathers of PP wave and PS wave, respectively; D _PP and D _PS denote the synthetic records of PP wave and PS wave, respectively; W _PP and W _PS denote the wavelet matrices corresponding to PP wave and PS wave, respectively; A _PP and A _PS denote the reflection coefficient parameter matrices of PP wave and PS wave, respectively; H is the matrix composed of odd and even components; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data, respectively; κ is the regularization parameter, and β is the adjustment factor; D _PP = W _PP ·A _PP ·H·F, D _PS = W _PS ·A _PS ·H·F;

S6.2. Use the alternating direction multiplier method to solve the objective function, update the model parameters, and obtain the inversion result F _INV of the odd and even component reflection coefficient pair vector, as shown in the following formula:

为求解目标函数反演得到的纵波速度反射率、横波速度反射率、密度反射率和各向异性参数EPS反射率、各向异性参数DEL反射率的奇分量反射系数向量；

为求解目标函数反演得到的纵波速度反射率、横波速度反射率、密度反射率和各向异性参数EPS反射率、各向异性参数DEL反射率的偶分量反射系数向量。in,

To solve the objective function, the odd-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropic parameter EPS reflectivity, and anisotropic parameter DEL reflectivity are obtained by inversion.

To solve the objective function, the even-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity are obtained by inversion.

S6.3. Using the even and odd component reflection coefficient pair results, the reflectivity result L _INV of the elastic parameter is calculated by inverse basis pursuit transformation, as follows:

分别表示反演得到的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS的反射率、各向异性参数DEL的反射率。in,

They respectively represent the reflectivity of the inverted longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DEL.

Furthermore, in the aforementioned step S7, the total error value is calculated using the PP wave prestack angle gather, the PP wave synthetic record, the PS wave prestack angle gather and the PS wave synthetic record according to the following formula:

表示PP波和PS波的叠前角度道集，

表示由反演结果计算的PP波、PS波的合成记录；α₁、α₂分别为PP波、PS波数据的权重参数；L_INV为步骤S6中得到的弹性参数的反射率结果；W_PP、W_PS表示PP波、PS波对应的子波矩阵，A_PP、A_PS分别表示PP波、PS波反射系数参数矩阵。Where misfit represents the total error, and the symbol ||*|| ₂ represents the L ₂ norm;

represents the pre-stack angle gathers of PP and PS waves,

represents the synthetic record of PP wave and PS wave calculated from the inversion results; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data respectively; L _INV is the reflectivity result of the elastic parameter obtained in step S6; W _PP and W _PS represent the wavelet matrices corresponding to PP wave and PS wave, and A _PP and A _PS represent the reflection coefficient parameter matrices of PP wave and PS wave respectively.

Compared with the prior art, the present invention has the following beneficial effects:

In the past, the inversion of VTI media mostly used Gaussian constraints to establish the objective function, and the resolution of the inversion results was insufficient. The present invention introduces basis pursuit decomposition into the anisotropic inversion of VTI media, and uses the sparse norm as a constraint to establish the objective function, which can effectively improve the resolution of the inversion results of velocity, density and anisotropy parameters, and provide data support for the subsequent fine exploration of shale reservoirs.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG1 is a flow chart of the high-resolution anisotropic multi-wave joint inversion method of the present invention.

Figure 2 shows the five-parameter results of conventional prestack inversion based on the VTI medium reflection coefficient approximation.

FIG. 3 is a graph showing five parameter results of multi-wave joint prestack inversion of VTI media based on basis pursuit decomposition according to the present invention.

FIG. 4 is a five-parameter result diagram of the VTI anisotropic multi-wave joint inversion of the present invention under noise test.

DETAILED DESCRIPTION

In order to better understand the technical content of the present invention, specific embodiments are given and described as follows in conjunction with the accompanying drawings.

Various aspects of the invention are described herein with reference to the accompanying drawings, in which many illustrative embodiments are shown. The embodiments of the invention are not limited to those described in the accompanying drawings. It should be understood that the invention is implemented by any of the various concepts and embodiments described above, as well as the concepts and embodiments described in detail below, because the concepts and embodiments disclosed in the invention are not limited to any implementation. In addition, some aspects disclosed in the invention may be used alone or in any appropriate combination with other aspects disclosed in the invention.

As shown in FIG1 , the high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media of the present invention comprises the following steps:

S1. Based on the PP wave prestack angle gathers and PS wave prestack angle gathers of the target study area, the PP seismic wavelet and PS seismic wavelet for inversion are extracted respectively, and the PP wave wavelet matrix and PS wave wavelet matrix are established respectively. Specifically:

Wherein, wPP(θ _i ) represents the wavelet sequence extracted from the PP wave data and varying with the angle θ _i (i=1, 2, ...., n).

Extract the PS wavelet sequence from the PS wave prestack seismic data and establish the PS wavelet matrix W _PS for inversion. The specific expression is:

Wherein, w _PS (θ _i ) represents the wavelet sequence extracted from the PS wave data and varying with the angle θ _i (i=1, 2, ..., n). S2. Obtain each elastic parameter at the wellhead position in the study area, calculate the reflectivity of each elastic parameter respectively, and construct the initial inversion model using the reflectivity of each elastic parameter at the wellhead position;

First, the elastic parameters at the wellhead position are obtained, including: longitudinal wave velocity V _P , shear wave velocity V _S , density DEN , anisotropy parameter EPS , and anisotropy parameter DEL .

Then the corresponding reflectivity is calculated based on the elastic parameters as follows:

l _VP =2(V _P2 -V _P1 )/(V _P2 +V _P1 ),

l _VS =2(V _S2 -V _S1 )/(V _S2 +V _S1 ),

l _DEN ＝2(DEN ₂ -DEN ₁ )/(DEN ₂ +DEN ₁ ), (3)

_EPS = EPS ₂ - EPS ₁ ,

l _DEL = DEL ₂ - DEL ₁

Among them, l _VP , l _VS , l _DEN , l _EPS and l _DEL represent the reflectivity of longitudinal wave velocity _VP , shear wave velocity _VS , density DEN, anisotropy parameter EPS and anisotropy parameter DEL, respectively; _VP1 , _VS1 , DEN ₁ , EPS ₁ and DEL ₁ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the upper medium of the reflection interface, respectively; _VP2 , _VS2 , DEN ₂ , EPS ₂ and DEL ₂ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the lower medium of the reflection interface, respectively.

Finally, the Backus average algorithm is used to extract the low-frequency information in each reflectivity result at the wellhead. The preset horizon interpretation data is used as the lateral constraint, and the Kriging interpolation method is used to obtain the initial inversion model used for inversion.

S3, according to the PP wave sub-wave matrix, the PS wave sub-wave matrix, and the inversion initial model, using the VTI medium reflection coefficient approximation formula, respectively obtain the PP wave synthesis record and the PS wave synthesis record; specifically including:

The PP wave synthesis record is obtained according to the following formula:

D _PP =W _PP ·R _PP =W _PP ·A _PP ·L, (4.1)

Among them, _DPP is the synthetic record of PP wave, _RPP is the reflection coefficient of PP wave, _WPP is the wavelet matrix of PP wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and _APP is the PP wave reflection coefficient parameter matrix.

The PS wave synthesis record is obtained according to the following formula:

D _PS =W _PS ·R _PS =W _PS ·A _PS ·L, (4.2)

Among them, D _PS is the synthetic record of PS wave, R _PS is the reflection coefficient of PP wave, W _PS is the wavelet matrix of PS wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and A _PS is the PS wave reflection coefficient parameter matrix.

In formula (4.1) and formula (4.2), the parameter vector L composed of reflectivity is calculated according to the following formula:

分别为第1层介质对应的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS反射率、各向异性参数DEL的反射率；

分别为第n层介质对应的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS的反射率、各向异性参数DET的反射率，T为转置符号。Among them, l _VP , l _VS , l _DEN , l _EPS , and l _DEL are parameter vectors composed of the reflectivity of longitudinal wave velocity, the reflectivity of shear wave velocity, the reflectivity of density, the reflectivity of anisotropy parameter EPS, and the reflectivity of anisotropy parameter DEL, respectively;

They are respectively the reflectivity of longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS, and anisotropy parameter DEL corresponding to the first layer of medium;

They are respectively the reflectivity of the longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DET corresponding to the nth layer of medium, and T is the transposition sign.

In formula (4.1), the PP wave reflection coefficient parameter matrix A _PP is obtained according to the following formula:

A _PP = [A _PP1 (i,θ _k ),A _PP2 (i,θ _k ),A _PP3 (i,θ _k ),A _PP4 (i,θ _k ),A _PP5 (i,θ _k )], (6)

Wherein, A _PP1 , A _PP2 , A _PP3 , A _PP4 , and A _PP5 represent the PP wave reflection coefficient parameters corresponding to the longitudinal wave velocity, the shear wave velocity, the density, the anisotropy parameter EPS, and the anisotropy parameter DEL, respectively; i represents the i-th reflection interface; and θ _k represents the k-th incident angle;

In formula (4.2), the PS wave reflection coefficient parameter matrix A _PS is obtained according to the following formula:

A _PS =[A _PS1 (i,θ _k ),A _PS2 (i,θ _k ),A _PS3 (i,θ _k ),A _PS4 (i,θ _k ),A _PS5 (i,θ _k )] ( 7)

Among them, A _PS1 , A _PS2 , A _PS3 , A _PS4 , and A _PS5 represent the PS wave reflection coefficient parameters corresponding to the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS, and anisotropy parameter DEL, respectively; i represents the i-th reflection interface, and θ _k represents the k-th incident angle.

In formula (6), A _PP1 , A _PP2 , A _PP3 , A _PP4 and A _PP5 are obtained according to the following formulas:

In formula (7), A _PS1 , A _PS2 , A _PS3 , A _PS4 and A _PS5 are obtained according to the following formulas:

Wherein, V _P1 and V _S1 represent the longitudinal wave velocity and transverse wave velocity of the upper layer medium of the reflection interface, respectively; V _P2 and V _S2 represent the longitudinal wave velocity and transverse wave velocity of the lower layer medium of the reflection interface, respectively; θ and φ represent the incident angles of the longitudinal wave and transverse wave, respectively. S4. For the initial inversion model, perform basis pursuit decomposition on it to obtain the vector of even and odd component reflection coefficients;

The constructed initial inversion model is characterized by a model parameter vector L composed of reflectivities:

L＝[l _VP ,l _VS ,l _DEN ,l _EPS ,l _DEL ] ^T (10)

Among them, l _VP , l _VS , l _DEN , l _EPS , and l _DEL represent parameter vectors composed of the reflectivity of longitudinal wave velocity, the reflectivity of shear wave velocity, the reflectivity of density, the reflectivity of anisotropy parameter EPS, and the reflectivity of anisotropy parameter DEL, respectively. Taking l _VP as an example, its specific form is:

分别代表第1、2和n个反射界面对应的纵波速度反射率。In the formula,

Represent the longitudinal wave velocity reflectivities corresponding to the 1st, 2nd and nth reflection interfaces respectively.

Perform basis pursuit decomposition on the reflectivity sequence to obtain the reflection coefficient pair vectors of the odd and even components:

L＝H·F (12)

in,

F＝[f _VP,odd ,f _VP,even ,f _VS,odd ,f _VS,even ,f _DEN,odd ,f _DEN,even ,f _EPS,odd ,f _EPS,even ,f _DEL,odd ,f _{DEL ,even} ] ^T (14)

In the formula, H is a matrix composed of odd and even components, h _odd and h _even represent the odd and even components of the decomposition dictionary; F represents the reflection coefficient pair vector composed of odd and even components; f _VP,odd , f _VS,odd , f _DEN,odd , f _EPS,odd , and f _DEL,odd are the odd component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively; f _VP,even , f _VS,even , f _DEN,even , f _EPS,even , and f _DEL,even are the even component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively.

Taking f _VP,odd and f _VP,even as examples, their specific forms are as follows:

分别表示在第1、2和n个反射界面纵波速度反射率的奇分量反射系数值；

分别表示在第1、2和n个反射界面处纵波速度反射率的偶分量反射系数值。in,

Respectively represent the odd component reflection coefficient values of the longitudinal wave velocity reflectivity at the 1st, 2nd and nth reflection interfaces;

They represent the even component reflection coefficient values of the longitudinal wave velocity reflectivity at the 1st, 2nd and nth reflection interfaces respectively.

S5. For the initial inversion model, the partial derivative of the reflectivity parameter is calculated using the VTI medium reflection coefficient approximation formula. Specifically, the partial derivative of the reflectivity parameter is calculated using the VTI medium reflection coefficient approximation formula, as shown in the following formula:

分别表示PP波反射系数和PS波反射系数对反射率参数向量L的偏导数；VTI介质反射系数近似式为：in,

They represent the partial derivatives of the PP wave reflection coefficient and the PS wave reflection coefficient with respect to the reflectivity parameter vector L respectively; the approximate formula of the VTI medium reflection coefficient is:

R _PP =A _PP ·L, R _PS =A _PS ·L, (18)

Among them, R _PP and R _PS represent the VTI medium reflection coefficient of PP wave and PS wave respectively; A _PP and A _PS represent the reflection coefficient parameter matrices of PP wave and PS wave respectively, and L is the parameter vector composed of reflectivity, that is, the initial inversion model.

S6. For the initial inversion model, establish a sparse constrained inversion objective function, use the alternating direction multiplier method to solve the objective function, update the model parameters, obtain the inversion results of the odd and even component reflection coefficient pair vectors, and then use the odd and even component reflection coefficient pair vectors to calculate the reflectivity results of the elastic parameters, specifically:

First, for the inversion initial model L, the L ₁ -L ₂ sparse norm is selected as the constraint to establish the inversion objective function, as follows:

分别表示PP波、PS波的叠前角度道集，D_PP、D_PS分别表示PP波、PS波的合成记录；W_PP、W_PS分别表示PP波、PS波对应的子波矩阵，A_PP、A_PS分别表示PP波、PS波反射系数参数矩阵；H为奇偶分量组成的矩阵；α₁、α₂分别为PP波、PS波数据的权重参数，κ为正则化参数，β为调节因子；D_PP＝W_PP·A_PP·H·F，D_PS＝W_PS·A_PS·H·F；Wherein, J(F) represents the inversion objective function; F represents the pair vector of odd and even component reflection coefficients, i.e., the inversion objective parameter; the symbol ||*|| ₁ represents the L ₁ norm, and the symbol ||*|| ₂ represents the L ₂ norm;

denote the prestack angle gathers of PP wave and PS wave, respectively; D _PP and D _PS denote the synthetic records of PP wave and PS wave, respectively; W _PP and W _PS denote the wavelet matrices corresponding to PP wave and PS wave, respectively; A _PP and A _PS denote the reflection coefficient parameter matrices of PP wave and PS wave, respectively; H is the matrix composed of odd and even components; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data, respectively; κ is the regularization parameter, and β is the adjustment factor; D _PP = W _PP ·A _PP ·H·F, D _PS = W _PS ·A _PS ·H·F;

In the second step, the alternating direction multiplier method is used to solve the objective function, update the model parameters, and obtain the inversion result F _INV of the odd and even component reflection coefficient pair vector, as shown in the following formula:

为求解目标函数反演得到的纵波速度反射率、横波速度反射率、密度反射率和各向异性参数EPS反射率、各向异性参数DEL反射率的奇分量反射系数向量；

为求解目标函数反演得到的纵波速度反射率、横波速度反射率、密度反射率和各向异性参数EPS反射率、各向异性参数DEL反射率的偶分量反射系数向量。in,

To solve the objective function, the odd-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropic parameter EPS reflectivity, and anisotropic parameter DEL reflectivity are obtained by inversion.

In order to solve the objective function, the even-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity and anisotropy parameter DEL reflectivity are obtained by inversion.

Finally, the reflectivity result L _INV of the elastic parameter is calculated by using the odd and even component reflection coefficient pair results through the inverse basis pursuit transformation, as shown in the following formula:

分别表示反演得到的纵波速度的反射率、横波速度的反射率、密度的反射率、各向异性参数EPS的反射率、各向异性参数DEL的反射率。S7、利用PP波叠前角度道集、PP波合成记录、PS波叠前角度道集与PS波合成记录计算获得总误差值；in,

They respectively represent the reflectivity of the P-wave velocity, the reflectivity of the S-wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DEL obtained by inversion. S7. The total error value is calculated using the PP wave prestack angle gather, the PP wave synthetic record, the PS wave prestack angle gather and the PS wave synthetic record;

The total error value is calculated according to the following formula:

表示PP波和PS波的叠前角度道集，

表示由反演结果计算的PP波、PS波的合成记录；α₁、α₂分别为PP波、PS波数据的权重参数；L_INV为步骤S6中得到的弹性参数的反射率结果；W_PP、W_PS表示PP波、PS波对应的子波矩阵，A_PP、A_PS分别表示PP波和PS波反射系数参数矩阵。Where misfit represents the total error, and the symbol ||*|| ₂ represents the L ₂ norm;

represents the pre-stack angle gathers of PP and PS waves,

represents the synthetic record of PP wave and PS wave calculated from the inversion results; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data respectively; L _INV is the reflectivity result of the elastic parameter obtained in step S6; W _PP and W _PS represent the wavelet matrices corresponding to PP wave and PS wave, and A _PP and A _PS represent the reflection coefficient parameter matrices of PP wave and PS wave respectively.

S8, determine whether the total error value is less than the preset error value, if yes, output the reflectivity result of the elastic parameter, then perform channel integration processing on the reflectivity result of the elastic parameter to obtain the final result of the elastic parameter, otherwise execute step S9, S9, return to execute steps S3 to S8, until the total error value is less than the preset error value, and obtain the reflectivity result of the elastic parameter corresponding to the total error value, then perform channel integration processing on the reflectivity result of the elastic parameter to obtain the final result of the elastic parameter, including longitudinal wave velocity _VP , shear wave velocity _VS , density DEN, anisotropy parameter EPS, and anisotropy parameter DEL.

Fig. 2 shows the five-parameter results of conventional prestack inversion based on the VTI medium reflection coefficient approximation, including P-wave velocity _VP , _S -wave velocity VS, density DEN, anisotropy parameter EPS and anisotropy parameter DEL. The black dotted line is the initial inversion model, the gray solid line is the real logging data, and the black solid line is the five-parameter inversion result. It can be seen that the conventional inversion result is not ideal, the overall trend is relatively smooth, the resolution is low, and the detailed information cannot be highlighted.

FIG3 is a five-parameter result of multi-wave joint prestack inversion of VTI media based on basis pursuit decomposition of the present invention, which also includes P-wave velocity _VP , _S -wave velocity VS, density DEN, anisotropy parameter EPS and anisotropy parameter DEL; wherein the black dotted line is the initial inversion model, the gray solid line is the real logging data, and the black solid line is the five-parameter inversion result; it can be seen that the result has a great improvement in resolution compared with the conventional method in FIG2, especially it can highlight many detail changes.

FIG4 is a five-parameter result of the VTI anisotropic multi-wave joint inversion of the present invention under the noise test; the black dotted line is the initial inversion model, the gray solid line is the real logging data, and the black solid line is the five-parameter inversion result; it can be seen that even when noise is added to the input data, the method of the present invention can still obtain an inversion result with a higher resolution. Although the accuracy is lower than that of FIG3 (no noise added), it is still significantly better than the conventional method in FIG2, which verifies that the present invention has good stability.

Although the present invention has been described above with preferred embodiments, it is not intended to limit the present invention. A person skilled in the art of the present invention may make various modifications and improvements without departing from the spirit and scope of the present invention. Therefore, the scope of protection of the present invention shall be determined by the claims.

Claims (9)

1. A high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media, characterized in that it comprises the following steps: S1. Based on the PP wave prestack angle gathers and PS wave prestack angle gathers of the target study area, the PP seismic wavelet and PS seismic wavelet for inversion are extracted respectively, and the PP wave wavelet matrix and PS wave wavelet matrix are established respectively; S2, obtaining various elastic parameters at the wellhead position of the study area, respectively calculating the reflectivity of each elastic parameter, and constructing an initial inversion model using the reflectivity of each elastic parameter at the wellhead position; S3, according to the PP wave sub-wave matrix, the PS wave sub-wave matrix, and the inversion initial model, using the VTI medium reflection coefficient approximation formula, respectively obtain the PP wave synthesis record and the PS wave synthesis record; S4, performing basis pursuit decomposition on the inversion initial model to obtain the odd and even component reflection coefficient vectors; S5. For the inversion initial model, the partial derivative of the reflectivity parameter is obtained using the VTI medium reflection coefficient approximation formula; S6. For the inversion initial model, establish a sparse constrained inversion objective function, use the alternating direction multiplier method to solve the objective function, update the model parameters, obtain the inversion results of the odd and even component reflection coefficient pair vectors, and then use the odd and even component reflection coefficient pair vectors to calculate the reflectivity results of the elastic parameters; S7, calculating and obtaining a total error value by using the PP wave prestack angle gather, the PP wave synthetic record, the PS wave prestack angle gather, and the PS wave synthetic record; S8, judging whether the total error value is less than a preset error value, if yes, outputting the reflectivity result of the elastic parameter, and then performing channel integration processing on the reflectivity result of the elastic parameter to obtain the final result of the elastic parameter, otherwise, executing step S9; S9, return to execute step S3 to step S8, until the total error value is less than the preset error value, and obtain the reflectivity result of the elastic parameter corresponding to the total error value, and then perform channel integration processing on the reflectivity result of the elastic parameter to obtain the final elastic parameter result. 2. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 1 is characterized in that step S2 specifically comprises: S2.1. Obtain various elastic parameters at the wellhead position in the study area: P-wave velocity V _P , S-wave velocity V _S , density DEN , anisotropy parameter EPS , anisotropy parameter DEL ; S2.2. Calculate the reflectivity of each elastic parameter according to the following formula: l _VP =2(V _P2 -V _P1 )/(V _P2 +V _P1 ), l _VS =2(V _S2 -V _S1 )/(V _S2 +V _S1 ), l _DEN = 2(DEN ₂ -DEN ₁ )/(DEN ₂ +DEN ₁ ), _EPS = EPS ₂ - EPS ₁ , l _DEL = DEL ₂ - DEL ₁ , Wherein, l _VP , l _VS , l _DEN , l _EPS and l _DEL represent the reflectivity of longitudinal wave velocity _VP , the reflectivity of shear wave velocity _VS , the reflectivity of density DEN, the reflectivity of anisotropy parameter EPS and the reflectivity of anisotropy parameter DEL respectively; _VP1 , _VS1 , DEN ₁ , EPS ₁ and DEL ₁ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the upper layer medium of the reflection interface respectively; _VP2 , _VS2 , DEN ₂ , EPS ₂ and DEL ₂ represent the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS and anisotropy parameter DEL of the lower layer medium of the reflection interface respectively; S2.3. Extract the low-frequency information from each reflectivity result at the wellhead, use the preset layer interpretation data as the lateral constraint, and use the reflectivity of each elastic parameter to build the initial inversion model. 3. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 1, characterized in that, in step S3, the PP wave synthetic record is obtained according to the following formula: D _PP =W _PP ·R _PP =W _PP ·A _PP ·L, Among them, D _PP is the synthetic record of PP wave, R _PP is the reflection coefficient of PP wave, W _PP is the wavelet matrix of PP wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and A _PP is the PP wave reflection coefficient parameter matrix; The PS wave synthesis record is obtained according to the following formula: D _PS =W _PS ·R _PS =W _PS ·A _PS ·L, Among them, D _PS is the synthetic record of PS wave, R _PS is the reflection coefficient of PP wave, W _PS is the wavelet matrix of PS wave, L is the parameter vector composed of reflectivity, that is, the inversion initial model, and A _PS is the PS wave reflection coefficient parameter matrix. 4. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 3, characterized in that in step S3, the parameter vector L composed of reflectivity is calculated according to the following formula:

Among them, l _VP , l _VS , l _DEN , l _EPS , and l _DEL are parameter vectors composed of the reflectivity of longitudinal wave velocity, the reflectivity of shear wave velocity, the reflectivity of density, the reflectivity of anisotropy parameter EPS, and the reflectivity of anisotropy parameter DEL, respectively;

They are respectively the reflectivity of the longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DEL corresponding to the first layer of the medium;

They are respectively the reflectivity of the longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DET corresponding to the nth layer of medium, and T is the transposition sign. 5. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 3, characterized in that, in step S3, the PP wave reflection coefficient parameter matrix A _PP is obtained according to the following formula: A _PP = [A _PP1 (i,θ _k ),A _PP2 (i,θ _k ),A _PP3 (i,θ _k ),A _PP4 (i,θ _k ),A _PP5 (i,θ _k )], Wherein, A _PP1 , A _PP2 , A _PP3 , A _PP4 , and A _PP5 represent the PP wave reflection coefficient parameters corresponding to the longitudinal wave velocity, the shear wave velocity, the density, the anisotropy parameter EPS, and the anisotropy parameter DEL, respectively; i represents the i-th reflection interface; and θ _k represents the k-th incident angle; The PS wave reflection coefficient parameter matrix A _PS is obtained according to the following formula: A _PS =[A _PS1 (i,θ _k ),A _PS2 (i,θ _k ),A _PS3 (i,θ _k ),A _PS4 (i,θ _k ),A _PS5 (i,θ _k )], Among them, A _PS1 , A _PS2 , A _PS3 , A _PS4 , and A _PS5 represent the PS wave reflection coefficient parameters corresponding to the longitudinal wave velocity, shear wave velocity, density, anisotropy parameter EPS, and anisotropy parameter DEL, respectively; i represents the i-th reflection interface, and θ _k represents the k-th incident angle. 6. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 4, characterized in that, in step S4, according to the inversion initial model L = [l _VP , l _VS , l _DEN , l _EPS , l _DEL ] ^T , basis pursuit decomposition is performed on it to obtain the odd and even component reflection coefficient pair vectors, as shown in the following formula: L = H·F

F＝[f _VP,odd ,f _VP,even ,f _VS,odd ,f _VS,even ,f _DEN,odd ,f _DEN,even ,f _EPS,odd ,f _EPS,even ,f _DEL,odd ,f _{DEL ,even} ] ^T , In the formula, H is a matrix composed of odd and even components, h _odd and h _even represent the odd and even components of the decomposition dictionary; F represents the reflection coefficient pair vector composed of odd and even components; f _VP,odd , f _VS,odd , f _DEN,odd , f _EPS,odd , and f _DEL,odd are the odd component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively; f _VP,even , f _VS,even , f _DEN,even , f _EPS,even , and f _DEL,even are the even component reflection coefficient vectors corresponding to the longitudinal wave velocity reflectivity, shear wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity, and anisotropy parameter DEL reflectivity, respectively. 7. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 6, characterized in that, in step S5, for the inversion initial model, the partial derivative of the VTI medium reflection coefficient approximation formula to the reflectivity parameter is used, as shown in the following formula:

in,

They represent the partial derivatives of the PP wave reflection coefficient and the PS wave reflection coefficient with respect to the reflectivity parameter vector L respectively; the approximate formula of the VTI medium reflection coefficient is: R _PP =A _PP ·L, R _PS =A _PS ·L, Among them, R _PP and R _PS represent the VTI medium reflection coefficient of PP wave and PS wave respectively; A _PP and A _PS represent the PP wave reflection coefficient parameter matrix and PS wave reflection coefficient parameter matrix respectively, and L is the parameter vector composed of reflectivity, that is, the inversion initial model. 8. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 7, characterized in that step S6 comprises the following sub-steps: S6.1. For the inversion initial model L, the L ₁ -L ₂ sparse norm is selected as the constraint to establish the inversion objective function, as follows:

Wherein, J(F) represents the inversion objective function; F represents the pair vector of odd and even component reflection coefficients, i.e., the inversion objective parameter; the symbol ||*|| ₁ represents the L ₁ norm, and the symbol ||*|| ₂ represents the L ₂ norm;

denote the prestack angle gathers of PP wave and PS wave, respectively; D _PP and D _PS denote the synthetic records of PP wave and PS wave, respectively; W _PP and W _PS denote the wavelet matrices corresponding to PP wave and PS wave, respectively; A _PP and A _PS denote the reflection coefficient parameter matrices of PP wave and PS wave, respectively; H is the matrix composed of odd and even components; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data, respectively; κ is the regularization parameter, and β is the adjustment factor; D _PP = W _PP ·A _PP ·H·F, D _PS = W _PS ·A _PS ·H·F; S6.2. Use the alternating direction multiplier method to solve the objective function, update the model parameters, and obtain the inversion result F _INV of the odd and even component reflection coefficient pair vector, as shown in the following formula:

in,

To solve the objective function, the odd-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropic parameter EPS reflectivity, and anisotropic parameter DEL reflectivity are obtained by inversion.

In order to solve the objective function, the even-component reflection coefficient vectors of the P-wave velocity reflectivity, S-wave velocity reflectivity, density reflectivity, anisotropy parameter EPS reflectivity and anisotropy parameter DEL reflectivity are obtained by inversion. S6.3. Using the even and odd component reflection coefficient pair results, the reflectivity result L _INV of the elastic parameter is calculated by inverse basis pursuit transformation, as follows:

in,

They respectively represent the reflectivity of the inverted longitudinal wave velocity, the reflectivity of the shear wave velocity, the reflectivity of the density, the reflectivity of the anisotropy parameter EPS, and the reflectivity of the anisotropy parameter DEL. 9. The high-resolution multi-wave joint prestack inversion method for shale reservoir VTI media according to claim 8, characterized in that in step S7, the total error value is calculated according to the following formula using PP wave prestack angle gathers, PP wave synthetic records, PS wave prestack angle gathers and PS wave synthetic records:

Where misfit represents the total error, and the symbol ||*|| ₂ represents the L ₂ norm;

Represents the pre-stack angle gathers of PP and PS waves,

represents the synthetic record of PP wave and PS wave calculated from the inversion results; α ₁ and α ₂ are the weight parameters of PP wave and PS wave data respectively; L _INV is the reflectivity result of the elastic parameter obtained in step S6; W _PP and W _PS represent the wavelet matrices corresponding to PP wave and PS wave, and A _PP and A _PS represent the reflection coefficient parameter matrices of PP wave and PS wave respectively.

Images