Linear noise suppression method and device based on rotation time shift FK transformation

Abstract

The invention relates to a linear noise suppression method and device based on rotation time shift FK transformation, and the method comprises the steps: leveling a linear event in a rotation manner based on the idea of multi-focus separation, and then carrying out the random time shift of each time point in the horizontal direction; frequency differences are formed between events similar to linear noise in apparent velocity due to random time shifting of data. And finally, carrying out linear noise suppression in the frequency-wave number FK domain by applying a compressed sensing reconstruction algorithm. According to the technical scheme, random time shifting is applied, energy near the dominant frequency of the FK domain is evenly distributed to all frequencies, the sparseness of a transform domain is enhanced, and therefore the algorithm can effectively suppress linear noise and reduce event energy loss with the apparent velocity and the frequency similar to those of the linear noise. The technical scheme is suitable for seismic data prestack linear noise suppression, powerful technical means support is provided for high-fidelity and high-signal-to-noise-ratio processing, and the problem that seismic data generally have a low signal-to-noise ratio is solved.

Description

Technical Field

[0001] The invention relates to the technical field of oil and gas exploration, and in particular to a linear noise suppression method and device based on rotational time-shift FK transformation. Background Art

[0002] Linear noise interference is a long-standing problem in seismic data processing and subsequent imaging, resulting in unsatisfactory continuity and quality of the event axis. Therefore, in recent years, many methods have been developed to suppress seismic linear and quasi-linear noise. Traditional methods can be roughly divided into two categories: filtering and wavefield conversion. The filtering method transforms the seismic data to the corresponding domain according to the difference between the signal and the noise in the transform domain, and designs the corresponding filter to attenuate the noise. For example, Yilmaz proposed fk filtering based on the difference in frequency and wavenumber between the signal and the noise to suppress linear noise. Turner used Radon transform to filter linear noise and also suppressed random noise. Douglas used correlation coefficient filtering to separate non-correlated noise from correlated signals on the basis of each frequency. Lu et al. proposed a local linear noise attenuation method based on local polynomial approximation. Henley proposed a coherent noise attenuation technique based on radial trace transformation, which can effectively attenuate noise while retaining the lateral details between seismic traces. Another type of method is based on wavefield transformation, flattens the linear noise event axis, and uses the corresponding transformation to fit the target linear noise. Then, the target linear noise is inversely transformed and then subtracted from the original seismic record to obtain a noise-free seismic record. For example, Yedlin et al. used the Karhunen-Loève (KL) transform to eliminate the correlation between the signal components and to separate the linear noise. Lu et al. proposed an adaptive singular value decomposition (SVD) to enhance the non-horizontal event axis and then flattened the linear noise by the estimated slope and data rotation. Huang et al. used the morphological difference between the signal and the coherent noise in the spatial direction to apply mathematical morphological filtering (MMF) to linear noise attenuation.

[0003] However, after suppressing seismic linear and quasi-linear noise, the above methods still have the problem of low signal-to-noise ratio of seismic data, and there is currently no method to suppress linear noise and improve the signal-to-noise ratio of seismic data by using rotational time-shift FK transform. Summary of the invention

[0004] The purpose of the present invention is to propose a linear noise suppression method and device based on rotational time-shift FK transform, aiming to solve the technical problem that the seismic data still has a low signal-to-noise ratio after suppressing seismic linear and quasi-linear noise in the prior art.

[0005] To achieve the above object, the present invention proposes a linear noise suppression method based on rotational time-shift FK transform, comprising the following steps:

[0006] S1: Identify linear noise and quasi-linear noise in seismic data and determine apparent velocity and distribution range;

[0007] S2: Use the event tracking method to track the shape of the quasi-linear noise and perform flattening processing;

[0008] S3: Rotate according to the apparent velocity of the linear noise event axis to make the linear noise level;

[0009] S4: Set random numbers and randomly rearrange the data containing linear noise horizontally according to the random numbers;

[0010] S5: Perform Fourier transform along the horizontal and vertical directions to identify the noise distribution range and characteristics in the FK domain;

[0011] S6: Using the compressed sensing algorithm, threshold processing is performed on the distribution area of ​​the nonlinear noise in the FK domain;

[0012] S7: Inversely transform the processed data to the time domain, and use the random number of step S4 to perform reverse rearrangement, and rotate in the opposite direction to the original angle, perform reverse flattening processing, and obtain linear noise;

