A seismic exploration data regularization processing method and device
By using a dip weighting function and least squares optimization in the frequency-wavenumber domain, seismic data is regularized, which solves the reconstruction artifacts and pseudo-Gibbs phenomenon in ALFT-type algorithms when dealing with irregular or missing data, and improves the imaging quality of seismic data.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- CHINA NAT PETROLEUM CORP
- Filing Date
- 2023-12-22
- Publication Date
- 2026-06-12
Smart Images
Figure CN120195732B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of seismic exploration data processing technology, and in particular to a method and apparatus for regularizing seismic exploration data processing. Background Technology
[0002] In current seismic exploration, the acquired data often fails to meet the requirements of subsequent processing and imaging for spatial regularity and spatial sampling density. Sparse or irregular spatial sampling attributes can affect the inversion process and migration imaging results of seismic data. Seismic data regularization processing techniques can improve the spatial sampling attributes of seismic data observation systems to a certain extent, thereby improving the inversion and imaging results of seismic data.
[0003] Currently, commonly used seismic data regularization methods include anti-leakage Fourier transform (ALFT) algorithms. The ALFT algorithm is a high-dimensional seismic data regularization algorithm that has better amplitude fidelity and can fill in data gaps to a certain extent. Due to its good application effect and ease of high-dimensional implementation, it has become the mainstream regularization method in the industry. Summary of the Invention
[0004] The inventors of this application have discovered that while traditional ALFT-type algorithms can fill in data gaps to some extent, they are generally only applicable to the regularization of uniformly and randomly sampled seismic data. When seismic data contains regularly missing data with large gaps, their processing effect is not ideal. On the one hand, the influence of periodic aliasing may produce reconstruction artifacts, resulting in inaccurate data reconstruction results; on the other hand, pseudo-Gibbs phenomena may occur at discontinuities in the phase axis. This leads to poor seismic data processing results, affecting the subsequent imaging quality and effect, and causing distortion in the regularized seismic data and imaging.
[0005] In view of the above problems, the present invention is proposed to provide a method and apparatus for regularizing seismic exploration data in order to overcome or at least partially solve the above problems.
[0006] This invention provides a method for regularizing seismic exploration data, comprising:
[0007] Based on the acquired time-space domain seismic data, the wavenumber spectrum data volume of the seismic data is obtained by conversion;
[0008] The following regularization process is performed on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume:
[0009] The maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice is determined based on the tilt angle weighting function; the tilt angle weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions.
[0010] Based on the Gram matrix, the wavenumber component corresponding to the maximum amplitude is optimized using least squares. Based on the least squares optimization data and the Gram matrix, the wavenumber spectrum of the current frequency slice is updated. The Gram matrix is pre-established based on the wavenumber spectrum.
[0011] Return to the step of determining the maximum amplitude until the preset number of iterations is reached, and obtain the wavenumber spectrum after the current frequency slice regularization process;
[0012] Based on the regularized wavenumber spectrum of all frequency slices, the regularized wavenumber spectrum data volume is obtained;
[0013] The regularized wavenumber spectrum data volume is transformed to obtain regularized time-space domain seismic data.
[0014] In some optional embodiments, the step of converting the acquired time-space domain seismic data into wavenumber spectrum data of the seismic data includes:
[0015] The acquired time-space domain seismic data is subjected to a Fourier transform along the time direction to obtain the frequency-space domain data volume;
[0016] A non-uniform Fourier transform is performed on the frequency slices in the frequency-space domain data volume along the spatial direction to obtain the wavenumber spectrum of the frequency slices, and a wavenumber spectrum data volume is generated based on the wavenumber spectra of all frequency slices.
[0017] In some optional embodiments, performing a non-uniform Fourier transform along the spatial direction on the frequency slices in the frequency-spatial domain data volume to obtain the wavenumber spectrum of the frequency slices includes:
[0018] For frequency slices f(x) in the frequency-spatial domain data volume m The wavenumber spectrum F(k) of the frequency slice is obtained by performing a non-uniform Fourier transform along the spatial direction using the following formula. n ):
[0019]
[0020] Where i is the virtual unit, M is the number of samples for spatial point x, and N is the number of samples for wave number k;
[0021] The superscripts 1, 2, 3, and 4 represent four different spatial dimensions, and the subscript m represents the sample number of spatial point x.
[0022] The superscripts 1, 2, 3, and 4 represent four different wavenumber dimensions, and the subscript n represents the sample number of wavenumber k.
[0023] In some optional embodiments, determining the maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice according to the tilt angle weighting function includes:
[0024] Obtain the amplitude spectrum of the wavenumber spectrum of the current frequency slice, and search for the position of the maximum amplitude based on the amplitude spectrum and the tilt weight function to obtain the wavenumber position corresponding to the maximum amplitude in the amplitude spectrum.
[0025] In some optional embodiments, the search for the location of maximum amplitude is performed based on the amplitude spectrum and the tilt weight function, including:
[0026] The maximum amplitude location is searched using the following formula:
[0027] |w(f,p j C j (p j )|≥|w(f,n)C j (n)|, 1≤n≤N
[0028] Wherein, the tilt angle weighting function is w(f,k)=(S -1 TS|D(f,k)|) q In the above formula, k in the tilt angle weighting function on the left and right sides are respectively p j and n; p j And n represents the wave number position; C j (p j ) indicates the wavenumber position p j wavenumber component at C; j (n) represents the wavenumber component at wavenumber position n;
[0029] In the tilt weighting function, D(f,k) is the wavenumber spectrum data volume, T is the smoothing operator, S is the amplitude energy operator, representing the sum of amplitude energy along the specified ray direction in the wavenumber spectrum data volume, and the exponent q represents the sharpness of the tilt weighting function: S -1 Let S be the inverse operator; where the expression for the amplitude energy operator S is as follows:
[0030] S(θ 1 ,θ 2 ,θ 3 ,θ 4 )=∫ r |D(θ 1 ,θ 2 ,θ 3 ,θ 4 ,r)|dr
[0031]
[0032] (θ 1 ,θ2 ,θ 3 ,θ 4 (f, k) represents the Cartesian coordinates in the frequency-wavenumber domain. 1 ,k 2 ,k 3 ,k 4 The corresponding polar coordinates, r is the radial length, θ j The angle representing the dimension corresponding to the j-th wavenumber of wavenumber k; k j This represents the dimension corresponding to the j-th wavenumber of wavenumber k.
[0033] In some optional embodiments, the wavenumber component corresponding to the maximum amplitude is optimized using least squares based on the Gram matrix, and the wavenumber spectrum of the current frequency slice is updated based on the least squares optimized data and the Gram matrix, including:
[0034] Based on the extracted Gram submatrix and the wavenumber spectrum subvector extracted from the wavenumber spectrum, the wavenumber component corresponding to the maximum amplitude is optimized by least squares to obtain least squares optimized data; the Gram submatrix is obtained by extracting specified elements from the Gram matrix according to the number of iterations of the current frequency slice;
[0035] The wavenumber spectrum of the current frequency slice is updated based on the least squares optimization data and the extracted Gram sub-matrix; the Gram sub-matrix is obtained by extracting the corresponding columns from the Gram matrix according to the number of iterations.
[0036] In some optional embodiments, the wavenumber component corresponding to the maximum amplitude is optimized by least squares based on the extracted Gram subarray and the wavenumber spectral subvector extracted from the wavenumber spectrum, to obtain least squares optimized data, including:
[0037] The preset wavenumber component matrix V0 is updated based on the wavenumber component data at the wavenumber position corresponding to the maximum amplitude.
[0038] Based on the extracted Gram submatrix A and wavenumber spectral vector B, the wavenumber component matrix V0 is processed using the least squares method to obtain the least squares optimized matrix V. Correspondingly, based on the least squares optimized data and the extracted Gram submatrix, the wavenumber spectrum of the current frequency slice is updated, including: calculating the product of the extracted Gram submatrix E and the least squares optimized matrix V, and updating the wavenumber spectrum of the current frequency slice based on the difference between the wavenumber spectrum of the current frequency slice and this product.
[0039] In some optional embodiments, based on the extracted Gram submatrix A and wavenumber spectral vector B, the wavenumber component matrix V0 is processed using the least squares method to obtain the least squares optimized matrix V, including:
[0040] The least squares optimization matrix V is obtained by using the following least squares equation:
[0041] (A+0.01I)V=(B+0.01V0)
[0042] Where I is the identity matrix;
[0043] A is according to (p) l ,p r The position p extracted from the Gram matrix j A square matrix of order l, ...,j, r = 1, ...,j
[0044] B is pressed p. l The position is extracted from the wavenumber spectrum F(k) p j A dimensional vector, l = 1, ..., j;
[0045] j represents the number of iterations;
[0046] Accordingly, the product of the extracted Gram submatrix E and the least squares optimization matrix V is calculated, including:
[0047] Get N rows of p j The gram submatrix E is arranged according to p. l The position is extracted from the Gram matrix p j The matrix formed by the columns, l = 1, ..., j, where j represents the iteration number;
[0048] Calculate the product of the Gram submatrix E and the least squares optimization matrix V.
[0049] In some optional embodiments, the regularized wavenumber spectrum data volume of the seismic data is transformed to obtain regularized time-space domain seismic data, including:
[0050] A non-uniform inverse Fourier transform is performed on the wavenumber spectrum of each frequency slice in the regularized wavenumber spectrum data volume of the seismic data to obtain the regularized frequency-spatial domain data volume.
[0051] An inverse Fourier transform is performed on the frequency-space domain data volume after regularization to obtain regularized time-space domain seismic data.
[0052] This invention also provides a seismic exploration data regularization processing device, comprising:
[0053] The first data acquisition module is used to convert the acquired time-space domain seismic data into wavenumber spectrum data volume of the seismic data.
[0054] The regularization processing module is used to perform the following regularization processing on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume:
[0055] The maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice is determined based on the tilt angle weighting function; the tilt angle weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions.
[0056] Based on the Gram matrix, the wavenumber component corresponding to the maximum amplitude is optimized using least squares. Based on the least squares optimization data and the Gram matrix, the wavenumber spectrum of the current frequency slice is updated. The Gram matrix is pre-established based on the wavenumber spectrum.
[0057] Return to the step with the maximum amplitude until the preset number of iterations is reached, and obtain the wavenumber spectrum after the current frequency slice is regularized.
[0058] Based on the regularized wavenumber spectrum of all frequency slices, the regularized wavenumber spectrum data volume is obtained;
[0059] The second data acquisition module is used to convert the regularized wavenumber spectrum data volume to obtain regularized time-space domain seismic data.
[0060] This invention also provides a computer storage medium storing computer-executable instructions, which, when executed by a processor, implement the above-described seismic exploration data regularization processing method.
[0061] This invention also provides a terminal device, including: a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the program to implement the above-described seismic exploration data regularization processing method.
[0062] The beneficial effects of the above-described technical solutions provided in the embodiments of the present invention include at least the following:
[0063] The frequency-wavenumber domain seismic data regularization method provided in this invention performs regularization processing on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume of seismic data converted to the frequency-wavenumber domain. During the processing, it improves the regularization of seismic data with large gaps and missing regularities. A dip weight function is generated during dip scanning in the frequency-wavenumber domain and used to constrain the energy picking process in each iteration, enhancing the regularization algorithm's resistance to spatial spurious frequencies and effectively avoiding potential reconstruction artifacts and inaccurate data reconstruction results. After selecting the wavenumber component with the highest energy at each step, a least-squares optimization is used to update all selected wavenumber components, improving the pseudo-Gibbs phenomenon at discontinuities on the phase axis. Data regularization improves the spatial sampling attributes of seismic exploration data, thereby enhancing the imaging quality of seismic exploration data.
[0064] Other features and advantages of the invention will be set forth in the description which follows, and will be apparent in part from the description, or may be learned by practicing the invention. The objects and other advantages of the invention may be realized and obtained by means of the structures particularly pointed out in the written description, claims, and drawings.
[0065] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description
[0066] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used in conjunction with embodiments of the invention to explain the invention and do not constitute a limitation thereof. In the drawings:
[0067] Figure 1 This is a flowchart of the seismic exploration data regularization processing method in an embodiment of the present invention;
[0068] Figure 2 This is a flowchart of a frequency slice regularization process in an embodiment of the present invention;
[0069] Figure 3 This is an example diagram of seismic data in an embodiment of the present invention;
[0070] Figure 4 This is an example diagram of the data regularization results in an embodiment of the present invention when there are no tilt angle weight constraints and no least squares optimization.
[0071] Figure 5 This is an example diagram of the data regularization results in an embodiment of the present invention when there are tilt angle weight constraints and no least squares optimization;
[0072] Figure 6 This is an example diagram of the data regularization result after introducing least squares optimization with tilt angle weight constraints in an embodiment of the present invention;
[0073] Figure 7 This is a schematic diagram of the structure of the seismic exploration data regularization processing device in an embodiment of the present invention. Detailed Implementation
[0074] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.
[0075] To address the issues in existing technologies where traditional ALFT-type algorithms, when processing seismic data with significant gaps and irregularities, may produce reconstruction artifacts due to periodic aliasing and pseudo-Gibbs phenomena at discontinuities in the same phase axis, resulting in suboptimal seismic data processing and impacting subsequent imaging quality, this invention provides a frequency-wavenumber domain seismic data regularization method. This method generates a dip weight function through dip scanning in the frequency-wavenumber domain, using this function to constrain the energy picking process during each iteration, thus enhancing the regularization algorithm's resistance to spatial aliasing. After selecting the wavenumber component with the highest energy at each iteration, a least-squares optimization is used to update all selected wavenumber components, improving the pseudo-Gibbs phenomenon at discontinuities in the same phase axis, thereby improving the spatial sampling properties of the seismic data and enhancing its imaging quality.
[0076] This invention provides a method for regularizing frequency-wavenumber domain seismic data, the process of which is as follows: Figure 1 As shown, it includes the following steps:
[0077] Step S11: Based on the acquired time-space domain seismic data, convert it to obtain the wavenumber spectrum data volume of the seismic data.
[0078] Step S12: Perform the following regularization process on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume to obtain the wavenumber spectrum after frequency slice regularization.
[0079] Regularization is performed on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume. For each frequency slice, the maximum amplitude in the amplitude spectrum of the current frequency slice's wavenumber spectrum is determined according to the tilt weight function. By setting the number of iterations, least squares optimization is performed on the wavenumber component corresponding to the maximum amplitude based on the Gram matrix. The wavenumber spectrum of the current frequency slice is updated based on the least squares optimization data and the Gram matrix. Through multiple iterations with a preset number of iterations, the regularized wavenumber spectrum of the current frequency slice is obtained. After processing all slices in the wavenumber spectrum data volume, the wavenumber spectra of all frequency slices are obtained.
[0080] Step S13: Based on the regularized wavenumber spectrum of all frequency slices, obtain the regularized wavenumber spectrum data volume.
[0081] Step S14: Convert the regularized wavenumber spectrum data volume to obtain regularized time-space domain seismic data.
[0082] In some optional embodiments, step S11 above, based on the acquired time-space domain seismic data, converts the seismic data into wavenumber spectrum data, including:
[0083] 1) Perform a fast Fourier transform on the acquired time-space domain seismic data along the time direction to obtain the frequency-space domain data volume.
[0084] Optionally, the acquired time-space domain seismic data can be subjected to a fast Fourier transform along the time direction to obtain a frequency-space domain data volume. Alternatively, other transformation methods can be used, as long as they can transform the time-space domain seismic data to the frequency-space domain.
[0085] 2) Perform a non-uniform Fourier transform on each frequency slice f(x) in the frequency-space domain data volume along the spatial direction to obtain the wavenumber spectrum F(k) of the frequency slice, and generate the wavenumber spectrum data volume D(f,k) based on the wavenumber spectra of all frequency slices.
[0086] For frequency slices f(x) in the frequency-spatial domain data volume m The wavenumber spectrum F(k) of the frequency slice is obtained by performing a non-uniform Fourier transform along the spatial direction using the following formula. n ):
[0087]
[0088] Where i is the virtual unit, M is the number of samples for spatial point x, and N is the number of samples for wave number k;
[0089] The superscripts 1, 2, 3, and 4 represent four different spatial dimensions, and the subscript m represents the sample number of the spatial point x; that is, the superscripts are used to distinguish the dimensions corresponding to the four spatial directions.
[0090] The superscripts 1, 2, 3, and 4 represent four different wavenumber dimensions, and the subscript n represents the sample number of wavenumber k.
[0091] Seismic data can be generated using artificially synthesized 3D seismic records. The following example uses four spatial directions: the x and y coordinates of the common center point, the source-receiver distance, and the source-receiver azimuth. For ease of explanation, in this embodiment, the number of spatial samples before and after regularization for the source-receiver distance and azimuth is set to 1. The sample number in the x-coordinate direction is cmp, and the sample number in the y-coordinate direction is line. Figure 3 The seismic traces of the data volume are displayed, showing the regular missing 3D seismic traces displayed in cmp-line order. There are 9 missing traces out of every 10 traces in the line direction of the data volume.
[0092] In step S12 above, the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume is regularized. The process of regularizing a single frequency slice is described in [link to documentation]. Figure 2 As shown, it includes:
[0093] Step S121: Determine the maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice according to the tilt angle weighting function; wherein, the tilt angle weighting function is determined according to the amplitude energy of the wavenumber spectrum data volume in different ray directions.
[0094] In this step, the amplitude spectrum of the wavenumber spectrum of the current frequency slice is obtained. Based on the amplitude spectrum and the tilt weight function, the maximum amplitude position is searched to obtain the wavenumber position corresponding to the maximum amplitude in the amplitude spectrum. This step uses the tilt weight function to constrain the energy picking process in each iteration to select the wavenumber component corresponding to the maximum energy, thereby enhancing the regularization algorithm's resistance to spatial aliasing.
[0095] Step S122: Based on the Gram matrix, perform least squares optimization on the wavenumber component corresponding to the maximum amplitude, and update the wavenumber spectrum of the current frequency slice based on the least squares optimization data and the Gram matrix; wherein, the Gram matrix is pre-established based on the wavenumber spectrum.
[0096] In this step, the wavenumber component corresponding to the maximum amplitude is optimized using least squares based on the extracted Gram submatrix and the wavenumber spectral subvector extracted from the wavenumber spectrum, resulting in least squares optimized data. The Gram submatrix is obtained by extracting specified elements from the Gram matrix based on the iteration number of the current frequency slice. The wavenumber spectrum of the current frequency slice is updated based on the least squares optimized data and the extracted Gram submatrix. The Gram submatrix is obtained by extracting the corresponding columns from the Gram matrix based on the iteration number. This step improves the pseudo-Gibbs phenomenon at the discontinuity point of the in-phase axis by updating all selected wavenumber components through least squares optimization after each selection of the wavenumber component corresponding to the maximum energy.
[0097] Step S123: Has the preset number of iterations been reached? If yes, proceed to step S124; otherwise, return to step S121 to determine the maximum amplitude.
[0098] Step S124: Obtain the wavenumber spectrum after the current frequency slice regularization process.
[0099] In some alternative embodiments, the Gram matrix and tilt weight function can be pre-established before the frequency slices are regularized.
[0100] The Gram matrix is based on the wavenumber spectrum and can be an N-order square matrix G. The N-order square matrix G = (g1, ..., g2) can be calculated from the wavenumber spectrum. N ), where g is the nth (n≤N) column vector of G. n as follows:
[0101]
[0102] The tilt weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions. Specifically, for the wavenumber spectrum data volume D(f,k), the sum of amplitude energies along different ray directions passing through the origin can be calculated, and then the tilt weighting function w(f,k) can be obtained based on the sum of amplitude energies.
[0103] w(f,k)=(S -1 TS|D(f,k)|) q
[0104] Where D(f,k) is the wavenumber spectrum data volume, T is the smoothing operator, S is the amplitude energy operator, representing the sum of amplitude energy along the specified ray direction in the wavenumber spectrum data volume, and the exponent q represents the sharpness of the tilt angle weighting function: S -1 Let S be the inverse operator; where the expression for the amplitude energy operator S is as follows:
[0105] S(θ 1 ,θ 2 ,θ 3 ,θ4 )=∫ r |D(θ 1 ,θ 2 ,θ 3 ,θ 4 ,r)|dr
[0106]
[0107] Where, (θ 1 ,θ 2 ,θ 3 ,θ 4 (f, k) represents the Cartesian coordinates in the frequency-wavenumber domain. 1 ,k 2 ,k 3 ,k 4 The corresponding polar coordinates, r is the radial length, θ j The angle representing the dimension corresponding to the j-th wavenumber of wavenumber k; k j Let S represent the dimension corresponding to the j-th wavenumber of wavenumber k. Operator S represents the sum of ray amplitude energies along the angular ray direction in the frequency-wavenumber domain; operator T represents smoothing and other processing; and the exponent q is used to control the sharpness of the tilt weight function w. -1 This represents the inverse process of transformation S, meaning that (θ) 1 ,θ 2 ,θ 3 ,θ 4 The energy value of the ray direction corresponding to the ray is repeatedly placed into the corresponding (k) of different frequencies f. 1 ,k 2 ,k 3 ,k 4 )Location.
[0108] In step S12 above, the regularization process of steps S121-124 is performed for each frequency slice. For a frequency slice, the number of iterations can be set, and regularization is achieved through multiple iterations. For example, let the iteration number j = 1, copy the frequency slice wavenumber spectrum F(k) in the data volume D(f, k) to C(k), and use it as the current frequency slice for regularization.
[0109] Optionally, in step S121 above: when determining the maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice according to the tilt angle weighting function, the amplitude spectrum of the wavenumber spectrum of the current frequency slice can be obtained. Based on the amplitude spectrum and the tilt angle weighting function, the maximum amplitude position is searched to obtain the wavenumber position corresponding to the maximum amplitude in the amplitude spectrum. The maximum amplitude position search can be performed using the following formula:
[0110] |w(f,p j C j (p j)|≥|w(f,n)C j (n)|, 1≤n≤N
[0111] Wherein, the tilt angle weighting function is w(f,k)=(S -1 TS|D(f,k)|) q In the above formula, k in the tilt angle weighting function on the left and right sides are respectively p j and n; p j And n represents the wave number position; C j (p j ) indicates the wavenumber position p j wavenumber component at C; j (n) represents the wavenumber component at wavenumber position n.
[0112] Continuing with the example above, for the amplitude spectrum of the current frequency slice C(k), the maximum amplitude position is searched using the following formula based on the amplitude spectrum and the tilt weight function, to obtain the corresponding maximum amplitude |C j (p j The wavenumber position p of | j , position the wavenumber p j wavenumber component C at the location j (p j Insert it into the corresponding position of vector V0.
[0113] Optionally, in step S122 above, the wavenumber component corresponding to the maximum amplitude is optimized using least squares based on the Gram matrix, and the wavenumber spectrum of the current frequency slice is updated based on the least squares optimized data and the Gram matrix, including the following implementation process:
[0114] 1) Extract specified elements from the Gram matrix based on the current frequency slice iteration number to obtain the Gram submatrix, and extract the wavenumber spectral subvector from the wavenumber spectrum.
[0115] Based on the current iteration number j of the frequency slice, a sub-matrix A is obtained by extracting specified elements from the Gram matrix G, and a sub-vector B is extracted from the wavenumber spectrum, where:
[0116] A is according to (p) l ,p r The position p extracted from the Gram matrix j A square matrix of order l, ..., j, r = 1, ..., j; j represents the iteration number;
[0117] B is pressed p. l The position is extracted from the wavenumber spectrum F(k) p j A dimensional vector, l = 1, ..., j.
[0118] 2) Based on the extracted Gram subarray and the wavenumber spectral vector extracted from the wavenumber spectrum, perform least squares optimization on the wavenumber component corresponding to the maximum amplitude to obtain least squares optimized data.
[0119] The preset wavenumber component matrix V0 is updated based on the wavenumber component data at the wavenumber position corresponding to the maximum amplitude.
[0120] Based on the extracted Gram submatrix A and wavenumber spectral vector B, the wavenumber component matrix V0 is processed using the least squares method to obtain the least squares optimized matrix V; the least squares optimized matrix V can be obtained using the following least squares equation:
[0121] (A+0.01I)V=(B+0.01V0)
[0122] Where I is the identity matrix.
[0123] After obtaining the least squares optimization matrix V, matrix V can be copied to V0, that is, V0 can be updated with matrix V.
[0124] 3) Extract the corresponding columns from the Gram matrix according to the number of iterations to obtain the Gram submatrix.
[0125] Get N rows of p j The gram submatrix E is arranged according to p. l The position is extracted from the Gram matrix p j The matrix formed by the columns, l = 1, ..., j, where j represents the iteration number; specifically, N rows p can be extracted from the Gram matrix based on the iteration number j of the current frequency slice. l The columnar matrix yields matrix E, where N is the number of wavenumber samples.
[0126] 4) Update the wavenumber spectrum of the current frequency slice based on the least squares optimized data and the extracted Gram submatrix.
[0127] Calculate the product of the extracted Gram submatrix E and the least-squares optimization matrix V, and update the wavenumber spectrum of the current frequency slice based on the difference between the current frequency slice's wavenumber spectrum and this product. That is, calculate N rows of p... j Column matrices E and p j The product of the vectors V and the result of subtracting the wavenumber spectrum F(k) of the current frequency slice from the product are used to update the wavenumber spectrum C(k).
[0128] If step S123 determines that the preset number of iterations has not been reached, let the iteration number j = j + 1, and repeat steps S121 and S122 until the preset number of iterations is reached, that is, all wavenumber components of the frequency slice are obtained, and the wavenumber spectrum after the current frequency slice is regularized is obtained. Then, the obtained least squares optimization matrix V is stored in the corresponding position of the output matrix H.
[0129] The output matrix H is a matrix that stores the regularized wavenumber spectra of all frequency slices. The matrix V of the regularized wavenumber spectrum of each frequency slice is stored in the corresponding position until all frequency slices are stored. The regularized wavenumber spectrum data volume can be obtained from the matrix H.
[0130] In some optional embodiments, step S14 above transforms the regularized wavenumber spectrum data volume to obtain regularized time-space domain seismic data, including: performing a non-uniform inverse Fourier transform on the wavenumber spectrum of each frequency slice in the regularized wavenumber spectrum data volume of the seismic data to obtain a regularized frequency-space domain data volume; and performing an inverse Fourier transform on the regularized frequency-space domain data volume to obtain regularized time-space domain seismic data. Continuing with the example above, after obtaining all wavenumber components of all frequency slices, the output matrix H is then inversely transformed back to the time-space domain as the data regularization result.
[0131] The frequency-wavenumber domain seismic data regularization method provided in this invention performs regularization processing on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume of the seismic data converted to the frequency-wavenumber domain. During the processing, the regularization processing of seismic data with regularity loss and large gaps is improved. The dip angle scan in the frequency-wavenumber domain generates a dip angle weight function, which is used to constrain the energy picking process in each iteration, thereby enhancing the anti-spatial aliasing ability of the regularization algorithm. After selecting the wavenumber component with the maximum energy at each step, all selected wavenumber components are updated through a least squares optimization, which improves the pseudo-Gibbs phenomenon at the discontinuity points of the phase axis.
[0132] The above method incorporates phase axis tilt constraints to enhance the regularization algorithm's resistance to spatial aliasing. The addition of least squares optimization improves the pseudo-Gibbs phenomenon at discontinuities in the phase axis, thereby improving the spatial sampling properties of seismic data and enhancing the imaging quality of seismic data.
[0133] Using the above-mentioned seismic exploration data regularization method to process artificially synthesized 3D seismic record data can achieve excellent results. The following analysis and explanation are based on the graphical results under different conditions, including the comparison of data regularization processing results with and without dip constraints and with and without least squares optimization.
[0134] Figure 4This represents the regularization result without tilt weight constraints and without least squares optimization. Without tilt weight constraints means that each element of the tilt weight function w in step S121 has a value of 1, i.e., no tilt weight constraints are applied; without least squares optimization means that the least squares optimization in step S122 is not performed, and V0 is directly copied to V.
[0135] Figure 5 This is the regularization result when there is an inclination constraint and no least squares optimization. That is, the least squares optimization in step S122 is not performed, and V0 is directly copied to V.
[0136] Figure 6 The result is the regularized result after least squares optimization with tilt angle constraints, i.e., the result after regularization using the processing steps S121-124.
[0137] right Figure 4 and Figure 5 It is known that for seismic data lacking rules, due to the interference of periodic pseudo-frequency energy, the regularization process will produce inaccurate reconstruction results due to incorrect energy search when there is no dip constraint. After adding the dip constraint of the phase axis, the regularization effect is significantly improved, but pseudo-Gibbs phenomenon will occur at the discontinuity point of the phase axis.
[0138] contrast Figure 5 and Figure 6 It can be seen that after introducing least squares optimization, the pseudo Gibbs phenomenon at the discontinuity point of the in-phase axis is suppressed, and the reconstruction effect is further improved.
[0139] The above practices demonstrate that the frequency-wavenumber domain seismic data regularization method provided in this invention can improve the spatial sampling attributes of seismic data and enhance the imaging quality of seismic data after processing data with missing rules. The regularization effect is significantly improved, and the reconstruction effect is further improved.
[0140] Based on the same inventive concept, embodiments of the present invention also provide a seismic exploration data regularization processing device. This device can be installed in a computer device with data processing capabilities, and its structure is as follows: Figure 7 As shown, it includes: a first data acquisition module 101, a rule-based processing module 102, and a second data acquisition module 103.
[0141] The first data acquisition module 101 is used to convert the acquired time-space domain seismic data into wavenumber spectrum data volume of the seismic data.
[0142] The regularization processing module 102 is used to perform the following regularization processing on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume:
[0143] The maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice is determined based on the tilt angle weighting function; the tilt angle weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions.
[0144] Based on the Gram matrix, the wavenumber component corresponding to the maximum amplitude is optimized by least squares. Based on the least squares optimization data and the Gram matrix, the wavenumber spectrum of the current frequency slice is updated. The Gram matrix is pre-established based on the wavenumber spectrum.
[0145] Return to the step with the maximum amplitude until the preset number of iterations is reached, and obtain the wavenumber spectrum after the current frequency slice is regularized.
[0146] Based on the regularized wavenumber spectrum of all frequency slices, the regularized wavenumber spectrum data volume is obtained;
[0147] The second data acquisition module 103 is used to convert the regularized wavenumber spectrum data volume to obtain regularized time-space domain seismic data.
[0148] This invention also provides a computer storage medium storing computer-executable instructions, which, when executed by a processor, implement the above-described seismic exploration data regularization processing method.
[0149] This invention also provides a terminal device, including: a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the program to implement the above-described seismic exploration data regularization processing method.
[0150] Regarding the seismic exploration data regularization processing device in the above embodiments, the specific methods by which each module performs its operations have been described in detail in the embodiments related to the method, and will not be elaborated here.
[0151] Unless otherwise specifically stated, terms such as processing, calculation, operation, determination, display, etc., may refer to the actions and / or processes of one or more processing or computing systems or similar devices that represent the manipulation and conversion of data representing physical (e.g., electronic) quantities within the registers or memory of the processing system into other data similarly representing physical quantities within the memory, registers, or other such information storage, transmission, or display devices of the processing system. Information and signals can be represented using any of a variety of different techniques and methods. For example, data, instructions, commands, information, signals, bits, symbols, and chips mentioned throughout the above description can be represented by voltage, current, electromagnetic waves, magnetic fields or particles, light fields or particles, or any combination thereof.
[0152] It should be understood that the specific order or hierarchy of steps in the disclosed process is an example of an exemplary method. Based on design preferences, it should be understood that the specific order or hierarchy of steps in the process may be rearranged without departing from the scope of this disclosure. The appended method claims provide elements of various steps in an exemplary order and are not intended to limit the scope to the specific order or hierarchy described.
[0153] In the detailed description above, various features are combined together in a single embodiment to simplify this disclosure. This approach to disclosure should not be construed as reflecting an intention that embodiments of the claimed subject matter require more features than are explicitly stated in each claim. Rather, as reflected in the appended claims, the invention is presented with fewer features than all of the features in a single disclosed embodiment. Therefore, the appended claims are hereby explicitly incorporated into the detailed description, with each claim representing a separate preferred embodiment of the invention.
[0154] Those skilled in the art will also understand that the various illustrative logic blocks, modules, circuits, and algorithm steps described in conjunction with the embodiments herein can be implemented as electronic hardware, computer software, or a combination thereof. To clearly illustrate the interchangeability between hardware and software, the various illustrative components, blocks, modules, circuits, and steps described above are generally described in terms of their functionality. Whether such functionality is implemented as hardware or software depends on the specific application and the design constraints imposed on the overall system. Those skilled in the art can implement the described functionality in alternative ways for each specific application; however, such implementation decisions should not be construed as departing from the scope of this disclosure.
[0155] The steps of the methods or algorithms described in conjunction with the embodiments herein can be directly embodied in hardware, software modules executed by a processor, or a combination thereof. The software modules can reside in RAM memory, flash memory, ROM memory, EPROM memory, EEPROM memory, registers, hard disks, removable disks, CD-ROMs, or any other form of storage medium well known in the art. An exemplary storage medium is connected to the processor, enabling the processor to read information from and write information to the storage medium. Of course, the storage medium can also be a component of the processor. The processor and storage medium can reside in an ASIC. The ASIC can reside in a user terminal. Alternatively, the processor and storage medium can exist as discrete components in the user terminal.
[0156] For software implementation, the techniques described in this application can be implemented using modules (e.g., procedures, functions, etc.) that perform the functions described in this application. This software code can be stored in memory units and executed by a processor. The memory units can be implemented within the processor or outside the processor; in the latter case, they are communicatively coupled to the processor via various means, as is well known in the art.
[0157] The foregoing description includes examples of one or more embodiments. It is certainly impossible to describe all possible combinations of components or methods in order to describe the above embodiments, but those skilled in the art will recognize that further combinations and arrangements of the various embodiments are possible. Therefore, the embodiments described herein are intended to cover all such changes, modifications, and variations that fall within the scope of the appended claims. Furthermore, the term "comprising" as used in the specification or claims is interpreted in a manner similar to the term "including," as interpreted when used as a conjunction in the claims. Additionally, the use of any term "or" in the specification of the claims is intended to mean "non-exclusive or."
Claims
1. A method for regularizing seismic exploration data, characterized in that, include: Based on the acquired time-space domain seismic data, the wavenumber spectrum data volume of the seismic data is obtained by conversion; The following regularization process is performed on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume: The maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice is determined based on the tilt angle weighting function; the tilt angle weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions. Based on the Gram matrix, the wavenumber components corresponding to the maximum amplitude are optimized using least squares. Specifically, this includes: updating the preset wavenumber component matrix based on the wavenumber component data at the wavenumber position corresponding to the maximum amplitude; processing the wavenumber component matrix using the least squares method based on the extracted Gram submatrix and wavenumber spectrum subvector to obtain the least squares optimized matrix; the Gram submatrix is obtained by extracting specified elements from the Gram matrix according to the iteration number of the current frequency slice, and the Gram matrix is pre-established based on the wavenumber spectrum; The wavenumber spectrum of the current frequency slice is updated based on the least squares optimization data and the Gram matrix. Specifically, this includes: calculating the product of the extracted Gram submatrix and the least squares optimization matrix, and updating the wavenumber spectrum of the current frequency slice based on the difference between the current wavenumber spectrum and the product; the Gram submatrix is obtained by extracting the corresponding columns from the Gram matrix according to the number of iterations. Return to the step of determining the maximum amplitude until the preset number of iterations is reached, and obtain the wavenumber spectrum after the current frequency slice regularization process; Based on the regularized wavenumber spectrum of all frequency slices, the regularized wavenumber spectrum data volume is obtained; The regularized wavenumber spectrum data volume is transformed to obtain regularized time-space domain seismic data.
2. The method as described in claim 1, characterized in that, The wavenumber spectrum data of the acquired time-space domain seismic data is obtained by converting the seismic data, including: The acquired time-space domain seismic data is subjected to a Fourier transform along the time direction to obtain the frequency-space domain data volume; A non-uniform Fourier transform is performed on the frequency slices in the frequency-space domain data volume along the spatial direction to obtain the wavenumber spectrum of the frequency slices, and a wavenumber spectrum data volume is generated based on the wavenumber spectra of all frequency slices.
3. The method as described in claim 2, characterized in that, The step of performing a non-uniform Fourier transform along the spatial direction on the frequency slices in the frequency-spatial domain data volume to obtain the wavenumber spectrum of the frequency slices includes: Frequency slices in frequency-spatial domain data volumes The wavenumber spectrum of the frequency slice is obtained by performing a non-uniform Fourier transform along the spatial direction using the following formula. : ； Where i is the virtual unit, M is the number of samples for spatial point x, and N is the number of samples for wave number k; The superscripts 1, 2, 3, and 4 represent four different spatial dimensions, and the subscript m represents the sample number of spatial point x. The superscripts 1, 2, 3, and 4 represent four different wavenumber dimensions, and the subscript n represents the sample number of wavenumber k.
4. The method as described in claim 1, characterized in that, Based on the tilt angle weighting function, determine the maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice, including: Obtain the amplitude spectrum of the wavenumber spectrum of the current frequency slice, and search for the position of the maximum amplitude based on the amplitude spectrum and the tilt weight function to obtain the wavenumber position corresponding to the maximum amplitude in the amplitude spectrum.
5. The method as described in claim 4, characterized in that, The search for the location of maximum amplitude is based on the amplitude spectrum and the tilt weight function, including: The maximum amplitude location is searched using the following formula: ； Wherein, the tilt angle weighting function is In the above formula, k in the tilt angle weighting function on the left and right sides are respectively and n; And n represents the wave number position; Indicates wavenumber position wavenumber components at that location; Let f represent the wavenumber component at wavenumber position n, and f represent the frequency. In the tilt angle weighting function, Let T be the wavenumber spectrum data volume, S be the smoothing operator, S be the amplitude energy operator, representing the sum of amplitude energy along a specified ray direction in the wavenumber spectrum data volume, and q be the exponent representing the sharpness of the tilt angle weighting function. Let S be the inverse operator; where the expression for the amplitude energy operator S is as follows: ； ； Cartesian coordinates in the frequency-wavenumber domain The corresponding polar coordinates, where r is the radial length. The angle representing the dimension corresponding to the j-th wavenumber of wavenumber k; Let f represent the dimension corresponding to the j-th wavenumber of wavenumber k, and let f represent the frequency.
6. The method as described in claim 1, characterized in that, Based on the extracted Gram submatrix A and wavenumber spectral vector B, the wavenumber component matrix V0 is processed using the least squares method to obtain the least squares optimized matrix V, which includes: The least squares optimization matrix V is obtained by using the following least squares equation: ； Where I is the identity matrix; A is according to The position is extracted from the Gram matrix. Square array , B is according to The position is extracted from the wavenumber spectrum F(k). dimensional vector, ; j represents the iteration number; Accordingly, the product of the extracted Gram submatrix E and the least-squares optimization matrix V is calculated, including: Get N rows The gram submatrix E is arranged according to... Positions are extracted from the Gram matrix The matrix formed by columns j represents the number of iterations; Calculate the product of the Gram submatrix E and the least squares optimization matrix V.
7. The method according to any one of claims 1-6, characterized in that, The regularized wavenumber spectrum data volume of seismic data is transformed to obtain regularized time-space domain seismic data, including: Non-uniform inverse Fourier transform is performed on the wavenumber spectrum of each frequency slice in the regularized wavenumber spectrum data volume of the seismic data to obtain the regularized frequency-spatial domain data volume. An inverse Fourier transform is performed on the regularized frequency-space domain data volume to obtain regularized time-space domain seismic data.
8. A seismic exploration data regularization processing device, characterized in that, include: The first data acquisition module is used to convert the acquired time-space domain seismic data into wavenumber spectrum data volume of the seismic data. The regularization processing module is used to perform the following regularization processing on the wavenumber spectrum of each frequency slice in the wavenumber spectrum data volume: The maximum amplitude in the amplitude spectrum of the wavenumber spectrum of the current frequency slice is determined based on the tilt angle weighting function; the tilt angle weighting function is determined based on the amplitude energy of the wavenumber spectrum data volume in different ray directions. Based on the Gram matrix, the wavenumber components corresponding to the maximum amplitude are optimized using least squares. Specifically, this includes: updating the preset wavenumber component matrix based on the wavenumber component data at the wavenumber position corresponding to the maximum amplitude; processing the wavenumber component matrix using the least squares method based on the extracted Gram submatrix and wavenumber spectrum subvector to obtain the least squares optimized matrix; the Gram submatrix is obtained by extracting specified elements from the Gram matrix according to the iteration number of the current frequency slice, and the Gram matrix is pre-established based on the wavenumber spectrum; The wavenumber spectrum of the current frequency slice is updated based on the least squares optimization data and the Gram matrix. The Gram matrix is pre-established based on the wavenumber spectrum and specifically includes: calculating the product of the extracted Gram sub-matrix and the least squares optimization matrix, and updating the wavenumber spectrum of the current frequency slice based on the difference between the current frequency slice's wavenumber spectrum and the product. The Gram sub-matrix is obtained by extracting the corresponding columns from the Gram matrix according to the number of iterations. Return to the step with the maximum amplitude until the preset number of iterations is reached, and obtain the wavenumber spectrum after the current frequency slice is regularized. Based on the regularized wavenumber spectrum of all frequency slices, the regularized wavenumber spectrum data volume is obtained; The second data acquisition module is used to convert the regularized wavenumber spectrum data volume to obtain regularized time-space domain seismic data.
9. A computer storage medium, characterized in that, The computer storage medium stores computer-executable instructions, which, when executed by a processor, implement the seismic exploration data regularization processing method according to any one of claims 1-7.
10. A terminal device, characterized in that, include: A memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor, when executing the program, implements the seismic exploration data regularization processing method according to any one of claims 1-7.