# Off-great-circle propagation of intermediate-period surface waves observed on a dense array in the French Alps

N. Cotte, H. A. Pedersen, M. Campillo, V. Farra

► To cite this version:

N. Cotte, H. A. Pedersen, M. Campillo, V. Farra. Off-great-circle propagation of intermediate-period surface waves observed on a dense array in the French Alps. Geophysical Journal International, 2000, 142, pp.825-840. 10.1046/j.1365-246X.2000.00187.x. insu-03606311

HAL Id: insu-03606311

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

Submitted on 11 Mar 2022

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
# Off-great-circle propagation of intermediate-period surface waves observed on a dense array in the French Alps

N. Cotte,<sup>1,2,*</sup> H. A. Pedersen,<sup>1</sup> M. Campillo,<sup>1</sup> V. Farra<sup>3</sup> and Y. Cansi<sup>2</sup>

Accepted 2000 April 3. Received 2000 January 17; in original form 1999 July 1

## SUMMARY

Array analysis is performed on surface waves recorded in the French Alps using a small-aperture (25 km) temporary array of six broad-band stations. The analysis shows that both Rayleigh and Love waves deviate relative to the great-circle path. The deviations are particularly strong, up to 30°, between 20 and 40 s period. To interpret these observations, we first study the effect of large-scale structures using ray tracing in a smooth, laterally heterogeneous model of the Earth. Second, we evaluate the local effect by considering a model for the French Alps including strong lateral heterogeneities around the array that were not taken into account in the ray tracing. By combining these two possible causes of the observed deviations, we propose an explanation for the general trend in the observed deviations. Finally, we show that by taking into account azimuthal deviations, phase velocities measured at a regional scale can be significantly improved.

**Key words:** array analysis, azimuthal deviation, French Alps, phase velocity, surface waves.

## INTRODUCTION

The lateral heterogeneities in the lithosphere induce refraction and reflection of intermediate-period surface waves due to changes in their phase velocities (e.g. Oliver 1962) and their eigenfunctions. Such effects are particularly strong at continental margins, across which the crustal thickness varies considerably. McGarr (1969) studied amplitude variations of 20 s Rayleigh waves and explained the amplitude variations as being caused by the focusing and defocusing due to lateral heterogeneities. Sobel & von Seggern (1978) observed the same type of anomalies at LASA for surface waves at 20 s period. They compared their observations to models calculated by ray tracing so as to explain off-great-circle propagation and phenomena of focusing and multipathing. Capon (1970) has also shown strong evidence of these lateral variations in the 20–40 s period range for Rayleigh waves, determining multipathing as observed at LASA.

Several studies have been dedicated to the measurement of off-great-circle propagation using polarization analysis. Among the most recent studies, Lerner-Lam & Park (1989) studied the frequency-dependent polarization of 10–100 s surface waves in the Western Pacific. They confirmed the presence of significant refraction and multipathing. On data recorded in Iberia,

Paulssen et al. (1990) performed a time- and frequency-dependent polarization analysis and suggested that the anomalous surface wave observations could be due to the interference of surface wave signals. More recently, Levshin et al. (1994) studied polarization anomalies observed in the former Soviet Union for surface waves crossing northern and central Eurasia, and showed how a regional structure showing lateral heterogeneities in the crust and upper mantle can cause large deviations from great-circle propagation. Using beam-forming analysis on data from a temporary array in the Netherlands, Germany and Belgium, Alsina & Snieder (1996) measured deviations from great-circle propagation of up to 30° at 40 s period and showed that waves that had propagated across the Tornquist–Teisseyre Zone had particularly large deviations.

With intermediate-period surface waves recorded on the Kyrgyzstan broad-band array, Pavlis & Mahdi (1996) found different characteristics in the propagation of surface waves for three different seismic events (two earthquakes and a nuclear explosion) located at different azimuths and epicentral distances. One path showed that the ray propagates in a rather laterally homogeneous medium and the two others showed evidence of complicated multipathing and scattering effects.

Due to poor long-period and broad-band station coverage, there are few results available concerning surface wave propagation in the Alps. Knopoff et al. (1966) showed deviations from the great-circle path in the Alps for intermediate-period Rayleigh waves. By performing a three-station analysis, they

<footnote n="1"><sup>1</sup>Laboratoire de Géophysique Interne et Tectonophysique, Observatoire de Grenoble, Université Joseph Fourier, Grenoble, France</footnote>

<footnote n="2"><sup>2</sup>Laboratoire de Détection et Géophysique, Bruyères-le-Châtel, France</footnote>

<footnote n="3"><sup>3</sup>Institut de Physique du Globe, Paris, France</footnote>

<footnote>*Now at: IGPP–UCSD 0225, SIO, La Jolla, CA 92093-0225, USA.
E-mail: ncotte@mahi.ucsd.edu</footnote>
showed that deviations can reach up to 12.5° at 15 s period for a wave incident from the North Atlantic Ocean. The stations they used were far from each other, separated by at least 200 km, and most of the stations were located out of the central Alps.

We collected data over the period August–December 1996, using six broad-band stations installed in the French Alps near the town of Briançon as part of the GéoFrance3D programme (e.g. Paul *et al.* 1998). The six stations formed a circle with a diameter of 25 km, much smaller than the wavelengths of intermediate-period surface waves. The surface waves recorded by such a mini-array are therefore sufficiently coherent to be analysed using methods assuming high coherency. The small aperture means that the wave fronts can be considered planar and the choice of using teleseismic events ensures that the fundamental modes are well separated from higher modes.

We aimed to measure the polarization from data recorded by these stations, either by the method described in Keilis-Borok (1989) or by that of Roberts & Christoffersson (1990), but our results showed large standard deviations for the measurements. We therefore chose to measure the slowness vector within the array to make it possible to deduce the angle of arrival of the incident waves. Large deviations were identified that we interpret as a combination of global and local effects. To study the global effect, we carried out ray tracing in the global phase velocity models. Local effects were studied by numerical simulations and a simplified ray approach, taking into account changes in Moho depth. Through a combination of the effects of global and local heterogeneities the observed deviations can be explained fairly well. Finally, we show that by taking into account azimuthal deviations, phase velocities measured at a regional scale can be significantly improved.

