Seismic inversion method and device, storage medium and electronic equipment

Abstract

The invention provides a seismic inversion method and device, a storage medium and electronic equipment, and the method comprises the steps: obtaining initial inversion data, and obtaining first azimuth seismic data and second azimuth seismic data; constructing at least one first azimuth angle positive operator, and constructing at least one second azimuth angle positive operator; constructing a first target inversion function based on the at least one first azimuth angle forward operator and the first azimuth angle seismic data; constructing a second target inversion function based on the at least one second azimuth angle forward operator and the second azimuth angle seismic data; determining a target inversion result of the to-be-inverted parameter based on the first target inversion function, the second target inversion function and the initial inversion data; and based on the target inversion result of the to-be-inverted parameter, calculating the target inversion result of the anisotropic parameter. According to the embodiment of the invention, the accuracy of seismic inversion can be improved.

Description

Technical Field

[0001] The present invention relates to the field of oil and gas seismic exploration, and in particular to a seismic inversion method, device, storage medium and electronic equipment. Background Art

[0002] Seismic waves usually show azimuthal anisotropy of different types and intensities in underground media, and the reflection and transmission responses also vary with azimuth. HTI (Horizontal Transverse Isotropy) medium, as a medium model that causes azimuthal anisotropy, is the most effective anisotropy model and plays an important role in describing seismic wave propagation. Anisotropy is important for describing the orientation and intensity of fractures and the type of filling fluid. However, the anisotropy induced by fractures contributes little to seismic reflection, and the prediction of multiple solutions directly from azimuthal seismic data is strong. In this regard, the related technology uses the inverted stiffness matrix to calculate the anisotropy parameters, which is ill-posed and leads to low accuracy of seismic inversion. Based on this, there is currently no good solution to improve the accuracy of seismic inversion. Summary of the invention

[0003] In view of this, the embodiments of the present invention provide a seismic inversion method, device, storage medium and electronic device to solve the problems such as the low accuracy of seismic inversion in related technologies; that is, the embodiments of the present invention can realize seismic inversion through a first target inversion function and a second target inversion function to predict the anisotropic parameters of the HTI medium, thereby obtaining a target inversion result of anisotropic parameters with higher accuracy, which can effectively improve the accuracy of seismic inversion.

[0004] According to one aspect of an embodiment of the present invention, a seismic inversion method is provided, the method comprising: Acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; constructing at least one first azimuth forward operator, and constructing at least one second azimuth forward operator; constructing a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data; and constructing a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data; Determining a target inversion result of a parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data; Based on the target inversion results of the parameters to be inverted, target inversion results of the anisotropic parameters are calculated.

[0005] According to another aspect of an embodiment of the present invention, a seismic inversion device is provided, the device comprising: An acquisition unit, used for acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; A processing unit, configured to construct at least one first azimuth forward operator and at least one second azimuth forward operator; The processing unit is further configured to construct a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data; and to construct a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data; The processing unit is further used to determine a target inversion result of the parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data; The processing unit is further used to calculate the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted.

[0006] According to another aspect of an embodiment of the present invention, an electronic device is provided, comprising a processor and a memory storing a program, wherein the program comprises instructions, which, when executed by the processor, enable the processor to perform the above-mentioned method.

[0007] According to another aspect of an embodiment of the present invention, a non-transitory computer-readable storage medium storing computer instructions is provided, wherein the computer instructions are used to enable a computer to execute the above-mentioned method.

[0008] According to another aspect of an embodiment of the present invention, a computer program product is provided, comprising a computer program, wherein when the computer program is executed by a processor, it is used to enable a computer to perform the above-mentioned method.

[0009] The embodiment of the present invention can construct at least one first azimuth forward operator and at least one second azimuth forward operator after obtaining the initial inversion data, and obtaining the first azimuth seismic data and the second azimuth seismic data. Then, a first target inversion function can be constructed based on the at least one first azimuth forward operator and the first azimuth seismic data; and a second target inversion function can be constructed based on the at least one second azimuth forward operator and the second azimuth seismic data. Further, the target inversion result of the parameter to be inverted can be determined based on the first target inversion function, the second target inversion function and the initial inversion data; and the target inversion result of the anisotropic parameter can be calculated based on the target inversion result of the parameter to be inverted. It can be seen that the embodiment of the present invention can realize seismic inversion through the first target inversion function and the second target inversion function to predict the anisotropic parameters of the HTI medium, thereby obtaining the target inversion result of the anisotropic parameters with higher accuracy, which can effectively improve the accuracy of seismic inversion. BRIEF DESCRIPTION OF THE DRAWINGS

[0010] Further details, features and advantages of the invention are disclosed in the following description of exemplary embodiments in conjunction with the accompanying drawings, in which: Figure 1 A schematic flow chart of a seismic inversion method according to an exemplary embodiment of the present invention is shown; Figure 2 A schematic flow chart of another seismic inversion method according to an exemplary embodiment of the present invention is shown; Figure 3 A schematic diagram showing a reflection coefficient according to an exemplary embodiment of the present invention; Figure 4 A schematic diagram showing another reflection coefficient according to an exemplary embodiment of the present invention; Figure 5 A method for C according to an exemplary embodiment of the present invention is shown. 11 Schematic diagram of contribution analysis; Figure 6 A method for C according to an exemplary embodiment of the present invention is shown. 13 Schematic diagram of contribution analysis; Figure 7 A method for C according to an exemplary embodiment of the present invention is shown. 33 Schematic diagram of contribution analysis; Figure 8 A method for C according to an exemplary embodiment of the present invention is shown. 44 Schematic diagram of contribution analysis; Fig. 9 A method for C according to an exemplary embodiment of the present invention is shown. 55 Schematic diagram of contribution analysis; Fig.10 A schematic diagram showing a contribution analysis for ρ according to an exemplary embodiment of the present invention is shown; Fig.11 A schematic diagram showing an inversion result according to an exemplary embodiment of the present invention; Fig.12 A schematic diagram showing another inversion result according to an exemplary embodiment of the present invention; Fig.13 An anisotropy parameter δ according to an exemplary embodiment of the present invention is shown. (v) Schematic diagram of Fig.14 Another anisotropy parameter δ according to an exemplary embodiment of the present invention is shown. (v) Schematic diagram of Fig.15 An anisotropy parameter ε according to an exemplary embodiment of the present invention is shown.(v) Schematic diagram of Fig.16 Another anisotropy parameter ε according to an exemplary embodiment of the present invention is shown. (v) Schematic diagram of Fig.17 An anisotropy parameter γ according to an exemplary embodiment of the present invention is shown. (v) Schematic diagram of Fig.18 Another anisotropy parameter γ according to an exemplary embodiment of the present invention is shown. (v) Schematic diagram of . DETAILED DESCRIPTION

[0011] Embodiments of the present invention will be described in more detail below with reference to the accompanying drawings. Although certain embodiments of the present invention are shown in the accompanying drawings, it should be understood that the present invention can be implemented in various forms and should not be construed as being limited to the embodiments described herein, which are instead provided for a more thorough and complete understanding of the present invention. It should be understood that the drawings and embodiments of the present invention are only for exemplary purposes and are not intended to limit the scope of protection of the present invention.

[0012] It should be understood that the various steps described in the method embodiments of the present invention may be performed in different orders and / or in parallel. In addition, the method embodiments may include additional steps and / or omit the steps shown. The scope of the present invention is not limited in this respect.

[0013] The term "including" and its variations used in this document are open inclusions, that is, "including but not limited to". The term "based on" means "based at least in part on". The term "one embodiment" means "at least one embodiment"; the term "another embodiment" means "at least one other embodiment"; the term "some embodiments" means "at least some embodiments". Relevant definitions of other terms will be given in the following description. It should be noted that the concepts of "first", "second", etc. mentioned in the present invention are only used to distinguish different devices, modules or units, and are not used to limit the order or interdependence of the functions performed by these devices, modules or units.

