CN112946732B

CN112946732B - Processing method and system for jointly suppressing multiple single cable of offshore stereo observation system

Info

Publication number: CN112946732B
Application number: CN202110140513.2A
Authority: CN
Inventors: 张进; 董博艺; 邢磊; 彭阳阳; 高俊杰; 王林飞; 郭绪兵; 尹燕欣
Original assignee: Ocean University of China
Current assignee: Ocean University of China
Priority date: 2021-02-02
Filing date: 2021-02-02
Publication date: 2022-09-30
Anticipated expiration: 2041-02-02
Also published as: CN112946732A

Abstract

The invention belongs to the technical field of seismic wave exploration, and discloses a processing method and a processing system for jointly pressing a multiple single cable by a marine stereoscopic observation system, wherein a seismic record received by a horizontal cable is extended to the seabed in the forward direction of the propagation direction of seismic waves; performing wave field separation on the vertical cable seismic record, wherein the separated downlink wave is just a multiple reflected wave received after being reflected by the sea surface; carrying out forward continuation on the separated downlink wave along different paths based on the positions of the detectors with different depths in an F-K domain to obtain a predicted multiple wave field; the horizontal cable wave field is matched with the vertical cable multiple wave field based on a single-channel least square matched filtering method, so that multiple information suppression in the horizontal cable wave field is realized, multiple is removed, and effective wave field information is reserved. The method for jointly pressing multiple waves can achieve good effects, and fully shows the adaptability and the effectiveness of the method in the field of marine multiple wave pressing treatment.

Description

Processing method and system for jointly suppressing multiple single cable of offshore stereo observation system

Technical Field

The invention belongs to the technical field of seismic wave exploration, and particularly relates to a processing method and a processing system for jointly suppressing multiple single cables by using an offshore three-dimensional observation system.

Background

Currently, the closest prior art: in the traditional seismic wave exploration method, a multiple is one of the most common interference noises, particularly in marine seismic exploration, due to the existence of stronger multiple, the signal-to-noise ratio of seismic data is reduced, which seriously interferes with the identification of effective waves by people, and makes velocity analysis, prestack and poststack migration in subsequent work very difficult, further influences the authenticity and reliability of seismic wave imaging, and leads to false imaging and serious influence on seismic interpretation work. Therefore, how to adopt a better processing means to explore the causative characteristics of the multiples and inhibit the attenuation multiple interference always has very important practical significance and development prospect.

At present, the multiple attenuation method can be mainly divided into two categories, one is to predict and attenuate multiple waves according to separability and periodicity between the effective waves and the multiple waves, time difference and the like. The second method is a wave equation prediction subtraction method, which is mainly based on the elastic wave equation theory, and is used for predicting multiples in original data by a forward simulation or inversion method, and then matching attenuation is performed on the predicted multiples from the original data. The amplitude fidelity of the seismic record is excellent by the prediction subtraction method, the algorithm is independent of the speed, an underground medium model is not needed or rarely needed, and the method is applied to the attenuation analysis of different kinds of multi-wave and has wide application prospect.

In the field of marine geophysical prospecting, free-surface multiple interference is a coherent interference that is particularly pronounced in seismic data. Due to the existence of multiple waves, the signal-to-noise ratio of seismic data is reduced, and the recognition rate of effective waves in seismic records is reduced. The method has great influence on subsequent velocity analysis and pre-stack and post-stack migration work, and causes great interference on the authenticity of seismic imaging. How to better suppress multiples on the basis of preserving the energy of the primary waves is very important for a subsequent series of work.

In recent years, the technology of multiple pressing treatment has been developed and is receiving more and more attention. Discussion and exploration of the multiple problem are always the prominent focus of the conference at the annual SEG annual meeting. With the development of science and technology, people pay more and more attention to the problem of multiple waves in seismic data. Since the 50 s of the 20 th century, more and more multiple attenuation processing methods have been proposed through continuous analysis and development. The predictive deconvolution method belongs to a processing means which is provided relatively early in time, and generally speaking, the method has a better processing result for 1D geological conditions. When the local underground medium has a larger positive velocity gradient, a parabolic Radon transformation method for transforming a parabolic in-phase axis into a lambda-F domain to form a straight line with a certain slope and filtering the straight line can be adopted, a hyperbolic-type seismic in-phase axis is subjected to Radon transformation to be used for hyperbolic-type Radon transformation of multiple attenuation, and a beam-focusing filtering method for separating multiples and primary waves in an F-X domain to eliminate the multiples is adopted to perform multiple attenuation suppression. These are all filtering methods based on time difference, and the effect is also obvious. However, if the time difference between the effective reflected wave and the multiple is small, especially when the offset distance is small, such a filtering method often causes a certain adverse effect on the effective wave; in addition, for the complex underground medium, the seismic response of the effective wavefield no longer strictly follows hyperbolic or parabolic laws. Due to the complex terrain conditions, the received effective reflected wavefield may be distorted. These factors all impose certain limitations on the use of filtering methods. Therefore, the conventional prediction deconvolution method and various filtering methods are only effective for some simple situations, such as the problem of multiples that the distinguishing features of primary waves and multiples are obvious, the periodicity is good, and the separability is strong.

In the 80 s of the 20 th century, the concept of free surface multiples and internal multiples was proposed in the field of multiple analysis, and further analysis on a multiple generation mechanism found that free surface multiples derived from back-and-forth reflections between a free surface and a reflection interface are the main multiple energy in seismic data. Later, wavefield extrapolation methods based on wavefield continuation theory appeared. The method is proposed for use in the multiple processing problem in marine seismic data. Loewenthal first predicted multiples using wavefield continuation methods and applied them to Common Depth Point (CDP) gathers. Kennett proposes a simulation method of free surface multiples in a one-dimensional space and an inversion method thereof, and the method needs prior estimation of underground information; riley and claerbb propose a 2D space free surface multiple simulation algorithm; berkhout proposes a multidimensional inversion algorithm for eliminating free surface multiples, which can process the situation that underground media are more complex; based on the minimum energy criterion, Verschuur proposes a method to suppress free-interface-related multiples, and Berkhout further gives predictions based on the wave equation and theory of free-surface multiple attenuation. In the process of predicting the multiples, the macroscopic velocity field does not need to be known in advance, so the adaptability of the prediction method is improved. Aiming at the condition of near migration data loss in seismic data, the multiple processing method given by Wang is a development based on a fluctuation theory method; morley deeply analyzes the derivation process of a wave equation operator for predicting the submarine multiples, and simplifies the prediction process of the multiples; wiggins perfected and further developed the method of Morley; the method of Berryhill is similar in principle to the Wiggins method, but it is implemented differently. The wave equation method is a trend in the development process of multiple wave processing technology by virtue of good application results.

The wave field continuation method is a prediction method based on wave field continuation theory, when wave field information received by a sea surface wave detector is transmitted downwards in a sea water layer again for a travel twice the depth of sea water, effective reflected wave field (primary wave) received for the first time can be converted into double-pass reflected wave, namely first-order multiple, the multiple series in the sea water layer is increased, and the received multiple wave field information is removed from the original wave field, so that the effect of suppressing multiple can be realized. Such a method has an advantage in that damage to the effective wave can be accurately minimized while attenuating multiples. However, the multiple wavefields predicted by the wavefield extrapolation differ from the multiple wavefields actually recorded in terms of their amplitude and phase and their arrival times. In this case, if the actual data is directly subtracted from the predicted multi-wave model, there is a high possibility that the primary reflected wave will be damaged. Therefore, a suitable matching algorithm needs to be selected for performing matched filtering processing in the later stage.

In recent decades, among various wave field continuation methods based on the wave equation of sound waves, the wave equation extrapolation method of a frequency-wave number domain (F-K domain) has the advantages of high calculation speed, simple form, small data storage capacity and the like. The wave field continuation method of the frequency wavenumber domain (F-K domain) proposed by Stolt under the premise of constant medium speed adds fast Fourier transform to the continuation process of the constant-speed frequency wavenumber domain. This is much faster than other methods. In addition, the method utilizes wave equation solution, so that the method is suitable for seismic wave fields under the condition of large inclination angle, the dispersion problem is small, and the calculation accuracy is high. Especially under the condition of the same time sampling interval, the calculation result after the wave field continuation is still very stable. However, the extended method of constant velocity media is only applicable to seismic wavefields where the seismic wave velocity is constant. Gazdag, 1978, therefore proposed phase shift wave equation extrapolation, achieved with zero offset, with the property of enabling precise shift excursions in the vertical direction. Later, Gazdag et al provided a phase shift interpolation method, which can handle the change in velocity of longitudinal and transverse seismic waves of the subsurface medium, improve the imaging quality of the seismic wavefield, and achieve the homing of the seismic wavefield at large dip angles. Lai and G.H.F.Gadner. then applied the phase shift interpolation method to seismic model forward modeling and succeeded. In 1985, Claerbout proposed a phase shift and phase shift interpolation wave field continuation method based on a common midpoint gather.

The phase shift method wave field extrapolation is realized by a synthetic seismic record method, and compared with the F-K method of Stolt, the method can obtain a more accurate wave field extrapolation result under the condition of wave equation vertical speed change. However, this method also has certain problems: in the process of downwardly extending a wave field, the event of the seismic wave event can generate a turn-back effect, and the existence of the turn-back effect can cause severe extension noise and interference on the event of the seismic migration profile, so that the signal-to-noise ratio of the seismic profile is reduced, and the imaging quality is poor. Although the boundary reflection effect can be eliminated by a method of zero filling on two sides of the seismic record, the data volume processed is doubled, the calculation efficiency is greatly reduced, and the memory occupation of a computer is large. In 1984-1985, shizhenhua and g.h.f.gardner adopt a method of space domain boundary absorption to effectively solve the continuation noise caused by the foldback effect, but because the condition of the absorption boundary is set to be realized in the space domain, and the realization process of wave field continuation is completed in the frequency and wave number domain, therefore, Fourier forward transformation and Fourier inverse transformation must be respectively performed along the space direction in the process of each continuation step, and when the continuation steps are more, the calculation speed is still slow.

Many scholars in China also do a great deal of analysis work on the problems of multiple prediction and attenuation, develop various different multiple processing methods and obtain certain results. The Zhang jin intensity adopts a self-adaptive method to attenuate free surface multiples, and the effect is obvious; the free surface multiple theory and the pressing method are comprehensively researched by Shenzhao and Nipponhua, forward simulation of the free surface multiple of a wave equation method and a pressed series sequence formula and an iterative formula are deduced, test processing is carried out on a theoretical model, the correctness of the formula is verified, and a certain application effect is obtained in practical application; the Tanshosha spring and the like apply a free surface multiple wave pressing method based on a wave equation in practical data and obtain better results; the method has the advantages that 1, due to systematic analysis, on the basis of data consistency prediction and theory, implementation process and theoretical model trial processing for suppressing free surface multiples, the effect of suppressing the multiples related to the free surface is very obvious; the method of combining the parabolic transformation and the wave equation, such as the Wangweihong method, suppresses the multiple waves, is based on Radon transformation and the wave equation theory, and has ideal effect; when the Zhang Xing Fang is used for carrying out free surface multiple attenuation, the Curvelet transformation and the feedback iteration method are combined to carry out the multiple removal work together, and the method has good self-adaptability on the multi-order multiple pressing the free surface. At present, the multiple suppression attenuation technology based on the wave equation principle is not developed deeply, and in general, the analysis of the multiple suppression technology and method cannot be mature.

In summary, the problems of the prior art are as follows:

at present, the multiple attenuation methods can be mainly divided into two categories, one category is to predict and attenuate multiple according to the separability and periodicity between the effective wave and the multiple, time difference and other characteristics, for example, f-k domain filtering, Radon filtering, tau-p domain transformation and other methods are used for suppressing multiple, and the method utilizes the difference between multiple and one time in time, so that the method is suitable for seismic data with obvious time difference between multiple and one time, but the effect of suppressing multiple of small offset data is not obvious, and effective signals are possibly damaged. The second method is a wave equation prediction subtraction method, which mainly predicts multiples in original data by a forward modeling method based on an elastic wave equation theory, and then performs matched attenuation on the predicted multiples from the original data, such as an SRME method. The prediction subtraction method is excellent in amplitude fidelity to seismic records, but has the problems of wave field phase distortion after multiple suppression, primary wave amplitude damage and the like.

In the two methods, the data of the towing cable is adopted for multiple suppression, and compared with the towing cable, the vertical cable can acquire abundant up-going waves and down-going waves, the multiple separation is simpler, and the data quality is higher. Therefore, a method for performing multiple combined suppression by using a stereo observation system formed by the marine streamer and the vertical cable is provided.

The difficulty of solving the technical problems is as follows:

how to separate multiple wave data from the vertical cable, how to extend the multiple wave data of the vertical cable to the horizontal towing cable, and performing matched attenuation to achieve the purpose of suppressing the multiple waves.

