# Phase-velocity dispersion curves and small-scale geophysics using noise correlation slantstack technique

Pierre Gouedard, C. Cornou, P. Roux

► To cite this version:

Pierre Gouedard, C. Cornou, P. Roux. Phase-velocity dispersion curves and small-scale geophysics using noise correlation slantstack technique. Geophysical Journal International, 2008, 172 (3), pp.971 à 981. 10.1111/j.1365-246X.2007.03654.x. insu-00333887

HAL Id: insu-00333887

https://insu.hal.science/insu-00333887

Submitted on 21 Aug 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Phase-velocity dispersion curves and small-scale geophysics using noise correlation slantstack technique

Pierre Gouédard<fnref n="1" />, Cécile Cornou<fnref n="1" />,<fnref n="2" /> and Philippe Roux<fnref n="1" />

Accepted 2007 October 10. Received 2007 October 10; in original from 2007 June 28

## SUMMARY

It has been demonstrated both theoretically and experimentally that the Green's function between two receivers can be retrieved from the cross-correlation of isotropic noise records. Since surface waves dominate noise records in geophysics, tomographic inversion using noise correlation techniques have been performed from Rayleigh waves so far. However, very few numerical studies implying surface waves have been conducted to confirm the extraction of the true dispersion curves from noise correlation in a complicated soil structure. In this paper, synthetic noise has been generated in a small-scale (<1 km) numerical realistic environment and classical processing techniques are applied to retrieve the phase velocity dispersion curves, first step toward an inversion. We compare results obtained from spatial autocorrelation method (SPAC), high-resolution frequency-wavenumber method (HRFK) and noise correlation slantstack techniques on a 10-sensor array. Two cases are presented in the (1–20 Hz) frequency band that corresponds to an isotropic or a directional noise wavefield. Results show that noise correlation slantstack provides very accurate phase velocity estimates of Rayleigh waves within a wider frequency band than classical techniques and is also suitable for accurately retrieving Love waves dispersion curves.

**Key words:** Surface waves and free oscillations; Computational seismology.

## 1 INTRODUCTION

The Green's function of a medium between two points A and B represents the record we would get at A if we put an impulse source at B. The use of random noise to reconstruct the Green's function has already been applied successfully in various fields of wave physics such as helioseismology, oceanography or geophysics at large scales and ultrasonics or non-destructive evaluation at small scales. Historically speaking, helioseismology was the first field where ambient noise cross-correlation performed from recordings of the sun surface random motion was used to produce coherent images of the sun subsurface (Duvall et al. 1993; Gilles et al. 1997). More recently, a seminal paper (Weaver & Lobkis 2001) has shown how diffuse thermal noise recorded and cross-correlated at two ultrasonic transducers fastened to one face of a duralumin sample provided the complete Green's function between these two points. Similar results were then obtained at a much larger scale in shallow underwater acoustics where both direct and reflected wave fronts were retrieved from ambient noise cross-correlation (Roux & Kuperman 2004). At the same time, the cross-correlation process was investigated from the perspective of multiple uncorrelated sources that were considered as noise sources (Campillo & Paul 2003; Derode et al. 2003). By summing the contribution of all sources to the correlation, it has been shown numerically that the correlation contains the causal and

acausal Green's function of the medium (Wapenaar 2004). At ultrasonic scales, the noise-extracted Green's function has been used to detect defects in duralumin samples (Larose et al. 2006a).

In seismology, Aki (1957) proposed fifty years ago to use seismic noise to retrieve the dispersion properties of the subsoil. For elastic waves, it has been theoretically showed that the convergence of noise correlation to the Green's function was bounded by the equipartition condition of the different components of the elastic field (Sánchez-Sesma et al. 2006; Sánchez-Sesma & Campillo 2006). In other words, the emergence of the Green's function is effective after a sufficient self-averaging process that is provided by random spatial distribution of the noise sources when considering long time series as well as scattering (Campillo 2006; Larose et al. 2006b). Shapiro & Campillo (2004) reconstructed the surface wave part of the Green's function by correlating seismic noise at stations separated by distance of hundreds of kilometres, and measured their dispersion curves at periods ranging from 5 to 20 s. This method led to the first application of passive seismic imaging in California (Sabra et al. 2005; Shapiro et al. 2005). Since then, several studies focused on seismic tomography at lithospheric scale by using such noise correlation techniques (Yao et al. 2006; Bensen et al. 2007; Lin et al. 2007; Pedersen et al. 2007). Some authors have also used this technique to address the issue of noise source spatial distribution in the cross-correlation process (Pedersen et al. 2007;

<footnote n="1">¹LGIT, CNRS 5559, Université Joseph Fourier, Grenoble, France. E-mail: pierre.gouedard@ujf-grenoble.fr</footnote>

<footnote n="2">²IRD, Institut de Recherche pour le Développement</footnote>
Stehly *et al.* 2007). However, despite convincing tomographic im-
ages of the subsurface shear velocity at large scales (Sabra *et al.*
2005; Shapiro *et al.* 2005) there is no study that confirms that the
noise correlation process yields the true phase velocity dispersion
curves in a realistic environment where several surface wave modes
are present. Such a work would imply (1) the computation of the
elastic field produced by multiple noise sources in a given subsurface
model, (2) the recording of the numerical noise traces on a seismic
array, (3) the computation of the noise-correlation function between
each receiver pair and (4) the extraction of the phase-velocity dis-
persion curves from the noise-correlation process for a comparison
to the actual dispersion curves directly obtained from the numerical
model. This is the objective of this paper in the frame of earthquake
engineering implying small-scale seismic networks (<500 m) and
high-frequency seismic noise (>1 Hz).

