Nonlinearly constrained time-domain high-resolution spectrum inversion method

Abstract

The invention provides a non-linear constrained time domain high-resolution spectrum inversion method, which relates to the field of seismic data broadband processing, and comprises the following steps: constructing a wedge-shaped model library composed of odd-even reflection coefficient pairs; in a time domain, carrying out convolution on the wedge-shaped model library and the wavelets to obtain a wedge-shaped model library responded by a simple layer; forming expected inversion seismic data through the combination of different layer responses; subtraction is carried out on the seismic data and an expected result, and an objective function is designed by applying L2-Cauchy constraint; solving a non-linear problem of regularization constraint by using a fast iterative shrinkage threshold algorithm to obtain a high-resolution reflection coefficient profile; and carrying out profile convolution shaping on the high-resolution reflection coefficient to obtain high-resolution wide-spectrum inversion seismic data. The method has the beneficial effects that the high-resolution reflection coefficient profile convolution shaping broadband wavelet is performed to obtain the high-resolution wide-spectrum inversion seismic data, the low frequency and the high frequency of the seismic data are compensated, and thin-layer information is effectively identified.

Description

Technical Field

[0001] This invention relates to the field of broadband seismic data processing, and in particular to a nonlinearly constrained time-domain high-resolution spectral inversion method. Background Technology

[0002] As my country's oil and gas exploration targets gradually shift towards thin-layer and unconventional complex oil and gas reservoirs, the resolution of existing conventional seismic data is low due to the influence of thin-layer structural wave group interference and formation energy absorption factors. In particular, in the identification of thin layers, the resolution of conventional seismic data attribute extraction and interpretation techniques is significantly insufficient, which has caused certain difficulties for thin-layer seismic identification and cannot meet the fine requirements of exploration and development.

[0003] To address the issue of low resolution in seismic data, spectral inversion techniques can be used to eliminate the influence of seismic wavelets and obtain high-resolution seismic wave reflection coefficient profiles, enabling thin-layer identification. However, conventional methods often suffer from multiple solutions, low resolution, or high computational costs during data processing, failing to meet the current industrial requirements for high-precision reservoir identification. Therefore, developing higher-resolution spectral inversion techniques to more precisely characterize subsurface geological structures, thereby more accurately predicting subsurface reservoir thickness and its spatial distribution, improving the accuracy of stratigraphic interpretation, and maintaining high computational efficiency is one of the core issues in improving seismic data resolution.

[0004] Chinese patent application CN201610020556.6, entitled "A Seismic Reflection Coefficient Inversion Method Based on Total Variation Minimization Constraint," mentions a similar reflection coefficient inversion method: a trace-by-trace, time-window-by-time inversion is performed on post-stack seismic data. For a single time window, a Fourier transform is first performed on the post-stack seismic data and pre-extracted seismic wavelets within the window to obtain the frequency domain expression of the reflection coefficient. Then, a Fourier transform is performed on the time-domain reflection coefficient. Based on the parity decomposition of the reflection coefficient, the real and imaginary parts are extracted, and corresponding solution equations are constructed. The equations are solved using a minimum total variation constraint based on the traditional conjugate gradient algorithm to obtain the odd and even reflection coefficients, and the original reflection coefficients within the time window are reconstructed. All time windows are inverted sequentially to obtain the reflection coefficients of a single seismic record. The inversion is then performed sequentially on all seismic traces to obtain the reflection coefficients of each seismic record. This invention's method can effectively invert seismic reflection coefficients, which is beneficial for improving seismic data resolution and reservoir prediction accuracy. However, this application employs a channel-by-channel, time-window-by-time inversion method, requiring frequency and time domain transformations of the seismic data within each time window and the construction of corresponding solution equations. This computational load is substantial, resulting in a lengthy inversion process. Furthermore, the processing results presented in the invention document indicate that this method exhibits poor lateral continuity.

