GB2524487A

GB2524487A - Removal of sea surface effects from seismic data

Info

Publication number: GB2524487A
Application number: GB1405182.5A
Authority: GB
Inventors: Lasse Amundsen; Johan Olof Anders Robertsson
Original assignee: Statoil Petroleum ASA
Current assignee: Equinor Energy AS
Priority date: 2014-03-24
Filing date: 2014-03-24
Publication date: 2015-09-30
Anticipated expiration: 2034-03-24
Also published as: WO2015144453A1; GB2524487B; GB201405182D0

Abstract

A method of removing sea surface effects from a gather of seismic data including the steps of: injecting a wavefield along an acquisition datum on a velocity model using a two-way wavefield propagation method to produce a modelled wavefield; and predicting multiples from the modelled wavefield. The step of predicting the multiples may comprise defining a velocity model without a free surface at the top and including in the model at least one reflecting interface associated with free-surface related multiples to be removed. The method may further comprise removing free-surface related multiples. Injecting the wavefield refers to the process of introducing a wavefield recorded in the field or in a previous simulation in a finite-difference grid.

Description

Removal of sea surface effects from seismic data

FIELD OF THE INVENTION

The invention relates to the removal of sea surface effects from seismic data.

BACKGROUND OF THE INVENTION

Historically, seismic data processing methods considered only primary reflections, i.e., energy which is emitted from a source, which reflects at a single interface in the sub-surface, and which is recorded at a receiver when reaching the earth's surface again.

In reality, waves follow many other different propagation paths, and sea-surface ghosts have a profound impact on marine seismic data. The source ghost corresponds to energy that is radiated upward from the source and immediately reflects downwards to interfere with the emitted wavetield that is directly radiated downwards. Similarly, the receiver ghost corresponds to sub-surface reflections that reflect from the sea surface downwards interfering with the recorded up-going wavefield. The main effect of the source-and receiver-ghosts is to severely reduce the bandwidth of the seismic data (Amundsen, 1993). However, in addition, if the sea-surface is rough, the ghosts will introduce statics shifts on the source-and receiver side as well as noise due to scattering of the seismic wavefield from the sea-surface (Laws and Kragh, 2002).

Noise due to rough-sea surface scattering can be of particular importance in time-lapse seismic applications (Laws and Kragh, 2002).

Multiples correspond to energy that reflects several times in opposite directions before being recorded at the surface and are classified according to how they reflect at interfaces before arriving at the receiver as illustrated in Figure 1. An important category of multiples in marine seismics is surface-related multiples (SRM5). These are multiples that reflect from the sea-surface at least once before being recorded (excluding the source and receiver ghosts). For surface-related multiples, the sea surface is referred to as the multiple-generating horizon. In Figure 1 we illustrate three different types of surface-related multiples which can be of different orders (first-, second-and third-order, etc.) depending on how many downwards reflections from the sea-surface they undergo.

An interbed multiple, on the other hand, is a multiple that has a sub-surface interface as the multiple-generating horizon instead of the sea surface (Figure 1). Multiples can of course also constitute both parts that are surface related and parts that are of interbed character.

A special case of a surface-related mutiples that is olten present in marine seismic data is so-called water-layer peglegs, which correspond to energy that reverberate several times within the water column either on the source side (source-side peglegs) or the receiver side (receiver-side peglegs).

The effect of multiples in data processing varies depending on the severity of the multiples and the sub-surface complexty. Sometimes the multiples can be clearly recognized as a time-delayed "copy" of an earlier primary. Other times, the multiple can be very difficult to recognize and lead to interpretation of false structure in the sub-surface. So-called scattered multiples are particularly difficult to deal with. These are multiples that arrive at the receiver (or emerge from the source) with a considerable cross sail-line propagation component. Since spatial sampling in the cross-line direction often is poor, multiple attenuation methods tend to be somewhat inefficient in handling these types of multiples. Finally, improperly attenuated multiples can also severely affect the signal-to-noise ratio after imaging.

Multiple removal methods can broadly be characterized as belonging in one of three different categories, namely methods based on separability, periodicity or methods based on modelling.