## THE GÉOFRANCE3D PROGRAMME AND THE BROAD-BAND STATIONS IN THE FRENCH ALPS

The first phase of the GéoFrance3D program (e.g. Fréchet *et al.* 1998; Paul *et al.* 1998) took place in 1996. The main objective of that phase was to study the mechanism of the collision that led to the formation of the Western Alps. A multidisciplinary approach was carried out using passive seismology as well as geodesy and geology.

### Description of the broad-band seismic experiment

The six broad-band stations formed a mini-array as shown in Fig. 1 (the station locations are given in Table 1). They were installed from August to December 1996. The stations were equipped with Güralp CMG3 100 s sensors and were located

<table><caption>**Table 1.** Locations of the temporary broad-band stations installed for the GéoFrance3D experiment shown in Fig. 1. The station OGAG is a permanent station of the Observatory of Grenoble.</caption><thead><tr><td>Station</td><td>Longitude E (°)</td><td>Latitude N (°)</td><td>Elevation (m)</td></tr></thead><tbody><tr><td>ARVD</td><td>6.751</td><td>44.764</td><td>1680</td></tr><tr><td>CERD</td><td>6.725</td><td>44.859</td><td>1740</td></tr><tr><td>FRED</td><td>6.536</td><td>44.761</td><td>1420</td></tr><tr><td>PUYD</td><td>6.614</td><td>44.886</td><td>1440</td></tr><tr><td>SCRD</td><td>6.618</td><td>44.716</td><td>1370</td></tr><tr><td>VALD</td><td>6.481</td><td>44.859</td><td>1630</td></tr><tr><td>OGAG</td><td>6.541</td><td>44.788</td><td>1250</td></tr></tbody></table>

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Map of the locations of the six broad-band stations installed in the French Alps for the GéoFrance3D programme from August to December 1996. The temporary mini-array (black squares) describes a circle with a diameter of approximately 25 km. The black triangle shows the permanent station OGAG of the Observatory of Grenoble.</figcaption></figure>
on a circle with a diameter of approximately 25 km. In addition, the array included the permanent station OGAG of the Observatory of Grenoble (Coutant *et al.* 1999). All stations recorded continuously with a sampling rate of 125 Hz.

The configuration of the array was chosen to be able to determine the direction of the incident surface waves using methods for which the coherency between waves recorded at different stations needs to be high. Such an analysis therefore requires the array to be of small aperture, which also has the advantage that no station is needed in the centre of the array as there is no spatial aliasing at the periods of interest.

## Data selection

36 teleseismic events were recorded with a high signal-to-noise ratio by at least three of the six stations of the array. The event locations are plotted in Fig. 2 using a polar representation, along with the phase velocity model calculated by Ekström *et al.* (1997) (referred as the ETL model) for the 35 s Rayleigh wave. The distance from the centre gives the epicentral distance directly, between 0° and 180°. The azimuthal coverage was sufficient to perform an analysis of surface wave deviation from the great-circle path as a function of both epicentral distance and theoretical backazimuth.

We deconvolved the data by instrument response so as to be able to use data from the permanent station (OGAG), which

is equipped with a different sensor (Streckeisen STS-2). Data were filtered between 5 and 100 s and decimated to obtain a sampling rate of 1.024 s.

In Fig. 3(a), we show the vertical component of seismograms for an event in Japan on the 1996 October 19 at 14:44:40.7 UT, which was recorded by five of the six temporary broad-band stations and by OGAG. The main phase is the Rayleigh wave. In Fig. 3(b), we show seismograms for an event in Peru on 1996 November 12 at 16:59:44.0 UT, recorded on the transverse component by four of the six temporary broad-band stations. The main phase is the Love wave. For each event the seismograms are very similar for different stations, in both amplitude and phase. The similarity confirms the validity of analysing the data using methods that are based on the assumption of a coherent field.

## ARRAY ANALYSIS: PROCESSING AND RESULTS

In this section we present the polarization and array analysis used to measure the deviations from a great-circle path in the French Alps for both Rayleigh and Love waves.

### Polarization

When there is only one dominant surface wave mode it is possible to measure the propagation directions of the waves by polarization analysis. This has the advantage that only one three-component station is needed. We used two different methods. First, we used the method of Keilis-Borok (1989) to analyse particle motion. This method uses the eigenvectors of the covariance matrix to estimate arrival angles of Rayleigh and Love waves in different time and frequency windows. The quality of the measurement in each time-frequency window is estimated by the energy and a quality factor. Results obtained in some continental areas offer very stable measurements (e.g. Levshin *et al.* 1994), but in our case the polarization measurements changed rapidly as a function of time and frequency, so it was very difficult—often impossible—to deduce any particular arrival angle. The main input parameter, which is the frequency filter width, was varied within a broad range, but no satisfactory results were obtained in terms of the stability of the polarization. The second method we used was that proposed by Roberts & Christoffersson (1990), in which polarization measurements are achieved either by a maximum likelihood analysis or by a transfer function analysis applied to the covariance matrix. The analysis was performed on a moving time window, and data were filtered in different frequency bands prior to analysis. The window length and the *a priori* frequency filtering play an important role in the analysis, but the best results (typically with a window length twice the maximum period of the frequency filter) still showed uncertainties of at least 15°–20° using the 90 per cent confidence interval, and in several cases no stable direction could be measured.

To further understand why polarization measurements were not successful when applied to data, we performed the same analysis on synthetic seismograms for Rayleigh waves propagating in the French Alps. The model and the method for calculating the synthetic seismograms are more explicitly stated in the subsection ‘Local effect’, where all types of wave coupling are taken into account. Even on the synthetic

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Locations of the 36 events (white circles) used for the array analysis in a polar representation. The greyscale in the background shows the phase velocity model (Ekström *et al.* 1997) for the 35 s period Rayleigh wave. The representation is centred on the centre of the array GEOF (see Fig. 1) and the distance from the centre gives the epicentral distance between 0 and 180° directly.</figcaption></figure>
961019 Japan - Vertical component

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

961112 Peru - Transverse component

