Wavefield Decomposition for Cross-Line Survey

Abstract

It is described a method for spatially filtering 3D seabed seismic data acquired in a cross-line geometry including methods of how to decompose the data using a filter with two spatial components.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1A illustrates an ocean bottom acquisition of a 3D seismic wavefield using an OBC and a source towed by a seismic vessel;

FIG. 1B illustrates the spread of data point after shooting a single shot line of the acquisition of FIG. 1A;

FIG. 1C illustrates the acquisition of 3D seismic wavefield using streamers towed by a first vessel and a source towed by a second seismic vessel;

FIG. 2 is a diagram illustrating steps in accordance with an example of the invention; and

FIG. 3 compares the performance of filter approximations s in accordance with examples of the invention.

EXAMPLES

In FIG. 1A there is illustrated an example of a seismic survey in cross-line geometry. The survey is a marine seismic survey in a body of water 10 between a seafloor 101 and a sea surface 102. A receiver cable 11 with a plurality of receivers 111 is laid out on the seafloor 101. The receivers 111 are preferably 4C sensors, though, as will be apparent from the following description, the sensors may comprise a 1C geophone and a hydrophone only, being thus capable of recording at least the vertical component of velocity and pressure at seafloor level.

A seismic vessel 12 tows a marine seismic source 13 close to the sea surface 102. The airgun 13 emits at precisely determined time intervals an impulse of acoustic energy referred to as “shot”. A dashed line 132 indicates the path of the towed airgun 13. The projection 133 of the dashed line 132 onto the seafloor 102 intersects the receiver line 112 at approximately 90 degrees. Though it is preferable to aim for a near-orthogonal orientation of receiver lines to shot lines, deviations are inevitable under real survey conditions. To facilitate the following description the receiver line or in-line direction is denoted as x direction, the shot-line or cross-line direction is marked as y direction and the vertical direction is taken as the z direction.

During a survey, the sources 131 are fired at intervals and the receivers 121 “listen” within a frequency and time window for acoustic signals such as reflected and/or refracted signals that are caused features in path of the emitted wavefield. After shooting a line, the vessel performs a u-turn in order to shoot a subsequent line usually with an offset in receiver line or x direction.

In the general practice, it is assumed that Green's functions which describe the wave propagation between source and receiver points are invariant for translation of the source and receiver in the cross-line direction. Hence, an offset between shooting lines can be regarded as an equal shift of the receiver line. As a result, data points obtained from a single source cross-line form a carpet on the seafloor, which is illustrated in FIG. 1B. In FIG. 1B, the triangles 122 denote the location of data points. As the wavefield is recorded in two spatial dimensions and in time, the resulting data are referred to as 3D wavefield.

In FIG. 1, there is shown a schematic cross-line marine survey with two vessels. A first vessel 15 tows five streamers 151 below the sea surface following a path 152. Simultaneously, a second vessel 16 on path 162 tows a seismic source 161 below the surface. As above in FIG. 1A, the resulting shot and receiver lines are essentially orthogonal to each other.

In the following description and the accompanying FIG. 2, steps are described leading to a decomposition of the 3D wavefield into up- and downgoing components.

After obtaining 20 the wavefield data set as acquired through the seismic receivers, the data are first preferably calibrated 21 to compensate for the differences between geophone and hydrophone recordings. Any suitable calibration may be used including for example the methods described in the International patent application PCT/GB03/04190. Following those methods, the calibration can be done using an in-line shot-line for the P, v_z and v_z components and using a cross-line shot-line for the v_y component.

Acoustic wavefield decomposition 22 is usually carried out on the pressure component P (involving spatial filtering of v_z). Instead in this example decomposition filters are applied to the vertical geophone component v_z (involving spatial filtering of P). The advantage of this example is that the spatial components of the filter only act on P as shown in PCT/GB03/04190.

Accordingly, acoustic wavefield decomposition into up- and down-going constituents above the seafloor can be achieved by solving the following equation:

$\begin{array}{cc}v\frac{+}{z}\left(f,k_x,k_y\right)=\frac{1}{2}a\left(f\right)v_z\left(f,k_x,k_y\right)\pm \frac{k_z\left(f,k_x,k_y\right)}{2\rho \omega }P\left(f,k_x,k_y\right)& \left[2\right]\end{array}$

