Earthquake reflection signal reconstruction method

Abstract

The invention discloses a seismic reflection signal reconstruction method, which relates to the technical field of oil exploration, and comprises the following steps: inputting a missing seismic signal; constructing a target function, solving the target function by adopting a self-adaptive accelerated sparse optimization algorithm, and transforming the three-dimensional seismic data volume containing the primary reflected waves and the multiples missing 50% into a hyperbolic sparse domain to obtain a three-dimensional sparse constraint hyperbolic Radon positive transformation model; and converting the three-dimensional sparse constraint hyperbolic Radon positive transformation model into a seismic data field, and obtaining a seismic reflection signal reconstruction result until all data volumes are reconstructed. According to the rapid three-dimensional seismic data reconstruction method, an effective signal is reconstructed by adopting a three-dimensional hyperbolic Radon transform method based on an adaptive acceleration sparse optimization algorithm, missing data in a three-dimensional seismic data volume is better reconstructed, and the calculation cost is remarkably reduced.

Description

Technical Field

[0001] This invention relates to the field of petroleum exploration technology, and in particular to a method for reconstructing seismic reflection signals. Background Technology

[0002] Current seismic acquisition methods are limited by high-density costs, resulting in actual data channel spacing far exceeding high-density standards. Furthermore, urban buildings and obstacles often cause missing channels, leading to insufficient subsequent imaging resolution. Conventional reflection signal reconstruction, based on sparse representation, low-rank completion, or five-dimensional interpolation, can partially repair missing channels, but it requires iteratively solving large matrices. The computational load increases cubically with the data volume, making it difficult to meet the needs of rapid on-site processing.

[0003] Therefore, there is an urgent need for a new method that balances accuracy and efficiency, reducing computational complexity to linear or quasi-linear levels while maintaining amplitude fidelity, in order to adapt to low-cost acquisition and real-time decision-making scenarios. Summary of the Invention

[0004] To address the aforementioned technical problems, this invention discloses a method for reconstructing seismic reflection signals. Compared to traditional seismic reflection signal reconstruction methods, this invention introduces a dynamic weight adjustment mechanism to enhance signal sparsity constraints, significantly improving inversion resolution. By combining spectral space momentum transfer technology, the algorithm's convergence rate is increased from sublinear O(1 / k) of the reweighted conjugate gradient method to O(1 / k²), and the computational cost of the Radon forward and inverse transforms in a single iteration is reduced to two iterations, improving computational efficiency by more than three times.

[0005] To achieve the above objectives, the present invention adopts the following technical solution:

[0006] A method for reconstructing seismic reflection signals includes the following steps:

[0007] Step 1. Input the missing seismic signal;

[0008] Step 2. Construct the objective function and solve the objective function using an adaptive accelerated sparse optimization algorithm. Transform the 3D seismic data volume containing 50% missing primary and secondary reflections into the hyperbolic sparse domain to obtain the 3D sparse-constrained hyperbolic Ladonite transform model.

[0009] Step 3. Convert the three-dimensional sparse-constrained hyperbolic Ladonite transform model obtained in Step 2 into the seismic data domain to obtain the seismic reflection signal reconstruction results;

[0010] Step 4. Determine if all data volumes have been reconstructed. If not, return to step 2; if completed, output the reconstructed 3D seismic data volume, and the data reconstruction is complete.

[0011] Optionally, in step 1, for the three-dimensional CMP record obtained by sequential excitation, the sparse-constrained random subsampling process is represented as follows:

[0012] (1);

[0013] (2);

[0014] (3);

[0015] In the formula, x and y are the vertical and horizontal coordinates in three-dimensional space, and t is time. For time in the new coordinate system, yes A positive transformation belongs to a positive instruction, i.e., another instruction. equal , For none A regular and complete high-density 3D CMP record with positive transformation. for The high-density 3D CMP record is regular and complete after positive transformation. N x With N y These represent the number of horizontal and vertical survey lines, respectively. for Incomplete 3D CMP records obtained by random subsampling after forward transformation. for The high-density 3D CMP record is regular and complete after positive transformation. It is a diagonal matrix used for random subsampling.

[0016] Optionally, step 2 specifically includes:

[0017] Step 21. The time-domain continuous form of the three-dimensional hyperbolic Radon transform is expressed as:

[0018] (4);