<figure><img src="image_5.png" /><figcaption>**Figure 3.** (a) Seismograms (vertical component) of an event in Japan (1996 October 19 at 14:44:40.7 UT) recorded by five of the six temporary stations and by the permanent station OGAG of the Observatory of Grenoble. The dominant phase is the Rayleigh wave. (b) Seismograms (transverse component) of an event in Peru (1996 November 12 at 16:59:44.0 UT) recorded by four of the six temporary stations. The dominant phase is the Love wave. For both (a) and (b) the amplitude scale is the same for all traces. Signals shown here are filtered between 5 and 100 s and are deconvolved by instrument response.</figcaption></figure>
seismograms, the polarization analysis was unstable. In the case of the synthetic seismograms, where a single-mode Rayleigh wave is incident on the model, coupling of the Rayleigh wave with other wave types, including reflected the Rayleigh wave, distorts the waveforms slightly, apparently enough to make the polarization measurement unstable. For data from the array, a further complication is that the Love wave coda and the Rayleigh wave overlap in time, and do not necessarily have perfectly perpendicular displacements in the horizontal plane. We concluded therefore that polarization methods are not well adapted to our data for the desired precision of no more than 2°. However, due to the small aperture of the array, waveforms were very coherent from one station to the other, making it possible to analyse data by array analysis based on the high coherence.

## Array analysis and processing

The method we used is adapted from that described by Barker *et al.* (1996) and originally proposed by Poupinet *et al.* (1984) for short-period non-dispersive signals. The idea was to measure the time delay between pairs of stations using a moving time window. As we were analysing dispersive signals we adapted the method for calculating the average time delay for different frequency bands, rather than using the slope of the phase of the cross-spectrum. Signals were filtered in narrow frequency bands and the phase of the cross-spectrum averaged within each frequency band, neglecting the dispersion. We thus obtained an average time delay corresponding to a narrow frequency band, using a weighting function determined by both the spectral amplitude and the coherency, as proposed by Poupinet *et al.* (1984). We bandpass filtered the data in three different frequency ranges—0.010–0.025, 0.025–0.050 and 0.050–0.2 Hz; narrower frequency filters would have caused the signal to become almost stationary. The maximum spectral amplitude is located in the frequency band 0.025–0.050 Hz. As suggested by Barker *et al.* (1996), we therefore chose time windows twice the length of the longest period for each frequency filter, except for the first one, for which the time window was 150 s long. Neighbouring time windows overlap by 75 per cent.

The time delay, $\Delta t_i$, measured for a pair of stations $i$ depends linearly on the slowness vector $\mathbf{p}(p_n, p_e)$ across the array, defined by the predicted time delay, $\partial t_i / \partial t_j = -p_n \Delta x_n^i - p_e \Delta x_e^i$, where $\Delta x_n^i$ and $\Delta x_e^i$ are the interstation distances along the north and east directions between the two stations. The slowness vector is estimated using a damped least-squares inversion. The measured backazimuth, $BAz$, of the incident wave in degrees from the north is given for each time window by $BAz = \tan^{-1}(p_e/p_n)$ (adding 180° if $p_e < 0$), and the phase velocity is $c = (\sqrt{p_n^2 + p_e^2})^{-1}$.

The propagation directions of Rayleigh waves are measured on the vertical component and those of the Love waves are measured on the transverse component. Deviations from great-circle propagation of Rayleigh waves may have contaminated the transverse component, but this did not induce significant errors in our measurements as the Love and Rayleigh waves were well separated in time for the epicentral distances considered. Deviations of Love waves from great-circle propagation showed that Love waves were present in both the radial and the transverse components. However, this rotation influenced only the amplitude and not the phase of these waves, so our method of analysis remains valid.

## Results for Rayleigh and Love waves

In Fig. 4 we show the results of the analysis for the signal bandpasssed between 20 and 40 s for the data presented in Fig. 3(a). The analysis of the Rayleigh waves is presented in Fig. 4(a), and that of the Love waves in Fig. 4(b). Above each figure the seismogram for the station PUYD (see Fig. 1) is shown. For each time window the phase velocity (top), the

<figure><img src="image_6.png" /><figcaption>961019 Japan - Transverse component - Bandpasssed 20-40s</figcaption></figure>

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

**Figure 4.** Results of the array processing between 20 and 40 s period for the vertical (a) and transverse (b) components. Records used are for an event in Japan. The equivalent seismograms recorded by the station PUYD (see Fig. 1) are shown at the top of each figure. The first part of each figure shows the phase velocity in km s⁻¹ as a function of the time. The second shows the observed backazimuth. The theoretical backazimuth is shown by the solid line at 44.3°. The third shows the rms error as a fraction of the data sampling interval, which is 1.024 s. The solid circles are plotted at the centre of each time window.
backazimuth (centre) and the rms error expressed as a fraction of the sampling rate (bottom) are shown as points. The points are located in the centre of the equivalent time windows, and because of the 75 per cent overlap of the time windows every four points are independent. The phase velocities that we obtain for individual events are somewhat unstable and parameter-dependent, even though they are reasonable for the period range corresponding to the *a priori* frequency filtering. However, the observed backazimuth is very stable within the time windows of high amplitude. As we show below (see the section entitled 'Consequences of azimuthal deviations on phase velocity measurements'), including several events makes it possible to measure phase velocities within the array.

For the example presented in Fig. 4, the theoretical backazimuth, which is shown as a continuous line, is 44.3°. For both Rayleigh and Love waves this analysis shows that the wave is incident with a smaller backazimuth than is expected. We estimate the backazimuth as the average backazimuth in the time interval of peak amplitude. The difference between the observed and theoretical values is thus -14° for the Rayleigh wave and -10° for the Love wave. The incident wave propagating through the array therefore shows a direction of propagation oriented more north-south than the direction given by the theoretical backazimuth of the source as seen at the station. Indeed, the ray deviates to propagate in the medium with the highest phase velocities. This means that for the example previously discussed, the surface waves coming from Japan propagate north of the great circle between the station and the source, avoiding the Tibetan plateau, which is a low-velocity structure for the period under consideration.