The significance of solving the technical problems is as follows:

the existing method only limits the suppression of the multiple waves to the data of the towing cable, combines the horizontal towing cable and the vertical cable, can fully utilize the advantage that the multiple waves on the vertical cable are easy to separate and extend, and improves the suppression effect of the multiple waves of the horizontal cable.

Disclosure of Invention

Aiming at the problems in the prior art, the invention provides a processing method and a processing system for a combined multiple pressing single cable of an offshore stereo observation system.

The invention is realized in this way, a processing method of a marine stereoscopic observation system combined pressed multiple single cable, which comprises the following steps:

firstly, extending the seismic record received by a horizontal cable to the seabed in the forward direction of the propagation direction of seismic waves;

secondly, performing wave field separation on the vertical cable seismic record, wherein the separated down-going wave is a multiple reflected wave received after being reflected by the sea surface;

thirdly, wave field forward continuation is carried out on the separated downlink waves in an F-K domain along different paths based on the positions of the detectors with different depths, the wave fields are extended to a seabed datum plane in a unified mode to be in one-to-one correspondence with horizontal cable wave fields, and a predicted multiple wave field is obtained;

and fourthly, matching the horizontal cable wave field with the vertical cable multiple wave field based on a single-channel least square matched filtering method, suppressing multiple information in the horizontal cable wave field, removing the multiple, reserving effective wave field information to obtain a result of removing the multiple, and performing reverse continuation on the obtained result of removing the multiple to return to the sea surface.

Further, after the fourth step, the following steps are also performed: verifying the application effect of reversely extending to the sea surface through different data models, which specifically comprises the following steps:

(1) firstly, predicting a multi-time wave model by using vertical cable data, and performing data matching correction on a vertical cable and a horizontal cable;

(2) then, a single-channel matching adaptive subtraction method is used for suppressing and attenuating the multiple waves, and the advantages and the disadvantages of the processing effects under different models are compared;

the different data models include:

the horizontal layered model adopts relatively ideal horizontal layered stratum media, the uppermost layer is a seawater layer, the lower layer is a layered horizontal stratum, five layers of media are counted, and the sea level is a good free reflection interface;

a stratigraphic dip marine geological model for verifying the reliability of the combined compaction process when lateral shifting is not large.

The invention provides a method for jointly suppressing multiple waves by a stereo observation system, which combines the marine vertical cable technology with the streamer technology, and jointly performs the prediction and the joint suppression process of the marine multiple waves by utilizing the characteristics that a marine vertical cable wave detector is placed along with the depth, can simultaneously receive an upper traveling wave and a lower traveling wave and is easy to separate a seismic wave field.

Further, in the third step, the F-K domain wave field separation method specifically comprises the following steps:

firstly, originally inputting data information by a vertical cable, carrying out F-K conversion on the original vertical cable seismic data, displaying the original vertical cable seismic data as wave field energy information of different interval positions in a frequency-wave number domain, reserving an upgoing wave field, removing the wave field information of a downgoing wave in the frequency-wave number domain, carrying out inverse F-K conversion on the remaining wave field information, returning to a t-x domain, and then obtaining a reflected wave field record which is separated out a downgoing wave and only reserves an upgoing wave;

and (4) performing wave field continuation on the downlink wave on the vertical cable to predict multiple waves, and performing matched pressing.

Further, the method for performing the wave field continuation comprises the following steps:

for the calculation form exp [ -i ω (t ± r/v) ], r in the formula represents the relative length of the seismic source position and the observation position, a one-way wave equation is utilized to carry out wave field continuation, only a wave field advancing towards a single direction is adopted, an uplink one-way wave is propagated upwards, a downlink one-way wave is propagated downwards, a negative sign in the formula represents a downlink one-way wave, and a positive sign in the formula represents an uplink one-way wave; the three-dimensional wave equation is as follows:

carrying out Fourier transformation to obtain:

and during wave field continuation, performing wave field extrapolation on the uplink wave and the downlink wave, and performing wave field extrapolation along the forward propagation mode of the actual seismic wave field from the seismic source position by adopting a forward continuation mode to realize a forward simulation process.

Further, the Fourier transform method specifically includes:

p(k _x ，k _z ，ω)＝∫∫∫p(x，z，t)exp(ik _x x+ik _z -iωt)dxdzdt；

p(x，z，t)＝∫∫∫p(k _x ，k _z ，ω)exp(-ik _x x-ik _z +iωt)dk _x dk _z dω；

solving the wave equation in the form of Fourier domain:

written in the form of

The upper type

Representing the wave dispersion form of a one-way wave equation, taking positive on the right of the equation as a downlink wave, and taking negative on the right of the equation as an uplink wave;

the corresponding expressions of the wave equation t-x domain and the f-k domain are as follows:

the wave dispersion form is changed into a one-way wave equation:

the method for extrapolating the phase-shifted wave field in the F-K domain comprises the following steps:

(1) first, the surface wavefield p (k) _x Z-0, ω) to the frequency-wavenumber domain p (k) _x ，z＝0，ω)；

(2) Calculating an extrapolation factor

(3) In the frequency-wavenumber domain, multiplying a seismic wave field by a factor C to obtain a wave field at a delta z position;

(4) repeating the steps (2) and (3) to gradually obtain frequency domain wave fields of all depths i delta z (i is 1,2,3 …) underground;

(5) wave field p (k) of frequency-wavenumber domain _x The z is i Δ z, ω) is inversely transformed back to the time-space domain p (x, z is i Δ z, t), and finally, the time-space domain seismic wave fields at different depths are obtained;

the method for extrapolating the F-K domain phase-shifted wave field comprises the following steps:

further, the method for F-K domain phase shift wavefield extrapolation further comprises: firstly, separating an upgoing wave and a downgoing wave from a total wave field of a vertical cable, and then respectively carrying out variable-depth wave field extrapolation of different paths and different distances on wave fields received by detectors at different depths based on a wave field continuation operator to obtain an extrapolated seismic wave field value;

and matching the multi-wave data obtained by extrapolation with the streamer data.

Further, in the fourth step, a horizontal cable wave field is matched with a vertical cable multiple wave field based on a single-channel least square matched filtering method, and multiple information suppression in the horizontal cable wave field is performed, wherein the least square matched filtering method comprises the following steps:

in the prediction process, firstly, vertical cable data predict multiple waves on a horizontal cable through wave field continuation; after the multiple waves are predicted, self-adaptive subtraction is carried out, and the multiple waves predicted to come out of the vertical cable are subtracted from the horizontal cable data;

converting the above formula into a time domain, and simultaneously limiting the length of a filter factor A in the time domain; then

The above equation is transformed into:

the time domain filtering factor a (t) in the formula is obtained, and the effective wave data P after the multiple wave is eliminated is obtained ₀ The energy of (2) is minimized as much as possible.

The invention also aims to provide the offshore stereo observation system combined pressing multiple single cable processing system of the offshore stereo observation system combined pressing multiple single cable processing method.

The invention also aims to provide an information data processing terminal for realizing the processing method of the combined multiple pressing single cable of the marine stereoscopic observation system.

It is another object of the present invention to provide a computer-readable storage medium, comprising instructions which, when executed on a computer, cause the computer to perform the method for processing a single cable of a combined multiple press for a marine stereoscopic vision system.

In summary, the advantages and positive effects of the invention are:

Because the horizontal cable and the vertical cable wave detector are arranged in different modes, the seismic wave fields of the horizontal cable and the vertical cable wave detector need to be calibrated to a uniform reference surface in a wave field continuation mode, and then matching pressing is carried out. The selection of the correction reference surface can select any plane in the water layer, so that the seabed surface is selected as the extension reference surface for convenience in processing and reducing extension noise as much as possible, and the matching pressing multiple process is uniformly calibrated to the seabed surface. Firstly, establishing an offshore stereo observation system for forward modeling, and respectively obtaining wave field information on an offshore towing cable and a vertical cable; then according to the relation between the total wave field on the vertical cable and the up-going wave and the down-going wave, adopting a proper filtering method to carry out wave field separation to obtain a separated multiple wave field of the vertical cable; respectively carrying out wave field extrapolation on seismic wave fields received at different depths on the vertical cable based on a wave field continuation theory so as to obtain wave field information of the multiple waves of the predicted towing cable; and finally, performing combined matching on the multiple waves predicted by the marine vertical cable and the streamer observation data based on a single-channel minimum quadratic matching filtering method, and finally realizing the process of jointly suppressing the multiple waves. The forward modeling data processing results of different geological models show that the method can achieve good multiple suppression attenuation effects, and the feasibility and effectiveness of the combined multiple suppression method under the marine three-dimensional observation system are fully demonstrated.

According to the method, aiming at the difference of the observation modes of the towing cable and the vertical cable, the horizontal correction is carried out on the data of the vertical cable based on the F-K domain phase shift wave field continuation method, the matching filtering is carried out on the data of the vertical cable and the horizontal cable, the calculation mode is simple, the speed is high, and the method has a good application effect.

Aiming at the problem of noise generated in the wave field continuation process, an F-K domain wave field absorption filtering method is adopted for noise suppression, so that the foldback effect can be better improved, and the processing and imaging quality is improved.

Drawings

Fig. 1 is a flowchart of a processing method for jointly suppressing multiple single cables by using a marine stereo observation system according to an embodiment of the present invention.

Fig. 2 is a full-scale multiple chart provided by the embodiment of the invention.

Fig. 3 is a diagram of short-range multiples provided by an embodiment of the present invention.

FIG. 4 is a diagram of a pegleg diagram according to an embodiment of the present invention.

Fig. 5 is a ghost map provided by an embodiment of the present invention.

FIG. 6 is a schematic diagram of the spatial location of acoustic wave field components provided by an embodiment of the present invention.

Fig. 7 is a diagram of a seawater-seabed two-layer forward modeling provided by the embodiment of the invention.

FIG. 8 is a forward modeling diagram of a streamer multiple provided by an embodiment of the invention.

Fig. 9 is a diagram of a marine vertical cable observation system provided by an embodiment of the invention.

FIG. 10 is a graph of the type of VSP seismic waves provided by an embodiment of the present invention.

Fig. 11 is a diagram of a vertical cable forward modeling seismic recording according to an embodiment of the present invention.

FIG. 12 is a graph of F-K domain spectral analysis of a vertical cable wavefield according to an embodiment of the present invention.

FIG. 13 is a diagram of a F-K domain wavefield separation method according to an embodiment of the present invention.

Fig. 14 is a diagram illustrating the effect of separating the up-going wave and the down-going wave by the F-K conversion method according to the embodiment of the present invention.

FIG. 15 is a schematic diagram of median filtered vertical cable wavefield separation according to an embodiment of the present invention.

Fig. 16 is a diagram of the effect of separating the up-going wave and the down-going wave by the median filtering method according to the embodiment of the present invention.

Fig. 17 is a graph illustrating the effect of the F-K domain filtering method compared with the median filtering method according to the embodiment of the present invention.

FIG. 18 is a diagram of a marine VSP reception mode provided by an embodiment of the present invention.

FIG. 19 is a flow chart of phase-shift wavefield extrapolation, according to an embodiment of the present invention.

FIG. 20 is a graph illustrating the foldback suppression effect of the F-K domain filtering method according to the embodiment of the present invention.

Fig. 21 is a schematic diagram of continuation correction of a stereo observation system according to an embodiment of the present invention.

Fig. 22 is a diagram of a two-cable continuation matching process provided by an embodiment of the invention.

Fig. 23 is a comparison graph of a two-wire continuation waveform provided by an embodiment of the present invention.

FIG. 24 is a comparison of streamer and vertical cable propagation waveshapes provided by embodiments of the present invention.

FIG. 25 is a graph of the primary wavefronts left after subtraction using matching as provided by embodiments of the present invention.

FIG. 26 is a diagram of predicted multiples provided by an embodiment of the present invention.

FIG. 27 is a diagram of a horizontal layer model provided by an embodiment of the present invention.

FIG. 28 is a diagram of streamer recording (a) and vertical cable seismic recording (b) as provided by an embodiment of the invention.

FIG. 29 is a diagram of wave field continuation predicted multiples provided by an embodiment of the present invention.

Fig. 30 is a least squares matched filter graph provided by an embodiment of the present invention.

FIG. 31 is a diagram of a model of a dipping formation provided by an embodiment of the invention.

FIG. 32 is a diagram of streamer recording (a) and vertical cable seismic recording (b) as provided by an embodiment of the invention.

FIG. 33 is a graph showing the combined multiple-suppressing results of the model of the inclined sea floor according to the embodiment of the present invention.

FIG. 34 is a rugged subsea model diagram provided by embodiments of the present invention.

FIG. 35 is a graph of the combined rough seafloor model hull multiples results provided by an embodiment of the present invention.

Detailed Description

