阅读说明：本技术 一种三维叠前地震数据层间多次波预测和压制方法 (Three-dimensional prestack seismic data interbed multiple prediction and suppression method ) 是由 胡天跃 黄建东 安圣培 谢飞 刘韬 朱成宏 于 2021-08-23 设计创作，主要内容包括：本发明公布了一种三维叠前地震数据层间多次波预测和压制方法,通过三维叠前地震数据层间多次波预测和匹配压制实现三维叠前地震数据层间多次波衰减；包括：对含有层间多次波的叠前数据进行预处理；在时间域根据产生层间多次波层位界面确定一次反射波同相轴；在频率-空间域构造虚一次反射波信号,将虚一次反射波信号与地震数据褶积构造层间多次波；对预测的层间多次波进行匹配处理,将层间多次波进行多次迭代预测和匹配压制,压制三维叠前地震数据层间多次波,突出有效波。本发明能够避免叠后层间多次波同相轴与一次反射波同相轴出现混叠而无法进行有效的层间多次波预测衰减,进一步提高地震资料处理的精度。(The invention discloses a method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data, which realizes interbed multiple attenuation of the three-dimensional prestack seismic data through interbed multiples prediction and matched suppression of the three-dimensional prestack seismic data; the method comprises the following steps: preprocessing prestack data containing interbed multiples; determining a primary reflection wave homophase axis according to the interlayer multiple layer position interface in a time domain; constructing a virtual primary reflected wave signal in a frequency-space domain, and performing convolution on the virtual primary reflected wave signal and seismic data to construct an interbed multiple; and performing matching processing on the predicted interbed multiples, performing multiple iterative prediction and matching suppression on the interbed multiples, suppressing the interbed multiples of the three-dimensional pre-stack seismic data, and highlighting the effective waves. The invention can avoid the aliasing of the in-phase axis of the interbedded multiples and the in-phase axis of the primary reflected wave after the stacking, so that the effective interbedded multiples can not be predicted and attenuated, and further improve the precision of seismic data processing.)

1. A three-dimensional prestack seismic data interbed multiple prediction and suppression method, realize three-dimensional prestack seismic data interbed multiple attenuation through three-dimensional prestack seismic data interbed multiple prediction and three-dimensional prestack seismic data interbed multiple match suppression; the method comprises the following steps: firstly, preprocessing prestack data containing interbed multiples to obtain regularized prestack data only containing primary reflected waves and interbed multiples; then determining a primary reflection wave homophase axis according to the interlayer multiple layer position interface in a time domain, and recording a time line of the homophase axis; constructing a virtual primary reflected wave signal in a frequency-space domain by using a virtual homophase axis method, and then performing convolution on the virtual primary reflected wave signal and original seismic data to construct an interlaminar multiple; furthermore, predicted interbed multiples are subjected to matching processing through an effective matching filtering method, and the purposes of suppressing the interbed multiples of the three-dimensional prestack seismic data and highlighting effective waves are achieved; wherein:

A. the method for predicting the interbed multiples of the three-dimensional prestack seismic data comprises the following steps:

A1) regularizing three-dimensional pre-stack seismic data; the method comprises the following steps:

defining a two-dimensional rectangular observation system received on the surface, determining the position relation between a shot point and a demodulator probe, providing each channel of seismic data according to field acquired channel head information, and carrying out regularization processing on the field acquired seismic data;

A2) preparing interbed multiple suppression data; the method comprises the following steps:

performing static correction on the seismic data; correcting the regularized seismic data to a datum plane; simultaneously picking up direct wave signals for cutting off, wherein the direct waves do not participate in prediction of interbed multiples of the prestack seismic data; the seismic data used are three-dimensional pre-stack seismic data which are subjected to preprocessing including surface multiple suppression and contain primary reflection waves and interbed multiples;

A3) reading the preprocessed shot gather data, and preparing data required by three-dimensional prestack seismic data interbed multiple prediction;

determining a primary reflected wave signal generating an interlayer multiple layer position, surrounding an effective signal by a time frame, and recording time frame information surrounding an in-phase axis of the effective signal; the method specifically comprises the following steps:

A31) determining a primary reflected wave generating an interlayer multiple layer position from top to bottom in shot gather data, adopting a variable time frame to surround a primary reflected wave homophase axis, recording time frame data information, traversing information of the reflection layer corresponding to the primary reflected wave on all shot gathers, and recording; recording all the layers of the shot gather generating the interbed multiples as input data of interbed multiple prediction;

A32) then regularizing the seismic data;

assembling the processed three-dimensional pre-stack seismic data and the primary reflection wave time line data according to a transmission-reflection-transmission WRW model to obtain regularized WRW model data, so that interlayer multiple prediction calculation is facilitated; the WRW model is a two-dimensional data matrix, wherein one dimension represents common shot point seismic data, and the other dimension represents common geophone point seismic data;

A4) predicting interlayer multiples of the three-dimensional pre-stack seismic data;

dividing shot gather data into an upper data body and a lower data body according to time line data, and directly and accurately predicting interlayer multiples of the three-dimensional pre-stack seismic data; the specific process of prediction comprises the following steps:

A41) dividing the regularized WRW model data into an upper part and a lower part according to the time line data of the primary reflected wave;

