Residual multiple suppression method and system

Abstract

The invention relates to the technical field of oil exploration, and provides a residual multiple suppression method and system. The method comprises the following steps: acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transformation processing on the original seismic data to obtain first seismic data; performing bunching filtering processing on the first seismic data to obtain second seismic data; and performing slowness filtering processing on the second seismic data to obtain target seismic data. According to the embodiment of the invention, the residual multiples in the seismic data can be effectively suppressed on the premise of not damaging effective signals.

Description

Technical Field

[0001] The present disclosure relates to the technical field of petroleum exploration, and in particular to a residual multiple wave suppression method and system. Background Art

[0002] When the reflected wave encounters an interface with a large difference in wave impedance during the upward propagation process, the reflected wave may propagate downward from this interface. When the reflected wave encounters the reflection interface, it is reflected again, resulting in the formation of multiple waves. As the focus and difficulty of seismic data processing, multiple waves have the characteristics of wide distribution range, strong energy, and many types, becoming the main bottleneck affecting the imaging of medium and deep geological structures. In terms of reflection mechanism, multiple waves include surface multiple waves, interlayer multiple waves, and micro-bend multiple waves; in terms of time difference, they are divided into long-period multiple waves and short-period multiple waves, among which the offshore seismic data also has the phenomenon of ringing. In conventional seismic data processing, multiple waves are noise signals that seriously interfere with the primary wave, and are often superimposed with the primary wave and difficult to separate, thereby reducing the signal-to-noise ratio of the seismic data, affecting the imaging effect of the seismic data, and hindering the correct interpretation of the seismic data. Therefore, studying the suppression method of multiple waves has considerable practical significance and application value.

[0003] According to the difference of residual time difference, the periodicity of multiple waves, and the predictability of multiple waves, the existing multiple wave suppression technology can be divided into two categories: wave theory prediction subtraction method and filtering method. Both methods have their own advantages and disadvantages. The wave theory prediction subtraction method is driven by velocity model or seismic data, predicts the multiple wave model from the seismic data, and then adaptively subtracts it from the seismic data. The prediction subtraction method can be divided into two types of methods: model-driven and data-driven. The advantage of this method is that it has high suppression accuracy and strong adaptability to complex structures, but due to its disadvantages such as dependence on velocity accuracy and large amount of calculation, its application scope is limited. The filtering method is based on the premise that multiple waves and primary waves have obvious characteristic differences in a specific transform domain. Various mathematical algorithms are used to calculate this difference to eliminate multiple waves, such as differences in period or time and frequency. Appropriate parameters are selected to transform data into different domains to eliminate multiple waves, such as predictive deconvolution, Radon transform, FK transform, etc. In addition, since the difference between multiple waves and primary waves at near offset is very small, after filtering and separation, there are still residual multiple waves concentrated at near offset, which has an adverse effect on seismic imaging. In order to reduce the adverse effects of residual multiple waves on seismic data stacking, the incision method is usually adopted to cut off some near-offset data, but this will affect the imaging effect of seismic data and have an adverse impact on interpretation. Most filtering methods in the prior art are based on simple geological assumptions. When the actual geological conditions are more complex, they often fail to achieve the desired effect and even damage harmful signals. Summary of the invention

[0004] In view of this, the embodiments of the present disclosure provide a residual multiple wave suppression method and system to solve the problem that most filtering methods in the prior art are based on simple geological assumptions and often fail to achieve ideal results or even damage harmful signals when the actual geological conditions are more complex.

[0005] A first aspect of an embodiment of the present disclosure provides a residual multiple wave suppression method, comprising:

[0006] Acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0007] Performing a beamforming filter process on the first seismic data to obtain second seismic data;

[0008] The second seismic data is subjected to slowness filtering to obtain target seismic data.

[0009] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0010] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0011] Generate a corresponding prediction identity matrix according to the initial random noise;

[0012] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0013] The second seismic data is calculated based on the original seismic data and the filter.

[0014] In some embodiments, the first calculation formula includes:

[0015] x=Bp+Cv+u

[0016] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0017] In some embodiments, the second calculation formula includes:

[0018] H p =G(A HQ -1 A) -1 A H Q -1

[0019] Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix;

[0020] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0021] In some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0022] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0023] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0024] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0025] In some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0026] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0027] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0028] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0029] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0030] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0031] In some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0032] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0033] Q(x,m)=P(x,t=m+T 1 (x))

[0034] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0035] In some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0036] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0037]

[0038] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0039] In some embodiments, the step of smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0040] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0041] U(x,m)=V(x,m)R(x,m)

[0042]

[0043] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0044] A second aspect of the embodiments of the present disclosure provides a residual multiple wave suppression system, comprising:

[0045] An acquisition module is used to acquire original seismic data, and perform pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0046] A first processing module, configured to perform a beamforming filter process on the first seismic data to obtain second seismic data;

[0047] The second processing module is used to perform slowness filtering on the second seismic data to obtain target seismic data.

