A method and device for directly suppressing multiple waves between onshore layers in the depth domain of a well control area

Abstract

The present disclosure relates to a method and apparatus for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area. The method comprises: determining a pre-stack CRP gather from which primary reflection wave information is removed based on acquired work area logging data and a stacked profile of pre-stack depth migration seismic results; calibrating the positions of the multiple waves; determining a stacked profile of the multiple waves based on the calibrated CRP gather; subtracting the stacked profile of the pre-stack depth migration seismic results from a stacked profile of multiple waves with a depth unit of m; and performing pulse deconvolution on the stacked profile after the subtraction to obtain a pre-stack depth migration result after direct suppression of multiple waves between onshore layers in the depth domain of the well-controlled area. This method solves the problem that in the past, when suppressing inter-layer multiple waves through Radon transform, the depth domain migration processing results need to be converted to the time domain through time-depth conversion. Distortion occurs during the time-depth conversion process, resulting in the obtained time domain gather having an inaccurate depiction of the onshore inter-layer multiple waves, which is not conducive to subsequent suppression processing.

Description

Technical Field

[0001] The present disclosure relates to the technical field of seismic data processing, and in particular to a method and device for directly suppressing multiple waves between onshore layers in a depth domain of a well-controlled area. Background Art

[0002] Improving the fidelity of results in seismic data processing is the primary goal. Onshore seismic data usually contain interlayer multiple waves, which seriously reduce the fidelity of processing results and even affect earthquake prediction results. There are many methods to remove multiple waves. The most commonly used is to apply Radon transform to suppress multiple waves in the τ-p domain. Radon transform was proposed by Austrian mathematician Radon in 1917. Claerbout was the first to introduce τ-p transform into seismic exploration; Hampson first proposed parabolic Radon transform to suppress multiple waves on CMP data sets; Kostov least squared the parabolic Radon transform so that its transform operator has a Toeplitz structure and can be solved quickly; Wang applied Radon transform domain adaptive filtering to improve the resolution of Radon transform; NG used Gaussian iteration method to calculate similarity coefficient weighting, and simultaneously sparsed the time axis and offset axis, which significantly improved the resolution of Radon transform; Sacchi and Ulrych introduced nonlinearity in the operator of least squares parabolic Radon transform. Regularized diagonal matrices and the proposed high-resolution Radon transform address the problem of data divergence in the Radon domain and improve the ability to distinguish between multiples and primaries in the Radon domain. Mauricio et al. proposed the conjugate gradient method to rapidly solve the high-resolution Radon transform. Herrmann et al. proposed the anti-aliasing high-resolution parabolic Radon transform to address both resolution and spatial aliasing issues. Xiong Deng et al. applied these theories to the attenuation of multiples, achieving good results. Liu Xiwu et al. used the least squares inversion method to study the parabolic Radon transform and the hyperbolic Radon transform, and proposed an algorithm for the sparse constrained preconditioned high-resolution Radon transform, which effectively suppressed multiples. Most of the above Radon transform-based suppression methods are related to seismic frequency, especially the high-precision Radon transform, which performs the Radon transform and inverse transform on a frequency-by-frequency basis in the frequency domain.

[0003] In processing a large amount of real-world onshore seismic data, we found that prestack depth migration is more capable of depicting complex structures than prestack time migration. Both the profiles and depth-domain CRP gathers produced by prestack depth migration have higher fidelity. However, the aforementioned methods for suppressing multiples are all implemented in the time domain. The depth migration results must be converted from time to depth using a velocity model. This time-depth conversion process introduces various distortions, resulting in time-domain gathers that are less accurate in depicting onshore interlayer multiples than depth-domain gathers. Furthermore, if the depth migration results are crudely converted from meters to milliseconds and then subjected to a Radon transform, frequency distortion and significant distortion in the inverse transform will occur. Therefore, simply processing depth migration results using time-domain methods is not feasible. Therefore, a method for suppressing onshore interlayer multiples in the depth domain is urgently needed. Summary of the Invention

[0004] The present disclosure proposes a method and device for directly suppressing onshore interlayer multiple waves in the depth domain of a well-controlled area, in order to solve the problem that, when suppressing interlayer multiple waves by Radon transform in the past, the depth domain offset processing results need to be converted to the time domain by time-depth conversion, and distortion occurs during the time-depth conversion process, resulting in the obtained time domain data set not being accurate enough in depicting the onshore interlayer multiple waves, which is not conducive to subsequent suppression processing.

[0005] According to one aspect of the present disclosure, a method for directly suppressing onshore interlayer multiple waves in a well-controlled area in a depth domain is provided, comprising:

[0006] Obtain the well logging data and pre-stack depth migration seismic stacked profiles of the work area. The depth units of the well logging data and pre-stack depth migration seismic stacked profiles are both in meters.

[0007] Determine typical wells in the work area, and based on the well logging curves of the typical wells, use the forward modeling formula considering the propagation effect to obtain pre-stack CRP gathers in the depth domain with the unit of m, removing the primary reflection wave information;