In order to make the objects, technical solutions and advantages of the present invention more apparent, the present invention is further described in detail with reference to the following embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not intended to limit the invention.

In the prior art, the data of the towline is adopted to carry out multiple suppression, the information of the types of the collected upgoing wave, the downgoing wave and the seismic wave is not rich, and the quality of the collected seismic wave data is poor.

1. The present invention will be described in detail below with reference to the accompanying drawings.

The system analyzes the forming mechanism and characteristics of the sea multi-wave, adopts a sea stereo observation system, obtains seismic data of a horizontal cable and a vertical cable by using forward simulation, realizes the extraction of the vertical cable multi-wave by using a wave field separation method, predicts the multi-wave by using a wave field continuation method, performs combined matching and pressing with the horizontal cable record by using a self-adaptive subtraction matching method, and finally obtains the seismic record processing effect of pressing the multi-wave. The processing and inspection results of different model data fully prove the practicability of the method, and have certain theoretical significance and application value.

In the invention, wave field characteristic analysis under a stereo observation system is carried out. And (3) establishing an offshore vertical cable and horizontal cable combined observation system, performing forward modeling on different seabed geological models, respectively obtaining seismic records containing effect waves and multiple waves in wave fields of the vertical cables and the horizontal cables, and analyzing the wave fields.

As shown in fig. 1, a method for processing a single cable by combining multiple suppression for a marine stereo observation system provided by an embodiment of the present invention includes:

a processing method for a single cable with multiple waves pressed by a marine stereoscopic observation system in a combined mode comprises the following steps:

in the first step, the seismic records received by the horizontal cable are extended to the seabed in the forward direction of the propagation direction of the seismic waves.

And secondly, performing wave field separation on the vertical cable seismic record, wherein the separated down-going wave is a multiple reflected wave received after being reflected by the sea surface.

And thirdly, performing wave field forward continuation on the separated downlink wave in an F-K domain along different paths based on the positions of the detectors with different depths, uniformly extending to a seabed reference surface to realize one-to-one correspondence with the horizontal cable wave field, and obtaining a predicted multiple wave field.

In the embodiment of the invention, wave field characteristic analysis under a stereo observation system is carried out. And (3) establishing an offshore vertical cable and horizontal cable combined observation system, performing forward modeling on different seabed geological models, respectively obtaining seismic records containing effect waves and multiple waves in wave fields of the vertical cables and the horizontal cables, and analyzing the wave fields.

And (5) performing offshore vertical cable wave field separation analysis. And (3) performing up-down wave separation on the vertical cable seismic wave field based on wave field separation methods such as an F-K domain filtering method, a median filtering method, a linear Radon transformation method and the like.

And (5) analyzing a double-cable wave field continuation matching method. Based on the F-K domain wave field continuation method, multiple prediction is carried out by utilizing vertical cable data, and combined matching correction is carried out with horizontal cable data, so that the process of jointly suppressing multiple by double cables is realized.

In the embodiment of the invention, aiming at the difference of the observation modes of the towing cable and the vertical cable, the horizontal correction is carried out on the data of the vertical cable based on the F-K domain phase shift wave field continuation method, and the matched filtering is carried out on the data of the vertical cable and the horizontal cable, so that the calculation mode is simple, the speed is high, and the application effect is better.

For the problem of insufficient illumination of a three-dimensional observation system under a complex seabed model, the observation system is changed based on different geological conditions, a multi-vertical cable system is introduced to carry out data acquisition and processing together, corresponding multiple matching suppression is carried out, and a better processing effect can be obtained compared with a simple observation system.

2. The present invention is further described below with reference to the wave field characteristics and processing method of the multiple waves of the marine stereo observation system.

2.1 Marine stereoscopic observation system

Nowadays, global marine oil and gas exploration faces the challenges of more complex structures such as deep water, deep layers, salt deposit and the like and lithologic stratum trap, and has higher requirements on the aspects of imaging precision, resolution, deep layer reflection strength and the like of marine seismic exploration data, and the requirements of the more complex geological problems cannot be met by the conventional marine streamer seismic exploration method. The industry always carries out technical research and development and innovation in aspects of equipment technology, data acquisition mode, navigation positioning, construction quality control, data processing method and the like, and the development of the current marine seismic exploration technology is supported by the continuous integration of science and technology.

In order to meet the requirement of accurate imaging of seismic data in a new exploration form, international marine petroleum seismic exploration has emerged a plurality of new data acquisition and processing technologies in recent years. In the excitation mode and the receiving mode of the seismic waves, the combined three-dimensional joint excitation and receiving of a multi-stage seismic source and a multi-cable are changed from a plane seismic source and a horizontal cable; in marine three-dimensional seismic exploration, narrow azimuth data acquisition of a straight line air route is developed to multi-ship multi-source multi-cable wide azimuth data acquisition and annular omnibearing data acquisition; three-dimensional seismic imaging also evolves from large to small planar elements. The new seismic acquisition technical method is applied, effectively overcomes the defects of seismic exploration of marine streamer cables, enriches azimuth information, improves the wave field illumination of steep dip stratum, complex fault and complex fault block, and also improves the signal-to-noise ratio of deep effective reflection signals and the imaging effect and precision, thereby meeting the requirements of exploration and development of complex oil and gas fields.

Marine seismic exploration is most commonly used with streamer operation. The streamers are sometimes laid on the sea surface and sometimes submerged below the sea surface. However, conventional streamer seismic methods at sea have serious drawbacks, as the trapping caused by the ocean multiples can severely affect the spectral content of the seismic data.

With the development of marine high-precision exploration, marine detection equipment is updated day by day, and exploration technologies including digital vertical cables, OBS, OBC and the like are gradually developed from original horizontal towing cables. Although the seismic exploration means of the conventional method has mature performance in suppressing multiple waves according to the data acquired by the horizontal towrope, the offshore vertical cable is used as a new data acquisition method for marine oil and gas exploration, a new visual field is provided for data acquisition, and the reception of uplink waves and downlink waves by the offshore vertical cable is different from that of the horizontal towrope, so that a new idea is provided for suppressing the multiple waves.

The marine stereo observation system combines the advantages of horizontal towing cable and vertical cable observation modes and performs acquisition and recording work of marine exploration seismic wave fields together. The marine survey vessel is provided with a proper observation system, seismic waves excited by a seismic source are diffused in a three-dimensional space formed by combining a horizontal cable and a vertical cable, wave field information is received by wave detectors positioned in different positions and different height spaces, and the characteristics of the seismic waves positioned at different wave detection positions are recorded. Due to the respective acquisition characteristics and complementarity of the different receiving cables, the effective wavefield reflected back through the sea bottom reflection interface is simultaneously received by the horizontal streamer and the vertical cable, and the multiples reflected back to the sea bottom through the sea surface free surface are also received by the horizontal streamer and the vertical cable. Due to the difference of the arrangement position and the receiving characteristics of the geophones, the wave field information received by the horizontal streamer on the seismic record is different from and related to the wave field information received by the vertical streamer. How to carry out reasonable wave field information joint correction and matching on wave field characteristics obtained by different receiving modes, and carry out recognition suppression and processing work on marine multiples together to realize joint acquisition processing of a stereo observation system is an important analysis content of the invention. The characteristics of the marine seismic wavefield and the wavefield continuation theory are explained in detail.

2.2 sea multiple wave field characteristics

2.2.1 Generation of multiple waves

In order to apply systematic analysis, discrimination and attenuation means to the multiples, firstly, the generation mechanism and self characteristics of the multiples need to be comprehensively analyzed, and the differences and the connections between the generation mechanism and the effective reflection are explored.

Multiple wave generation is caused by the repeated reflection of seismic waves between two or more interfaces in the subsurface. When the reflected wave propagates to the sea surface (ground), the interface between the sea surface and the air is an interface with a significant wave impedance difference. The reflected wave may reflect downward from the interface, and when propagating to the seafloor reflection interface, the reflection again propagates back to the surface, thus forming multiples.

The interface with better reflection coefficient can generate multiple waves, the seismic wave field can generate reflection energy when encountering the reflection interface in the downward transmission process, but if the reflection coefficient of the reflection interface is small, the wave field energy after being reflected by the interface is very small, if multiple reflections are continuously carried out between stratums, only the multiple waves with extremely small energy and extremely unobvious energy can be generated, and only multiple reflected waves with strong energy generated on the reflection interface with larger reflection coefficient can be obviously identified on the seismic wave receiving record, so that the interface can generate stronger and obvious multiple waves compared with strong reflection interfaces such as sea water surface, seabed, bedrock surface, limestone, strong unconformity surface and the like.

Type of 2.2.2 multiples

Generally speaking, the reflection coefficient of sea level is-1, the wave field will generate stronger reflection when propagating to sea level, similarly, when sea bottom also has stronger reflection coefficient, after seismic wave is excited by seismic source between sea level and sea bottom surface, the wave field will propagate between sea water layer many times, and when meeting reflection interface, it will generate multiple reflection wave. When the depth of the seawater is small, multiple waves generated in the seawater layer are called as ringing; when the depth of the seawater is large, the multiple waves generated in the seawater layer are called as free surface multiple waves. And based on the wave field propagation characteristics of multiples, the multiples can be classified mainly as follows:

(1) the whole process is carried out for multiple times: the whole-course multiple refers to a reflected wave generated at a deep reflection interface, which is reflected after propagating back to the sea surface (ground surface), and then propagates downward to the same interface to be reflected again, and is reflected back and forth between the two reflection interfaces for a plurality of times, which is called a whole-course multiple, and may also be called a simple multiple, as shown in fig. 2. According to the geometrical principle of seismic wave propagation, the multiple time-distance curve equation of the horizontal seabed stratum can be expressed as

In the formula: t is t ₀ 2 h/v; h is the depth of the reflection interface seabed; v is the propagation speed of seismic waves in seawater; x is the distance from the demodulator probe to the shot point; n represents the nth reflected wave of the sea bottom, when n>1, it is called a multiple reflection wave.

(2) Short-range multiples: sometimes, the seismic wave propagates to a strong reflecting interface in the subsurface and then reflects from the deep interface back to the sea surface, where it again propagates down the sea surface and then reflects at a shallower interface, called a short-range multiple, also called a local multiple, as shown in fig. 3.

(3) Micro-bending multiples: multiple reflections are carried out between certain reflection interfaces, and the propagation path of the multiple reflections is asymmetric; or multiple reflections occur within a thin layer. In general, the short range multiples and the peglegs do not differ in a strict sense, as shown in fig. 4.

(4) Ghost (ghost wave): the explosive source may form a virtual reflection in the presence of a wave impedance interface, as shown in fig. 5. When the seismic source is excited in the sea, a part of the excitation energy is transmitted upwards, and after being transmitted to the sea surface, the excitation energy is reflected downwards from the sea surface. This is known as ghost, and ghost wavefield records will have a delay time t, the value of which is the time of the two-way travel of the wave from the source to the surface of the water, as compared to seismic waves traveling directly down from the excitation location.

Free surface multiples: there is a special multiple at sea, called the free surface multiple, also known as the ringing, which is caused by multiple reflections in the water layer. This is because sea level and sea bottom surface are strong reflection interfaces, which easily causes wave field to propagate back and forth, causing multiple interference. If the sea floor has variations in elevation, this noise is called reverberation because the wavefield scatters and interferes with multiples. If the sea floor is relatively flat and belongs to an interface where the reflection coefficient is relatively stable, the energy entering the water layer will cause multiple reflections in the water layer, which is known as ringing. In fact, reverberation and ringing are not substantially different. They have a stable sinusoidal-like waveform, last a long time and mask the active wave. In the processing of marine seismic data, the interference wave is the most obvious interference wave.

2.2.3 identification of multiples

The generation of multiples is closely related to the underground geological conditions, and because the geological conditions of the regions are different and the forming principles and types of the multiples are different, the characteristics of the multiples displayed on the data records are also different. How to identify and eliminate multiple waves in seismic data processing and recover the original characteristics of seismic data is an urgent problem to be solved in marine seismic data processing. In general, various seismic data can be used to identify multiples from several aspects:

(1) direct recognition of multiples on shot gathers, CMP gather, stack profiles

In some areas, the multiples are reflected very strongly on shot gathers or CMP gather, and can be judged as multiples with certainty.

(2) Resolving multiples by primary folding profile and offset profile

The preliminary stacking section can be used for further analyzing the distribution condition of the multiple waves, the performance difference of the multiple waves on the stacking sections with different offset distances is large, and the multiple waves can be identified accordingly. Based on the morphological expression of the seismic section, the multiple waves are presented on the seismic section with smaller offset distance, and the primary waves are presented on the seismic section with larger offset distance.

(3) Identifying multiples using autocorrelation and spectral analysis

Since the multiple waves, particularly the ringing, have obvious periodicity and repeatability in waveform, the correlation is relatively good. Therefore, autocorrelation can be used to further determine the periodicity of multiples and select the parameters for suppressing multiples based on the autocorrelation.

