Seismic profile denoising method based on improved laplace-filtered reverse time migration

Abstract

A seismic profile denoising method based on improved Laplace-filtered reverse time migration, comprising: collecting seismic data, and obtaining an imaging result after reverse time migration processing; differentiating and expressing second-order derivatives of a Laplacian operator in an x direction and a z direction as an eight-neighborhood second-order differential form, and then splitting the eight-neighborhood second-order differential form into four Laplacian filters comprising two opposite directions to obtain improved Laplacian filters; and filtering and denoising the imaging result after reverse time migration processing by using an improved Laplacian filter that is only related to a second-order differential in a z-axis direction. A Laplacian filtering operator that retains approximate normal direction components of second-order differentials while omitting approximate parallel direction components of the second-order differentials is selected for filtering processing, so that the in-phase axis continuity and amplitude characteristics of effective signals are maintained while better eliminating imaging noise generated by back-reflected waves.

Description

Seismic profile denoising method based on improved Laplace filter reverse time migration Technical Field

[0001] The present invention relates to a denoising method for seismic sections based on improved Laplace filtering reverse time migration, and is a method for solving the problems of event discontinuity and residual noise still existing on imaging sections after Laplace filtering is applied after migration, and belongs to the field of exploration geophysics. Background Art

[0002] Prestack depth migration algorithms can be categorized as ray-based and wave-equation-based. One-way wave-equation migration cannot accurately image steeply dipping interfaces, while the Kirchhoff integral method (Zhu and Lines, 1998) can only describe wave propagation in smooth media, making it difficult to obtain accurate imaging profiles of complex structures. Reverse time migration is a seismic data processing method that uses the two-way wave equation to migrate seismic reflection data, effectively capturing subsurface images that describe geological structures (Baysal et al., 1983; Loewenthal et al., 1983). It offers numerous advantages, including the ability to utilize various waveforms, including reflections, multiples, and reversals, and its independence from formation dip. Cross-correlation imaging is the most widely used and simple method for reverse time migration (RTM). It provides accurate dynamic characteristics for imaging (Claerbout, 1971). However, compared to single-way wave migration, cross-correlation imaging results in the presence of a significant amount of high-amplitude, low-frequency noise due to backreflections (Liu et al., 2010; Du et al., 2013). This noise is unavoidable for acoustic, elastic (Yan and Sava, 2008; Du et al., 2012), and anisotropic RTM (Zhang et al., 2009). This noise significantly impacts the quality of the imaged sections and further processing and interpretation. Backreflections are waves generated when the source or receiver wavefields are extended and encounter strong impedance interfaces. Cross-correlation with the normally propagating waves produces high-amplitude, low-frequency noise in the imaged sections.

[0003] Since the 1980s, numerous researchers at home and abroad have proposed various methods for suppressing low-frequency noise in reverse time migration. In 1983, Baysal proposed suppressing low-frequency noise in post-stack RTM by matching wave impedances to suppress reflections at interfaces, thereby obtaining a two-way reflectionless wave equation. Loewenthal et al. (1987) used a smoothed velocity model to reduce backreflections, but these two methods have not been very effective in pre-stack applications. More recently, an imaging method based on traveling wave separation has been proposed (Yoon and Marfurt, 2006; Liu et al., 2011). This method fully utilizes all wavefield types, separating upgoing and downgoing waves based on the geometric relationship between the incident and reflected waves. The imaging results are then cross-correlated between components in different vertical directions in the wavefield, effectively suppressing the low-frequency noise generated by backreflections. In 2015, Fei et al. (2015) proposed retaining only the upgoing wave of the source wavefield and the downgoing wave of the receiver wavefield, applying these two components to cross-correlation imaging to suppress low-frequency noise and migration artifacts. Wang et al. (2016) further proposed separating the wavefield into left-going, right-going, upgoing, and downgoing waves, and cross-correlating these components with the source and receiver wavefields, thereby effectively separating low-frequency noise and imaging crosstalk. Least-squares migration (Nemeth et al., 1999; Dai et al., 2011; Hu et al., 2016) is a development direction of conventional migration. Its core concept is to find an exact solution in the model space under linear inversion theory. This is equivalent to considering amplitude compensation and waveform correction on the basis of conventional migration structural imaging, thereby obtaining imaging results with better amplitude preservation and higher resolution. Later, the least-squares concept was combined with the reverse-time migration method to produce least-squares reverse-time migration (LSRM) (Dong et al., 2012; Zhang et al., 2015; Liu et al., 2016, 2017). Compared with traditional migration methods, LSRM provides reflection wave imaging results with more balanced amplitudes and higher resolution, and can eliminate migration noise. However, LSRM typically requires multiple iterations to obtain a converged final image, which places high demands on computational cost and speed. Accordingly, techniques such as Poynting vectors have also been introduced into LSRM (Yang and Zhang, 2018; Wang et al., 2021), improving the efficiency of the least-squares algorithm and achieving good results. Technical issues