[0048] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0049] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0050] Generate a corresponding prediction identity matrix according to the initial random noise;

[0051] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0052] The second seismic data is calculated based on the original seismic data and the filter.

[0053] In some embodiments, the first calculation formula includes:

[0054] x=Bp+Cv+u

[0055] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0056] In some embodiments, the second calculation formula includes:

[0057] H p =G(A H Q -1 A) -1 A H Q -1

[0058] Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu]H , E[] represents the mathematical expectation value of the matrix;

[0059] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0060] In some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0061] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0062] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0063] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0064] In some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0065] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0066] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0067] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0068] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0069] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0070] In some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0071] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0072] Q(x,m)=P(x,t=m+T 1 (x))

[0073] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0074] In some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0075] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0076]

[0077] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0078] In some embodiments, the step of smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0079] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0080] U(x,m)=V(x,m)R(x,m)

[0081]

[0082] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0083] According to a third aspect of the embodiments of the present disclosure, an electronic device is provided, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the steps implemented when the processor executes the computer program include:

[0084] Acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0085] Performing a beamforming filter process on the first seismic data to obtain second seismic data;

[0086] The second seismic data is subjected to slowness filtering to obtain target seismic data.

[0087] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0088] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0089] Generate a corresponding prediction identity matrix according to the initial random noise;

[0090] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0091] The second seismic data is calculated based on the original seismic data and the filter.

[0092] In some embodiments, the first calculation formula includes:

[0093] x=Bp+Cv+u

[0094] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0095] In some embodiments, the second calculation formula includes:

[0096] H p =G(A H Q -1 A) -1 A H Q -1

[0097] Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix;

[0098] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0099] In some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0100] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0101] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0102] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0103] In some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0104] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0105] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0106] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0107] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0108] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0109] In some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0110] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0111] Q(x,m)=P(x,t=m+T 1 (x))

[0112] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0113] In some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0114] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0115]

[0116] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0117] In some embodiments, the step of smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0118] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0119] U(x,m)=V(x,m)R(x,m)

[0120]

[0121] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0122] In a fourth aspect of the embodiments of the present disclosure, a computer-readable storage medium is provided, wherein the computer-readable storage medium stores a computer program, and the steps implemented when the computer program is executed by a processor include:

[0123] Acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0124] Performing a beamforming filter process on the first seismic data to obtain second seismic data;

[0125] The second seismic data is subjected to slowness filtering to obtain target seismic data.

[0126] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0127] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0128] Generate a corresponding prediction identity matrix according to the initial random noise;

[0129] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0130] The second seismic data is calculated based on the original seismic data and the filter.

[0131] In some embodiments, the first calculation formula includes:

[0132] x=Bp+Cv+u

[0133] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0134] In some embodiments, the second calculation formula includes:

[0135] H p =G(A H Q -1 A) -1 A H Q -1

[0136] Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix;

[0137] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0138] In some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0139] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0140] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0141] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0142] In some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0143] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0144] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0145] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0146] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0147] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0148] In some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0149] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0150] Q(x,m)=P(x,t=m+T 1 (x))

[0151] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0152] In some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0153] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0154]

[0155] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0156] In some embodiments, the step of smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0157] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0158] U(x,m)=V(x,m)R(x,m)

[0159]

[0160] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0161] In a fifth aspect of the embodiments of the present disclosure, a computer program product is provided. The computer program product includes a computer program or an instruction. When the computer program or the instruction is executed by a processor, the steps implemented include:

[0162] Acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0163] Performing a beamforming filter process on the first seismic data to obtain second seismic data;

[0164] The second seismic data is subjected to slowness filtering to obtain target seismic data.

[0165] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0166] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0167] Generate a corresponding prediction identity matrix according to the initial random noise;

[0168] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0169] The second seismic data is calculated based on the original seismic data and the filter.

[0170] In some embodiments, the first calculation formula includes:

[0171] x=Bp+Cv+u

[0172] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0173] In some embodiments, the second calculation formula includes:

[0174] H p =G(A H Q -1 A) -1 A H Q -1

[0175] Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix;

[0176] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0177] In some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0178] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0179] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0180] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0181] In some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0182] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0183] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0184] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0185] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0186] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0187] In some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0188] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0189] Q(x,m)=P(x,t=m+T 1 (x)

[0190] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0191] In some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0192] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0193]

[0194] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0195] In some embodiments, the step of smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0196] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0197] U(x,m)=V(x,m)R(x,m)

[0198]

[0199] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0200] The remaining multiple waves are suppressed by using the beamforming filter and slowness filter technology in series combination, and the residual multiple waves in the post-stack data are suppressed by layer flattening and instantaneous slowness filtering after stacking, so as to achieve a better balance between retaining effective waves and suppressing multiple waves, so as to achieve the purpose of suppressing multiple waves and protecting effective signals, and also have the advantages of better amplitude preservation, better retention of low frequencies, and more conducive to observing geological structures such as volcanic structures contained in seismic data. It can effectively suppress the remaining multiple waves in seismic data without damaging effective signals. BRIEF DESCRIPTION OF THE DRAWINGS