(4) Identification of multiples by velocity spectra or velocity sweeps

In the velocity spectrum, the velocity of the multiple wave appears very clearly as a series of multiple reflections corresponding to the second intense energy bolus. Alternatively, the multiple can be identified by scanning the superimposed section at a constant speed using a series of velocities, the shape and imaging of which will be very different at different velocities for the same location on the section.

The identification of multiples sometimes requires only one method, and sometimes requires a combination of methods to identify.

2.3 conventional treatment method for offshore multiples

Conventional marine streamer wavelet flattening methods typically use wave equation inversion to flatten the multiples to predict a multiple model that is equivalent to the multiples contained in the original seismic record. The predicted multiples model is not exactly the same as the multiples in the original recording, they differ in amplitude phase and arrival time, and if directly subtracted, they do not serve the purpose of suppressing the multiples in the original recording. The multiple wave model predicted is consistent with the multiple wave in the original record through the adaptive matched filtering method, and then the process of attenuating the multiple wave is realized by subtracting the predicted multiple wave after the adaptive matched filtering from the actual original record. The method of first predicting multiples and then subtracting multiples is called a multiple prediction subtraction method, which is called a prediction subtraction method for short.

The prediction subtraction method is generally divided into two steps, first predicting the multi-wave model, and then adaptive matched filtering. For the prediction of the multi-wave model, prior information does not need to be known in advance according to the needs, and the method is divided into a model driving method and a data driving method. The principle of the model driving method is that multiple waves are simulated based on an acoustic wave equation according to the dynamic characteristics of a seismic wave field according to some prior information, such as the depth of a sea water layer, the sea bottom interface reflection coefficient, the seismic wave velocity and the like. The data driving method for predicting multiples does not need to use prior and posterior information, but actually records and inverts multiples, and iteration feedback methods such as Berkhout and Verschuur, backscattering series (ISS) methods such as Weglein and constant interpolation methods such as Sen belong to the data driving method. The adaptive matched filtering method may be classified into a least square method (2-norm method), a 1-norm method, a logarithmic-norm method, and a pattern recognition method according to a matching criterion.

The prediction subtraction method can be applied to complex underground structures, needs little or no prior information of the underground structures, can make better prediction on multiples in various forms, does not cause interference on effective reflection signals, inhibits the multiples, provides good seismic data for subsequent prestack data processing, and has much larger calculation amount than a direct filtering method.

The traditional method for suppressing multiples by using conventional streamer data is based on a wave field continuation method, and is used for predicting and suppressing multiples by combining forward continuation and reverse continuation of wave field information. For wave field information received on a sea surface streamer, firstly adopting a certain wave field continuation means to carry out forward continuation of a wave field towards the seabed direction, wherein the process is equivalent to that when a section of wave field is artificially added in a sea water layer, the process is equivalent to that a seismic wave field forwards advances for a distance along a time axis on a time section; then, wave field information of sea level is extended downwards to a seabed interface in a reverse extension mode, the process is equivalent to that when a section of travel is reduced, a distance is 'backed' on a time section of a wave field, which is equivalent to that of a seismic wave field, and the seismic wave field is returned to a certain previous moment; and finally, filtering wave field information generated by two times of wave field continuation by adopting a certain self-adaptive matched filtering method on the seabed, wherein the time difference between the wave field information and the seabed is just the time difference expression of multiples, so that the effective wave field information recorded by suppressing the multiples is obtained after the filtering treatment. And then extending the effective wave field to the sea surface in a reverse extension mode, thereby realizing the pressing process of the multiple waves. Where the forward continuation may be less accurate when the seafloor topography is complex, this disadvantage may be reduced accordingly if the location of the matching subtraction is made on the seafloor. The forward continuation uses the convolution of the simulated single-shot seismic wave field record and the green function corresponding to different positions in the sea bottom and the sea water layer, and the reverse continuation is related processing.

2.4 marine multiple forward modeling method

On the basis of understanding the mechanism, classification, characteristics and identification method of multiple generation, the part uses a forward simulation mode to carry out simulation analysis on multiple problems encountered under the receiving condition of the marine stereoscopic observation system, so as to obtain a simulated seismic data record. Through forward modeling of the received seismic records, the multiple wave problem in the actual operation process on the sea can be deeply known and understood, and a good foundation is laid for seismic data processing. The forward modeling method adopted by the invention is a finite difference method based on an acoustic wave equation. The multiple forward modeling analysis in this section can provide a good data base for the separation and joint suppression of multiples in later sections.

The seismic wave field forward modeling technology is helpful for people to better understand the geometrical characteristics and the dynamic characteristics of wave field information in geological conditions, and is used for determining appropriate field acquisition parameters and verifying correct interpretation of seismic data. The seismic wave field forward modeling methods are very many, and the forward modeling methods used in actual work are different according to different simulated media. The geophysical forward modeling method mainly comprises physical modeling and numerical modeling, wherein the physical modeling is usually implemented under excessively complex conditions and is not convenient to be widely applied, and the numerical modeling is simple, high in calculation speed, easy in parameter setting and high in practicability. The commonly used numerical simulation methods mainly include ray tracing, integral equation and wave equation. The ray tracing method reflects the kinematic and geometrical characteristics of the seismic waves; the wave equation method is based on the principle and technology of seismic wave mechanics, and analyzes the propagation rule of seismic wave in corresponding physical parameters of the underground medium model, and mainly reflects the dynamic characteristics of the seismic wave. The seismic wave field information coverage of the seismic wave field is wider, and the dynamic characteristics of the seismic waves can be explained more reasonably and more perfectly.

Longitudinal waves in conventional seismic wave propagation media are expansion and contraction strain waves that have the same characteristics as acoustic waves in fluids. The acoustic wave equation can be used approximately to analyze the actual propagation problem of the seismic wave. While the acoustic medium is the simplest approximation in the earth's medium, it is the most widely used in practical seismic exploration. The two-dimensional acoustic wave equation for isotropic heterogeneous media can be written as (Aki and Richards,1980):

in the formula (I), the compound is shown in the specification,

is a gradient operator; p is pressure; c is the propagation speed of the wave in the medium; p is density; f is the bulk stress. This expression is a second order partial differential equation in the time domain. Similarly, the second order partial differential equation can be expressed as a first order system of velocity and stress components. The first order stress-velocity equation form of the two-dimensional acoustic wave equation can be expressed as:

in the formula, v _x And v _z Is the particle velocity; p is the normal stress.

With the interleaved method grid, the spatial location of the corresponding acoustic wavefield component is shown in FIG. 6.

A2N-order spatial difference precision and a second-order time difference precision staggered grid high-order finite difference format are given below, and

respectively stress p, velocity v _x And v _z Discrete values at (i, j, k). The difference format of equation (2-3) is:

in the formula, a _n Is a finite difference coefficient; Δ x and Δ z are the grid spacing in the x and z directions, respectively; Δ t is a time step.

Fig. 7 and 8 show the forward modeling results of multiple waves on the free surface under the simple model of seawater-seabed two-layer medium, and it can be clearly seen that the multiple waves are distributed in a hyperbolic curve form and have better periodicity. The method can carry out clear forward modeling on multiple waves of different orders, and forward results are used as an important basis for subsequent seismic data processing.

3 the invention is further described below in connection with the offshore vertical cable multiple wavefield separation and continuation method.

The offshore vertical cable technology is that the cable which is originally horizontally placed is vertically placed below the sea surface to be in a vertically placed state, a remotely controllable buoy is placed on the upper part of the cable to enable the cable to have upward floating tension, and a base is arranged at the bottom of the cable and can be fixed on the sea bottom and is provided with an angle sensor. The acquisition process comprises seismic waves excited by an air gun seismic source on the marine survey vessel, and hydrophones are carried on the vertical cables to receive information from the bottom layer of the sea bottom and are transmitted back to the marine survey vessel through radio. The method is mainly used for geophysical survey for ocean stereo observation. As shown in fig. 9:

the principle of marine vertical cable observation is the same as that of land based Vertical Seismic Profiling (VSP). VSP is a method of placing receivers at different depths down the surface excitation, i.e. looking for the wavefield in the vertical direction. Since the geophones are located below the surface, it is a significant feature of the VSP method to receive both down-going waves (direct + down-multiples) and up-going waves (effective reflected + up-multiples) from above the geophones. The special position relation of the seismic source and the geophone breaks through the limitation of the traditional observation mode in the conventional surface seismic exploration, and more types of seismic waves can be observed than the conventional surface seismic. Therefore, there is a unique set of processing flow and method for processing VSP data.

The VSP seismic wave field mainly comprises direct waves, primary reflected waves, multiples and the like, wherein the VSP seismic wave field can be divided into uplink multiples and downlink multiples according to different propagation directions of the multiples.

(1) Direct wave. It is the seismic wavefield where the wavefield signal is directly propagated from the seismic source to the downhole detector, where a in fig. 10 is the direct wave.

(2) A primary reflected wave. The wave is a wave which is transmitted from a seismic source to a stratum below a detector and then reflected to the detector to be received after passing through a reflection interface, and b in the graph 10 is a reflected wave. The reflected wave is the effective wave required by seismic data processing, and reflects the reflection coefficient and the structural form of the underground reflection interface.

(3) Multiples. It is the wave received by the detector after multiple reflection is generated at the interface, c and d in fig. 10 show the multiple in the VSP wave field, the multiple includes the down multiple and the up multiple, c is the down multiple, it is the wave field that the multiple is transmitted to the detector after the surface reflection; d is an upgoing multiple belonging to the wave field of the multiple received by upward propagation after being reflected by the reflecting interface.

According to the analysis, the VSP records usually contain several different wave fields, and the different wave fields present different forms and distribution rules in the seismic records and are mixed together in an interleaving mode, so that the distinguishing and processing work of effective wave field information is seriously influenced. These several different wavefields contain subsurface medium formation and lithology information, and only if these wavefields are effectively separated, can the subsurface formation be accurately imaged. Therefore, it is critical to perform good VSP wavefield separation.

3.1 offshore vertical cable wavefield separation method

In conventional marine vertical cable seismic data processing, the original seismic records include both the upward wave that is reflected by the subsurface reflection interface and then travels upward, the downward wave that is transmitted to the ground surface after being reflected by the ground surface and then travels downward, and the direct wave that is directly transmitted to the detector position after being excited by the seismic source, as shown in fig. 11. The seismic wavefields recorded by the vertical cables are much richer than those recorded by the streamers, which is an advantage of the vertical cable observation. Therefore, for the data information of the vertical cable, a first task is to separate the up-going wave and the down-going wave of the vertical cable. Because there is downlink wave information such as direct wave in addition to uplink wave information in the actual vertical cable data record, and the uplink wave and downlink wave can be received and recorded simultaneously under the ground, but the uplink wave and downlink wave are usually mixed together, the conventional seismic processing work often uses the uplink wave field information or downlink wave field information alone for data interpretation, so the uplink wave field and downlink wave field are separated first before the subsequent processing work. In the vertical cable data, the up-going wave and the down-going wave have distinct apparent velocities, so the separation of the up-going wave and the down-going wave of the vertical cable is mainly based on the difference of the apparent velocities of the two waves. In vertical cable data, the travel time of a downlink wave increases with the distance of propagation depth, and the apparent velocity value is positive; the upstream wave travel decreases with increasing depth of propagation and is therefore negative in apparent velocity. Thus, the wave field separation can be performed according to the difference of the apparent velocities of the up-going wave and the down-going wave. The commonly used wave field separation methods mainly include F-K filtering, median filtering, Radon transformation, etc.

3.1.1F-K domain wave field separation method

The F-K domain wave field separation method is to convert the wave field information of t-x domain into F-K domain, i.e. frequency-wave number domain, by means of Fourier transform of the original seismic record. Therefore, the up-going wave and the down-going wave which are originally mixed together can be well separated from the F-K domain, and have different characteristics and position information, so that if the filtering processing is realized in the F-K domain, the up-going wave and the down-going wave can be well separated without damaging the characteristics of respective wave fields. And then the wave field of the F-K domain is inversely transformed back to the t-x domain, and the separation processing of the up-going wave and the down-going wave is completed. The F-K domain wave field separation filtering mode mainly has the following advantages:

(1) this method is very easy to implement because the up-and-down waves in the vertical cable records are at different interval positions after being transformed into the F-K domain, and therefore, the wavefields that distinguish and separate each other can be very easily processed, as shown in fig. 12;

(2) the original records with different apparent velocity wave field information, such as P-waves and S-waves, can be obviously identified in the F-K domain due to the different apparent velocities and processed through wave field separation.