For all 36 events, we performed the analysis for Rayleigh and Love waves in the three frequency bands with the condition that the signal-to-noise ratio was greater than 4. All results for the three frequency bands are shown in Fig. 5 from long (top) to short (bottom) periods using a polar representation. Results for Rayleigh waves are given in the left column and those for Love waves in the right column. The distance from the centre is the epicentral distance between 0° and 180°. Open circles correspond to an observed backazimuth smaller than the theoretical backazimuth (i.e. an anticlockwise deviation) and plus signs indicate larger backazimuths (i.e. a clockwise deviation). The sizes of the signs increase with increasing deviation from the theoretical backazimuth.

Our results show the existence of large angular differences between the observed and theoretical backazimuths. In some cases deviations are up to 30°. Note some differences between the three frequency bands, and also between Rayleigh and Love waves. In each case some strong geographical correlations exist; for example, all measurements corresponding to the events from southeast Asia show negative deviations for a backazimuth between 40° and 50° and an epicentral distance of approximately 80°–90°. However, the deviations are also a function of epicentral distance. For example, for periods longer than 40 s for the Love wave (top right in Fig. 5) and a backazimuth of 45° for different distances (85° and 135°), the measured deviations are of opposite sign.

## INTERPRETATIONS OF THE OBSERVED AZIMUTHAL DEVIATIONS

The interpretation of the observed deviations should ideally take into account the full wave propagation between the source

and the receiver. As the heterogeneities are of varying size as compared to the wavelength, an approximate method cannot be justified. However, as the global models are very smooth, ray tracing provides a reliable estimation of the effect of the large-scale structures for the off-great-circle path deviations. Despite the strong heterogeneity of the Alps, the detailed base of knowledge about their structure allows us to evaluate the effect of the local heterogeneities using full-wavefield modelling and simple considerations based on Snell's law. Finally, combining the two effects gives a good correspondence to the experimental data.

## The global effect: ray tracing in earth models

For periods shorter than 100 s, there are several isotropic models for Love and Rayleigh waves available; the most recent are those of Ekström *et al.* (1997) (see Fig. 2), Laske & Masters (1996) and Trampert & Woodhouse (1995). Others models such as CRUST 5.1 (Mooney *et al.* 1998) and 3SMAC (Nataf & Ricard 1996) are determined by *a priori* constraints on the earth structure.

First, we studied the multipathing effect using a ray tracer based on a shooting approach combined with a 2-D finite element ray tracer (Farra 1990) in Cartesian coordinates. The code is able to find all the different rays between the given source and receiver; as a result, waves propagating along different paths can be evaluated. The amplitudes of each of these waves were calculated. We show in Fig. 6(a) the rays calculated by this method within the 3SMAC model for a 35 s Rayleigh wave. Three solutions were found between the source (S) located in Taiwan (21.9°N; 121.5°E) and the receiver (R) located in the French Alps (44.8°N; 6.6°E). The model used is a Mercator transformation (Jobert & Jobert 1983) of the spherical velocity model where the great circle is rotated onto the equator to minimize distortion around the path caused by the transformation. The off-great-circle deviation measured on data for this station-source couple is -25° (with a large uncertainty) for periods between 20 and 40 s. In this example of ray tracing, off-great-circle deviations calculated at the station at the 35 s period are -11.9°, -9.7° and +1.8°. The three rays arrive simultaneously at the station (the differences in arrival times are smaller than 10 s and the total propagation time is 2640 s) and their amplitudes are similar. It is, however, possible to verify that the main part of the wave train is not composed of multiple arrivals. For this purpose we applied singular value decomposition and spectral matrix filtering to the data. Both of these methods showed that in our data only one wave was identified, so possible multipathing must be very small. This is not particularly surprising as ray theory can be successfully applied when the heterogeneities are much bigger than the wavelength (or represented by sharp boundaries), and this condition is not respected in the 3SMAC model (Fig. 6a). Phenomena such as wave front healing are not taken into account in the ray tracing.

Second, we also used the ray tracer presented by Laske & Masters (1996), working with spherical coordinates. It is based on the work of Woodhouse & Wong (1986), who described the formulation for solving the ray equation on a sphere. The authors proposed a linearized theoretical framework for interpreting off-great-circle arrival angles and amplitude anomalies. For a given frequency, the linear perturbation theory gives the relation between the tangent of the off-great-circle azimuth
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Observed deviations obtained for the three different bandpass filters from long periods (top) to short periods (bottom), with results for Rayleigh waves in the left column and those for Love waves in the right column. Open circles are for measured angles smaller than the expected ones and plus signs for larger angles.</figcaption></figure>

v(Δ) and the path integrals along the great circle:

$$v(\Delta) = - \frac{1}{\sin \Delta} \int_{0}^{\Delta} \sin \phi \frac{\partial}{\partial \theta} \left( \frac{\delta c(\frac{\pi}{2}, \phi)}{c_o} \right) d\phi, \quad (1)$$

where Δ is the epicentral distance, θ the colatitude, φ the longitude and c(θ, φ) the phase velocity. A shooting method is used to find the solution around the initial take-off angle whose tangent is v(Δ). The linear perturbation theory ignores multipathing; however, all the solutions due to multipathing

can be determined by the ray tracer by shooting in all possible directions.

In Fig. 6(b), we compare these two ray-tracing results within the 3SMAC model (Rayleigh wave, T=35 s), filtered using a maximum harmonic degree of 20 so that the structures are very smooth and multipathing does not occur. Small circles indicate the results using the method of Farra (1990) and crosses are those using the method of Laske & Masters (1996). In such a smooth model the two methods give very similar results with calculated deviations of −5.4° and −6.2°, respectively.
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Ray tracing in global models. (a) Ray tracing within the 3SMAC model (35 s period Rayleigh wave) by the 2-D ray tracer of Farra (1990). The source (S) is located in Taiwan and the receiver (R) in the Alps. Results clearly show multipathing as the three rays arrive simultaneously at the station (the difference in the arrival times is smaller than 10 s and the total propagation time is 2640 s). Off-great-circle path deviations are −11.9°, −9.7° and +1.8°. (b) Ray tracing within the 3SMAC model filtered for degrees smaller than 20. Superimposed we show the results provided by the method proposed by Laske & Masters (1996) as crosses. The two methods give respectively −5.4° and −6.2° of deviation from the great circle. (c) Ray tracing within the ETL model (Ekström *et al.* 1997) for the 35 s period Rayleigh wave. The two methods give respectively −4.2° and −4.0°.</figcaption></figure>

