阅读说明：本技术 一种深层叠前地震数据噪声压制方法及系统 (Deep pre-stack seismic data noise suppression method and system ) 是由 王晓凯 赵晨 陈文超 师振盛 于 2021-06-17 设计创作，主要内容包括：本发明公开了一种深层叠前地震数据噪声压制方法及系统,对原始叠前数据进行分块；对分块后的其中一块原始叠前地震数据利用中值估计算子估算该块数据的噪声方差等级；对估计过噪声方差等级的原始叠前地震数据分块进行高维连续小波变换；依据相干干扰的角度确定分块需要保留的连续小波变换变换结果角度范围；根据分块的噪声标准差确定其需要保留的子带内的阈值；利用确定好的每个保留子带的阈值对需要保留子带中的小波系数进行阈值处理；对阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号；重复以上步骤直至所有分块处理完成,将得到的所有高信噪比分块重新拼接为完整的高信噪比数据作为完成噪声压制的深层叠前地震数据。(The invention discloses a deep prestack seismic data noise suppression method and a system, which are used for partitioning original prestack data; estimating the noise variance level of one block of original pre-stack seismic data after being partitioned by using a median estimation operator; performing high-dimensional continuous wavelet transformation on the original pre-stack seismic data blocks with the estimated noise variance level; determining the angle range of continuous wavelet transform results required to be reserved for partitioning according to the angle of coherent interference; determining a threshold value in a subband needing to be reserved according to the standard deviation of the noise of the block; performing threshold processing on the wavelet coefficients in the sub-bands needing to be reserved by using the determined threshold of each reserved sub-band; carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band subjected to threshold processing to obtain a reconstructed high signal-to-noise ratio signal; and repeating the steps until all the blocks are processed, and splicing all the obtained high signal-to-noise ratio blocks into complete high signal-to-noise ratio data serving as deep pre-stack seismic data for completing noise suppression.)

1. A deep prestack seismic data noise suppression method is characterized by comprising the following steps:

s1, acquiring original pre-stack data and partitioning;

s2, estimating the noise variance grade of any original pre-stack seismic data after the block division in the step S1 by using a median estimation operator;

s3, performing high-dimensional continuous wavelet transform on the original pre-stack seismic data estimated in the step S2 to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

s4, determining the high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved in the coefficient sub-band obtained in the step S3 according to the angle of coherent interference;

s5, determining the threshold value in the sub-band of the high-dimensional continuous wavelet transform coefficient required to be reserved in the step S4 according to the noise variance level obtained in the step S2

S6, LiThreshold value determined in step S5Carrying out threshold processing on wavelet coefficients in high-dimensional continuous wavelet transform coefficient sub-bands needing to be reserved;

s7, carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band subjected to threshold processing in the step S6 to obtain a reconstructed signal with high signal-to-noise ratio

S8, repeating the steps S2-S7 until all the blocks acquired in the step S1 are processed completely, and then obtaining all the high signal-to-noise ratio signalsAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

2. The method of claim 1, wherein in step S1, the original pre-stack seismic data is divided into a plurality of nxnxnxn data blocks with overlaps disposed between the data blocks, the number of blocks being:

wherein n is the amount of overlap,denotes rounding up, N_x×N_y×N_zIs a three-dimensional prestack data volume.

3. The method of claim 1, wherein in step S2, the noise standard deviation levelThe method specifically comprises the following steps:

where MEDIAN is the MEDIAN operator, F (k) is the result of three-dimensional Fourier transform on the block, [ k ]_min,k_max]A high frequency range of F (k).

4. The method of claim 1, wherein in step S3, the high-dimensional continuous wavelet transform is performed on an arbitrarily selected unfinished block f (x) to obtain the high-frequency sub-band CWT (f; b, a, θ, φ) as follows:

wherein the content of the first and second substances,representing wavenumber domain data obtained by high-dimensional Fourier transform of an original data block,is the wavenumber domain obtained by the wavelet atom under the parameter group (a, theta, phi; k) through the high-dimensional Fourier transform,is composed ofIFFT represents an inverse fourier transform.

