Seismic reflection coefficient inversion method and system based on L0 norm constraint and medium

Abstract

The invention provides a seismic reflection coefficient inversion method and system based on L0 norm constraint and a medium, and belongs to the field of oil and gas geophysical exploration. The method comprises the following steps: step 100, acquiring post-stack seismic data; step 200, constructing an inversion equation based on L0 norm constraint; and step 300, solving an optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain a reflection coefficient. According to the method, the constraint of the L0 norm is added during inversion of the reflection coefficient, namely, the reflection coefficient is made to be 0 as far as possible, the wavelet can be reconstructed at the same time, the signal-to-noise ratio of the reflection coefficient is high, transverse continuity is good, the resolution capability of an earthquake is improved, the frequency band of the earthquake is improved, and the precision of seismic exploration is improved.

Description

Technical Field

[0001] This invention belongs to the field of oil and gas geophysical exploration, specifically relating to a seismic reflection coefficient inversion method, system, and medium based on L0 norm constraints. Background Technology

[0002] Seismic exploration resolution is a crucial issue throughout the entire process of seismic acquisition, processing, and interpretation. Traditional deconvolution methods based on convolution models are ill-suited to the practical needs of high-precision seismic exploration. To address these practical problems encountered in exploration production, researchers have drawn upon and improved upon theories and achievements from other disciplines to develop seismic nonlinear inversion methods and related technologies that can enhance seismic exploration resolution.

[0003] Tang Yuyuan et al. proposed a seismic reflection coefficient method based on L1-2 norm constraints (2013, 5th Annual Conference on Oil and Gas Geophysics). They introduced L1-2 norm sparsity regularization into seismic deconvolution and used the pDCA algorithm to solve the objective function, thus realizing the inversion of seismic reflection coefficients. Through numerical simulation data testing, the advantages of L1-2 norm were verified. It performed better when the seismic data contained noise, had good noise resistance, and could well protect the relative position and amplitude of the reflection coefficients.

[0004] Zhu Xiangyu et al. proposed a method to improve seismic data resolution based on reflection coefficient inversion (2022, Petroleum Geophysical Exploration). First, based on compressed sensing theory, they constructed an objective function under the constraint of a sparse regular operator by utilizing the sparsity of reflection coefficients and the compressibility of seismic data. Then, they inverted the reflection coefficients using a fast soft thresholding iterative method and used the low-frequency and high-frequency components in the reflection coefficient spectrum to compensate for the seismic data spectrum, thereby reconstructing the seismic data spectrum. Ultimately, this method achieves the goal of bidirectional bandwidth broadening to improve seismic data resolution.

[0005] Chinese patent publication CN115327624A discloses an inversion method and system for seismic wavelet and reflection coefficients. It assumes that the seismic wavelet is compactly supported and smooth, and that the reflection coefficients are relatively sparse, constructing a corresponding optimization problem for inverting the seismic wavelet and reflection coefficient sequences. By using alternating iterations, the joint inversion problem of seismic wavelet and reflection coefficients based on compact smoothness and relative sparsity is decomposed into a wavelet inversion subproblem and a reflection coefficient inversion subproblem, which are then solved using a near-end algorithm. Compared to existing seismic wavelet and reflection coefficient inversion methods, its advantages include easier selection of optimal parameters and better lateral continuity of the inverted reflection coefficients.

[0006] The reflection coefficients obtained by existing methods are relatively messy in the horizontal direction, with low signal-to-noise ratio and poor horizontal continuity, which is somewhat different from the actual geological conditions. Summary of the Invention

[0007] The purpose of this invention is to solve the problems existing in the prior art and provide a seismic reflection coefficient inversion method, system and medium based on L0 norm constraints, improve the accuracy of reflection coefficient inversion of post-stack seismic data, improve the lateral continuity of the profile, expand the seismic frequency band and enhance the accuracy of seismic exploration.

[0008] This invention is achieved through the following technical solution:

[0009] A first aspect of the present invention provides a method for inverting seismic reflection coefficients based on L0 norm constraints, comprising:

[0010] Step 100: Obtain post-stack seismic data;

[0011] Step 200: Construct the inversion equation based on the L0 norm constraint;

[0012] Step 300: Find the optimal solution for the constructed inversion equation based on L0 norm constraints to obtain the reflection coefficient.

[0013] A further improvement of the present invention is that:

[0014] The post-stack seismic data includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i}, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data.

[0015] A further improvement of the present invention is that:

[0016] In step 200, the inversion equation based on the L0 norm constraint is constructed as follows:

[0017]

[0018] Where β is another constant; H is a binary function H(0) = 0, and the others are 1; λ is a constant parameter, usually on the order of 0.005; wr = d + n, w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data; R is an auxiliary parameter.

[0019] A further improvement of the present invention is that:

[0020] Step 300 involves finding the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient. Specific operations include:

[0021] Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include:

[0022] Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as:

[0023]

[0024] Separate equation (4) into two parts:

[0025]

[0026] Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include:

[0027] Given the initial reflection coefficient, the optimal solution for the auxiliary parameters is obtained using the following formula:

[0028]

[0029] Where i = 1, 2, ..., N;

[0030] Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include:

[0031] Using formula (5), for r i Taking the partial derivative, we obtain the following formula:

[0032]

[0033] Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient.

[0034] A further improvement of the present invention is that:

[0035] Step 300 involves finding the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient. Specific operations also include:

[0036] Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include:

[0037] Using the reflection coefficient calculated in step 303 as input, repeat steps 302 and 303, iterate and optimize N times to obtain the final reflection coefficient.

[0038] A second aspect of the present invention provides a seismic reflection coefficient inversion system based on L0 norm constraints, comprising:

[0039] Acquisition unit, used to acquire post-stack seismic data;

[0040] The building unit is used to construct inversion equations based on L0 norm constraints;

[0041] The solver unit is used to find the optimal solution to the constructed inversion equation based on the L0 norm constraint, and obtain the reflection coefficient.

[0042] A further improvement of the present invention is that:

[0043] The post-stack seismic data acquired by the acquisition unit includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i}, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data.

[0044] A further improvement of the present invention is that:

[0045] The inversion equation based on the L0 norm constraint is constructed in the building unit as follows:

[0046]

[0047] Where β is another constant; H is a binary function H(0) = 0, and the others are 1; λ is a constant parameter, usually on the order of 0.005; wr = d + n, w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data; R is an auxiliary parameter.

[0048] A further improvement of the present invention is that:

[0049] The solving unit finds the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient, specifically by performing the following operations:

[0050] Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include:

[0051] Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as:

[0052]

[0053] Separate equation (4) into two parts:

[0054]

[0055] Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include:

[0056] Given the initial reflection coefficient, the optimal solution for the auxiliary parameters is obtained using the following formula:

[0057]

[0058] Where i = 1, 2, ..., N;

[0059] Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include:

[0060] Using formula (5), for r i Taking the partial derivative, we obtain the following formula:

[0061]

[0062] Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient;

[0063] Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include:

[0064] Using the reflection coefficient calculated in step 303 as input, repeat steps 302 and 303, iterate and optimize N times to obtain the final reflection coefficient.

[0065] A third aspect of the present invention provides a computer-readable storage medium storing at least one computer-executable program, which, when executed by the computer, causes the computer to perform the steps in the L0-norm-constrained seismic reflection coefficient inversion method.

[0066] Compared with the prior art, the beneficial effects of the present invention are:

[0067] This invention obtains seismic reflection coefficients through the inversion of post-stack seismic data. Because the L0 norm constraint is added during the reflection coefficient inversion, that is, the reflection coefficient is made as close to 0 as possible, but the wavelet can be reconstructed at the same time, resulting in a high signal-to-noise ratio and good lateral continuity of the reflection coefficient, which helps to improve the resolution of earthquakes and enhance the accuracy of seismic exploration by increasing the seismic frequency band. Attached Figure Description