In earthquake engineering, microseisms and microtremors (noise produced by anthropogenic sources) have been used for more than 30 yr for estimating subsurface S-wave velocities and for site effect purposes. Microtremor studies indeed originated in the pioneering work of Japanese authors (Kanai et al. 1954; Aki 1957; Nogoshi & Igarashi 1971; Nakamura 1989). In recent decades, the use of microtremor array measurements [mainly frequency-wavenumber based methods (Burg 1964; Capon 1969) and the spatial autocorrelation (SPAC) technique (Aki 1957)] applied to ambient-vibration noise wavefields have spread throughout the world as a substitute of borehole measures or active seismic methods, which are costly and even prohibitive in urban areas. Recently, new methods have emerged: Bettig et al. (2001), Ohori et al. (2002), Asten et al. (2004) and Cho et al. (2004, 2006) proposed modification and/or extension of the SPAC technique, while Louie (2001) presented the linear slantstack method for ambient noise recordings. Although noise correlation technique gave rise to numerous applications and studies in large-scale seismology since Shapiro & Campillo (2004), this technique was only applied by Chavez-Garcia et al. (2005) to subsurface structure imaging.

In summary, the goal of microtremor array measurements and noise correlation techniques is to measure dispersion curves from which subsurface shear wave profiles are then extracted. It has to be noted however that even if microtremor array analysis is performed from 2-D arrays, the interpretation of the resulting dispersion curves are mostly performed assuming wave propagation in 1-D structures.

Within the framework of the last symposium on the Effects of Surface Geology on seismic motion (ESG2006), an international noise blind test (Cornou et al. 2006) was proposed in order to compare surface wave dispersion characteristics derived from competing analysis approaches and to make a clear assessment regarding the potential of microtremor array studies for site effect estimation. This blind test involved both synthetic and real data sets and was opened to a large scientific community, with no restrictions regarding the choice of analysis approaches. Synthetic data provided the opportunity to perform a benchmark test where the site structure and the wavefield situation are fully controlled. For SPAC and FK based techniques—which were the most frequently used techniques by participants—this exercise pinpointed that surface wave dispersion characteristics were well retrieved over frequency bands controlled by the resolution capabilities of seismic array and the energy content of ambient noise. Interestingly, the only participating group that used noise correlation (Gouédard et al. 2006) has shown very good phase velocities estimates for both Rayleigh and Love fundamental and higher modes outside the classical array proxy-capabilities at high frequency, estimates at low frequency being as good as the ones provided by classical SPAC and FK techniques.

The scope of this paper is twofold. First, we numerically demon-
strate that the noise-correlation process leads to the true phase–
velocity dispersion curves even when several surface wave modes
are present. Second, we clarify reasons why noise correlation en-
ables to provide phase velocity estimates of both Rayleigh and Love
waves at higher frequencies compared to classical processing ap-
proaches. This issue was beyond the scope of the ESG2006 noise
blind test and could not be fully addressed since participants used
different array layouts. We simulate in this paper seismic ambient
noise for one of the models proposed within the ESG2006 noise
blind test for an array configuration that is suitable for correlation
techniques. As a main issue in the convergence of the correlation
technique to the Green’s function is the influence of the temporal and
spatial distribution of noise sources, we computed noise synthetics
in Section 2 for both isotropic and unidirectional spatial distribution
of sources. In Sections 3 and 4, FK and SPAC analysis as well as
noise correlation analysis are then applied on the vertical compo-
nent. Results obtained by different methods are compared to the true
dispersion curves and to one another. Finally, we conclude in Sec-
tion 5 on practical considerations for the implementation of noise
correlation techniques for subsurface passive imaging.

2 AMBIENT NOISE SIMULATION
AND ARRAY CONFIGURATION

## 2.1 Ambient noise simulation

Ambient noise was simulated in a complex shallow structure with
strong impedance contrast and complex layering including low-
velocity zones. This ground profile is similar to model N102 pro-
posed within the framework of ESG2006 Noise Blind test (Cornou
et al. 2006). The compressional-(Vp) and shear wave (Vs) profiles of
this model as well as corresponding dispersion curves are displayed
in Fig. 1.

Regarding noise synthetics generation, noise sources were ap-
proximated by subsurface forces located at 0.5-m depth with ran-
dom force orientation and amplitude (Moczo & Kristek 2002). Dis-
tribution of sources is random in time. In this study, two distinct
distributions of source locations were considered: an isotropic and
a unidirectional distribution as depicted in Fig. 2. The source time
function employed at each point location is a delta-like signal with
a frequency Fourier amplitude spectrum flat from 0.1 to 20 Hz.
Computation up to 20 Hz of the associated wave field has then been
performed using the wavenumber-based technique of Hisada (1994,
1995) for 1-D horizontally layered structures. In case of isotropic
source distribution, noise data set (hereafter called isotropic noise) is
composed of three sets of 30 min of noise synthetics, each of them
implying source shots at 74–268 different locations excited at 60
different random times, while for the directional source distribution
one data set of 30 min was considered (hereafter called directional
noise). Such duration of noise recordings is similar to duration usu-
ally considered in real world experiment. Fig. 3 displays example of
noise synthetics.

## 2.2 Choice of an optimized seismic array layout

This paper goal is to compare correlation technique, SPAC and
FK techniques. To do so, the array was designed to be suitable for
correlation technique as well as for SPAC and FK techniques and
to mimic classical ambient noise studies. Following criteria were
followed.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** (Left-hand panel) P- and S-wave velocity profiles (grey and black curves, respectively); (right-hand panel) corresponding Rayleigh and Love waves dispersion curves for the fundamental and first higher modes.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** (a) Source locations (red dots) for isotropic noise. (b) Source locations (red dots) for directional noise. Receiver locations are indicated by black squares. (c) Fourier amplitude spectra of the delta-like signal used as source time function.</figcaption></figure>

(1) The number of receivers is restricted to 10.

(2) As the noise is not necessarily isotropic, a large number of azimuths have to be scanned by the different pairs of receivers to average over angle distribution.