5. The method according to claim 1, wherein in step S4, the noise tilt θ in the data volume is determined by coherent interference angle_{coherent_noise}(ii) a Keeping the angle theta smaller than theta in the continuous wavelet transform result CWT (f; b, a, theta, phi)_{coherent_noise}Sub-band ofI.e. retaining only [ -theta ]_{coherent_noise}+Δθ,θ_{coherent_noise}-Δθ]Sub-bands within an angular range.

6. The method according to claim 1, wherein in step S5, the threshold value thre (a) of each high frequency subband is:

wherein the content of the first and second substances,is a threshold coefficient of each sub-band, which takes the value of 1/a³Is in direct proportion; c is a constant controlling the degree of noise suppression, typically taken to be 3,is the standard deviation.

7. The method according to claim 1, wherein in step S6, the hard clipping filtering is performed according to the threshold coefficient d determined in step S5, specifically:

wherein, CWT_thr(f; b, a, theta, phi) are wavelet coefficients that have been threshold filtered, and CWT (f; b, a, theta, phi) are wavelet coefficients that have not been filtered.

8. The method of claim 1, wherein in step S7, the signal with high signal-to-noise ratio is used

Wherein, CWT_thrAnd (f, b, a, theta and phi) are wavelet coefficients after threshold filtering, a is a scale factor of the wavelet atoms, phi is the angle of the polarization direction of the wavelet atoms on the azimuth angle, and theta is the angle of the planned direction of the wavelet atoms on the dip angle.

9. The method according to claim 1, wherein step S8 is specifically: after all the data blocks divided in the step S1 are subjected to the operations from the step S2 to the step S7, the data after noise suppression are obtained by splicing the data blocks again according to the arrangement during division; before splicing, two adjacent data blocks are respectively removed from the corresponding directionsN is the overlap amount.

10. A deep prestack seismic data noise suppression system, comprising:

the blocking module is used for acquiring original pre-stack data and blocking the original pre-stack data;

the transformation module estimates the noise variance grade of any original pre-stack seismic data after the blocking of the blocking module by using a median estimation operator;

the estimation module is used for carrying out high-dimensional continuous wavelet transform on the original pre-stack seismic data estimated by the transformation module to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

the selection module determines a high-dimensional continuous wavelet transform coefficient sub-band which needs to be reserved in the coefficient sub-band obtained by the estimation module according to the angle of the coherent interference;

a threshold module for determining the threshold value in the sub-band of the high-dimensional continuous wavelet transform coefficient to be reserved in the selection module according to the noise variance grade obtained in the transformation module

Treatment ofModule, threshold value determined by threshold value moduleCarrying out threshold processing on wavelet coefficients in high-dimensional continuous wavelet transform coefficient sub-bands needing to be reserved;

a reconstruction module for performing high-dimensional wavelet inverse transformation on the wavelet transformation sub-band processed by the threshold value processing module to obtain a reconstructed high signal-to-noise ratio signal

A pressing module for processing all the blocks obtained by the block module and then obtaining all the high signal-to-noise ratio signalsAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

Technical Field

The invention belongs to the technical field of exploration geophysics, and particularly relates to a deep prestack seismic data noise suppression method and system.

Background

In recent years, as shallow hydrocarbon resources are being explored in large quantities, hydrocarbon resource exploration is gradually shifted to the deep subterranean zones. The seismic waves are attenuated faster in the stratum, and the deep reflection energy is weak, so that the signal-to-noise ratio of deep seismic data is lower. In addition, the surface wave will have a larger influence in the deep layer, further reducing the signal-to-noise ratio of the seismic data, and although the post-stack data with higher signal-to-noise ratio can be obtained by stacking the seismic data (pre-stack data) acquired for many times, many pre-stack processes and analyses still require pre-stack data with high signal-to-noise ratio.

Prior art 1: prediction based filtering methods (e.g., FX, TX, polynomial fitting, etc.). Such methods utilize a linear predictive relationship between adjacent traces to estimate a reflection signal from a noisy seismic signal. The method has good performance when the target stratum signal is a linear in-phase axis, and has poor denoising result when the target stratum signal is a nonlinear in-phase axis; furthermore, such methods do not have the resolution capability to distinguish between in-phase axis angles.