[0068] Figure 1 This is a flowchart of a seismic reflection coefficient inversion method based on L0 norm constraints in an embodiment of the present invention;

[0069] Figure 2 This is a flowchart illustrating the iterative optimization process of the reflection coefficient in an embodiment of the present invention.

[0070] Figure 3 The post-stack seismic data obtained in the embodiments of the present invention;

[0071] Figure 4 This is a sequence diagram of reflection coefficients obtained by inversion using conventional methods;

[0072] Figure 5 This is a sequence diagram of reflection coefficients obtained by inversion using the method of this invention. Detailed Implementation

[0073] The present invention will now be described in further detail with reference to the accompanying drawings:

[0074]

Example 1

[0075] This invention provides a method for inverting seismic reflection coefficients based on L0 norm constraints, such as... Figure 1As shown, the main steps include:

[0076] Step 100: Obtain post-stack seismic data; wherein, the post-stack seismic data includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i}, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data;

[0077] Step 200: Construct the inversion equation based on the L0 norm constraint;

[0078] Step 300: Find the optimal solution for the constructed inversion equation based on L0 norm constraints to obtain the reflection coefficient.

[0079] This invention adds L0 norm constraints during the inversion of reflection coefficients, that is, makes the reflection coefficients as close to 0 as possible, while still being able to reconstruct the wavelet, resulting in a high signal-to-noise ratio and good lateral continuity of the reflection coefficients. This helps to improve the resolution of earthquakes and enhance the accuracy of seismic exploration by increasing the seismic frequency band.

[0080]

Example 2

[0081] Step 200 involves constructing the inversion equation based on the L0 norm constraint. The specific operations include:

[0082] Equations were established using a convolution model on the post-stack seismic data:

[0083] wr=d+n (1)

[0084] Where w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data, n is random noise, and n ~ (0.σ 2 The variance is generally on the order of 0.001. In practice, it is necessary to solve for the unknown quantity r. Generally, r is assumed to be sparse, which is consistent with actual geophysical problems. In this invention, the constraint of the L0 norm is given, that is, the reflection coefficient r is made as sparse as possible. Therefore, the solution of equation (1) is:

[0085] r = argmin[(wr - d)] T (wr-d)+λ|r|0] (2)

[0086] Where λ is a constant parameter, typically on the order of 0.005; |r|0 represents the L0 norm, which is an integer representing the number of non-zero values. Note that r is a vector containing multiple components.

[0087] The physical meaning of the entire expression is that the reflection coefficients should be as sparse as possible while still representing the original signal. Formula (2) is difficult to find an optimal solution for; it is generally transformed into the L1 norm, etc. Therefore, formula (2) is modified to obtain the inversion equation based on the L0 norm constraint:

[0088]

[0089] Where β is another constant; H is a binary function H(0) = 0, and the others are 1; λ is a constant parameter, usually on the order of 0.005; wr = d + n, w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data; R is an auxiliary parameter.

[0090]

Example 3

[0091] Step 300 involves finding the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient. Specific operations include:

[0092] Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include:

[0093] Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as:

[0094]

[0095] The optimized strategy selects a distribution splitting method, dividing equation (4) into two parts:

[0096]

[0097] Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include:

[0098] The objective function of formula (6) may have two minimum values, when R i When r = 0, the objective function value is r. i 2 In other cases, the objective function value is λ / β, and R0 i =r i Then, the minimum value among the two objective function values ​​is the optimal solution for the auxiliary parameters, as shown in formula (7):

[0099]

[0100] Input the initial reflection coefficient and substitute it into formula (7) to obtain the optimal solution of the auxiliary parameters; the number of sample points of the reflection coefficient sequence is N, i = 1, 2, ..., N;

[0101] The initial reflection coefficient can be the reflection coefficient obtained by existing methods, or it can be a set of random values ​​between (-1, 1).