[0019] In the formula, D h (t,x,y) represents the value of the input seismic signal at time t in coordinate (x,y), where x and y are along the longitudinal and transverse survey lines, respectively; For HRT model, Describing the Diclave function, A basis function of HRT; variables , and These represent the intercept and the slowness along the longitudinal and transverse survey lines, respectively.

[0020] Step 22. Perform a three-dimensional hyperbolic Radon transform in the frequency domain, and use equation (5) to transform equation (4) to obtain equation (6):

[0021] (5);

[0022] (6);

[0023] In the formula, The intercept in the new coordinate system;

[0024] From equation (4) to equation (6) The forward transformation process is abbreviated as:

[0025] (7);

[0026] Defined according to equation (5) The inverse transformation process is as follows:

[0027] (8);

[0028] Therefore, equation (6) can be transformed into equation (4) using equation (7). The inverse transformation process is abbreviated as:

[0029] (9);

[0030] In the formula, yes The inverse transformation is a reverse order, that is, setting t to equal... ;

[0031] Step 23. By transforming the coordinates of the time axis, the three-dimensional hyperbolic Radon transform is transformed into a three-dimensional parabolic Radon transform in the new coordinate system. Applying a Fourier transform to equation (6) yields:

[0032] (10);

[0033] In the formula, i is the imaginary unit. It is a frequency domain seismic signal. For frequency, It is a frequency domain HRT model;

[0034] Equation (10) is discretized into the discrete form shown in Equation (11):

[0035] (11);

[0036] In the formula, and These are the number of seismic traces along the longitudinal and transverse directions, respectively. and These represent the number of slow sampling points in the longitudinal and transverse survey lines, respectively.

[0037] Equation (11) is based on frequency For a complex linear system with variables, the matrix form can be used as follows:

[0038] (12);

[0039] In the formula, and Both are vectors and matrices with frequency as the variable; for each frequency, and These are the frequency-domain HRT operators along the longitudinal and transverse survey lines, respectively. For HRT model; For the frequency domain three-dimensional hyperbolic Radon transform operator; matrix and They respectively have the forms shown in equation (13) and equation (14):

[0040] (13);

[0041] (14);

[0042] Step 24. Obtain the three-dimensional hybrid domain hyperbolic Radon transform according to equation (12), which is used to balance the computational efficiency in the frequency domain and the sparsity in the time domain, as shown in the following equation:

[0043] (15);

[0044] In the formula, and These are the forward and inverse Fourier transforms, To obtain a three-dimensional CMP record using a three-dimensional hyperbolic Radon transform;

[0045] set up After positive transformation, the regular and complete high-density 3D CMP record can be sparsely represented in the 3D hyperbolic sparse domain. Combining equation (2) and equation (15), we get:

[0046] (16);

[0047] Step 25. Substitute equation (16) into equation (2) and simultaneously construct the objective function based on sparse constraint theory as follows:

[0048] (17);

[0049] In the formula, Let be the objective function. It is a trade-off parameter used to adjust the weight between the L2 norm and the L1 norm;

[0050] Step 26. Solve equation (17) using an adaptive accelerated sparse optimization algorithm to obtain the sparse-constrained three-dimensional hyperbolic Lardon transform model. The algorithm iteration formula is as follows:

[0051] (18);

[0052] (19);

[0053] (20);

[0054] (twenty one);

[0055] In the formula, r represents the residual, and k represents the number of iterations. Let c be the step size, and c be the iteration step size. This represents the current solution of the algorithm in the k-th iteration. It is an auxiliary variable of the algorithm in the k-th iteration, used to accelerate the convergence of the algorithm. It is obtained by updating the momentum term. The soft threshold function is expressed as follows:

[0056] (twenty two);

[0057] In the formula, For symbolic functions, The threshold parameter has the following empirical expression:

[0058] (twenty three).

[0059] Optionally, in step 3, the three-dimensional hyperbolic Radon inverse transform operator is used with... Inverse transformation to recover the complete high-density 3D CMP record The process is expressed as follows:

[0060] (twenty four);

[0061] (25).

[0062] The beneficial effects of this invention are as follows: This invention is based on a sparse-constrained seismic reflection signal reconstruction method using a fast three-dimensional hyperbolic Radon transform. It employs an adaptive accelerated sparse optimization algorithm to solve for the sparse representation basis in the seismic reflection signal reconstruction method. Compared with the conventional reweighted conjugate gradient method, this invention introduces a dynamic weight adjustment mechanism to enhance the sparsity constraint of the signal, significantly improving the inversion resolution. Combined with spectral space momentum transfer technology, the convergence rate of the algorithm is increased from sublinear O(1 / k) of the reweighted conjugate gradient method to O(1 / k²), and the computational cost of a single iteration of Radon forward and inverse transforms is reduced to 2, improving computational efficiency by more than 3 times. This results in a more accurate reconstruction of seismic reflection signals in three-dimensional CMP records. Attached Figure Description