For comparison we performed the same calculations on
the model ETL (Ekström *et al.* 1997) for the 35 s Rayleigh
wave, for which the maximum harmonic degree is 40 but which
corresponds according to the authors to structures determined
by a harmonic degree of 20. The results are shown in Fig. 6(c).
The results from the two methods are again almost equal,
−4.2° and −4.0°, respectively. Discrepancies between the results

calculated within the two models are small and may be due to
interpolation of the model in Cartesian coordinates. We there-
fore chose to use the ETL model as it is a model calculated from
data and is rather smooth, so multipathing effects do not occur.
Ray tracing showed that the arrival angle given by the linear
perturbation theory (Woodhouse & Wong 1986) is sufficient
for calculating the arrival angle of the incident ray at the station.
Due to the simplicity of the calculation, we therefore used the linear perturbation theory to estimate the influence of the global structure of the observed off-great-circle deviation.

We then calculated the predicted backazimuth for a large number of source locations. Fig. 7 shows predicted deviations from the great-circle path for the 35 s period Rayleigh wave using a polar representation, where the greyscale indicates the value and the sign of the deviation. Note that for a source located close to the station (centre of the circle), the calculated deviation is small (±2°). Note also that the deviation depends on both epicentral distance and backazimuth. For example, for a source located with a backazimuth of 225°, the deviation is positive for epicentral distances smaller than 140° and negative for larger epicentral distances. The results obtained with the 3SMAC model for the same period are very similar. The only significant difference is for a source located with a backazimuth of 300°–330° and an epicentral distance between 40° and 140° (this area corresponds to the north of Canada, the west coast of the USA and the Pacific Ocean), for which deviations are positive in the 3SMAC model and negative in the ETL model.

In summary, ray tracing shows that surface waves recorded in the Alps can be strongly deviated by large-scale structures. Some calculated deviations are greater than 20° (see Fig. 7), but ray tracing alone fails to explain the observed deviations for some events.

<figure><img src="image_10.png" /><figcaption>**Figure 7.** Results of the estimation of the angular deviation of the incident Rayleigh wave from the great circle using the linear perturbation theory (Woodhouse & Wong 1986) for a large number of possible source locations. The model used is that of Ekström *et al.* (1997) for a 35 s period Rayleigh wave. The greyscale indicates the sign and value of the deviation calculated at the station.</figcaption></figure>

## Local effect: the influence of sharp lateral variations in the French Alps

The array was located in the Briançonnais zone, east of the Penninic Frontal Thrust and west of the Ivrea zone (e.g. Debelmas & Kerckhove 1980). As shown by Ménard & Thouvenot (1984), Bouguer anomalies present significant lateral variations through the Briançonnais zone, implying that the structure is laterally very heterogeneous. Ménard & Thouvenot (1984) interpreted the Bouguer anomalies by lithospheric flaking in the crust beneath the Briançonnais zone and by crustal thickening, because the Moho reaches a depth of approximately 50 km under this zone (Perrier 1979). Such a complex medium could have a strong influence on intermediate-period surface waves propagating through the Alps. To confirm this we first used numerical simulations to study surface wave refraction and diffraction effects induced by these lateral heterogeneities and then predicted refraction effects for the different events.

Fig. 8(a) shows the crustal model proposed by Ménard, adapted from Ménard & Thouvenot (1984). The structure is locally oriented N30°E. The crust is 30 km thick towards the northwest and 35 km thick towards the southeast, and increases to 50 km under the Briançonnais zone. Lithospheric flaking in this area is proposed based on gravimetric profiles (Ménard & Thouvenot 1984) and reflection seismics (ECORS-CROP Deep Seismic Sounding Group 1989). We performed forward modelling in this model using the indirect boundary element method (IBEM) in the formulation of Pedersen *et al.* (1996). In these simulations, a plane and single-mode surface wave is obliquely incident upon a multilayered and laterally heterogeneous 2-D structure. The computation is made in three dimensions and the coupling of all wave types is taken into account. The calculations were performed in the frequency domain and synthetic seismograms were obtained by multiplication with a source function followed by an inverse Fourier transform. The velocities used in our model were $V_P=6.4$ km s⁻¹ and $V_S=3.7$ km s⁻¹ in the crust, and $V_P=7.9$ km s⁻¹ and $V_S=4.55$ km s⁻¹ in the upper mantle. Densities were respectively 2.8 and 3.1 g cm⁻³.

We calculated synthetic seismograms for a Rayleigh wave incident on the structure from the north (i.e. 60° incidence). By adding the phase due to the propagation along strike to the linear profile used for numerical simulations, we reconstituted the array configuration that we have in the Alps, that is, six stations describing a circle (see Fig. 1). This was possible for this particular ('2.5-D') modelling because the wavefield is described by a phase delay that depends on the input parameters, that is, the phase velocity of the incident wave and its incidence angle. The array is located above the lithospheric flaking and slightly to the west of it. Synthetic seismograms were convolved with a Ricker wavelet with a dominant period of 30 s, corresponding to the period with maximum amplitude in the data.

The incident wave was set to be from the north, that is, with a backazimuth of zero. We performed the same array analysis on these synthetic seismograms as that described in the previous section. The results are shown in Fig. 8(b) for the 20–40 s band-pass. The analysis shows that the waves are slightly deviated to the west because the deviation is approximately −8° for a period of 30 s, and can reach up to −15° for shorter periods.

As the deviations are negative, this suggests that the average velocity under the array is smaller than that in the medium
<figure><img src="image_11.png" /><figcaption>**Figure 8.** (a) The crustal model used in the numerical simulations (after Ménard & Thouvenot 1984). The 2-D structure is oriented N30°E, with the Briançonnais zone characterized by the Moho plunge and lithospheric flaking. Velocities are $V_p=6.4$ km s$^{-1}$ and $V_s=3.7$ km s$^{-1}$ in the crust, and $V_p=7.9$ km s$^{-1}$ and $V_s=4.55$ km s$^{-1}$ in the upper mantle. Densities are respectively 2.8 and $3.1$ g cm$^{-3}$. Synthetic seismograms are calculated for a linear profile shown by triangles. (b) Example of array analysis that we performed on synthetic seismograms between 20 and 40 s period.</figcaption></figure>