[0004] The purpose of this invention is to provide a reverse time migration method based on an improved Laplace filter, aiming to solve the problems of event discontinuity and residual noise existing in Laplace filtering after migration.

[0005] The Laplace operator effectively eliminates imaging noise, but it also alters the amplitude and phase characteristics of the seismic profile events, necessitating additional correction. To obtain higher-quality imaging profiles, the Laplace operator can be modified to suppress low-frequency noise. Given that the Laplace operator cannot effectively preserve the seismic profile events and phase characteristics, this paper proposes an improved Laplace filtering algorithm that can suppress imaging noise while producing high-quality seismic profiles. Technical Solutions

[0006] A seismic profile denoising method based on improved Laplace filter reverse time migration is characterized by comprising the following steps:

[0007] 1) Collect seismic data using geophones;

[0008] 2) Use the reverse time migration method to obtain the source wavefield and the receiver wavefield, and apply the zero-delay cross-correlation imaging condition to obtain the imaging results after RTM processing:

[0009] The zero-delay cross-correlation imaging condition is expressed as

[0010] Where (x, z) represents the position of the imaging point, T max represents the maximum recording time, S and R are the source wavefield and receiver wavefield, and I(x,z) represents the imaging result after RTM processing;

[0011] 3) Filtering:

[0012] The imaging result is subjected to a Laplace filter to suppress the low-frequency noise after imaging. The Laplace filter is expressed as

[0013] Where I lap represents the imaging result after Laplace filtering, (x,z) represents the imaging point position, I(x,z) represents the imaging result after RTM processing, and Laplace operator is a second-order differential operator in n-dimensional Euclidean space. The Laplace operator is defined as

[0014] Substituting formula (1) into formula (2), we can get the imaging result after Laplace filtering, which is expressed as

[0015] The derivative of the discrete function is calculated in the form of difference. The Laplace operator filter of the discrete function is obtained by differentiating the second-order derivatives in the x and z directions, which is expressed as

[0016] I lap (x,z) = [I(x - 1,z - 1) + I(x - 1,z + 1)(five)

[0017] + I(x + 1,z - 1) + I(x + 1,z + 1)

[0018] + I(x + 1,z) + I(x - 1,z)

[0019] + I(x,z - 1) + I(x,z + 1)]

[0020] - 8I(x,z)

[0021] Split the eight - neighborhood second - order differential Laplace filter. Since the two opposite directions in the eight - neighborhood have the same mathematical meaning, it is split into four Laplace filters to obtain an improved Laplace filter:

[0022] Use any one of the improved Laplace filters shown in equations (six), (seven), (eight), and (nine) to operate on all the imaging results after RTM, so as to denoise the seismic profile.

[0023] In step 2), the cross - correlation imaging condition is modified using the normalization of the source wavefield, and the obtained imaging result is expressed as:

[0024] In step 2), the RTM processing includes conventional acoustic wave RTM or elastic wave RTM.

[0025] In step 3), since I1 lap and I3 lap are symmetric in mathematical meaning and have little significance for denoising the seismic profile, it is preferred to use the improved Laplace filter represented by I2 lap and I4 lap

[0026] In step 3), it is preferred to use the improved Laplace filter represented by I2 lap Beneficial effects

[0027] ​​Based on its application in reverse time migration imaging, this paper proposes and implements an improved Laplace filtering method. Reverse time migration is a two-way wave-based migration imaging method that boasts high precision imaging of complex structures and is widely used. Zero-delay cross-correlation imaging is a simple and effective imaging condition in reverse time migration. However, the cross-correlation between the backreflected wave and the detection wavefield in this method generates a large amount of strong-amplitude low-frequency noise. Among current low-frequency noise suppression methods, Laplace filtering works well and can better protect the effective signal. However, it suffers from issues such as discontinuous events and uneven amplitudes.