[0008] Calibrate the positions of the multiple waves on the prestack CRP gather from which the primary reflection wave information is removed, and replace the depth unit m of the calibrated multiple wave prestack CRP gather with ms;

[0009] For the calibrated multiple wave CRP gathers with depth unit of ms, NMO flattening and removal of the far-track stretching part are performed through velocity analysis to form a multiple wave superposition profile;

[0010] Replacing the depth unit of the multiple wave stack section with ms by m, subtracting the pre-stack depth migration seismic result stack section from the multiple wave stack section with a depth unit of m to obtain a subtracted stack section, and replacing the depth unit of m in the subtracted stack section with ms;

[0011] Pulse deconvolution processing is performed on the subtracted stacked sections with a depth unit of ms to obtain pre-stack depth migration results after direct suppression of onshore interlayer multiple waves in the depth domain of the well-controlled area.

[0012] Preferably, the method for determining typical wells in a work area comprises:

[0013] Determine that the well with the deepest vertical depth among all the wells in the work area is the typical well.

[0014] Preferably, the method of obtaining a prestack CRP gather with a depth unit of m and primary reflection wave information removed based on the well logging curve of the typical well by using a forward modeling formula considering propagation effect comprises:

[0015] Substituting the P-wave and S-wave velocities of the logging curve of the typical well into the forward modeling formula considering the propagation effect, a depth domain wave equation forward modeling full wavefield CRP gather in units of m is obtained;

[0016] Substituting the typical well acoustic wave curve into the forward modeling formula considering propagation effect, ignoring all propagation effects except the primary reflection, the depth domain primary reflection wave CRP gather in units of m is obtained by forward modeling;

[0017] The full-wavefield CRP gathers are forward modeled using the depth-domain wave equation, and the depth-domain primary reflection wave CRP gathers are subtracted to obtain the pre-stack CRP gathers without the primary reflection wave information.

[0018] Preferably, the method for calibrating the positions of multiple waves on a pre-stack CRP gather from which primary reflection wave information has been removed comprises:

[0019] If the event on the pre-stack CRP gather with the primary reflection wave information removed is bent, the depth position of the event corresponding to the offset distance of 0 is marked as the multiple wave existence position, thereby obtaining a calibrated multiple wave pre-stack CRP gather.

[0020] Preferably, the method of performing NMO flattening and removing the far-track stretched portion on the calibrated multiple wave CRP gathers with a depth unit of ms through velocity analysis to form a multiple wave superposition profile comprises:

[0021] determining weighted stacked seismic traces based on the calibrated multiple wave CRP gathers;

[0022] Performing velocity analysis on the weighted stacked seismic traces to find the NMO velocity;

[0023] flattening the calibrated multiple pre-stack CRP gathers using the NMO velocity to obtain a flattened CRP gather;

[0024] The far-path stretching part is removed from the flattened CRP gather to obtain the removed CRP gather;

[0025] The CRP gathers after the excision are horizontally stacked to form a multiple wave stacking section.

[0026] Preferably, the method of performing a process of removing the far-track stretching portion on the flattened CRP gather to obtain the removed CRP gather comprises:

[0027] For the flattened CRP gather, the first predetermined percentage and the second predetermined percentage of the offset are cut off to obtain the multiple wave superposition profile.

[0028] Preferably, before replacing the depth unit m of the subtracted stacked section with ms, the method further includes: performing quality control on the subtracted stacked section, the method comprising:

[0029] A cross-correlation calculation is performed on the superposition profile after the subtraction and the multiple wave superposition profile to determine whether the cross-correlation result is less than a third predetermined percentage. If so, the quality control passes; otherwise, the quality control fails.

[0030] Preferably, the length m of the inverse filter operator of the pulse deconvolution processing is: if the total trace length of the seismic records in the work area is less than or equal to 10,000 ms, the length of the inverse filter operator is 400 ms; if the total trace length of the seismic records in the work area is greater than 10,000 ms, the length of the inverse filter operator is 600 ms.

[0031] Preferably, the method further comprises: performing quality control on the obtained prestack depth migration results, wherein the method comprises:

[0032] The depth unit ms in the prestack depth migration result is replaced by m;

[0033] A cross-correlation calculation is performed on the pre-stack depth migration result with a depth unit of m and the multiple wave stacking section to determine whether the cross-correlation result is less than a fourth predetermined percentage. If so, the quality control passes; otherwise, a new round of iterative pulse deconvolution processing is performed on the stacking section after the subtraction until the cross-correlation result between the pre-stack depth migration result and the multiple wave stacking section is less than the fourth predetermined percentage.

[0034] According to one aspect of the present disclosure, a device for directly suppressing multiple waves between onshore layers in a well-controlled area in a depth domain is provided, comprising:

[0035] The acquisition unit is used to obtain the well logging data and the stacked section of the pre-stack depth migration seismic results in the work area. The depth unit of the well logging data and the unit of the stacked section data of the pre-stack depth migration seismic results are both in meters.