In equation [2], a(f) denotes the optional frequency-dependent calibration filter that corrects for imperfections in the recording of the geophones, v_z⁻ denotes the up-going constituent of the vertical component of particle velocity, and v_z⁺ denotes the down-going constituent of the vertical component of particle velocity. The velocity v_z is the recorded or estimated vertical component of particle velocity in the frequency f—wavenumber domain, P(f, k_x, k_y) is the recorded pressure, and ρ is the density in the recording medium.

The term k_z, which can be expressed as

k _z(f, k_x, k_y)=√{square root over ((2πf/c)²−k_x²−k_y²)}  [3]

is the absolute value of the vertical wavenumber expressed in terms of horizontal wavenumbers in the in-line direction k_x and the cross-line direction k_y, and the velocity c of the recording medium. It should be noted that the decomposition could also be achieved by computing the up-going component of the recorded pressure P⁻ using equations [1], leading expression which include terms of 1/k_z. Such terms can be approximated using similar expansions as described below.

In known decomposition methods using any of the above equations [1, 2], the cross-line or y-directions is mostly ignored or a radial symmetry is assumed, with the vertical wavenumber then being computed using an approximation based exclusively on a one-dimensional direction, i.e. the in-line wavenumber k_x or the radial wavenumber k_r. When 3D effects of the sub-surface or acquisition geometry are significant, such approximations are no longer valid.

Equations as those described herein can be implemented in the common mid-point domain, which is proposed by Vermeer (1994) and Thomas (2000). It is however fruitful to rewrite or approximate equation [3] into a form constituting a cascade (sum or product) of one-dimensional (1D) spatial filters acting in the x- or y-directions only. This represents a computationally attractive way of filtering the data (both in terms of CPU and resorting data between different domains). One way to obtain filters of this form is to make suitable Taylor expansions of the horizontal wavenumbers in the square-root term around zero wavenumbers.

This approximation remains valid for data corresponding to propagating waves at k_x²+k_y²<(2πf/c)²).

The expression for the vertical wavenumber can be rewritten and expanded in k_x and ky to produce a few different alternative expansions that can be implemented using a cascade of filters that only act in the cross-line or in-line direction one at a time:

$\begin{array}{cc}{k}_z\left(f,k_x,k_y\right)\approx \sqrt{{\left(2\pi f/c\right)}^2-k_x^2}\left(\begin{array}{c}\begin{array}{c}1-\frac{k_y^2}{2\left({\left(2\pi f/c\right)}^2-k_x^2\right)}-\\ \frac{k_y^4}{8{\left({\left(2\pi f/c\right)}^2-k_x^2\right)}^2}+\end{array}\\ 0\left(k_x^6,k_y^6\right)\end{array}\right).& \left[4a\right]\\ \begin{array}{c}{k}_z\left(f,k_x,k_y\right)\approx 2\pi f/c\\ \left(1-\frac{k_x^2+k_y^2}{2{\left(2\pi f/c\right)}^2}-\frac{\begin{array}{c}{k}_x^4+k_y^4+\\ 2k_x^2+k_y^2\end{array}}{8{\left(2\pi f/c\right)}^4}+0\left(k_x^6,k_y^6\right)\right)\end{array}& \left[4b\right]\\ \begin{array}{c}{k}_z\left(f,k_x,k_y\right)\approx \sqrt{{\left(2\pi f/c\right)}^2-k_x^2}+\sqrt{{\left(2\pi f/c\right)}^2-k_y^2}-\\ 2\pi f/c\left(\begin{array}{c}2+\frac{2k_x^2k_y^2}{8{\left(2\pi f/c\right)}^4}+\\ 0\left(k_x^6,k_y^6\right)\end{array}\right)\end{array}& \left[4c\right]\end{array}$

Equations [4a-4c] represent different ways to proceed with an implementation of the filter with two spatial components that only rely on being able to filter the data along two perpendicular spatial directions one at a time which is exactly what can be achieved using the method described above in a cross-line geometry. Note that this does not mean that only a “cross” of mid-points are used in filtering the data. All cross-terms of multiplications of terms with different horizontal wavenumbers will result in a “virtual” carpet of data being used of dimension of the length of the spatial filters in both directions.