A42) calculating to obtain virtual homophase axis data;

A43) performing convolution operation on the virtual event axis data and the lower part of seismic data, and predicting all orders of interbed multiples corresponding to the layer; wherein, the dip angle of the top-bottom interface of the stratum generating the interbed multiples is not required;

B. performing interlayer multiple matching suppression on the three-dimensional pre-stack seismic data, and performing multiple iterative prediction and matching suppression on the interlayer multiple in the three-dimensional pre-stack seismic data; the method comprises the following steps:

B1) designing a matched filter by using a first norm and a second norm matching method;

the matched filter includes:

the single-channel filter module is used for performing norm matching on only one channel in the seismic data to obtain a single-channel filtering factor;

the multi-channel filter module is used for obtaining a multi-channel matched filtering factor by taking a certain channel In the seismic data as a central channel and taking an adjacent seismic channel In the In-line direction as an auxiliary channel;

the square filter module takes one of the seismic data as a central channel and takes the seismic channels around as auxiliary data to obtain a filter factor; wherein the peripheral seismic channels comprise seismic channels In-line and Cross-line directions;

B2) matching the predicted data with the actual data by using a matched filter to generate multi-wave data between pre-stack layers after matching;

the matching specifically comprises the following steps: calculating matched interbed multiples by multiplying the predicted prestack interbed multiples according to the obtained filter factors,

B3) finally, pre-stack seismic data after interlayer multiple wave attenuation is obtained, and seismic data after interlayer multiple wave suppression is obtained; carrying out multiple iterative cycles on the prediction and matching pressing process of the interbed multiples, and finally obtaining effective wave data after all interbed multiples are pressed;

the iterative loop comprises an outer loop and an inner loop; wherein, the iterative process generated by moving the primary reflected wave is an external loop and is used for suppressing the interbed multiples generated at different layers; multiple prediction and matching pressing iteration processes in the outer loop calculation process are inner loops, and the inner loops are used for pressing all order interbed multiples generated by a specific layer;

through the steps, the attenuation of the multiple between the layers of the seismic data before the dimension stack is realized.

2. The method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data as recited in claim 1, wherein the step a4) of performing three-dimensional prestack interbed multiples prediction includes constructing a virtual event and predicting interbed multiples; the method comprises the following specific steps:

the virtual event and the interbed multiples constructed in the corresponding frequency domain using the virtual event method are respectively expressed as:

D_I(x_s,y_s,ω,x_r,y_r)＝∫∫D_V(x_s,y_s,ω,x,y)D′₀(x,y,ω,x_r,y_r) dxdy type (2)

Wherein D is_V(x_s,y_s,ω,x_r,y_r) Representing constructed virtual homophase axes, D_I(x_s,y_s,ω,x_r,y_r) Representing constructed interbed multiples, D₀(x,y,ω,x_r,y_r) And D'₀(x,y,ω,x_r,y_r) Respectively represent an upper part and a lower part of original data which are not overlapped with each other,represents D₀(x) is a conjugate of (C)_s,y_s) Representing source coordinates, (x)_r,y_r) Representing the coordinates of the receiving point;

2) distinguishing data to realize accurate prediction of three-dimensional interbed multiples;

when data are distinguished, the upper part and the lower part do not contain the same phase axis, otherwise, the inaccuracy of prediction is increased; meanwhile, the lower part data body D 'in the seismic data'₀Not only has a reflected wave, but also includes an interbed multiple, and the predicted result includes interbed multiples; the seismic data are data of which the surface multiples are completely suppressed;

according to the formulas (1) and (2), accurate prediction of interbed multiples can be realized by using prestack seismic data without prior information and assumption on an underground medium.

3. The method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data as recited in claim 2, wherein the discretizing of the equations (1) and (2) facilitates the processing of discretized seismic data signals;

the discretization and iterative computation process in the iteration step is represented as:

wherein, P_m(x_r,y_r,z₀；x_s,y_s,z₀(ii) a ω) represents a shot point (x) from the surface_r,y_r) To the surface inspection point (x)_s,y_s) The multiple between the layers of (a) and (b),representing seismic data from which all interbed multiples waves associated with the first layer have been removed; s^-1(ω) represents the inverse matrix of the source wavelet;represents the subsurface impulse response without interbed multiples; r (z)₁) Representing a wave field downward scattering operator positioned at the first layer interface, wherein the downward scattering operator of the earth surface is a negative identity matrix; the formula (8) represents: the interbed multiples associated with the first layer interface are equivalent to convolution of the wavefield pulse reflection from which the interbed primary reflection is removed, the layer down-scatter operator, and the wavefield that does not contain interbed multiples;

rewrite equation (8) so that all wavefields to the right of the equation are surface excitations and surface receptions, expressed as:

recording:

substituting the virtual pulse reflection without multiple information into an interlayer multiple model to obtain formula (11):

the imaginary pulse reflection without multiple information in equation (11) is the modified imaginary iso-axis, which is equivalent to the convolution of the primary pulse reflection associated with the layer interface and the wavefield pulse without the interbed multiple associated with the layer;