[0201] In order to more clearly illustrate the technical solutions in the embodiments of the present disclosure, 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 for the embodiments, and ordinary technicians in this field can obtain other drawings based on these drawings without paying any creative work.

[0202] Figure 1 is a schematic diagram of an application scenario of the residual multiple wave suppression method provided according to an embodiment of the present disclosure;

[0203] Figure 2 is a flow chart of some embodiments of a residual multiple wave suppression method provided according to an embodiment of the present disclosure;

[0204] Figure 3 is a flow chart of some other embodiments of another residual multiple wave suppression method provided according to an embodiment of the present disclosure;

[0205] Figure 4 It is a simplified structural diagram of a residual multiple wave suppression system provided according to an embodiment of the present disclosure;

[0206] Figure 5 is a schematic diagram of an electronic device provided according to an embodiment of the present disclosure;

[0207] Figure 6 is a schematic diagram comparing the effects before and after beamforming filtering according to an embodiment of the present disclosure;

[0208] Figure 7 is a schematic diagram comparing the effects before and after slowness filtering according to an embodiment of the present disclosure;

[0209] Figure 8 It is a schematic diagram of the effect of a residual multiple wave suppression method provided in accordance with an embodiment of the present disclosure;

[0210] Fig. 9 is a schematic diagram of a preliminary stacked cross section of region A provided according to an embodiment of the present disclosure;

[0211] Fig.10 It is a comparative diagram of the processing effects of pre-stack deconvolution and high-precision Radon transform for region A provided in an embodiment of the present disclosure;

[0212] Fig.11 is a schematic diagram of the effect of performing a residual multiple wave suppression method on area A according to an embodiment of the present disclosure;

[0213] Fig.12 It is a comparative diagram of the traditional method for removing residual multiple waves and the residual multiple wave suppression method provided by the embodiment of the present disclosure. DETAILED DESCRIPTION

[0214] Embodiments of the present disclosure will be described in more detail below with reference to the accompanying drawings. Although certain embodiments of the present disclosure are shown in the accompanying drawings, it should be understood that the present disclosure can be implemented in various forms and should not be construed as being limited to the embodiments set forth herein. On the contrary, these embodiments are provided to provide a more thorough and complete understanding of the present disclosure. It should be understood that the drawings and embodiments of the present disclosure are only for exemplary purposes and are not intended to limit the scope of protection of the present disclosure.

[0215] It should also be noted that, for ease of description, only the parts related to the present disclosure are shown in the drawings. In the absence of conflict, the embodiments and features in the embodiments of the present disclosure can be combined with each other.

[0216] It should be noted that the concepts such as "first" and "second" mentioned in the present disclosure are only used to distinguish different systems, devices, modules or units, and are not used to limit the order or interdependence of the functions performed by these systems, devices, modules or units.

[0217] It should be noted that the modifications of "one" and "plurality" mentioned in the present disclosure are illustrative rather than restrictive, and those skilled in the art should understand that unless otherwise clearly indicated in the context, it should be understood as "one or more".

[0218] The names of the messages or information exchanged between multiple devices in the embodiments of the present disclosure are only used for illustrative purposes and are not used to limit the scope of these messages or information.

[0219] The present disclosure will be described in detail below with reference to the accompanying drawings and in conjunction with embodiments.

[0220] Embodiment 1:

[0221] Figure 1 It is a schematic diagram of an application scenario of the residual multiple wave suppression method according to the first embodiment of the present disclosure.

[0222] exist Figure 1 In the application scenario, first, the computing device 101 can obtain original seismic data 102, perform pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data 102, and obtain first seismic data 103;

[0223] Secondly, the computing device 101 may perform a beamforming filter process on the first seismic data 103 to obtain second seismic data 104;

[0224] Finally, the computing device 101 may perform slowness filtering on the second seismic data 104 to obtain target seismic data 105 .

[0225] It should be noted that the computing device 101 can be hardware or software. When the computing device is hardware, it can be implemented as a distributed cluster consisting of multiple servers or terminal devices, or as a single server or a single terminal device. When the computing device is embodied as software, it can be installed in the hardware devices listed above. It can be implemented as multiple software or software modules for providing distributed services, or as a single software or software module. No specific limitation is made here.

[0226] It should be understood that Figure 1 The number of computing devices in the embodiment is only illustrative. Any number of computing devices may be provided according to implementation requirements.

[0227] Embodiment 2:

[0228] Continue to refer Figure 2 , shows a process 200 of the second embodiment of the residual multiple wave suppression method according to the present disclosure. The method can be Figure 1 The method for suppressing residual multiple waves comprises the following steps:

[0229] Step 201, acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data.

