Surface wave retrieval in layered media using seismic interferometry by multidimensional deconvolution

Geophysical Journal International Geophys. J. Int. (2014) 196, 230–242 Advance Access publication 2013 November 15 doi: 10.1093/gji/ggt389 Surface wave retrieval in layered media using seismic interferometry by multidimensional deconvolution Karel N. van Dalen,1 Kees Wapenaar1 and David F. Halliday2 1 Department Accepted 2013 September 23. Received 2013 September 17; in original form 2013 May 3 SUMMARY Virtual-source surface wave responses can be retrieved using the crosscorrelation (CC) of wavefields observed at two receivers. Higher mode surface waves cannot be properly retrieved when there is a lack of subsurface sources that excite these wavefields, as is often the case. In this paper, we present a multidimensional-deconvolution (MDD) scheme that is based on an approximate convolution theorem. The scheme introduces an additional processing step in which the CC result is deconvolved by a so-called point-spread tensor. The involved point-spread functions capture the imprint of the lack of subsurface sources and possible anelastic effects, and quantify the associated spatial and temporal smearing of the virtualsource components that leads to the poor surface wave retrieval. The functions can be calculated from the same wavefields as used in the CC method. For a 2-D example that is representative of the envisaged applications, we show that the deconvolution partially corrects for the smearing. The retrieved virtual-source response only has some amplitude error in the ideal situation of having the depth of the required vertical array equal to the depth penetration of the surface waves. The error is due to ignored cross-mode terms in the approximate convolution theorem. Shorter arrays are also possible. In the limit case of only a single surface receiver, the retrieved virtual-source response is still more accurate than the CC result. The MDD scheme is valid for horizontally layered media that are laterally invariant, and includes exclusively multicomponent point-force responses (rather than their spatial derivatives) and multicomponent observations. The improved retrieval of multimode surface waves can facilitate dispersion analyses in shallow-subsurface inversion problems and monitoring, and surface wave removal algorithms. Key words: Interferometry; Surface waves and free oscillations; Interface waves. 1 I N T RO D U C T I O N Seismic interferometry is a technique that uses observations of a wavefield to create virtual seismic sources at locations where only receivers are present. More specifically, by crosscorrelating observations at two different seismic receivers, one retrieves or reconstructs an approximation to the Green’s function (i.e. the impulse response or virtual-source response) as if one of the receivers were a source (e.g. Campillo & Paul 2003; Larose et al. 2006; Wapenaar & Fokkema 2006; Schuster 2009; Snieder et al. 2009). This technique has been used in many different applications. For example, virtual subsurface sources were created using real surface sources (Bakulin & Calvert 2006), and virtual-source reflected waves were retrieved from background noise recordings (Draganov et al. 2007); both examples show that body waves can be retrieved successfully. Surface wave retrieval has received relatively wide attention, which can be explained by the fact that wavefield recordings are often dominated by surface waves, particularly when the sources are located close to the Earth’s surface. Applications exist on dif- 230  C ferent scales. In regional seismology, virtual-source surface wave responses can be retrieved by applying interferometry to so-called passive wavefields that are excited by ambient noise sources (such as ocean storms; Shapiro et al. 2005) or by earthquakes (by exploiting the seismic coda; Campillo & Paul 2003). Using such ambient noise sources, one can retrieve the virtual-source response at frequencies that are difficult to generate using active sources (Wathelet et al. 2004; Park et al. 2007). The results are used to determine group velocity images, which is often done based on the fundamental-mode Rayleigh wave only (e.g. Shapiro et al. 2005; Gerstoft et al. 2006; Bensen et al. 2007); in general, the surface wave response in a layered medium consists of a superposition of modes. In exploration seismology, virtual-source surface waves can be retrieved using active sources at the surface. The retrieved surface wave responses are used to guide filters designed to suppress surface waves in seismic data as they overshadow the much weaker body waves reflected at deeper targets of interest (e.g. Halliday et al. 2007, 2010). In near-surface seismology and geotechnical engineering, retrieved surface waves, obtained from either active or passive sources, are The Authors 2013. Published by Oxford University Press on behalf of The Royal Astronomical Society. GJI Marine geosciences and applied geophysics of Geoscience and Engineering, Delft University of Technology, Stevinweg 1, NL-2628 CN Delft, The Netherlands. E-mail: 2 Schlumberger Gould Research, High Cross, Madingley Road, Cambridge, CB3 0EL, UK Surface wave retrieval in layered media analysis to determine the modal scale factors of the proposed approximate convolution theorem in Section 3. In Section 4, we derive the MDD scheme; numerical examples are shown in Section 5. The envisaged applications of the MDD scheme in a geophysical context are discussed in Section 6. Finally, we summarize our conclusions in Section 7. 2 S TAT I O N A RY- P H A S E A N A LY S I S O F C O N VO LU T I O N T H E O R E M In this section, we apply the convolution-type reciprocity theorem to surface wave responses and evaluate the associated integral using the stationary-phase method. Snieder (2004) and Halliday & Curtis (2008) did a similar derivation for the correlation-type reciprocity theorem. In the next section, we use the results to derive an approximate convolution theorem for surface waves. Throughout this paper, we adopt the following Fourier transform over time for an arbitrary function f (x, t):  ∞ f (x, t) exp(−iωt) dt, (1) fˆ(x, ω) = −∞ where ω denotes angular frequency, t denotes time and x = [x, y, z]T is a vector containing spatial coordinates. The hat refers to the (x, ω) domain. From here onwards, the ω dependence is left out for brevity [i.e. fˆ(x) = fˆ(x, ω)]. All considered functions are real-valued in the space–time domain, and it is therefore sufficient to consider ω ≥ 0. The space–frequency domain elastodynamic convolution-type reciprocity theorem reads (de Hoop 1995; Aki & Richards 2002)   Ĝ in (x R , x)n j cn jkl (x)∂k Ĝ lm (x, x S ) Ĝ im (xR , x S ) = S  − n j cn jkl (x)∂k Ĝ il (x R , x)Ĝ nm (x, x S ) dS, (2) where S denotes the enclosing boundary of a volume V, which has outward-pointing unit normal nj (see Fig. 1, where we consider the specific case of a cylindrical volume in view of the analysis below); summation is invoked over repeated indices, but n (...truncated)