[0013] S8: Subtract the noisy data from the obtained linear noise to obtain seismic data after suppressing the linear noise.

[0014] Furthermore, in step S1, the seismic data is pre-stack seismic data, and linear noise is identified through the interactive function of GeoEast software. The linear noise is linearly or quasi-linearly distributed on the pre-stack data, which is significantly different from the hyperbolic event axis.

[0015] Furthermore, in step S2, the tracking direction of the event axis is determined according to the similarity between adjacent seismic records.

[0016] Further, in step S3, the rotation method is:

[0017]

[0018] S0 is the original seismic data containing linear noise, and S is defined as part of the data of S0;

[0019]

[0020] S contains linear noise. After testing, the more accurate S is, the better the linear noise suppression effect is. The rotation flattening operator is defined as Different from plane wave prediction (PWD), the rotation flattening operator used in this paper is The overall seismic data is rotated to make the linear noise horizontal while keeping the relative position of the linear noise and the other events unchanged. This can better protect the other events while suppressing the linear noise. The rotation and flattening operator corresponds to a mapping, namely:

[0021]

[0022]

[0023] Assume that the linear noise in S has a constant apparent velocity v a , the linear noise inclination angle can be calculated as θ = arctan(l / v a ), then the data size after rotation is updated to

[0024]

[0025] The rotated coordinates can be described as:

[0026]

[0027] For the reverse rotation, replace θ with θ, i.e.

[0028]

[0029] The rotation and flattening operators shown in formula (3) and formula (4) can be obtained: The positive and negative transformation.

[0030] Furthermore, the principle of step 4 is that after the linear noise is rotated and flattened, the non-horizontal event axis is treated as noise. An effective idea is that the linear noise is approximately horizontal, and the main energy is concentrated near the zero wavenumber in the FK domain. The non-horizontal event axis energy is distributed at each wavenumber and frequency. Sparse promotion transformation can be applied to further deteriorate the energy concentration of the remaining non-horizontal event axes. According to the revelation that the secondary sources of aliased data are randomly distributed in the common detection point domain, multiple sources can be separated according to the sparsity of the main source energy in the sparse transformation domain. Similar to this idea, we apply a horizontal random time-shift operator to shuffle the data at each time point, shuffle the non-horizontal event axes of some frequency bands into the full frequency band, reduce the sparsity of the data, and then apply the compressed sensing principle in the FK domain to extract linear noise;

[0031] Finally, the data after linear noise suppression can be obtained by subtracting the original data s. Based on the above, in step S4, the lateral time shift operator P = FΓF -1 , where F and F -1 represents the forward and reverse Fourier transform, and FF -1 =I. The block matrix Γ represents the time shift and can be described as

[0032]

[0033] Among them, Γ i (i=1, ..., m′) represents the phase shift operator of the i-th channel of seismic data, in the form of a diagonal matrix, where ω i (i=1, ..., n′) represents the circular frequency τ i represents the random time shift. The flattened data S′ needs to be vectorized horizontally, denoted by S′ r = Row(S′), and Similarly, the inverse transverse vectorization is denoted as S′=IRow(S′ r ), then the data after horizontal random time shift can be expressed as Ignoring horizontal vectorization and inverse quantization, the rotation time shift operation (RTS) of the data S can be written as:

[0034] Further, the purpose of the random time shift of the distribution range of step 5 can be used Figure 2 Give a description. Figure 2 There are three event axes on a, S1 and S2 are horizontally distributed, S3 is in a hyperbolic shape, and its single-channel waveform is shown by the black arrow; it is considered that S2 is the linear noise to be extracted, and S1 is an effective event axis with the same slope as S2. Figure 2 a is subjected to horizontal random time shift. Figure 2 c and the black arrows, S2 is still linear, while S1 and S3 have become random noises similar to “co-detection point secondary sources”. Figure 2 b and Figure 2 d), the green and red boxes are the gathering places of S2 and part of S1 energy (horizontal phase axis, wave number is approximately 0). It can be observed that the red box ( Figure 2 The energy focusing of d) is much higher than that of the green box ( Figure 2 b), which increases the sparsity of linear noise in the transform domain and is beneficial to the extraction of linear noise.

[0035] Furthermore, in step S6, the compressed sensing algorithm extracts linear noise based on a random noise suppression algorithm of compressed sensing theory, and the extraction equation is as follows:

[0036]

[0037] In order to simplify the equation, the lateral vectorization and inverse quantization steps are ignored. B represents a transformation with sparse capability. Since the linear noise after flattening has an infinite apparent velocity and the wave number is approximately 0, B is selected as the FK transformation in this paper. The corresponding sparse promotion transformation is It is called RTS-FK transformation.

[0038] Furthermore, the inverse transformation operation of step 7 needs to be determined according to steps 2 to 5. For example, if a 30° clockwise rotation is performed in step 3, a 30° counterclockwise rotation is required in step 7.

[0039] Furthermore, the target equation is solved by using an iterative soft threshold algorithm, that is, repeatedly iterating the following equation:

[0040]

[0041] Among them, α represents the iteration step size, T λ represents the soft threshold function in the time domain, θ is a constant, and x represents the input value:

[0042]

[0043] Furthermore, the step S8 can be performed by subtracting the remaining linear noise S from the noisy seismic record S. linear to obtain a noise-free signal.

[0044] The present invention also proposes a linear noise suppression device based on rotation time-shift FK transformation, comprising:

[0045] A memory, a processor, and a linear noise suppression program stored in the memory and executable on the processor, wherein the linear noise suppression program is executed by the processor to perform the steps of the linear noise suppression method.

[0046] The present invention also proposes a computer-readable storage medium, on which a linear noise suppression program is stored. When the linear noise suppression program is executed by a processor, the steps of the above-mentioned linear noise suppression method are implemented.

[0047] Compared with the prior art, the present invention has the following beneficial effects:

[0048] This technical solution uses random time shift to evenly distribute the energy near the main frequency in the FK domain to each frequency, enhancing the sparsity of the transform domain. Therefore, this algorithm can effectively suppress linear noise and reduce the energy loss of the event axis whose apparent velocity and frequency are similar to linear noise. This technical solution is suitable for pre-stack linear noise suppression of seismic data, providing a powerful technical means to support high-fidelity and high signal-to-noise ratio processing, and solving the problem of low signal-to-noise ratio commonly found in seismic data. BRIEF DESCRIPTION OF THE DRAWINGS

[0049] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on the structures shown in these drawings without paying creative work.

[0050] Figure 1 It is a schematic diagram of the process of the present invention;

[0051] Figure 2 Schematic diagram of the lateral random time shift of the present invention, wherein (a) original data; (b) Figure 2 a) FK spectrum; (c) lateral random time shift of data; (d) Figure 2 FK spectrum of c;

[0052] Figure 3 This is a synthetic data test of the present invention, wherein: (a) data without linear noise; (b) synthetic data with linear noise (S1-S5);

[0053] Figure 4 The synthetic data denoising test diagram of the present invention, wherein (a) is the Radon transform denoising test; (b) is the FK filtering algorithm after rotation; (c) is the RTS-FK-CS algorithm of the present invention; (d) (e) and (f) are respectively Figure 4 a, 4b, and 4c correspond to the linear noise suppressed by the algorithm;

[0054] Figure 5 The data of the present invention after rotation at a speed of S1 and random time shift;

[0055] Figure 6 To test the noise resistance of the algorithm of the present invention;

[0056] Figure 7 Actual seismic data for testing the present invention;

[0057] Figure 8 Seismic data after linear noise is suppressed by the algorithm of the present invention;

[0058] Fig. 9 The linear noise suppressed by the algorithm of the present invention and the local enlarged view of the linear noise; wherein, (a) the linear noise suppressed by the algorithm; (b) the local enlarged view of the linear noise. DETAILED DESCRIPTION

[0059] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0060] It should be noted that if the embodiments of the present invention involve directional indications (such as up, down, left, right, front, back, etc.), the directional indications are only used to explain the relative position relationship, movement status, etc. between the components in a certain specific posture. If the specific posture changes, the directional indication will also change accordingly.

[0061] In addition, if there are descriptions involving "first", "second", etc. in the embodiments of the present invention, the descriptions of "first", "second", etc. are only used for descriptive purposes and cannot be understood as indicating or implying their relative importance or implicitly indicating the number of technical features indicated. Therefore, the features limited to "first" and "second" may explicitly or implicitly include at least one of the features. In addition, if "and / or" or "and / or" appears in the full text, its meaning includes three parallel solutions. Taking "A and / or B" as an example, it includes solution A, solution B, or solutions that satisfy both A and B. In addition, the technical solutions between the various embodiments can be combined with each other, but it must be based on the ability of ordinary technicians in this field to implement. When the combination of technical solutions is contradictory or cannot be implemented, it should be deemed that such combination of technical solutions does not exist and is not within the scope of protection required by the present invention.