[0005] Chinese patent application CN201811031301.5, entitled "A Reflection Coefficient Inversion Method for Improving Seismic Data Resolution," mentions a similar reflection coefficient inversion method. This method is based on the L1 norm and is used to invert seismic reflection coefficients. It first introduces a low-frequency model to control the amplitude variation range of the reflection coefficients, then introduces a smoothing matrix, and compares the changes in seismic reflection coefficients before and after introducing a low-pass filter. When constructing the low-pass filter, its spectrum is first built in the frequency domain, then transformed to the time domain to obtain the filter function. Finally, the filter matrix is ​​constructed, and the linear equation after introducing the low-frequency model to control the reflection coefficients is solved to obtain the reflection coefficient array. This invention improves the bandwidth of reflection coefficient inversion by using an improved algorithm with a low-pass filter, thus improving the seismic data resolution. However, the introduction of a low-frequency model to control the amplitude variation range of the inverted coefficients may impose certain limitations on the inversion results. If the variations in the subsurface medium exceed the range of the low-frequency model, it may lead to distortion of the inversion results, especially in non-target areas. Furthermore, introducing the low-frequency model and smoothing matrix requires a more complex mathematical modeling and solution process for the inversion problem. This may increase the computational complexity and time consumption.

[0006] Chinese patent application CN202010761385.9, entitled "An Intelligent Method and System for Inverting Reflection Coefficients of Seismic Data," mentions a similar reflection coefficient inversion method: based on a velocity model and the constant density assumption, a reflection coefficient model is obtained; a wavelet is designed and convolved with the reflection coefficient model to obtain post-stack seismic data; the post-stack seismic data and corresponding reflection coefficients are padded into N channels to form labeled pseudo-2D seismic data training pairs; an inversion network is designed; the inversion network is trained using the labeled pseudo-2D seismic data training pairs; the input of the inversion network is the post-stack seismic data, and the output of the inversion network is the reflection coefficients; the post-stack seismic data to be tested is padded into N channels and input into the inversion network, outputting N channels of reflection coefficients; the N channels of reflection coefficients are weighted and superimposed to obtain the reflection coefficients of the post-stack seismic data to be tested. Compared with the prior art, this invention has good stability in intelligent reflection coefficient inversion. However, this application requires a large amount of labeled pseudo-2D seismic data for training the inversion network to obtain efficient and accurate inversion results. Therefore, this method may not be applicable or may require additional data processing and annotation work in situations where a large amount of labeled data is lacking, and the processing results show that its lateral continuity is poor.

[0007] To address the aforementioned technical issues, we have developed a time-domain high-resolution spectral inversion method with nonlinear constraints. We construct an odd-even wedge model library and use a fast iterative shrinkage threshold algorithm to solve the nonlinear problem of regularization constraints, thereby obtaining high-resolution wide-spectrum inversion seismic data and solving the above technical problems. Summary of the Invention

[0008] This invention provides a high-resolution spectral inversion method in the time domain with nonlinear constraints. It constructs an odd-even wedge model library and uses a fast iterative shrinkage threshold algorithm to solve the nonlinear problem of regularization constraints, thereby obtaining high-resolution wide-spectrum inversion seismic data. Both low and high frequencies of the seismic data are compensated, the dominant frequency band is significantly increased, deep energy is well compensated, and thin-layer information is effectively identified.

[0009] Firstly, this invention provides a high-resolution time-domain spectral inversion method with nonlinear constraints, comprising the following steps:

[0010] We study the theory of odd and even reflection coefficient pairs and construct a wedge model library composed of odd and even reflection coefficient pairs.

[0011] In the time domain, the wedge model library and wavelet convolution are used to obtain a wedge model library with simple layer response;

[0012] The expected inversion seismic data is constructed by combining responses from different layers;

[0013] The difference between the seismic data and the expected results is calculated, and the objective function is designed using L2-Cauchy constraints.

[0014] The nonlinear problem with regularization constraints is solved using a fast iterative shrinkage threshold algorithm to obtain a high-resolution reflection coefficient profile.

[0015] High-resolution spectral inversion seismic data is obtained by convolutional shaping wavelet of high-resolution reflection coefficient profile.

[0016] Furthermore, the wedge-shaped model library composed of odd and even reflection coefficient pairs is constructed using the following steps:

[0017] The total reflection coefficients are expressed separately as even-numbered reflection coefficient pairs and odd-numbered reflection coefficient pairs:

[0018]

[0019] Where n refers to the number of samples between the two reflectors; t refers to the time point where the reflective layer is located; r o This refers to the odd-numbered reflection coefficient pair component; Δt refers to the layer thickness; r e This refers to the even-numbered reflection coefficient component.

[0020] Scaling the even-numbered and odd-numbered reflection coefficient pairs and combining these two parts results in:

[0021] ar e +br0=(a+b)δ(t)+(ab)δ(t+nΔt) (2)