Prior art 2: a filtering method based on a high-dimensional discrete wavelet domain. The discrete wavelet transform is suitable for sparse representation of point singular targets, and can realize better suppression of random noise by methods such as threshold processing and the like. However, the discrete wavelet atom has the characteristic of isotropy, which makes it impossible to represent the direction information of the signal, and thus it is difficult to process structural noise having a certain direction.

Disclosure of Invention

The technical problem to be solved by the present invention is to provide a method and a system for suppressing deep prestack seismic data noise, which perform direction and threshold filtering on deep prestack seismic signals in a high-dimensional wavelet domain, so as to obtain seismic data with high signal-to-noise ratio.

The invention adopts the following technical scheme:

a deep prestack seismic data noise suppression method comprises the following steps:

s1, acquiring original pre-stack data and partitioning;

s2, estimating the noise variance grade of any original pre-stack seismic data after the block division in the step S1 by using a median estimation operator;

s3, performing high-dimensional continuous wavelet transform on the original pre-stack seismic data estimated in the step S2 to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

s4, determining the high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved in the coefficient sub-band obtained in the step S3 according to the angle of coherent interference;

s5, determining the threshold value in the sub-band of the high-dimensional continuous wavelet transform coefficient required to be reserved in the step S4 according to the noise variance level obtained in the step S2

S6, utilizing the threshold determined in step S5Carrying out threshold processing on wavelet coefficients in high-dimensional continuous wavelet transform coefficient sub-bands needing to be reserved;

s7, carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band subjected to threshold processing in the step S6 to obtain a reconstructed signal with high signal-to-noise ratio

S8, repeating the steps S2-S7 until all the blocks acquired in the step S1 are processed completely, and then obtaining all the high signal-to-noise ratio signalsAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

Specifically, in step S1, the original pre-stack seismic data is divided into a plurality of nxnxnxnxn data blocks, and an overlap is provided between the data blocks, where the number of the divided blocks is:

wherein n is the amount of overlap,denotes rounding up, N_x×N_y×N_zIs a three-dimensional prestack data volume.

In particularIn step S2, the noise standard deviation levelThe method specifically comprises the following steps:

where MEDIAN is the MEDIAN operator, F (k) is the result of three-dimensional Fourier transform on the block, [ k ]_min,k_max]A high frequency range of F (k).

Specifically, in step S3, a high-dimensional continuous wavelet transform is performed on an arbitrarily selected block f (x) that is not processed, so as to obtain a high-frequency subband CWT (f; b, a, θ, Φ) as follows:

wherein the content of the first and second substances,representing wavenumber domain data obtained by high-dimensional Fourier transform of an original data block,is the wavenumber domain obtained by the wavelet atom under the parameter group (a, theta, phi; k) through the high-dimensional Fourier transform,is composed ofIFFT represents an inverse fourier transform.

Specifically, in step S4, the noise tilt angle θ in the data volume is determined by the coherent interference angle_{coherent_noise}(ii) a Keeping the angle theta smaller than theta in the continuous wavelet transform result CWT (f; b, a, theta, phi)_{coherent_noise}Sub-bands of, i.e. retaining only [ -theta [ - ]_{coherent_noise}+Δθ,θ_{coherent_noise}-Δθ]Sub-bands within an angular range.

Specifically, in step S5, the threshold thre (a) of each high-frequency subband is:

wherein the content of the first and second substances,is a threshold coefficient of each sub-band, which takes the value of 1/a³Is in direct proportion; c is a constant controlling the degree of noise suppression, typically taken to be 3,is the standard deviation.

Specifically, in step S6, the hard clipping filtering is performed according to the threshold coefficient d determined in step S5, specifically:

wherein, CWT_thr(f; b, a, theta, phi) are wavelet coefficients that have been threshold filtered, and CWT (f; b, a, theta, phi) are wavelet coefficients that have not been filtered.

Specifically, in step S7, the signal with high snr is obtained

Wherein, CWT_thrAnd (f, b, a, theta and phi) are wavelet coefficients after threshold filtering, a is a scale factor of the wavelet atoms, phi is the angle of the polarization direction of the wavelet atoms on the azimuth angle, and theta is the angle of the planned direction of the wavelet atoms on the dip angle.

Specifically, step S8 specifically includes: the symbol obtained in step S1After the data blocks are subjected to the operations from the step S2 to the step S7, splicing the data blocks again according to the arrangement during the division to obtain the data after noise suppression; before splicing, two adjacent data blocks are respectively removed from the corresponding directionsN is the overlap amount.