[0230] In some embodiments, the execution body of the residual multiple wave suppression method (such as Figure 1 The computing device 101 shown in the figure can be connected to the target device through a wired connection or a wireless connection, and then, the original seismic data is obtained, and the original seismic data is processed by pre-stack prediction deconvolution and high-precision Radon transform to obtain the first seismic data. The original seismic data may refer to seismic data that has not been processed by multiple wave suppression. After the original seismic data is processed by pre-stack prediction deconvolution and high-precision Radon transform (that is, conventional multiple wave suppression is performed), the first seismic data with residual multiple waves remaining after multiple wave suppression is obtained.

[0231] It should be noted that the above-mentioned wireless connection methods may include but are not limited to 3G / 4G / 5G connection, WiFi connection, Bluetooth connection, WiMAX connection, Zigbee connection, UWB (ultra wideband) connection, and other wireless connection methods currently known or developed in the future.

[0232] Step 202: Perform beamforming filtering on the first seismic data to obtain second seismic data.

[0233] In some embodiments, the execution subject may perform beamforming filtering on the first seismic data to obtain second seismic data. The second seismic data may refer to seismic data obtained after beamforming filtering. Beamforming filtering may refer to a multi-channel filtering method that can extract distortion-free signals under the premise of minimizing output noise energy. Figure 6 As shown, Figure 6 The left diagram in shows the seismic signal before beamforming filtering, and the right diagram shows the seismic signal after beamforming filtering.

[0234] In some optional implementations, the execution subject may perform beamforming filtering on the first seismic data to obtain second seismic data through the following steps:

[0235] The first step is to decompose the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information containing multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model. The initial effective signal may refer to the incompletely determined target seismic data (i.e., data without multiple waves) decomposed according to the first calculation formula. Coherent noise may refer to relevant noise data related to the formation that is inevitably generated during the operation. Random noise may refer to noise that occurs randomly during the operation, such as noise caused by passing trains, cars or other sudden situations. Initial coherent noise information may refer to coherent noise decomposed and predicted from the first seismic data. Initial random noise may refer to random noise decomposed and predicted from the first seismic data. The effective signal filter operator may refer to parameters related to the initial effective signal. The coherent noise filter operator may refer to parameters related to coherent noise (i.e., parameters related to multiple waves).

[0236] In some optional implementations, the first calculation formula may include:

[0237] x=Bp+Cv+u

[0238] Among them, x represents the first seismic data, B represents the effective signal filter operator, p represents the predicted effective signal, C represents the coherent noise filter operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format. It should be pointed out that the predicted effective signal may refer to the predicted data interval corresponding to the second seismic data preliminarily decomposed according to x, which is not the same as the second seismic data.

[0239] The second step is to generate a corresponding prediction identity matrix according to the initial random noise. The prediction identity matrix may refer to an identity matrix obtained by multiplying the matrix composed of the initial random noise obtained by decomposition with itself.

[0240] The third step is to calculate the target filter of the multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator and the prediction unit matrix, and the second calculation formula constructed based on the Lagrange multiplier method.

[0241] In some optional implementations, the second calculation formula includes:

[0242] H p =G(A H Q -1 A) -1 A H Q -1

[0243] Among them, H prepresents the target filter, G = (I, 0), I represents the unit matrix when the random noise is normally distributed and each channel has the same energy, A = (B, C), Q = E[uu] H , E[] represents the mathematical expectation value of the matrix;

[0244] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0245] The fourth step is to calculate the second seismic data according to the original seismic data and the filter. The second seismic data can be obtained by multiplying the original seismic data and the filter.

[0246] Step 203: Perform slowness filtering on the second seismic data to obtain target seismic data.

[0247] In some embodiments, the above-mentioned execution entity may perform slowness filtering on the second seismic data to obtain target seismic data. The target seismic data may refer to the target seismic data (signal) obtained after further eliminating the multiple waves after the beamforming filtering and slowness filtering. Slowness filtering may refer to the filtering method of flattening the interlayer multiple waves within a given range below the strong reflection layer according to the time of the picked up overlying strong reflection layer, defining the sync axis with instantaneous slowness close to 0 as the interlayer multiple waves, using Radon filtering to remove the flattened interlayer multiple waves, and then flattening the formation to obtain the section after suppressing the interlayer multiple waves. Figure 7 and Figure 8 As shown, Figure 7 The left diagram in the figure shows the seismic signal before slowness filtering, and the right diagram shows the seismic signal after slowness filtering. Figure 8 The left diagram in FIG. 1 shows a seismic signal with multiple wavefronts removed according to the method disclosed in the present invention. Figure 8 The middle diagram in FIG. 1 shows the seismic signal after multiple waves are removed according to the method disclosed in the present invention. Figure 8 The right side diagram in shows the multiple waves removed according to the method of the present disclosure.

[0248] In some optional implementations, the execution subject may perform slowness filtering on the second seismic data based on the following steps to obtain target seismic data:

[0249] In the first step, the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface are extracted from the second seismic data. The first layer bottom interface may refer to the layer bottom interface of the stratum corresponding to the generation of the interlayer multiple waves. The second layer bottom interface may refer to the layer bottom interface preset (e.g., determined by experiment) for setting the flattening interval.

