Prestack AVO inversion optimization method based on DTW matching technology

Abstract

The present invention provides a pre-stack AVO inversion optimization algorithm based on DTW matching technology. Seismic traces with minimal error are selected as reference traces, a Euclidean distance matrix is ​​established, and the cumulative Euclidean distance matrix is ​​calculated to characterize the global similarity between the sequence reference trace and the trace to be corrected. The cumulative error matrix is ​​calculated to determine the optimal path for correcting the trace. This path stores the amount of stretching and compression at each point, representing the optimal time-domain stretching and compression scheme for the trace to be corrected. By traversing all traces to be corrected, an optimized pre-stack trace gather is obtained. Using this optimized trace gather for pre-stack AVO inversion improves inversion accuracy. The present invention only performs local stretching and compression on the time axis, ensuring the accuracy of subsequent AVO analysis and inversion. Furthermore, the stretching and compression amount of each trace to be corrected can be automatically calculated as the correction amount. In situations where the effective reflection signal is weak and background noise interference is strong, a control layer can be added to make the matching results more reliable.

Description

Technical Field

[0001] The present invention relates to the technical field of pre-stack AVO inversion optimization, and in particular to a pre-stack AVO inversion optimization method based on DTW matching technology. Background Art

[0002] Prestack AVO inversion is one of the most commonly used methods in seismic exploration. This method primarily exploits the fact that the reflection coefficient of seismic reflection signals at elastic reflection interfaces varies with the incident angle (offset), thereby obtaining lithologic and physical properties that reflect complex subsurface oil and gas reservoirs. The theory is based on the Zoeppritz equation, which describes the energy distribution of seismic waves at interfaces. Specifically, when a seismic wave ray is incident at a certain angle, the reflection coefficient at that point can be calculated using the incident angle and the P- and S-wave velocities and densities of the layers above and below the interface.

[0003] The conventional inversion process typically requires a given initial elastic parameter model, forward-derived seismic gathers using the Zoeppritz equation, and then, by comparing the model with actual observational data, determining model corrections to invert elastic parameters sufficient to characterize the actual seismic gathers. However, actual seismic data often contain both effective waves and noise, and coupled with the influence of certain acquisition conditions, data from certain regions, especially onshore, can exhibit poor event continuity and uneven gathers. This poses significant challenges and interference to prestack inversion. Summary of the Invention

[0004] The present invention addresses the problem of low pre-stack inversion accuracy caused by poor event continuity and uneven gathers in seismic data from certain regions. A pre-stack AVO inversion optimization algorithm based on DTW matching technology is provided. To achieve the above objectives, the present invention adopts the following technical solutions:

[0005] A pre-stack AVO inversion optimization algorithm based on DTW matching technology has the following specific steps:

[0006] The basic idea of ​​DTW is to calculate the similarity between each point in two discrete sequences, and then use dynamic programming to adaptively stretch and compress the time series according to the similarity between each point, ultimately achieving the best match between the two discrete sequences.

[0007] Suppose there is a seismic gather with a time length of t. The seismic gather with the smallest error is selected as the reference gather S0(t). The other gathers to be corrected are S n (t), where n = 1, 2...ntrace, ntrace is the number of gathers to be corrected. According to the actual situation of the gathers, the stretching and compression time window L is given, and the Euclidean distance matrix is ​​established:

[0008] E i，l =|S0(i)-Sn (i+l)| 2

[0009] Where i = 1, 2, ..., t, 1 = -L, ..., -1, 0, 1, ...L, i is the number of sampling points, 1 is the amount of tension and compression, S0(i) represents the value assigned to the i-th sampling point of the reference track, S n (i+l) represents the amplitude value of the nth seismic trace to be corrected after the current sampling point is shifted back by one sampling point. This is used to establish the Euclidean distance matrix of the model, which is also the error matrix of the model, denoted as E. Matrix E reflects the point-to-point similarity between the reference trace and the trace to be corrected. However, a cumulative Euclidean distance matrix D is still required to represent the global similarity between the sequence reference trace and the trace to be corrected. First, the cumulative initial matrix is ​​initialized so that D[0, l] = E[0, l]. The remaining cumulative error value is calculated as follows:

[0010] D[i, l]=E[i, l]+min{D[i-1, l-1], D[i-1, l+1], D[i-1, l]}

[0011] The meaning of this formula is: find the minimum D value from the three directions of the previous sampling point, add the error value E of the current sampling, calculate layer by layer until i=t-1, and complete the calculation of the entire cumulative error matrix.

[0012] By calculating the cumulative error matrix, we've actually found many paths (that is, multiple stretching and compression schemes). Now we need to find the path with the lowest total cost. Therefore, starting from the last sample point i=t-1 in the cumulative error matrix, we iterate through each l and determine the l that minimizes D[t-1, l], and store the current l value in u[t-1]. Now that the path for the last sample point has been determined, the optimal path for the previous sample point t-2 can only come from three directions: D[t-2, l+1], D[t-2, l-1], and D[t-2, l]. We select the one with the smallest value and store the current l value in u[t-2]. Using this method backwards, we ultimately obtain a sequence u that stores the stretching and compression amount for each point, which is the optimal stretching and compression scheme for the corrected path in the time domain.

[0013] The problem of low prestack inversion accuracy caused by poor event continuity and uneven gathers can be solved by using the corrected gathers for prestack inversion.

[0014] The present invention improves DTW and applies it to the matching and correction of the actual observation gather events in pre-stack gather inversion, which can effectively improve the continuity of the events, reduce the low inversion resolution caused by gather unevenness, and thus improve the resolution and accuracy of pre-stack inversion.

[0015] Compared with the prior art, the advantages of the present invention are:

[0016] Only local stretching and compression are performed on the time axis, without compromising the amplitude information of the gathers. This ensures the accuracy of subsequent AVO analysis and lays a solid foundation for subsequent pre-stack inversion. The stretching and compression amount for each trace to be corrected can be automatically calculated. In situations where the effective reflection signal is weak and background noise interference is strong, a control layer can be added to make the matching results more reliable. BRIEF DESCRIPTION OF THE DRAWINGS

[0017] Figure 1 is a flow chart of an embodiment of the present invention;

[0018] Figure 2 is a pre-stack gather to be corrected according to an embodiment of the present invention;

[0019] Figure 3 This is a corrected pre-stack gather from an embodiment of the present invention. DETAILED DESCRIPTION

[0020] In order to make the objectives, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below with reference to the accompanying drawings and examples.

[0021] like Figure 1 As shown in FIG, a prestack AVO inversion optimization method based on DTW matching technology includes the following steps:

[0022] Step 1: Input the pre-stack gather S to be processed and the initial inversion model M.

[0023] Step 2: Enter Figure 2 In the pre-stack gather shown, a seismic trace with a smaller error is selected as the reference trace S0. Usually, a small-angle gather stacking gather is selected, such as a stacking gather with an incident angle of 0 to 5°.

[0024] Step 3: Calculate the error matrix E between the track to be corrected and the reference track:

[0025] E i，l =|S0(i)-S n (i+l)| 2

[0026] Where i = 1, 2, ..., t, l = -L, ..., -1, 0, 1, ...L, i is the number of sampling points, 1 is the amount of tension and compression, S0(i) represents the value assigned to the i-th sampling point of the reference track, S n (i+1) represents the amplitude value of the nth seismic trace to be corrected after the current sampling point moves back one sampling point.

[0027] Step 4: Calculate the cumulative error matrix D between the track to be corrected and the reference track:

[0028] Let: D[0, l] = E[0, l], then the cumulative error value of the remaining part can be expressed as:

[0029] D[i, l]=E[i, l]+min{D[i-1, l-1], D[i-1, l+1], D[i-1, l]}