Another technical solution of the present invention is a deep prestack seismic data noise suppression system, including:

the blocking module is used for acquiring original pre-stack data and blocking the original pre-stack data;

the transformation module estimates the noise variance grade of any original pre-stack seismic data after the blocking of the blocking module by using a median estimation operator;

the estimation module is used for carrying out high-dimensional continuous wavelet transform on the original pre-stack seismic data estimated by the transformation module to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

the selection module determines a high-dimensional continuous wavelet transform coefficient sub-band which needs to be reserved in the coefficient sub-band obtained by the estimation module according to the angle of the coherent interference;

a threshold module for determining the threshold value in the sub-band of the high-dimensional continuous wavelet transform coefficient to be reserved in the selection module according to the noise variance grade obtained in the transformation module

A processing module for determining the threshold value by using the threshold value moduleCarrying out threshold processing on wavelet coefficients in high-dimensional continuous wavelet transform coefficient sub-bands needing to be reserved;

a reconstruction module for performing high-dimensional wavelet inverse transformation on the wavelet transformation sub-band processed by the threshold value processing module to obtain a reconstructed high signal-to-noise ratio signal

A pressing module for processing all the blocks obtained by the block module and then obtaining all the high signal-to-noise ratio signalsAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

Compared with the prior art, the invention has at least the following beneficial effects:

the invention relates to a deep prestack seismic data noise suppression method, which uses continuous high-dimensional wavelet transformation to perform noise suppression processing on deep prestack seismic data; firstly, data is blocked according to factors such as the size of a computer memory and the like to be suitable for processing, then the noise standard deviation of each block is estimated according to the median of the absolute value of each block, a new fast continuous wavelet transform algorithm based on fast Fourier transform is used for carrying out high-dimensional continuous wavelet transform on the blocks of the original high-dimensional seismic data to obtain a series of high-dimensional continuous wavelet transform coefficient sub-bands, then the high-dimensional continuous wavelet transform coefficient sub-bands needing to be stored are determined according to the direction of coherent noise in the seismic data, the threshold value of each sub-band is determined according to the estimated noise standard deviation, threshold processing is carried out on the reserved sub-bands, wavelet inverse transform is carried out on the sub-bands subjected to the threshold processing to obtain blocked seismic data subjected to noise suppression, and the data blocks are spliced to obtain the complete noise suppression data after all the blocks are processed. Compared with the conventional method for suppressing the noise of the pre-stack seismic data, the method for suppressing the noise of the pre-stack seismic data uses the high-dimensional continuous wavelet method, can simultaneously process random noise and coherent noise, and can obtain a good noise suppression result.

Furthermore, the problems of computer memory overflow and overlarge operation time expense caused by overlarge data are solved by partitioning the data. The overlapping is arranged between the data blocks, so that the subsequent processing can solve the problem of edge degradation caused by discrete Fourier transform in a cutting mode.

Furthermore, the new rapid continuous wavelet transform algorithm is used, so that the time complexity of calculation is favorably reduced, and the high-dimensional continuous wavelet transform has practical value in seismic data noise suppression.

Further, the standard deviation of the noise in the data is estimatedThe standard deviation is proportional to a constant c controlling the filter strength and a constant cIs the threshold d (a, θ, φ) for that subband.

Furthermore, directional noise is filtered by removing the sub-band, so that coherent noise in the direction different from that of the reflected signal can be filtered.

Furthermore, the threshold processing of the sub-band to be reserved is beneficial to removing most random noise while reserving the effective signal, because the random noise is distributed more uniformly in the same sub-band, and the effective signal is distributed more intensively, which shows that only a small part of wavelet coefficients in the same sub-band have larger absolute values and contain most effective signals, and the other part with smaller wavelet coefficients contains most noise, therefore, most random noise can be suppressed and the effective signal can be reserved through the threshold processing.

Further, the rapid continuous wavelet inverse transformation is used in step S7, which is beneficial to accelerating the continuous wavelet inverse transformation process and further improving the processing efficiency.