[0028] Therefore, we further proposed an improved Laplace filtering method. Unlike the Laplace filter, in the second-order differential approximation, we further deduced that only the components related to the second-order differential in the z-axis direction are considered, while those related to the x-axis direction are ignored. The components related to the x-axis second-order differential are more sensitive to tilted structures, while the components related to the z-axis second-order differential are more responsive to both parallel and tilted structures, further removing the offset artifact. Theoretical analysis in the wavenumber domain also demonstrated the feasibility of this method. BRIEF DESCRIPTION OF THE DRAWINGS

[0029] Figure 1 Schematic diagram of source wave field vector and receiver wave field vector;

[0030] Figure 2: Reverse time migration low-frequency noise generation mechanism;

[0031] Figure 3. Reverse time migration results under full wave equation conditions;

[0032] Figure 4 Improved Laplace filter I1 lap the result;

[0033] Figure 5 Improved Laplace filter I2 lap the result;

[0034] Figure 6 Improved Laplace filter I3 lap the result;

[0035] Figure 7 Improved Laplace filter I4 lap the result;

[0036] FIG8 is a flowchart of seismic profile denoising based on improved Laplace filter reverse time migration shown in Example 1;

[0037] Figure 9 Schematic diagram of the Sigsbee2 model;

[0038] FIG10 shows the denoised imaging result of the seismic profile based on the improved Laplace filter reverse time migration shown in Example 1. Best Mode for Carrying Out the Invention

[0039] The seismic profile denoising method based on improved Laplace filter reverse time migration includes the following steps:

[0040] 1) Collect seismic data using geophones;

[0041] 2) Use the reverse time migration method to obtain the source wavefield and the receiver wavefield. The conventional reverse time migration algorithm applies the zero-delay cross-correlation imaging condition to obtain the imaging (i.e., the imaging result after RTM):

[0042] The zero-delay cross-correlation imaging condition can be expressed as

[0043] Among them, (x, z) represents the position of the imaging point, T max represents the maximum recording time, S and R are the source wavefield and receiver wavefield, and I represents the imaging result after RTM. This imaging condition is simple to apply and provides stable structural images. However, the image amplitude of cross-correlation imaging has no clear physical relationship with the reflection coefficient and is arbitrarily scaled according to the source strength and receiver coverage. The cross-correlation imaging condition can be modified by normalizing the source wavefield, and the resulting imaging result can be expressed as:

[0044] During reverse time migration, strong-amplitude, low-frequency noise often contaminates the migration results, especially in shallow layers and above strong reflecting interfaces. When back reflections are developed, the cross-correlation between the source wavefield and the receiver wavefield generates noise due to the cross-correlation between the forward-propagating source wavefield and the back-reflected receiver wavefield, and vice versa. This noise exhibits low-frequency characteristics in the frequency domain. Low-frequency noise strongly interferes with the migration profile, and the greater the difference in wave impedance between interfaces, the stronger the back-reflected wave energy, and the stronger the low-frequency noise. Compared with one-way wave migration, the one-way wave migration method restricts the propagation direction of the wavefield. The wavefield of the forward extension of the shot point only considers downward propagation when encountering a reflecting surface. Therefore, imaging only occurs in the direction opposite to the propagation direction of the shot point and the receiver point, resulting in the absence of strong energy noise.

[0045] 3) Filtering:

[0046] After conventional acoustic wave RTM or elastic wave RTM is performed on the imaging result, the low-frequency noise after imaging is suppressed using a Laplace filter. The Laplace filter is expressed as

[0047] Where I lap represents the imaging result after Laplace filtering, (x, z) represents the position of the imaging point, I(x, z) represents the imaging result after RTM processing (i.e., the seismic profile after RTM), and the Laplace operator is a second-order differential operator in n-dimensional Euclidean space. The Laplace operator is defined as

[0048] Substituting formula (2) into formula (3), we can get the imaging result after Laplace filtering, which is expressed as

[0049] The derivative of a discrete function can be calculated in the form of a difference. The Laplace operator filter of the discrete function is obtained by differentiating the second-order derivatives of the Laplace operator in the x and z directions, which is expressed as

[0050] I lap (x,z)=[I(x-1,z-1)+I(x-1,z+1) (6)

[0051] +I(x+1,z-1)+I(x+1,z+1)

[0052] +I(x+1,z)+I(x-1,z)

[0053] +I(x,z-1)+I(x,z+1)]

[0054] -8I(x,z)

[0055] The improved Laplace filter is characterized by splitting the eight-neighborhood second-order differential Laplace filter. Since the two relative directions in the eight-neighborhood have the same mathematical meaning, it is split into four Laplace filters to obtain:

[0056] Any of the improved Laplace filters shown in equations (7), (8), (9), and (10) is used to operate on all RTM imaging results, thereby denoising the seismic profile.

[0057] Furthermore, if we analyze in the two-dimensional wavenumber domain, we can get the following through Fourier transform:

[0058] Where k is the imaging domain wave number vector, k x is the wave number along the x direction at the reflection point, k z is the wave number along the z direction at the reflection point;

[0059] Figure 1 is a schematic diagram of the source wavefield vector S and the receiver wavefield vector R of the wave, as well as the obtained imaging wavefield vector n, where θ is the angle between the source wavefield vector S and the imaging wavefield vector R.

[0060] The imaging domain wave number vector k is expressed as: k = k R -k S(12)

[0061] Where k R is the wave number vector of the detector wave field, k S is the source wave field wave number vector; using the cosine theorem, we can obtain:

[0062] Among them, v is the seismic wave velocity and ω is the angular frequency. It can be seen from the above formula that the Laplace filter after Fourier transform is an angle domain filter.

[0063] Figure 2 shows the receiver wavefield and source wavefield. The angle between the source wavefield vector and the receiver wavefield vector is 2θ. In conventional reverse time migration, images are taken at all angles. After Laplace filtering, an angle-dependent weight is added to the imaging result. The weight is large within the range of 0 to 60°, which effectively preserves the imaging result. When the angle θ between the source wavefield vector and the receiver wavefield vector is between 60 and 90°, the weight is relatively small, thereby filtering out large-angle imaging wavefields, i.e., noise.

[0064] From the above examples, we can see that the Laplace filter as an angle domain filter has the following characteristics:

[0065] In the step 2), due to I1 lap and I3 lap It is symmetrical in mathematical sense, but has little significance for denoising seismic profiles. lap and I4 lap The resulting modified Laplace filter,

[0066] I2 lap Converting to the wavenumber domain for analysis yields:

[0067] I4 lap Converting to the wavenumber domain for analysis yields:

[0068] Where α is the angle between the imaging wave field vector and the z direction, which is equivalent to adding the weighting factor cos to the wave number vector wave field after ordinary Laplace filtering. 2 α and sin 2 α, for complex underground strata, the angle α is generally small, usually in the range of 0 to 30°, and the weighting factor cos 2 α can well protect the formation imaging wave field vector, while sin 2 α will suppress the formation imaging wave field vector; on the other hand, the noise angle α is large, and the weighting factor cos 2 α can effectively suppress the noise imaging wave field vector, and correspondingly sin 2α will play a protective role against noise imaging wave field vector, so the weighted factor cos 2 α, that is, the filter corresponding to Equation (14), so the Laplace operator I2 lap The imaging noise can be further suppressed.

[0069] The improved Laplace filter reverse time migration method, that is, the four Laplace operators shown in equations (7) to (10), are applied to the imaging results. Figure 3 is filtered respectively, and the imaging results shown in Figures 4-7 are obtained:

[0070] All Laplace operators suppress the low-frequency noise generated by back reflection, but each Laplace operator responds differently to low-frequency noise, which is related to the second-order differential in the x-axis direction. lap (Figure 7) shows the migration image that is sensitive to the tilted structure, and I1 lap and I3 lap (Figures 4 and 6) show that the migration images are sensitive to the structures in two opposite directions. Compared with the ideal migration results, there is also serious crosstalk noise, which causes destructive interference in the final stacked image. I2 is only related to the second-order differential in the z-axis direction. lap (As shown in Figure 5) It has good imaging effects for both inclined and horizontal structures, and the low-frequency noise generated by back reflection is also well suppressed.

[0071] The use of the modified Laplace operator to suppress low-frequency noise relies on the second-order differential approximation. Unlike the Laplace operator, which consists of both parallel and normal terms, the modified Laplace operator selects only the parallel terms to form the Laplace operator. We performed reverse time migration of the fault model using cross-correlation imaging (Figure 3) and then used the modified Laplace operator (Figure 5). Comparing the two, the reverse time migration results using the modified Laplace operator clearly demonstrate the suppression of low-frequency noise generated by backreflections. Modes for Carrying Out the Invention

[0072] Referring to FIG8 , this example includes the following steps:

[0073] (1) The source wave field is simulated based on the wave equation and the boundary value is maintained at the imaging time point.

[0074] (2) Using the wave equation to reversely extend the detector wave field;

[0075] Steps (1) and (2) are both simulated using the second-order time and eighth-order space finite difference wave equations, and the obtained source wave fields are normalized.