converting primary waves into interbed multiples and converting n interbed multiples into n +1 interbed multiples when self-adapting virtual homophase axes and seismic data are convolved; the expansion according to the Noelman series can be obtained:

wherein:

by the iterative method represented by the formulas (12) to (14), effective suppression of interbed multiples of the three-dimensional prestack seismic data is achieved.

4. The method as claimed in claim 1, wherein step B is specifically performed by matching interbed multiples using a plurality of modules, including a single-channel norm matching module, a multi-channel norm matching module, and a square norm matching module, to perform interbed multiples matching and suppressing on the three-dimensional prestack seismic data, and comprises:

utilizing a two-norm matching method; based on the principle of minimum energy, the energy is minimum after the predicted three-dimensional interlayer multi-wave data is subtracted from the input three-dimensional interlayer multi-wave data;

for the single pass two-norm method, the minimum objective function is:

the filter operator a is subjected to partial derivation on both sides of the formula (3) and is made equal to 0, and the formula (3) is minimized, and then converted into a linear equation:

M^TMa＝M^Td (4)

the coefficient matrix on the left side of the formula (4) is a Topelitz matrix, and the calculation is simple;

solving a filtering operator a through a Levinson recursive algorithm or SVD decomposition, and directly performing division operation under the condition of stable equation;

one norm matching method is as follows:

e＝|d-Ma₁| (5)

the objective function is a singular function, transforming equation (3) into a weighted two-norm minimization, represented as equation (6):

wherein, W represents a weighting matrix,ε＝max(d)/100；r_j＝d-Ma。

5. the method for predicting and suppressing the interbed multiples of the three-dimensional prestack seismic data as claimed in claim 4, wherein the matching effect is improved by changing the size of epsilon to adjust the weight of the first norm and the second norm calculation, and the single-channel first norm calculation equation is obtained by:

M^TW^TWMa＝M^TW^TWd (7)

all seismic data after the three-dimensional interbed multiples are attenuated can be output.

6. The method of predicting and suppressing interbed multiples of three-dimensional prestack seismic data according to claim 4, wherein the length of the filter operator a is selected; the length is more stable when the length is shorter, and the length is more sensitive when the length is longer.

Technical Field

The invention belongs to the technical field of geophysical exploration, relates to a prediction and suppression method of stratum interbed multiples, and particularly relates to a prediction and matching suppression method of three-dimensional prestack seismic data interbed multiples.

Background

In land seismic exploration, the problem of interbed multiples is one of the major factors affecting seismic data processing. Especially in complex deep regions, the multiple waves generated by the overlying high-low speed interlayer can cover the primary reflected wave signals of the underlying stratum, which is not beneficial to the processing and interpretation of seismic data. Western regions of china have a large number of carbonate reservoirs, but underground structures are very complex, and faults and holes develop, which can generate strong interbed multiples. The existence of the interbed multiples makes the seismic wave field of the area more complex, influences the interpretation work of geophysicists and geologists on the seismic data processing result, and greatly increases the uncertainty of reservoir interpretation.

The generation mechanism of the interbed multiples is very complex, so that the prediction difficulty is high. The theory of the interbed multiples correlation is not perfect either with the theory of the surface multiples correlation. Current methods for attenuating interbed multiples are mainly classified into three major categories: the first is a filtering method, such as a prediction deconvolution method, a Radon transform method and the like, and the method suppresses the interbed multiples by the difference of the frequency, the period and the normal time difference of prestack data of the primary waves and the interbed multiples; the second category is based on model-driven wave field continuation predictive subtraction, such as confocal methods [1 ]; the third type is data-driven prediction subtraction method based on wave equation, which includes backscattering series method, virtual homophase axis method [2], etc. The filtering method is mainly used for suppressing the interbed multiples by means of obvious difference between the effective waves and the interbed multiples, and the application effect on actual seismic data is poor; the model-driven interbed multiple attenuation method is too dependent on a macroscopic model, and the accurate velocity model cannot be obtained by actual seismic data; the pure data-driven interbed multiple pressing method does not depend on a model or obvious characteristic difference between primary waves and multiple waves, and has practical significance on attenuation of interbed multiple waves.

The virtual event axis method is proposed by the professor Ikelle in 2006 [3], and like the surface multiple removal method proposed by the professor Berkout, the interlayer multiple prediction part adopts a complete data driving mode, does not need any underground prior information, and can accurately predict the interlayer multiple in the seismic record. The virtual event method is therefore a pure data driven method of interbed multiple attenuation. Wujing et al effectively suppress interbed multiples in seismic data in the south yellow sea area based on a virtual homophase axis method. Liujiahui et al (2018) further provides an iterative virtual primary wave method for suppressing interbed multiples based on a virtual event method, and the method is applied to actual seismic data processing in a Tarim area, so that the seismic imaging precision of the area is greatly improved [4 ].