Further, edge cropping is performed on the processed data block to avoid edge degradation caused by fourier transform, and since an overlap amount is already left when the data block is partitioned in step S1, it is ensured that cropping does not cause data loss. And then complete deep pre-stack seismic data subjected to noise suppression is obtained through data splicing.

Furthermore, the directional noise of the data collected in the actual work area is often similar to the distribution of the formation data or has a certain rule, so that the data of the whole work area can be processed by using the same parameters, and the processing efficiency can be improved. [ xk1]

In summary, the invention can effectively and rapidly realize noise suppression on deep pre-stack seismic data, and the adopted method is rapid high-dimensional continuous wavelet transformation, direction subband selection and hard amplitude limiting threshold processing, and has the characteristic of better processing directional noise and can remove random noise.

The technical solution of the present invention is further described in detail by the accompanying drawings and embodiments.

Drawings

FIG. 1 is a flow chart of the present invention;

FIG. 2 is a graph of a signal for an actual CMP gather;

FIG. 3 is a graph of the signal of FIG. 2 after noise suppression by the method of the present patent;

FIG. 4 is a difference section of FIGS. 2 and 3;

FIG. 5 is a line graph of a single coverage inline before an actual deflection;

FIG. 6 is a graph of the noise suppressed result of the signal of FIG. 5 by the method of the present invention;

fig. 7 is a difference section of fig. 5 and 6.

Detailed Description

The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are some, not all, embodiments of the present invention. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.

It will be understood that the terms "comprises" and/or "comprising," when used in this specification and the appended claims, specify the presence of stated features, integers, steps, operations, elements, and/or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components, and/or groups thereof.

It is also to be understood that the terminology used in the description of the invention herein is for the purpose of describing particular embodiments only and is not intended to be limiting of the invention. As used in the specification of the present invention and the appended claims, the singular forms "a," "an," and "the" are intended to include the plural forms as well, unless the context clearly indicates otherwise.

It should be further understood that the term "and/or" as used in this specification and the appended claims refers to and includes any and all possible combinations of one or more of the associated listed items.

Various structural schematics according to the disclosed embodiments of the invention are shown in the drawings. The figures are not drawn to scale, wherein certain details are exaggerated and possibly omitted for clarity of presentation. The shapes of various regions, layers and their relative sizes and positional relationships shown in the drawings are merely exemplary, and deviations may occur in practice due to manufacturing tolerances or technical limitations, and a person skilled in the art may additionally design regions/layers having different shapes, sizes, relative positions, according to actual needs.

The invention provides a deep prestack seismic data noise suppression method, which uses high-dimensional continuous wavelet transformation to perform noise suppression processing on deep prestack seismic data. The invention firstly blocks the data according to factors such as the size of a computer memory and the like to be suitable for processing, then estimates the noise standard deviation according to the median of the absolute value of each block, uses a new fast continuous wavelet transform algorithm based on fast Fourier transform to carry out high-dimensional continuous wavelet transform on the blocks of the original high-dimensional seismic data to obtain coefficient sub-bands, then determines the sub-bands needing to be stored according to the direction of coherent noise in the seismic data, determines a threshold value according to the estimated noise standard deviation, carries out threshold value processing on the reserved sub-bands, carries out wavelet inverse transform on the sub-bands subjected to threshold value processing to obtain the block seismic data after noise suppression, and splices the data blocks after all block processing to obtain the complete noise suppression data.

Referring to fig. 1, the deep prestack seismic data noise suppression method of the present invention includes the following steps:

s1, partitioning the original pre-stack data;

take a three-dimensional prestack data volume as an example, with specification N_x×N_y×N_z. Selecting proper block size according to the memory and parallel computing capability of a computer, and dividing the original pre-stack seismic data into a plurality of NXNXN data blocks. In order to solve the problem that the boundary of the blocks is not smooth due to the selection of thresholds of different subsequent data blocks and numerical calculation errors, partial overlapping needs to be carried out among the data blocks so as to ensure that the discontinuity of the boundary of the blocks can be solved through corresponding processing when splicing is carried out after the processing is finished.

If the overlapping amount in the x direction, the y direction and the z direction is n points, the number of the blocks will be nWherein the symbolsIndicating rounding up.

S2, estimating the noise variance level of one block of original pre-stack seismic data after being partitioned in the step S1 by using a median estimation operator;