[0014] It should be noted that the modifications of "one" and "plurality" mentioned in the present invention are illustrative rather than restrictive, and those skilled in the art should understand that, unless otherwise clearly indicated in the context, it should be understood as "one or more".

[0015] The names of the messages or information exchanged between multiple devices in the embodiments of the present invention are only used for illustrative purposes, and are not used to limit the scope of these messages or information.

[0016] It should be noted that the execution subject of the seismic inversion method provided in the embodiment of the present invention may be one or more electronic devices, which is not limited by the present invention; wherein, the electronic device may be a terminal (i.e., a client) or a server, then when the execution subject includes multiple electronic devices, and the multiple electronic devices include at least one terminal and at least one server, the seismic inversion method provided in the embodiment of the present invention may be jointly executed by the terminal and the server. Accordingly, the terminal mentioned here may include but is not limited to: smart phones, laptops, desktop computers, intelligent voice interaction devices, and the like. The server mentioned here may be an independent physical server, or a server cluster or distributed system composed of multiple physical servers, or a cloud server that provides cloud services, cloud databases, cloud computing (cloud computing), cloud functions, cloud storage, network services, cloud communications, middleware services, domain name services, security services, CDN (Content Delivery Network), and basic cloud computing services such as big data and artificial intelligence platforms, and the like.

[0017] Based on the above description, an embodiment of the present invention proposes a seismic inversion method, which can be executed by the electronic device (terminal or server) mentioned above; or, the seismic inversion method can be executed by the terminal and the server together. For the convenience of explanation, the following description will be given by taking the electronic device executing the seismic inversion method as an example; Figure 1 As shown, the seismic inversion method may include the following steps S101-S105: S101, obtaining initial inversion data, and obtaining first azimuth seismic data and second azimuth seismic data.

[0018] Optionally, the initial inversion data may be the initial inversion data of any area; correspondingly, the first azimuth seismic data and the second azimuth seismic data may be the first azimuth seismic data and the second azimuth seismic data of any area. Optionally, the initial inversion data may include the initial value of each benchmark parameter in each sampling medium layer in at least one sampling medium layer among multiple benchmark parameters, and may also include the initial inversion parameter data (i.e., initial value) of the parameter to be inverted, and so on; the embodiment of the present invention is not limited to this. For the sake of ease of explanation, the following description will be based on the example that the initial inversion data may include the initial value of each benchmark parameter in each sampling medium layer. Optionally, multiple benchmark parameters may include, but are not limited to: shear wave velocity (V S ), longitudinal wave velocity (V P ), density (ρ), and anisotropy parameter (δ (v) , ε (v) , γ (v) )wait.

[0019] Optionally, the initial inversion data, the first azimuth seismic data, and the second azimuth seismic data may be data in a well logging model (i.e., data generated in a model test), or data in actual seismic data, etc.; this is not limited in the embodiment of the present invention. Optionally, the initial inversion data, the first azimuth seismic data, and the second azimuth seismic data may be obtained in the following ways, but are not limited to: The first acquisition method: the electronic device's own storage space can store initial inversion data, first azimuth seismic data and second azimuth seismic data of each of multiple areas. In this case, the initial inversion data, first azimuth seismic data and second azimuth seismic data of any area can be acquired from its own storage space.

[0020] The second acquisition method: the electronic device can acquire the well logging data of any area, so as to acquire the initial inversion data, the first azimuth seismic data and the second azimuth seismic data through the well logging data.

[0021] The third acquisition method: the electronic device can display a data input interface. In this case, the user can perform a data input operation on the data input interface to input initial inversion data, first azimuth seismic data, and second azimuth seismic data through the data input operation. In this case, the electronic device can respond to the data input operation and acquire the initial inversion data, first azimuth seismic data, second azimuth seismic data, and so on.

[0022] S102: construct at least one first azimuth forward operator, and construct at least one second azimuth forward operator.