A common multiple attenuation method based on separability involves transformation of the data using the so-called parabolic radon transform (Sheriff and Geldart, 1995).

The transform maps data from offset-time domain into a domain of parabolic curvature vs. intercept time of events at zero offset. Since multiples tend to have a greater degree of curvature compared to primary reflections at similar arrival times (due to the fact that the multiples usually have propagated in shallower parts of the sub-surface where propagation velocities are lower), they will often separate well from primary reflections after the parabolic radon-transform. The effect can also be enhanced by applying a normal-rnoveout correction with primary propagation velocities before the radon transform. After surgically isolating multiples from signal, the isolated multiples can be inversely transformed back to the offset-time domain and adaptively subtracted from the original data.

Multiple attenuation methods based on periodicity involve the linear radon-transform where data in the offset-time domain are mapped into the slowness-zero-offset-intercept-time domain (often referred as the tau-p domain; Sheriff and Geldart, 1995).

After such a transformation, multiples can be removed through predictive deconvolution of data traces corresponding to different slowness (or) values.

In multiple attenuation methods based on modelling, a model of the multiples is first produced before (adaptively) subtracting it from the data. The most common and highly successful method for surface-related multiple attenuation is referred to as SRME (Surface Related Multiple Attenuation) and was originally developed by researchers at Technical University of Delft in the Netherlands (Verschuur et al., 1992; Dragoset and Jericevic, 1998). The method is based on the observation that any surface-related multiple can be regarded as being composed of separate legs corresponding to individual primary reflections. Through a spatial convolution of the recorded primaries with themselves (over both source and receiver coordinates), a model of surface-related multiples can be produced without any knowledge of the sub-surface. The method can be highly effective but relies on adequate sampling of data on the surface, particularly for short source-receiver offsets. In its full 3D formulation the method can become highly computationally intensive due to the large number of convolutions and interpolation between different source/receiver locations that is required.

Another multiple method based on modelling is the wave-equation modelling (WEM) method based on the so-called one-way wave equation where multiples are modelled using a synthetic model of background velocities and reflectors (Stork et al., 2006).

Amundsen (2001) introduced yet another multiple elimination method that relies on decomposition of the recorded wavefield into its up-and down-going constituents. The down-going wavefield can be deconvolved from the up-going wavefield removing both surface-related multiples as well as the source signature. Decomposition of the recorded data is closely analogous to deghosting the data and will typically require multicomponent recordings. Until recently, this was only available for seabed recordings where the method by Amundsen (2001) has been applied with some success (Brunelliere et al., 2004).

In seabed seismic data where multicomponent recordings are commonly made, it is possible to decompose the recorded data into up-and down-going components both as they would be perceived just above the seafloor (in the water) or just below inside the seabed (Amundsen, 1993). If the data are decomposed into their up-and down-going components in the seabed (below the water layer), the receiver-side water layer peglegs will automatically be excluded from the up-going part of the wavefield. This has often proved to be a very effective means for removing significant multiple energy in seabed recordings (Stewart et al., 2007).

The level of the sea surface is affected by for instance tides altering the length of the propagation leg of surface-related multiples in the water. This aspect was exploited by Calvert (2005), Robertsson and Kostov (2010) and Ikelle (2004) to develop multiple attenuation methods using seismic data acquired in two different tidal states.

Finally, multiples do not necessarily have to be regarded as noise. Provided that the multiples can be processed properly they convey information of sub-surface structures and can even help to image structures that have not been illuminated by primary reflections (e.g., when undershooting surface obstructions such as oil platforms or in the vicinity of complex salt structures). This aspect has been highlighted recently in a similar context to simultaneous source acquisition where it has been proposed to fundamentally change processing of seismic data moving away from the view of acquisition and separation of sequential shot records. Instead an alternative approach would be to consider illumination of the sub-surface achieved by any means of acoustic energy radiated from different locations on the surface be it from simultaneous sources or surface-related multiple reflections (Muijs et al., 2007a; Muijs et al., 2007b; Berkhout and Verschuur, 2012; Wapenaar et al., 2012).

Wavefield injection