and (f) (x) performing three-dimensional Fourier transform on an arbitrarily selected unfinished block f (x) to obtain a result F (k). Selecting its high frequency range [ k ]_min,k_max]Within this range, the median of the absolute value | F (k) | of the Fourier transform coefficient is obtainedWhere MEDIAN is the MEDIAN estimator for | F (k) |. The estimated noise standard deviation level in the dataMedian of | F (k) | and | F (k) |The following steps are involved:

s3, performing high-dimensional continuous wavelet transform on the original pre-stack seismic data selected in the step S2 to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

taking the unfinished block f (x) in step S2 as an example, the result CWT (f; b, a, θ, Φ) of the fast continuous wavelet transform is:

wherein, the parameter b represents the offset vector of the wavelet atom (a three-dimensional vector under the three-dimensional data volume), ψ is a wavelet function with direction selectivity selected in advance, a represents the scale factor of wavelet atom scaling, θ represents the rotation angle of the wavelet atom around the inclination angle, and φ represents the rotation angle of the wavelet atom around the azimuth angle. r is_θ,φ(x-b) represents rotation of the vector x-b about the inclination angle, theta, and about the azimuth angle, phi. The result CWT (f; b, a, theta, phi) of a three-dimensional original pre-stack seismic data after continuous wavelet transform is six-dimensional, and the complexity of the calculation time is O (N)⁶N_aN_θN_φ) N is the edge length of the data volume block, N_aIs the number of scales, N_θIs the number of the inclination angles, N_φIs the number of azimuths.

Since the above calculation amount is extremely large and the calculation cost is extremely high, the inverse fourier transform is used to rewrite equation (1) to accelerate:

wherein, the IFFT refers to inverse Fourier transform,representing wavenumber domain data obtained by high-dimensional Fourier transform of an original data block,refers to the wave number domain form obtained by the wavelet atom under the parameter group (b, a, theta, phi) through the high-dimensional Fourier transform,is composed ofConjugation of (1). The time complexity is O (N)³N_aN_θN_φ) (ii) a The time consumption is reduced to the original value by the acceleration methodUnder the condition of a fixed scale a, an inclination angle theta and an azimuth angle phi, the continuous wavelet transform can be rapidly calculated by using the formula (2), and the obtained continuous wavelet transform result is called a sub-band.

S4, determining the high-dimensional continuous wavelet transform coefficient sub-band required to be reserved in the step S3 according to the angle of coherent interference;

if the dip angle of coherent interference in the pre-stack seismic data is theta_{coherent_noise}Considering that the visual velocity of coherent interference is small, only the angle smaller than theta in CWT (f; b, a, theta, phi) is reserved_{coherent_noise}Sub-bands of, i.e. retaining only [ -theta [ - ]_{coherent_noise}+Δθ,θ_{coherent_noise}-Δθ]And (3) a high-dimensional continuous wavelet transform coefficient sub-band in an angle range (wherein delta theta is a smaller angle value for ensuring the angle filtering effect, and the value is generally 5 degrees).

S5, determining the threshold value in each subband needing to be reserved in the step S4 according to the noise standard deviation obtained in the step S2;

the noise standard deviation is obtained in the step S3 equation (2)On the basis, the threshold value thre (a) within each subband to be retained is determined:

wherein the content of the first and second substances,is a threshold coefficient of each sub-band, which takes the value of 1/a³Is in direct proportion; c is a constant that controls the degree of noise suppression, typically taken to be 3.

S6, using the threshold value of each sub-band needing to be reserved determined in the step S5Carrying out threshold processing on the wavelet coefficient in the subband needing to be reserved;

in the sub-band determined in step S4 to be reserved, according to the threshold determined in S5And (3) carrying out hard limiting filtering:

s7, carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band subjected to threshold processing in the step S6 to obtain a reconstructed signal with high signal-to-noise ratioFinishing noise reduction processing;

thresholding a set of wavelet coefficients CWT_thr(f, b, a, theta, phi) to obtain reconstructed high signal-to-noise ratio signal by fast inverse wavelet transform

Wherein, C_ψIs a constant associated with the wavelet function.

S8, repeating the steps S2-S7 until all the blocks are processed completely, and then obtaining the high signal-to-noise ratio signalAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