(3) The spatial extension of the array has to be adapted to the targeted wavelength interval. Expecting a 500 m s⁻¹ surface wave velocity at 5 Hz, the wavelength of interest is centred around 100 m.

(4) The distance between neighbouring stations must sample a range interval with minimal gaps and over sampling, covering be-

tween 1 and 5 wavelengths to provide enough resolution in the F-K domain.

The spiral shape of the array shown in Fig. 4 fulfils the above criteria. This array roughly presents an even interreceiver distance range sampling from 72 to 637 m.

## 2.3 Array response

Any seismic array configuration can be considered as a discrete spatial sampling of the continuous seismic wavefield in space and

<figure><img src="image_4.png" /><figcaption>**Figure 3.** A typical signal (velocity) for the vertical component.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Figure 4.** (a) array layout; (b) azimuthal distribution of receiver pairs; (c) Range sampling for receiver pairs obtained using the array in (a) (circles). Triangles correspond to the projected distribution of ranges according to a directional noise with a backazimuth of 61° from the North.</figcaption></figure>

time. The sampling theorem then holds and the short-wavelength part of the wavefield cannot be recovered uniquely (spatial aliasing). For linear equidistantly spaced sensor arrays, the relation between the interstation distance $d_{min}$ of neighbouring stations and the spatial Nyquist frequency $\lambda_{Nyq}$ is such that each wavelength needs to be sampled (equidistantly) by at least two discrete sampling locations:

$$ \lambda_{Nyq} = \lambda_{min} = 2d_{min} \qquad (1) $$

The resolution capability of a seismic array is defined as the ability to separate two waves propagating at closely spaced wavenumbers. This resolution is related to the maximum interstation distance $d_{max}$, that is, to the array aperture:

$$ \lambda_{max} = d_{max} \qquad (2) $$

Eq. (2) is in a strict sense only valid for 1-D array layouts and the conventional beamforming algorithm (Burg 1964; Lacoss et al. 1969). The sensor geometry having in general a 2-D irregular shape, both the aliasing condition and the resolution capability depend on the direction of the incident wavefield and the effective smallest/largest interstation distance along the wave propagation direction (Henstridge 1979; Asten & Henstridge 1984; Gaffet 1998; Ohnberger 2005; Asten 2006; Okada 2006; Wathelet et al. 2007). Therefore, for arbitrary 2-D array geometries there is no simple analytic expression which relates the array shape with its spatial aliasing or resolution limits. Despite this fact, observation and numerical evaluation of the array response function have shown that the wavelength limits ($\lambda_{min}$, $\lambda_{max}$) as derived above can serve as a first order proxy to specify the resolution capabilities of a seismic array. However, both high-resolution frequency-wavenumber (hereafter called HRFK) techniques (Capon 1969) as well as the SPAC technique (Aki 1957) show improved resolution capability compared to the conventional beamformer. While for the high-resolution FK method it has been repeatedly reported (e.g. Woods & Lintz 1973; Asten & Henstridge 1984; Tokimatsu 1997; Satoh et al.

2001; Okada 2003) that the longest resolvable wavelength $\lambda_{max}$ is around three to six times longer than for the conventional technique ($\lambda_{max} \sim 3-6 d_{max}$), many authors report a superior performance for the SPAC technique compared to HRFK for the longer wavelength limit (e.g. Okada 2003; Cornou et al. 2006). According to many studies (Horike 1985; Miyakoshi 1996; Asten et al. 2004) the longest resolvable wavelength $\lambda_{max}$ is in the order of 10–15 times the radius of the array configuration ($\lambda_{max} \sim 5-7.5 d_{max}$).

Since we use SPAC, HRFK and correlation techniques in this paper, we define for the array resolution capability the following wavelength limits ($\lambda_{min} = 2d_{min}$, $\lambda_{max} = 3d_{max}$).

# 3 PHASE VELOCITY DISPERSION CURVES OBTAINED FROM HRFK AND SPAC ANALYSIS

## 3.1 High-resolution frequency-wavenumber (HRFK) analysis