By wavefield injection we refer to the process of introducing a wavefield recorded in the field or in a previous simulation on a finite-difference (FD) grid. The wavefield needs to contain both pressure and particle velocity recordings and should be injected along a closed surface in the finite-difference solution. Of key importance is the fact that wavefield can be reconstructed perfectly preserving its directionality inside the injection surface. Several authors have described how this can be done using for instance the Kirchhoff integral (Amundsen et al., 1993, Mittet, 1994) or using the so-called FD-injection technique (Robertsson and Chapman, 2000). Because of its conceptual simplicity we will be describing the wavefield injection process using the FD-injection method for a second-order accurate acoustic staggered finite-difference scheme (Virieux, 1984). However, the method applies equally well to other methods of injection or higher-order finite-difference methods.

Let us consider a 2D staggered FD grid in the vicinity of a closed injection surface S The region outside the surface is referred to as the external region Y whereas the region inside the surface is referred to as the internal region. The process of introducing the source wavefield is straightforward and involves manipulating the update of the wavedfield only in the vicinity of S for the update of points where the spatial extent of the FD-stencil intersects S (where the wavefield in the grid is discontnuous).

The wavefield to be injected is known (either from recorded data or from a previous simulaton) in the vicinity of the surface S for all times. The wavefield wil interact with the model throughout the FD grid within Y and and generate a scattered secondary wavefield. In both the injected wavefield and a scattered secondary wavefield is present whereas in only the secondary wavefield is present. Since the injected wavefield is known this can be added and subtracted as appropriate for the update of points along the injection surface for the parts of the spatial stencils that intersect (subtract the source wavefield at the appropriate points when updating points in the immediate exterior of the injection surface and adding the injected wavefield at the appropriate points when updating points in the immediate interior of the injection surface).

In each FD time step the pressure is first updated in the entire grid. When the update is complete, we go back and correct the update at the points where the spatial FD stencil intersected the surface S Inside S in Y, the wavefield is updated as if the injected wavefield was propagating through the entire grid. Therefore we must add the injected wavefield to the particle velocity components corresponding to the parts of the stencil that are outside S in. For a second-order accurate scheme, we therefore need to know the injected wavefield along the closest grid point outside S. Outside S in the wavefield is updated as if no wavefield was injected. Therefore we must subtract the injected wavefield from the particle velocity components corresponding to the parts of the stencil that are inside S in Y. Next, we advance the calculation by half a FD time step and update the particle velocities in the entire grid. The wavefield injection is performed using the same procedure by adding and subtracting the pressure of the injection wavefield at the grid point around. In total, we therefore need to know the values at all times of pressure and particle velocities along a single layer of grid-points around S, staggered appropriately both in time and space. By iterating these two steps of the update, the entire FD-simulation is stepped through and the wavefield is injected along the surface

S

Robertsson and Chapman (2000) described how to use FD-injection to recompute the response of finite-difference solutions after local model alterations. The method will help us to understand the mechanics of the injection process that will be useful to describe our methods for imaging and multiple elimination.

In the method by Robertsson and Chapman (2000), an FD-simulation is initially carried out on an unperturbed model. During this simulation the wavefield is recorded along a closed surface around a target of interest. In a second step, the model is altered inside the closed surface. Next we wish recompute the wavefield after the model alteration.

This can be done by injecting the wavefield from the first simulation on the same surface but now on the perturbed model. Part of the wavefield will now leak through the injection surface. This wavefield corresponds to the difference in wavefields before and after the perturbation and is exact as long as long-range interactions between the perturbed and unperturbed model are accounted for. Robertsson and Chapman (2000) proposed to truncate the grid of the perturbed model in order to make substantial computational savings. This can be done at the expense of neglecting second-and higher-order interactions of the perturbed wavefield with the unperturbed model that propagate back inside the region of change again.

From the description in the previous paragraph it should be clear that the FD-injection method will result in an apparent discontinutity in wavefields outside and inside the injection surface. Through the manipulation of wavefield quantities during the update of the FD-simulation, this discontinuity is not experienced by the two superimposed wavefields: the injected wavefield and the secondary wavefield. Note that if we would have injected the wavefield on an unperturbed model, inside the injection surface we would have observed a complex wavefield exactly corresponding to the wavefield in the first simulation in that part of the model, whereas outside the injection surface, the solution would have been zero (only the secondary wavefield is present there which in such a case is zero).