In recent years, with the rapid development of three-dimensional seismic data acquisition and computer performance, three-dimensional multiple suppression algorithms have been further developed and perfected, such as a three-dimensional surface multiple removal method [5] [6], a three-dimensional radon transformation method [7], a three-dimensional sparse inversion primary wave method [8], and the like. These algorithms all achieve some effect in three-dimensional seismic data. The method for pressing interbed multiples in seismic data by using the virtual event method has obvious advantages, and although Liujiahui and the like (2018) derive the theoretical solution of the three-dimensional virtual event method, pressing interbed multiples in three-dimensional actual seismic data has not made substantial progress. This also limits the application of the three-dimensional virtual event method to the compression of interbed multiples in three-dimensional seismic data. In theory, the three-dimensional algorithm avoids various assumptions of the two-dimensional algorithm, is consistent with the dimension of actual seismic data, and has higher prediction precision than the two-dimensional algorithm. Therefore, the three-dimensional multiple theory better conforms to the actual situation and has better application value. However, at present, theoretical research based on a virtual event axis method stays in a two-dimensional stage, and is not directly applied to three-dimensional prestack interbed multiple pressing. With the normalization of the acquisition and processing of three-dimensional seismic data in the future, the research of a three-dimensional interbed multiple pressing method with higher precision is urgent.

Therefore, the technical scheme of establishing the three-dimensional virtual homodyne axis method for suppressing the interbed multiples has very important significance for improving the processing effect of the three-dimensional seismic data and the seismic imaging precision.

Reference documents:

[1]Berkhout A J,and Verschuur D J.2005.Removal of internal multiples with the common-focus-point(CFP)approach:Part 1-Explanation of the theory.Geophysics 70(3):V45-V60.

[2] wujing, Wuzhiqiang, Hutianyun, He Yuhua, Wanpu, Yan Guijing, Lilin, 2013. pressing multiple waves between earthquake layers based on the construction of virtual homophase axis, geophysical reports 56(3):985, 994).

[3]Ikelle L T.2006.A construct of internal multiples from surface data only:the concept of virtual seismic events.Geophysical Journal International 164(2):383-393.

[4]Liu J H,Hu T Y,Peng G X and Cui Y F.2018.Removal of internal multiples by iterative construction of virtual primaries,Geophysical Journal International 215(1):81-101.

[5]Berkhout A J,and Verschuur D J.1997.Estimation of multiple scattering by iterative inversion,Part I:Theoretical considerations.Geophysics 62(5):1586-1595.

[6] Stone, royal vitamin red, plum blossom, well flood 2013, three-dimensional surface multiple prediction method based on wave equation, geophysical report 56(2), 2023 and 2032.

[7] Tanghua, maoweiwei 2014.3D high-order parabolic Radon transform seismic data amplitude-preserving reconstruction geophysical newspaper 57(9):2918-2927

[8]van Dedem E J,and Verschuur D J.2002.3D surface-related multiple prediction using sparse inversion:Experience with field data.SEG Technical Program Expanded Abstracts:2094-2097.

Disclosure of Invention

Aiming at the problems in the prior art and aiming at overcoming the bottleneck of the three-dimensional prestack seismic data interbed multiple suppression application technology, the invention provides a three-dimensional prestack seismic data interbed multiple attenuation method based on the virtual homomorphic axis principle. The invention discloses a method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data. The method comprises the steps of accurately picking up a primary reflected wave signal of a certain layer, constructing a virtual primary reflected wave signal of the layer (which cannot be recorded in actual seismic records), constructing all orders of interbed multiples related to the layer by utilizing a convolution mode, and suppressing the interbed multiples by utilizing an adaptive matched filter.

The core of the invention is: the invention fully considers the problem of bottleneck of prediction and suppression of interbed multiples in actual seismic data, summarizes and finds that the interbed multiples denoising effect is better than that of poststack interbed multiples denoising effect, accurately predicts and attenuates interbed multiples directly on three-dimensional shot gather data, avoids aliasing of a poststack interbed multiples homophase axis and a primary reflected wave homophase axis, cannot perform effective interbed multiples prediction attenuation, and further improves the accuracy of seismic data processing.

In order to solve the technical problems, the technical scheme provided by the invention is as follows:

a three-dimensional prestack seismic data interbed multiple attenuation method comprises a three-dimensional prestack seismic data interbed multiple prediction method and a three-dimensional prestack seismic data interbed multiple matching suppression method.

In a first aspect, the present invention provides a method for predicting an interbed multiple of three-dimensional prestack seismic data, including:

regularizing three-dimensional pre-stack seismic data, defining a two-dimensional rectangular observation system received on the surface, determining the position relationship between a shot point and a demodulator probe, providing each channel of seismic data according to field acquired channel head information, and regularizing the field acquired seismic data to obtain regularized three-dimensional pre-stack seismic data.

Preparing interbed multiple suppressed data, wherein static correction needs to be carried out on seismic data before interbed multiple prediction is carried out, and the regularized three-dimensional prestack seismic data is corrected to a datum plane, so that subsequent seismic data can be conveniently processed; meanwhile, direct wave signals need to be picked up for cutting, and the direct waves cannot participate in prediction of interbed multiples of the prestack seismic data. The seismic data used in the invention has been preprocessed, including surface multiple suppression, to obtain three-dimensional pre-stack seismic data containing primary reflections and interbed multiples.

And reading the three-dimensional pre-stack seismic data (shot gather data) processed by the steps, and preparing data required by interlayer multiple prediction of the three-dimensional pre-stack seismic data. Determining a primary reflected wave signal generating an interlayer multiple layer position, surrounding the effective signal by a time frame, and recording time frame information surrounding the in-phase axis of the effective signal. The method specifically comprises the following steps:

determining primary reflection waves generating interlayer multiple wave positions from top to bottom in shot gather data, surrounding a primary reflection wave homophase axis by a variable time frame, recording time frame data information, traversing information of the reflection layer corresponding to the primary reflection waves on all shot gathers, and recording. And recording all the layers of the shot gather generating the interbed multiples as input data of interbed multiple prediction.

Then, regularizing the seismic data, wherein the conventional format seismic data are inconvenient to use for interlayer multiple prediction calculation, and the three-dimensional pre-stack seismic data and the primary reflected wave time line data need to be assembled according to a propagation-reflection-propagation (WRW) model to obtain regularized WRW model data; wherein the WRW model is a two-dimensional data matrix, wherein one dimension represents common shot point seismic data and the other dimension represents common geophone point seismic data.

Predicting interlayer multiples of the three-dimensional pre-stack seismic data, dividing shot gather data into an upper data body and a lower data body according to time line data, and directly and accurately predicting the interlayer multiples of the three-dimensional pre-stack seismic data; the specific process of prediction comprises the following steps:

dividing the regularized WRW model data into an upper part and a lower part according to the time line data of the primary reflection waves, firstly carrying out correlation calculation to obtain virtual homophase axis data, then carrying out convolution operation on the virtual homophase axis data and the seismic data of the lower part, and predicting all orders of interbed multiples corresponding to the horizon. Wherein, no requirement is made on the dip angle of the top-bottom interface of the stratum generating the interbed multiples.

In a second aspect, the present invention provides a method for matching and suppressing interbed multiples in three-dimensional prestack seismic data, which performs iterative prediction and matching suppression on interbed multiples in three-dimensional prestack seismic data for multiple times, including:

and designing a matched filter by using a first norm and second norm matching method to match the predicted data with the actual data, generating multi-time wave data between pre-stack layers after matching, and finally obtaining pre-stack seismic data after the multi-time wave between the layers is attenuated.

In one embodiment, a matched filter includes:

the single-channel filter module is used for performing norm matching on only one channel in the seismic data to obtain a single-channel filtering factor;

the multi-channel filter module is used for obtaining a multi-channel matched filtering factor by taking a certain channel In the seismic data as a central channel and taking an adjacent seismic channel In the In-line direction as an auxiliary channel;

and the square filter module also takes a certain channel In the seismic data as a central channel, and takes the seismic channels (In-line and Cross-line direction seismic channels) around the central channel as auxiliary data to obtain a filter factor.

According to the obtained filtering factor, multiplying the predicted prestack interbed multiples to calculate matched interbed multiples, and obtaining seismic data after suppressing the interbed multiples; and carrying out multiple iterative cycles on the prediction and matching pressing process of the interbed multiples, and finally obtaining the effective wave data after all the interbed multiples are pressed.

The whole process of the invention relates to two iteration loops, and the iteration process generated by moving a reflected wave is called as an outer loop and mainly aims to suppress the multiple waves among layers generated at different layers; in the calculation process of the outer loop, a plurality of times of prediction and matching suppression iterative processes are called as the inner loop, and the inner loop is mainly used for suppressing all orders of interbed multiples generated by a specific horizon.

From the above description, it can be seen that the present invention provides a method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data. Firstly, preprocessing prestack data containing interbed multiples to obtain regularized prestack data only containing primary reflected waves and interbed multiples, then determining a primary reflected wave homophasic axis according to a layer position interface generating the interbed multiples in a time domain, recording a time line of the homophasic axis, constructing a virtual primary reflected wave signal in a frequency-space domain by using a virtual homophasic axis principle, and then performing convolution on the virtual primary reflected wave signal and original seismic data to construct the interbed multiples; furthermore, predicted interbed multiples are matched through an effective matched filtering method, so that the aims of suppressing interbed multiples of the three-dimensional prestack seismic data and highlighting effective waves are fulfilled. The method directly performs prediction attenuation on the interbed multiples of the three-dimensional prestack seismic data by using a virtual homophase axis method. The method can directly denoise the seismic data before the stack without depending on the obvious difference between the primary reflection wave and the interbed multiples and an underground accurate velocity model, and avoids various problems of hypothesis and difficult processing of the seismic data after the stack. The invention is a pure data driving method based on wave equation, has good practicability for the interlayer multiple pressing of the actual three-dimensional prestack seismic data, can be widely applied to the industry, greatly and effectively presses the interlayer multiple, and highlights the proportion of the effective wave in the prestack seismic data. By adopting the three-dimensional prestack seismic data interbed multiple prediction method and the matching suppression method, the attenuation of the three-dimensional prestack seismic data interbed multiple can be realized.

In modern three-dimensional seismic exploration, hundreds of thousands of traces of seismic data are involved. The method provided by the invention can accurately predict the multiple between the three-dimensional prestack layers only by two data volumes, reduces the requirement on the memory of a computer and greatly improves the calculation efficiency.

Compared with the prior art, the invention has the beneficial effects that:

the method realizes the three-dimensional prestack seismic data interbed multiple suppression, accurately predicts the three-dimensional prestack interbed multiples based on the virtual event method, and greatly improves the calculation efficiency. The invention has the advantages that:

firstly, accurately predicting the interbed multiples of all orders of an interbed multiple layer by picking up the information of primary reflected waves of the interbed multiple layer;

secondly, the pressing is carried out without depending on the obvious difference between the multiple and the primary wave between layers, so that the requirement on data is reduced;

(III) the method does not need any prior information underground and is a pure data-driven method;

fourthly, calculating the multiple between the prestack layers based on a correlation and convolution mode, thereby reducing the calculation amount;

(V) utilizing an effective filter to effectively suppress the interbed multiples;

the method for predicting and suppressing the interbed multiples of the three-dimensional prestack seismic data improves the proportion of effective signals in the prestack data and indirectly improves the resolution of the post-stack effective signals;

the invention has good industrial value and can be widely used for processing actual three-dimensional seismic data.

Drawings

In order to more clearly describe the technical scheme of the invention, the following description is made by simply describing some drawings in the embodiment.

Fig. 1 is a flow chart of a method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data according to an embodiment of the present invention.

FIG. 2 is a schematic diagram of a two-dimensional rectangular observation system defined in an embodiment of the present invention;

wherein, (a) is common shot point seismic data; (b) is common receiver seismic data.

Fig. 3 is a schematic diagram illustrating an implementation process of the three-dimensional prestack interbed multiple prediction method according to an embodiment of the present invention.

Fig. 4 is a detailed flowchart of an iteration method in the embodiment of the present invention.

FIG. 5 is a velocity model diagram of an embodiment of the present invention.

FIG. 6 is pre-stack seismic data obtained from the model of FIG. 5, in accordance with an embodiment of the present invention.

FIG. 7 is a diagram of a shot gather data (seismic data regularized according to the WRW model) obtained according to the model of FIG. 5 in accordance with an exemplary embodiment of the present invention.

FIG. 8 is a diagram of predicted data for a shot virtual event axis according to an embodiment of the present invention.

FIG. 9 shows a diagram of a predicted interval multiple before a shot stack according to an embodiment of the present invention.

FIG. 10 shows a diagram of a shot prestack interbed multiples for matching in accordance with an embodiment of the present invention.

FIG. 11 is a diagram of prestack seismic data after multiple interbed wave passes in accordance with an exemplary embodiment of the present invention.

FIG. 12 is a diagram of prestack seismic data after multiple interbed wave passes with iterative methods in accordance with an exemplary embodiment of the present invention.

Detailed Description

The technical embodiments of the present invention will be more clearly described below by way of examples with reference to the accompanying drawings, without limiting the scope of the invention in any way. The embodiments of the present invention are provided so that this disclosure will be fully conveyed to those skilled in the art.

Referring to fig. 1, an embodiment of the present invention provides a flowchart of an embodiment of a method for predicting and suppressing interbed multiples of three-dimensional prestack seismic data. The method comprises the following steps 101, 102, 103, 104 and 105.

Step 101, acquiring time domain three-dimensional shot gather seismic data containing interbed multiples and preprocessing the data.

In one possible implementation, acquiring information corresponding to three-dimensional pre-stack seismic data, including the abscissa and the ordinate of a shot point and the abscissa and the ordinate of a demodulator probe, defining observation system information (fig. 2), where a triangle represents a geophone, a quadrangle star represents a seismic source, a diagram a represents a certain shot excitation, all demodulator probes receive (representing common shot point seismic data), a diagram b represents a certain demodulator probe reception, and other points excite the seismic source (representing common demodulator probe seismic data); the seismic data of each cannon is preprocessed, the shot point data and the demodulation point data are corrected to the same horizontal reference surface, subsequent regularization data can be conveniently conducted, and direct wave signals are picked up at the same time. The existence of direct waves can affect the structure of interbed multiples, so the direct waves in the three-dimensional prestack seismic data need to be cut off before the next step is carried out.

Step 102, time frame data of the primary reflected wave is produced.

The method of utilizing the virtual same-phase axis to predict the interbed multiples needs to construct virtual primary signals. Only one reflected wave signal is contained in the time frame, and the virtual in-phase axis signal can be accurately constructed. And meanwhile, the seismic data are regularized again according to the WRW model, namely, the common shot point seismic data or the common receiving point seismic data are all subjected to one-dimensional (arranged into a one-dimensional array), the regularized seismic data are arranged according to the abscissa, namely the common shot point seismic data, and the ordinate is the common geophone point data.

Returning to FIG. 1, step 103, three-dimensional prestack interbed multiples prediction, which directly predicts interbed multiples on three-dimensional shot gather data. This is the core technology of the present invention.

The invention discloses a method for predicting three-dimensional prestack interlayer multiples, which comprises the following two steps: the first is to construct virtual homophase axes; the second is predicting the interbed multiples. The method specifically comprises the following steps:

1) constructing virtual homophase axes and interbed multiples;

the virtual event and the interbed multiples constructed in the corresponding frequency domain using the virtual event method are respectively expressed as:

D_I(x_s，y_s，ω，x_r，y_r)＝∫∫D_V(x_s，y_s，ω，x，y)D′₀(x，y，ω，x_r，y_r) dxdy type (2)