We used the HRFK technique (Capon 1969) as implemented in the sesarray software (http://www.geopsy.org; Wathelet et al. 2007). Operating with sliding time windows and narrow frequency bands, this method provides the wave propagation parameters (azimuth and slowness as a function of frequency) of the most coherent plane wave arrivals. The wave propagation on the vertical component was measured by using 200 frequency bands between 0.5 and 20 Hz. The central frequency of each band was selected to be equally spaced in logarithm scale. A fraction of the central frequency $f_c$ defined the frequency bandwidth (0.97 $f_c$–1.03 $f_c$). We selected the time window length as 300 times the central period corresponding to the analysed frequency band $f_c$. Fig. 5 displays the normalized histograms of the phase velocity estimates as a function of frequency for both noise data sets (isotropic and directional noise). Consistently
<figure><img src="image_6.png" /><figcaption>**Figure 5.** (Top) HRFK frequency–velocity normalized histograms for isotropic and directional noise; (Bottom) Normalized histograms derived from SPAC analysis for both isotropic and directional noise. Limits of array response given the following wavelength limits [$\lambda_{min} = 2d_{min}$, $\lambda_{max} = 3d_{max}$] are indicated by dashed lines. Full lines correspond to theoretical Rayleigh fundamental (R0) and first higher mode (R1).</figcaption></figure>

with the resolution capabilities of the array, measurable phase velocities lie between resolution ($\lambda_{max} = 3d_{max}$) and aliasing limits of the array ($\lambda_{min} = 2d_{min}$). Consistently also with the known capability of FK technique that performs better for directional wave propagation, the first higher mode branch of Rayleigh waves dispersion is well retrieved for directional noise.

## 3.2 Spatial autocorrelation (MSPAC) analysis

Spatially average autocorrelation coefficients were computed on the vertical component following the modified vertical SPAC method (MSPAC) developed by Bettig *et al.* (2001) and implemented by

Wathelet *et al.* (2004). This method is particularly suitable here since it allows arbitrary array layouts. Autocorrelation coefficients were estimated by using sliding time windows as being three hundred times the central period corresponding to the analysed frequency band $f_c$. As for FK analysis we use 200 frequency bands between 0.5 and 20 Hz. Histograms of all possible phase velocities corresponding to measured autocorrelation coefficients are displayed in Fig. 5 for both isotropic and directional noise data sets. In case of isotropic noise, SPAC provides good estimates of phase velocities in between the array limits. For directional noise, SPAC estimates are deteriorated compared to results derived from the isotropic noise case due to the sparse distribution of interstation azimuths for such
directional sources. We also observe that SPAC fails in identifying higher modes which is an inherent limitation of classical SPAC technique (e.g. Okada 2003).

# 4 PHASE VELOCITY DISPERSION CURVES OBTAINED FROM NOISE CORRELATIONS SLANTSTACK (NCSS)

## 4.1 Principle

It has been demonstrated both theoretically and experimentally that the Green's function between two receivers can be retrieved from the cross-correlation of isotropic noise records. This technique is applied here to synthetic noise to obtain the Green's function for each receiver pair of a seismic array. For a displacement field $u(\vec{r}, t)$ measured simultaneously at two receiver locations $\vec{r}_1$ and $\vec{r}_2$, we define the normalized cross correlation over a recording time window $T$ as:

$$C_{1,2}(\tau) = \frac{\int_{0}^{T} u(\vec{r}_1, \tau) u(\vec{r}_2, t+\tau) dt}{\sqrt{\int_{0}^{T} u^2(\vec{r}_1, t) dt \int_{0}^{T} u^2(\vec{r}_2, t) dt}} \quad (3)$$

The denominator in (3) is a normalization factor that practically mitigates local station effects and helps in determining the actual coherence of the noise field between the stations. Indeed, the autocorrelation $C_{1,1}(\tau)$ has a maximum of 1 while strongly decorrelated noise signals (i.e. when stations are too far away or when the incoming noise is incoherent) do not exhibit any peak in the normalized noise correlation function but a random temporal signal with a variance of $\frac{1}{2T\Delta\omega} \ll 1$, when dealing with frequency-limited noise sources of bandwidth $\Delta\omega$ (Sabra et al. 2005; Weaver & Lobkis 2005).

For each pair of the three-component receivers, the nine-component correlation tensor is constructed by using each combination of the receiver components in the cross-correlation process. Transverse and radial components are relative to the azimuthal direction of receiver pair and are computed for each pair separately by projecting the North and East components towards transverse and radial directions. In the following only the vertical–vertical (Z–Z) component of the correlation tensor will be studied for a fair comparison with SPAC and HRFK methods. Rayleigh waves are expected to be mainly observed on Z–Z and R–R components. Transverse–

transverse (T–T) components should also provide useful information on Love waves, as mentioned later on. Other components of the correlation tensor do not provide direct information to measure dispersion curves.

As the medium is horizontally stratified, the extracted Green's function only depends on the distance between receivers used for computing correlation. Correlation functions obtained for all receiver pairs are then plotted versus distance to construct a seismic section, as if a classical active seismic experiment was carried out with a linear array. Phase velocity dispersion curves are subsequently obtained by applying classical FK transform to the reconstructed seismic section. For sake of conciseness we will refer hereafter this above-mentioned algorithm as noise correlation slantstack technique (NCSS).

For the analysis of both isotropic and directional noise data sets that will be presented in the following sections, synthetics are pre-processed to whiten their frequency spectrum between 0.5 and 20 Hz in order to limit the impact of high-frequencies damping in the medium. Amplitude of synthetics is also made one bit (i.e. only the sign of the time-domain whitened signal is kept, see Larose et al. 2004) to get rid off temporal variation of noise sources amplitude.

## 4.2 Analysis of isotropic noise

A seismic section obtained for each receiver pair from the cross-correlation of the Z-component noise synthetics is plotted in Fig. 6(a). Contrary to classical seismic sections, the obtained seismic section does not exhibit equally spaced Green's functions in spatial domain, the minimum distance between two Green's functions being $d_{min} = 0.3$ m and the maximum distance being $d_{max} = 637.4$ m. Fig. 7(a) displays the frequency–velocity diagrams derived by applying classical FK analysis to the Z–Z seismic section. Limits of the array response are computed according to minimum and maximum distances ($d_{min}$, $d_{max}$) previously defined. Rayleigh waves' dispersion is nicely retrieved on this diagram.

## 4.3 Analysis of directional noise

In the case of directional noise, the correlation tensor is no longer connected to the Green's tensor since the wave arrival between each receiver pair depends on the noise directivity angle. However, assuming a 1-D medium, the correlation for each receiver pair is the Green's function along the direction of incident wave. If the

<figure><img src="image_7.png" /><figcaption>**Figure 6.** (a) Zoom between 0 and 300 m of the seismic section obtained for the Z–Z component of the correlation tensor accumulating all receiver pairs. (b) The same section obtained in the case of directional noise. The signal-to-noise ratio is higher despite shorter records in time, as all noise sources coherently contribute to the Green's function reconstruction.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Frequency–velocity diagrams for NCSS applied to the Z–Z component of the correlation tensor in the case of isotropic (a) and directional (b) noise. The two first theoretical Rayleigh modes (R0 and R1 modes) are indicated by full lines. Limits of array response given the following wavelength limits [$\lambda_{min} = 2d_{min}$, $\lambda_{max} = 3d_{max}$] are indicated by dashed lines (low velocity limits corresponding to $\lambda_{min}$ is very close to the frequency axis).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 8.** Matched field processing indicates impinging sources (red-colour region) propagating with a backazimuth of 61° (Y positive direction indicates North and X positive direction indicates East). Purple stars indicate station locations.</figcaption></figure>

direction of incident wave is known, it is then possible to account for such apparent Green's function by defining, for each receiver pair, an effective distance as being the projection of distance between receivers along the noise direction (or the difference in source–receiver distances if the source is in the near field). The noise correlation then converges to the Green's function related to this effective distance. The projection along the noise direction modifies each receiver pair distance, and thereafter resolution capabilities of the array. In this paper, direction of noise is first determined by applying a matched field processing (Baggeroer et al. 1993). The average source backazimuth is found to be 61° (Fig. 8) and effective distances are computed together with Green's function. As illustrated in Fig. 4(c), the minimum and maximum distances between receiver pair are reduced to 6 and 436 m, respectively. The seismic section obtained from the Z–Z correlation tensor is shown in Fig. 6(b). Compared to the isotropic

noise data set, the signal-to-noise ratio is strongly increased as each source is now coherently contributing to the Green's function reconstruction when distances are projected along the noise direction. Frequency–velocity diagram derived from the FK processing is displayed in Fig. 7(b). Since the signal-to-noise ratio of seismic sections has been increased, dispersion curve of Rayleigh waves are better retrieved than for the isotropic noise data set.

## 5 RESULTS COMPARISON

Fig. 9 displays phase velocity estimates as a function of frequency obtained for both isotropic and directive noise data sets and derived by HRFK and NCSS. In order to better quantify the ability of each technique in estimating phase velocities, we have displayed in Fig. 10 the relative deviation of measured phase velocities from the true velocities (see Cornou et al. 2006 for details about the plotting procedure). The most striking feature is the ability of NCSS technique in estimating Rayleigh wave phase velocities at higher frequency than HRFK technique for both directional and isotropic noise. Indeed, HRFK provides phase velocity estimates between 1 and 3 Hz, while NCSS technique provides estimates between 1 and 12 Hz. For the array layout used in this paper, the superiority of the NCSS technique mainly lies in the reconstruction of a 1-D seismic section with spatial sampling between reconstructed Green's functions smaller than the minimum distance between neighbourhood receivers. As a consequence, the spatial Nyquist frequency is reduced and shorter wavelength parts of the wavefield can thus be recovered compared to classical array approaches. On the opposite, the previously defined array proxy-capability to resolve the longest wavelength parts of the wavefield is similar for both approaches: the maximum distance between receivers for both HRFK and SPAC analysis or between Green's functions for the reconstructed seismic profile are indeed similar, except when specific processing is performed to account for directive sources. In such a case, projection of distances between receivers along the direction of the incident wave field leads to effective distances smaller than the initial distance sampling, which deteriorates array capability for longest wavelengths.
<figure><img src="image_10.png" /><figcaption>**Figure 9.** (Top) Phase velocity estimates of fundamental and higher modes of Rayleigh waves (black squares) derived from frequency–velocity normalized histograms (Fig. 5) using HRFK analysis. (Bottom) Phase velocity estimates of fundamental and higher modes of Rayleigh (black squares) and Love waves (black crosses) derived from frequency–velocity diagrams (Fig. 7) using NCSS. Limits of array response given the following wavelength limits [$\lambda_{min} = 2d_{min}$, $\lambda_{max} = 3d_{max}$] are indicated by dashed lines. Theoretical dispersion curves of fundamental and higher modes of Rayleigh and Love waves are indicated by red and blue lines, respectively.</figcaption></figure>

Table 1 summarizes the minimum and maximum distances for both HFRK/SPAC and NCSS approaches that control the array proxy-capabilities as defined in this study. One has to mention however that, besides considerations based on the theoretical array proxy-capabilities, the actual performance of techniques to resolve lowest frequencies (i.e. below 1 Hz) can not be fully addressed here due to the high pass filtering effect of the sediment layer leading to the vanishing of vertical spectral energy below 1 Hz (see Cornou et al. 2006).

Regarding effects of noise source spatial distribution, phase velocities of Rayleigh first mode (R1) are better retrieved for both techniques when noise is directional (Figs 9 and 10). This is due to the better ability of FK analysis for a directive source and to the improved signal-to-noise ratio for NCSS technique as mentioned in the previous section. For HRFK analysis and both directional and isotropic noise, phase velocities estimates are very close to the theoretical dispersion curves, the median value of absolute relative deviation ranging between 3 and 8 per cent. The estimates are however slightly biased to higher velocities, especially at longest wavelengths (below 1.2 Hz), as a consequence of insufficient array resolution (see Cornou et al. 2006). For correlation technique and directional noise, fundamental Rayleigh wave phase velocities are particularly well retrieved above 1.2 Hz: the median value of absolute relative deviation is less than 1 per cent for Rayleigh fundamental and first modes), while deviation to higher velocities can be observed below 1.2 Hz. For isotropic noise, phase velocity estimates are more scattered, especially for the fundamental Rayleigh mode below 2.7 Hz. However, median value of absolute relative velocity deviation is 2 per cent above 2.7 Hz and 7.6 per cent be-