[0250] The second step is to perform a slowness smoothing process on the third seismic data to obtain fourth seismic data. The third seismic data may refer to the interlayer data between the bottom interface of the first layer and the bottom interface of the second layer. The slowness smoothing process may refer to the data obtained after the interlayer data is processed by slowness processing, smoothing and other steps.

[0251] In some optional implementations, the third seismic data is subjected to slowness smoothing processing to obtain fourth seismic data, including: flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data; performing an event-strengthening transformation on the flattened seismic data to obtain event-strengthening seismic data; calculating the instantaneous slowness of each data point in the event-strengthening seismic data according to a difference calculation method, and smoothing each instantaneous slowness to obtain smoothed seismic data; performing an event-weakening transformation on the smoothed seismic data to obtain event-weakened seismic data; and performing an anti-layer flattening process on the event-weakened seismic data to obtain the fourth seismic data.

[0252] In some optional implementations, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0253] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0254] Q(x,m)=P(x,t=m+T 1 (x)

[0255] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0256] In some optional implementations, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0257] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0258]

[0259] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0260] In some optional implementations, smoothing each instantaneous slowness to obtain smoothed seismic data includes:

[0261] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0262] U(x,m)=V(x,m)R(x,m)

[0263]

[0264] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0265] It should be pointed out that the calculation formula used to transform the smoothed seismic data into weakened event and the calculation formula used to perform anti-layer flattening processing on the weakened event seismic data can be inversely calculated from the above calculation formulas, and will not be repeated here.

[0266] The third step is to replace the third seismic data in the second seismic data with the fourth seismic data to obtain the target seismic data.

[0267] The first to fifth calculation formulas mentioned above are all used to distinguish different calculation formulas in the embodiments of the present disclosure and can be set as needed without any limitation.

[0268] The beneficial effects of one of the implementation methods in the above-mentioned embodiments of the present disclosure include at least: suppressing the remaining multiple waves by using the beamforming filter and slowness filter technology in series combination, suppressing the residual multiple waves in the post-stack data by means of layer flattening and instantaneous slowness filtering after stacking, thereby achieving a better balance between retaining effective waves and suppressing multiple waves, achieving the purpose of suppressing multiple waves and protecting effective signals, and also having the advantages of better amplitude preservation, better retention of low frequencies, and being more conducive to observing geological structures such as volcanic structures contained in seismic data. It is possible to effectively suppress the remaining multiple waves in seismic data without damaging effective signals.

[0269] The present disclosure is further described below through a specific embodiment:

[0270] By studying the seismic profile of area A, the following relevant diagrams and simple explanations can be obtained for different processing steps:

[0271] like Fig. 9 , showing the initial stacking section of the area.

[0272] After pre-stack deconvolution and high-precision Radon transform of the seismic data in this area, the first processed related data are obtained, such as Fig.10 As shown, Fig.10 The left diagram in the figure shows the seismic data with multiple wavefronts removed. Fig.10 The middle diagram in the figure shows the seismic data after deconvolution and high-precision Radon transformation. Fig.10 The right figure in the figure shows the multiple waves removed by deconvolution and high-precision Radon transform.

[0273] After the seismic data processed for the first time is processed by the present disclosure, the comparison effect is as follows: Fig.11 As shown, Fig.11 The left diagram in FIG. 1 shows seismic data with multiple wavefronts removed according to the disclosed method. Fig.11 The middle diagram in FIG. 1 represents the seismic data processed according to the method disclosed in the present invention. Fig.11 The right side of FIG. 3 illustrates the multiple waves removed according to the method of the present disclosure.

[0274] The cross section after the residual multiple waves are removed by the traditional method and the residual multiple waves are removed by the present disclosure is as follows: Fig.12 As shown, the left side represents the seismic data after the residual multiple waves are removed by the traditional method, and the right side represents the seismic data after the residual multiple waves are removed by using the present disclosure.

[0275] Embodiment three:

[0276] Continue to refer Figure 3 , shows a process 300 of the third embodiment of the residual multiple wave suppression method according to the present disclosure, which can be Figure 1 The residual multiple wave suppression method includes:

[0277] Step 301, acquiring original seismic data, performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data, and obtaining first seismic data.

[0278] Step 302, decompose the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information containing multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model.

[0279] Step 303: Generate a corresponding prediction unit matrix according to the initial random noise.

[0280] Step 304 , calculate the target filter of the multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator and the prediction unit matrix, and the second calculation formula constructed based on the Lagrange multiplier method.

[0281] Step 305, calculating second seismic data based on the original seismic data and the filter.

[0282] Step 306: extract the first layer bottom interface generating the interlayer multiple waves and the third seismic data between the preset second layer bottom interfaces from the second seismic data.

[0283] Step 307, flatten the third seismic data along the bottom interface of the first layer to obtain flattened seismic data.

[0284] Step 308, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data.