[0062] In this embodiment, the present invention synthesized Figure 3 The pre-stack data shown in a contains 256 channels, 256 time sampling points, a sampling interval of 2ms, and a trace interval of 6.25m. The area where the linear noise exists in step 1 is clearly visible, and the apparent velocity is known: the apparent velocity of the linear noise S1 is close to S0 (S1: 3906m / s, S0: 4464m / s), and the apparent velocity of S2-S4 gradually decreases (S2: 3125m / s, S3: 1563m / s, S4: 1093m / s). S5 is used to simulate the linear noise in the opposite direction. The traditional FK-based method (FK filtering after rotation and flattening, abbreviated as R-FK-Filtering), the RTS-FK-CS proposed in the present invention and the linear Radon transform are used to attenuate the linear noise, and the structural similarity (SSIM) is selected to measure the denoising performance. The denoising results are shown in the figure. Figure 4 shown.

[0063] The step 2 is invalid in this embodiment because the noise in this embodiment is purely linear and does not require the tracking of the event axis.

[0064] Apply steps 3 and 4 to rotate and randomly shift the linear noise, corresponding to Figure 5 As shown, it can be seen that the linear noise is still linear, while the effective signal and random noise both become random noise.

[0065] exist Figure 4 In a and 4d, transformation noise due to the inverse Radon transform is generated. Due to the poor energy concentration in the Radon domain, many valid signals are contained in the noise. In addition, S5 in the opposite direction cannot be completely removed using the linear Radon transform. Denoising results ( Figure 4 a) and no linear noise data ( Figure 3 a) The SSIM between them is 0.48. The test is repeated 10 times with a denoising step of 0.3 seconds.

[0066] like Figure 4 As shown in b and 4e, R-FK-Filtering can attenuate most of the linear noise. However, without random time shift, it is difficult to distinguish the linear noise S1 in the FK domain from the first arrival S0, resulting in a serious loss of the first arrival wave of the FK transform. Denoising result ( Figure 4 b) and non-linear noise input data ( Figure 3 a) The SSIM between them is 0.61. The test is repeated 10 times and the denoising step takes 1.5 seconds.

[0067] exist Figure 4 c and 4f, RTS-FK-CS improves the transform domain sparsity by random time shift, which helps to identify S1 and the first arrival S0 ( Figure 5 , S1 is in the red box, and S0 is randomized by random time shift). RTS-FK-CS can better protect the first arrival wave S0 with a speed close to S1, while attenuating the linear noise S1. Denoising result ( Figure 4 c) and noise-free input data ( Figure 3 a) has an SSIM of 0.77. Repeating the test 10 times, the denoising step takes 3.3 seconds. The key parameter is the threshold in the hard threshold function. Since the linear noise energy using the RTS-FK transform is concentrated near the zero wavenumber, we usually set the threshold to 0.1 times the maximum value.

[0068] In order to test the robustness of the algorithm, in this example, random noise is added to the model (after low-pass filtering), and the amplitude of the noise is 1%-100% of the maximum amplitude of the noise-free data. As the amplitude of the random noise increases, the similarity (SSIM) between the linear noise suppressed by the three methods and the actual linear noise decreases. However, the SSIM of RTS-FK-CS is still higher than that of RTS-FK-CS and linear Radon transform, which shows that RTS-FK-CS is also quite robust to random noise ( Figure 6 ).

[0069] This embodiment also uses actual seismic data to perform linear noise suppression tests, such as Figure 7 The figure shows selected seismic data containing linear noise. The deep layer in the data is covered by linear noise. After being processed by the method of the present invention, the deep linear noise is suppressed ( Figure 8 ), and there is almost no effective signal in the suppressed random noise ( Fig. 9 ).

[0070] The invention processes the time domain of the signal and combines it with petroleum geological geophysical exploration, belonging to a method for suppressing noise in petroleum geophysical exploration.

[0071] Through this scenario example, it is verified that the linear noise suppression method and device based on rotational time-shift FK transformation provided by the implementation mode of the present application can indeed solve the problem of current data being contaminated by linear noise.