After all the data blocks are processed according to the flow, splicing the data body of each block with noise suppressed to obtain the noise suppressed data; due to edge degradation caused by Fourier transformation, a part of invalid data of the edge needs to be cut off; because different data blocks have certain differences in noise suppression effect due to different threshold selection, numerical calculation errors and the like, the data can generate obvious boundary lines during splicing; therefore, when the data are divided in step S1, there is a partial overlap between the data, and this partial overlap needs to be cut off by cutting in the splicing process in this step.

Taking the x direction as an example, assume a data block f₁And f₂If n points overlap for the data blocks adjacent in the x direction, the data block f₁After removal in the x-directionPoint of (2), data block f₂Before removal in the x directionAnd (4) splicing. The value of n needs to be determined by small batch experiments.

In another embodiment of the present invention, a deep pre-stack seismic data noise suppression system is provided, which can be used to implement the deep pre-stack seismic data noise suppression method described above, and specifically, the deep pre-stack seismic data noise suppression system includes a blocking module, a transformation module, an estimation module, a selection module, a threshold module, a processing module, a reconstruction module, and a suppression module.

The blocking module is used for acquiring original pre-stack data and blocking the original pre-stack data;

the transformation module estimates the noise variance grade of any original pre-stack seismic data after the blocking of the blocking module by using a median estimation operator;

the estimation module is used for carrying out high-dimensional continuous wavelet transform on the original pre-stack seismic data estimated by the transformation module to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform;

the selection module determines a high-dimensional continuous wavelet transform coefficient sub-band which needs to be reserved in the coefficient sub-band obtained by the estimation module according to the angle of the coherent interference;

a threshold module for determining the threshold value in the sub-band of the high-dimensional continuous wavelet transform coefficient to be reserved in the selection module according to the noise variance grade obtained in the transformation module

A processing module for determining the threshold value by using the threshold value moduleCarrying out threshold processing on wavelet coefficients in high-dimensional continuous wavelet transform coefficient sub-bands needing to be reserved;

a reconstruction module for performing high-dimensional wavelet inverse transformation on the wavelet transformation sub-band processed by the threshold value processing module to obtain a reconstructed high signal-to-noise ratio signal

A pressing module for processing all the blocks obtained by the block module and then obtaining all the high signal-to-noise ratio signalsAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

In yet another embodiment of the present invention, a terminal device is provided that includes a processor and a memory for storing a computer program comprising program instructions, the processor being configured to execute the program instructions stored by the computer storage medium. The Processor may be a Central Processing Unit (CPU), or may be other general purpose Processor, a Digital Signal Processor (DSP), an Application Specific Integrated Circuit (ASIC), an off-the-shelf Programmable gate array (FPGA) or other Programmable logic device, a discrete gate or transistor logic device, a discrete hardware component, etc., which is a computing core and a control core of the terminal, and is adapted to implement one or more instructions, and is specifically adapted to load and execute one or more instructions to implement a corresponding method flow or a corresponding function; the processor of the embodiment of the invention can be used for the operation of the deep prestack seismic data noise suppression method, and comprises the following steps:

acquiring original pre-stack data and partitioning; estimating the noise variance grade of any original pre-stack seismic data after the block division by using a median estimation operator; performing high-dimensional continuous wavelet transform on the estimated original pre-stack seismic data to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform; determining a high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved in a series of high-dimensional continuous wavelet transform coefficient sub-bands according to the angle of coherent interference; determining threshold value in high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved according to noise variance levelUsing thresholds within subbands of high-dimensional continuous wavelet transform coefficients that need to be retainedCarrying out threshold processing on the wavelet coefficient in the subband needing to be reserved; carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band after threshold processing to obtain a reconstructed high signal-to-noise ratio signalRepeating the above steps until all the blocks are processed completely, and then obtaining the high signal-to-noise ratio signalAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

In still another embodiment of the present invention, the present invention further provides a storage medium, specifically a computer-readable storage medium (Memory), which is a Memory device in a terminal device and is used for storing programs and data. It is understood that the computer readable storage medium herein may include a built-in storage medium in the terminal device, and may also include an extended storage medium supported by the terminal device. The computer-readable storage medium provides a storage space storing an operating system of the terminal. Also, one or more instructions, which may be one or more computer programs (including program code), are stored in the memory space and are adapted to be loaded and executed by the processor. It should be noted that the computer-readable storage medium may be a high-speed RAM memory, or may be a non-volatile memory (non-volatile memory), such as at least one disk memory.