[0036] The multiple wave pre-stack CRP gather determination unit is used to determine the typical wells in the work area, and obtain the pre-stack CRP gathers with the primary reflection wave information removed in the depth domain unit of m based on the logging curves of the typical wells using the forward modeling formula considering the propagation effect;

[0037] a multiple wave calibration unit, configured to calibrate the positions of the multiple waves on the prestack CRP gather from which the primary reflection wave information has been removed, and replace the depth unit m of the calibrated multiple wave prestack CRP gather with ms;

[0038] The multiple wave superposition profile determination unit is used to perform NMO flattening and remove the far-track stretching part of the calibrated multiple wave CRP gather with a depth unit of ms through velocity analysis to form a multiple wave superposition profile;

[0039] a primary reflection wave stacking profile determining unit, configured to replace ms in the depth unit of the multiple wave stacking profile with m, subtract the pre-stack depth migration seismic result stacking profile from the multiple wave stacking profile with a depth unit of m to obtain a subtracted stacking profile, and replace m in the depth unit of the subtracted stacking profile with ms;

[0040] The pre-stack depth migration result determination unit is used to perform pulse deconvolution processing on the subtracted stacked sections with a depth unit of ms to obtain the pre-stack depth migration result after direct suppression of multiple waves between onshore layers in the depth domain of the well control area.

[0041] The present invention has at least the following beneficial effects:

[0042] The present disclosure proposes a method and apparatus for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area. The method uses logging data to separate multiple wave pre-stack CRP gathers through wave equation forward modeling. The NMO velocity is found by velocity analysis of the multiple wave CRP gathers. The NMO velocity is used to flatten and remove the far-path stretching portion, and the multiple wave stacking profile is formed after stacking. The multiple wave profile is subtracted from the original PSDM profile before pulse deconvolution processing is performed. The present disclosure does not perform time-to-depth conversion throughout the process. Instead, the multiple waves of the pre-stack depth migration seismic stacking profile, which depicts the multiple waves more clearly, are directly suppressed. This avoids distortion caused by conversions between the depth domain and the time domain, thereby obtaining a higher-fidelity pre-stack depth migration result after multiple wave suppression. BRIEF DESCRIPTION OF THE DRAWINGS

[0043] The accompanying drawings herein are incorporated into and constitute a part of the specification. These drawings illustrate embodiments consistent with the present disclosure and, together with the specification, are used to explain the technical solutions of the present disclosure.

[0044] Figure 1 A flow chart illustrating a method for directly suppressing onshore interlayer multiple waves in the depth domain of a well-controlled area according to an embodiment of the present disclosure is shown. DETAILED DESCRIPTION

[0045] Various exemplary embodiments, features, and aspects of the present disclosure will be described in detail below with reference to the accompanying drawings. The same reference numerals in the accompanying drawings represent elements with the same or similar functions. Although various aspects of the embodiments are shown in the accompanying drawings, the drawings are not necessarily drawn to scale unless otherwise indicated.

[0046] The word “exemplary” is used exclusively herein to mean “serving as an example, example, or illustration.” Any embodiment described herein as “exemplary” is not necessarily to be construed as preferred or advantageous over other embodiments.

[0047] The term "and / or" herein simply describes an association relationship between associated objects, indicating that three relationships can exist. For example, "A and / or B" can represent the existence of three situations: A alone, A and B simultaneously, and B alone. Furthermore, the term "at least one" herein refers to any combination of at least two of any one or more of a plurality of items. For example, "at least one of A, B, and C" can represent any one or more elements selected from the set consisting of A, B, and C.

[0048] In addition, numerous specific details are provided in the following detailed description to better illustrate the present disclosure. Those skilled in the art will appreciate that the present disclosure can be practiced without certain specific details. In some instances, methods, means, components, and circuits well known to those skilled in the art are not described in detail in order to highlight the main points of the present disclosure.

[0049] Figure 1 A flow chart showing a method for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area according to an embodiment of the present disclosure is shown. Figure 1As shown, a method for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area comprises: step S01: obtaining well logging data and a stacked profile of pre-stack depth migration seismic results in a work area, wherein the depth unit of the well logging data and the unit of the stacked profile data of the pre-stack depth migration seismic results are both m; step S02: determining a typical well in the work area, and obtaining a pre-stack CRP gather with a depth unit of m after removing primary reflection wave information based on the well logging curve of the typical well by using a forward modeling formula considering propagation effect; step S03: calibrating the positions of the multiple waves on the pre-stack CRP gather after removing the primary reflection wave information, and replacing the depth unit m of the calibrated pre-stack CRP gather with m s; Step S04: for the calibrated multiple wave CRP gather with a depth unit of ms, NMO flattening and removal of the far-path stretching part are performed through velocity analysis to form a multiple wave stacking section; Step S05: replacing the depth unit of ms in the multiple wave stacking section with m, subtracting the pre-stack depth migration seismic result stacking section from the multiple wave stacking section with a depth unit of m to obtain the subtracted stacking section, and replacing the depth unit of m in the subtracted stacking section with ms; Step S06: performing pulse deconvolution processing on the subtracted stacking section with a depth unit of ms to obtain the pre-stack depth migration result after direct suppression of the multiple waves between onshore layers in the depth domain of the well control area.