low 2.7 Hz. Such scattering can be explained by the low signal-to-noise ratio observed on the reconstructed seismic profile compared to the one derived from directional noise (Fig. 6). Considering longer time series would certainly improve the signal-to-noise ratio and ensure that correlation converges towards more robust Green's functions.

Phase velocity dispersion curves may also be extracted from *R*–*R* and *T*–*T* components of the cross-correlation tensor. This is a clear advantage of the NCSS technique with respect to HRFK. Rayleigh waves are expected on the *R*–*R* component and Love waves on the *T*–*T* component. Using *R*–*R* component is useful to improve estimate of Rayleigh dispersion curve by providing an independent measure. *T*–*T* component of the correlation tensor offers a way to measure Love wave dispersion curves with a high accuracy as they are separated from Rayleigh waves. Figs 9 and 11 show Love wave dispersion curves obtained using NCSS analysis for both isotropic and directive noise source data sets. Phase velocity estimates are very close to the theoretical dispersion curves outlining the ability of such technique to also retrieve information on Love waves propagation.

## 6 CONCLUSION

In this paper, we have simulated ambient seismic noise for both isotropic and directional spatial distribution of noise sources on a given array layout. Noise synthetics were then analysed using HRFK, SPAC technique and NCSS technique. Results outlined that NCSS provides phase velocity estimates of Rayleigh waves within a wider frequency band (between 1 and 12 Hz) than classical
<figure><img src="image_11.png" /><figcaption>**Figure 10.** Relative deviation of phase estimates (squares) from the true dispersion curves obtained for both isotropic (right-hand panel) and directional (left-hand panel) noise data sets. Methods used to measure phase velocities are indicated (CORR stands for correlation technique, while HRFK stands for High-Resolution frequency-wavenumber analysis) together with the surface wave modes (R0 and R1 stands for fundamental and first higher Rayleigh modes, respectively; L0 stands for fundamental Love mode). The relative deviation is given in percent; grey bars correspond to the following deviation range: [−50–50 per cent], [−100–50 per cent], [50–100 per cent].</figcaption></figure>