Wherein D is_V(x_s，y_s，ω，x_r，y_r) Representing constructed virtual homophase axes, D_I(x_s，y_s，ω，x_r，y_r) Representing constructed interbed multiples, D₀(x，y，ω，x_r，y_r) And D'₀(x，y，ω，x_r，y_r) Respectively represent an upper part and a lower part of original data which are not overlapped with each other,represents D₀(x) is a conjugate of (C)_s，y_s) Representing source coordinates, (x)_r，y_r) Representing the coordinates of the receiving point.

2) Distinguishing data to realize accurate prediction of three-dimensional interbed multiples;

when data are distinguished, the upper part and the lower part cannot contain the same phase axis, otherwise, the inaccuracy of prediction is increased. At the same time, due to the 'lower' partial data volume D 'in the raw seismic data'₀Not only the primary reflection but also the interbed multiples are included, so the predicted result includes the interbed multiples. Seismic data used in the present invention refers to data in which surface multiples have been completely suppressed. To better describe the computational principle and physical meaning of the method, fig. 3 shows the computational process of the simple three-dimensional virtual event method, wherein fig. 3a shows the process of constructing a three-dimensional virtual event and fig. 3b shows the process of constructing an interbed multiple. The specific calculation method is as follows: first, using equation (1) to calculate virtual event data (i.e., as shown in fig. 3 a), the surface-excited and received seismic signals are converted into surface-excited and subsurface-received virtual event signals; next, the lower part data volume and the virtual event axis data are convolved by equation (2) to construct the interbed multiples of the corresponding horizon (i.e. as shown in fig. 3 b). The seismic data and the data are subjected to convolution in a time-space domain by using the formula (1) and the formula (2), and accurate prediction of three-dimensional interbed multiples can be realized. Meanwhile, according to the formulas (1) and (2), the method does not need any prior information such as an underground macroscopic velocity model and any assumption on an underground medium, can accurately predict the interbed multiples by using the pre-stack seismic data of the method, and is a pure data driving method.

Returning to FIG. 1, step 104, the three-dimensional interbed multiples are suppressed.

The predicted interbed multiples obtained in the last step have certain differences with interbed multiples corresponding to actual seismic data in amplitude, frequency domain and phase properties. The predicted three-dimensional interbed multiples need to be matched to the actual interbed multiples using a filter factor. As shown in fig. 1, multiple modules may be used to match the interbed multiples, including a single-pass norm matching module, a multi-pass norm matching module, and a square norm matching module. The method comprises the following specific steps:

in the invention, the three-dimensional interlayer multiple matching method adopts a first norm matching method and a second norm matching method. Firstly, a two-norm matching method is specifically described, and the energy is minimum after the input three-dimensional interlayer multiple data is subtracted by the predicted three-dimensional interlayer multiple data based on the principle of minimum energy. For the single pass two-norm method, the minimum objective function is:

the filter operator a is subjected to partial derivation on both sides of the formula (3) and is made equal to 0, and the formula (3) is minimized, and then converted into a linear equation:

M^TMa＝M^Td (4)

because the coefficient matrix on the left side of the formula (4) is a Topelitz matrix, the calculation is simple, the surface operator a can be solved through a Levinson recursive algorithm or SVD decomposition, and the division operation can be directly carried out under the condition that the equation is relatively stable. It should be noted that the length of the operator a is selected, and is generally stable when the length is short and sensitive when the length is long.

The norm matching method is specifically as follows:

e＝|d-Ma|₁ (5)

since the objective function is a singular function, equation (3) is transformed into a weighted two-norm minimization problem:

wherein, W represents a weighting matrix,ε＝max(d)/100；r_jd-Ma. Can be used forThe weight calculated by the first norm and the second norm is adjusted by changing the size of epsilon, so that the matching effect is improved. Further, a norm calculation equation of a single channel can be obtained as follows:

M^TW^TWMa＝M^TW^TWd (7)

and 105, outputting all seismic data after the three-dimensional interbed multiples are attenuated.

The iterative process between step 103 and step 104 in fig. 1 represents a process of distinguishing data according to an artificial interface, and according to actual seismic data, a main layer position generated by interbed multiples is judged, and a primary reflection wave signal of the layer position is determined, so that the number of artificial interfaces, namely the number of iterations in fig. 1, is determined. The process between them can be described in detail by fig. 4.

Fig. 4 provides a detailed flowchart of the iterative process based on step 103 to step 104 in fig. 1 according to the embodiment of the present invention. In step 401, regularization data designed according to the WRW model is input. In step 402, a human interface is selected based on the preprocessed timeline data, the three-dimensional data volume is divided into an 'up' portion and a 'down' portion, and the two data volumes are simultaneously fourier transformed into the frequency-space domain. In step 403, three-dimensional virtual event data is constructed using the virtual event method of FIG. 3. In step 404, the imaginary in-phase axis constructed in the previous step is convolved with the 'down' portion data to predict the interbed multiples of the corresponding horizon, and the seismic data is returned to the time-space domain using an inverse fourier transform. In step 405, the predicted three-dimensional interbed multiples are adaptively matched, a filter operator is solved, and the true interbed multiples in the shot gather are subtracted. And setting certain iteration times, and iterating the steps 403 to 405 to attenuate all order interbed multiples corresponding to the layer. For those skilled in the art, the iteration from step 403 to step 405 may be understood as a feedback iteration process of the surface multiples, and the surface multiples are continuously removed. Generally, the number of iterations set to 3-5 here can be considered to be well suppressed for all orders of interbed multiples corresponding to the horizon. In step 406, the position of the artificial interface is judged, if the position meets the condition, iteration is carried out (the artificial interface is moved downwards), the seismic data is divided into an 'upper' part and a 'lower' part again, and the processes from step 403 to step 406 are executed again; and if the condition is not met, outputting a calculation result.