[0285] Step 309: The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data.

[0286] Step 310, performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data.

[0287] Step 311, performing anti-layer flattening processing on the event-weakening seismic data to obtain fourth seismic data.

[0288] Step 312, replacing the third seismic data in the second seismic data with the fourth seismic data to obtain target seismic data.

[0289] In some embodiments, the specific implementation of steps 301-312 and the technical effects brought about can be referred to Figure 2 The corresponding steps in the second embodiment are not described again here.

[0290] All the above optional technical solutions can be combined in any way to form optional implementation methods of the present application, and will not be described in detail here.

[0291] Embodiment 4:

[0292] The following are system embodiments of the present disclosure, which can be used to execute the method embodiments of the present disclosure. For details not disclosed in the system embodiments of the present disclosure, please refer to the method embodiments of the present disclosure.

[0293] Further references Figure 4 As an implementation of the above methods in the above figures, the present disclosure provides an embodiment of a residual multiple wave suppression system. The system embodiment is similar to Figure 2 The above-mentioned embodiment 2 corresponds.

[0294] like Figure 4 As shown, the residual multiple wave suppression system 400 of this embodiment includes:

[0295] An acquisition module is used to acquire original seismic data, and perform pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data;

[0296] A first processing module, configured to perform a beamforming filter process on the first seismic data to obtain second seismic data;

[0297] The second processing module is used to perform slowness filtering on the second seismic data to obtain target seismic data.

[0298] In some embodiments, performing a beamforming filter process on the first seismic data to obtain second seismic data includes:

[0299] Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model;

[0300] Generate a corresponding prediction identity matrix according to the initial random noise;

[0301] Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method;

[0302] The second seismic data is calculated based on the original seismic data and the filter.

[0303] In some optional implementations of some embodiments, the first calculation formula includes:

[0304] x=Bp+Cv+u

[0305] Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

[0306] In some optional implementations of some embodiments, the second calculation formula includes:

[0307] H p =G(A H Q -1 A) -1 A H Q -1

[0308] Among them, H prepresents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix;

[0309] In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

[0310] In some optional implementations of some embodiments, performing slowness filtering on the second seismic data to obtain target seismic data includes:

[0311] Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data;

[0312] Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data;

[0313] The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

[0314] In some optional implementations of some embodiments, performing slowness smoothing on the third seismic data to obtain fourth seismic data includes:

[0315] Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data;

[0316] Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data;

[0317] The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data;

[0318] Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data;

[0319] The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

[0320] In some optional implementations of some embodiments, flattening the third seismic data along the first layer bottom interface to obtain the flattened seismic data includes:

[0321] The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes:

[0322] Q(x,m)=P(x,t=m+T 1 (x))

[0323] Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, T 1 Indicates the bottom interface of the first layer that generates interlayer multiple waves.

[0324] In some optional implementations of some embodiments, performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes:

[0325] According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes:

[0326]

[0327] Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

[0328] In some optional implementations of some embodiments, the smoothing of each instantaneous slowness to obtain smoothed seismic data includes:

[0329] According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes:

[0330] U(x,m)=V(x,m)R(x,m)

[0331]

[0332] Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

[0333] It is understood that the modules described in the system 400 are similar to those described in the reference Figure 2 Therefore, the operations, features and beneficial effects described above for the method are also applicable to the system 400 and the modules included therein, and will not be described in detail here.

[0334] Embodiment five:

[0335] like Figure 5 As shown, the electronic device 500 may include a processing device (e.g., a central processing unit, a graphics processing unit, etc.) 501, which can perform various appropriate actions and processes according to a program stored in a read-only memory (ROM) 502 or a program loaded from a storage device 508 into a random access memory (RAM) 503. In the RAM 503, various programs and data required for the operation of the electronic device 500 are also stored. The processing device 501, the ROM 502, and the RAM 503 are connected to each other via a bus 504. An input / output (I / O) interface 505 is also connected to the bus 504.

[0336] Typically, the following devices may be connected to the I / O interface 505: an input device 506 including, for example, a touch screen, a touch pad, a keyboard, a mouse, a camera, a microphone, an accelerometer, a gyroscope, etc.; an output device 507 including, for example, a liquid crystal display (LCD), a speaker, a vibrator, etc.; a storage device 508 including, for example, a magnetic tape, a hard disk, etc.; and a communication device 509. The communication device 509 may allow the electronic device 500 to communicate with other devices wirelessly or by wire to exchange data. Although Figure 5 The electronic device 500 is shown with various devices, but it should be understood that it is not required to implement or possess all the devices shown. More or fewer devices may be implemented or possessed instead. Figure 5 Each block shown in the figure may represent one device, or may represent multiple devices as required.