[0063] Figure 1 This is a schematic flowchart of a seismic reflection signal reconstruction method according to the present invention;

[0064] Figure 2 This is a 3D seismic data volume containing 50% missing primary and secondary reflections;

[0065] Figure 3 To utilize the three-dimensional hyperbolic sparse domain zero-intercept time slice obtained by this invention;

[0066] Figure 4 To obtain a zero-intercept time slice for a three-dimensional hyperbolic sparse domain using the traditional reweighted conjugate gradient method;

[0067] Figure 5 The reconstructed three-dimensional data volume using this invention;

[0068] Figure 6 This is a 3D data volume reconstructed using the traditional reweighted conjugate gradient method. Detailed Implementation

[0069] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, not all of them. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention. Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the claimed invention, but merely to represent selected embodiments of the invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0070] A method for reconstructing seismic reflection signals, such as Figure 1 As shown, it includes the following steps:

[0071] Step 1. Input the missing seismic signal; for the 3D CMP record obtained by sequential excitation, the sparse-constrained random subsampling process is represented as follows:

[0072] (1);

[0073] (2);

[0074] (3);

[0075] In the formula, x and y are the vertical and horizontal coordinates in three-dimensional space, and t is time. For time in the new coordinate system, yes Forward transform, For none A regular and complete high-density 3D CMP record with positive transformation. for The high-density 3D CMP record is regular and complete after positive transformation. N x With N y These represent the number of horizontal and vertical survey lines, for Incomplete 3D CMP records obtained by random subsampling after forward transformation. for The high-density 3D CMP record is regular and complete after positive transformation. It is a diagonal matrix used for random subsampling.

[0076] Step 2. Construct the objective function and solve it using an adaptive accelerated sparse optimization algorithm. Transform the 3D seismic data volume containing 50% missing primary and secondary reflections into a hyperbolic sparse domain to obtain a 3D sparse-constrained hyperbolic Ladonite transform model; specifically including:

[0077] Step 21. The time-domain continuous form of the three-dimensional hyperbolic Radon transform is expressed as:

[0078] (4);

[0079] In the formula, D h (t,x,y) represents the value of the input seismic signal at time t in coordinate (x,y), where x and y are along the longitudinal and transverse survey lines, respectively; For HRT model, This represents Delta (Dikra function). A basis function of HRT; variables , and These represent the intercept and the slowness along the longitudinal and transverse survey lines, respectively.

[0080] Step 22. The frequency domain three-dimensional hyperbolic Radon transform has the advantage of high computational efficiency. In order to perform the three-dimensional hyperbolic Radon transform in the frequency domain, it is necessary to eliminate the radical in equation (4) and transform equation (4) using equation (5) to obtain equation (6):

[0081] (5);

[0082] (6);

[0083] In the formula, The intercept in the new coordinate system;

[0084] From equation (4) to equation (6) The forward transformation process is abbreviated as:

[0085] (7);

[0086] Accordingly, according to the definition in equation (5) The inverse transformation process is as follows:

[0087] (8);

[0088] Therefore, equation (6) can be transformed into equation (4) using equation (7). The inverse transformation process is abbreviated as:

[0089] (9);

[0090] In the formula, yes Inverse transform;

[0091] Step 23. By transforming the coordinates of the time axis, the three-dimensional hyperbolic Radon transform is transformed into a three-dimensional parabolic Radon transform in the new coordinate system. Applying a Fourier transform to equation (6) yields:

[0092] (10);

[0093] In the formula, i is the imaginary unit. It is a frequency domain seismic signal. For frequency, It is a frequency domain HRT model;

[0094] Equation (10) is discretized into the discrete form shown in Equation (11):

[0095] (11);

[0096] In the formula, and These are the number of seismic traces along the longitudinal and transverse directions, respectively. and These represent the number of slow sampling points in the longitudinal and transverse survey lines, respectively.

[0097] Equation (11) is based on frequency For a complex linear system with variables, the matrix form can be used as follows:

[0098] (12);