[0022] Here, a and b refer to the scaling parameters of the odd and even reflection coefficients with respect to the components, and their magnitudes are not necessarily equal.

[0023] The theory of even and odd reflection coefficient pairs allows us to decompose any pair of reflection coefficients into a sum of even and odd reflection coefficient pairs in a certain proportion. The even reflection coefficient pairs have the same amplitude and polarity, while the odd reflection coefficient pairs have the same amplitude and opposite polarity. These pairs can be combined using a and b to form any pair of reflection coefficients, as shown below:

[0024] r(t,n,Δt)=cδ(t)+dδ(t+nΔt)=ar e +br0 (3)

[0025] Where c = (a + b); d = (ab).

[0026] By introducing a wedge matrix, formula (3) is expressed using matrix notation, and the matrix representation is as follows:

[0027]

[0028] For simple layer models, there is only one pair of scaling parameters a. i and b i Since the value is non-zero, all other scaling parameters should be zero. This ensures that only the target layer is non-zero. Formula (4) can be rewritten using a matrix as follows:

[0029]

[0030] Where U is the even reflection coefficient wedge matrix, a is the scaling parameter sequence of the even reflection coefficient with respect to the component, V is the odd reflection coefficient wedge matrix, and b is the scaling parameter sequence of the odd reflection coefficient with respect to the component.

[0031] Furthermore, the construction steps of the wedge model library for simple layer responses are as follows:

[0032] After convolving each column of formula (5) with wavelet, each column contains wavelets and is mostly non-zero. Since the convolution introduces additional samples, its length will be slightly longer. The convolution with wavelet is:

[0033]

[0034] Representing formula (6) using a matrix is ​​as follows:

[0035]

[0036] Where s is the constructed layer response; It is the convolution of wavelet and U; It is the convolution of wavelet and V.

[0037] Furthermore, the expected inverted seismic data is constructed by combining responses from different layers. The construction steps are as follows:

[0038] Since the sampling rate is Δt, a pair of equal impulse functions with an interval of nΔt constitutes the model for each odd-even wedge reflection coefficient pair. The reflector kernel matrix of the reflection coefficient pair is constructed by shifting the reflection coefficient pair along the time axis by mΔt, where m ranges from the first data point to the total number of samples in the seismic trace. Therefore, each odd-even wedge reflection coefficient pair can be written as:

[0039]

[0040] Where m is the total number of samples from the first data point to all data points in the seismic trace; n is the number of samples between two reflectors; t is the time point where the reflecting layer is located; r o This refers to the odd-numbered reflection coefficient pair component; Δt refers to the layer thickness; r e This refers to the even-numbered reflection coefficient component.

[0041] Put these two parts together:

[0042]

[0043] Where c = (a + b); d = (ab); It is r e Shifted data through resampling; r 0m It is the shifted data of r0 through resampling.

[0044] By introducing wedge matrices, formula (9) can be represented using matrix notation, where each pair of wedge matrices contains only one pair of scaling parameters a. i and b i Since the first value is non-zero, all other scaling parameters should be zero. Therefore, the number of non-zero scaling parameters equals the number of wedge matrices. The matrix representation is as follows:

[0045]

[0046] For each wedge matrix, rearranging the matrix equation of formula (10) transforms it into:

[0047]

[0048] Formula (11) can be rewritten using a matrix method as follows:

[0049]

[0050] At this stage, U and V can be interpreted as cyclic matrices. Therefore, when the computational load is high, matrix multiplication can be implemented as convolution in the frequency domain, which can significantly improve computational efficiency.

[0051] The result obtained by wavelet convolution is:

[0052]

[0053] Then rearrange:

[0054]

[0055] in, The length of each is N+M+1 times N; a M b M Each has N lengths.

[0056] Then the matrix notation is represented as:

[0057]

[0058] in, The length is N+M+1 times N×M; m a m b The length of is N×M; the length of s is N+M+1.

[0059] Furthermore, the objective function designed using L2-Cauchy constraints is constructed through the following steps:

[0060] min(||s-Gm||2+λ|m|1) (16)

[0061] Where s is the input seismic data; G is the constructed model library; λ is the regularization parameter; and m is the result to be obtained, m = [m a m b ] T =[a0 a1 … a M b0 b1 … b M ] T .

[0062] Furthermore, the steps for solving the nonlinear problem with regularization constraints using the fast iterative shrinkage threshold algorithm are as follows:

[0063] By inverting the time-domain high-resolution spectral inversion results of the nonlinear constraints, the inversion first involves converting the iterative solution x... k y k Initialize by setting step size γ0 = 1 and t0 = 1, then iteratively calculate x. k+1 :

[0064]

[0065] Among them, t k and a k This is an intermediate variable for the iteration step size, with the aim of providing the optimal iteration step size for each iteration to achieve the best iteration efficiency.

[0066] y k =x k +a k (x k -x k-1 (18)

[0067]

[0068] Here, sgn is the sign function, which returns an integer variable indicating the sign of the parameter; max is the function to find the maximum value; and λ is the regularization parameter.

[0069] If ||x k+1 -x k If ||<ε, then output the iterative solution x. k+1 Otherwise, repeat iteration steps (17)-(19).

[0070] In a second aspect, embodiments of the present invention provide an electronic device, characterized in that it includes: a processor and a memory, wherein the processor is used to execute a program for a nonlinearly constrained time-domain high-resolution spectral inversion method to implement the nonlinearly constrained time-domain high-resolution spectral inversion method.

[0071] Fourthly, embodiments of the present invention provide a storage medium, characterized in that the storage medium stores one or more programs, which can be executed by one or more processors to implement the nonlinearly constrained time-domain high-resolution spectral inversion method.

[0072] The beneficial effects of the technical solution provided by the embodiments of the present invention are as follows: The beneficial effects of the nonlinear constraint time-domain high-resolution spectral inversion method of the present invention are that the method is based on the standard least squares inversion principle, adopts rigorous theoretical derivation, and has great accuracy; the wedge model library of odd and even reflection coefficient pairs, the wedge model library of simple layer response, the objective function of L2-Cauchy constraint, and the fast iterative shrinkage threshold algorithm involved in the theory all have a rigorous theoretical basis, and the inversion results are consistent with the characteristics of seismic data.

[0073] The invention has a clear physical meaning, and the described nonlinear constraint mechanism is readily apparent, involving broadband processing of seismic data. This invention can be used to efficiently process seismic data to obtain high-resolution reflection coefficient profiles, and can also be used to reconvolve and shape broadband wavelets to obtain high-resolution broadband seismic data. Attached Figure Description

[0074] To more clearly illustrate the technical solution of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. Obviously, the drawings listed below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0075] Figure 1 This is a flowchart illustrating the time-domain high-resolution spectral inversion method with nonlinear constraints according to the present invention.

[0076] Figure 2 The synthetic seismic signal types provided in this embodiment of the present invention are as follows: a) 30Hz Ricker wavelet; b) synthetic seismic signal with Gaussian white noise added at a signal-to-noise ratio of 30dB; c) synthetic seismic signal with Gaussian white noise added at a signal-to-noise ratio of 10dB; d) synthetic seismic signal with Gaussian white noise added at a signal-to-noise ratio of 5dB.

[0077] Figure 3 The following are the spectral inversion results of the synthetic seismic signal provided in this embodiment of the present invention: a) is the true reflection coefficient sequence; b) is the spectral inversion reflection coefficient sequence of the synthetic seismic signal with a signal-to-noise ratio of 30dB; c) is the spectral inversion reflection coefficient sequence of the synthetic seismic signal with a signal-to-noise ratio of 10dB; d) is the spectral inversion reflection coefficient sequence of the synthetic seismic signal with a signal-to-noise ratio of 5dB.

[0078] Figure 4 This is the original seismic data;

[0079] Figure 5 The results are from the high-resolution processing of the original seismic data via spectral inversion.

[0080] Figure 6 The left and right images show the time spectrum of the original seismic data, the middle image shows the time spectrum after conventional deconvolution, and the right image shows the time spectrum after high-resolution spectral inversion.

[0081] Figure 7 This is the spectrum of the original seismic data;

[0082] Figure 8 This is the spectrum after high-resolution spectral inversion. Detailed Implementation

[0083] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to embodiments. Of course, the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0084] Example 1:

[0085] like Figure 1 As shown, Figure 1 This is a flowchart of a nonlinearly constrained time-domain high-resolution spectral inversion method according to the present invention.

