Estimation of primaries by sparse inversion incuding the ghost
Eric Verschuur, Delft University of Technology SUMMARY
Today, the problem of surface-related multiples, especially in shallow water, is not fully solved. Although surface-related multiple elimination (SRME) method has proved to be suc-cessful on a large number of data cases, the involved adaptive subtraction acts as a weak link in this methodology, where pri-maries can be distorted due to their interference with multiples. Therefore, recently, SRME has been redefined as a large-scale inversion process, called estimation of primaries by sparse in-version (EPSI). In this process the multi-dimensional primary impulse responses are considered as the unknowns in a large-scale inversion process. By parameterizing these impulse re-sponses as spikes in the space-time domain, and using a spar-sity constraint in the update step, the algorithm looks for those primaries that, together with their associated multiples, explain the total input data. As the objective function in this mini-mization process truly goes to zero, the tendency for distorting primaries is greatly reduced. An additional advantage is that imperfections in the data can be included in the forward model and resolved simultaneously, such as the missing near offsets. In this paper it is demonstrated that the ghost effect can also be included in the EPSI formulation after which a ghost-free primary estimate can be obtained, even in the case the ghost notch is within the desired spectrum.
INTRODUCTION
The surface-related multiple elimination (SRME) method (a.o. Berkhout, 1982; Verschuur et al., 1992; Berkhout and Ver-schuur, 1997; Weglein et al., 1997) has developed itself as one of the standard multiple removal tools in today’s seismic data processing sequences. During the last decade, the interest in a full 3D implementation of the SRME method has grown sig-nificantly (Biersteker, 2001; Lin et al., 2004; van Dedem and Verschuur, 2005; Moore and Bisley, 2005; Baumstein et al., 2005; Dragoset et al., 2008; Aaron et al., 2008).
The attractive feature of the SRME method is that it in theory it can predict multiples without any knowledge on the subsur-face. The theory of SRME requires an impulse response of the subsurface, whereas in practice the measured data is used as this operator. Therefore, the predicted multiples exhibit a wrong observed wavelet, which needs to be corrected for in adaptive subtraction. The latter is usually based on a minimum energy criterion, which is known to be not always optimal for SRME (Nekut and Verschuur, 1998; Guitton and Verschuur, 2004; van Groenestijn and Verschuur, 2008).
Therefore, recently a new approach to multiple removal was developed by van van Groenestijn and Verschuur (2009a): esti-mation of primaries by sparse inversion (EPSI). The main dif-ference with SRME is that the two-stage processing method, being prediction and adaptive subtraction, is replaced by a full waveform inversion process: the primary reflection events are
the unknowns in this algorithm and are parameterized in a suit-able way. In van van Groenestijn and Verschuur (2009a) the adopted parameterization consists of band-limited spikes and an effective source wavelet. Baardman et al. (2010) discussed a refinement, where the wavelet was made time-variant in or-der to include the change of the observed seismic wavelet in case of complex propagation effects (fine layering, dispersion and absorption). Lin and Herrmann (2009, 2010) redefined EPSI in the curvelet domain, and van Groenestijn and Ver-schuur (2009b) and Savels et al. (2011) have shown various applications to complex synthetic and field datasets.
One advantage of writing multiple removal or – with better words – primary estimation as a large-scale inversion problem is the fact that imperfections in the input data, that usually have a distorting effect on the SRME-output, within EPSI can be made part of the inversion process. This was already shown in van Groenestijn and Verschuur (2009a) for the missing near offset data, particularly for shallow water, where EPSI could recover them in order to get optimum multiple prediction. Another effect that need be accounted for in field data is the ghost. If not removed properly, the traditional SRME result will be sub-optimum. Weglein et al. (1997) already noted that input data for SRME should be properly deghosted in order to satisfy the physical relationships. The standard procedure to handle this is to apply (advanced) deghosting algorithms, such as described in Wang et al. (2012), before SRME is carried out. However, with current acquisition techniques aiming at broad-band data and/or deeper towed streamers, the ghost notch usu-ally appears inside the desired frequency spectrum. Therefore, we will include the ghost within the EPSI formalism in order to invert for it simultaneously, such that we are estimating the primary impulse responses without the ghost influence.
SURFACE-RELATED MULTIPLE ELIMINATION
In Berkhout and Verschuur (1997) it has been proposed to rewrite the surface-related multiple removal scheme of Ver-schuur et al. (1992) as an iterative procedure:
P(i+1)0 = P − P(i)0 AP, (1)
P(i)0 representing the pre-stack data containing the estimated upgoing primaries and internal multiples in iteration i, P be-ing the total upgobe-ing data (primaries and all multiples) and A representing the so-called surface operator:
A= S−1R∩, (2)
in which the inverse source properties (S−1) are combined with the reflectivity at the free surface (R∩). The matrix notation is taken from Berkhout (1982). Each column of a data matrix, e.g. P, contains a wave field (or a shot record) for one fre-quency. Thus, the full matrix contains a pre-stack data volume for one frequency component. The upgoing primary data P0
can be written as the source matrix times the primary impulse response matrix:
P0= X0S. (3)
Each column of the source matrix S contains the effective down-going wavefield for one shot record.
In practice, the directivity effects are often neglected or taken into account in a separate preprocessing step (a deghosting process), such that matrix A can be written as a frequency de-pendent scalar A(ω) (Verschuur and Berkhout, 1997). Thus, the prediction of the surface-related multiples can be written as:
ˆ
M(i+1)= P(i)0 P, (4)
after which A(ω) is determined in an adaptive subtraction pro-cess, often based on minimum energy in the primaries:
P(i+1)0 = P − A(ω) ˆM(i+1). (5) Usually, this adaptive subtraction is applied in a localized fash-ion (Verschuur and Berkhout, 1997), such that (residual) di-rectivity effects and prediction imperfections due to 3D effects can be (partly) compensated for.
ESTIMATION PRIMARIES BY SPARSE INVERSION
It has been demonstrated that the subtraction of predicted mul-tiples is the weak link in the SRME process, because it allows multiples to locally match to strong primary energy, yielding distortions of the primaries and, as a consequence, leaving residual multiple energy behind (see e.g. Nekut and Verschuur, 1998; Guitton and Verschuur, 2004; Abma et al., 2005; van Groenestijn and Verschuur, 2008). Therefore, the EPSI algo-rithm was designed to avoid this subtraction process by making the primaries the unknowns in a large-scale inversion process. To describe the EPSI algorithm (van Groenestijn and Verschuur, 2009a) we should again consider equations 1 and 3. If we take
S(ω) = S(ω)I (meaning assuming a constant source wavelet for all shots) and we assume the surface reflectivity to be a scalar R∩(being approximately−1) we get:
P= X0S+ X0R∩P. (6)
Through full waveform inversion we try to estimate the un-known, multidimensional primary impulse responses X0 and
source wavelet S such that the primaries X0S together with the
surface multiples X0R∩P can explain the observed total
up-going data P. The unknown dataset X0 is parameterized in
the time domain with spikes. The difference between the to-tal upgoing data P and the estimated primaries and multiples,
X0S− X0R∩P, is the residual V:
V= P − X0S− X0R∩P. (7)
The EPSI algorithm drives the residual V to zero, i.e. it is minimizing the following objective function:
J=P− ˆX0S− ˆX0ˆ R∩P2, (8) where ˆX0and ˆS represent the estimate of the primary impulse
responses and the estimate of the source wavelet, respectively.
This is done in an iterative way where the primary impulse re-sponse data volume X0is built up slowly during the iteration
process in the time domain. In this way the adaptive subtrac-tion is avoided and interference between primaries and multi-ples is better handled.
INCLUDING THE GHOST IN EPSI
When the ghost effect is included in EPSI, where we consider the receiver side ghost, the forward model changes to:
Pd= DX0S+ DX0R∩[D]−1Pd. (9)
where the detector operator D contains the ghost effect at the receiver side and Pdrepresents the measured data including the
detector ghost. We can see in the term at the right hand side that we first have to remove the ghost effect from the measure-ments Pd, creating the upgoing wavefield at the surface, after
which it is convolved with the surface reflectivity and the im-pulse response X0 in order to predict the multiples. Finally,
the ghost effect has to be included in the predicted multiples in order to match it with the observed data. However, assum-ing that the subsurface structures are moderate, such that the arrival angles are not too different from the incident angles for each event, the two ghost response matrices on the right hand side cancel, yielding:
Pd≈ DX0S+ X0R∩Pd. (10)
Thus, the forward model of seismic data only contains one ex-tra ghost response matrix D in order to match the estimated primary impulse responses with the observed primaries, which include the ghost effect. Knowing the receiver depth, the ghost operator is deterministic and can easily be included in the EPSI algorithm. It turns out that in the calculation of the update for
X0the ghost operator needs be involved in a conjugate mode,
but otherwise the EPSI algorithm remains largely the same.
500 1000 1500 2000 2500 3000 Depth (m) 1500 2500 3500 4500 Velocity (m/s) 500 1000 1500 2000 2500 3000 Depth (m) 1000 2000 3000 4000 Density (kg/m^3)
a) velocity profile b) density profile
Figure 1: Subsurface model based on the velocity (a) and den-sity (b) profile.
0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m) 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m)
a) shot with multiples b) shot without multiples
Figure 2: Synthetic dataset from a subsurface model based on the profiles of Figure 1. a) Shot gather with all multiples in-cluded. b) Shot gather with primaries and internal multiples.
EXAMPLES
We will demonstrate the effect of ghosts on SRME and EPSI with an example for a horizontally layered model based on the velocity and density profile as shown in Figure 1. The water bottom depth is moderate, being 400 m.
First, we forward model a seismic response that does not in-clude a ghost effect at the receiver side. This resembles the situation of perfect deghosting and is shown in Figure 2a. For comparison we also display the modeled multiple-free data in Figure 2b.
Next, we make the data more realistic and bring in the ghost ef-fect at the receiver side. We have done this quite dramatically, with a streamer depth of 30m. This results in a spectral notch at 25 Hz, being in the middle of the spectrum (the data is mod-eled up to 60 Hz). The result is shown in Figure 3b. The effect of the receiver ghost is clearly visible in the observed signature of the data, when compared to the data without ghost (Figure 2a). Standard practice is to apply a deghosting process via a deconvolution in the wavenumber-frequency domain. How-ever, due to the notch in the spectrum, this deconvolution has to be stabilized quite strongly. The deghosted data is displayed in Figure 3b. Note that the effect of the notch is visible as residual ringing.
Now, we will investigate the effect of applying the traditional SRME algorithm to the ideal ghost-free data and the data af-ter deghosting.We use three iaf-terations and employ typical pa-rameters for the adaptive subtraction. The SRME outputs are displayed in Figure 4. Note that even with the ghost-free data, the SRME output suffers from adaptive subtraction issues, pro-ducing a leakage of multiples into our primaries (Figure 4a). If
0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m) 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m)
a) shot with 30 m ghost b) shot after deghosting
Figure 3: a) Modeled dataset with a 30m receiver ghost effect. b) Result of stabilized deghosting.
we use the data after the deghosting process, SRME provides the result in Figure 4b. Note that this result looks quite accept-able, although the deghosting artifacts are still visible and the subtraction problems are similar as with the ghost-free result (Figure 4a).
Because SRME does not give optimum results, even for the perfect ghost-free input data, we consider the EPSI method for estimation of the primaries. First, EPSI is applied to the ghost-free input data, resulting in the output plotted in Figure 5a, where we see the band-limited version of the estimated primary impulse response data, i.e. the estimated X0. Note
that the primaries are very well recovered and can be compared well with the forward modeled primaries in Figure 2b. For the sake of comparison to SRME, the data with the 30 m receiver ghost after deghosting (Figure 3b) is also fed into the EPSI algorithm, with the same parameter settings as the pre-ceding example. The result is shown in Figure 5b, again being the estimated band-limited impulse response data X0.
Surpris-ingly enough, this is quite a good result. It has the favorable properties of the EPSI result on the ghost-free data in terms of multiple suppression and correct primary recovery. However, the ringing effects form the deghosting are still visible. Note that this EPSI result is much better than the SRME result for the deghosted data (Figure 4b).
However, the best result is expected when the ghost effect is included in the EPSI procedure. Therefore, this extended ver-sion of EPSI, described in the previous section, is applied to the data with ghost, shown in Figure 3b. The result of this EPSI process is displayed in Figure 6a. Note that the estimated im-pulse response does not show the ringing effects as visible in the deghosted EPSI result (Figure 5b). The suppression of the ghost effect in the EPSI result is confirmed by comparing the
0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m) 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m)
a) SRME output, ghost-free b) SRME after deghosting
Figure 4: a) SRME output for the data modeled without ghost. b) SRME output for the data after stabilized deghosting.
FK spectrum of the input data with ghost, displayed in Figure 6b showing a clear ghost notch effect, with the FK-spectrum of the estimated impulse response, shown in Figure 6c. Note that the impulse response provides a broad-band primary result.
CONCLUSIONS
Although SRME has provided good results in numerous field data cases, it often has some problems in perfectly subtract-ing the predicted multiples from the input data, resultsubtract-ing in primary distortion and residual multiple energy in the SRME output. For this reason the EPSI process has been developed, in which primaries are estimated via a large-scale, full waveform inversion process, where the estimated primaries, together with their corresponding surface multiples, should explain the ob-served seismic data. Thus, the residual data is really driven to zero without the trend of multiples distorting primary energy. One complicating factor in practice is the ghosts that are present in the data. SRME is based on ghost-free data and, therefore, a deghosting process should be applied in advance. For low-frequency enhancing deep water tows, the ghost will have a notch inside the desired data spectrum. This complicates the deghosting process, yielding artifacts in the SRME result. EPSI appears to be much less sensitive to these notch effects, due to the non-linear optimization that is involved. When EPSI is directly applied to the deghosted data, still an acceptable primary estimate can be obtained.
However, the best solution is making the ghost effect part of the EPSI procedure, such that EPSI works directly on the data with the ghost still included. In this way primary estimation and ghost removal is optimally combined in one process and the obtained output has the desired broadband characteristics.
0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m) 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m)
a) X0from ghost-free input b) X0from deghosted input
Figure 5: Estimated primary impulse responses, after band-limitation for display purpose, from EPSI for a) the ghost-free input data of Figure 2a and b) the deghosted input data of Fig-ure 3b.
ACKNOWLEDGMENTS
The author thanks the sponsoring companies of the Delphi consortium for their support. Furthermore, he thanks prof. Berkhout for useful discussions on this topic.
0 20 40 60 80 -0.04 -0.02 0 0.02 0 20 40 60 80 -0.04 -0.02 0 0.02 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 Time (s) -2000 -1000 0 1000 2000 Offset (m)
a) X0from modified EPSI
b) FK spectrum of input
c) FK spectrum of X0
Figure 6: Result from EPSI for the input data with the 30 m receiver ghost effect (Figure 3b). a) Estimated primary im-pulse responses, after band-limitation for display purpose. b) FK spectrum of the input data with ghost. c) FK spectrum of the result of the modified EPSI method.
REFERENCES
Aaron, P., R. O’Toole, S. Barnes, R. F. Hegge, and R. G. van Borselen, 2008, True-azimuth versus zero-azimuth 3-D multiple prediction in WATS processing: 78th Ann. Inter-nat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 2431–2435.
Abma, R. L., M. M. N. Kabir, K. H. Matson, S. Mitchell, S. A. Shaw, and W. McLain, 2005, Comparisons of adaptive subtraction methods for multiple attenuation: The Leading Edge, 24, 277–280.
Baardman, R. H., D. J. Verschuur, R. G. van Borselen, M. O. Frijlink, and R. F. Hegge, 2010, Estimation of primaries by sparse inversion using dual sensor data: 80th Ann. Inter-nat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 3468–3472.
Baumstein, A., M. T. Hadidi, D. L. Hinkley, and W. S. Ross, 2005, A practical procedure for application of 3D SRME to conventional marine data: The Leading Edge, 24, 254–258. Berkhout, A. J., 1982, Seismic migration, imaging of acoustic energy by wave field extrapolation, A: theoretical aspects: Elsevier (second edition).
Berkhout, A. J., and D. J. Verschuur, 1997, Estimation of mul-tiple scattering by iterative inversion, part I: theoretical con-siderations: Geophysics, 62, 1586–1595.
Biersteker, J., 2001, Shell’s surface multiple attenuation tech-nique: 71st Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 1301–1304. Dragoset, W. H., I. Moore, M. Yu, and W. Zhao, 2008, 3D
gen-eral surface multiple prediction: An algorithm for all sur-veys: 78th Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 2426–2430. Guitton, A., and D. J. Verschuur, 2004, Adaptive subtraction of
multiples using the l1-norm: Geophys. Prosp., 52, 27–38. Lin, D., J. Young, and Y. Huang, 2004, 3D SRME application
in the Gulf of Mexico: 66th Ann. Internat. Mtg., EAGE, Ex-panded Abstracts, Eur. Ass. of Expl. Geophys., ExEx-panded abstracts, G035.
Lin, T. Y., and F. J. Herrmann, 2009, Unified compressive sensing framework for simultaneous acquisition with pri-mary estimation: 79th Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 3113– 3117.
——–, 2010, Sparsity-promoting migration from surface-related multiples: 80th Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 3333– 3337.
Moore, I., and R. Bisley, 2005, 3D surface-related multiple prediction (SMP): A case history: The Leading Edge, 24, 270–284.
Nekut, A. G., and D. J. Verschuur, 1998, Minimum energy adaptive subtraction in surface-related multiple elimination: 68th Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 1507–1510.
Savels, T., K. de Vos, and J. W. de Maag, 2011, Surface mul-tiple attenuation through sparse inversion; results for com-plex synthetics and real data: First Break, 29, 55–64. van Dedem, E. J., and D. J. Verschuur, 2005, 3D
surface-related multiple prediction: A sparse inversion approach:
Geophysics, 70, V31–V43.
van Groenestijn, G. J. A., and D. J. Verschuur, 2008, Towards a new approach for primary estimation: 78th Ann. Inter-nat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Expanded abstracts, 2487–2491.
——–, 2009a, Estimating primaries by sparse inversion and application to near-offset data reconstruction: Geophysics,
74, A23–A28.
——–, 2009b, Estimation of primaries and near offsets by sparse inversion: Marine data applications: Geophysics, 74, R119–R128.
Verschuur, D. J., and A. J. Berkhout, 1997, Estimation of mul-tiple scattering by iterative inversion, part II: practical as-pects and examples: Geophysics, 62, 1596–1611.
Verschuur, D. J., A. J. Berkhout, and C. P. A. Wapenaar, 1992, Adaptive surface-related multiple elimination: Geophysics,
57, 1166–1177.
Wang, Z., A. B. Weglein, J. D. Mayhan, P. Terenghi, and C. Rivera, 2012, Greens theorem derived deghosting: fun-damental analysis, numerical test results, and impact on iss free-surface multiple elimination: 82nd Ann. Internat. Mtg., SEG, Expanded Abstracts, Soc. Expl. Geophys., Ex-panded abstracts, 5p.
Weglein, A. B., F. A. Gasparotto, P. M. Carvalho, and R. H. Stolt, 1997, An inverse scattering series method for attenu-ating multiples in seismic reflection data: Geophysics, 62, 1975–1989.
EDITED REFERENCES
Note: This reference list is a copy-edited version of the reference list submitted by the author. Reference lists for the 2013 SEG Technical Program Expanded Abstracts have been copy edited so that references provided with the online metadata for each paper will achieve a high degree of linking to cited sources that appear on the Web.
REFERENCES