<table><caption>**Table 1.** Minimum ($d_{min}$) and maximum ($d_{max}$) distances serving as array-resolution proxies.</caption><thead><tr><th rowspan="2">Sources</th><th colspan="2">$d_{min}$</th><th colspan="2">$d_{max}$</th></tr><tr><th>SPAC/HRFK</th><th>NCSS</th><th>SPAC/HRFK</th><th>NCSS</th></tr></thead><tbody><tr><td>Random (m)</td><td>72.3</td><td>0.3</td><td>637.4</td><td>637.4</td></tr><tr><td>Directional (m)</td><td>72.3</td><td>0.18</td><td>637.4</td><td>436.2</td></tr></tbody></table>

For SPAC and HRFK techniques, $d_{min}$ stands for the minimum distance between neighbourhood stations and $d_{max}$ stands for the array aperture. For noise correlation slantstack (NCSS) technique, $d_{min}$ and $d_{max}$ stands for the minimum and maximum distances between Green's functions observed on the reconstructed seismic profile.

techniques (between 1 and 3 Hz). We explain this superiority of NCSS technique by the processing step that consists in deriving a 1-D seismic section composed of reconstructed Green's functions with an effective spatial sampling smaller than the absolute minimum interstation distance. As a result, the effective spatial Nyquist frequency is smaller for correlation techniques than for classical approaches. Applying a frequency–wavenumber transform to the reconstructed seismic section allows the NCSS technique to recover shorter wavelengths of the wavefield than classical array approaches.

NCSS technique is a very promising tool for subsurface structure imaging as the estimates only present here 1–2 per cent error, compared to 3–8 per cent for HRFK. Furthermore, Love wave dispersion curves can also be accurately obtained with no extra effort. Here, the choice of the array geometry is obviously a key factor for retrieving dispersion curves at higher frequencies than classical

techniques. The use of typical array geometries (circles, triangles, L-shape arrays) commonly used in ambient noise array technique would certainly not lead to such good results at high frequency since the spatial sampling of such arrays do not span a wide range of distances. Further studies should be performed in order to derive the typical most appropriate array layouts suitable for NCSS techniques.

In this paper, 30 min of noise recording is enough to measure accurately Rayleigh dispersion curves from NCSS technique in the targeted frequency band. Convergence and sensitivity studies have however to be carried out to adjust the recording time to practical cases. Besides, NCSS technique performs better in the case of directional noise since it takes advantage of noise directivity to increase signal-to-noise ratio. In case of several non-isotropic noise directions, NCSS technique is no more suitable as several possible projected ranges exist for one receiver pair. The cross-correlation is then a mix between different Green's functions.

Finally, we used a here FK method, as one among many processing techniques used in classical active surveys, to measure dispersion curves from the reconstructed seismic profile. The possibility to retrieve other information than dispersion curves (e.g. ellipticity, attenuation, etc.) from cross-correlation has to be explored.

## ACKNOWLEDGMENTS

Noise synthetics were performed at the Service Commun de Calcul Intensif (SCCI) of Grenoble, France observatory (OSUG). We thank the two anonymous reviewers for their constructive comments.
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Frequency–velocity diagrams for the *T–T* component of the correlation tensor in the case of isotropic (a) and directional (b) noise. The fundamental theoretical mode of Love waves (L0) is indicated in full line. Limits of array response given the following wavelength limits [$\lambda_{min} = 2d_{min}$, $\lambda_{max} = 3d_{max}$] are indicated by dashed lines (low velocity limits corresponding to $\lambda_{min}$ is very close to the frequency axis).</figcaption></figure>

## REFERENCES

Aki, K., 1957. Space and time spectra of stationary stochastic waves with special reference to microtremors, *Bull. Earthq. Res. Inst.*, **35**, 415–456.

Asten, M.W., 2006. On bias and noise in passive seismic data from finite circular array data processed using SPAC methods, *Geophysics*, **71**(6), V153–V162.

Asten, M.W. & Henstridge, J.D., 1984. Array estimators and use of microseisms for reconnaissance of sedimentary basins, *Geophysics*, **49**, 1828–1837.

Asten, M.W., Dhu, T. & Lam, N., 2004. Optimised array design for microtremor array studies applied to site classification; observations, results and future use, in *Proceedings of the 13th World Conference of Earthquake Engineering*, Vancouver.

Baggeroer, A.B., Kuperman, W.A. & Mikhalevsky, P.N., 1993. An overview of matched fields methods in ocean acoustics, *IEEE J. Ocean. Eng.*, **18**(4), 401–424.

Bensen, G.D., Ritzwoller, M.H., Barmin, M.P., Levshin, A.L., Lin, F., Moschetti, M.P., Shapiro, N.M. & Yang, Y., 2007. Processing seismic ambient noise data to obtain reliable broad-band surface wave dispersion measurements, *Geophys. J. Int.*, **169**, 1239–1260.