[0099] In the formula, and Both are vectors and matrices with frequency as the variable; for each frequency, and These are the frequency-domain HRT operators along the longitudinal and transverse survey lines, respectively. For HRT model; For the frequency domain three-dimensional hyperbolic Radon transform operator; matrix and They respectively have the forms shown in equation (13) and equation (14):

[0100] (13);

[0101] (14);

[0102] Step 24. Obtain the three-dimensional hybrid domain hyperbolic Radon transform according to equation (12), which is used to balance the computational efficiency in the frequency domain and the sparsity in the time domain, as shown in the following equation:

[0103] (15);

[0104] In the formula, and These are the forward and inverse Fourier transforms, To obtain a three-dimensional CMP record using a three-dimensional hyperbolic Radon transform;

[0105] set up After positive transformation, the regular and complete high-density 3D CMP record can be sparsely represented in the 3D hyperbolic sparse domain. Combining equation (2) and equation (15), we get:

[0106] (16);

[0107] Step 25. Substitute equation (16) into equation (2) and simultaneously construct the objective function based on sparse constraint theory as follows:

[0108] (17);

[0109] In the formula, Let be the objective function. It is a trade-off parameter used to adjust the weighting between the L2 norm and the L1 norm. In this embodiment... =0.5;

[0110] Step 26. Solve equation (17) using an adaptive accelerated sparse optimization algorithm to obtain the sparse-constrained three-dimensional hyperbolic Lardon transform model. The algorithm iteration formula is as follows:

[0111] (18);

[0112] (19);

[0113] (20);

[0114] (twenty one);

[0115] In the formula, r represents the residual, and k represents the number of iterations. Let c be the step size, and c be the iteration step size. This represents the current solution of the algorithm in the k-th iteration. It is an auxiliary variable of the algorithm in the k-th iteration, used to accelerate the convergence of the algorithm. It is obtained by updating the momentum term. The soft threshold function is expressed as follows:

[0116] (twenty two);

[0117] In the formula, For symbolic functions, The threshold parameter has the following empirical expression:

[0118] (twenty three).

[0119] Step 3. Convert the three-dimensional sparse-constrained hyperbolic Radon forward transform model obtained in Step 2 into the seismic data domain to obtain the seismic reflection signal reconstruction results; use the three-dimensional hyperbolic Radon inverse transform operator and... Inverse transformation to recover the complete high-density 3D CMP record The process is expressed as follows:

[0120] (twenty four);

[0121] (25).

[0122] Step 4. Determine if all data volumes have been reconstructed. If not, return to step 2; if completed, output the reconstructed 3D seismic data volume, and the data reconstruction is complete.

[0123] The present invention provides a rapid three-dimensional seismic data reconstruction method, which, when applied to three-dimensional seismic data volumes, achieves ideal computational results.

[0124] Figure 2 The three-dimensional seismic data volume containing 50% missing primary and secondary reflections shows a severe lack of effective signals. Figure 3 This is a zero-intercept time slice obtained using the 3D seismic data reconstruction optimization algorithm of this invention; for comparison, Figure 4The diagram shows a zero-intercept time slice obtained by the conventional reweighted conjugate gradient method. A comparison reveals that the method of this invention has higher resolution than the conventional method. Therefore, the 3D seismic data reconstruction method of this application has higher resolution in the slow-motion direction. Figure 5 The image shows the reconstructed 3D data volume obtained using the 3D seismic data reconstruction optimization algorithm of this invention; for comparison, Figure 6 The reconstructed 3D data volume obtained using the traditional reweighted conjugate gradient method is shown. By comparison, it is demonstrated that the method of the present invention achieves better results in terms of reconstruction effect and computational efficiency.

[0125] 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 method for reconstructing seismic reflection signals, characterized in that, Includes the following steps: Step 1. Input the missing seismic signal; Step 2. Construct the objective function and solve the objective function using an adaptive accelerated sparse optimization algorithm. Transform the 3D seismic data volume containing 50% missing primary and secondary reflections into the hyperbolic sparse domain to obtain the 3D sparse-constrained hyperbolic Ladonite transform model. Step 3. Convert the three-dimensional sparse-constrained hyperbolic Ladonite transform model obtained in Step 2 into the seismic data domain to obtain the seismic reflection signal reconstruction results; Step 4. Determine if all data volumes have been reconstructed. If not, return to step 2; if completed, output the reconstructed 3D seismic data volume, and the data reconstruction is complete.