[0023] It should be noted that the embodiments of the present invention can establish the precise elastic impedance tensor of HTI medium P wave (i.e. longitudinal wave), SV wave (i.e. wave in which particle vibration occurs in a plane perpendicular to the wave propagation plane) and SH wave (i.e. wave in which particle vibration occurs in a plane parallel to the wave propagation plane) based on the relationship between stress and strain (i.e. the relationship between rock stress and strain); the elastic impedance tensor is defined by the stress field and the velocity field, and can fully describe the kinematic characteristics of the rock, and can accurately determine the reflection and transmission coefficient at any incident angle and any elastic contrast interface. Among them, in the coordinate system with the z axis vertically downward, the relationship between stress and strain can be expressed as shown in Formula 1.1: Formula 1.1 Where σ and u represent stress and displacement (σ xz represents the stress in the x-direction and the z-direction, σ yz represents the stress in the y and z directions, σ zz represents the stress in the z direction and the z direction, u x represents the displacement in the x direction, uy represents the displacement in the y direction, u z represents the displacement in the z direction), C represents the stiffness matrix of the medium, which can be defined by the stress-strain relationship, C ab It can represent the elements in the stiffness matrix (which can be called stiffness coefficients, a and b are positive integers). A stiffness matrix can include the values ​​of each stiffness coefficient in multiple stiffness coefficients. In the embodiment of the present invention, the HTI medium requires 5 independent elastic parameters to characterize, and its stiffness coefficient can be expressed as shown in formula 1.2: Formula 1.2 Among them, V S represents the shear wave velocity, V P represents the longitudinal wave velocity, ρ represents the density, δ (v) , ε (v) , γ (v) Represents the anisotropy parameter. Accordingly, according to the displacement equation of seismic wave propagation, Formula 1.1 can be written in the form of a matrix-vector, as shown in Formula 1.3: Formula 1.3 Where i is an imaginary unit, w is the angular frequency, and τ is the phase function whose gradient is the normal direction of the wavefront. The velocity field can also be expressed in terms of displacement as shown in Formula 1.4: Formula 1.4 Among them, v is the propagation speed of seismic waves in the medium, and t represents time. Based on this, formula 1.3 is the stress field related to displacement, formula 1.4 is the characteristic of the velocity field, and the relationship between the stress field, velocity field and elastic impedance tensor is shown in formula 1.5: Formula 1.5 Where E represents the elastic impedance tensor of the medium. Based on this, by combining Formula 1.3-Formula 1.5, the specific expression of the elastic impedance tensor of the HTI medium can be obtained as shown in Formula 1.6: Formula 1.6 The subscript m (m=P, SV, SH) can represent P waves, SV waves or SH waves, and “+” and “-” represent the downward and upward propagation of seismic waves, i.e., downgoing and upgoing waves; x and p y They represent the horizontal slowness along the x and y axes respectively, and q represents the vertical slowness, which can be specifically expressed as shown in formula 1.7: Formula 1.7 From the above, it can be seen that Formula 1.6 gives the precise elastic impedance tensor of P waves, SV waves and SH waves of HTI media, which can be expressed by stiffness coefficients and slowness; based on this, the elastic impedance tensor is defined by the stress field and the velocity field, which can fully describe the kinematic characteristics of the rock, and can accurately determine the reflection and transmission coefficients at any incident angle and any elastic contrast interface. Furthermore, based on the precise elastic impedance tensor of the HTI medium, the embodiment of the present invention derives the precise and weak contrast approximate expressions of the reflection coefficient represented by the stiffness matrix, constructs the propagation matrix of seismic wave propagation in the HTI medium, simplifies the precise reflection coefficient equation by using the relationship between the uplink and downlink propagation matrices, and derives the approximate formula of the reflection coefficient of the HTI medium under the weak contrast assumption, thereby improving the accuracy of the reflection coefficient at large angles. Correspondingly, the polarization vector of the HTI medium can be expressed as shown in Formula 1.8: Formula 1.8 Where U is the polarization vector, G P and G SV is the coefficient of the polarization vector, which can be specifically expressed as shown in formula 1.9: Formula 1.9 Wherein, the subscript P(SV) indicates that the subscript is P (indicating P waves) or SV (indicating SV waves). Correspondingly, the propagation matrix is ​​an important function that describes the propagation characteristics of seismic waves in the medium, which can be expressed by the elastic impedance tensor and the polarization vector as shown in Formula 1.10: Formula 1.10 Wherein K is the HTI medium propagation matrix. In the embodiment of the present invention, it can be considered that the relationship between the uplink and downlink elastic impedance tensors and the polarization vector of the HTI medium is as shown in Formula 1.11: Formula 1.11 Among them, the matrix F is a diagonal matrix. Substituting formula 1.11 into formula 1.10, the relationship between the uplink and downlink wave propagation matrices can be obtained as shown in formula 1.12: Formula 1.12 Based on this, the propagation matrix and polarization vector can be used to characterize the precise reflection coefficient equation of the HTI medium, as shown in formula 1.13: Formula 1.13 Wherein, the superscripts (1) and (2) represent the upper and lower media, and R is the reflection coefficient composed of qP wave (i.e. P wave), qSV wave (i.e. SV wave) and qSH wave (i.e. SH wave), which can be specifically expressed as shown in Formula 1.14: Formula 1.14 Among them, AB is used to indicate AB waves (A and B are P, SV or SH respectively), that is, the incident wave is A wave and the reflected wave is B wave, such as PP is used to indicate PP wave (that is, the incident wave is P wave and the reflected wave is also P wave), PSV is used to indicate PSV wave (that is, the incident wave is P wave and the reflected wave is SV wave), and so on; R AB is the reflection coefficient of the AB wave. Based on this, substituting formula 1.12 into formula 1.13, the exact reflection coefficient equation can be simplified to formula 1.15: Formula 1.15 Among them, the matrices S and Q can be expressed by the propagation matrix and polarization vector as shown in formula 1.16: Formula 1.16 Wherein, the superscript n (n=1,2) indicates the upper or lower layer medium. Assuming that the elastic and anisotropic parameters of the upper and lower layers vary slightly, that is, Q=(Q (1) +Q (2) ) / 2, S=(S (1) +S (2) ) / 2, ∆Q=Q (2) -Q (1) , ∆S=S (2) -S (1) If true, formula 1.15 can be further simplified to formula 1.17: Formula 1.17 Based on this, according to the elastic dynamics equation of seismic wave propagation, the matrix S and Q have a relationship Q T S=S T Q=Λ, Λ is the scale factor. Therefore, the PP wave reflection coefficient of HTI medium can be approximately expressed as shown in formula 1.18: Formula 1.18 Among them, Q P , S P ,Λ P It can be shown as formula 1.19: Formula 1.19 Based on this, by substituting formula 1.19 into formula 1.18, the display form of the PP wave reflection coefficient of the HTI medium can be obtained as shown in formula 1.20: Formula 1.20 Among them, v 11 、v 13 、v 33 、v 44 and v 55 It can be shown as formula 1.21: Formula 1.21 Among them, the first coefficient term Λ 1 and the second coefficient term Λ 2 It can be shown as formula 1.22: Formula 1.22 Among them, the above Γ 11 =C 11 / C 33 , Γ 13 =C 13 / C 33 , Γ 44 =C 44 / C 33， Γ 55 =C 55 / C 33 , Γ ρ =C ρ / C 33 It should be understood that Formula 1.20 and Formula 1.21 are weak contrast approximations of the PP wave reflection coefficient of the HTI medium. When the anisotropy parameter is 0, the derived reflection coefficient is equivalent to the Aki approximation equation in the isotropic medium.

[0024] Furthermore, an embodiment of the present invention establishes a forward model of a fracture plane (i.e., a first azimuth angle of 90°) and a symmetry axis plane (i.e., a second azimuth angle of 0°), constructs a target functional (i.e., a first target inversion function and a second target inversion function) that considers the mapping relationship between the stiffness matrix and the anisotropic parameters; and analyzes the variation law of the contribution of each parameter to the reflection coefficient with the azimuth, determines the dominant azimuth for seismic inversion in the fracture plane and the symmetry axis plane, and introduces a mapping relationship between the stiffness matrix and the anisotropic parameters, thereby overcoming the problem of ill-posedness in calculating the anisotropic parameters using the stiffness matrix.

[0025] It should be noted that according to the derived HTI medium PP wave reflection coefficient equation (Formula 1.20 and Formula 1.21), when the azimuth (i.e., azimuth angle) is 90°, the slowness p x =0, then the coefficient term is as shown in formula 1.23: Formula 1.23 Based on this, parameter C 11 , C 13 and C 55 There is no contribution to the reflection coefficient in the crack plane, which means that only C 33 , C 44 and ρ, the number of parameters to be inverted is reduced from 6 to 3. 11 , C 13 and C 55The contribution to the reflection coefficient is the largest in the plane of the symmetry axis. Therefore, the embodiment of the present invention can use a step-by-step inversion method to predict the stiffness matrix and density: the first step is to predict C in the crack plane. 33 , C 44 and ρ; the second step is to predict C in the plane of the symmetry axis 11 , C 13 and C 55 .

[0026] Furthermore, the disturbance terms of each parameter in the reflection coefficient can be written in the form of logarithmic difference, as shown in formula 1.24: Formula 1.24 Then, according to the reflection coefficient equation, the forward model in the crack plane is established as shown in formula 1.25: Formula 1.25 Among them, the superscript 1 =90° is the first azimuth angle, Represents seismic data with an azimuth of 90° (also called first azimuth seismic data), which is composed of h incident angles θ under the first azimuth 1 ,θ 2 …θ h The gathers are composed of X 1 =[X 33 ,X 44 ,X ρ ] T It can be the vector form of the first parameter to be inverted, which is represented by the parameter C 33 , C 44 and the logarithm of ρ; Optionally, the number of the first forward operators in the at least one first forward operator may be one; Optionally, the first forward operator may be expressed as shown in Formula 1.26: Equation 1.26 Wherein, W is a seismic wavelet kernel matrix (also referred to herein as a first seismic wavelet kernel matrix), which can be obtained by the seismic wavelet time shift, such as a wavelet extracted from seismic data at a first azimuth; L can be a difference matrix, which can be set according to experience or actual needs, and the embodiment of the present invention does not limit this; all parameters before WL can be referred to as first azimuth forward operator perturbation data, that is, the first azimuth forward operator perturbation data may include each first perturbation parameter v 33 、v 44 and v ρDisturbance values ​​at various incident angles and a first azimuth; optionally, the disturbance value of a disturbance parameter at one incident angle and one azimuth may be affected only by the corresponding incident angle and azimuth, that is, other parameters may be eliminated by the above formula to determine the disturbance value of the disturbance parameter, or, the value of a disturbance parameter at one incident angle and one azimuth may be set according to experience or actual needs.

[0027] Based on this, when constructing at least one first azimuth forward operator, the electronic device can determine the first azimuth forward operator perturbation data, and determine the first seismic wavelet kernel matrix; and construct at least one first azimuth forward operator based on the first azimuth forward operator perturbation data and the first seismic wavelet kernel matrix. Optionally, when determining the first azimuth forward operator perturbation data, h incident angles under the first azimuth can be determined, and the first azimuth forward operator perturbation data can be determined based on the first azimuth and the h incident angles under the first azimuth; the first azimuth forward operator perturbation data includes each first parameter (i.e., X 33 , X 44 and X ρ ) corresponding to the disturbance parameters (i.e., each first disturbance parameter, such as X 33 , X 44 and X ρ Respectively with v 33 、v 44 and v ρ Correspondingly) the disturbance values ​​at the first azimuth and h incident angles under the first azimuth, where h is a positive integer. Then correspondingly, the electronic device can use formula 1.26 to construct at least one first azimuth forward operator.