Bettig, B., Bard, P.-Y., Scherbaum, F., Riepl, J., Cotton, F., Cornou, C. & Hatzfeld, D., 2001. Analysis of dense array noise measurements using the modified spatial auto-correlation method MSPAC—application to the Grenoble area, *Bollettino di Geofisica Teorica e Applicata*, **42**, 281–304.

Burg, J.P., 1964. Three-dimensional filtering with an array of seismometers, *Geophysics*, **29**, 693–713.

Campillo, M., 2006. Phase and correlation in ‘random’ seismic fields and the reconstruction of the Green function, *Pure Appl. Geophys.*, **163**, 475–502.

Campillo, M. & Paul, A., 2003. Long-range correlations in the diffuse seismic coda, *Science*, **299**, 547–549.

Capon, J., 1969. High-Resolution frequency-wavenumber spectrum analysis, *Proc. IEEE*, **57**(8), 1408–1418.

Chavez-Garcia, F.J., Rodriguez, M. & Stephenson, W.R., 2005. An alternative approach to the SPAC analysis of microtremors: exploiting stationarity of noise, *Bull. seism. Soc. Am.*, **95**(1), 277–293.

Cho, I., Tada, T. & Shinozaki, Y., 2004. A new method to determine phase velocities of Rayleigh waves from microseisms, *Geophysics*, **69**(6), 1535–1551.

Cho, I., Tada, T. & Shinozaki, Y., 2006. A generic formulation for microtremor exploration methods using three-component records from a circular array, *Geophys. J. Int.*, **165**, 236–258.

Cornou, C., Ohrnberger, M., Boore, D.M., Kudo, K. & Bard P.-Y., 2006. Using ambient noise array techniques for site characterisation: results from an international benchmark, in *Proc. 3rd Int. Symp. on the Effects of Surface Geology on Seismic Motion, Grenoble, 29 August – 01 September, 2006*, eds Bard, P.Y., Chaljub, E., Cornou, C., Cotton, F. and Guéguen, P., LCPC Editions, Paper NBT.

Derode, A., Larose, E., Campillo, M. & Fink, M., 2003. How to estimate the Green’s function of a heterogeneous medium between two passive sensors? Application to acoustic waves, *Appl. Phys. Lett.*, **83**(15), 3054–3056.

Di Giulio, G., Cornou, C., Ohrnberger, M., Wathelet, M. & Rovelli, A., 2006. 2-D small aperture arrays for velocity profiles estimation using ambient seismic noise in a small-size alluvial basin (Colfiorito, Italy), *BSSA*, **96**(5), 1915–1933.

Duvall, T.L., Jefferies, S.M., Harvey, J.W. & Pomerantz, M.A., 1993. Time distance helioseismology, *Nature*, **362**, 430–432.

Fäh, D., Stamm, G. & Havenith, H.-B., 2006. Analysis of ambient vibration array data in a blind test; Love and Rayleigh wave phase velocities and ellipticity, in *Proc. 3rd Int. Symp. on the Effects of Surface Geology on Seismic Motion, Grenoble, 29 August – 01 September, 2006*, eds Bard, P.Y., Chaljub, E., Cornou, C., Cotton, F. and Guéguen, P., LCPC Editions, Paper N07.

Gilles, P.M., Duvall, T.L., Scherrer, P.H. & Bogart, R.S., 1997. Subsurface flow of material from the Sun equator’s to its poles, *Nature*, **390**, 63–64.

Gaffet, S., Larroque, C., Deschamps, A. & Tressols, F., 1998. A dense array experiment for the observation of waveform perturbations, *Soil Dyn. Earthq. Eng.*, **17**, 475–484.

Gouédard, P., Roux, P. & Campillo, M., 2006. Dispersion curves measurement using noise cross-correlation technique, in *Proc. 3rd Int. Symp. on the Effects of Surface Geology on Seismic Motion, Grenoble, 29 August – 01 September, 2006*, eds Bard, P.Y., Chaljub, E., Cornou, C., Cotton, F. and Guéguen, P., LCPC Editions, Paper N06.

Henstridge, J.D., 1979. A signal processing method for circular arrays, *Geophysics*, **44**, 179–184.

Hisada, Y., 1994. An efficient method for computing Green’s functions for a layered halfspace with sources and receivers at close depths, *Bull. seism. Soc. Am.*, **84**, 1456–1472.
Hisada, Y., 1995. An efficient method for computing Green's functions for a layered halfspace with sources and receivers at close depths (Part 2), *Bull. seism. Soc. Am.*, **85*(4), 1080–1093.

Horike, M., 1985. Inversion of phase velocity of long-period microtremors to the S-wave-velocity structure down to the basement in urbanized areas, *J. Phys. Earth*, **33**(2), 59–96.

Kanai, K. & Tanaka, T., 1954. Measurement of microtremor, *Bull. Earth. Res. Inst.*, **32**, 199–209.

Köhler, A., Ohrnberger, M., Scherbaum, F., Wathelet, M. & Cornou, C., 2007. Assessing the reliability of the modified three-component spatial autocorrelation technique, *Geophys. J. Int.*, **168**(2), 779–796.

Lacoss, R.T., Kelly, E.J. & Toksöz, M.N., 1969. Estimation of seismic noise structure using arrays, *Geophysics*, **34**(1), 21–38.

Larose, E., Derode, A., Campillo, M. & Fink, M., 2004. Imaging from one-bit correlation of wide-band diffuse wavefield, *J. Appl. Phys.*, **95**, 8393–8399.

Larose, E., Lobkis, O.I. & Weaver, R.L., 2006a. Passive correlation imaging of a buried scatterer, *J. Acoust. Soc. Am.*, **119**, 3549–3552.