[0086] In step 101, a wedge-shaped model library composed of even and odd reflection coefficient pairs is constructed. The total reflection coefficient pairs are represented separately as even reflection coefficient pairs and odd reflection coefficient pairs. The theory of even and odd reflection coefficient pairs allows any reflection coefficient pair to be decomposed into a sum of even and odd reflection coefficient pairs in a certain proportion. The even reflection coefficient pairs have the same amplitude and the same polarity, while the odd reflection coefficient pairs have the same amplitude and opposite polarity. Any reflection coefficient pair can be formed by combining odd and even reflection coefficient pairs.

[0087] In step 102, a wedge model library of simple layer responses is constructed. Each column of the wedge model library, which consists of pairs of odd and even reflection coefficients, is convolved with a wavelet. Each column contains a wavelet and is mostly non-zero. Since the convolution introduces additional samples, its length will be slightly longer.

[0088] In step 103, the expected inverted seismic data is constructed by combining responses from different layers. The reflector kernel matrix of the reflection coefficient pairs is constructed by shifting the reflection coefficient pairs along the time axis, where the shift range is from the first data point to all samples in the seismic trace. By introducing wedge matrices, where only one pair of scaling parameters in every two wedge matrices has a non-zero sum, and the remaining scaling parameters should be zero, the number of non-zero scaling parameters equals the number of wedge matrices. This stage can be interpreted as a cyclic matrix. Therefore, when the computational load is large, matrix multiplication can be implemented as convolution in the frequency domain, which can significantly improve computational efficiency.

[0089] In step 104, the objective function of the L2-Cauchy constraint is constructed.

[0090] In step 105, the fast iterative shrinkage threshold algorithm is used to solve the nonlinear problem of regularization constraints, and the time-domain high-resolution spectral inversion result of the nonlinear constraints is obtained through inversion.

[0091] In step 106, high-resolution spectral inversion seismic data is obtained by convolving and shaping the high-resolution reflection coefficient profile wavelet. This invention addresses the problems of low seismic data resolution, difficulty in thin-layer identification, and low computational efficiency by constructing a wedge model library of odd-even reflection coefficient pairs and a wedge model library of layer responses, designing an objective function with L2-Cauchy constraints, and using a fast iterative shrinkage threshold algorithm to solve the least-squares problem, thus efficiently obtaining high-resolution broadband inversion seismic data.

[0092] As an optimization, in step 101, a wedge-shaped model library composed of odd and even reflection coefficient pairs is constructed;

[0093] The total reflection coefficients are expressed separately as even-numbered reflection coefficient pairs and odd-numbered reflection coefficient pairs:

[0094]

[0095] Where n refers to the number of samples between the two reflectors; t refers to the time point where the reflective layer is located; r o This refers to the odd-numbered reflection coefficient pair component; Δt refers to the layer thickness; r e This refers to the even-numbered reflection coefficient component.

[0096] Scaling the even-numbered and odd-numbered reflection coefficient pairs and combining these two parts results in:

[0097] ar e +br0=(a+b)δ(t)+(ab)δ(t+nΔt) (2)

[0098] Here, a and b refer to the scaling parameters of the odd and even reflection coefficients with respect to the components, and their magnitudes are not necessarily equal.

[0099] The theory of even and odd reflection coefficient pairs allows any pair of reflection coefficients to be decomposed into a sum of even and odd reflection coefficient pairs in a certain proportion. The even reflection coefficient pairs have the same amplitude and polarity, while the odd reflection coefficient pairs have the same amplitude and opposite polarity. This can be further explained by a... i and b i They can be combined into any pair of reflection coefficients, as shown below:

[0100] r(t,n,Δt)=cδ(t)+dδ(t+nΔt)=ar e +br0 (3)

[0101] Where c = (a + b); d = (ab).

[0102] By introducing a wedge matrix, formula (3) is expressed using matrix notation, and the matrix representation is as follows:

[0103]

[0104] For simple layer models, there is only one pair of scaling parameters a. i and b i Since the value is non-zero, all other scaling parameters should be zero. This ensures that only the target layer is non-zero. Formula (4) can be rewritten using a matrix as follows:

[0105]

[0106] Where U is the even reflection coefficient wedge matrix, a is the scaling parameter sequence of the even reflection coefficient with respect to the component, V is the odd reflection coefficient wedge matrix, and b is the scaling parameter sequence of the odd reflection coefficient with respect to the component.

[0107] As an optimization, in step 102, a wedge model library of simple layer responses is constructed, as follows;