[0028] Similarly, the forward model in the plane of the symmetry axis is shown in formula 1.27: Formula 1.27 Among them, the superscript 0 =0° is the second azimuth, Seismic data with an azimuth of 0° (also called second azimuth seismic data); X 2 =[X 11 ,X 13 ,X 55 ] T It can be the vector form of the second parameter to be inverted, which is represented by the parameter C 11 , C 13 and C 55 Optionally, the number of the second forward operators in at least one second forward operator may be 2, that is, the second forward operator may be included. and Optionally, the second forward operator may be as shown in formula 1.28: Formula 1.28 Wherein, W here may also be referred to as the second seismic wavelet kernel matrix, such as a wavelet extracted from the second azimuth seismic data, or a wavelet extracted from the first azimuth seismic data, and so on; Optionally, the second seismic wavelet kernel matrix may be the same as or different from the first seismic wavelet kernel matrix, and this is not limited in the embodiment of the present invention; Optionally, the second azimuth forward operator disturbance data may include each second disturbance parameter v 11 、v 13 and v 55 The disturbance value at each incident angle and the second azimuth; optionally, the h incident angles at the first azimuth and the h incident angles at the second azimuth may be the same or different, and the embodiment of the present invention is not limited to this; for the sake of ease of explanation, the following explanation is taken as an example that the h incident angles at the first azimuth and the h incident angles at the second azimuth are the same.

[0029] Based on this, when constructing at least one second azimuth forward operator, the electronic device can determine the second azimuth forward operator perturbation data, and determine the second seismic wavelet kernel matrix; and construct at least one second azimuth forward operator based on the second azimuth forward operator perturbation data and the second seismic wavelet kernel matrix. Optionally, when determining the second azimuth forward operator perturbation data, h incident angles at the second azimuth are determined, and based on the second azimuth and the h incident angles at the second azimuth, the second azimuth forward operator perturbation data is determined; the second azimuth forward operator perturbation data includes the perturbation values ​​of the perturbation parameters corresponding to each first parameter at the h incident angles at the second azimuth and the second azimuth, and each second parameter in the second parameter to be inverted (i.e., X 11 , X 13 and X 55 ) corresponding to the disturbance parameters (i.e., each second disturbance parameter, such as X 11 , X 13 and X 55 Respectively with v 11 、v 13 and v 55 Correspondingly) the disturbance values ​​at the second azimuth and the h incident angles under the second azimuth. Then correspondingly, the electronic device can use formula 1.28 to construct at least one second azimuth forward operator.

[0030] S103, constructing a first target inversion function based on at least one first azimuth forward operator and first azimuth seismic data; and constructing a second target inversion function based on at least one second azimuth forward operator and second azimuth seismic data.

[0031] S104, determining a target inversion result of the parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data.

[0032] Optionally, when the number of sampling medium layers in at least one sampling medium layer is multiple, the target inversion result of the parameter to be inverted (i.e., X) may include the target inversion results of each parameter in the parameter to be inverted (i.e., each first parameter in the first parameter to be inverted and each second parameter in the second parameter to be inverted) under each sampling medium layer.

[0033] S105, calculating the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted.

[0034] Optionally, the target inversion result of the anisotropic parameter may include the target inversion result of the anisotropic parameter under each sampling medium layer.

[0035] The embodiment of the present invention can construct at least one first azimuth forward operator and at least one second azimuth forward operator after obtaining the initial inversion data, and obtaining the first azimuth seismic data and the second azimuth seismic data. Then, a first target inversion function can be constructed based on the at least one first azimuth forward operator and the first azimuth seismic data; and a second target inversion function can be constructed based on the at least one second azimuth forward operator and the second azimuth seismic data. Further, the target inversion result of the parameter to be inverted can be determined based on the first target inversion function, the second target inversion function and the initial inversion data; and the target inversion result of the anisotropic parameter can be calculated based on the target inversion result of the parameter to be inverted. It can be seen that the embodiment of the present invention can realize seismic inversion through the first target inversion function and the second target inversion function to predict the anisotropic parameters of the HTI medium, thereby obtaining the target inversion result of the anisotropic parameters with higher accuracy, which can effectively improve the accuracy of seismic inversion.

[0036] Based on the above description, the embodiment of the present invention further proposes a more specific seismic inversion method. Accordingly, the seismic inversion method can be executed by the electronic device (terminal or server) mentioned above; or, the seismic inversion method can be executed by the terminal and the server together. For the convenience of explanation, the following description is based on the example of an electronic device executing the seismic inversion method; please refer to Figure 2 , the seismic inversion method may include the following steps S201-S207: S201, obtaining initial inversion data, and obtaining first azimuth seismic data and second azimuth seismic data.

[0037] S202: construct at least one first azimuth forward operator, and construct at least one second azimuth forward operator.

[0038] S203, constructing a first target inversion function based on at least one first azimuth forward operator and first azimuth seismic data; and constructing a second target inversion function based on at least one second azimuth forward operator and second azimuth seismic data.

[0039] In the embodiment of the present invention, formula 1.25 and formula 1.27 show the forward model in the fracture plane and the symmetry axis plane of the HTI medium, and the forward model can predict the stiffness matrix by the step-by-step inversion method. However, the calculation of anisotropy using the stiffness matrix is ​​ill-posed, that is, a small perturbation of the stiffness matrix will seriously affect the prediction of the anisotropic parameters. Therefore, the present invention introduces a mapping relationship between the stiffness matrix and the anisotropic parameters in the inversion, and the mapping relationship can be expressed by a function as shown in formula 2.1: Formula 2.1 Among them, the function g can represent the mapping relationship between the stiffness matrix and the anisotropic parameter (i.e., the anisotropic parameter of the HTI medium). Optionally, the parameters to be inverted may include a first parameter to be inverted and a second parameter to be inverted; based on this, when constructing a first target inversion function based on at least one first azimuth forward operator and first azimuth seismic data, the electronic device can construct a first target inversion function based on at least one first azimuth forward operator, the first azimuth seismic data, and the second parameter to be inverted. Optionally, the electronic device can use formula 2.2 to construct the first target inversion function: Formula 2.2 Wherein, k represents the number of iterations (i.e., in the kth iteration, the inversion parameter data of the second parameter to be inverted at the kth iteration is used as the value of the second parameter to be inverted in the first target inversion function, so as to predict the inversion parameter data of the first parameter to be inverted at the k+1th iteration), λ 1 is the first regularization parameter, which is used to control the smoothness of the stiffness matrix and density inversion results; optionally, the first regularization parameter and parameter α can be set according to experience or actual needs, which is not limited in the embodiment of the present invention; optionally, the first objective inversion function can represent the stiffness matrix prediction in the crack plane.

[0040] Accordingly, when constructing a second target inversion function based on at least one second azimuth forward operator and second azimuth seismic data, the electronic device may construct the second target inversion function based on at least one second azimuth forward operator, the second azimuth seismic data and the first parameter to be inverted. Optionally, the electronic device may use formula 2.3 to construct the second target inversion function: Formula 2.3 Among them, λ 2is the second regularization parameter; optionally, the second regularization parameter can be set according to experience or actual needs, and the embodiment of the present invention does not limit this; optionally, the second target inversion function can represent the stiffness matrix prediction in the plane of the axis. It can be seen that the embodiment of the present invention can use the inversion parameter data of the first parameter to be inverted at the k+1th iteration as the value of the first parameter to be inverted in the second target inversion function, so as to predict the inversion parameter data of the second parameter to be inverted at the k+1th iteration. Based on this, the inversion parameter data of the first parameter to be inverted at the k+1th iteration can be inverted by formula 2.2, and then the inversion parameter data of the first parameter to be inverted at the k+1th iteration is substituted into formula 2.3 to invert the inversion parameter data of the second parameter to be inverted at the k+1th iteration, and then the inversion parameter data of the second parameter to be inverted at the k+1th iteration is substituted into formula 2.2 to invert the inversion parameter data of the first parameter to be inverted at the k+2th iteration, so as to invert the first parameter to be inverted and the second parameter to be inverted by multiple alternating inversions until convergence.