The wavefield injection method injects wavefields constituents in the ingoing and outgoing directions. The wavefield injected in the ingoing direction will recreate the injected wavefield in its interior. It will constitute of direct waves as well as multiple interactions between various parts inside and outside the injection surface in the unperturbed model. The wavefield injected in the outgoing direction is injected to destructively interfere with the previousy injected wavefield propagating through the injection surface. In the unperturbed model, the destructive interference will be perfect.

In a perturbed model, however, the secondary wavefield will leak through.

SUMMARY OF THE INVENTION

The invention provides methods, apparatus, and a computer-readable medium as set out in the accompanying claims.

Embodments of the invention will now be described, by way of example only, with reference to the accompanying drawings.

BRIEF DESCRIPTION OF THE FIGURES

Figure 1 illustrates different categories of multiples; Figure 2 shows a model used in the wavefield injection of data recorded along a horizontal surface; and Figure 3 shows a computing device suitable for carrying out methods described below.

DESCRIPTION OF PREFERRED EMBODIMENTS

Each feature disclosed or illustrated in the present specification may be incorporated in the invention, whether alone or in any appropriate combination with any other feature disclosed or illustrated herein.

We describe new methods for the prediction and removal of sea surface effects from marine seismic data. This includes the removal of source-and receiver-side ghosts as well as surface-related multiples. In one of the embodiments, the method requires a velocity model and a reflectivity model similar to what is required in WesternOeco's WEMM (Wave Equation Migration Modelling) multiple elimination method (Stork et al., 2006). In another embodiment we employ the method of exact boundary conditions (van Manen et al., 2007; Vasmel et al., 2013; Vasmel et a., 2014) to predict multiples

from wavefield injected data.

We describe a method which can be referred to as Full Wave Wave Equation Migration Modelling, or FW-WEMM, and which differs from WEMM in that it uses the so-called FD-injection method and a two-way wave equation solver. It also relies on the acquisition of multicomponent data recordings and is therefore particularly applicable to towed marine multicomponent streamer data, seabed cable or node data, or VSP data.

It can be generalized to work on data such as slanted cable or over/under geometries that enable wavefield decomposition. It can also be generalized to conventional towed marine seismic recordings using hydrophones only provided that the data can be separated into its up-and down-going constituents.

Figure 2 shows the wavefield injection of data recorded along a horizontal surface 10 where the source 8 was located below the surface 10. A multiple-generating horizon (ie reflector) 24 has been included in the otherwise smooth velocity model.

Wavefield injection simulating marine seismic experiments with and without a free surface Let us consider a marine recording geometry where a source B is excited below a streamer 10 where pressure and particle velocities are recorded.

First, let us consider what happens if we inject the recorded data in a model with the free surface (ie sea surface) in its correct location. This is analogous to the experiment carried out by Mittet (2002). In practice, a marine recording geometry does not correspond to a closed surface but as we shall see that is not relevant for our application. We will however be subject to edge ditfractions at the front and at the tail of the streamer but will ignore those effects for now. The first question that arises is what is inside and outside the injection surface in our case. We can either choose to reconstruct the wavefield above the streamer or below the streamer. Note that if the model exactly corresponds to the real world where the data were recorded we will observe that the wavefield is zero on ether side of the streamer depending on what direction we choose to reconstruct the wavefield in. Mittet's (2002) point was that injection cannot be used to image both primaries and multiples in Reverse Time Migration (RTM). This statement is incorrect but we need to perturb the model first.

The perturbation that we will introduce is to remove the sea surface and replace it with an absorbing boundary 12. Absorbing boundaries are also used on the other three sides, as shown in Figure 2. We will regard the interior to the injection surface 10 to be above the streamer 10 where both injected and secondary wavefield will be reconstructed. Below the streamer 10 only the secondary wavefield will be reconstructed. The secondary wavefield corresponds to the part of the wavefield that should have had an interaction with the free surface that now is missing. Let us consider the impact on a few different arrivals of interest, which are illustrated in Figure 2.