The frequency-wavenumber domain wavefield separation principle steps are shown in fig. 13. Firstly, the original input data information of the vertical cable, as can be seen from the figure, in the original vertical cable wave field, the upgoing wave and the downgoing wave are relatively easy to identify but are mixed together, and the wave field separation work is difficult to carry out, so that the processing directly in the t-x domain is difficult. After the original vertical cable seismic data are subjected to F-K conversion, wave field energy information of different interval positions can be displayed in a frequency-wave number domain, because the upper traveling wave and the lower traveling wave are positioned at different plane positions, the wave field separation is very easy, if the energy at the position of the lower traveling wave is directly filtered and suppressed by reserving the upper traveling wave, the wave field information of the lower traveling wave can be obtained in the same way. The up-wave field is reserved in the figure, so that the wave field information of the down-wave is removed in the frequency-wave number domain, the remaining wave field information is subjected to inverse F-K conversion, and after the wave field information returns to the t-x domain, the wave field record which is used for separating the down-wave and only reserving the up-wave is obtained. It can be seen that after the F-K transformation and the inverse F-K transformation, the wave field separation effect is very good, and a good data base is laid for the subsequent migration and continuation work.

In order to improve the data processing precision, the invention adopts an F-K domain wave field separation method to separate the vertical cable up-going wave and the vertical cable down-going wave. And in the wave field continuation process, continuation noise and foldback effect are suppressed by using an F-K domain filtering method. And (4) performing wave field continuation on the downlink wave on the vertical cable to predict multiple waves in subsequent processing, and performing matching pressing. As shown in fig. 14.

3.1.2 median filtering separation method

The median filtering method has the advantages that the time for introducing digital signal processing is very early, the method has wide application in data processing of geophysical data, the method also achieves better effect in vertical cable wave field separation, and the improved median filtering method has more ideal effect on the vertical cable wave field separation.

The median filtering method comprises the following steps: when filtering the digital signal sequence x (j) (∞ < j < ∞), an odd-length L-long window is first defined, where L is 2N +1 and N is a positive integer. At some point, the signal samples in the window are x (i-N), …, x (i), …, x (i + N), where x (i) is the signal sample value at the center of the window. After the L signal sample values are arranged from small to large, the median value x (i) is taken as an output value.

The process of median filtering is essentially to perform a smoothing process on a series of data and then output a median value, and is an optimal filtering method based on the 'minimum absolute error' criterion. Table 1 describes the procedure of the median filtering method.

TABLE 1 median Filter Algorithm

The main principles of median filtering of vertical cable wavefield data can be summarized as: as the apparent velocity signs of the up-going wave and the down-going wave in the wave field data are opposite and accord with the difference of positive and negative signs, the slopes of the in-phase axes displayed in the seismic section are different, and the vertical cable wave field separation median filtering method is used for carrying out VSP wave field separation by utilizing the characteristic.

For vertical cable seismic data, the median filtering wave field separation method comprises the following steps (taking the retention of traveling waves as an example):

(1) aligning and leveling the recorded downlink waves;

(2) performing median filtering operation on the leveled record to obtain a median result;

(3) performing reverse time shifting on the obtained result;

(4) and subtracting the filtered result from the original data to obtain the uplink wave in the original record.

The basic process of the median filtered wavefield separation method is shown in fig. 15. The key point is that the non-horizontal homophase axis is filtered by selecting a proper span interval, so that the seismic record with useless wave fields eliminated can be obtained, and the effective seismic wave homophase axis is kept without interference.

The median filtering method can achieve very good results for vertical cable wavefield separation, as shown in fig. 16. The median filtering method is relatively simple to implement and is relatively practical. However, the data it can process must satisfy the assumption that the waveforms are substantially consistent. Data received by complex regions or different lithological strata are complex, and the median filtering effect in the regions is not ideal.

Taking the vertical cable down-going wave field as an example, the processing effect of the F-K domain wave field separation method and the median filtering method on the vertical cable seismic data is compared, as shown in fig. 17. The left graph represents a downlink wave field obtained after being processed by an F-K domain wave field separation method, and the right graph represents a downlink wave field obtained after being processed by a median filtering method. Therefore, no matter which filtering method is adopted, the requirements on the separation of the up-going wave and the down-going wave of the offshore vertical cable can be met. However, in the actual processing process, the filtering span value of the median filtering method is sometimes difficult to select, and the span value is selected too much, which may damage the signal energy of the effective wave field and affect the retention of effective information; if the span value is too small, the filtering effect may not be complete and the wavefield separation may not be clean. Therefore, the invention adopts the F-K domain wave field separation method, reasonably filters out useless wave field energy through the interval analysis of the F-K frequency spectrum, and can more effectively separate the up-going wave and the down-going wave of the vertical cable wave field.

3.2 offshore vertical cable wave field continuation method

The method is used for processing vertical cable seismic data based on a wave equation method, and a very key step is how to carry out wave field continuation. The wave field values at different space points can be known by a wave field continuation method, and the simulation processing work of the seismic wave field is carried out. For the relatively complex stratum condition, the wave equation method numerical simulation has better processing effect than the ray tracing method, but the operation speed is slower. When wave field calculation is performed by using a wave equation method, a one-way up-and-down traveling wave is generally adopted to perform wave field continuation and extrapolation, and compared with a full wave field extrapolation mode, the wave field calculation method has the advantages of simpler calculation form and higher calculation speed. As shown in figure 18 marine VSP reception mode.

When the wave equation carries out wave field continuation, the most key calculation form is exp [ -i ω (t ± r/v) ], and r in the formula represents the relative length of the seismic source position and the observation position. Generally, the one-way wave equation is utilized to carry out wave field continuation, namely, the wave field advancing towards a single direction only, the uplink one-way wave is propagated upwards, and the downlink one-way wave is propagated downwards, so that the calculation form can be greatly simplified. The sign in the formula represents the downlink one-way wave if the sign is negative, and represents the uplink one-way wave if the sign is positive. In order to calculate as accurately as possible, and to reduce other influencing factors, it is generally assumed artificially that the propagation medium of the wave field is a homogeneous medium.

The three-dimensional wave equation is written as follows, as shown in equation (3-1):

fourier transformation is carried out on the formula (3-1) to obtain:

when the wave field continuation is carried out, the up-down traveling wave can carry out the wave field extrapolation process, and the forward simulation process can be realized by adopting a forward continuation mode, namely performing wave field extrapolation along a forward propagation mode of an actual seismic wave field from a seismic source position; the wave field information obtained from the position of the wave detector is transmitted against the propagation direction of the seismic waves in a reverse continuation mode, and the final result cannot be different due to the selection of the continuation mode. But in the actual process, which extrapolation method is more suitable to be selected can be emphasized for the actual target problem.

The wave equation continuation method used by the invention is a one-way wave extrapolation forward modeling method based on the Huygens principle, and has the advantages of simple and easy operation form and high operation speed in the data calculation and processing process, and the seismic data profile obtained by forward modeling is very clear.

3.2.1 one-way wave equation

Wave equation migration methods for seismic data are based on wave field continuation extrapolation, and are generally based on one-way wave equation wave field extrapolation methods in seismic migration processing. Seismic data is typically transformed into the F-K domain by a 2D Fourier transform, while the velocity function is constant, time frequency is converted into vertical wavenumbers in the frequency-wavenumber domain (Stolt, 1978), and an inverse Fourier transform is performed to obtain seismic data images. The mapping is for seismic data post-stack migration. In general, the 2D Fourier transform is applicable under a constant velocity condition, that is, the seismic wave velocity does not change with the change of the spatial coordinate position. Similarly, the phase shift migration used by Gazdag is also performed in the F-K domain, but the frequency transformation is not required, and the seismic wavefield only needs to be phase-shifted once in the frequency-wavenumber domain, so as to realize the process of extrapolating the wavefield downward by one step. The seismic wave velocities of different depths are respectively used during each phase shift, so that the condition that the seismic wave velocity changes along with the depth can be adapted.

However, neither the Stolt method nor the Gazdag phase shift method can solve the problem of the transverse change of the seismic wave velocity. The phase shift plus interpolation method proposed by Gazdag (1984) is a good improvement over conventional phase shift method shifts, which can be applied approximately to simple lateral shift situations. The time domain finite difference method based on the delay coordinate system can well adapt to the longitudinal and transverse changes of the seismic wave velocity, but the finite difference method has very obvious dip angle limiting conditions, so that the method is very difficult when facing the problem of the high-order time domain wave equation. The splitting step Fourier method proposed by Stoffa (1990) can be seen as a further improvement of the phase shift method. It can handle moderate shifting conditions of formation conditions.

Starting from a two-dimensional scalar wave equation, one can describe the form of wavefield propagation in a medium with constant density and longitudinal velocity v (x, z):

where x is the horizontal axis of space, z is the depth axis (positive in the vertical direction), t is time, p is the wave function, and v represents the velocity of longitudinal waves. In seismic exploration, it is generally considered that the velocity of a seismic wave in the time for acquiring seismic data does not change significantly with time, so equation (3-3) can be subjected to Fourier transform over time, and Fourier transform can be used in the x and z directions provided that the velocity does not change spatially. The invention adopts the following Fourier transformation form:

p(k _x ，k _z ，ω)＝∫∫∫p(x，z，t)exp(ik _x x+ik _z -iωt)dxdzdt (3-4)

p(x，z，t)＝∫∫∫p(k _x ，k _z ，ω)exp(-ik _x x-ik _z +iωt)dk _x dk _z dω (3-5)

replacing the Fourier forward transform (3-4) into equation (3-3) to solve the wave equation in the form of Fourier domain:

write (3-6) to the form

Equation (3-7) represents the wave dispersion form of the one-way wave equation, with positive going waves on the right of the equation and negative going waves on the right of the equation.

The corresponding expressions of the wave equation t-x field and f-k field are proposed by Claerbout:

bringing (3-8) into (3-7), the wave dispersion form of which can be changed into a one-way wave equation:

equation (3-10) is the theoretical basis of the phase shift method wave field continuation and finite difference wave field continuation method.

3.2.2F-K Domain phase shifted wavefield extrapolation

In 1978 Gazdag proposed phase shift method shifts. It can be verified that the harmonic solution of equation (3-4) is:

the solution of the up-going wave and the down-going wave represented by the above equation is:

based on equations (3-11) and (3-12), the wavefield information at the depth z + Δ z position can be derived from the wavefield information at the depth z position, for the up-going wave:

for down-going waves

Suppose that the wavefield information p (k) for the earth's surface has been obtained in advance _x And z is 0, ω), then the seismic wavefields for the different spatial location points in the subsurface may be derived step by step based on either (3-13) or (3-14).

To summarize, the computation steps for the phase shift method to achieve wavefield extrapolation are as follows (above-going wave-down extrapolation is an example):

(1) first, the surface wavefield p (k) _x Z-0, ω) to the frequency-wavenumber domain p (k) _x ，z＝0，ω)。

(2) Calculating an extrapolation factor

(3) In the frequency-wavenumber domain, the seismic wavefield is multiplied by a factor C to obtain the wavefield at Δ z.

(4) And (3) repeating the steps (2) and (3) to gradually obtain frequency domain wave fields of all underground depths i delta z (i is 1,2 and 3 …).

(5) Wave field p (k) of frequency-wavenumber domain _x And z is equal to i delta z and omega, and is inversely transformed back to a time-space domain p (x, z is equal to i delta z and omega), and finally, time-space domain seismic wave fields at different depths are obtained. As shown in the phase shift wavefield extrapolation flowchart of fig. 19.

3.2.3F-K Domain continuation noise suppression

Due to the fast and accurate algorithm, the phase shift wave equation extrapolation is a very important migration processing method in the field of seismic exploration. However, when the method is actually applied to the seismic data migration processing process, the phase shift extrapolation method is prone to have a serious foldback effect at the boundary, each continuation can cause strong continuation noise, and the identification of effective waves is greatly interfered. In order to reduce the effect of the fold-back effect, the corresponding fold-back response can be suppressed by additionally adding one time of blank seismic traces on both sides of the seismic section and eliminating zero values after migration, but the input section is enlarged, the data volume is artificially doubled, and simultaneously, the calculation workload and the requirement on a computer memory are increased, so that the calculation efficiency is greatly reduced. And a certain absorption method can be adopted at the boundary of the space domain to suppress the folding back effect. However, since the absorption is done in the spatial domain and the wavefield is extrapolated in the wavenumber domain, the forward and inverse Fourier transforms must be performed in sequence in each extrapolation step. The operation speed is slower when the continuation step length is more. The invention improves the frequency-wavenumber domain (F-K domain) attenuation method, performs spectrum analysis on the continuation process, suppresses linear interference by using the F-K filtering method, has simple algorithm and high calculation speed, does not need to additionally increase the burden on a calculation system, and can keep the effective wave frequency band unaffected. Simulation results show that the filtering method can effectively overcome the boundary effect in the frequency-wavenumber domain. Because the phase shift offset and the wave number filtering are carried out in the wave number domain at the same time, the forward and backward Fourier transformation does not need to be realized in the process of each continuation step length. Therefore, the processing method can greatly improve the calculation efficiency of the phase shift method offset processing.