2. The seismic reflection signal reconstruction method as described in claim 1, characterized in that, In step 1, for the three-dimensional CMP record obtained by sequential excitation, the sparse-constrained random subsampling process is represented as follows: （1）； （2）； （3）； In the formula, x and y are the vertical and horizontal coordinates in three-dimensional space, and t is time. For time in the new coordinate system, yes Forward transform, For none Regular and complete high-density 3D CMP records with positive transformation. for The high-density 3D CMP record is regular and complete after positive transformation. N x With N y These represent the number of horizontal and vertical survey lines, for Incomplete 3D CMP records obtained by random subsampling after forward transformation. for The high-density 3D CMP record is regular and complete after positive transformation. It is a diagonal matrix used for random subsampling.

3. The seismic reflection signal reconstruction method as described in claim 2, characterized in that, Step 2 specifically includes: Step 21. The time-domain continuous form of the three-dimensional hyperbolic Radon transform is expressed as: （4）； In the formula, D h (t,x,y) represents the value of the input seismic signal at time t in coordinate (x,y), where x and y are along the longitudinal and transverse survey lines, respectively; For HRT model, Describing the Diclave function, A basis function of HRT; variables , and These represent the intercept and the slowness along the longitudinal and transverse survey lines, respectively. Step 22. Perform a three-dimensional hyperbolic Radon transform in the frequency domain, and use equation (5) to transform equation (4) to obtain equation (6): （5）； （6）； In the formula, The intercept in the new coordinate system; From equation (4) to equation (6) The forward transformation process is abbreviated as: （7）； Defined according to equation (5) The inverse transformation process is as follows: （8）； Therefore, equation (6) can be transformed into equation (4) using equation (7). The inverse transformation process is abbreviated as: （9）； In the formula, yes Inverse transform; Step 23. By transforming the coordinates of the time axis, the three-dimensional hyperbolic Radon transform is transformed into a three-dimensional parabolic Radon transform in the new coordinate system. Applying a Fourier transform to equation (6) yields: （10）； In the formula, i is the imaginary unit. It is a frequency domain seismic signal. For frequency, It is a frequency domain HRT model; Equation (10) is discretized into the discrete form shown in Equation (11): （11）； In the formula, and These are the number of seismic traces along the longitudinal and transverse directions, respectively. and These represent the number of slow sampling points in the longitudinal and transverse survey lines, respectively. Equation (11) is based on frequency For a complex linear system with variables, the matrix form can be used as follows: （12）； In the formula, and Both are vectors and matrices with frequency as the variable; for each frequency, and These are the frequency-domain HRT operators along the longitudinal and transverse survey lines, respectively. For HRT model; For the frequency domain three-dimensional hyperbolic Radon transform operator; matrix and They respectively have the forms shown in equation (13) and equation (14): （13）； （14）； Step 24. Obtain the three-dimensional hybrid domain hyperbolic Radon transform according to equation (12), which is used to balance the computational efficiency in the frequency domain and the sparsity in the time domain, as shown in the following equation: （15）； In the formula, and These are the forward and inverse Fourier transforms, To obtain a three-dimensional CMP record using a three-dimensional hyperbolic Radon transform; set up After positive transformation, the regular and complete high-density 3D CMP record can be sparsely represented in the 3D hyperbolic sparse domain. Combining equation (2) and equation (15), we get: （16）； Step 25. Substitute equation (16) into equation (2) and simultaneously construct the objective function based on sparse constraint theory as follows: （17）； In the formula, Let be the objective function. It is a trade-off parameter used to adjust the weight between the L2 norm and the L1 norm; Step 26. Solve equation (17) using an adaptive accelerated sparse optimization algorithm to obtain the sparse-constrained three-dimensional hyperbolic Lardon transform model. The algorithm iteration formula is as follows: （18）； （19）； （20）； （21）； In the formula, r represents the residual, and k represents the number of iterations. Let c be the step size, and c be the iteration step size. This represents the current solution of the algorithm in the k-th iteration. It is an auxiliary variable of the algorithm in the k-th iteration. The soft threshold function is expressed as follows: （22）； In the formula, For symbolic functions, The threshold parameter has the following empirical expression: （23）。 4. The seismic reflection signal reconstruction method as described in claim 3, characterized in that, In step 3, the three-dimensional hyperbolic Radon inverse transform operator is used with... Inverse transformation to recover the complete high-density 3D CMP record The process is expressed as follows: （24）； （25）。

Images