[0102] Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include:

[0103] Using formula (5), for r i Taking the partial derivative, we obtain the following formula:

[0104]

[0105] Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient;

[0106] Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include:

[0107] like Figure 2 As shown, using the reflection coefficient calculated in step 303 as input, steps 302 and 303 are repeated, and the optimization is performed N times to obtain the final reflection coefficient.

[0108]

Example 4

[0109] This embodiment uses actual post-stack earthquake data as an example to perform reflection coefficient inversion using the method of the present invention, in order to further explain the technical effects of the method of the present invention.

[0110] Figure 3 The post-stack seismic data actually acquired in this embodiment of the invention were inverted using both conventional methods and the method of this invention, and the resulting reflection coefficient sequence diagram is shown below. Figure 4 and Figure 5 As shown, there are Figure 4 and Figure 5 It can be seen that, compared with the existing technology, the reflection coefficient obtained by the method of the present invention has a higher signal-to-noise ratio and higher inversion accuracy.

[0111]

Example 5

[0112] This invention provides a seismic reflection coefficient inversion system based on L0 norm constraints, comprising:

[0113] The acquisition unit is used to acquire post-stack seismic data; the post-stack seismic data includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i}, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data.

[0114] The building unit is used to construct inversion equations based on L0 norm constraints;

[0115] The inversion equation based on the L0 norm constraint is constructed as follows:

[0116]

[0117] Where β is another constant; H is a binary function H(0) = 0, and the others are 1; λ is a constant parameter, usually on the order of 0.005; wr = d + n, w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data; R is an auxiliary parameter.

[0118] The solving unit is used to find the optimal solution to the constructed inversion equation based on the L0 norm constraint, and obtain the reflection coefficient. Specifically, it performs the following operations:

[0119] Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include:

[0120] Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as:

[0121]

[0122] The optimized strategy selects a distribution splitting method, dividing equation (4) into two parts:

[0123]

[0124] Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include:

[0125] The objective function of formula (6) may have two minimum values, when R i When r = 0, the objective function value is r. i 2 In other cases, the objective function value is λ / β, and R0 i =r i Then, the minimum value among the two objective function values ​​is the optimal solution for the auxiliary parameters, as shown in formula (7):

[0126]

[0127] Input the initial reflection coefficient and substitute it into formula (7) to obtain the optimal solution of the auxiliary parameters; the number of sample points of the reflection coefficient sequence is N, i = 1, 2, ..., N;

[0128] The initial reflection coefficient can be the reflection coefficient obtained by existing methods, or it can be a set of random values ​​between (-1, 1).

[0129] Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include:

[0130] Using formula (5), for r i Taking the partial derivative, we obtain the following formula:

[0131]

[0132] Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient;

[0133] Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include:

[0134] Using the reflection coefficient calculated in step 303 as input, repeat steps 302 and 303, iterate and optimize N times to obtain the final reflection coefficient.

[0135]

Example 6

[0136] This invention provides a computer-readable storage medium storing at least one computer-executable program, which, when executed by the computer, causes the computer to perform the steps in the L0 norm-constrained seismic reflection coefficient inversion method.

[0137] Those skilled in the art will understand that all or part of the processes in the methods of the above embodiments can be implemented by a computer program instructing related hardware. The computer program can be stored in a non-volatile computer-readable storage medium, and when executed, it can include the processes of the embodiments of the above methods. Any references to memory, storage, databases, or other media used in the embodiments provided in this application can include non-volatile and / or volatile memory. Non-volatile memory can include read-only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM), or flash memory. Volatile memory can include random access memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in various forms, such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), dual data rate SDRAM (DDRSDRAM), enhanced SDRAM (ESDRAM), synchronous link DRAM (SLDRAM), Rambus direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM), etc.

[0138] The above technical solution is only one embodiment of the present invention. For those skilled in the art, based on the principles disclosed in the present invention, it is easy to make various types of improvements or modifications, and not limited to the technical solutions described in the specific embodiments of the present invention. Therefore, the foregoing description is only a preferred option and is not restrictive.