The wave field continuation process of the wave equation migration of the phase shift method is mainly carried out in a frequency wave number domain, and the main formula is as follows:

p(k _x ，z，ω)＝p(k _x ，z＝0，ω)exp[-ik _z z] (3-15)

according to the above formula, p (k) _x Z, ω) and a phase shift factor exp [ -ik) _z z]The seismic wavefield received at the surface can be extended backward to any depth z in the subsurface. Of course, equation (3-15) is a solution obtained when velocity v is constant, let z be n × Δ z, v (z) be v (n Δ z), and if Δ z is small enough that velocity is constant within Δ z, then extrapolate Δ z step by step from z 0 position to the subsurface depth as a step, using different velocities v (n Δ z) in each continuation step, then the vertical shifting of the wavefield continuation formula can be achieved using the constant velocity formula:

p(k _x ，z _i+1 ，ω)＝p(k _x ，z _i ，ω)exp[-ik _z Δz] (3-16)

in the formula: Δ z ═ z _i+1 -z _i 。

Since the in-phase axis on the time profile always moves toward the tilt-up direction after the offset processing is performed, it is an essential step to perform the field extension processing for the offset. For the wavefield information at a certain position on the event of the seismic section, the data point will move a distance from the original position in the upward-inclined direction when the wavefield descent delay is performed, in other words, if the data point is not originally at the edge position of the section, after the wavefield descent delay step, it may reach or approach the edge position of the seismic section, which is called energy edge clustering. Because wavefield continuation proceeds step-wise downward, energy clustering is a continuous dynamic process, and although the energy at the initial boundary position is attenuated, at the time of downward continuation, the newly added unsuppressed energy occurs again at the boundary position, which still causes the boundary reentry effect. A feasible approach is to add boundary absorption conditions for each continuation step and all the energy reaching the boundary is absorbed in time. Since the whole process is energy absorption in the space domain and extension step is implemented in the wavenumber domain, for p (k) in the formula (3-16) _x ，z _i ω) is k _x By inverse Fourier transform of (1), can obtain

p(k _x ，z _i+1 ，ω)＝p1(x，z _i ，ω) (3-17)

Then adding absorption conditions B (x) to the formula (3-17) to obtain p ₂ (x，z _i ω), B (x) are selected according to the following principle: the decay curve continues as narrowly as possible to minimize the impact factor on the seismic traces.

By adopting a recursion method, order

K for the absorbed wavefield _x Is transformed and multiplied by a phase shift factor exp < -ik _z Δz]To obtain

p(k _x ，z _i+1 ，ω)＝p(k _x ，z _i+1 ，ω)exp[-ik _z Δz) (3-19)

Because the space domain absorption method has the defects of low calculation speed, low efficiency and the like, the invention adopts a frequency wave number domain (F-K domain) absorption method to suppress the foldback effect. Because the phase shift extrapolation is realized in the frequency wavenumber domain, the positive and negative Fourier transform is not needed for the wavenumber domain absorption, and the method can improve the calculation efficiency of the phase shift offset. Wave number domain absorption does not need to convert k into _x The wavenumber is converted into the spatial domain, but k for wavenumber _x Multiplying a certain attenuation coefficient corresponds to a filter factor, so that this is also true when the number of waves is subjected to filter processing, that is:

changing the weighting factor k in the equation (3-20) _x The filtering processing of the wave number can be realized.

Because the absorption and continuation of the wave field information are carried out in the F-K domain at the same time, the wave number domain absorption can weaken the K _x The high wave number signal energy in the direction can obviously inhibit the reentry effect, weaken the dispersion phenomenon and improve the quality of the offset processing.

When the vertical cable wave field continuation is carried out, although the rapid F-K phase shift algorithm adopted by the invention has the advantages of simple calculation form, high calculation speed and the like, due to the defects of the phase shift method, the seismic wave field energy can be gradually gathered towards two sides and slowly expanded to the boundary when the phase shift method is continued downwards by one step, so that stronger boundary reflection energy is caused. The presence of the reentry effect can seriously interfere with the identification and analysis of the active waves in the seismic wavefield, and the linear interference in the form of "cross lines" not only causes the damage of the extended imaging quality, but also always mixes with the effective energy, greatly reduces the signal-to-noise ratio of the wavefield information and must be removed. The invention adopts the F-K domain noise absorption method, and carries out absorption attenuation treatment on continuation noise caused by each wave field continuation in the frequency-wave number domain, thereby better suppressing the interference energy on the seismic section and laying a good foundation for matching correction and filtering algorithm in the subsequent work.

FIG. 20 shows the processing effect of the F-K domain noise absorption method. The graph A is a wave field continuation result without filtering processing, and strong linear noise interference can be seen, which is a boundary reentry effect generated by the defects of a phase shift method, and obvious energy aggregation can be seen near end points on two sides of an earthquake homophase axis. This has a great influence on the post-processing. And the B picture is the result of performing absorption filtering in a frequency wave number domain by using an F-K filtering method, and the B picture shows that the suppression effect is very good, the linear interference is effectively inhibited, the boundary foldback effect energy is weakened, and the end effect is not obvious any more. Moreover, the processing result can also show that the dispersion phenomenon of the signals after wave field extension is attenuated to some extent compared with the wave field which is not absorbed, which shows that the wave number domain not only can well absorb the boundary reflection, but also can have the filtering function of reducing dispersion and enhancing the energy of effective signals.

4 the present invention is further described below in conjunction with a method for processing multiples by combining a marine stereo observation system.

After the seismic records obtained by the vertical cables are subjected to wave field separation by using the method, the separated upgoing wave and downgoing wave seismic records are obtained. Next, multiple wave prediction is needed and matched with the horizontal cable data, wherein the invention mainly analyzes the downward wave of the vertical cable. Since the seawater surface is a strong free reflection interface with a reflection coefficient of-1, this is equivalent to the propagation mechanism of multiple waves on the free surface. When extra travel is added between the sea surface and the sea bottom, the seismic wave field which is reflected by the sea surface and then transmitted back to the sea bottom and received by the wave detector is just the wave field of free surface multiples. However, because the vertical cable detectors are arranged along the depth direction, the received multiple wave field cannot be directly used in the matching pressing process of the multiple waves. Therefore, a horizontal correction process is required for the received multiple wavefield, and after a certain transformation, the wavefield of the free surface multiple received by the vertical cable is corrected to the multiple wavefield received by the horizontal cable. Because the wave fields received by the vertical cables comprise the up-going wave field and the down-going wave field, the multiple waves are easier to identify, and the wave field separation is easier to realize according to the apparent velocity difference, the multiple wave field separated and predicted by the vertical cables is used for matching and suppressing the multiple wave interference received on the horizontal cables, so that the method has better precision, simple processing method and higher data processing quality.

In order to apply the vertical cable multiple information to the aspect of horizontal cable multiple suppression, the seismic wave fields of the vertical cable multiple information and the horizontal cable multiple information need to be calibrated to a uniform reference surface in a wave field continuation mode, and then matched suppression is carried out. The selection of the correction reference surface can select any plane in the water layer, and when continuation correction is carried out, the key steps are that the wave field of seismic waves at a certain depth is subjected to wave field continuation based on a wave equation theory, and the wave field with known depth is extrapolated to an unknown depth position to obtain the wave field characteristics at the position. The vertical cable wave detectors are arranged along the depth gradually, so that wave field information received by the wave detectors at different depths has a certain calculation form and quantity relation in a wave field continuation equation. Through a certain wave field continuation mode, the wave field information received by the detectors with different depths can be subjected to wave field extrapolation along different paths, for example, the wave field value at a certain detector position is continued to be subjected to forward continuation of the wave field, and the wave field information at the corresponding position can be obtained by continuing to be continued to the sea bottom or the sea surface. Similar to the elevation static correction problem in land seismic exploration, the basic method is to set a horizontal correction plane with a certain height, obtain the vertical elevation between each detection point and a reference surface, and divide the vertical elevation by the seismic wave velocity to obtain the corrected time difference, so that the time difference correction can be carried out on the seismic records. However, the acquisition modes of the vertical cables are different, and particularly when the vertical cables and the detectors of the horizontal towline jointly receive seismic reflected wave signals in a stereo observation system, primary effective wave and multiple wave fields received by the detectors on the vertical cables are not in a simple time difference relationship, and the reflected waves are not vertically transmitted in the transmission process, so that errors of wave field information positions are inevitably caused if only simple time difference correction is carried out, and vertical cable wave field correction cannot be carried out according to an elevation correction method of land seismic exploration. Therefore, the correction means adopted by the method is frequency-wavenumber domain (F-K) wave field continuation cable depth correction, 2D Fourier conversion is carried out on seismic wave field data to be in a frequency-wavenumber domain, and continuation correction is carried out on vertical cable data along the wave field propagation direction according to a one-way wave propagation theory. As shown in fig. 21 for extended correction of the stereo observation system.

4.1 continuation correction method of vertical cable data to streamer data

Since the vertical cable receiving mode is different from the horizontal cable, the detector is arranged along the vertical direction, and the reflected wave has different travel time difference when the wave propagates along the vertical direction. Thus at different depths, the received wavefield information is different. The time variation quantity delta t of the seismic wave field is linearly increased along with the increase of the offset distance when the seismic records are received on the horizontal cable, but the seismic data received by the vertical cable do not meet the relation, because the vertical cable seismic data and the horizontal cable seismic data lead the time of the upgoing wave (primary reflected wave) received by the vertical detector to be earlier than the primary reflected wave received by the sea surface detector by delta t, and the distance in the horizontal direction is correspondingly reduced by delta x; the downward ghost (ghost) received by the vertical cable receivers is delayed in time by Δ t relative to the primary reflections received by the surface receivers, and the distance in the horizontal direction is increased by a distance of Δ x. Aiming at the difference, a wave field continuation method is needed to be adopted to respectively extend the seismic wave field information received by the detectors at different depths to corresponding positions of the horizontal cable detectors so as to realize reasonable matching with the wave field received by the horizontal cable. In the invention, free surface multiples, i.e. ghost wavefields that propagate to the sea floor again after reflection at the sea level, need to be processed and suppressed. This appears as a down-going wavefield above the vertical cable receive record. Therefore, wave field continuation is carried out by utilizing a downward wave field of the vertical cable, seismic wave fields received by detectors at different depths are respectively pushed out to the sea bottom surface, and simultaneously, seismic records containing free surface multi-wave data received by the horizontal cable are also prolonged to the sea bottom surface, and the data of the two are matched. Because the data of the vertical cable continuation is the multiple wave data reflected back from the sea level, namely the ghost reflection data, the prediction work of the multiple waves between the sea surface and the sea bottom is completed in the wave field continuation process. As shown in the two-cable continuation matching process of fig. 22.

Therefore, based on the wave equation theory, the invention adopts a frequency-wavenumber domain wave field continuation method to extract the upgoing wave and the downgoing wave from the data of the vertical cable detectors with different depths, namely, firstly, the upgoing wave and the downgoing wave are separated from the total wave field of the vertical cable, and then based on a wave field continuation operator, the wave fields received by the detectors with different depths are respectively subjected to variable-depth wave field extrapolation with different paths and different distances, so as to obtain the seismic wave field value after extrapolation. And matching the multi-wave data obtained by extrapolation with the streamer data. The specific implementation method of the two-cable data matching will be described in detail.

In the marine seismic exploration process, as the sea water-air surface is a good reflection interface, for the receiving situation of a plurality of detectors, the surface is sunk to the depth z _i The wave field p (x, z) received by the detector _i T) the main component of the information includes the up-going wave field u (x, z) propagating towards the sea surface after being reflected by a certain reflection interface at the sea bottom _i T), and a down-going wavefield d (x, z) that travels down-going by reflecting again after traveling to the surface of the seawater _i ，t)。

p(x，z _i ，t)＝u(x，z _i ，t)+d(x，z _i ，t) (4-1)

Will (4-1) make 2D Fourier transform for offset x and time t, then

P(k _x ，z _i ，ω)＝U(k _x ，z _i ，ω)+D(k _x ，z _i ，ω) (4-2)

Wave field d (x, z) therein _i T) represents the seismic upgoing wavefield u (x, z) _i T) a downward wave field propagating down to the sea surface after reflection, so that when the reflection index at sea level is assumed to be-1, a wave field can be obtained

D(k _x ，z _i ，ω)＝U(k _x ，z _i ，ω)exp(-j2k _z z _i ) (4-3)

Substituting the formula (4-3) into the formula (4-2) to obtain

P(k _x ，z _i ，ω)＝U(k _x ，z _i ，ω)G(k _x ，z _i ，ω) (4-4)