[0050] The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area provided by the embodiment of the present invention specifically includes the following steps:

[0051] Step S01: Obtain the well logging data and pre-stack depth migration seismic results stacked profile Stk in the work area R , where the depth unit of the well logging data and the unit of the stacked profile data of the pre-stack depth migration PSDM seismic results are both in meters.

[0052] Step S02: Determine a typical well in the work area, and based on the well logging curves of the typical well, use a forward modeling formula that considers propagation effects to obtain a pre-stack CRP gather with a depth unit of m, removing the primary reflection wave information.

[0053] In the present disclosure, the method for determining a typical well in a work area includes: determining the well with the deepest vertical depth among all the wells in the work area as the typical well.

[0054] In the present disclosure, the method for obtaining a pre-stack CRP gather with a depth unit of m and the primary reflection wave information removed based on the well logging curve of the typical well by using a forward modeling formula considering propagation effects includes: substituting the P-wave and S-wave velocities of the well logging curve of the typical well into the forward modeling formula considering propagation effects to obtain a depth-domain wave equation forward modeling full-wavefield CRP gather with a unit of m; substituting the acoustic wave curve of the typical well into the forward modeling formula considering propagation effects, ignoring all propagation effects except the primary reflection, and forward modeling to obtain a depth-domain primary reflection wave CRP gather with a unit of m; and forward modeling the full-wavefield CRP gather using the depth-domain wave equation to subtract the depth-domain primary reflection wave CRP gather to obtain a pre-stack CRP gather with the primary reflection wave information removed.

[0055] In the embodiment of the present disclosure, the forward modeling formula considering the propagation effect includes:

[0056]

[0057] Where, (R d ) new The matrix wave equation depth domain seismic forward modeling results considering the effects of interface reflection and transmission, R d and R u is the reflection propagation matrix of the downgoing wave field and the upgoing wave field, T d and T u is the transmission propagation matrix of the downgoing wave field and the upgoing wave field, E d and E u The momentum equation of the elastic system representing the downgoing wave field and the upgoing wave field, r d represents the downlink wave reflection coefficient sequence.

[0058] Substituting the P- and S-wave velocities from typical well logging curves into the forward modeling formulas that account for propagation effects yields a full-wavefield CRP gather A in the depth domain, representing a unit of m. Substituting the primary reflection wave from a typical well acoustic wave curve into the forward modeling formulas that account for propagation effects yields a full-wavefield CRP gather B in the depth domain, representing a unit of m. Subtracting the primary reflection wave CRP gather B from the full-wavefield CRP gather A yields a prestack CRP gather C in the depth domain, representing a unit of m. This gather is primarily composed of multiple reflections and includes other wavefield information, such as converted waves. This gather is essentially a prestack CRP gather C with the primary reflection wave information removed.

[0059] Step S03: calibrating the positions of the multiple waves on the pre-stack CRP gather from which the primary reflection wave information is removed, and replacing the depth unit m of the calibrated multiple wave pre-stack CRP gather with ms.

[0060] In the present disclosure, the method for calibrating the positions of multiple waves on a pre-stack CRP gather from which primary reflection wave information has been removed includes: if a phase axis on the pre-stack CRP gather from which primary reflection wave information has been removed is bent, marking the depth position of the phase axis corresponding to an offset of 0 as the position where the multiple waves exist, thereby obtaining a calibrated multiple wave pre-stack CRP gather.

[0061] In the disclosed embodiment, the obtained pre-stack CRP gather C with the primary reflection wave information removed is analyzed and calibrated. When obvious coaxial bending occurs on the gather, the position with an offset of 0 depth is marked as the location where multiple waves exist. After all markings are completed, a calibrated multiple wave pre-stack CRP gather is obtained.

[0062] The depth unit m (meter) of the calibrated multiple wave pre-stack CRP gather is directly replaced by ms (milliseconds). No conversion is required during the replacement.

[0063] Step S04: For the calibrated multiple wave CRP gather with a depth unit of ms, NMO flattening and removal of the far-track stretched portion are performed through velocity analysis to form a multiple wave superposition profile.

[0064] In the present disclosure, the method for forming a multiple wave stacking profile by performing NMO flattening and removing the far-track stretched portion on a calibrated multiple wave CRP gather with a depth unit of ms through velocity analysis includes:

[0065] According to the calibrated multiple wave CRP gather, a weighted stacked seismic trace is determined; a velocity analysis is performed on the weighted stacked seismic trace to find the NMO velocity; the calibrated multiple wave pre-stack CRP gather is flattened using the NMO velocity to obtain a flattened CRP gather; the flattened CRP gather is subjected to a far-path stretched portion removal process to obtain a removed CRP gather; the removed CRP gather is horizontally stacked to form a multiple wave stacking profile.