[0041] S204: Determine a target inversion result of the parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data.

[0042] In the embodiment of the present invention, a nonlinear inversion optimization algorithm for HTI media based on azimuthal amplitude difference is developed. Due to the introduction of the mapping relationship function between the stiffness matrix and the anisotropic parameters, the target functional exhibits nonlinear characteristics. The embodiment of the present invention can use the Newton method to iteratively solve the target functional. According to the Newton method, the update rules of Formula 2.2 and Formula 2.3 can be shown as Formula 2.4: Formula 2.4 Among them, the first step length μ 1 and the second step length μ 2 It can be set according to experience or according to actual needs, which is not limited in the embodiment of the present invention; that is, the values ​​of the first step length and the second step length at the kth iteration can be set according to experience or actual needs; optionally, the values ​​of the first step length and the second step length at different iterations can be the same or different, which is not limited in the embodiment of the present invention. and It can represent the first-order and second-order partial derivative operators. Optionally, the first-order partial derivative and the second-order partial derivative of the first target inversion function (i.e., the partial derivative of the first target inversion function under the first parameter to be inverted, which can also be called the partial derivative of the first target inversion function with respect to the first parameter to be inverted) can be expressed as shown in Formula 2.5: Formula 2.5 Among them, J 1 Let g be the function of vector X. 1 The Jacobian matrix, H1 Let g be the function of vector X. 1 Correspondingly, the first-order partial derivative and the second-order partial derivative of the second objective inversion function (i.e., the partial derivative of the second objective inversion function under the second parameter to be inverted, also known as the partial derivative of the second objective inversion function for the second parameter to be inverted) can be expressed as shown in formula 2.6: Formula 2.6 Among them, J 2 Let g be the function of vector X. 2 The Jacobian matrix, H 2 Let g be the function of vector X. 2 The Hessian matrix of . Optional, J 1 , J 2 , H 1 , H 2 It can be shown as formula 2.7: Formula 2.7 Based on this, the Jacobian matrix and Hessian matrix of the function can be obtained according to formula 2.7, so as to calculate the first and second derivatives of the target functional (i.e., the first target inversion function and the second target inversion function), and iteratively invert the first parameter to be inverted and the second parameter to be inverted based on the update rule.

[0043] In an embodiment of the present invention, the electronic device may determine the inversion parameter data of the first parameter to be inverted at the 0th iteration (i.e., the initial inversion parameter data of the first parameter to be inverted) and the inversion parameter data of the second parameter to be inverted at the 0th iteration (i.e., the initial inversion parameter data of the second parameter to be inverted) based on the initial inversion data. Exemplarily, when the initial inversion data include the initial values ​​of each reference parameter under each sampling medium layer, the initial stiffness matrix may be determined based on the initial inversion data (e.g., determined using formula 1.2, which may include the initial stiffness values ​​of each stiffness coefficient under each sampling medium layer), so as to determine the initial inversion parameter data of the first parameter to be inverted and the initial inversion parameter data of the second parameter to be inverted based on the initial stiffness matrix (e.g., performing logarithmic operations on the stiffness values ​​of the stiffness coefficients, etc.), a stiffness matrix may include stiffness data of each stiffness coefficient among a plurality of stiffness coefficients, and the stiffness data of a stiffness coefficient may include the stiffness values ​​of the corresponding stiffness coefficient under each sampling medium layer. Based on this, one inversion parameter data may include the inversion parameter value of the first parameter to be inverted or the second parameter to be inverted in each sampling medium layer at one iteration number.

[0044] Correspondingly, the electronic device can determine the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration based on the first target inversion function, the second target inversion function and the initial inversion data. That is to say, the initial inversion parameter data of the first parameter to be inverted and the initial inversion parameter data of the second parameter to be inverted can be determined based on the initial inversion data, so as to determine the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted based on the first target inversion function, the second target inversion function, the initial inversion parameter data of the first parameter to be inverted and the initial inversion parameter data of the second parameter to be inverted, where k is a positive integer.

[0045] Based on this, the partial derivative of the first target inversion function at the first parameter to be inverted (i.e., the partial derivative of the first target inversion function with respect to the first parameter to be inverted) can be calculated based on the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration. For example, the partial derivative of the first target inversion function at the first parameter to be inverted (which may include first-order partial derivatives and second-order partial derivatives) can be calculated through formula 2.5; and based on the partial derivative of the first target inversion function at the first parameter to be inverted, the inversion parameter data of the first parameter to be inverted at the k+1th iteration is determined. For example, the first parameter to be inverted is updated through formula 2.4 to obtain the inversion parameter data of the first parameter to be inverted at the k+1th iteration.

[0046] Correspondingly, the partial derivative of the second target inversion function at the second parameter to be inverted at the kth iteration (i.e., the partial derivative of the second target inversion function with respect to the second parameter to be inverted) can be calculated based on the inversion parameter data of the second parameter to be inverted at the kth iteration and the inversion parameter data of the first parameter to be inverted at the k+1th iteration. For example, the partial derivative of the second target inversion function at the second parameter to be inverted (which may include first-order partial derivatives and second-order partial derivatives) can be calculated through formula 2.6; and based on the partial derivative of the second target inversion function at the second parameter to be inverted, the inversion parameter data of the second parameter to be inverted at the k+1th iteration is determined. For example, the second parameter to be inverted is updated through formula 2.4 to obtain the inversion parameter data of the second parameter to be inverted at the k+1th iteration.

[0047] Further, if the inversion convergence condition is not met, k is increased by 1 so that k is equal to k+1, and k=k+1 can be executed; and the iterative execution is continued based on the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration, and the partial derivative of the first target inversion function at the first parameter to be inverted is calculated, and then the first parameter to be inverted and the second parameter to be inverted are iteratively updated until the inversion convergence condition is met; if the inversion convergence condition is met, the inversion parameter data of the first parameter to be inverted at the k+1th iteration is used as the target inversion result of the first parameter to be inverted, and the inversion parameter data of the second parameter to be inverted at the k+1th iteration is used as the target inversion result of the second parameter to be inverted, so as to achieve the target inversion result of the parameter to be inverted.

[0048] Optionally, the electronic device may determine the second parameter update difference based on the inversion parameter data of the second parameter to be inverted at the k+1th iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration, and judge whether the inversion convergence condition is met based on the second parameter update difference; Optionally, the electronic device may use formula 2.8 to judge whether the inversion convergence condition is met: Formula 2.8 The update difference threshold η may be set according to experience or actual needs, and the embodiment of the present invention does not limit this. Based on this, when the second parameter update difference is greater than the update difference threshold, it can be determined that the inversion convergence condition is not met; when the second parameter update difference is less than or equal to the update difference threshold, it can be determined that the inversion convergence condition is met.

[0049] Optionally, the electronic device may also determine the first parameter update difference based on the inversion parameter data of the first parameter to be inverted at the k+1th iteration and the inversion parameter data of the first parameter to be inverted at the kth iteration, and judge whether the inversion convergence condition is met based on the first parameter update difference, etc.; the embodiment of the present invention is not limited to this.