[0076] The Sigsbee2 model has a depth of 2.6 km and a width of 7.6 km (see Figure 9). It uses a Ricker wavelet with a main frequency of 30 Hz and a time sampling interval of 1 ms. Shots are placed from 0 to 7500 m above the surface, with a shot spacing of 500 m, for a total of 16 shots. Receivers are evenly spaced on the surface at 10 m intervals, for a total of 750, with a total recording length of 4 s.

[0077] (3) The zero-delay cross-correlation imaging condition is Where, T max represents the maximum recording time, S and R are the source wavefield and the receiver wavefield.

[0078] (4) Set Laplace filter Where I represents the imaging result, and (x, z) represents the position of the imaging point.

[0079] (5) Differentiate the second-order derivatives of the Laplace operator in the x and z directions to obtain the discrete form of the Laplace filter I lap (x,z)=[I(x-1,z-1)+I(x-1,z+1)+I(x+1,z-1)+I(x+1,z+1)+I(x+1,z)+I(x-1,z)+I(x,z-1)+I(x,z+1)]-8I(x,z).

[0080] (6) Split the eight-neighborhood second-order differential Laplace operator and select the Laplace operator related to the horizontal second-order differential

[0081] (7) I2 lap The denoising result obtained by applying it to the RTM imaging results is shown in Figure 10. Industrial Applicability

[0082] Computational results using Marmousi model data demonstrate the effectiveness of the improved Laplace filter. While maintaining event continuity, it suppresses low-frequency noise generated by backreflections, resulting in imaging results with good amplitude balance and a higher signal-to-noise ratio. Applying the improved Laplace filter after reverse time migration is a simple and direct method. It filters the migrated data volume and is dependent only on the data volume's grid size. It offers the advantages of high computational efficiency and practicality.

Claims

1. A seismic profile denoising method based on improved Laplace filter reverse time migration, characterized in that The following steps are involved: 1) Collect seismic data using geophones; 2) Use the reverse time migration method to obtain the source wave field and the receiver wave field, and apply the zero-delay cross-correlation imaging condition to obtain the imaging results after RTM processing: The zero-delay cross-correlation imaging condition is expressed as Where (x, z) represents the position of the imaging point, T max represents the maximum recording time, S and R are the source wave field and the receiver wave field, and I(x,z) represents the imaging result after RTM processing; 3) Filtering: The imaging result is subjected to a Laplace filter to suppress the low-frequency noise after imaging. The Laplace filter is expressed as Where I lap represents the imaging result after Laplace filtering, (x,z) represents the imaging point position, I(x,z) represents the imaging result after RTM processing, and Laplace operator is a second-order differential operator in dimensional Euclidean space. The Laplace operator is defined as Substituting formula (1) into formula (2), we get the imaging result after Laplace filtering, which is expressed as The derivative of the discrete function is calculated in the form of difference. The Laplace operator filter of the discrete function is obtained by differentiating the second-order derivatives of the Laplace operator in the x and z directions, which is expressed as I lap (x,z) = [I(x - 1,z - 1)+I(x - 1,z + 1) (Five) +I(x+1,z-1)+I(x+1,z+1) +I(x+1,z)+I(x-1,z)+I(x,z-1) +I(x,z+1)]-8I(x,z) The eight-neighborhood second-order differential Laplace filter is split into four Laplace filters to obtain an improved Laplace filter because the two relative directions in the eight-neighborhood have the same mathematical meaning: The imaging result after the RTM processing is operated by using any one of the improved Laplace filters shown in equations (six), (seven), (eight) and (nine), so as to achieve denoising of the seismic profile.

2. The seismic profile denoising method based on improved Laplace filter reverse time migration as claimed in claim 1, characterized in that In step 2), the normalization of the source wave field is used to modify the cross-correlation imaging condition, and the obtained imaging result is expressed as:

3. The seismic profile denoising method based on improved Laplace filter reverse time migration as claimed in claim 1, characterized in that In the step 2), the RTM process includes conventional sonic wave RTM or elastic wave RTM.

4. The seismic profile denoising method based on improved Laplace filter reverse time migration as claimed in claim 1, characterized in that In the step 3), due to I1 lap and I3 lap It is symmetrical in mathematical sense, but has little significance for denoising seismic profiles. lap and I4 lap Represents the improved Laplace filter.

5. The seismic profile denoising method based on improved Laplace filter reverse time migration as claimed in claim 4, characterized in that In the step 3), preferably, I2 lap Represents the improved Laplace filter.

Images