[0108] After convolving each column of formula (5) with wavelet, each column contains wavelets and is mostly non-zero. Since the convolution introduces additional samples, its length will be slightly longer. The convolution with wavelet is:

[0109]

[0110] Representing formula (6) using a matrix is ​​as follows:

[0111]

[0112] Where s is the constructed layer response; It is the convolution of wavelet and U; It is the convolution of wavelet and V.

[0113] In step 103, the expected inversion seismic data is constructed by combining the responses of different layers;

[0114] Since the sampling rate is Δt, a pair of equal impulse functions with an interval of nΔt constitutes the model for each odd-even wedge reflection coefficient pair. The reflector kernel matrix of the reflection coefficient pair is constructed by shifting the reflection coefficient pair along the time axis by mΔt, where m ranges from the first data point to the total number of samples in the seismic trace. Therefore, each odd-even wedge reflection coefficient pair can be written as:

[0115]

[0116] Where m is the total number of samples from the first data point to all data points in the seismic trace; n is the number of samples between two reflectors; t is the time point where the reflecting layer is located; r o This refers to the odd-numbered reflection coefficient pair component; Δt refers to the layer thickness; r e This refers to the even-numbered reflection coefficient component.

[0117] Put these two parts together:

[0118]

[0119] Where c = (a + b); d = (ab); It is r e Shifted data through resampling; r 0mIt is the shifted data of r0 through resampling.

[0120] By introducing wedge matrices, formula (9) can be represented using matrix notation, where each pair of wedge matrices contains only one pair of scaling parameters a. i and b i Since the first value is non-zero, all other scaling parameters should be zero. Therefore, the number of non-zero scaling parameters equals the number of wedge matrices. The matrix representation is as follows:

[0121]

[0122] For each wedge matrix, rearranging the matrix equation of formula (10) transforms it into:

[0123]

[0124] Formula (11) can be rewritten using a matrix method as follows:

[0125]

[0126] At this stage, U and V can be interpreted as cyclic matrices. Therefore, when the computational load is high, matrix multiplication can be implemented as convolution in the frequency domain, which can significantly improve computational efficiency.

[0127] The result obtained by wavelet convolution is:

[0128]

[0129] Then rearrange:

[0130]

[0131] in, The length of each is N+M+1 times N; a M b M Each has N lengths.

[0132] Then the matrix notation is represented as:

[0133]

[0134] in, The length is N+M+1 times N×M; m a m b The length of is N×M; the length of s is N+M+1.

[0135] As an optimization, in step 104, an objective function for the L2-Cauchy constraint is constructed;

[0136] min(||s-Gm||2+λ|m|1) (16)

[0137] Where s is the input seismic data; G is the constructed model library; λ is the regularization parameter; and m is the result to be obtained, m = [m a m b ] T =[a0 a1 … a M b0 b1 … b M ] T .

[0138] As an optimization, in step 105, the fast iterative shrinkage threshold algorithm is used to solve the nonlinear problem of regularization constraints;

[0139] The high-resolution spectral inversion results in the time domain of the nonlinear constraint are obtained through inversion. The inversion first involves converting the iterative solution x... k y k Initialize by setting step size γ0 = 1 and t0 = 1, then iteratively calculate x. k+1 :

[0140]

[0141] Among them, t k and a k This is an intermediate variable for the iteration step size, with the aim of providing the optimal iteration step size for each iteration to achieve the best iteration efficiency.