After the decomposition 22, multiples could be removed 23 from the data set. Further processing steps and/or filtering steps 24 could be performed on the decomposed data set. Using what is commonly referred to as imaging or migration 25 the data set can be further processed to yield an image of subterranean layers. These images are used for hydrocarbon exploration and reservoir characterization. The optional steps 21 and 23-25 are indicated in FIG. 2 as dashed blocks.

FIG. 3 shows a panel of the exact wavefield decomposition filter using equation [3] in the top left, difference between equation [3] and the filter approximation [4a] in the top right, difference between equation [3] and the filter approximation [4b] in the bottom left, and difference between equation [3] and the filter approximation [4c] in the bottom right. Note that these plots only assess the accuracy of the filter approximations and do not include the error due to their discretization. In other words, the plots do not show inaccuracies related to how the different terms in the spatial filter approximations [4a-4c] are implemented (e.g., using 3-point or 5-point derivative approximation). This would of course introduce a dependence on frequency as well. However, this is of secondary importance as appropriate approximations that are sufficiently accurate are straightforward to find.

From the top right of FIG. 3, it can be seen that equation [4a] which was obtained by making a Taylor expansion in the y-direction only results in the best approximation of the three examples for azimuths close to the in-line cable direction. Equation [4b] which is used in the difference plot at the bottom left of FIG. 3 results in an approximation which is equally good along all azimuths. An advantage with this filter is that it can be fully implemented as a compact filter. Equation [4c] which is used in the difference plot at the bottom right of FIG. 3 results in a fully accurate approximation both along the in-line and cross-line azimuths. The filter can be implemented using a compact filter approximation for the cross-term only. Exactly which of the alternative implementation [4a-4c] that is most attractive may vary depending on different combinations of requirements in terms of computational cost (CPU and data access) and accuracy.

While the invention has been described in conjunction with the exemplary embodiments described above, many equivalent modifications and variations will be apparent to those skilled in the art when given this disclosure. Accordingly, the exemplary embodiments of the invention set forth above are considered to be illustrative and not limiting. Various changes to the described embodiments may be made without departing from the spirit and scope of the invention.

The above approximations or similar approximations can be used for example with the separations developed by Amundsen et al. in the above cited U.S. Pat. No. 6,101,408 or in: “Multiple attenuation and P/S splitting of multicomponent OBC data at a heterogeneous sea floor”, Wave Motion 32 (2000), 67-78. In the latter document, demultiple or decomposition equations are found for elastic decomposition (particle velocity, traction) or P/S wave splitting below the sea floor.

Claims

1. A method of decomposing a seismic wavefield, wherein a 3D wavefield is obtained by a cross-line acquisition and filtered applying a decomposition filter having two spatial directions to obtain a decomposed wavefield.

2. The method of claim 1 wherein the decomposition is for at least one of a group consisting of up- /down going decomposition, P/S decomposition, elastic decomposition and acoustic decomposition.

3. The method of claim 1 wherein the filter comprises in-line (kx) and cross-line components (ky) or a spatial representation of the in-line (kx) and cross-line components (ky).

4. The method of claim 1 wherein the filter is applied as a cascaded filter.

5. The method of claim 1 wherein the filter is a compact filter.

6. The method of claim 1 wherein the filter filters an obtained pressure wavefield.

7. The method of claim 1 wherein the filter exclusively filters an obtained pressure wavefield.

8. The method of claim 1 wherein the step of applying the filter is preceded by a calibration step to match geophone recordings with hydrophone recordings.

9. The method of claim 1 wherein the step of applying the filter is followed by a step of removing multiples from a component of the decomposed wavefield.

10. The method of claim 1 wherein the step of applying the filter is followed by a step of imaging or migrating the filtered wavefield to generated an image of subterranean formations.

11. The method of claim 1 wherein the wavefield is obtained through receivers located on the sea floor.

12. The method of claim 1 wherein the wavefield is obtained through receivers towed by a vessel.