a) The direct wave propagating downwards (not shown) into the sub-surface from the source 8. This has not crossed the injection boundary 10 and is not generated either above or below the injection surface 10.

b) The direct wave 16 propagating upwards towards the sea surface. This will be injected above the injection surface 10 and will be absorbed by the absorbing boundary 12.

c) The downward-propagating ghost 18 of the upward-propagating direct wave 16 will be injected downwards with opposite polarity to destructively interfere with its downwardly propagating counterpart (that would have been reflected from the upwardly propagating direct wave reflected at the sea surface). However, since we have replaced the sea surface by an absorbing boundary 12 it will propagate below the injection surface 10 with opposite polarity.

d) Primary events (ie primary reflections) 19 will be injected above the injection surface only and will be absorbed by the absorbing boundary 12.

e) Ghosts 20 of primary reflection events 19 will be injected downwards with opposite polarity to destructively interfere with their downwardly propagating counterparts (generated from upwardly propagating energy that would have reflected off the sea surface if it would have been included in the model). However, since we have replaced the sea surface by an absorbing boundary 12 they will propagate below the injection surface 10 with opposite polarity.

f) Surface-related multiples 22 are generated from the ghosts 20 of the primary or lower-order surface-related multiple reflection events and are therefore present below the injection surface 10 with opposite polarity. The injection surface 10 will also inject surface-related multiples with the correct polarity upwards only. They will interfere above the injection surface 10 with the corresponding propagating events 22 and should in principle cancel even though the absorbing boundary 12 would compensate for anything remaining.

g) Ghosts 20 of surface-related multiples, just as for the ghosts of primary events, these will be injected downwards with opposite polarity to destructively interfere with their downwardly propagating counterparts (that would have been generated after reflections from the sea surface which is not included in the model). However, since we have replaced the sea surface by an absorbing boundary 12 they will propagate below the injection surface 10 with opposite polarity.

In summary, if we record the wavefield below the injection surface 10, the method will result in a recording of the source ghost 18, the receiver ghost and all surface-related multiples. Note that all ghosts (including those of surface-related multiples) are directly recorded and properly injected from the data. The last leg of the surface-related multiples have been modeled in the FD model. This is important since the error of modelling surface-related multiples of higher order does not accumulate. All that remains to be done is to adaptively subtract these data from the raw data resulting in a data set consisting of ghost-free primaries only.

Example embodiments

Figure 2 shows a preferred embodiment of our application of injection of a wavefield on a smooth velocity model without a sea surface for surface-related multiple (SRM) removal. As opposed to methods such as the multi-dimensional deconvolution based multiple suppression by Amundsen (2001), this preferred embodiment bears more resemblance with the wave-equation migration modelling method for multiple suppression by Stork et al. (2006) where individual SRM's are targeted and modelled.

Super-imposed on the velocity model shown in Figure 2 we have included a reflector 24 that generates SRM's. The reflector 24 can either be modeled by having two piece-wise smooth velocity models with a discontinuity at the reflector 24 or as a jump in density contrast so that reflections are generated when propagating a wavefield in the model.

If we inject the recorded wavefield on the model in Figure 2 we will observe the following. Below the receiver level 10 we will observe what usually is considered as noise in reflection seismic processing: the ghost 18 of the direct wave 16, all receiver ghosts of the data 20, the source ghost 30 of the primary reflection from a chosen reflector, and all SRM's associated with the chosen reflector 22 and the part of 20 not corresponding to the receiver ghost of the primaries. The wavefield observed above the receiver level 10 has been deghosted on the receiver side, has had the direct wave and the primary from the selected target source-side deghosted as well as all surface-related multiples associated with the target of interest removed.

In principle the wavefield modeled above the receiver level is the desired output.

However, modeling the reflection coefficient accurately from the target of interest is difficult in practice and we will likely need to rely on adaptive subtraction of a noise estimate from the data instead. For a practical workflow we propose to proceed as follows: 1. Remove the direct wave 16 from the data (using for instance a wavefield separation method for the water layer model bounded by a sea surface as described by Amundsen et al., 2014).