[0050] Optionally, the electronic device may also determine that the inversion convergence condition is met when the number of iterations reaches a preset iteration number threshold; and / or determine that the inversion convergence condition is not met when the number of iterations does not reach the preset iteration number threshold, etc.; the embodiment of the present invention does not limit this. The optional preset iteration number threshold may be set according to experience or according to actual needs, and the embodiment of the present invention does not limit this.

[0051] S205, determining a target stiffness matrix based on the target inversion result of the parameter to be inverted, where the target stiffness matrix includes stiffness data of each stiffness coefficient among the plurality of stiffness coefficients.

[0052] In an embodiment of the present invention, the electronic device can perform an exponential operation on the target inversion results of the inversion parameters to be inverted (i.e., perform an exponential operation on the target inversion results of the inversion parameters to be inverted through the exp function) to obtain a target stiffness matrix, thereby obtaining the stiffness values ​​of each stiffness coefficient in each sampling medium layer (here, the target stiffness value).

[0053] S206, based on the mapping relationship between each stiffness coefficient and a plurality of reference parameters, respectively determine a calculation method corresponding to each reference parameter in the plurality of reference parameters, the plurality of reference parameters including anisotropic parameters.

[0054] Among them, the above function g can be used to express the mapping relationship between each stiffness coefficient and multiple reference parameters. That is, the items on the right side of equation 2.1 can be used to calculate V in turn. P 、V S ,ρ,δ (v) , ε (v) , γ (v) , which can represent the mapping relationship between each stiffness coefficient and multiple reference parameters, so that the calculation method corresponding to each reference parameter in the multiple reference parameters can be determined respectively. Based on this, the electronic device can use formula 2.1 to determine the calculation method corresponding to each reference parameter in the multiple reference parameters.

[0055] S207, determining a target inversion result of anisotropic parameters according to the calculation method corresponding to each benchmark parameter and the target stiffness matrix.

[0056] In the embodiment of the present invention, the electronic device can calculate the target inversion result of each reference parameter according to the calculation method and target stiffness matrix corresponding to each reference parameter, thereby obtaining the target inversion result of the anisotropic parameter. Based on this, the embodiment of the present invention can obtain the target inversion result of the anisotropic parameter under each sampling medium layer.

[0057] Optionally, the electronic device may also determine a PP wave reflection coefficient term; and determine a target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix, and so on. Based on this, when the number of sampling medium layers in at least one sampling medium layer is multiple, for any two sampling medium layers among the multiple sampling medium layers (i.e., any two HTI medium layers, one is an upper medium and the other is a lower medium), Formula 1.21 may be used to determine the PP wave reflection coefficient term corresponding to any two sampling medium layers based on the target inversion results of each benchmark parameter under any two sampling medium layers, such as determining the inversion mean of each benchmark parameter based on the target inversion results of each benchmark parameter under any two sampling medium layers, and determining the PP wave reflection coefficient term corresponding to any two sampling medium layers based on the inversion mean of each benchmark parameter using Formula 1.21, and so on; Optionally, the electronic device may use Formula 1.20 to determine the target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix, and the difference of any stiffness coefficient (i.e., Value) can be the difference between the target stiffness value of any stiffness coefficient of the lower medium in any two sampling medium layers and the target stiffness value of any stiffness coefficient of the upper medium. The value of any stiffness coefficient substituted into formula 1.20 can be the average of the target stiffness value of any stiffness coefficient of the lower medium and the target stiffness value of any stiffness coefficient of the upper medium, and so on.

[0058] In an embodiment of the present invention, in order to further verify the feasibility of the anisotropic parameter prediction method proposed in the embodiment of the present invention, on the one hand, the accuracy of the proposed HTI medium reflection coefficient approximation is tested on a double-layer medium, where both the upper and lower media are HTI media, and the differences in various parameters are small. By comparing with the precise value of the reflection coefficient, the validity and accuracy of the proposed reflection coefficient approximation equation are tested, and the variation characteristics of the contribution of each parameter to the reflection coefficient with the azimuth are analyzed, and the dominant azimuth for the inversion of the stiffness matrix and the anisotropic parameters are selected; on the other hand, the feasibility of the proposed inversion method (i.e., the seismic inversion method) is tested on the model and actual seismic data. In the model test, noise is added to the synthetic seismic record to test the stability and noise resistance of the inversion method; in the actual seismic data test, by comparing with the conventional anisotropic seismic inversion method, the reliability and accuracy of the proposed inversion in the anisotropy prediction of fractured reservoirs are studied.

[0059] For example, when testing the accuracy of the reflection coefficient approximation equation on a double-layer medium, it is assumed that the stiffness matrix and density of the upper layer are C 11 =12.19Gpa (gigapascal), C 13 =4.88Gpa, C 33 =13.54Gpa, C 44 =4.97Gpa, C 55=3.98Gpa, ρ=2.64Gpa, anisotropy parameter is δ (v) =-0.05, ε (v) =-0.05, γ (v) =-0.1; the stiffness matrix and density of the lower medium are C 11 =15.19Gpa, C 13 =7.50Gpa, C 33 =16.88Gpa, C 44 =6.08Gpa, C 55 =4.25Gpa, ρ=2.70Gpa, anisotropy parameter is δ (v) =-0.1, ε (v) =-0.1, γ (v) =-0.15. Based on this, Figure 3 shows different orientations (i.e., azimuth, i.e., The reflection coefficient of different values ​​of ) varies with the incident angle (i.e., θ), where the black line is the exact reflection coefficient, the short dashed line is the classic Rüger approximation (a weak anisotropic approximation), and the double dashed line (i.e., the long dashed line) is the reflection coefficient proposed in the embodiment of the present invention, R PP represents the PP wave reflection coefficient; Figure 4 The variation of the reflection coefficient with azimuth at different incident angles is shown; it can be seen that the reflection coefficient proposed in the embodiment of the present invention is in better agreement with the precise value, indicating that the HTI medium reflection coefficient proposed in the present invention is effective, that is, the reflection coefficient calculated by the target inversion result obtained by the inversion method proposed in the embodiment of the present invention is effective, that is, the accuracy is higher.

[0060] For another example, the contribution of each parameter to the reflection coefficient was analyzed on the double-layer model. Figure 5-Figure 10 The stiffness matrix C is shown 11 , C 13 , C 33 , C 44 , C 55 The contribution of the density to the reflection coefficient in different directions varies with the incident angle. It can be seen that the stiffness coefficient C 11 , C 13 and C 55 In the crack plane ( =90°) has no contribution to the reflection coefficient, and the stiffness coefficient C 11 , C 13 and C 55 In the plane of the symmetry axis ( =0°) contributes more to the reflection coefficient, and so on.

[0061] For another example, the feasibility of the proposed anisotropic parameter prediction method is tested on a model. The model is a well logging curve obtained from a reservoir with vertical fractures. It is worth noting that its anisotropic parameters are calculated by a rock physics model, and the seismic records are synthesized using an accurate reflection coefficient equation. By adding a certain amount of noise, noise-free seismic data with a signal-to-noise ratio of 8 is generated. Fig.11 are elastic and anisotropic parameters predicted from noise-free seismic data; km / s represents kilometers per second, g / cm 3 It can represent grams per cubic centimeter; the gray line can be the initial model (i.e., the initial inversion data), which is obtained by low-pass filtering the true model; the black line is the true model, the dotted line is the inversion result of the conventional linear method, and the dotted line is the target inversion result predicted by the method proposed in the embodiment of the present invention. In the absence of noise, the target inversion result predicted by the method proposed in the embodiment of the present invention is consistent with the true model. Correspondingly, Fig.12 This is the inversion result with a signal-to-noise ratio of 8. Compared with the noise-free inversion result, the signal-to-noise ratio and accuracy of the anisotropic parameters are slightly reduced, but the method proposed in the embodiment of the present invention is still better than the conventional linear method, which proves that the inversion method proposed in the present invention is effective.