In the formula G (k) _x ，z _i ω) is an imaginary reflection operator, also called imaginary reflection filter operator, and has

G(k _x ，z _i ，ω)＝1-exp(-j2k _z z _i ) (4-5)

Is set to a depth z ₁ And z ₂ The seismic wave fields received by the two detectors are p (x, z) respectively ₁ Y) and p (x, z) ₂ T) where the upgoing wavefield is u (x, z), respectively ₁ T) and u (x, z) ₂ T), for the frequency-wavenumber domain, the wavefield information for two depth positions can be derived based on the equation (4-1):

P(k _x ，z ₁ ，ω)＝U(k _x ，z ₁ ，ω)G(k _x ，z ₁ ，ω) (4-6)

P(k _x ，z ₂ ，ω)＝U(k _x ，z ₂ ，ω)G(k _x ，z ₂ ，ω) (4-7)

the up-going wavefield U (k) of the top detector can be known from the wave field continuation principle _x ，z ₁ ω) can be considered as the upgoing wavefield U (k) at the lower detector _x ，z ₂ ω) to the upper detector, there is

U(k _x ，z ₁ ，ω)＝U(k _x ，z ₂ ，ω)exp(-jk _z Δz) (4-8)

That is to say

U(k _x ，z ₂ ，ω)＝U(k _x ，z ₁ ，ω)exp(jk _z Δz) (4-9)

Δz＝z ₂ -z ₁ ；

Substituting the formula (4-8) into the formula (4-7) to obtain

P(k _x ，z ₂ ，ω)＝U(k _x ，z ₁ ，ω)G(k _x ，z ₂ ，ω)×exp(jk _z Δz) (4-10)

Order to

G ₁ (k _x ，ω)＝1-exp(-j2k _z z ₁ ) (4-11)

G ₂ (k _x ，ω)＝[1-exp(-j2k _z z ₂ )]exp[jk _z (z ₂ -z ₁ )] (4-12)

Then equations (4-8) and (4-9) can be written as

P(k _x ，z ₁ ，ω)＝U(k _x ，z ₁ ，ω)G ₁ (k _x ，ω) (4-13)

P(k _x ，z ₂ ，ω)＝U(k _x ，z ₂ ，ω)G ₂ (k _x ，ω) (4-14)

Therefore, (4-13) and (4-14) can be converted into seismic wave fields obtained by an upper depth detector and a lower depth detector

P(k _x ，ω)＝G(k _x ，ω)，U(k _x ，z ₁ ，ω) (4-15)

In the formula

P(k _x ，ω)＝[P(k _x ，z ₁ ，ω)，P(k _x ，z ₂ ，ω)] ^T (4-16)

G(k _x ，ω)＝[G ₁ (k _x ，ω)，G ₂ (k _x ，ω)] ^T (4-17)

Down-going wavefield U (k) using least squares optimal solution equation (4-8) _x ，z ₁ ω), even if minimized by the following objective function, where H is the Hermitian operator.

min：J＝[P(k _x ，ω)-G(k _x ，ω)U(k _x ，z ₁ ，ω)] ^H (4-18)

Simplified down-going wavefield U (k) _x ，z ₁ ω) of the measured value

Therefore, based on the wave field extrapolation algorithm, the wave field signals received by the vertical cable detectors with different depths can be gradually corrected into the wave field signals in the horizontal direction, and the prediction process of the free surface multiples is realized. The continuation algorithms used by the invention are all wave field continuation equations obtained based on a phase shift method, a predicted multiple wave field is obtained from the vertical cable, and then a proper matched filtering algorithm needs to be selected to be combined with the towing cable data for matched suppression.

4.2 multiple self-adaptive matching filtering method

After the multiples are predicted, the next operation is to subtract the predicted multiples from the original data. However, the predicted multiples are not true multiples. The method is not completely consistent with the actual multiples in the original data, and has larger difference compared with the actual multiples, if the predicted multiples are directly subtracted from the actual data, an error result can be obtained, and the method is not scientific. Aiming at the difference between the amplitude phase and the arrival time, a better filtering method needs to be selected for carrying out matched filtering processing. Therefore, the quality of the selected adaptive subtraction algorithm will directly affect the result of the multiple attenuation process. In recent years, different scholars at home and abroad do significant work on the adaptive algorithm, and a plurality of effective adaptive algorithms are analyzed.

FIG. 23 shows a comparison of waveforms after vertical and horizontal cable seismic data, respectively, have been extended to the seafloor. It can be seen that after the seismic data obtained by different receiving modes is subjected to wave field continuation, the seismic waveforms of the seismic data change, so that certain differences exist between the two waveforms in terms of amplitude, phase and arrival time, and a corresponding processing result cannot be obtained by direct subtraction. Therefore, selecting a suitable matched filtering method is an indispensable important step.

4.2.1 least squares Filtering

In the prediction process, the vertical cable data is used for predicting multiple waves on the horizontal cable through wave field continuation. After the multiples are predicted, the next step is to perform adaptive subtraction, and subtract the multiples predicted to come out of the vertical cable from the horizontal cable data.

If the above equation is converted into the time domain, a limit is imposed on the length of the filtering factor a. Then the

The above equation is transformed into:

in order to obtain the time domain filtering factor a (t) in the above formula, the effective wave data P after the multiple wave is eliminated is tried to be obtained ₀ The energy of (2) is minimized as much as possible.

Least squares filtering may also be referred to as wiener filtering. This filtering theory was first adopted by Wiener. And (4) analyzing the equation (4-20), if the energy of the output effective wave field is to be minimized, minimizing the output at the left side of the equation, and generally adopting a method of eliminating zero values. In this case, the equations (4-20) can be simplified as:

therefore, the filter factor a (f) can be directly obtained from the above equation, and the complete expression is:

in equations (4-23), for single frequency multiple data

If the absolute value is smaller, instability of the algorithm may increase.

To solve the above problems, the prior art proposes a solution, and a specific expression is as follows:

epsilon in the equation ² As a stabilizing factor, M ^* Is the complex conjugate of M. In this equation, if the amplitude of the multiple waves is relatively small, approaching 0, the value of the denominator of equations (4-24) also gradually approaches the magnitude of the stability factor. In other words, the filter factor value is solved according to the stability factor value. The value of the filter factor is generally chosen to be larger when the signal-to-noise ratio of the signal is relatively low.

Assuming that the input signal is a multiple m (t), the filter coefficient is a (t), the actual output is y (t) ═ a (t) × m (t), the desired output is p (t), and the output error is e (t) ═ p (t) -y (t), then the output energy error is q ═ Σ ₁ e ² (t)＝∑ ₁ (p(t)-y(t)) ² Based on the least squares method, the output energy error q is minimized as much as possible. The partial derivative is calculated for a based on the output energy error equation and changed to 0.

And (5) simplifying again to obtain:

∑ _τ a(τ)∑ ₁ m(t-τ)m(t-s)＝∑ ₁ p(t)m(t-s) (4-26)

where s represents the length of the filter factor. According to the definitions of auto-correlation and cross-correlation, wherein

Bringing the above equation into (4-26), then (4-26) can be rewritten as:

rewriting the above equation into a matrix form:

by solving this matrix, the filter factor can be obtained.

However, if the amplitude of the predicted multiples may be relatively small, there may be instability with respect to the filter factor calculation. To address this drawback, a stability factor may be added after the equation that calculates the energy error. Then the error energy expression becomes:

epsilon in the equation ² Is the stability factor coefficient. The formula (4-30) can be changed into the following formula according to the above calculation form:

the corresponding matrix becomes:

4.2.2 Single pass adaptive matched Filtering

The purpose of matched filtering is to minimize the difference between the predicted multiples and the multiples in actual seismic data, and then to implement more effective multiple attenuation after adaptive multiple attenuation, and the most important step of the method is to select the optimal matched filtering operator. The basic theory of single pass matched filtering is explained by a certain mathematical model.

x (t) is the seismic input trace to be processed; h (t) is a matched filter operator; y (t) is the actual output of the earthquake after filtering treatment, and the formula is (4-33); z (t) is the desired output; e (t) is the output error, and the formula is (4-34); q is the total mean square error, and when the input x (t) is a function of time, the expression of Q is formula (4-35).

y(t)＝h(t)*x(t)＝∑ _τ h(τ)x(t-τ) (4-33)

e(t)＝y(t)-z(t) (4-34)

Q＝∑ _t [y(t)-z(t)] ² (4-35)

It can be seen that for the matched filter factor h (t), there is always a corresponding mean square error Q, and the main role is to calculate the best matching operator to minimize the mean square error Q of the actual output y (t) and the desired output z (t) according to the input x (t) and the desired output z (t). I.e. a least squares based criterion for minimizing the residual energy.

It is assumed that the matched filter operator h (t) is finite, such as:

h (t) ═ (h (0), h (1), … …, h (s)), where (t) ═ 0,1, … … s)

Substituting the formula (4-32) into the formula (4-34) to obtain the formula (4-35):

according to the criterion of minimum residual energy, each value h (λ) of the matched filter operator h (t) should comply with the following constraint:

the equation for solving the single-channel matched filter operator can be obtained as follows

h(t)*r _xx (t)＝r _zx (t)(t＝0，1，...，s) (4-38)

First, a minimum mean square error relation is defined, which is obtained by the following formula (4-36):

both sides are divided by gamma _zx (0)，

Order to

Wherein gamma is _zz (0)≠0， (4-40)

Then the following results are obtained:

epsilon is the normalized mean square error that needs to be obtained, also called the total error energy. If e is 0, the actual output at that time is identical to the desired output, i.e. in the case of the most perfect match. Therefore, the method of calculating the matched filter operator h (t) according to the minimum residual energy criterion is substantially equal to finding the best matched filter operator h (t) to minimize ε.

Therefore, the single-pass adaptive matched filtering process mainly comprises the following two steps:

(1) when the seismic data is processed, x (t) represents a predicted multiple model; the desired output is raw seismic data containing multiples represented by z (t). By solving equations (4-33), the matched filter operator h (t) can be found.

Assuming that the length of the matched filter operator is s +1, equations (4-33) are written in matrix form:

equations (4-42) are Toeplitz equations. By calculating the Toeplitz equation, the single-pass matched filter operator h (t) can be calculated.

(2) After the matching operator is solved, the multiple can be reduced by formula (4-30)

p(t)＝z(t)-x(t)*h(t) (4-43)

Where p (t) is the seismic data after subtraction of multiples, z (t) is the original seismic data containing multiples, x (t) is the predicted multiples model, h (t) is the matched filter operator, and denotes convolution. In fact, x (t) h (t) is the single matched multiple, i.e., the multiple that is actually subtracted.

FIG. 24 shows a comparison of the extended wavelet waveforms of a streamer and a vertical cable after a wavefield extension process. It can be clearly seen that the two wavelet waveforms have a certain difference, and are not completely overlapped, and the direct subtraction cannot obtain a correct processing result, so that the matched filtering processing is required. Fig. 25 shows the waveform of the primary effective wave left after the adaptive matching subtraction algorithm, and it can be seen that the multiples can be well suppressed after the matching subtraction. The predicted multiples waveform is shown in FIG. 26, which has a waveform with good similarity to the streamer true multiples waveform.

5 the invention is further described below in connection with the application of the effect analysis.

And establishing a proper calculation model based on the multi-wave prediction subtraction method for inspection. The seismic model adopts a rapid F-K phase shift method to carry out wave field continuation, and realizes matched filtering based on single-channel self-adaptive minimum two-multiplication. In the wave field continuation process, due to the fact that the arrangement modes of the horizontal cables and the vertical cables are different, the space positions of the received seismic wave fields are different, and therefore the seismic wave fields of the horizontal cables and the vertical cables need to be calibrated to a uniform reference surface in a wave field continuation mode, and then matching pressing is conducted. The selection of the correction reference surface can select any plane in the water layer, but in order to facilitate processing and reduce extension errors as much as possible, the sea bottom surface is selected as the extension reference surface, and the matching pressing multiple process is uniformly calibrated to the sea bottom surface.

Firstly, extending the seismic record received by a horizontal cable to the seabed in the forward direction of the propagation direction of seismic waves; secondly, performing wave field separation on the vertical cable seismic record, wherein the separated down waves are just multiple reflected waves received after being reflected by the sea surface; thirdly, wave field forward continuation is carried out on the separated downlink waves along different paths based on the positions of the detectors with different depths in an F-K domain, the separated downlink waves are extended to a seabed datum plane in a unified mode to be in one-to-one correspondence with the horizontal cable wave field, and a predicted multiple wave field is obtained; and the fourth step is to match the horizontal cable wave field with the vertical cable multiple wave field based on a single-channel least square matched filtering method, so as to realize the multiple information suppression in the horizontal cable wave field, remove the multiple and keep the effective wave field information. Therefore, after the series of processing, the result of removing the multiple is finally obtained, but the result is the result of carrying out matching subtraction on the sea bottom, so that the reverse continuation is finally carried out to the sea surface.