Larose, E. et al., 2006b. Correlation of random wavefields: an interdisciplinary review, *Geophysics*, **71**(4), S11–S121.

Lin, F., Ritzwoller, M.H., Townend, J., Bannister, S. & Savage, M., 2007. Ambient noise Rayleigh wave tomography of New Zealand, *Geophys. J. Int.*, **170**, 649–666.

Louie, J., 2001. Faster, better: shear-wave velocity to 100 meters depth from refraction microtremor arrays, *Bull. seism. Soc. Am.*, **91**, 347–364.

Miyakoshi K., 1996. A range of wavelengths possible to estimate phase velocities of surface waves in microtremors, in *Proc. of the 94th SEGJ Conf., Society of Exploration Geophysicists*, Japan, pp. 178–182 (in Japanese).

Moczo, P. & Kristek, J., 2002. "FD code to generate noise synthetics", *Sesame report* D02.09, available on http://SESAME-FP5.obs.ujf-grenoble.fr.

Nakamura, Y., 1989. A method for dynamic characteristics estimation of subsurface using microtremors on the ground surface, *Quart. Rept. Railway Tech. Res. Inst. Tokyo*, **30**, 25–33.

Nogoshi, M. & Igarashi, T., 1971. On the amplitude characteristics of microtremor (part 2) (in Japanese with english abstract), *J. Seismol. Soc. Jpn.*, **24**, 26–40.

Ohori, M., Nobata, A. & Wakamatsu, K., 2002. A Comparison of ESAC and FK methods of estimating phase velocity using arbitrarily shaped microtremor arrays, *Bull. seism. Soc. Am.*, **92**(6), 2323–2332.

Okada, H., 2003. "*The Microseismic Survey Method*", *Society of Exploration Geophysicists of Japan*, Geophysical Monographs Vol 12, translated by Koya Suto, Society of Exploration Geophysicists, Tulsa.

Okada, H., 2006. Theory of efficient array observations of microtremors with special reference to the SPAC method, *Explor. Geophys.*, **37**, 73–85.

Ohrnberger, M., 2005. Report on the FK/SPAC capabilities and limitations, *SESAME Deliverable D19.06*, pp. 43, http://sesame-fp5.obs.ujf-grenoble.fr/Delivrables/Del-D19-Wp06.pdf

Pedersen, H., Krüger, F. & the SVEKALAPKO Seismic Tomography Working Group, 2007. Influence of the seismic noise characteristics

on noise correlations in the Baltic shield, *Geophys. J. Int.*, **168**, 197–210.

Roux, P. & Kuperman, W.A., 2004. Extracting coherent wavefronts from acoustic ambient noise in the ocean, *J. Acoust. Soc. Am.*, **116**(4), 1995–2003.

Sabra, K.G., Gerstoft, P., Roux, P. & Kuperman, W.A., 2005. Surface wave tomography from microseisms in southern California, *Geophys. Res. Lett.*, **32**, L14311.

Sánchez-Sesma, F.J. & Campillo, M., 2006. Retrieval of the Green function from cross correlation: the canonical elastic problem, *Bull. seism. Soc. Am.*, **85**, 269–284.

Sánchez-Sesma, F.J., Perez-Ruiz, J., Campillo, M. & Luzón, F., 2006. Elastodynamic 2-D Green function retrieval from cross-correlation: canonical inclusion problem, *Geophys. Res. Lett.*, **33**, L13305.

Satoh, T., Kawase, H., Iwata, T., Higashi, S., Sato, T., Irikura, K. & Huang, H.-C., 2001. A-wave velocity structure of the Taichung basin, Taiwan, estimated from array and single-station records of microtremors, *Bull. seism. Soc. Am.*, **91**(5), 1267–1282.

Shapiro, N.M. & Campillo, M., 2004. Emergence of broadband Rayleigh waves from correlations of the ambient seismic noise, *Geophys. Res. Lett.*, **31**, L07614.

Shapiro, N.M., Campillo, M., Stehly, L. & Ritzwoller, M.H., 2005. High-resolution surface wave tomography from ambient seismic noise, *Science*, **307**, 1615–1618.

Stehly, L., Campillo, M. & Shapiro, N.M., 2007. Travel time measurements from noise correlation: stability and detection of instrumental time-shifts, *Geophys. J. Int.*, **171**, 223–230.

Tokimatsu, K., 1997. Geotechnical site characterization using surface waves, in *Proc. 1st Intl. Conf. Earthquake Geotechnical Engineering*, ed. Ishihara, Balkema, pp. 1333–1368.

Wapenaar, K., 2004. Retrieving the elastodynamic Green's Function of an arbitrary inhomogeneous medium by cross-correlation, *Phys. Rev. Lett.*, **93**, 254301.

Wathelet, M., Jongmans, D. & Ohrnberger, M., 2004. Surface wave inversion using a direct search algorithm and its application to ambient vibration measurements, *Near Surface Geophys.*, **2**, 211–221.

Wathelet, M., Jongmans, D., Ohrnberger, M. & Bonnefoy-Claudet, S., 2007. Array performances for ambient vibrations on a shallow structure and consequences over vs inversion, *J. Seismol.*, in press.

Weaver, R.L. & Lobkis, O.I., 2001. Ultrasonics without a source: thermal fluctuation correlations at MHz frequencies, *Phys. Rev. Lett.*, **87**(13), 134301.

Weaver, R.L. & Lobkis, O.I., 2005. The mean and variance of diffuse field correlations in finite bodies, *J. Acoust. Soc. Am.*, **118**, 3447–3456.

Woods, J.W. & Lintz, P.L., 1973. Plane waves at small arrays, *Geophysics*, **38**, 1023–1041.

Yao, H., Van Der Hilst, R.D. & de Hoop, M.V., 2006. Surface-wave array tomography in SE Tibet from ambient seismic noise and two-station analysis – I. Phase velocity maps, *Geophys. J. Int.*, **166**, 732–744.