Equations (1) and (2) are continuous integration processes and cannot be applied to discrete data signal processing. The actual seismic data is discrete, and equations (1) and (2) need to be discretized to facilitate processing of discrete seismic data signals. The discretization and iterative computation process involved in the iteration step of fig. 4 is as follows:

assuming that the surface multiples of the current seismic data have been removed and that there are now interbed multiples associated with the first layer interface, the iterative feedback prediction model for the interbed multiples can be expressed as:

wherein, P_m(x_r，y_r，z₀；x_s，y_s，z₀(ii) a ω) represents a shot point (x) from the surface_r，y_r) To the surface inspection point (x)_s，y_s) The multiple between the layers of (a) and (b),representing seismic data from which all interbed multiples waves associated with the first layer have been removed. Furthermore, S^-1(ω) represents the inverse matrix of the source wavelet, then S^-1 Representing the subsurface impulse response without interbed multiples. R (z)₁) The operator representing the downward scattering of the wave field at the interface of the first layer, generally taking the complex conjugate Δ P of the reflection of a pulse at the interface of this layer₁ ^*(z₁，z₁)S^-1(ω), then according to this definition, the downward scattering operator of the surface is a negative identity matrix. The physical meaning of formula (8) is: the interbed multiples associated with the first layer interface are equivalent to removing the layer primary reflectionThe wavefield pulse reflection of the wave, the layer down scatter operator, and the wavefield without the interbed multiples are convolved. However, the three wave fields of the right side of equation (8) are surface-excited subsurface reception, subsurface-excited subsurface reception, and subsurface-excited surface reception, which are not easily obtained in general practical seismic exploration, so that the equation (8) is rewritten such that all the wave fields of the right side of the equation are surface-excited and surface-received, i.e., the equation (8) is rewritten such that all the wave fields of the right side are surface-excited and surface-received

Note the book

The reflection of virtual pulse without multiple information can be obtained by substituting into interlayer multiple wave model

The imaginary pulse reflection without multiples information in equation (11) is the modified imaginary iso-axis, which is equivalent to the convolution of the primary pulse reflection associated with the layer interface and the wavefield pulse without the interbed multiples associated with the layer.

When self-adaptive virtual event and seismic data are convolved, primary waves are converted into interbed multiples, and n interbed multiples are converted into n +1 interbed multiples. The expansion according to the Noelman series can be obtained:

wherein:

equations (12) through (14) are the methods used in the iteration of fig. 4, and effective suppression of interbed multiples of three-dimensional prestack seismic data is achieved based on the method of the present invention.

Returning to fig. 1, all the steps in fig. 1 are executed, and all the prestack seismic data after three-dimensional interbed multiple wave waves are suppressed can be obtained. For further details of the present embodiment, please refer to the above description.

Application example

This example uses a three-dimensional 4-layer layered model with a size of 1km x 1km and a grid spacing of 5 m. The speed of the model is 2km/s, 1.5km/s, 2.5km/s and 1.5km/s from top to bottom in sequence. Both the shot spacing and the geophone point spacing were set to 0.04 km. The specific speed configuration is shown in fig. 5. And carrying out three-dimensional numerical simulation by using a finite difference method. FIG. 6 shows shot gather data (after the direct wave is excised) received from a seismic source located at the center of the earth's surface on a model. FIG. 7 is the re-regularized shot gather data according to the WRW model. By observing the single shot data of fig. 6 and 7, the development of interbed multiples is known.

By applying the method, a virtual event model can be constructed as shown in fig. 8, and by comparing the virtual event model with the real shot gather data shown in fig. 7, the virtual event model can be known to be unobservable in actual seismic records and only predict the intermediate process of interbed multiples. FIG. 9 is a three-dimensional interbed multiples predicted using a virtual event. By comparing with fig. 7, it can be seen that the predicted interbed multiples are accurate. Fig. 10 shows real interbed multiples matched by filtering using an adaptive matching method. FIG. 11 is pre-stack seismic data after flattening a plurality of wavefronts between three-dimensional layers. It can be seen that the interbed multiples have been significantly attenuated, but some residue remains. FIG. 12 is a seismic data utilizing an iterative method to suppress multiple wave-fronts between layers. Numerical results indicate that interbed multiples are effectively suppressed.

It is noted that the disclosed embodiments are intended to aid in further understanding of the invention, but those skilled in the art will appreciate that: various substitutions and modifications are possible without departing from the spirit and scope of the invention and appended claims. Therefore, the invention should not be limited to the embodiments disclosed, but the scope of the invention is defined by the appended claims.