One or more instructions stored in a computer-readable storage medium may be loaded and executed by a processor to implement the corresponding steps of the method for noise suppression of seismic data prior to deep prestack in the embodiments described above; one or more instructions in the computer-readable storage medium are loaded by the processor and perform the steps of:

acquiring original pre-stack data and partitioning; estimating the noise variance grade of any original pre-stack seismic data after the block division by using a median estimation operator; performing high-dimensional continuous wavelet transform on the estimated original pre-stack seismic data to obtain a series of coefficient sub-bands of high-dimensional continuous wavelet transform; determining a high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved in a series of high-dimensional continuous wavelet transform coefficient sub-bands according to the angle of coherent interference; determining threshold value in high-dimensional continuous wavelet transform coefficient sub-band needing to be reserved according to noise variance levelUsing thresholds within subbands of high-dimensional continuous wavelet transform coefficients that need to be retainedThresholding wavelet coefficients in subbands requiring preservationValue processing; carrying out high-dimensional wavelet inverse transformation on the wavelet transformation sub-band after threshold processing to obtain a reconstructed high signal-to-noise ratio signalRepeating the above steps until all the blocks are processed completely, and then obtaining the high signal-to-noise ratio signalAnd re-splicing into complete high signal-to-noise ratio pre-stack seismic data.

In order to make the objects, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are some, but not all, embodiments of the present invention. The components of the embodiments of the present invention generally described and illustrated in the figures herein may be arranged and designed in a wide variety of different configurations. Thus, the following detailed description of the embodiments of the present invention, presented in the figures, is not intended to limit the scope of the invention, as claimed, but is merely representative of selected embodiments of the invention. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.

Take an actual seismic data as an example.

Referring to FIG. 2, FIG. 2 shows pre-stack CMP gather data for a physical area where noise buries most of the available signal, the noise including a significant amount of random noise and a significant amount of directional structural noise. Referring to fig. 3 and 4, fig. 3 is data of noise suppressed by the method of the present invention, in which most of random noise and structural noise are removed and effective signals are revealed, which indicates that the method has high noise suppression efficiency for deep prestack seismic data, and fig. 4 is a difference section of fig. 3 and 2, in which no effective data structure is visible to the naked eye, which indicates that the method does not damage the effective signals. Referring to fig. 5, noise includes a large amount of random noise and a small amount of structural noise for a single overlay inline line before an actual offset. Referring to fig. 6 and 7, fig. 6 shows the data after the noise suppression according to the present invention, most of the noise is removed, and the signal structure originally buried can be seen clearly. Fig. 7 is a differential profile of fig. 5 and 6, excluding the effective signal structure observable to the naked eye.

In summary, the deep pre-stack seismic data noise suppression method and system provided by the invention have good suppression effects on random noise and directional structural noise in deep pre-stack seismic data, and do not damage effective signals. Compared with the prior art, the method can remove both random noise and structural noise.

As will be appreciated by one skilled in the art, embodiments of the present application may be provided as a method, system, or computer program product. Accordingly, the present application may take the form of an entirely hardware embodiment, an entirely software embodiment or an embodiment combining software and hardware aspects. Furthermore, the present application may take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, and the like) having computer-usable program code embodied therein.

The present application is described with reference to flowchart illustrations and/or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the application. It will be understood that each flow and/or block of the flow diagrams and/or block diagrams, and combinations of flows and/or blocks in the flow diagrams and/or block diagrams, can be implemented by computer program instructions. These computer program instructions may be provided to a processor of a general purpose computer, special purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, create means for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.

These computer program instructions may also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to function in a particular manner, such that the instructions stored in the computer-readable memory produce an article of manufacture including instruction means which implement the function specified in the flowchart flow or flows and/or block diagram block or blocks.

These computer program instructions may also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.

The above-mentioned contents are only for illustrating the technical idea of the present invention, and the protection scope of the present invention is not limited thereby, and any modification made on the basis of the technical idea of the present invention falls within the protection scope of the claims of the present invention.