Claims

1. A method for inverting seismic reflection coefficients based on L0 norm constraints, characterized in that, include: Step 100: Obtain post-stack seismic data; Step 200: Construct the inversion equation based on the L0 norm constraint; Step 300: Find the optimal solution for the constructed inversion equation based on L0 norm constraints to obtain the reflection coefficient.

2. The method according to claim 1, characterized in that, The post-stack seismic data includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i }, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data.

3. The method according to claim 1, characterized in that, In step 200, the inversion equation based on the L0 norm constraint is constructed as follows: Where β is a constant; H is a binary function H(0) = 0, and all others are 1; λ is a constant parameter, typically on the order of 0.005; wr = d + n, where w is the seismic wavelet matrix, r is the reflection coefficient, and d is the post-stack seismic data; R is an auxiliary parameter, R = {R i } 4. The method according to claim 3, characterized in that, Step 300 involves finding the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient. Specific operations include: Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include: Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as: Separate equation (4) into two parts: Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include: Given the initial reflection coefficient, the optimal solution for the auxiliary parameters is obtained using the following formula: Where i = 1, 2, ..., N; Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include: Using formula (5), for r i Taking the partial derivative, we obtain the following formula: Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient.

5. The method according to claim 4, characterized in that, Step 300 involves finding the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient. Specific operations also include: Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include: Using the reflection coefficient calculated in step 303 as input, repeat steps 302 and 303, iterate and optimize N times to obtain the final reflection coefficient.

6. A seismic reflection coefficient inversion system based on L0 norm constraints, characterized in that, include: Acquisition unit, used to acquire post-stack seismic data; The building unit is used to construct inversion equations based on L0 norm constraints; The solver unit is used to find the optimal solution to the constructed inversion equation based on the L0 norm constraint, and obtain the reflection coefficient.

7. The system according to claim 6, characterized in that, The post-stack seismic data acquired by the acquisition unit includes multiple post-stack seismic sub-data, and the post-stack seismic record data is denoted as {d}. i }, i = 1, 2, ..., N, where N represents the number of traces in the post-stack seismic data.

8. The system according to claim 6, characterized in that, The inversion equation based on the L0 norm constraint constructed in the building unit is as follows: Where β is another constant; H is a binary function H(0) = 0, and the others are 1; λ is a constant parameter, usually on the order of 0.005; wr = d + n, w is the seismic wavelet matrix, r is the reflection coefficient, d is the post-stack seismic data; R is an auxiliary parameter.

9. The system according to claim 8, characterized in that, The solving unit finds the optimal solution for the constructed inversion equation based on the L0 norm constraint to obtain the reflection coefficient, specifically by performing the following operations: Step 301 involves splitting the constructed inversion equations based on L0 norm constraints. Specific operations include: Since each component of the reflection coefficient r can be optimized individually, the inversion equation can be written as: Separate equation (4) into two parts: Step 302: Obtain the optimal solution for the auxiliary parameters. Specific operations include: Given the initial reflection coefficient, the optimal solution for the auxiliary parameters is obtained using the following formula: Where i = 1, 2, ..., N; Step 303: Using the optimal solution of auxiliary parameters, the reflection coefficient is calculated. Specific operations include: Using formula (5), for r i Taking the partial derivative, we obtain the following formula: Substitute the optimal solution of the obtained auxiliary parameters into formula (9) to calculate the reflection coefficient; Step 304, iterative optimization, to obtain the final reflection coefficient. Specific operations include: Using the reflection coefficient calculated in step 303 as input, repeat steps 302 and 303, iterate and optimize N times to obtain the final reflection coefficient.

10. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores at least one computer-executable program, which, when executed by the computer, causes the computer to perform the steps in the seismic reflection coefficient inversion method based on L0 norm constraints as described in any one of claims 1-5.

Images