[0142] y k =x k +a k (x k -x k-1 (18)

[0143]

[0144] Here, sgn is the sign function, which returns an integer variable indicating the sign of the parameter; max is the function to find the maximum value; and λ is the regularization parameter.

[0145] If ||x k+1 -x k If ||<ε, then output the iterative solution x. k+1 Otherwise, repeat iteration steps (17)-(19).

[0146] Constructing synthetic seismic signals such as Figure 2 As shown, the spectral inversion results are as follows: Figure 3 As shown;

[0147] Next, actual seismic data was tested. Figure 4 This is the original seismic data; Figure 5 The results are from the high-resolution processing of the original seismic data via spectral inversion. Figure 6The left and right images show the time spectrum of the original seismic data, the middle image shows the time spectrum after conventional deconvolution, and the right image shows the time spectrum after high-resolution spectral inversion. Figure 7 This is the spectrum of the original seismic data; Figure 8 This is the spectrum after high-resolution spectral inversion. After high-resolution spectral inversion, both the low and high frequencies of the amplitude spectrum are compensated, the dominant frequency band is significantly increased, and deep energy is well compensated, resulting in better resolution and readability of the processed seismic data. The processed seismic data can effectively identify thin-layer information and shallow structures because the low-frequency portion of the amplitude spectrum is enhanced during processing, allowing seismic waves to propagate better at depths. Furthermore, by enhancing the high-frequency portion, the resolution of the seismic data is improved, making subtle geological structures and anomalies more apparent. High-resolution spectral inversion optimizes the imaging of deep structures, and the processed seismic data can better reflect velocity changes and complex structures in the subsurface medium, thus improving the imaging quality of seismic data.

[0148] Example 2:

[0149] This invention provides an electronic device, characterized in that it includes a processor and a memory, wherein the processor is used to execute a program for a nonlinearly constrained time-domain high-resolution spectral inversion method to implement the nonlinearly constrained time-domain high-resolution spectral inversion method.

[0150] An electronic device includes at least one processor, memory, at least one network interface, and other user interfaces. The various components of the electronic device are coupled together via a bus system. It is understood that the bus system is used to enable communication between these components. In addition to a data bus, the bus system also includes a power bus, a control bus, and a status signal bus.

[0151] The user interface may include a display, keyboard, or clicking device (e.g., mouse, trackball, touchpad, or touchscreen). It is understood that the memory in this embodiment may be volatile memory or non-volatile memory, or may include both.

[0152] In this embodiment of the invention, the processor executes the method steps provided in each method embodiment by calling a program or instruction stored in the memory, specifically a program or instruction stored in an application program.

[0153] In some implementations, the memory stores elements such as executable units or data structures, or subsets thereof, or extended sets thereof: operating systems and applications.

[0154] The operating system includes various system programs, such as the framework layer, core library layer, and driver layer, used to implement various basic business functions and handle hardware-based tasks. The application programs include various applications, such as media players and browsers, used to implement various application functions. The program implementing the method of this invention can be included in the application programs.

[0155] Example 3:

[0156] This invention provides a storage medium, characterized in that the storage medium stores one or more programs, which can be executed by one or more processors to implement the nonlinearly constrained time-domain high-resolution spectral inversion method.

[0157] The method steps described in conjunction with Embodiment 1 disclosed herein can be implemented using hardware, a software module executed by a processor, or a combination of both. The software module can be located in random access memory (RAM), main memory, read-only memory (ROM), electrically programmable ROM, electrically erasable programmable ROM, registers, hard disk, removable disk, CD-ROM, or any other form of storage medium known in the art.

[0158] Of course, the above description is not intended to limit the present invention, and the present invention is not limited to the examples given above. Any changes, modifications, additions or substitutions made by those skilled in the art within the scope of the present invention should also fall within the protection scope of the present invention.

Claims

1. A nonlinearly constrained time-domain high-resolution spectral inversion method, characterized in that, include: Construct a wedge model library consisting of odd and even reflection coefficient pairs; In the time domain, the odd-even wedge model library and wavelet convolution are used to obtain a wedge model library with simple layer response; The expected inversion seismic data is constructed by combining responses from different layers; The objective function is designed by subtracting the seismic data from the expected results and applying regularization constraints. The nonlinear problem with regularization constraints is solved using a fast iterative shrinkage threshold algorithm to obtain a high-resolution reflection coefficient profile. High-resolution broadband inversion seismic data is obtained by convolving high-resolution reflection coefficient profiles with broadband wavelets.

2. The nonlinearly constrained time-domain high-resolution spectral inversion method according to claim 1, characterized in that, The wedge-shaped model library composed of odd and even reflection coefficient pairs is constructed using the following steps: The total reflection coefficients are expressed separately as even-numbered reflection coefficient pairs and odd-numbered reflection coefficient pairs: In the formula, n refers to the number of samples between the two reflectors, t refers to the time point where the reflective layer is located, and r o This refers to the odd-numbered reflection coefficient pair component, where Δt refers to the layer thickness, and r... e This refers to the even-numbered reflection coefficient for the component; Scaling is applied to even-numbered and odd-numbered reflection coefficient pairs, and these two parts are then combined for representation: ar e +br0=(a+b)δ(t)+(a-b)δ(t+nΔt) (2) Where a and b are the scaling parameters of the odd and even reflection coefficients for the components; Any pair of reflection coefficients can be decomposed into a sum of even and odd reflection coefficient pairs in a certain proportion. These pairs can then be combined using a and b to form any pair of reflection coefficients, as shown below: r(t,n,Δt)=cδ(t)+dδ(t+nΔt)=ar e +br0 (3) Where c = (a + b); d = (ab); By introducing a wedge matrix, formula (3) is expressed using matrix notation, and the matrix representation is as follows:

3. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 2, characterized in that, For simple layer models, there is only one pair of scaling parameters a. i and b i Since the value is non-zero, all other scaling parameters should be zero. Formula (4) can be rewritten using a matrix as follows: Where U is the even reflection coefficient wedge matrix, a is the scaling parameter sequence of the even reflection coefficient with respect to the component, V is the odd reflection coefficient wedge matrix, and b is the scaling parameter sequence of the odd reflection coefficient with respect to the component.

4. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 3, characterized in that, The construction steps for the simple layer response parity wedge model library are as follows: The wavelet convolution of each column of formula (5) is: Formula (6) can be rewritten using a matrix as follows: Where s is the constructed layer response; It is the convolution of wavelet and U; It is the convolution of wavelet and V.

5. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 1, characterized in that, The expected inverted seismic data is constructed by combining responses from different layers. The construction steps are as follows: Each pair of odd and even wedge-shaped reflection coefficients can be written as: Where m is the total number of samples from the first data point to all data points in the seismic trace; n is the number of samples between two reflectors; t is the time point where the reflecting layer is located; r o This refers to the odd-numbered reflection coefficient pair component; Δt refers to the layer thickness; r e This refers to the even-numbered reflection coefficient for the component; Put these two parts together: Where c = (a + b) and d = (ab), It is r e Through resampled shift data, r 0m It is the shifted data of r0 through resampling; By introducing wedge matrices, formula (9) can be represented using matrix notation, where each pair of wedge matrices contains only one pair of scaling parameters a. i and b i Since the first value is non-zero, all other scaling parameters should be zero. Therefore, the number of non-zero scaling parameters equals the number of wedge matrices. The matrix representation is as follows: For each wedge matrix, rearranging the matrix equation of formula (10) transforms it into: Formula (11) can be rewritten using a matrix method as follows: At this stage, U and V are interpreted as cyclic matrices; The result obtained by wavelet convolution is: Then rearrange: in, The length of each is N+M+1 times N; a M b M Each has N lengths; Then the matrix notation is represented as: in, The length is N+M+1 times N×M; m a m b The length of is N×M; the length of s is N+M+1.

6. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 1, characterized in that, The regularization constraint design objective function is constructed through the following steps: min(||s-Gm||2+λ|m|1) (16) Where s is the input seismic data; G is the constructed model library; λ is the regularization parameter; and m is the result to be obtained, m = [m a m b ] T =[a0 a1…a M b0 b1…b M ] T .

7. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 1, characterized in that, The nonlinear problem with regularization constraints is solved using a fast iterative shrinkage threshold algorithm to obtain a high-resolution reflection coefficient profile, including: The high-resolution spectral inversion results in the time domain of the nonlinear constraint are obtained through inversion. The inversion first involves converting the iterative solution x... k y k Initialize by setting step size γ0 = 1 and t0 = 1, then iteratively calculate x. k+1 : Among them, t k and a k This is an intermediate variable for the iteration step size, with the aim of providing the optimal iteration step size for each iteration to achieve the best iteration efficiency. y k =x k +a k (x k -x k-1 ) (18) Here, sgn is the sign function, which returns an integer variable indicating the sign of the parameter; max is the function to find the maximum value; and λ is the regularization parameter.

8. The time-domain high-resolution spectral inversion method with nonlinear constraints according to claim 7, characterized in that, If ||x k+1 -x k If ||<ε, then output the iterative solution x. k+1 Otherwise, repeat iterations (17)-(19).

9. An electronic device, characterized in that, include: A processor and a memory, the processor being used to execute a program for a nonlinearly constrained time-domain high-resolution spectral inversion method to implement the high-resolution spectral inversion processing method based on strong axis stripping as described in any one of claims 1 to 8.

10. A storage medium, characterized in that, The storage medium stores one or more programs, which can be executed by one or more processors to implement the nonlinear constraint time-domain high-resolution spectral inversion method according to any one of claims 1 to 8.

Images