outside the laterally heterogeneous model, as predicted by Snell's law. This indicates that the most important effect on wave deviation is the Moho depth, which is more important than the lithospheric flaking. Note that Moho depth is poorly taken into account in the previous ray tracing as the global model wavelengths are greater than the lateral variations in the Alps.

As it is not computationally practical to perform the numerical simulations for all of the events, we instead applied the simple approximation given by Snell's law based on the assumptions that the medium where the incident wave propagates corresponds to a 30 km thick crust and the medium in which the transmitted wave propagates corresponds to a 50 km thick crust, with velocities and densities as for the numerical simulations. For a wave incident from the north, the deviation is approxi-

mately -6°. This is slightly different from the deviation measured on the synthetic seismograms (-8°), but we consider that this first approximation is satisfactory, especially taking into account the uncertainties on the predicted angle in the global models.

We calculated the deviation given by the difference between the incident and the refracted angles for all incidence angles on the structure that correspond to recorded seismic events. To take into account the shape of the Alps, we approximated the structure by four linear parts, as shown in Fig. 9. The orientation chosen depended on the backazimuth of the incident wave. It was N75 for a backazimuth between 5° and 55°, N0 between 55° and 210°, N30 between 210° and 240° and N45 between 240° and 5°.
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Topographic map of the Alps. The orientation of the curved structure, used for calculating the refraction by Snell's law, depends on the backazimuth of the incident wave: N75 for a backazimuth between 5° and 55°, N0 between 55° and 210°, N30 between 210° and 240° and N45 between 240° and 5°.</figcaption></figure>

In summary, the local structure can induce deviations of up to 20°, but fails to explain even the main trends in the observed deviations.

## Combining global propagation and local effects

To combine the effects of lateral refractions of global propagation with those of the local structure, we first took the backazimuths calculated by ray tracing and then calculated the effect of the local structure by Snell's law. For example, for the source in Japan (event of 1996 October 19) presented in Figs 3(a) and 4(a), the theoretical backazimuth is 44.3°. The ray tracing gives a deviation of −3.2° for a station in the Alps. Thus, the incident wave on the local structure has a backazimuth of 41.1°. Finally, the angle calculated by Snell's law is equal to −8.0°, so the predicted backazimuth is 33.1°, i.e. a deviation of −11°.

This analysis was carried out for all events. Fig. 10 shows the results for Rayleigh and Love waves. Solid circles are for *T*=35 s and crosses are for *T*=50 s. Predicted deviations

are shown for the global model (Figs 10a and b), the local structure (Figs 10c and d) and the superposition of the two (Figs 10e and f). We see that using either of the two methods alone does not explain the observed deviations, as observations and predictions are not correlated and the range of predicted deviations is significantly smaller than the observed deviation. However, by combining global and local effects we predict deviations that are as large as the observed deviations and show a better fit to the data. The disadvantage of this type of representation is that all geographical information is lost.

To remedy this, Fig. 11 shows a comparison between predicted and observed deviations in polar plots. Observed deviations are shown in the left column and predicted deviations in the right column for Rayleigh waves (top) and Love waves (bottom). Observed deviations are given for periods between 20 and 40 s, with predictions for 35 s period waves. Grey shading highlights points for which the sign of the deviation is the same in the predictions and the observations.

For Rayleigh waves only four points of the 26 show a significant discrepancy. There is no geographical relationship
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Results of the deviations calculated for Rayleigh (left column) and Love waves (right column) as a function of the observed deviations. The solid circles are for *T*=35 s and the crosses for *T*=50 s. (a), (b) Predictions calculated in global models. (c), (d) Predictions calculated in the local structure. (e), (f) Total predictions as the superposition of the global and local structures.</figcaption></figure>

between these points and the difference may be due to the
structure around the source. For Love waves the fit is not as
good as for the Rayleigh waves, with six points out of 20 for
which the discrepancy is large. Note that four of these points
are located with a backazimuth between -40° and 40°. This

indicates that either the global model or the local structure
does not correct for sources located north of the array, as both
of them yield negative deviations for these events. For the
other two points there does not seem to be any geographical
relationship.
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Comparison between measured (left column) and predicted (right column) deviations at a station located in the French Alps. Top: Rayleigh waves; bottom: Love waves. Measurements are for signals between 20 and 40 s, and predictions are for 35 s period. Grey shading corresponds to points with a reasonable fit between observations and predictions.</figcaption></figure>

Note that in most cases we can explain the deviations for Rayleigh and Love waves as observed by the mini-array, with a greater similarity between the two types of waves than in the data. Even if we cannot explain the exact value of all the observed deviations, we can at least predict in most cases the sign of the deviation.

## CONSEQUENCES OF AZIMUTHAL DEVIATIONS FOR PHASE VELOCITY MEASUREMENTS

The array analysis demonstrates the existence of large deviations between the direction of propagation of the surface waves through the array and the theoretical backazimuth determined from the source and the station. Such large off-great-circle deviations may induce significant errors in phase velocity

measurements. We therefore measured phase velocities assuming that the waves have propagated along the great circle between the source and the receivers, so as to study the changes in phase velocity induced by correcting for the observed deviations.

### Phase velocity measurement within the mini-array

To obtain phase velocities, the phase difference is usually measured between stations that are separated by at least the wavelength of interest. This condition is not theoretically necessary as phase differences can be measured for small distances, but the influence of lateral refraction, multiple wave trains and noise can make measurements unstable.

Assuming that the wave front of the surface waves is plane through the array (Fig. 1), we calculated the distance between each pair of stations projected onto the propagation direction
of the wave. Two distances were calculated: that with the propagation direction along the great-circle direction, and that with the propagation direction as observed with the array analysis.