2. Forward propagate the data (by wavefield injection) using the model in Figure 2 outputting the wavefield immediately below the receiver level 10.

3. Model the receiver ghosts of the original data using a method for wavefield separation (for instance the finite-difference injection method with a model with water layer only without sea surface as described by Amundsen et al., 2014).

4. Subtract the receiver ghosts in step (3) from the noise estimate in step (2) which now only contains the up-going surface-related multiples.

5. Reghost the noise estimate for instance by using the method described by Amundsen et al. (2014) whereby the noise estimate is injected on a mode with a water layer only bounded by the sea surface on the top. Reghosting is the opposite of deghosting, where a ghost is added to its corresponding up-going event.

6. Adaptively subtract the surface-related multiples estimate from step (5) from the original data.

In an alternative preferred embodiment we can predict all surface related multiples in one step without the knowledge of the sub-surface or velocity model. The method resembles that of Vasmel et al. (2014). However, the main difference is that we use the method of wavefield injection in combination with the exact boundary conditions used by Vasmel et al. (2014) to model surface-related multiples only [Vasmel et al. (2014) directly remove the multiples without first predicting them similar to Amundsen (2001)]. An example of a preferred workflow is as follows: 1. Remove the direct wave 16 from the data (using for instance a wavefield separation method for the water layer model bounded by a sea surface as described by Amundsen et al., 2014).

2. Forward propagate the data (by wavefield injection) using a homogenous water layer model with absorbing boundaries all sides except for the bottom edge where exact boundary conditions are used as in Vasmel et al. (2014). These exact boundary conditions uses real recorded data as Green's functions. Note that in this embodiment a model such as the one in Figure 2 is not used. Instead we use a model of a water layer only with absorbing boundaries on all sides except for the bottom where we used exact boundary conditions.

4. Subtract the receiver ghosts in step (3) from the noise estimate in step (2) which now only contains (all) up-going surface-related multiples.

5. Reghost the noise estimate for instance by using the method described by Amundsen et al. (2014) whereby the noise estimate is injected on a mode with a water layer only bounded by the sea surface on the top. As noted above, reghosting is the opposite of deghosting, where a ghost is added to its corresponding up-going event.

In the section on wavefield injection we did not discuss the fact that the injection surface is staggered and that recorded data is not staggered in space. Also, to deal with higher order operators the method may use Kirchhoff extrapolation to predict the wavefield constituents at the required number of levels and staggered locations.

Instead of using FD-injection, it is possible to simply inject the wavefield along monopole and dipole sources. However, this may generate some slight noise in FD-modelling due to modelling of discontinuities near point sources.

The method has been described in the context of a towed multicomponent streamer.

The method can also be used with seabed cables and nodes, towed slanted streamers and over/under configurations where wavefield decomposition can also be achieved using hydrophone data only.

As described here, the source should be located below the streamer. This is not always practical and in any case, the wavefield in the vicinity of the source will likely have to be manipulated due to acquisition geometries and possibly clipping that might occur for near offsets in the recordings. An alternative is to simply mute the direct arrival n which case it does not matter if the source is located above or below the streamer. Muting will work particularly well in deep-water environments where the direct wave is well separated from the reflections. Note that if the direct wave is muted the source source-ghosts will remain on the primaries after applying the method.

In one of the embodiments, just like in WEMM, our FW-WEMM technique provides the ability to focus efforts and attention on a particular multiple generator. A complete model is not required but rather only the velocity model above the multiple generator(s) and the geometry of the multiple generator(s).

FW-WEMM is more expensive compared to WEMM since depth and time become the same dimension in WEMM. However, the cost of FW-WEMM is still small compared to Reverse Time Migration (RTM). RTM needs 3 FD simulations to produce and image whereas FW-WEMM needs one.