[0030] Step 5: Find the path u with the minimum error: Starting from the last sampling point i=t-1 in the cumulative error matrix, traverse each l and determine the l that minimizes the value of D[t-1, l], and place the current l value in u[t-1]. At this point, the path of the last sampling point has been determined. Therefore, the optimal path to the previous sampling point t-2 can only come from three directions: D[t-2, l+1], D[t-2, l-1], and D[t-2, l]. Select the one with the smallest value and place the current l value in u[t-2]. Use this method to push back forward. The final sequence u stores the stretch and compression amount of each point.

[0031] Step 6: Apply the tensile and compressive values ​​of each point calculated in step 5 to the trace to be corrected to obtain the corrected seismic trace.

[0032] Step 7: Repeat steps 3, 4, 5, and 6 to obtain the complete corrected seismic gather S new like Figure 3 shown.

[0033] Step 8: Calculate synthetic seismic data using the Zoeppritz equation using the model parameters.

[0034] Step 9: Calculate the difference between the synthetic seismic data calculated in Step 8 and the corrected observation data obtained in Step 7. Determine whether to terminate the inversion by checking whether the difference is less than the termination condition. If so, terminate the inversion immediately and output the model parameters. Otherwise, proceed to Step 10.

[0035] Step 10: Use the difference obtained in step 9 and the partial derivative of the Zoeppritz equation with respect to the model parameters, i.e., the Jacobian matrix, to calculate the model update amount. This update amount is used to update the model parameters, and steps 8 and 9 are repeated until the termination condition is met.

[0036] Those skilled in the art will appreciate that the embodiments described herein are intended to help readers understand the implementation methods of the present invention, and it should be understood that the scope of protection of the present invention is not limited to such specific descriptions and embodiments. Those skilled in the art can make various other specific variations and combinations based on the technical teachings disclosed in the present invention without departing from the essence of the present invention, and such variations and combinations are still within the scope of protection of the present invention.

Claims

1. The pre-stack AVO inversion optimization algorithm based on DTW matching technology is characterized by: The following steps are involved: The seismic gather has a time length of t. The seismic gather with the smallest error is selected as the reference gather S0(t). The other gathers to be corrected are S n (t), where n = 1, 2...ntrace, ntrace is the number of gathers to be corrected; the stretching and compression time window L is given according to the actual situation of the gathers, and the Euclidean distance matrix is ​​established: E i，l =|S0(i)-S n (i+l)| 2 Where i = 1, 2, ..., t, 1 = -L, ..., -1, 0, 1, ...L, i is the number of sampling points, 1 is the amount of tension and compression, S0(i) represents the value assigned to the i-th sampling point of the reference track, S n (i+l) represents the amplitude value of the current sampling point of the nth seismic trace to be corrected after moving back one sampling point, and the Euclidean distance matrix of the model is established based on this, which is also the error matrix of the model, denoted as E; The cumulative Euclidean distance matrix D is defined to characterize the global similarity between the sequence reference trace and the trace to be corrected. First, the cumulative initial matrix is ​​initialized so that D[0, l] = E[0, l]. The remaining cumulative error value is calculated by the following formula: D[i, l]=E[i, l]+min{D[i-1, l-1], D[i-1, l+1], D[i-1, l]} The meaning of this formula is: find the minimum D value from the three directions of the previous sampling point, add the error value E of the current sampling, calculate layer by layer until i = t-1, and complete the calculation of the entire cumulative error matrix; Starting from the last sampling point i=t-1 of the cumulative error matrix, we traverse each l and determine the l that minimizes the value of D[t-1, l], and put the value of l at that time into u[t-1]; Determine the path of the last sampling point. The optimal path for the previous sampling point t-2 comes from three directions: D[t-2, l+1], D[t-2, l-1], and D[t-2, l]. Select the one with the smallest value and put the l value at that time into u[t-2]. Push it back in sequence. The obtained sequence u stores the stretching and compression amount of each point, which is the optimal stretching and compression scheme for the correction channel in the time domain. The corrected gathers are used for prestack inversion.

Images