[0337] In particular, according to the present embodiment, the process described above with reference to the flowchart can be implemented as a computer software program. For example, the present embodiment includes a computer program product, which includes a computer program or instruction carried on a computer-readable medium, and the computer program or instruction contains a program code for executing the method shown in the flowchart. In the present embodiment, the computer program can be downloaded and installed from the network through the communication device 509, or installed from the storage device 508, or installed from the ROM 502. When the computer program is executed by the processing device 501, the following steps can be performed: acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data; performing beamforming filtering processing on the first seismic data to obtain second seismic data; performing slowness filtering processing on the second seismic data to obtain target seismic data.

[0338] It should be noted that the computer-readable medium in the present embodiment may be a computer-readable signal medium or a computer-readable storage medium or any combination of the two. The computer-readable storage medium may be, for example, but not limited to, an electrical, magnetic, optical, electromagnetic, infrared, or semiconductor system, device or device, or any combination of the above. More specific examples of computer-readable storage media may include, but are not limited to: an electrical connection with one or more wires, a portable computer disk, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber, a portable compact disk read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the above. In the present embodiment, the computer-readable storage medium may be any tangible medium containing or storing a program that may be used by or in combination with an instruction execution system, device or device. In the present embodiment, the computer-readable signal medium may include a data signal propagated in a baseband or as part of a carrier wave, which carries a computer-readable program code. This propagated data signal may take a variety of forms, including but not limited to electromagnetic signals, optical signals, or any suitable combination of the above. The computer readable signal medium may also be any computer readable medium other than a computer readable storage medium, which may send, propagate or transmit a program for use by or in conjunction with an instruction execution system, apparatus or device. The program code contained on the computer readable medium may be transmitted using any suitable medium, including but not limited to: wires, optical cables, RF (radio frequency), etc., or any suitable combination of the above.

[0339] In some embodiments, the client and the server may communicate using any currently known or future developed network protocol such as HTTP (HyperText Transfer Protocol), and may be interconnected with any form or medium of digital data communication (e.g., a communication network). Examples of communication networks include a local area network ("LAN"), a wide area network ("WAN"), an internet (e.g., the Internet), and a peer-to-peer network (e.g., an ad hoc peer-to-peer network), as well as any currently known or future developed network.

[0340] The computer-readable medium may be included in the device, or may exist independently without being installed in the electronic device. The computer-readable medium carries one or more programs. When the one or more programs are executed by the electronic device, the electronic device performs the following steps: acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data; performing beamforming filtering processing on the first seismic data to obtain second seismic data; and performing slowness filtering processing on the second seismic data to obtain target seismic data.

[0341] The computer program code for performing the operation of the present embodiment can be written in one or more programming languages ​​or a combination thereof, including object-oriented programming languages, such as Java, Smalltalk, C++, and conventional procedural programming languages, such as "C" language or similar programming languages. The program code can be executed entirely on the user's computer, partially on the user's computer, as an independent software package, partially on the user's computer and partially on a remote computer, or entirely on a remote computer or server. In the case of a remote computer, the remote computer can be connected to the user's computer through any type of network, including a local area network (LAN) or a wide area network (WAN), or can be connected to an external computer (e.g., using an Internet service provider to connect through the Internet).

[0342] The flow chart and block diagram in the accompanying drawings illustrate the possible architecture, function and operation of the system, method and computer program product according to various embodiments of the present disclosure. In this regard, each square box in the flow chart or block diagram can represent a module, a program segment or a part of a code, and the module, the program segment or a part of the code contains one or more executable instructions for realizing the specified logical function. It should also be noted that in some implementations as replacements, the functions marked in the square box can also occur in a sequence different from that marked in the accompanying drawings. For example, two square boxes represented in succession can actually be executed substantially in parallel, and they can sometimes be executed in the opposite order, depending on the functions involved. It should also be noted that each square box in the block diagram and / or flow chart, and the combination of the square boxes in the block diagram and / or flow chart can be implemented with a dedicated hardware-based system that performs a specified function or operation, or can be implemented with a combination of dedicated hardware and computer instructions.

[0343] The modules described in this embodiment may be implemented by software or hardware. The modules described may also be set in a processor, for example, may be described as:

[0344] Acquisition module, first processing module and second processing module. For example, the acquisition module can also be described as "a module for acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data".

[0345] The functions described above herein may be performed at least in part by one or more hardware logic components. For example, without limitation, exemplary types of hardware logic components that may be used include: field programmable gate arrays (FPGAs), application specific integrated circuits (ASICs), application specific standard products (ASSPs), systems on chips (SOCs), complex programmable logic devices (CPLDs), and the like.

[0346] The above descriptions are only some preferred embodiments of the present disclosure and an explanation of the technical principles used. Those skilled in the art should understand that the scope of the invention involved in the embodiments of the present disclosure is not limited to the technical solutions formed by a specific combination of the above-mentioned technical features, but should also cover other technical solutions formed by any combination of the above-mentioned technical features or their equivalent features without departing from the above-mentioned inventive concept. For example, the above-mentioned features are replaced with the technical features with similar functions disclosed in the embodiments of the present disclosure (but not limited to) and the technical solutions formed.

Claims