The invention verifies the application effect of the single-channel matching adaptive subtraction method for jointly suppressing multiples under a stereo observation system through different data models. The method comprises the steps of firstly predicting a multiple wave model by using vertical cable data, realizing data matching correction of a vertical cable and a horizontal cable, then performing suppression attenuation on multiple waves by using a single-channel matching adaptive subtraction method, and comparing the advantages and the disadvantages of processing effects under different models. The model adopts relatively ideal horizontal layered stratum media, the uppermost layer is a seawater layer, the lower layer is a layered horizontal stratum, the total of five layers of media are provided, and the sea level is a good free reflection interface; the model II is a stratum-inclined marine geological model, the seismic wave speed of the stratum model has slow horizontal transverse change, and the reliability of the combined pressing method is verified when the transverse speed change is not large; the model III is a relatively complex rugged seabed stratum model, is relatively consistent with seabed stratum conditions in a real environment, and can well verify the rationality and applicability of the multiple-wave suppression attenuation method; finally, aiming at a special model with large fluctuation and severe transverse speed change under a specific condition, a multi-cable combined observation and data merging strategy is provided, and the purposes of improving the suppression effect and enhancing the signal to noise ratio are achieved by reasonably changing the arrangement method of an observation system, improving the illumination range and the data precision of a stereo observation system and realizing the purposes of improving the suppression effect and enhancing the signal to noise ratio. Model trial calculation shows that the method for jointly suppressing multiples can achieve good effects on different seabed stratum models, and the adaptability and the effectiveness of the method in the field of offshore multiple suppression are fully demonstrated.

5.1 horizontal stratum seabed model

The model is composed of five layers of horizontal layered media, the uppermost layer is a seawater medium, and the finite difference method of an acoustic wave equation is adopted for synthesizing the seismic record. The method for processing the absorption boundary and the reflection boundary of the model is adopted respectively, the sea level is used as a good free reflection interface with a reflection coefficient of-1, shot gather records containing free surface multiples can be simulated by a numerical simulation algorithm, and seismic wave fields received by a horizontal streamer and a vertical cable are recorded simultaneously by a three-dimensional system.

In forward modeling of different geological models, the following three-dimensional observation systems are adopted, Ricker wavelets are used for seismic wavelets, the number of seismic record tracks is set to 512 tracks, the distance between 5m tracks is adopted for a towline, the distance between 1m tracks is adopted for a vertical cable, the number of sampling points is 2048, the sampling interval is 0.002s, and the depth of seawater is set to 500 m. The source firing position is 1250m, and the midpoint is shot. After the seismic source is activated, the receivers on the streamer and the vertical cable simultaneously receive seismic wavefield information.

FIG. 27 shows a horizontal layer model created by forward modeling, in which the depth of the sea water layer is 500m and the seismic velocity of the sea water layer is 1500 m/s. The seismic wave velocity of the seabed stratum model is sequentially increased from top to bottom as marked in the figure. The sea surface is set as a free surface, the seismic source is located in the middle of the horizontal cable on the sea surface, and the seismic forward simulation recordings received by the streamer and the vertical cable are respectively shown in fig. 28(a) and (b).

Fig. 28(a) is a seismic record received by a streamer, and it can be seen that wave field components on a wave field seismic section are very complex due to interference of free surface multiples, multiple hyperbolic co-phase axes are mixed with each other and are not easily distinguished, and great difficulty is brought to subsequent processing of effective wave reflection waves, so effective multiple separation must be performed; fig. 28(b) is a vertical cable seismic record, and since the up and down traveling waves are both straight lines in the case of zero offset, the up and down traveling waves can be separated by a conventional F-K filtering method. The later wave field continuation needs to utilize the down-going wave, so that the traveling wave is reserved for multiple prediction.

FIG. 29(a) shows the seismic wavefield after the continuation of the streamer data to the seafloor, and FIG. 29(b) shows the predicted multiples wavefield for the vertical cable data. Aiming at the difference of double-cable acquisition modes of a stereo observation system, a wave field continuation method is needed to be adopted to respectively extend seismic wave field information received by detectors at different depths to corresponding positions of detectors of a horizontal cable, so that reasonable matching of the seismic wave field information received by the detectors of the horizontal cable is realized. In the invention, free surface multiples, i.e. ghost wavefields that propagate to the sea floor again after reflection at the sea level, need to be processed and suppressed. This appears as a down-going wavefield above the vertical cable receive record. Therefore, the wave field continuation is carried out by utilizing the downward wave field of the vertical cable, the seismic wave fields received by the geophones at different depths are respectively pushed out to the sea bottom surface, meanwhile, the seismic record containing the free surface multi-time wave data received by the horizontal cable is also prolonged to the sea bottom surface, and the two data are matched. Because the data of the vertical cable continuation is multi-time wave data reflected back by the sea level, namely ghost reflection data, the prediction work of the multi-time waves between the sea surface and the seabed is completed in the wave field continuation process.

Fig. 30 is a least squares matched filter graph provided by an embodiment of the present invention. The predicted multiples model after wave field extension is not completely matched with the multiples in the actual seismic record, that is, there is a great difference in the phases and amplitudes of the multiples and the jump time of the waveform, and the more perfect the matching between the predicted multiples and the actual multiples is, the better the matching is to subtract the multiples from the original seismic data. The method adopts a single-channel adaptive matching algorithm based on a least square method, and effectively suppresses the multiple waves in the seismic data after filtering processing, thereby obtaining good effect.

5.2 inclined stratum seabed model

Fig. 31 is a tilt model established by forward modeling. Because the stratum inclines, the seismic wave velocity has slow change in the transverse direction, so that in the course of wave field continuation the seismic wave homophase axis and horizontal stratum have different difference, and the amplitude phase position of seismic wave field also has obvious difference. Based on the matching suppression method, a good suppression effect can still be obtained, and the method for jointly suppressing multiple waves under the stereo observation system is feasible.

FIG. 32(a) is a seismic recording received by a streamer and FIG. 32(b) is a seismic recording received by a vertical cable. Under the condition of stratum inclination, the same-phase axis is not a strict hyperbolic curve any more, but the multiples in the red dotted line frame still appear periodically compared with the effective wave, so that the method can be continuously adopted for pressing treatment by utilizing the periodic difference of the multiples. FIG. 33 shows the results after multiple-pass compression, and it can be seen that the combined matched compression still achieves good results in the case of formation dip.

5.3 rugged stratum seabed model

Fig. 34 shows a rugged seabed model, the vertical distribution of the seismic wave velocity is similar to that of a horizontal lamellar model, but due to the fluctuation influence of the actual seabed terrain, the real seabed model is not always horizontal, but has a certain fluctuation change. At this time, the sea bottom reflected wave no longer passes through the horizontal interface, the seismic record receives the standard hyperbola, the distribution of the same phase axis of the seismic record is influenced to a certain extent under the condition that the sea bottom is rugged, and the characteristics of the multiple waves are also different. After the wave field continuation is done, the multiple wave field in the horizontal cable seismic record is shown in the dashed frame of fig. 35, based on the multiple wave suppression method of the present invention, after data processing, it is obvious from the figure that still a good suppression effect can be obtained, the multiple waves in the dashed frame are suppressed clean basically, and only the primary reflected waves corresponding to different reflection interfaces with strong effective energy are left. The method shows that the combined suppression effect of multiple waves under the stereoscopic observation system can be achieved under the conditions that the transverse fluctuation of the seabed stratum is not large and the transverse speed change of seismic waves is not severe.

The present invention will be further described with reference to effects.

The invention is based on the wave field continuation theory, realizes the wave field extraction of vertical cable seismic wave data and the continuation correction of horizontal cable multiple wave data under an offshore three-dimensional observation system, and adopts a self-adaptive matching subtraction method to suppress multiple waves. Through the processing and inspection of different model data, a better effect can be obtained, and the following conclusion can be obtained:

(1) in marine seismic exploration, multiples are effectively removed as interference waves. Aiming at the limitation of the conventional horizontal towline multiple suppression method, the invention fully considers the wave field characteristics of multiple under a stereo observation system, realizes the towline-vertical cable combined suppression method, and can obtain better suppression effect for different geological models.

(2) The invention extrapolates the seismic wave field based on the F-K domain wave field continuation theory. Aiming at the problems of continuation noise and boundary reflection in the extrapolation process, the invention combines the F-K domain filtering and pressing method to improve the wave field continuation attenuation operator so as to greatly reduce the continuation noise.

The above description is only for the purpose of illustrating the preferred embodiments of the present invention and is not to be construed as limiting the invention, and any modifications, equivalents and improvements made within the spirit and principle of the present invention are intended to be included within the scope of the present invention.

Claims

1. A processing method of a single cable for pressing multiple waves by combining an offshore stereo observation system is characterized by comprising the following steps:

fourthly, matching a horizontal cable wave field with a vertical cable multiple wave field based on a single-channel least square matched filtering method, suppressing multiple information in the horizontal cable wave field, removing the multiple, reserving effective wave field information to obtain a result of removing the multiple, and performing reverse continuation on the obtained result of removing the multiple to return to the sea surface;

the F-K domain wave field separation method in the third step specifically comprises the following steps:

firstly, inputting original vertical cable seismic data, carrying out F-K conversion on the original vertical cable seismic data, displaying the original vertical cable seismic data as wave field energy information of different interval positions in a frequency-wave number domain, reserving an upgoing wave field, removing the wave field information of a downgoing wave in the frequency-wave number domain, carrying out inverse F-K conversion on the remaining wave field information, returning to a t-x domain, and obtaining a reflected wave field record which is separated into the downgoing wave and only reserves the upgoing wave;

predicting multiple waves by performing wave field continuation on the basis of the downlink waves on the vertical cable, and performing matched pressing;

the method for wave field continuation comprises the following steps:

for computing forms

The method comprises the following steps that (1) r in a formula represents the relative length of a seismic source position and an observation position, a one-way wave equation is utilized to carry out wave field continuation, only wave fields advancing in a single direction are used, an uplink one-way wave is propagated upwards, a downlink one-way wave is propagated downwards, a negative sign in the formula represents the downlink one-way wave, and a positive sign in the formula represents the uplink one-way wave; the three-dimensional wave equation is as follows:

；

carrying out Fourier transformation to obtain:

；

during wave field continuation, wave field extrapolation is carried out on the uplink wave and the downlink wave, and the wave field extrapolation is carried out from the position of a seismic source in a forward continuation mode along the forward propagation mode of the actual seismic wave field to realize a forward simulation process;

the Fourier transformation method specifically comprises the following steps:

；

；

solving the wave equation in the form of Fourier domain:

；

written in the form of

；

The upper type

Representing the wave dispersion form of the one-way wave equation, taking the positive on the right of the equation as a downstream wave, and taking the negative on the right of the equation as an upstream wave;

；

；

the wave dispersion form is changed into a one-way wave equation:

；

(1) first, the earth surface wave field

Conversion into frequency-wavenumber domain

；

(2) Calculating an extrapolation factor

；

(3) In the frequency-wavenumber domain, the seismic wave field is multiplied by a factor C to obtain

A wave field of (c);

(4) repeating the steps (2) and (3) to gradually obtain all underground depths

A frequency domain wavefield of; (ii) a Wherein i =1,2,3 …;

(5) wave field of frequency-wavenumber domain

Inverse transform back to the time-space domain

Finally, obtaining time-space domain seismic wave fields at different depths;

the method for F-K domain phase-shifted wavefield extrapolation further comprises: firstly, separating an upgoing wave and a downgoing wave from a total wave field of a vertical cable, and then respectively carrying out variable-depth wave field extrapolation of different paths and different distances on wave fields received by detectors at different depths based on a wave field continuation operator to obtain an extrapolated seismic wave field value;

2. The offshore stereo observation system joint suppressed multiple single-cable processing method of claim 1, wherein in the fourth step, matching the horizontal cable wave field with the vertical cable multiple wave field based on a single-pass least squares matched filtering method, for suppressing multiple information in the horizontal cable wave field, the least squares matched filtering method comprises:

；

converting the above formula into a time domain, and simultaneously limiting the length of a filter factor A; then

The above equation is transformed into:

；

the time domain filtering factor a (t) in the formula is obtained, and the effective wave data after the multiple wave is eliminated is obtained

The energy of (2) is reduced to a minimum.

3. The offshore stereo observation system combined pressing multiple single cable processing system of the offshore stereo observation system combined pressing multiple single cable processing method according to claim 1.

4. An information data processing terminal for realizing the processing method of the combined suppression multiple single cable of the marine stereoscopic observation system according to any one of claims 1-2.

5. A computer-readable storage medium comprising instructions that, when executed on a computer, cause the computer to perform the method for joint multiple press single cable processing of a marine stereo vision system as claimed in any one of claims 1-2.