[0066] In the embodiment of the present disclosure, before determining the weighted stacked seismic traces based on the calibrated multiple wave CRP gathers, the method further includes: performing random noise suppression processing on the calibrated multiple wave CRP gathers based on prediction error analysis to obtain a noise-suppressed CRP gather x j (i),.

[0067] In the embodiment of the present disclosure, the method for determining weighted stacked seismic traces based on the calibrated multiple wave CRP gathers includes: j(i) Determine a standard trace; determine a weighting coefficient based on the standard trace and the noise-suppressed CRP trace gather; determine a weighted seismic trace based on the weighting coefficient and the noise-suppressed CRP trace gather; and stack the weighted seismic traces to obtain a weighted stacked seismic trace.

[0068] CRP gathers after noise suppression x j (i) (i=1,2,…,M,j=1,2,…,N), where M is the number of samples and N is the number of channels in the channel set. j In (i), a standard track y(i) (i=1, 2, ..., M) is determined so that the difference between the standard track and each recording track is minimized.

[0069] The method for determining the standard path includes: obtaining the standard path by using the least squares principle calculation formula; wherein the least squares principle calculation formula for the standard path includes:

[0070]

[0071] In the formula, i = 1, 2, ..., M, and the standard trace is the average of N traces (N is the effective seismic trace) after superposition.

[0072] After calculating the standard trace, the jth seismic trace x j After multiplying the weighting coefficient, (i) should be closest to the standard track y(i) in the least square sense.

[0073] Among them, determine the weighting coefficient w j Method (i) includes: using formula (3), calculating the weighting coefficient w of point i0 j (i) is:

[0074]

[0075] Where T is the maximum recording track length.

[0076] According to the obtained weighting coefficient w j (i), using formula (4) to calculate x j (i) is weighted to obtain the weighted seismic trace r j (i);

[0077] r j (i) = x j (i)w j (i) (4)

[0078] The weighted seismic trace r j (i), use formula (5) to perform superposition and obtain the weighted stacked seismic trace s(i);

[0079]

[0080] In the disclosed embodiment, velocity analysis is performed on the obtained weighted stack seismic trace s(i) to find the NMO velocity. The calibrated multiple pre-stack CRP gathers are flattened using the NMO velocity to obtain the flattened CRP gathers.

[0081] In the present disclosure, the method of removing the far-path stretched portion of the flattened CRP gather to obtain a removed CRP gather comprises: removing a first predetermined percentage and a second predetermined percentage of the offset of the flattened CRP gather to obtain a multiple wave superposition profile.

[0082] In the embodiment of the present disclosure, the first predetermined percentage is 20% and the second predetermined percentage is 30%. For example, if the total offset of the leveled CRP gather is 5000m, the seismic data with offsets from 0m to 1000m and from 3500m to 1500m on the gather are removed.

[0083] In the embodiment of the present disclosure, the obtained CRP gathers after excision are horizontally stacked to obtain the multiple wave stacking profile Stk M .

[0084] Step S05: Replace the depth unit ms of the multiple wave stacking section with m, subtract the pre-stack depth migration seismic result stacking section from the multiple wave stacking section with a depth unit of m to obtain the subtracted stacking section, and replace the depth unit m of the subtracted stacking section with ms.

[0085] In the present disclosure, before replacing the depth unit m of the superimposed section after subtraction with ms, it also includes: performing quality control on the superimposed section after subtraction, and the method includes: performing cross-correlation calculation on the superimposed section after subtraction and the multiple wave superimposed section, and judging whether the cross-correlation result is less than a third predetermined percentage; if so, the quality control passes; otherwise, the quality control fails.

[0086] In the embodiment of the present disclosure, the pre-stack depth migration seismic result obtained in step S01 is stacked into a section Stk R , and the obtained multiple wave superposition profile Stk M Subtract and get the superimposed section Stk after subtraction (R-M) .

[0087] In the embodiment of the present disclosure, the third predetermined percentage is: 20%. (R-M) , and the multiple wave superposition profile Stk M Perform cross-correlation calculation to determine whether the cross-correlation result is less than 20%. If so, the quality control passes.

[0088] The superimposed profile Stk after subtraction of the quality control (R-M) The depth unit m is directly replaced by ms to get Stk (R-M)’ (t).

[0089] Step S06: performing pulse deconvolution processing on the stacked sections after subtraction with the depth unit being ms, and obtaining a pre-stack depth migration result after direct suppression of multiple waves between onshore layers in the depth domain of the well-controlled area.

[0090] In the embodiment of the present disclosure, for the seismic wavelet ω(t) = {ω(0), ω(1), L, ω(n)}, the expected output d(t) is a sharp pulse: d(t) = δ(t) = (1, 0, 0, L). The involved inverse filtering factor a(t) should theoretically be infinite, but in practice it is finite, and its main part is:

[0091] a(t)={a(-m0),a(-m0+1),L,a(-m0+m)} (6)

[0092] Where -m0 is the unknown starting time of the main part of the filter factor a(t), m is the total length of the filter factor a(t), and both m and m0 are positive integers.

[0093] When the seismic wavelet is at the minimum delay, m0=0, m=m, and the inverse filtering factor is a(t)={a(0), a(1), L, a(m)}. When t<0, ω(t)=0, and the matrix equation is:

[0094]

[0095] In formula (7), r ωω is the inverse filter factor corresponding to each sample point in the matrix. Solving this matrix equation can obtain the inverse filter operator a(t) with the expected output being a sharp pulse δ(t).

[0096] Use a(t) to input earthquake record Stk (R-M)’ (t) Deconvolution is performed to obtain the output seismic data after pulse deconvolution processing, that is, the pre-stack depth migration result Stk after direct suppression of multiple waves between onshore layers in the depth domain of the well control area (R-M)’DCN :

[0097]

[0098] In the present disclosure, the length m of the inverse filter operator of the pulse deconvolution processing is: if the total trace length of the seismic records in the work area is less than or equal to 10,000 ms, the length of the inverse filter operator is 400 ms; if the total trace length of the seismic records in the work area is greater than 10,000 ms, the length of the inverse filter operator is 600 ms.

[0099] The present disclosure further includes: performing quality control on the obtained pre-stack depth migration result, wherein the method includes: replacing ms of the depth unit of the pre-stack depth migration result with m;

[0100] A cross-correlation calculation is performed on the pre-stack depth migration result with a depth unit of m and the multiple wave stacking section to determine whether the cross-correlation result is less than a fourth predetermined percentage. If so, the quality control passes; otherwise, a new round of iterative pulse deconvolution processing is performed on the stacking section after the subtraction until the cross-correlation result between the pre-stack depth migration result and the multiple wave stacking section is less than the fourth predetermined percentage.

[0101] In the embodiment of the present disclosure, the pre-stack depth migration result Stk after the direct suppression of the multiple waves between the onshore layers in the depth domain of the well control area is obtained. (R-M)’DCN The depth unit is changed from ms to m, for Stk (R-M)’DCN Conduct quality control.

[0102] The fourth predetermined percentage is 10%. (R-M)’DCN Multiple wave stacking section Stk M Perform cross-correlation calculations, when Stk (R-M)’DCN and Stk M If the cross-correlation result is less than 10%, the quality control passes and the final result is output; if the cross-correlation result is greater than or equal to 10%, the pulse deconvolution process of step S06 is iteratively repeated until Stk (R-M)’DCN and Stk M When the cross-correlation result is less than 10%, the result Stk is output. (R-M)’DCN′ That is, it is the result of pre-stack depth migration after direct suppression of multiple waves between onshore layers in the depth domain of the well-controlled area.

[0103] Among them, due to Stk M From the forward modeling results of the well logging curve, so usually Stk M The record length is less than Stk (R-M)’DCN The length of the channel, here the cross-correlation calculation channel length used for quality control is Stk M The final output is Stk (R-M)’DCN All the Taoist priests.

[0104] It can be understood that the above-mentioned various method embodiments mentioned in the present disclosure can be combined with each other to form combined embodiments without violating the principle logic. Due to space limitations, the present disclosure will not elaborate on them.

[0105] The execution entity of the method for directly suppressing multiple waves between onshore layers in the depth domain of the well-controlled area may be a processing device for the method. For example, the method for directly suppressing multiple waves between onshore layers in the depth domain of the well-controlled area may be executed by a terminal device, a server, or other processing device. The terminal device may be a user equipment (UE), a mobile device, a user terminal, a terminal, a cellular phone, a cordless phone, a personal digital assistant (PDA), a handheld device, a computing device, an in-vehicle device, a wearable device, etc. In some possible implementations, the method for directly suppressing multiple waves between onshore layers in the depth domain of the well-controlled area may be implemented by a processor invoking computer-readable instructions stored in a memory.

[0106] Those skilled in the art will understand that in the above-mentioned method of the specific implementation method, the writing order of each step does not mean a strict execution order and does not constitute any limitation on the implementation process. The specific execution order of each step should be determined by its function and possible internal logic.