1. A residual multiple wave suppression method, characterized in that: include: Acquiring original seismic data, and performing pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data; Performing a beamforming filter process on the first seismic data to obtain second seismic data; The second seismic data is subjected to slowness filtering to obtain target seismic data.

2. The method according to claim 1, characterized in that Performing a beamforming filter process on the first seismic data to obtain second seismic data includes: Decomposing the first seismic data based on a preset first calculation formula to obtain an initial effective signal, an effective signal filter operator, initial coherent noise information including multiple waves, a coherent noise filter operator and initial random noise information, wherein the decomposition model is obtained by combining a fitting signal model and a regular coherent noise model; Generate a corresponding prediction identity matrix according to the initial random noise; Calculating a target filter of a multi-constraint problem according to the effective signal filter operator, the coherent noise filter operator, the prediction unit matrix, and a second calculation formula constructed based on a Lagrange multiplier method; The second seismic data is calculated based on the original seismic data and the filter.

3. The method according to claim 2, characterized in that The first calculation formula includes: x=Bp+Cv+u Among them, x represents the first seismic data, B represents the effective signal filtering operator, p represents the predicted effective signal, C represents the coherent noise filtering operator, v represents the initial coherent noise information, and u represents the initial random noise in matrix format.

4. The method according to claim 2, characterized in that: The second calculation formula includes: H p =G(A H Q -1 A) -1 A H Q -1 Among them, H p represents the target filter, G = (I, 0), I represents the unit matrix, A = (B, C), Q = E [uu] H , E[] represents the mathematical expectation value of the matrix; In addition, the constraint condition for the coherent noise in the second calculation formula is that the coherent noise has zero or minimum response, and the constraint condition for the predicted unit matrix is ​​that the random noise is normally distributed and each channel has the same energy.

5. The method according to claim 1, characterized in that Performing slowness filtering on the second seismic data to obtain target seismic data includes: Extracting the third seismic data between the first layer bottom interface generating the interlayer multiple waves and the preset second layer bottom interface from the second seismic data; Performing slowness smoothing processing on the third seismic data to obtain fourth seismic data; The third seismic data in the second seismic data is replaced by the fourth seismic data to obtain the target seismic data.

6. The method according to claim 5, characterized in that The step of performing slowness smoothing on the third seismic data to obtain fourth seismic data comprises: Flattening the third seismic data along the bottom interface of the first layer to obtain flattened seismic data; Performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data; The instantaneous slowness of each data point in the event-enhanced seismic data is enhanced according to the difference calculation method, and each instantaneous slowness is smoothed to obtain smoothed seismic data; Performing event-weakening transformation on the smoothed seismic data to obtain event-weakening seismic data; The event-weakening seismic data is subjected to anti-layer flattening processing to obtain the fourth seismic data.

7. The method according to claim 6, characterized in that The step of flattening the third seismic data along the first layer bottom interface to obtain flattened seismic data comprises: The third seismic data is flattened along the first layer bottom interface according to a preset third calculation formula to obtain flattened seismic data, wherein the third calculation formula includes: Q(x,m)=P(x,t=m+T1(x)) Among them, P() represents the third seismic data before leveling, Q(x,m) represents the leveled seismic data after leveling, t represents the time coordinate before leveling, m represents the time coordinate after leveling, x represents the coordinate of the common reflection point in the direction of the comprehensive survey line on the profile, and T1 represents the bottom interface of the first layer that produces interlayer multiple waves.

8. The method according to claim 7, characterized in that The step of performing event-enhanced transformation on the flattened seismic data to obtain event-enhanced seismic data includes: According to a preset fourth calculation formula, the flattened seismic data is transformed into event-enhanced seismic data to obtain event-enhanced seismic data, wherein the fourth calculation formula includes: Among them, R(x,m) represents the data after the phase axis is enhanced, Q(x,m) represents the flattened seismic data, and Q H (x,m) represents the Hilbert transform of the data Q(x,m).

9. The method according to claim 8, characterized in that The step of smoothing each instantaneous slowness to obtain smoothed seismic data includes: According to the preset fifth calculation formula, each instantaneous slowness is smoothed to obtain smoothed seismic data, wherein the fifth calculation formula includes: U(x,m)=V(x,m)R(x,m) Among them, U(x,m) represents the data after smoothing based on instantaneous slowness, J(x,m) represents instantaneous slowness, V() represents the smoothing processor based on instantaneous slowness, J represents instantaneous slowness, and △J represents the preset edging amount during smoothing based on instantaneous slowness.

10. A residual multiple wave suppression system, characterized in that: include: An acquisition module is used to acquire original seismic data, and perform pre-stack prediction deconvolution processing and high-precision Radon transform processing on the original seismic data to obtain first seismic data; A first processing module, configured to perform a beamforming filter process on the first seismic data to obtain second seismic data; The second processing module is used to perform slowness filtering on the second seismic data to obtain target seismic data.

11. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the steps of the method according to any one of claims 1 to 9 are implemented.

12. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 9 are implemented.

Images