FigureS shows a computing device 60, which may for example be a personal computer (PC), on which methods described herein can be carried out. The computng device 60 comprises a display 62 for displaying information, a processor 64, a memory 68 and an input device 70 for allowing information to be input to the computing device. The input device 70 may for example include a connection to other computers or to computer readable media, and may also include a mouse or keyboard for allowing a user to enter information. These elements are connected by a bus 72 via which information is exchanged between the components.

It should be appreciated that any of the methods described herein may also include the step of acquiring seismic or electromagnetic data, which may then be processed in accordance with the method.

References Amundsen L (1993) Wavenumber-based filtering of marine point source data.

Geophysics 58: 1335-1348.

Amundsen, L., 2001, Elimination of free-surface related multiples without need of the source wavelet: Geophysics, 66, 327-341.

Amundsen, L., B. Arntsen and R. Mittet, 1993, Depth imaging of offset vertical seismic profile data. Geophysical Prospecting, 41, 1009-1031.

Amundsen, L., Robertsson, J. 0. A., Pedersen, 0., 2014, Wave equation processing and imaging of marine multicomponent data beyond traditional RTM: 76th EAGE Conference, Amsterdam, June 2014.

Berkhout AJ and Verschuur DJ (2012) Full wavefield migration -utilization of multiples in seismic migration. 74th EAGE Conference & Technical Exhibition, Extended Abstracts: B041.

Brunellière J, Caprioli P, Orion 5, Tilling D, Amundsen L and Aronsen H (2004) Surface Multiple Attenuation By Up-down Wavefeld Deconvolution: an OBC Case Study: 74th Annual International Meeting, SEG, Expanded Abstracts: 849-852.

Calvert R (2005) 4D technology: where are we, and where are we going? Geophysical Prospecting 53: 161-171.

Dragoset WH and Jericevic Z (1998) Some remarks on surface multiple attenuation.

Geophysics 63: 772-789.

Ikelle L (2004) The Phenomenon of Low And High Tides In the Demultiple Process.

74th Annual International Meeting, SEG, Expanded Abstracts: 1858-1 860.

Laws RM and Kragh E (2002) Rough seas and time-lapse seismic. Geophysical Prospecting 50: 195-208.

Mittet, R., 1994, Implementation of the Krchhoff integral for elastic waves in staggered-grid modelling schemes: Geophysics, 59, 1894-1901.

Mittet, R., 2002, Multiple suppression by prestack reverse time migration: A nail in the coffin: Expanded abstract at the 64th EAGE Annual Meeting.

Muijs R, Robertsson JOA and Holliger K (2007a) Prestack depth migration of primary and surface-related multiple reflections: Part I -Imaging. Geophysics 72: S59-S69.

Muijs R, Robertsson JOA and Holliger K (2007b) Prestack depth migration of primary and surface-related multiple reflections: Part II -Identification and removal of residual multiples. Geophysics 72: S71 -S76.

Robertsson, J. 0. A., and C. H. Chapman, 2000, An efficient method for calculating finite-difference seismograms after model alterations: Geophysics, 65, 907-918.

Robertsson JOA and Kostov C (2010) Multiple attenuation method. US Patent: 7,710,821.

Sheriff RE and Geldart LP (1995) Exploration Seismology. Cambridge University Press.

Stewart J, Shatilo A, Jing C, Rape T, Duren R, Lewallen K and Szurek G (2007) A comparison of streamer and OBC seismic data at Beryl Alpha field, U. K. North Sea.

Geophysics 72: B69-B80.

Stork, C., Kapoor, J., Zhao, W., Dragoset B., and K. Dingwall, 2006, Predicting and removing complex 3D surface multiples with WEM modeling -an alternative to 3D SRME for wide azimuth surveys?: Expanded abstract at the 76th SEQ Annual Meeting.

van Manen, D. J., J. 0. A. Robertsson and C. Curtis, 2007, Exact wave field simulation for finite-volume scattering problems: Journal of Acoustical Society of America, 122(4), EL115-E Li 21 Vasmel, M., J. 0. A., Robertsson, D. J. van Manen and A. Curtis, 2013, Immersive experimentation in a wave propagation laboratory: J. Acoust. Soc. Am., 134, EL492-EL498.