[0072] The above descriptions are only optional embodiments of the present invention, and are not intended to limit the patent scope of the present invention. All equivalent structural changes made using the contents of the present invention's specification and drawings, or directly / indirectly applied in other related technical fields, are included in the patent protection scope of the present invention.

Claims

1. A linear noise suppression method based on rotational time-shift FK transform, characterized in that: The steps include: S1: Identify linear noise and quasi-linear noise in seismic data and determine apparent velocity and distribution range; S2: Use the event tracking method to track the shape of the quasi-linear noise and perform flattening processing; S3: Rotate according to the apparent velocity of the linear noise event axis to make the linear noise level; S4: Set random numbers and randomly rearrange the data containing linear noise horizontally according to the random numbers; S5: Perform Fourier transform along the horizontal and vertical directions to identify the noise distribution range and characteristics in the FK domain; S6: Using the compressed sensing algorithm, threshold processing is performed on the distribution area of ​​the nonlinear noise in the FK domain; S7: Inversely transform the processed data to the time domain, and use the random number of step S4 to perform reverse rearrangement, and rotate in the opposite direction to the original angle, perform reverse flattening processing, and obtain linear noise; S8: Subtract the noisy data from the obtained linear noise to obtain seismic data after suppressing the linear noise.

2. The linear noise suppression method based on rotational time-shift FK transform according to claim 1 is characterized in that: In the step S1, the seismic data is pre-stack seismic data.

3. The linear noise suppression method based on rotational time-shift FK transform according to claim 1, characterized in that: In step S2, the tracking direction of the event axis is determined according to the similarity between adjacent seismic records.

4. The linear noise suppression method based on rotational time-shift FK transform according to claim 1 is characterized in that: In step S3, the rotation method is: S0 is the original seismic data containing linear noise, and S is defined as part of the data of S0; S contains linear noise, and the rotation and flattening operator is defined as The rotation and flattening operator corresponds to a mapping, namely: Assume that the linear noise in S has a stable apparent velocity υ a , the linear noise inclination angle can be calculated as θ = arctan(1 / υ a ), then the data size after rotation is updated to The rotated coordinates can be described as: For the reverse rotation, replace θ with -θ, i.e. The rotation and flattening operators shown in formula (3) and formula (4) can be obtained: The positive and negative transformation.

5. The linear noise suppression method based on rotational time-shift FK transform according to claim 4 is characterized in that: In step S4, define the transverse time shift operator P = FΓF -1 , where F and F -1 represents the forward and reverse Fourier transform, and FF -1 =I, the block matrix Γ represents the time shift, which can be described as Among them, Γ i (i=1, ..., m′) represents the phase shift operator of the i-th channel of seismic data, in the form of a diagonal matrix, where ω i (i=1, ..., n′) represents the circular frequency τ i represents the random time shift, because The flattened data S′ needs to be vectorized horizontally, denoted by S′ r = Row(S′), and Similarly, the inverse transverse vectorization is denoted as S′=IRow(S′ r ), then the data after horizontal random time shift can be expressed as 6. The linear noise suppression method based on rotational time-shift FK transform according to claim 5 is characterized in that: In step S6, the compressed sensing algorithm extracts linear noise based on a random noise suppression algorithm of compressed sensing theory, and the extraction equation is as follows:

7. The linear noise suppression method based on rotational time-shift FK transform according to claim 6 is characterized in that: The target equation is solved by using the iterative soft threshold algorithm, that is, repeatedly iterating the following equation: Among them, α represents the iteration step size, T λ represents the soft threshold function in the time domain, λ is a constant, and x represents the input value:

8. The linear noise suppression method based on rotational time-shift FK transform according to claim 7 is characterized in that: The step S8 can be performed by subtracting the remaining linear noise S from the noisy seismic record S linear to obtain a noise-free signal.

9. A linear noise suppression device based on rotational time-shift FK transform, characterized in that: The linear noise suppression device includes: a memory, a processor, and a linear noise suppression program stored in the memory and executable on the processor. When the linear noise suppression program is executed by the processor, the steps of the linear noise suppression method based on rotational time-shift FK transformation as described in any one of claims 1 to 8 are implemented.

10. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a linear noise suppression program, which, when executed by a processor, implements the steps of the linear noise suppression method based on rotational time-shift FK transformation as claimed in any one of claims 1 to 8.

Images