[0062] For another example, the embodiment of the present invention may apply the proposed inversion method to a certain actual work area. Fig.13 and Fig.14 The anisotropic parameters predicted by the inversion method proposed in the embodiment of the present invention and the conventional linear method. The high-value area of ​​the profile indicates the reservoir location where vertical fractures are developed. The black curve (i.e., the longitudinal black curve) is the corresponding anisotropic parameter curve (used to indicate that the curve on the right side (i.e., the right side of the dotted line) indicates a higher value, and the curve on the left side indicates a lower value, i.e., the more the black curve deviates to the right, the larger the parameter value indicated). Compared with the conventional inversion method, the results predicted by the method proposed in the embodiment of the present invention have advantages in resolution and accuracy. Accordingly, Fig.15 and Fig.16 is the anisotropic parameter predicted by the inversion method proposed in the embodiment of the present invention and the conventional linear method, and Fig.17 and Fig.18 The anisotropic parameters predicted by the inversion method proposed in the embodiment of the present invention and the conventional linear method. By comparing the prediction results of the conventional method and the method proposed in the embodiment of the present invention, it can be seen that the inversion method proposed in the embodiment of the present invention effectively improves the stability and accuracy of the anisotropic parameter prediction.

[0063] The embodiment of the present invention can construct at least one first azimuth forward operator and at least one second azimuth forward operator after obtaining the initial inversion data, and obtaining the first azimuth seismic data and the second azimuth seismic data. Then, a first target inversion function can be constructed based on the at least one first azimuth forward operator and the first azimuth seismic data; and a second target inversion function can be constructed based on the at least one second azimuth forward operator and the second azimuth seismic data. Based on this, the target inversion result of the parameter to be inverted can be determined based on the first target inversion function, the second target inversion function and the initial inversion data; and the target stiffness matrix can be determined based on the target inversion result of the parameter to be inverted, and the target stiffness matrix includes the stiffness data of each stiffness coefficient in a plurality of stiffness coefficients. Further, based on the mapping relationship between each stiffness coefficient and a plurality of benchmark parameters, the calculation method corresponding to each benchmark parameter in a plurality of benchmark parameters can be determined respectively, and the plurality of benchmark parameters include anisotropic parameters; and the target inversion result of the anisotropic parameter can be determined according to the calculation method corresponding to each benchmark parameter and the target stiffness matrix. It can be seen that the embodiments of the present invention can effectively improve the accuracy of the target inversion results of the anisotropy parameters, that is, can improve the reliability and accuracy of the inversion in the anisotropy prediction of the fractured reservoir.

[0064] Based on the description of the above-mentioned related embodiments of the seismic inversion method, the embodiment of the present invention further proposes a seismic inversion device, which can be a computer program (including program code) running in an electronic device; the seismic inversion device can include an acquisition unit and a processing unit. The seismic inversion device can execute Figure 1 or Figure 2 The seismic inversion method shown, that is, the seismic inversion device can run the above units: An acquisition unit, used for acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; A processing unit, configured to construct at least one first azimuth forward operator and at least one second azimuth forward operator; The processing unit is further configured to construct a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data; and to construct a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data; The processing unit is further used to determine a target inversion result of the parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data; The processing unit is further used to calculate the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted.

[0065] In one embodiment, when constructing at least one first azimuth forward operator, the processing unit may be specifically used to: determine first azimuth forward operator disturbance data, and determine a first seismic wavelet kernel matrix; construct at least one first azimuth forward operator based on the first azimuth forward operator disturbance data and the first seismic wavelet kernel matrix; When constructing at least one second azimuth forward operator, the processing unit can be specifically used to: determine the second azimuth forward operator disturbance data, and determine the second seismic wavelet kernel matrix; and construct at least one second azimuth forward operator based on the second azimuth forward operator disturbance data and the second seismic wavelet kernel matrix.

[0066] In another embodiment, when determining the first azimuth forward operator disturbance data, the processing unit can be specifically used to: determine h incident angles under a first azimuth, and determine the first azimuth forward operator disturbance data based on the first azimuth and the h incident angles under the first azimuth; the first azimuth forward operator disturbance data includes disturbance values ​​of disturbance parameters corresponding to each first parameter in the parameters to be inverted at the first azimuth and the h incident angles under the first azimuth, where h is a positive integer; When determining the second azimuth forward operator perturbation data, the processing unit can be specifically used to: determine h incident angles under the second azimuth, and determine the second azimuth forward operator perturbation data based on the second azimuth and the h incident angles under the second azimuth; the second azimuth forward operator perturbation data includes the perturbation values ​​of the perturbation parameters corresponding to the respective first parameters at the second azimuth and the h incident angles under the second azimuth, and the perturbation values ​​of the perturbation parameters corresponding to the respective second parameters in the parameters to be inverted at the second azimuth and the h incident angles under the second azimuth.

[0067] In another embodiment, the parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted; when the processing unit constructs the first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data, it can be specifically used to: construct the first target inversion function based on the at least one first azimuth forward operator, the first azimuth seismic data and the second parameter to be inverted; When the processing unit constructs the second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data, it can be specifically used to: construct the second target inversion function based on the at least one second azimuth forward operator, the second azimuth seismic data and the first parameter to be inverted.

[0068] In another embodiment, the parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted; when the processing unit determines the target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data, it can be specifically used to: determine the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration based on the first target inversion function, the second target inversion function and the initial inversion data, and calculate the partial derivative of the first target inversion function at the first parameter to be inverted based on the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration, where k is a positive integer; and determine the inversion parameter data of the first parameter to be inverted at the k+1th iteration based on the partial derivative of the first target inversion function at the first parameter to be inverted; calculate the inversion parameter data of the second parameter to be inverted at the kth iteration based on the inversion parameter data of the first parameter to be inverted at the k+1th iteration; calculate the inversion parameter data of the first parameter to be inverted at the k+1th iteration based on the inversion parameter data of the second parameter to be inverted at the kth iteration and the inversion parameter data of the first parameter to be inverted The partial derivative of the second target inversion function under the second parameter to be inverted is calculated based on the inversion parameter data under the k+1th iteration; and the inversion parameter data of the second parameter to be inverted under the k+1th iteration is determined based on the partial derivative of the second target inversion function under the second parameter to be inverted; if the inversion convergence condition is not met, k is increased by 1 so that k is equal to k+1; and the inversion parameter data of the first parameter to be inverted under the kth iteration and the inversion parameter data of the second parameter to be inverted under the kth iteration are continued to be iteratively executed to calculate the partial derivative of the first target inversion function under the first parameter to be inverted until the inversion convergence condition is met; if the inversion convergence condition is met, the inversion parameter data of the first parameter to be inverted under the k+1th iteration is used as the target inversion result of the first parameter to be inverted, and the inversion parameter data of the second parameter to be inverted under the k+1th iteration is used as the target inversion result of the second parameter to be inverted, so as to achieve the determination of the target inversion result of the parameter to be inverted.

[0069] In another embodiment, when the processing unit calculates the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted, it can be specifically used to: determine a target stiffness matrix based on the target inversion result of the parameter to be inverted, the target stiffness matrix including stiffness data of each stiffness coefficient among multiple stiffness coefficients; determine a calculation method corresponding to each benchmark parameter among the multiple benchmark parameters based on a mapping relationship between each stiffness coefficient and multiple benchmark parameters, the multiple benchmark parameters including the anisotropic parameter; determine the target inversion result of the anisotropic parameter according to the calculation method corresponding to each benchmark parameter and the target stiffness matrix.