Vasmel, M., J. 0. A. Robertsson and L. Amundsen, 2014, A new solution to eliminate free surface related multiples inmulticomponent streamer recordings: 76th EAGE Conference, Submitted.

Verschuur DJ, Berkhout AJ and Wapenaar CPA (1992) Adaptive surface-related multiple elimination. Geophysics 57: 1166-1177.

Virireux, J., 1984, SH-wave propagation in heterogeneous media: Velocity-stress finite-difference method: Geophysics, 49, 1933-1 957.

Wapenaar CPA, van der Neut J and Thorbecke J (2012) Deblending by direct inversion. Geophysics 77: 9-12.

Claims

CLAIMS: 1. A method of removing sea surface effects from a gather of seismic data including the steps of: injecting a wavetield along an acquisition datum on a velocity model using a two-way wavefield propagation method, to produce a modelled wavefield; andpredicting multiples from the modelled wavefield.
2. The method of claim 1, wheren the step of predicting multiples from themodelled wavefield includes:defining a velocity model without a free surface at the top; and including in said velocity model at least one reflecting interface associated with free-surface related multiples to be removed.
3. The method of claim 2, where said velocity model is a kinematically good velocity model.
4. The method of claim 2 or 3, wherein said reflecting interface is modelled through a density contrast.
5. The method of claim 2 or 3, wherein said reflecting interface is modelled by having two piece-wise smooth velocity models with a discontinuity at the reflecting interface.
6. A method as claimed in any preceding claim, which further comprises removing free-surface related multiples.
7. The method of claim 1, wherein the step of predicting multiples from the modelledwavefield includes:defining a velocity model of a water layer without a free surface at the top; and using the method of exact boundary conditions at the bottom of the velocity model with Green's functions obtained from acquired seismic data.
8. The method of claim 1, wherein the step of predicting multiples from the modelledwavefield includes:convolving at least two different recorded data traces with each other.
9. The method of any preceding claim, wherein said seismic gather is a common shot gather.
10. The method of any one of claims 1 to 8, wherein said seismic gather is a common receiver gather.
ii. The method of any preceding claim, wherein said seismic data was produced at recording locations using a source located below the recording locations.
12. The method of any preceding claim, wherein said seismic data was produced using a muted direct wave.
13. The method of any preceding claim, wherein said two-way wavefield propagation method is a finite-difference method.
14. The method of any preceding claim, wherein said step of injecting a wavefield is carried out using the method of multiple point sources.
15. The method of any one of claims ito 13, wherein said step of injecting a wavefield is carried out using the FD-injection method.
16. The method of any of the prececEng claim, wherein the resulting wavefield is recorded at desired locations during propagation of the modelled wavefield.
17. The method of claim 16, wherein said step of recording the resulting wavefield at desired locations produces a recorded wavefield, and the method further includes adaptively subtracting the recorded wavefield from said seismic data to remove at least one surface-related multiple.
18. The method of claim 16, wherein said step of recording the resulting wavefield at desired locations produces a recorded wavefield, and the method further includes adaptively subtracting the recorded wavefield from said seismic data to remove receiver ghosts.
19. The method of any preceding claim, wherein said data were acquired using towed marine multicomponent seismic data.
20. The method of any one of claims 1 to 18, wherein said seismic data were acquired using multicomponent seabed seismic data.
21. The method of any one of claims 1 to 18, wherein said seismic data were acquired using towed marine over/under seismic data.
22. The method of any one of claims 1 to 18, wherein said seismic data were acquired using multicomponent seabed node data.
23. The method of any one of claims 1 to 18, wherein said seismic data were acquired using VSP data.
24. The method of any one of claims 1 to 18, wherein pressure and vertical component of particle velocity in said seismic data were acquired or estimated.
25. The method of any one of claims 1 to 18, wherein said seismic data were acquired using towed marine hydrophone seismic data where the recorded data can be decomposed into its up-and down-going constituents.
26. An apparatus comprising at least a processor, a memory, and an input device, wherein said apparatus is programmed to carry out the method of any preceding claim.
27. A computer-readable medium containing computer-readable instructions for performing a method as claimed in any of claims 1 to 25.