[0107] The present disclosure also provides a device for directly suppressing multiple waves between onshore layers in a well-controlled area in the depth domain, comprising: an acquisition unit, configured to acquire well logging data and a stacked profile of pre-stack depth migration seismic results in a work area, wherein the depth unit of the well logging data and the unit of the stacked profile data of the pre-stack depth migration seismic results are both in meters; a multiple wave pre-stack CRP gather determination unit, configured to determine a typical well in the work area, and obtain a pre-stack CRP gather with a depth domain unit of meters after removing primary reflection wave information based on the well logging curves of the typical well using a forward modeling formula considering propagation effects; a multiple wave calibration unit, configured to perform multiple wave calibration on the positions of multiple waves on the pre-stack CRP gather after removing primary reflection wave information, and replace the depth unit m of the calibrated multiple wave pre-stack CRP gather with ms; The multiple wave stacking profile determination unit is used to perform NMO flattening and remove the far-path stretching part of the calibrated multiple wave CRP data set with a depth unit of ms through velocity analysis to form a multiple wave stacking profile; the primary reflection wave stacking profile determination unit is used to replace the ms of the depth unit of the multiple wave stacking profile with m, subtract the pre-stack depth migration seismic result stacking profile from the multiple wave stacking profile with a depth unit of m to obtain the subtracted stacking profile, and replace the m of the depth unit of the subtracted stacking profile with ms; the pre-stack depth migration result determination unit is used to perform pulse deconvolution processing on the subtracted stacking profile with a depth unit of ms to obtain the pre-stack depth migration result after direct suppression of the multiple waves between onshore layers in the depth domain of the well control area.

[0108] In some embodiments, the functions or modules and units included in the device provided by the embodiments of the present disclosure can be used to execute the method described in the above method embodiments. The specific implementation can refer to the description of the above method embodiments. For the sake of brevity, it will not be repeated here.

[0109] The core of this disclosure lies in using well logging data to separate multiple channel gathers through wave equation forward modeling. After analyzing the multiple velocity, NMO flattening is performed, and the mid-50% offset information is superimposed to form a multiple profile. This multiple profile is then subtracted from the original PSDM profile before pulse deconvolution. The depth units, m and ms, are swapped twice to quickly complete the necessary processing without affecting the frequency. Finally, the depth-domain migration results after multiple suppression are output through quality control, enabling direct and rapid suppression of multiples between onshore layers without time-to-depth conversion.

[0110] There are two explanations on the technical applicability of the present disclosure. First, changing the unit of the depth domain processing result from m to ms does not affect the velocity analysis of the separated multiple wave channel set, and the velocity analysis result is independent of the frequency. Second, if the unit of the depth deviation result is directly changed to ms for Radon domain multiple wave suppression processing, many problems will arise. This type of method has a certain relationship with frequency, and the effective frequency band range of the result cannot be accurately estimated in the depth domain. The basis of high-precision Radon transform is to perform frequency one by one in the frequency domain, which will cause serious frequency distortion and inverse transformation distortion problems. The method disclosed in the present disclosure does not use frequency-related processing methods after the depth unit is replaced, thereby avoiding the occurrence of the above problems.

[0111] The present disclosure realizes the horizontal well platform area on the basis of not performing time-depth conversion in the whole process, which is usually defined as 10km 2 Directly suppressing multiple waves in the target area of ​​a horizontal well with a relatively flat structure in the depth domain can avoid the distortion of depth domain processing results when converting back and forth between the depth domain and the time domain. Good results have been seen in actual data processing in multiple working areas, and it has practicality and broad prospects for promotion and application.

[0112] While various embodiments of the present disclosure have been described above, the foregoing description is intended to be illustrative, non-exhaustive, and not limited to the disclosed embodiments. Many modifications and variations will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments. The terminology used herein is selected to best explain the principles of the embodiments, their practical applications, or technological improvements in the marketplace, or to enable others skilled in the art to understand the embodiments disclosed herein.

Claims

1. A method for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area, characterized in that: include: Obtain the well logging data and pre-stack depth migration seismic stacked profiles of the work area. The depth units of the well logging data and pre-stack depth migration seismic stacked profiles are both in meters. Determine typical wells in the work area, and based on the well logging curves of the typical wells, use the forward modeling formula considering the propagation effect to obtain pre-stack CRP gathers in the depth domain with the unit of m, removing the primary reflection wave information; Calibrate the positions of the multiple waves on the prestack CRP gather from which the primary reflection wave information is removed, and replace the depth unit m of the calibrated multiple wave prestack CRP gather with ms; For the calibrated multiple wave CRP gathers with depth unit of ms, NMO flattening and removal of the far-track stretching part are performed through velocity analysis to form a multiple wave superposition profile; Replacing the depth unit of the multiple wave stack section with ms by m, subtracting the pre-stack depth migration seismic result stack section from the multiple wave stack section with a depth unit of m to obtain a subtracted stack section, and replacing the depth unit of m in the subtracted stack section with ms; Pulse deconvolution processing is performed on the subtracted stacked sections with a depth unit of ms to obtain pre-stack depth migration results after direct suppression of onshore interlayer multiple waves in the depth domain of the well-controlled area.

2. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 1, characterized in that: The method for determining typical wells in a work area comprises: Determine that the well with the deepest vertical depth among all the wells in the work area is the typical well.

3. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 1, characterized in that: The method of obtaining a pre-stack CRP gather with a depth unit of m and primary reflection wave information removed based on the well logging curve of the typical well by using a forward modeling formula considering propagation effect includes: Substituting the P-wave and S-wave velocities of the logging curve of the typical well into the forward modeling formula considering the propagation effect, a depth domain wave equation forward modeling full wavefield CRP gather in units of m is obtained; Substituting the typical well acoustic wave curve into the forward modeling formula considering propagation effect, ignoring all propagation effects except the primary reflection, the depth domain primary reflection wave CRP gather in units of m is obtained by forward modeling; The full-wavefield CRP gathers are forward modeled using the depth-domain wave equation, and the depth-domain primary reflection wave CRP gathers are subtracted to obtain the pre-stack CRP gathers without the primary reflection wave information.

4. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 1, characterized in that: The method for calibrating the positions of multiple waves on the pre-stack CRP gather after removing the primary reflection wave information comprises: If the event on the pre-stack CRP gather with the primary reflection wave information removed is bent, the depth position of the event corresponding to the offset distance of 0 is marked as the multiple wave existence position, thereby obtaining a calibrated multiple wave pre-stack CRP gather.

5. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 1, characterized in that: The method for forming a multiple wave superposition profile by performing NMO flattening and removing a far-track stretched portion on a calibrated multiple wave CRP gather having a depth unit of ms through velocity analysis includes: determining weighted stacked seismic traces based on the calibrated multiple wave CRP gathers; Performing velocity analysis on the weighted stacked seismic traces to find the NMO velocity; flattening the calibrated multiple pre-stack CRP gathers using the NMO velocity to obtain a flattened CRP gather; The far-path stretching part is removed from the flattened CRP gather to obtain the removed CRP gather; The CRP gathers after the excision are horizontally stacked to form a multiple wave stacking section.

6. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 5, characterized in that: The method of performing a process of removing the far-track stretching portion on the flattened CRP gather to obtain the removed CRP gather comprises: For the flattened CRP gather, the first predetermined percentage and the second predetermined percentage of the offset are cut off to obtain the multiple wave superposition profile.

7. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to claim 1, characterized in that: Before replacing the depth unit m of the subtracted stacked section with ms, the method further includes: performing quality control on the subtracted stacked section, the method comprising: A cross-correlation calculation is performed on the superposition profile after the subtraction and the multiple wave superposition profile to determine whether the cross-correlation result is less than a third predetermined percentage. If so, the quality control passes; otherwise, the quality control fails.

8. The method for directly suppressing onshore interlayer multiple waves in the depth domain of the well control area according to claim 1, characterized in that: The length m of the inverse filter operator of the pulse deconvolution process is: if the total trace length of the seismic records in the work area is less than or equal to 10,000 ms, the inverse filter operator length is 400 ms; if the total trace length of the seismic records in the work area is greater than 10,000 ms, the inverse filter operator length is 600 ms.

9. The method for directly suppressing multiple waves between onshore layers in the depth domain of the well control area according to any one of claims 1 to 8, characterized in that: The method further includes: performing quality control on the obtained prestack depth migration results, wherein the method includes: The depth unit ms in the prestack depth migration result is replaced by m; A cross-correlation calculation is performed on the pre-stack depth migration result with a depth unit of m and the multiple wave stacking section to determine whether the cross-correlation result is less than a fourth predetermined percentage. If so, the quality control passes; otherwise, a new round of iterative pulse deconvolution processing is performed on the stacking section after the subtraction until the cross-correlation result between the pre-stack depth migration result and the multiple wave stacking section is less than the fourth predetermined percentage.

10. A device for directly suppressing multiple waves between onshore layers in the depth domain of a well-controlled area, characterized in that: include: The acquisition unit is used to obtain the well logging data and the stacked section of the pre-stack depth migration seismic results in the work area. The depth unit of the well logging data and the unit of the stacked section data of the pre-stack depth migration seismic results are both in meters. The multiple wave pre-stack CRP gather determination unit is used to determine the typical wells in the work area, and obtain the pre-stack CRP gathers with the primary reflection wave information removed in the depth domain unit of m based on the logging curves of the typical wells using the forward modeling formula considering the propagation effect; a multiple wave calibration unit, configured to calibrate the positions of the multiple waves on the prestack CRP gather from which the primary reflection wave information has been removed, and replace the depth unit m of the calibrated multiple wave prestack CRP gather with ms; The multiple wave superposition profile determination unit is used to perform NMO flattening and remove the far-track stretching part of the calibrated multiple wave CRP gather with a depth unit of ms through velocity analysis to form a multiple wave superposition profile; a primary reflection wave stacking profile determining unit, configured to replace ms in the depth unit of the multiple wave stacking profile with m, subtract the pre-stack depth migration seismic result stacking profile from the multiple wave stacking profile with a depth unit of m to obtain a subtracted stacking profile, and replace m in the depth unit of the subtracted stacking profile with ms; The pre-stack depth migration result determination unit is used to perform pulse deconvolution processing on the subtracted stacked sections with a depth unit of ms to obtain the pre-stack depth migration result after direct suppression of multiple waves between onshore layers in the depth domain of the well control area.

Images