For each incident surface wave on the array, we calculated for each distance the cross-spectrum for each pair of stations with distances larger than 5 km, as well as its phases. The phase was slightly smoothed to remove small fluctuations, but sufficiently little to keep major fluctuations. The phase $\phi$ was then converted into a time delay, $\Delta t = \phi/(2\pi f)$. For each frequency we thus obtained several measurements of time delays corresponding to different distances. Taking those frequencies for which the coherence was greater than 0.98, we performed a linear regression of the (time-delay/distance) points, constraining the line to pass through the origin. The slope of the line indicates the phase velocity and the variance on the slope gives its uncertainty.

## Comparison between measured phase velocities without and with azimuth corrections

A very strict data selection had to be applied to obtain stable measurements. Only four events for the Rayleigh waves and six events for the Love waves were used in the analysis. Using the procedure described above, we obtained 21 interstation measurements for Rayleigh waves and 35 for Love waves. The resulting phase velocities are presented in Fig. 12. The left-hand figure shows the results for incident waves using distances not corrected for the off-great-circle propagation. Open circles represent Love waves and solid circles Rayleigh waves. The final dispersion curves were not smoothed so the velocity at each frequency is independent of its neighbour. The grey bars are the uncertainties on the phase velocity. These standard deviations are large for long periods for which the spectral amplitude is low and where very small errors in the phase measurement induce large errors in the time delay. In the 20–40 s range, standard deviations are smaller, especially for

Rayleigh waves, and their smoothness illustrates the stability of the phase velocity measurement. Measurements for periods shorter than 15 s show very large instabilities and could not be used for further analysis.

The right-hand side of Fig. 12 shows the phase velocities obtained using distances corrected by observed deviations. Note that without any azimuth correction, phase velocities are generally higher than after correction, except for Rayleigh waves at periods longer than 50 s. Without the azimuth correction, the calculated phase velocities fluctuate from one period to the next, while such fluctuations are strongly diminished by the correction, especially between 20 and 40 s period. The error bars are smaller when the phase velocity measurements are corrected for azimuth deviations. Note that even by applying the correction, the average value of the phase velocity for Love waves for periods greater than 40 s is high for a continent and must be considered with caution. All these observations lead us to conclude that phase velocities measured at a regional scale must take into account off-great-circle deviations.

## DISCUSSION AND CONCLUSIONS

Our results show significant off-great-circle deviations for intermediate-period fundamental mode surface waves propagating in the French Alps. Polarization analysis was not possible, probably due to the strong crustal heterogeneity in the Alps, so the deviations were obtained by array analysis. Our observations show the following.

1. For intermediate-period surface waves, off-great-circle deviations can reach up to 30° in the French Alps. These deviations vary strongly as a function of both backazimuth and epicentral distance.

2. The deviations can be very different for Love and Rayleigh waves for the same source, even if in most cases the deviations are of the same sign. Even if velocity maps clearly predict this observation, the values measured for the two waves can be important.

<figure><img src="image_15.png" /><figcaption>**Figure 12.** Phase velocities presented as a function of the period. The left-hand figure shows the phase velocities obtained without the azimuth correction, and the right-hand figure shows the phase velocities for incident waves corrected for the off-great-circle deviations. Open and solid circles are respectively for Love and Rayleigh waves. Grey bars show the standard deviation of the phase velocity measurements.</figcaption></figure>
Correcting for azimuth deviations improves the phase velocity measurements. Dispersion curves were smoother and had smaller error bars when the azimuth correction was applied. Furthermore, this correction can significantly change the phase velocities.

To interpret the observed deviations, we determined the effect of propagation in a global phase velocity model and the effect of the local structure of the Alpine arc. For the ray tracing, the choice of the models proposed by Ekström *et al.* (1997) for Love and Rayleigh waves for 35 and 50 s periods did not introduce big differences compared to other models (e.g. 3SMAC). The ray tracing in a global model shows that deviations of up to 20° can be expected, even if the lateral variations in the models are very smooth compared to the wavelength (Fig. 7). For example, the rays tend to avoid large structures with low phase velocities such as northern and central Eurasia (e.g. Levshin *et al.* 1994; Pavlis & Mahdi 1996).

We also calculated the refraction and diffraction effects generated by a local model proposed for the Briançonnais zone (Ménard & Thouvenot 1984). The proposed model shows very large lateral heterogeneities a little smaller than the wavelengths in the data, including variations in Moho depth as well as lithospheric flaking in the crust (Fig. 8a). Array analysis performed on synthetic seismograms shows that for a Rayleigh wave incident at 60° on the structure, the waves can be refracted up to 15° (Fig. 8b). In spite of complex diffraction effects, Snell’s law provides a good approximation of the local effects. The model proposed for calculating the deviations only shows lateral variations in Moho depth, from 30 km outside to 50 km within the Alps. As the Alpine structure has a significant curvature, we adapted the strike of the structure according to the azimuth of the incident wave (Fig. 9). The locally produced refractions are up to 20° for Rayleigh waves at 30 s period.

It is necessary to combine the global propagation and the local effects to obtain a reasonable fit between predicted and observed deviations, as shown in Fig. 11. Discrepancies are due to the neglect of structure around the source, to the oversimplified model used for the Alps, and to the smoothing of the global models, which annihilates diffractions on the structures and multipathing. It is difficult to evaluate the degree to which these cause the discrepancies between predictions and observations.

Our observations show that one should be careful in the measurement of phase velocities on a regional scale between pairs of stations. Although it is mostly assumed that the surface wave propagates along the great circle, we obtained off-great-circle deviations up to 30°. Thus, even if stations are properly aligned on the great circle, the error in the phase velocity can be significant.

## ACKNOWLEDGMENTS

Numerical simulations presented in this paper were performed at the Service Commun de Calcul Intensif de l'Observatoire de Grenoble (SCCI). We thank G. Laske and G. Masters for providing the ray tracer and Y. Capdeville, J.-P. Montagner and S. Chevrot for communicating the global models. We also wish to thank O. Coutant for providing the data from the Observatoire de Grenoble, H.-C. Nataf and P. Cardin for their help, and R. Roberts for providing the polarization program. We are indebted to G. Ménard for fruitful discussions

and elaboration of a simplified model for the simulations, to R. Bark for a critical reading of the manuscript, and to J. Fréchet for encouragement with broad-band experiments in the Alps and participation in fieldwork. We also thank two anonymous reviewers for constructive comments. NC is supported by a Centre National de la Recherche Scientifique of France and Laboratoire de Détection Géophysique (Bruyères-le-Châtel) contract No. 72B 087/00. The field experiment was part of the GéoFrance3D programme and was supported by the InterReg programme of the European Union. All those who helped are gratefully thanked.