[0070] In another embodiment, the processing unit may also be used to: determine a PP wave reflection coefficient term; and determine a target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix.

[0071] According to one embodiment of the present invention, each unit in the seismic inversion device can be separately or completely combined into one or several other units to form a structure, or one (or some) of the units can be further divided into multiple functionally smaller units to form a structure, which can achieve the same operation without affecting the realization of the technical effects of the embodiments of the present invention. The above-mentioned units are divided based on logical functions. In practical applications, the function of one unit can also be realized by multiple units, or the functions of multiple units can be realized by one unit. In other embodiments of the present invention, any seismic inversion device can also include other units. In practical applications, these functions can also be implemented with the assistance of other units, and can be implemented by the collaboration of multiple units.

[0072] According to another embodiment of the present invention, the program can be executed by running a program on a general electronic device such as a computer including a central processing unit (CPU), a random access memory medium (RAM), a read-only memory medium (ROM), and other processing elements and storage elements. Figure 1 or Figure 2 The computer program (including program code) of each step involved in the corresponding method shown in the embodiment of the present invention is used to construct a seismic inversion device and implement the seismic inversion method of the embodiment of the present invention. The computer program can be recorded on a computer storage medium, for example, and loaded into the above-mentioned electronic device through the computer storage medium and run therein.

[0073] Based on the description of the above method embodiment and device embodiment, the exemplary embodiment of the present invention further provides an electronic device, including: at least one processor; and a memory connected to the at least one processor in communication. The memory stores a computer program that can be executed by the at least one processor, and the computer program is used to enable the electronic device to perform the method according to the embodiment of the present invention when executed by the at least one processor.

[0074] Exemplary embodiments of the present invention also provide a non-transitory computer-readable storage medium storing a computer program, wherein the computer program, when executed by a processor of a computer, is used to cause the computer to perform a method according to an embodiment of the present invention.

[0075] An exemplary embodiment of the present invention further provides a computer program product, comprising a computer program, wherein when the computer program is executed by a processor of a computer, the computer is used to enable the computer to perform a method according to an embodiment of the present invention.

[0076] Furthermore, it should be understood that what is disclosed above is only a preferred embodiment of the present invention, and certainly cannot be used to limit the scope of rights of the present invention. Therefore, equivalent changes made according to the claims of the present invention are still within the scope of the present invention.

Claims

1. A seismic inversion method, characterized in that: include: Acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; constructing at least one first azimuth forward operator, and constructing at least one second azimuth forward operator; constructing a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data; and constructing a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data; Determining a target inversion result of a parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data; Based on the target inversion results of the parameters to be inverted, target inversion results of the anisotropic parameters are calculated.

2. The method according to claim 1, characterized in that The constructing of at least one first azimuth forward operator comprises: Determining first azimuth forward operator disturbance data and determining first seismic wavelet kernel matrix; constructing at least one first azimuth forward operator based on the first azimuth forward operator disturbance data and the first seismic wavelet kernel matrix; The constructing of at least one second azimuth forward operator comprises: determining second azimuth forward operator disturbance data and determining second seismic wavelet kernel matrix; At least one second azimuth forward operator is constructed based on the second azimuth forward operator disturbance data and the second seismic wavelet kernel matrix.

3. The method according to claim 2, characterized in that The determining of the first azimuth forward operator disturbance data comprises: Determine h incident angles under a first azimuth, and determine first azimuth forward operator perturbation data based on the first azimuth and the h incident angles under the first azimuth; the first azimuth forward operator perturbation data includes perturbation values ​​of perturbation parameters corresponding to each first parameter in the parameters to be inverted at the first azimuth and the h incident angles under the first azimuth, where h is a positive integer; The determining of the second azimuth forward operator disturbance data comprises: Determine h incident angles under a second azimuth, and determine second azimuth forward operator perturbation data based on the second azimuth and the h incident angles under the second azimuth; the second azimuth forward operator perturbation data include perturbation values ​​of perturbation parameters corresponding to each of the first parameters at the second azimuth and the h incident angles under the second azimuth, and perturbation values ​​of perturbation parameters corresponding to each second parameter in the parameters to be inverted at the second azimuth and the h incident angles under the second azimuth.

4. The method according to any one of claims 1 to 3, characterized in that: The parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted; and constructing a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data includes: constructing a first target inversion function based on the at least one first azimuth forward operator, the first azimuth seismic data, and the second parameter to be inverted; The constructing a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data comprises: A second target inversion function is constructed based on the at least one second azimuth forward operator, the second azimuth seismic data, and the first parameter to be inverted.

5. The method according to any one of claims 1 to 3, characterized in that: The parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted; and determining a target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data includes: Based on the first target inversion function, the second target inversion function and the initial inversion data, determining the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration, and calculating the partial derivative of the first target inversion function at the first parameter to be inverted based on the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration, where k is a positive integer; and determining the inversion parameter data of the first parameter to be inverted at the k+1th iteration based on the partial derivative of the first target inversion function at the first parameter to be inverted; calculating the partial derivative of the second target inversion function under the second parameter to be inverted based on the inversion parameter data of the second parameter to be inverted at the kth iteration and the inversion parameter data of the first parameter to be inverted at the k+1th iteration; and determining the inversion parameter data of the second parameter to be inverted at the k+1th iteration based on the partial derivative of the second target inversion function under the second parameter to be inverted; If the inversion convergence condition is not met, k is increased by 1 so that k is equal to k+1; and the inversion parameter data of the first parameter to be inverted at the kth iteration and the inversion parameter data of the second parameter to be inverted at the kth iteration are continuously iterated to calculate the partial derivative of the first target inversion function under the first parameter to be inverted until the inversion convergence condition is met; If the inversion convergence condition is met, the inversion parameter data of the first parameter to be inverted at the k+1th iteration is used as the target inversion result of the first parameter to be inverted, and the inversion parameter data of the second parameter to be inverted at the k+1th iteration is used as the target inversion result of the second parameter to be inverted, so as to determine the target inversion result of the parameter to be inverted.

6. The method according to any one of claims 1 to 3, characterized in that: The target inversion result of the anisotropic parameter is calculated based on the target inversion result of the parameter to be inverted, including: Determining a target stiffness matrix based on a target inversion result of the parameter to be inverted, wherein the target stiffness matrix includes stiffness data of each stiffness coefficient among a plurality of stiffness coefficients; Based on the mapping relationship between each stiffness coefficient and a plurality of reference parameters, respectively determining a calculation method corresponding to each reference parameter in the plurality of reference parameters, the plurality of reference parameters including the anisotropy parameter; According to the calculation method corresponding to each benchmark parameter and the target stiffness matrix, the target inversion result of the anisotropic parameter is determined.

7. The method according to claim 6, characterized in that The method further comprises: Determine the PP wave reflection coefficient term; A target PP wave reflection coefficient is determined based on the PP wave reflection coefficient term and the target stiffness matrix.

8. A seismic inversion device, characterized in that: The device comprises: An acquisition unit, used for acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; A processing unit, configured to construct at least one first azimuth forward operator and at least one second azimuth forward operator; The processing unit is further configured to construct a first target inversion function based on the at least one first azimuth forward operator and the first azimuth seismic data; and to construct a second target inversion function based on the at least one second azimuth forward operator and the second azimuth seismic data; The processing unit is further used to determine a target inversion result of the parameter to be inverted based on the first target inversion function, the second target inversion function and the initial inversion data; The processing unit is further used to calculate the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted.

9. An electronic device, characterized in that: include: processor; as well as Memory for storing programs, The program includes instructions, which, when executed by the processor, cause the processor to perform the method according to any one of claims 1 to 7.

10. A non-transitory computer-readable storage medium storing computer instructions, wherein: The computer instructions are used to make a computer execute the method according to any one of claims 1-7.

Images