## REFERENCES

Alsina, D. & Snieder, R., 1996. Constraints on the velocity structure beneath the Tornquist–Teisseyre Zone from beam-forming analysis, *Geophys. J. Int.*, **126**, 205–218.

Barker, J.S., Campillo, M., Sánchez-Sesma, F.J., Jongmans, D. & Singh, S.K., 1996. Analysis of wave propagation in the valley of Mexico from a dense broad-band array of seismometers, *Bull. seism. Soc. Am.*, **86**, 1667–1680.

Capon, J., 1970. Analysis of Rayleigh-wave multipath propagation at Lasa, *Bull. seism. Soc. Am.*, **60**, 1701–1731.

Coutant, O., Deschamps, A. & Gaffet, S., 2000. Permanent broad-band networks in the French Alps, *J. Seism.*, submitted.

Debelmas, J. & Kerckhove, C., 1980. Les Alpes franco-italiennes, *Géol. alpine*, **56**, 21–58.

ECORS-CROP Deep Seismic Sounding Group, 1989. A new picture of the Moho under the Western Alps, *Nature*, **337**, 249–251.

Ekström, G., Tromp, J. & Larson, E., 1997. Measurements and global models of surface wave propagation, *J. geophys. Res.*, **102**, 8137–8157.

Farra, V., 1990. Amplitude computation in heterogeneous media by ray perturbation theory: a finite element method approach, *Geophys. J. Int.*, **103**, 341–354.

Fréchet, J. *et al.*, 1998. A dense temporary seismic network in the W. Alps, *Ann. Geophys.*, EGS meeting, Nice, 88.

Jobert, N. & Jobert, G., 1983. An application of ray theory to the propagation of waves along a laterally heterogeneous spherical surface, *Geophys. Res. Lett.*, **10**, 1148–1157.

Keilis-Borok, V.I., 1989. *Seismic Surface Waves in a Laterally Inhomogeneous Earth*, Kluwer Academic, Dordrecht.

Knopoff, L., Mueller, S. & Pilant, W.L., 1966. Structure of the crust and upper mantle in the Alps from the phase velocity of Rayleigh waves, *Bull. seism. Soc. Am.*, **56**, 1009–1044.

Laske, G. & Masters, G., 1996. Constraints on global phase velocity maps from long-period polarization data, *J. geophys. Res.*, **101**, 16 059–16 075.

Lerner-Lam, A.L. & Park, J.J., 1989. Frequency-dependent refraction and multipathing of 10-100 second surface waves in the western Pacific, *Geophys. Res. Lett.*, **16**, 527–530.

Levshin, A.L., Ritzwoller, M.H. & Ratnikova, L.I., 1994. The nature and cause of polarization anomalies of surface waves crossing northern and central Eurasia, *Geophys. J. Int.*, **117**, 577–590.

McGarr, A., 1969. Amplitude variations of Rayleigh waves—horizontal refraction, *Bull. seism. Soc. Am.*, **59**, 1307–1334.

Ménard, G. & Thouvenot, F., 1984. Ecaillage de la lithosphère européenne sous les Alpes occidentales: arguments gravimétriques et sismiques liés à l'anomalie d'Ivréa, *Bull. Soc. géol. France*, **7**, 875–884.

Mooney, W.D., Laske, G. & Masters, G., 1998. CRUST 5.1: a global crustal model at 5° × 5°, *J. geophys. Res.*, **103**, 727–747.

Nataf, H.-C. & Ricard, Y., 1996. 3SMAC: an a priori tomographic model of the upper mantle based on geophysical modeling, *Phys. Earth. planet. Inter.*, **95**, 101–122.

Oliver, J., 1962. A summary of observed seismic surface wave dispersion, *Bull. seism. Soc. Am.*, **52**, 81–86.
Paul, A., Thouvenot, F., Fréchet, J., Cattaneo, M., Spallarossa, D. & Béthoux, N., 1998. Local earthquake tomography of the southwestern Alps (GéoFrance3D 1996 experiment), *Ann. Geophys.*, EGS meeting, Nice, 90.

Paulssen, H., Levshin, A.L., Lander, A.V. & Snieder, R., 1990. Time-and frequency-dependent polarization analysis: anomalous surface wave observations in Iberia, *Geophys. J. Int.*, **103**, 483-496.

Pavlis, G.L. & Mahdi, H., 1996. Surface wave propagation in central Asia: observations of scattering and multipathing with the Kyrgyzstan broadband array, *J. geophys. Res.*, **101**, 8437-8455.

Pedersen, H., Maupin, V. & Campillo, M., 1996. Wave diffraction in multilayered media with the indirect boundary element method: application to 3-D diffraction of long-period surface waves by 2-D lithospheric structures, *Geophys. J. Int.*, **125**, 545-558.

Perrier, G., 1979. La structure des Alpes occidentales déduite des données géophysiques, *Ecl. geol. Helv.*, **73**, 407-424.

Poupinet, G., Ellsworth, W.L. & Fréchet, J., 1984. Monitoring velocity variations in the crust using earthquake doublets: an application to the Calaveras Fault, *J. geophys. Res.*, **89**, 5719-5731.

Roberts, R.G. & Christoffersson, A., 1990. Decomposition of complex single-station three-component seismograms, *Geophys. J. Int.*, **103**, 55-74.

Sobel, P.A. & von Seggern, D.H., 1978. Applications of surface-wave ray tracing, *Bull. seism. Soc. Am.*, **68**, 1359-1380.

Trampert, J. & Woodhouse, J.H., 1995. Global phase velocity maps of Love and Rayleigh waves between 40 and 150 seconds, *Geophys. J. Int.*, **122**, 675-690.

Vidale, J.E., 1986. Complex polarization analysis of particle motion, *Bull. seism. Soc. Am.*, **76**, 1393-1405.

Woodhouse, J.H. & Wong, Y.K., 1986. Amplitude, *Geophys. J. R. astr. Soc.*, **87**, 753-773.