# The shear wave velocity structure of the heterogeneous alluvial plain of Beirut (Lebanon): combined analysis of geophysical and geotechnical data

Nancy Salloum, Denis Jongmans, Cécile Cornou, Dalia Youssef Abdel Massih, Fadi Hage Chehade, Armand Mariscal, Christophe Voisin

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

## ► To cite this version:

Nancy Salloum, Denis Jongmans, Cécile Cornou, Dalia Youssef Abdel Massih, Fadi Hage Chehade, et al.. The shear wave velocity structure of the heterogeneous alluvial plain of Beirut (Lebanon): combined analysis of geophysical and geotechnical data. Geophysical Journal International, 2014, 199 (2), pp.894-913. 10.1093/gji/ggu294. hal-04756450

HAL Id: hal-04756450

https://hal.science/hal-04756450v1

Submitted on 28 Oct 2024

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_2.png" /></figure>
# The shear wave velocity structure of the heterogeneous alluvial plain of Beirut (Lebanon): combined analysis of geophysical and geotechnical data

Nancy Salloum,<fnref n="1" />,<fnref n="2" /> Denis Jongmans,<fnref n="1" /> Cécile Cornou,<fnref n="1" /> Dalia Youssef Abdel Massih,<fnref n="3" /> Fadi Hage Chehade,<fnref n="2" /> Christophe Voisin<fnref n="1" /> and Armand Mariscal<fnref n="1" />

Accepted 2014 July 22. Received 2014 July 21; in original form 2014 March 26

## SUMMARY

A geotechnical and geophysical campaign was performed at sites located in the alluvial plain of the river of Beirut (Lebanon), which is characterized by a significant lateral and vertical geological variability, along with anthropogenic disturbances in the first metres. The method combination has allowed detecting a shallow conductive low velocity layer of varying depth and thickness, corresponding to a soft clay layer embedded in coarser formations. This layer was found to exert strong control on the experimental dispersion curves (estimated from both active and passive experiments) characterized by a continuous mode superposition at high frequency, associated with an increase in phase velocity. Vs profiles in boreholes turned out to be of prime importance for adequately defining the parametrization before inversion and for ensuring the reliability of the inversion dispersive estimates at low frequency. A major output of this study is also to show that this low velocity layer, along with the strong shear velocity contrast at its bottom, significantly contributes to the site seismic response, and could make it difficult to use the measured H/V peak frequency as a proxy for the soil thickness over bedrock.

**Key words:** Guided waves; Site effects; Wave propagation.

## 1 INTRODUCTION

The shallow geology of most sites is characterized by considerable spatial variability in geotechnical properties, as a result of the natural processes of erosion, weathering and deposition that continuously reshape the landscape. In the last centuries, human activities (building, mining and pumping) have also contributed to modify the near surface, in particular in urbanized areas which can be covered by anthropogenic deposits locally reaching 10–20 m (among others, Jongmans & Campillo 1990; Pagliaroli et al. 2014). This spatial heterogeneity existing at several scales in the superficial formations could be of prime importance for seismic and gravity hazard assessments. Indeed, numerous observations in earthquake engineering have shown that vertical heterogeneity in a basin-like structure (mainly seismic velocity and density) could significantly modify ground motions, a phenomenon which was dramatically evidenced during the large 1985 Michoacan and 1989 Loma Prieta earthquakes (Chin & Aki 1991; Chavez-Garcia & Bard 1994). In the damaged areas—the Mexico basin and the San Francisco basin—it was also shown that ground motions are variable over short distances, as a result of the lateral variability in the su-

picial formations (Bard et al. 1988; Hanks & Brady 1991). Characterization of the spatial variability in soil seismic properties is then of main importance for assessing seismic ground motion using deterministic or statistical techniques (Assimaki et al. 2003; Thompson et al. 2009). Previous studies have also showed that small-scale heterogeneity greatly influences the liquefaction potential of soil deposits (Popescu et al. 2005) and can affect the failure surface geometry in landslides (Flageollet et al. 2000; Bièvre et al. 2011).

In geotechnical engineering, shear wave velocity (Vs) has increasingly emerged as a key geophysical parameter for characterizing soil layer compactness. It is now widely used for evaluating site effects in engineering codes or for accurate ground motion modelling (e.g. Chavez-Garcia & Faccioli 2000; Komatitsch et al. 2004; Makra et al. 2005; Chaljub et al. 2010), for liquefaction assessment (Finn 2000) or evaluating landslide-related deconsolidation (Jongmans et al. 2009). In situ Vs measurements can be performed in various ways including borehole tests, shear-wave refraction, shear-wave reflection and surface wave techniques (see Hunter et al. 2002 and Boore 2006, for a review). Most borehole techniques offer the advantages of accuracy and of constant resolution with depth but the

<footnote n="1">¹ISTerre, Univ. Grenoble Alpes, CNRS, IRD, IFSTTAR, BP53X F-38041 Grenoble Cedex 9, France. E-mail: denis.Jongmans@ujf-grenoble.fr</footnote>

<footnote n="2">²Department of Civil Engineering, Lebanese University, Hadath Beirut, Lebanon</footnote>

<footnote n="3">³Department of Civil Engineering, Notre Dame University, NDU, Louaize, Lebanon</footnote>
disadvantage of being invasive and of providing isolated spatial in-
formation. The downhole test, for which a SH source is operated at
the surface and one or two three-component geophones are lowered
in the hole (Jongmans 1992), requires just one borehole, is then
relatively economic but has a resolution which slightly decreases
with source-geophone distance (wave attenuation and lower signal
to noise ratio). In contrast, surface seismic methods are non-invasive
and investigate a large volume of ground, with the limit of resolu-
tion decreasing with depth. In recent years, surface wave methods
have been increasingly applied to measure dispersion curves which
are then inverted for deriving Vs vertical profiles (for a review, see
Socco & Strobbia 2004; Foti et al. 2011). Although some recent
advances have been done to include the three components of the
ground motion (e.g. Poggi & Fäh 2010), the majority of ambient
vibration studies focus only on vertical components and on the
Rayleigh wave modes. Methods are classically divided into those
using active or passive sources (Boore 2006). A main limitation
of active sources lies in the difficulty to generate low frequency
waves, which limits the penetration to a few tens of metres (Park
et al. 1999; Bard et al. 2010). To overcome this problem, the low
frequency part of the dispersion curve can be obtained from the
analysis of array measurements of ambient vibrations (e.g. Park
et al. 2005; Bard et al. 2010; Zor et al. 2010), which are rich in low
frequencies. Different techniques (frequency-wavenumber, spatial
autocorrelation) can be used to extract dispersion characteristics of
ambient vibrations in order to derive vertical Vs profiles (for a re-
view, see Wathelet et al. 2008; Foti et al. 2011). Active and passive
techniques then appear complementary to each other for deriving
the surface wave dispersion curve (e.g. Wathelet et al. 2004; Park
et al. 2007; Zor et al. 2010; Di Giulio et al. 2012). Some additional
information, like the peak frequency in the H/V spectrum, the shape
of the H/V curve or the ellipticity of Rayleigh waves, can also be
introduced in the inversion process, in order to better constrain the
Vs profile (Scherbaum et al. 2003; Arai & Tokimatsu 2005; Picozzi
et al. 2005; Poggi et al. 2012; Hobiger et al. 2013).

The motivation behind this work is to benefit from an exceptional
set of geotechnical data gathered in downtown Beirut (Lebanon),
a town which was severely damaged and destroyed several times
by earthquakes (e.g. the 551 and 1202 AD earthquakes, Daëron
et al. 2007). The study site is located in the alluvial plain of the
river Beirut, characterized by soil deposits exhibiting strong lateral
and vertical facies variations. In particular, geotechnical tests have
showed the presence of a shallow soft clay layer of variable thick-
ness, the recognition of which is of crucial importance in geotech-
nical engineering. Also, the first few metres below the surface has
been significantly modified or affected by human activity.

The objective of this work is to test the capabilities of commonly
used geophysical methods [electrical resistivity tomography (ERT)
profiles, active seismic profiles, down-hole tests, single-station and
array ambient vibration measurements] to assess in a comprehen-
sive way the spatial variability of geotechnical characteristics in
these geological and urban conditions. Particular emphasis will be
placed on the soft clay layer geometry and on the required prior
knowledge and strategy to derive reliable Vs profiles using surface
wave inversion methods.

2 STUDY SITE: GEOLOGICAL CONTEXT AND INVESTIGATION

The study site is located in the centre of Beirut city (Lebanon), in
the flat alluvial plain of the river of Beirut (Fig. 1). The Quaternary
alluvium filling of variable thickness (a few tens of metres) consists
of interbedded layers of pebble, gravel, sand and clay, overlying
marly limestone of Tertiary age (Dubertret 1944). This carbon-
ate formation outcrops in the hills surrounding the alluvial plain
(Fig. 1).

The study area includes two building sites (labeled *Si*1 and *Si*2) and one site (*Si*3) where dozens of old buses have been abandoned (location in Fig. 2). Detailed geotechnical campaigns were

<figure><img src="image_3.png" /><figcaption>**Figure 1.** (a) Location of Beirut in Lebanon. (b) Geological map of Beirut's centre, simplified from Dubertret (1944) and location of the site study (red rectangle). Coordinates are in Universal Transverse Mercator (UTM).</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Location of the three sites *Si1* to *Si3* on the geological map, along with the borehole implantation. Coordinates are in Universal Transverse Mercator (UTM).</figcaption></figure>

conducted at the two building sites in 2008 (ACTS 2009). Fig. 2 shows the location of the 13 boreholes (50 m deep) and the 6 boreholes (20 m deep) drilled in sites *Si1* and *Si2*, respectively. In most of the boreholes, standard penetration tests (SPTs) were performed with a spacing of 1.5 m. Additionally, particle size distribution tests and Atterberg limits were made on core samples to classify the soils encountered in the boreholes (USCS classification, Cornforth 2005). At these two sites, geotechnical investigation has shown the presence of a shallow soft clay layer (with low penetration resistance) interbedded in coarser formations, the thickness of which varies from 2 to 7 m. The newly investigated bus site (Si3, Fig. 2) is characterized by a strong spatial superficial heterogeneity, highlighted by the presence of concrete pavement remains, Earth heaps, half-filled trenches and subsurface sewer concrete pipes. All

the previous data gathered on the three sites are synthesized and interpreted jointly with the results obtained in this study.

## 3 METHODS

Four 20 m deep boreholes (labelled B20–B23 in Fig. 3) were drilled at site *Si3* in 2011 November. SPTs were conducted in these boreholes during the drilling, with the measurement of the penetration resistance (number of blows: *N*) every 1.5 m. After casing, shear wave down-hole tests were performed in the four boreholes. Identification laboratory tests were performed on 21 core samples. A geophysical campaign was performed in 2011, including electrical resistivity tomography (ERT) profiles, active seismic profiles,
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Location of geotechnical and geophysical tests. B1–B23: Boreholes; A1–A8: passive seismic arrays; E1–E3: electrical profiles. H/V measurements are located by crosses.</figcaption></figure>

single-site ambient vibration measurements and the deployment of passive seismic arrays. Location of the experiments is indicated in Fig. 3.

Four down-hole tests were carried out in boreholes B20–B23 (Fig. 3). Shear waves were generated at the surface (from 1 to 2 m far from the borehole head), using a 5 kg sledge hammer hitting a wood plank sidewise, and were recorded in the borehole by two 3C 10 Hz geophones spaced 2.625 m apart and mechanically coupled to the PVC casing. The upper geophone position was 2 m depth in boreholes B21, B22 and B23, and 4 m in borehole B20, due to the presence of a cavity between 1.4 and 2.9 m. Measurements were performed with a spacing of 2 m in the boreholes. S-wave traveltimes were manually picked and shear wave velocities were determined by assuming straight ray paths between the source and geophones. Three ERT profiles (labelled E1–E3, Fig. 3) were performed using a Wenner–Schlumberger configuration with 64 electrodes 1.5–2 m apart. Apparent resistivity data were inverted using the software RES2DINV with the L1 norm (Loke & Barker 1996). Active surface seismic experiments were performed along 13 profiles

(S1–S13, Figs 3 and 4), deploying 24 vertical 4.5 Hz Mark Products geophones spaced at an interval of 1.5 or 2 m and connected to a 24-channel Geometries Geode acquisition unit. The signals generated by 5 kg vertical hammer drops were acquired with a sampling rate of 4000 Hz and duration of 2 s, allowing both P waves and Rayleigh waves to be recorded. The source was fired 10 times at both ends of the profile and individual signals were stacked to enhance signal-to-noise ratio. Dispersion maps were computed applying the frequency wave-number method (FK; Lacoss et al. 1969), using time-series corresponding to source-to-receiver distances larger than four times the source offset distance (typically 2 m) in order to reduce near-field effects (Socco & Strobbia 2004). Phase velocities related to the most energetic frequencies were then manually picked. The two dispersion curves obtained for each source offset were then averaged to get a mean dispersion curve and the uncertainty.

Single-station ambient noise measurements were performed at 151 sites (Fig. 3), using either Cityshark™ acquisition unit (Chatelain et al. 2000) connected to 5s-Lennartz sensors (cut-off period of 5 s) or Taurus (Nanometrics) connected to Güralp CMG40T
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Location of the active seismic profiles (S1–S13) performed at site Si3. Profiles are alternatively in light blue and brown, for clarity purpose.</figcaption></figure>

sensors (cut-off period of 30 s). 40 additional measurements located outside the geographical area displayed in Fig. 3 were also measured with the same instruments. These 30 min to 1 hr duration records were used to estimate the fundamental resonance frequency of the soil cover using the *H/V* technique, originally proposed by Nogoshi & Igarashi (1971) and wide-spread by Nakamura (1989). It consists in estimating the spectral ratio between the horizontal and vertical components of ambient noise recorded at the surface. The *H/V* ratio has been experimentally and numerically proved to provide the fundamental resonance frequency (or period) of the site (e.g. Bonnefoy-Claudet *et al.* 2006; Haghshenas *et al.* 2008; Albarello & Lunedei 2010). In this study, *H/V* spectral ratios were computed on 50-s stationary time windows, selected by applying a STA/LTA anti-triggering algorithm with LTA = 1 s and STA = 30 s. Computed Fourier amplitude spectra were smoothed using the

procedure of Konno & Ohmachi (1998), using a *b*-value of 40. The two horizontal components were then combined by computing the quadratic mean and the *H/V* ratios obtained for individual windows were averaged.

At eight single-station sites, ambient vibrations were also recorded by seismic arrays (A1–A8, Fig. 3) composed of eight three-component Güralp CMG40T connected to Nanometric Taurus digitizers. Each array had a circular shape (one station in the centre and seven on a circle) with apertures of 4, 10 and 20 m, depending on the site. Seismic array recordings were used to derive dispersion estimates (autocorrelation and dispersion curves) of Rayleigh waves by using *FK* (Lacoss *et al.* 1969), high-resolution *FK* (*HRFK*, Capon 1969) and Modified Spatial AutoCorrelation *MSPAC* techniques (Aki 1957; Bettig *et al.* 2001), as implemented in the software *Geopsy* (http://www.geopsy.org; Wathelet *et al.* 2008). For all
techniques, phase velocities were determined in 100 frequency bands logarithmically spaced between 0.5 and 50 Hz. In each frequency band, the time window length was defined as 100–300 times the central period, depending on the array. Average FK dispersion curves and related standard deviation were derived from the phase velocity histograms (probability density functions) for wavenumbers above the classical array resolution (kmin), defined as the mid-height width of the main lobe peak in the array response function. For HRFK processing, the lower wavenumber limit is half the classical array resolution (kmin/2; Asten & Henstridge 1984; Wathelet et al. 2008). MSPAC autocorrelation ratios were estimated for each semicircular subarray (ring), using the set of all possible combination of azimuths and distances between all sensor pairs within the array (Bettig et al. 2001). Dispersion estimates were inverted in order to derive Vs profiles, using the conditional Neighbourhood Algorithm (Wathelet 2008).

## 4 BOREHOLE TESTS

Fig. 5 presents the SPT profiles and borehole logs at the three sites Si1, Si2 and Si3, along with Vs values in the four boreholes of the site Si3. The layered structure at site Si1 is shown by the 50 m deep borehole B6 (Fig. 5a), which encountered, from top to bottom, (1) a 11-m-thick relatively stiff and coarse layer of gravel (with pebbles), (2) a 9-m-thick layer of soft ($N < 20$) clayey sand and (3) a thick layer of sand, the penetration strength of which increases with depth. The other SPT tests performed at site Si1 (Figs 5b–d) exhibit a similar shape, with some variations which will be studied further. No boreholes reached the carbonate bedrock in this area. Further north (site Si2; Figs 5e–h), the 20 m deep boreholes and N values corroborated the presence of a soft clay layer of variable thickness (a few metres to 7 m) interbedded between a compact coarse upper layer and a lower medium-dense sandy layer. At site Si3 (Figs 5i–l), this lower layer disappears and the clayey layer directly overlies the marly bedrock ($N > 100$) at a depth ranging between 14 m (B23) and 18 m (B20, B22), illustrating lateral facies variations. To the north of this site, the clay thickens to more than 16 m at B21 (Fig. 5j) where the bedrock was not reached. Core sample analysis at site Si3 showed that the clayey layer is made of low plasticity clay (CL) to sandy clay (SC), according to USCS classification.

We made some statistics on the N values measured in the clay layer. At sites Si2 and Si3, consistent mean values of $14 \pm 7$ were found, while N is higher and more variable ($21 \pm 9$) at site Si1. This suggests that the clay layer is slightly stiffer at site Si1, in agreement with the higher content in sand (SC) found in boreholes. In three boreholes (B2–B4), the clay layer is even barely distinguishable from the surrounding layers. Considering the upper limit of N values for defining the clay layer at each site, we studied its spatial variation across the study area. The maps of the layer thickness and top depth are shown in Figs 6(a) and (b), respectively. The clay layer is thin (2–4 m) at site Si2 and dramatically thickens northward to reach more than 16 m in borehole B21. At site Si1, the layer thickness varies between 0 and 10 m over short distances, highlighting the strong spatial variability in the alluvial plain. However, as mentioned before, the clay layer is locally stiffer and less well defined in this area, making its limits more uncertain. The top of the clay layer (Fig. 6b) is located at a relatively similar depth (10–14 m) across the whole study area, slightly diminishing in depth (6–8 m) at the northwestern part of site Si3, where the bedrock is suboutcropping.

Water measurements in the boreholes (from a few m deep to 1 m above ground level), suggest the existence of two aquifers: a superficial one located in the upper permeable layer of gravel overlying the clay layer, and a confined aquifer located in the sand layer and connected to the carbonate bedrock.

Vs profiles derived from four down-hole tests at site Si3 are shown in Figs 5(i)–(l) (black dots). Velocities are found to be very stable in the clay layer (average value around $150 \text{ m s}^{-1}$), while they are higher and variable in the upper gravel layer, consistently with the spatial heterogeneity also witnessed by SPT measurements. At the top the marly limestone, Vs is about $400 \text{ m s}^{-1}$, indicating highly weathered rocks.

## 5 ELECTRICAL TOMOGRAPHY RESULTS

The electrical images obtained after 7–10 iterations with a rms less than 7 per cent are shown in Fig. 7. For the chosen configuration, the penetration depth is between 24 m (E2, E3) and 17 m (E1). Electrical resistivity ($\rho$) values are spread out over a wide range from 1 to $1100 \text{ }\Omega\text{m}$, indicating a large variability in the electrical properties of the alluvial layers. The electrical tomography E1 (Fig. 7a), performed at site Si1, has evidenced a layering in resistivity (from 5 to more than $1000 \text{ }\Omega\text{m}$). Below a few m thick highly resistive ($\rho > 430 \text{ }\Omega\text{m}$) shallow layer, corresponding to the superficial gravelly soil comprising pebbles, the resistivity image indicates the presence of interbedded finer sediment lenses ($12 < \rho < 430 \text{ }\Omega\text{m}$). A homogeneous very conductive layer ($\rho < 12 \text{ }\Omega\text{m}$) was found at $10–11 \text{ m}$ depth, which matches the soft clayey layer shown by SPT and borehole data. The profile E2 exhibits higher resistivity values (generally over $90 \text{ }\Omega\text{m}$), with lower resistivity zones appearing in the central part and in the northern end, indicating the development of sandy layers and the significant thinning of the clay layer, which was also shown by the geotechnical data (Fig. 7b). Finally, the profile E3 performed at the site Si3 is characterized by a strongly heterogeneous superficial layer and by the presence of a very conductive body (resistivity lower than $10 \text{ }\Omega\text{m}$), the wideness of which increases with depth (Fig. 7c). This conductive plume pattern probably results from the discharge of pollutant substances in a hole dug at the surface, which was observed in the middle of the profile (Fig. 7c). This pollution strongly affects the electrical image and contributes to mask the bedding. This effect is discussed in Appendix A.

## 6 RESONANCE FREQUENCY MAPPING

The H/V spectral ratio curves measured in the study area are highly variable and three main different shapes can be identified, depending on site location. This difference is illustrated in Fig. 8(a), which shows curves extracted at three different sites (labelled H1–H3 in Fig. 8a). At site H1, the H/V curve is almost flat with a very small peak at relatively high frequency (4 Hz), suggesting the presence of outcropping bedrock in this area. At site H3 (Fig. 8a), the H/V curve exhibits a sharp peak with an amplitude of about 6, which may indicate a high impedance contrast between the sedimentary cover and the seismic bedrock (Albarello & Lunedei 2011). At site H2, the H/V curve shows a more complicated shape, in between a broad peak and plateau like-shape. Such a feature was previously observed by numerous authors and was interpreted as the result of a 2-D irregular basement interface (e.g. Uebayashi 2003; Bonnefoy-Claudet et al. 2006; Özalaybey et al. 2011; Le Roux et al. 2012).
<figure><img src="image_7.png" /><figcaption>**Figure 5.** SPT profiles and borehole logs at the sites *Si1* (top line: a, b, c and d), *Si2* (middle line: e, f, g and h) and *Si3* (bottom line: i, j, k and l). Vs values measured during down-hole tests at site *Si3* are plotted for boreholes *B20–B23* (i–l). The two-letter symbols of USCS soil classification are indicated when core samples are analysed (CL, low plasticity clay; SC, clayey sand; SM, silty sand; SW, well-graded sand; GW, well-graded gravel; GP, poorly graded gravel).</figcaption></figure>

In this study, this type of curves is observed at the edge of the alluvial plain. The map of the fundamental resonance frequency $f_{HV}$ (and corresponding period) values is displayed over the geological map in Fig. 8(b). Sites located in the alluvial plain exhibit sharp *H/V*

peak with period (frequency) ranging from 0.2 to 0.8 s (1.25–5 Hz), the largest resonance periods being observed in the valley center. In contrast, records on Tertiary and Cretaceous marly limestone exhibit no *H/V* peaks (site H1, Fig. 8a).
<figure><img src="image_8.png" /><figcaption>Figure 6. (a) Map of the clay layer thickness in the study area over the geological map; (b) map of the clay layer top.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 7. Electrical tomography sections. (a) Electrical profile E1, five iterations, Abs error = 4.5 per cent; (b) electrical profile E2, 7 iterations, Abs error = 4.7 per cent; (c) electrical profile E3, five iterations, Abs error = 6.7 per cent. The profile E1 was performed on a platform excavated 4 m below the ground level during civil engineering works.</figcaption></figure>

The theoretical transfer functions for vertically incident SH waves were computed using the method proposed by Kennett (1983) for simplified Vs profiles derived from the downhole tests in boreholes B20–B23 (Fig. 9), assuming quality factors of 50 and 25 for P and S waves, respectively. These transfer functions are compared to the observed H/V curves in Fig. 9. The modelled resonance frequencies (1.7–3 Hz) well agree with $f_{HV}$ measured at boreholes B21, B22 and B23. This outlines that the main seismic impedance contrast in these boreholes could be the interface between the clay layer and the weathered marly limestone (Figs 5j–l). In contrast, the SH resonance frequency simulated at borehole B20 (2.5 Hz) significantly overestimates the measured resonance frequency (1.4 Hz; Fig. 9). This suggests the presence of a main seismic impedance contrast at a depth larger than the clay/marly limestone interface.

# 7 SURFACE WAVE INVERSION

## 7.1 Dispersion estimates

Active and passive seismic experiments were locally combined to derive Vs profiles below the eight seismic arrays. Fig. 10 shows the dispersion estimates computed for the passive seismic array A7, also using the data from active seismic profiles S7 and S10 (see location in Fig. 4). Figs 10(a)–(c) display the seismograms recorded (at 1.5 m spacing) along the active S10 profile, the Fourier amplitude spectra and the FK dispersion map, respectively, while the results of noise array analysis are shown in Figs 10(d)–(e). Most of the seismic energy released by the shots (Fig. 10b) is between 20 and 80 Hz, limiting thus the frequency band for the extraction of the dispersion curve (Fig. 10c). Analysis of the two nearby S7 and S10
<figure><img src="image_10.png" /><figcaption>**Figure 8.** (a) Average H/V ratios (thick black line) ± standard deviation (dashed black line) at *H1*, *H2* and *H3*, the location of which is indicated in (b). (b) Map of fundamental resonance periods derived from microtremor recordings with the geological units. in the background.</figcaption></figure>

seismic profiles led to the average dispersion curve (green dots in Fig. 10c). Phase velocity histograms derived from *HRFK* analysis of seismic noise (Fig. 10d) clearly indicate a rise in velocity (blue dots) with increasing frequency between 10 and 22 Hz, and a decrease in velocity above 30 Hz. Similar phase velocity increase with frequency for the same spectral band is observed in the histograms derived from *MSPAC* analysis (Fig. 10e). In the low frequency band, the *MSPAC* histogram and the autocorrelation curve (Fig. 10f) are disturbed by a band kink around 4 Hz, which only allows getting dispersion estimates in the 2.8–4 Hz range (red dots; Fig. 10e). This kink is neither caused by a lack of energy in the Fourier amplitude spectrum for the vertical component (shown in Fig. 10f) nor by a poor azimuthal distribution of sensor pairs. Our interpretation is that this kink results from a jump to a higher mode of propagation, as already observed by Claprood *et al.* (2011) and Ikeda *et al.* (2012). Other smaller mode jumps are also visible on the phase velocity histograms (Fig. 10d) at 7, 13 and 18 Hz.

Fig. 11 displays the combined dispersion curves obtained at all sites. They were gathered in Figs 11(a)–(c) for the arrays located

in site *Si2* and in the eastern part of site *Si3*, in the western part of *Si3* and, and in site *Si1*, respectively. The combination of active and passive experiments allowed deriving the phase velocity over a wide frequency band (3–4.5 Hz to 60–100 Hz), except at site *Si1* where dispersion estimates were extracted up to 22 Hz only, because of the absence of active surface experiments.

In a stratified medium where *Vs* increases with depth, the phase velocity dispersion curve monotonically decreases with frequency, as observed for arrays *A1* and *A2* (Fig. 11c). All other dispersion curves exhibit an abnormal shape, with an increase in phase velocity for frequencies rising from 5–10 to 22–30 Hz (Figs 11a and b). Over 30 Hz, the curves show variable shapes, depending on the array location. This variability probably results from strong natural and anthropogenic heterogeneity in the very shallow layers. Comparing Figs 11(a) and (b), phase velocities for arrays *A6* and *A7* are lower by about 50 m s⁻¹ than those for arrays *A3*, *A4* and *A5*. This difference is probably linked to the thickness of the clay layer, which is larger in the northwestern part of site *Si3* than in its eastern part (Fig. 6a). The dispersion curve obtained at array *A8* shows a large increase
<figure><img src="image_11.png" /><figcaption>**Figure 9.** From left to right for each panel: measured shear wave velocities (black squares) overlaying simplified velocity structure (thin black line), observed average *H/V* curves (thick black lines) with related standard deviation (thin black lines) and theoretical 1-D SH transfer functions (grey lines). Upper left-hand panel: borehole B20 and seismic array A5; Upper right-hand panel: borehole B21 and seismic array A6; Lower left-hand panel: borehole B22 and seismic array A7; Lower right-hand panel: borehole B23 and seismic array A8.</figcaption></figure>

in phase velocity at low frequency (<15 Hz), consistently with the thinning of the alluvial layers in this area. Another striking feature is that the phase velocities computed below 8 Hz at arrays A1 and A2 (Fig. 11c) are far larger than those obtained at the other sites (Figs 10a and b). This issue will be discussed further.

The increase in phase velocity with frequency in some parts of experimental dispersion curves has already been reported in literature (among others Socco & Strobbia 2004; Maraschini & Foti 2010; Socco et al. 2010; Bergamo et al. 2011; Ikeda et al. 2012). The conventional explanation is the presence at depth of a low velocity layer, causing higher frequencies to propagate at higher velocities in a limited band of the dispersion curve (fundamental mode). Another effect of the soft layer, strengthening the first one, is that higher modes propagate with more energy than the fundamental one on given frequency bands and make the apparent dispersion curve jumping to higher modes at high frequencies. In order to understand the modal pattern of the experimental dispersion curves, the shear-wave velocity profiles at the borehole sites B20–B23 were used to compute the theoretical dispersion curves for the 10 first modes. They are compared to the array experimental curves in Fig. 12. At boreholes B20, B21 and B22 (arrays A5, A6 and A7, respectively), the agreement between the curves is very good over 4 Hz, highlighting the continuous energy transfer to successive higher modes for increasing frequencies. Although the general shape of the measured dispersion curve is also reproduced at borehole B23 (array A8), some discrepancies with theoretical phase velocities are

observed, which may be explained by complex wave propagation effects in this area where the bedrock is suboutcropping (Fig. 2). Fig. 12 also shows that the simplified ground models underestimate the phase velocities below 3 Hz and 4 Hz at arrays A5 and A6, respectively. This suggests that Vs tends to increase with depth in the marly limestone at these sites.

The continuous jump of modes to high frequencies and the resulting observed modal superposition has already been reported in Maraschini & Foti (2010) from surface wave real data and, numerically, in Zhang & Lu (2003a,b). As mentioned in Bergamo et al. (2011), this modal superposition makes the inversion of dispersion curve especially challenging, owing to the lack of any evidence of mode jumps. Classical inversion schemes would obviously fail in retrieving the correct shear velocity structure and more sophisticated inversion schemes, excluding a priori mode identification, would be required: full waveform inversion (Forbriger 2003), multimodal inversion accounting for the distribution of maximum displacement of each mode at the surface (Lu & Zhang 2006; Lu et al. 2007), or minimization of the stiffness matrix determinant (Maraschini et al. 2010). Although the last method is the most appealing for its simplicity and its computational effectiveness, the used linearized inversion scheme requires selecting an initial sufficiently realistic model. Testing and/or developing the above-mentioned inversion techniques are however beyond the scope of this paper. Alternatively, we chose to make the best use of all available geotechnical, geological data and resonance frequencies to define the
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Surface wave analysis at the vicinity of borehole B22, combining active (a–c) and passive (d–f) seismic experiments. (a) Seismograms (vertical component) recorded along the *S*10 profile for a shot located at –36.5 m offset. The profile location is indicated in Fig. 3. (b) Corresponding Fourier amplitude spectra. (c) Corresponding *FK* dispersion map. The continuous and dashed black lines are the spread length resolution limit after O’Neill (2004) and the aliasing limit, respectively. (d) Histogram of Rayleigh wave phase velocities obtained by stacking histograms observed on individual arrays by applying *HRFK* analysis to seismic noise records. The darkest colour indicates highest value of the probability density function. The dashed and continuous black lines are the theoretical array resolution limits, namely *kmin* and *kmin/2* (see text for details). (e) Rayleigh wave phase velocity histograms derived from *MSPAC* average autocorrelation ratios and average dispersion curve along with their related standard deviation (dots). The darkest colour indicates highest value of the probability density function. In (c), (d) and (e), green dots indicate the average dispersion curve (± standard deviation) computed for the *S*10 and *S*7 profiles (hammer records), while blue and red dots show the Rayleigh wave dispersion curves (± standard deviation) obtained on individual noise array records using *HRFK* and *MSPAC* analysis, respectively. (f) Average *MSPAC* curves (black dots) obtained from the analysis of the rings having inner and outer radius of 15.2 and 16.1 m, respectively. Continuous grey line indicates Fourier amplitude spectrum of the vertical component measured at the centre of the array.</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Dispersion curves (average phase velocity ± standard deviation) obtained by combining *FK*, *SPAC* and *MASW* results: (a) site *Si2* and *Si3*, (b) site *Si3* and (c) site *Si1*.</figcaption></figure>

parametrization and constrain the inversion of the part of the dispersion curve identified as the fundamental mode of Rayleigh waves.

## 7.2 Parametrization

The parametrization choice is of prime importance for the reliability of the dispersion curve inversion process and has to be driven by the *a priori* knowledge of the site (Renalier et al. 2010). In the case of the alluvial plain of Beirut, characterized by the presence of a low velocity layer, this prior information is crucial for the inversion reliability. The chosen ground model is made of three soil layers overlying a homogeneous half-space (sound bedrock), in agreement with the borehole data. The first two upper layers (sandy gravel and clay) are encountered at all sites (*Si1*, *Si2* and *Si3*), while the third one consists of a sand layer found at sites *Si1* and *Si2*, pinching out northward and disappearing in boreholes *B20–B23* (Site *Si3*; Fig. 5).

At this site, the third layer is made of weathered marly limestone of unknown depth. Table B1 shows the chosen parameter ranges in the different layers at the array sites *A2–A8*. The dispersion estimates at array *A1* and *A4* are not inverted, as they are very similar to those determined at arrays *A2* and *A5*, respectively. In the two upper layers, *Vs* and layers bottom depth ranges were mainly driven by the down-hole tests (site *Si3*) and *SPT* tests (all sites), considering the minimum and maximum values. In the shallow heterogeneous gravelly layer (*GL*), *Vs* values are allowed to vary from 300 to 1000 m s⁻¹ with a thickness range changing with the array location,

accordingly with the borehole data (see Table B1). In the clay layer (CL), *Vs* varies between 100 and 200 m s⁻¹ at site *Si3*. The same *Vs* range is considered for site *Si2*, *N* (SPT) values in the clay (14 ± 7) being similar at the two sites. In contrast, *N* values in the clay at site *Si1* are larger (21 ± 9) and the allowed *Vs* range has been expanded (100–300 m s⁻¹). In the lower weathered marly limestone layer (*WM*; arrays *A5–A8*), *Vs* values were allowed to vary between 400 and 2000 m s⁻¹, assuming a linear increase with depth. The same approach was used for the lower sand layer (SL) found below arrays *A2* and *A3*, in which *Vs* values as low as 200 m s⁻¹ are allowed at the layer top. The bottom depth range for the third layer is very large (15–200 m) for most arrays, as this information is not constrained by geotechnical data. The Poisson’s ratio varies from 0.2 to 0.5 for all layers, in order to limit the ratio *Vp/Vs* within an acceptable range, while density is fixed to 1800, 1500 and 2000 kg m⁻³ in the three layers overlying bedrock.

## 7.3 Shear wave velocity profile determination

Inversion results were obtained after more than 100 000 theoretical model computations, using from two to five seed numbers in order to be independent of the random sampling. An ensemble of shear wave velocity profiles that explain the data within their uncertainty bounds was extracted, applying the ‘acceptable solution’ concept (Lomax & Snieder 1994; Souriau et al. 2011): all calculated dispersion curves that fall within the data uncertainty bounds are given
<figure><img src="image_14.png" /><figcaption>**Figure 12.** Dispersion curve (black dots) measured at arrays A5–A8, and theoretical Rayleigh wave dispersion curves for the first ten modes, computed from the simplified velocity structure at boreholes B20–B23 (Fig. 9), respectively.</figcaption></figure>

the threshold misfit value of 1. From this statistically acceptable ensemble of Vs profiles, we selected those that provide SH peak frequency consistent with the measured H/V peak frequency (values lying within the interval of one standard deviation on either side of the mean, the values of which are given in Table B1), assuming that H/V peak frequency provides the resonance frequency, as previously discussed. Such a selection of Vs profiles was applied to arrays A2, A3 and A5 only, since very good agreement was already found between SH peak frequencies inferred from down-hole Vs profiles and H/V peak frequencies in arrays A6, A7 and A8 (Fig. 9). For these arrays, the maximum inverted depth was fixed to half the maximum measured wavelength. For all arrays, a reasonable number of acceptable solutions (>600, Table B1) explaining the data and respecting the parametrization constraints was found, with however significant differences from one array to the other. The inverted shear wave profiles are shown at the same scale in Fig. 13. We first consider the arrays A3–A8, for which only the low frequency part of the measured dispersion curve belonging to the fundamental mode of Rayleigh waves was inverted (see Table B1). In the first gravel layer, Vs values are between 300 and 400 m s⁻¹, while Vs in the clay is usually in the range 150–200 m s⁻¹. At arrays A2, A3 and A5, shear-wave velocities below the clay layer slightly increase to

a minimum depth of about 50–60 m, where the seismic substratum ($Vs > 500 \text{ m s}^{-1}$) is encountered.

As outlined before, phase velocities computed below 8 Hz at array A2 are significantly higher than those obtained at the other sites (Fig. 11). Inversion performed with the assumption that the experimental dispersion curve is the Rayleigh fundamental mode did not succeed in finding models explaining the steep part of the curve below 8 Hz. These results suggest that the first higher mode could be dominant at low frequency, as already shown by different authors (e.g. Socco & Strobia 2004). For array A2, it was then assumed that the measured phase velocities belong to the first higher mode from 3.2 to 8 Hz and to the fundamental mode of Rayleigh waves for higher frequencies. Inversion found acceptable models with velocities compatible with the other sites. We tried to understand the reasons for this difference in modal repartition at array A2. During array noise measurements at site Si1, a geotechnical company was injecting concrete for pile foundation in a 50 m deep borehole located between boreholes B9 and B13 (Fig. 3). The recorded noise time-series (Fig. 14a) show unusual long duration transients, probably generated by the concrete injection. This hypothesis is supported by the azimuth-frequency map (Fig. 14b), which shows a strong dominance of waves propagating to N40° between 4 and
<figure><img src="image_15.png" /><figcaption>**Figure 13.** Ensemble of inverted shear wave velocity profiles that explain the observed dispersion data within their uncertainty bounds for arrays A2–A8. The Vs profile displayed in red for A2 is used for ambient noise modelling (Fig. 14).</figcaption></figure>

10 Hz. These results suggest that the modal energy distribution of Rayleigh waves could have been influenced by local deep sources at this site. This hypothesis was tested by numerically simulating four minutes of ambient noise at site Si1, following the procedure described in Bonnefoy-Claudet et al. (2006). The Vs model used in this modelling is shown in Fig. 13 (A2, red line). Noise sources were distributed at discrete depths from 1.5 to 49.5 m depth in a borehole located 60 m far from the centre of array A2. Noise sources were modelled as a vertical force and Green functions were convolved with a delta-like source time function having a constant Fourier amplitude spectrum between 0.1 and 40 Hz. The Fourier spectra (vertical component), obtained after normalizing each synthetic time-series by its maximum amplitude at the array centre, are displayed for various source depths in Fig. 14(c). For frequencies below 10 Hz, noise amplitude spectra corresponding to sources located below the clay layer (21.5 m) are almost one order of magnitude lower than those recorded for shallower sources. Interestingly, the Fourier amplitude spectrum observed at the A2 array centre (after normalizing the time-series by its maximum value) is very similar to the synthetic Fourier amplitude spectra for sources located at about 25 m (Fig. 14c). Corresponding dispersion curves (Fig. 14d) clearly indicate that sources located in the upper sand and clay layers (i.e. source depths ranging from 1.5 to 21.5 m)

excite the fundamental mode of Rayleigh wave between 2 and 10 Hz. In the opposite, sources located below the clay layer lead to an apparent dispersion curve that is ranging between the first and the second higher modes of Rayleigh waves in the band 6–10 Hz. The proximity of the dispersion curve to the first or second higher mode is clearly controlled by the source depth. The eigenfunctions for the first four modes at 6 Hz (Fig. 14e) show that sources located at the top of the deepest sand layer (source depth of 25.5 m) excite more efficiently the first higher mode. Despite the model simplicity, these results support the interpretation that the low-frequency phase velocities measured at arrays A2 (and A1) correspond to the first higher Rayleigh wave mode (Fig. 14d).

## 8 2-D GROUND STRUCTURE INTERPRETATION

Two cross-sections (PP' and QQ', see location in Fig. 3) synthesising the geological, seismic and geotechnical results are presented in Fig. 15. The first cross-section, oriented NS and 560 m long, shows the presence of four geological layers (Fig. 15a). The upper layer, made of horizons of sand and gravel with pebbles, thins out northward (from 14 to 4 m). It is characterized by Vs
<figure><img src="image_16.png" /><figcaption>**Figure 14.** (a) Ambient noise signals recorded at the centre of array A2. (b) Histogram of azimuths of ambient noise waves as a function of frequency. The reddest colour indicates highest value of the probability density function. (c) Fourier amplitude spectra of the noise synthetics computed at various source depths (coloured lines) and of the measured noise (black line), at the center of array A2. Each seismic trace was normalized by its maximum before computing Fourier amplitude spectrum. (d) Rayleigh wave dispersion curves derived from noise synthetics for different source depths (coloured lines) and theoretical dispersion curves (black lines, first 10 modes) for the ground model shown in Fig. 13 (A2, red line). The black dots indicate the measured dispersion curve in A2 array; (e) Rayleigh wave mode shape (eigenfunctions) for the vertical component at 6 Hz.</figcaption></figure>

values usually ranging from 300 to 400 m s⁻¹. This coarse layer overlies the 2–12-m-thick soft clay layer, which exhibits a significant thinning in the middle of the profile (site Si2). Vs in the clay decreases to 150–200 m s⁻¹. Below sites Si1 and Si2, the clay layer is underlain with a thick sand layer (30–40 m), the velocity of which increases gradually with depth, from 300 to 500 m s⁻¹. A seismic bedrock (Vs > 500 m s⁻¹) is found at a minimum depth of 50 m. This sand layer pinches out northward and disappears below site Si3, where the clay layer overlies the marly limestone at about 18–20 m. Vs values in this layer is low (about 400 m s⁻¹), indicating that the rock is highly weathered. The seismic bed-rock (Vs > 500 m s⁻¹) was found below array A5 at a minimum depth of 60 m. The water levels measured in the boreholes and deduced from refraction tests suggest the presence of a confined aquifer below the clay layer and of a superficial aquifer in the upper coarse gravel layer. The second EW cross-section (Fig. 15b) illustrates the same three-layer structure (gravel, clay and marly limestone) encountered below site Si3, with a slight thinning of the alluvial layers at the eastern section end, where the marly bedrock was found at 14 m.

Again, the good comparison between the resonance frequencies computed for the three upper layers and the measured *H/V* peak

frequencies (Fig. 9) prohibits to estimate the depth of the seismic bedrock (Vs > 500 m s⁻¹) below the most northeastern arrays (A6–A8), unlike the results below array A5. As shown in Fig. 13, this problem results from the presence of the low velocity clay layer and of the significant resulting seismic velocity contrast at its base. The latter generates a resonance frequency very close to the measured f<sub>HV</sub> values below some arrays, which cannot bring any constrain on the bedrock depth in the inversion process. In order to better understand the relation between the different resonance frequencies, Fig. 16 compares the experimental *H/V* curves with the theoretical 1-D SH transfer functions computed for all the acceptable velocity profiles obtained at arrays A2 (site Si1), A3 (site Si2) and A5 (site Si3), considering and not considering the seismic bedrock. It turned out that, depending on the velocity models, the two transfer functions could be very similar both in terms of frequency and shape, like at array A3. In this case, the seismic bedrock depth is hard to define (Fig. 13b). In contrast, when the two transfer functions are different (Figs 16a and c), a significant contrast in velocity at depth is found during inversion (Figs 13a and c). These results evidence that the existence of a relatively thick low velocity layer could strongly control the seismic response of the site. In such velocity structure, f<sub>HV</sub> values could not be considered as a reliable proxy
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Interpretative ground structure along the PP' and QQ' cross-sections indicated in Fig. 3, with the Vs values obtained from downhole tests and surface wave inversion.</figcaption></figure>

of the sediment thickness, as usually observed in sites exhibiting increasing shear wave velocities with depth (e.g. Hizen *et al.* 2004; Nguyen *et al.* 2004; Le Roux *et al.* 2012).

## 9 CONCLUSIONS

A geotechnical and geophysical campaign was carried out at three sites in the alluvial plain of the river of Beirut (Lebanon). The geotechnical structure of this urban area is characterized by a strong horizontal and vertical variability in the alluvial layers. In particular, both electrical methods and geotechnical tests have detected the presence of a shallow soft and conductive clay layer embedded in coarser formations, the thickness of which varies from 2 to 12 m. In the computed dispersion curves, the presence of this low velocity layer was shown to generate at high frequency a continuous transfer of energy from one higher mode to another, resulting in a continuous increase in phase velocity with frequency. Such modal superposition prohibits use of classical inversion schemes based on mode numbering. In this case study, borehole data (geological log and Vs values in the upper layers) were of prime importance to adequately define the parametrization and to ensure the inversion reliability of the dispersion curve branch belonging to the fundamental mode at low frequency. It was also found that this low velocity layer generates a resonance phenomenon that largely contributes to the seismic response of the sites. In such seismic structures, the mea-

sured H/V peak frequency cannot be reliably used as an estimator of the thickness of soft layers overlying bedrock. Human activity was also found to have a significant influence on the application of geophysical prospecting at this urban site. In addition to the numerous subsoil modifications (pavement, excavations, concrete pipes, trenches and fillings), which have influenced all shallow geophysical measurements, the discharge of pollutant substances at the surface was found to have significantly affected some electrical images, masking the bedding and the presence of the clay layer. A last, anthropogenic effect has been the activation of continuous seismic sources related to building activity during the survey, which have excited higher modes of surface waves at low frequency. This case study illustrates the need of combining investigation methods in order to understand the geophysical measurements in a complex medium and to reach a robust 2-D/3-D model.

## ACKNOWLEDGEMENTS

This work was partially funded by project LIBRIS ANR Risknat-2009-006, the Institut de Recherche pour le Développement (IRD) and the CEDRE French-Lebanese collaboration program. Part of seismic instruments used in this study belong to the French national pool of portable Seismic instruments Sismob (INSU-CNRS). ISTerre is part of Labex OSUG@2020 (ANR10 LABX56). We thank the Associate Editor and two anonymous reviewers for their constructive comments.
<figure><img src="image_18.png" /><figcaption>**Figure 16.** Measured average *H/V* curves (thick black lines) with related standard deviation (thin black lines) and, theoretical 1D SH transfer functions computed from the ensemble of shear wave velocity profiles shown in Fig. 13. The dark and light grey curves stand for the SH transfer functions computed from the *Vs* profiles with and without considering the seismic bedrock.</figcaption></figure>

## REFERENCES

ACTS, 2009. 4748 Towers project-Korniche Al-Naher-Beirut, internal geotechnical engineering report, 23 p.

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

Albarello, D. & Lunedei, E., 2010. Alternative interpretations of horizontal to vertical spectral ratios of ambient vibrations: new insights from theoretical modeling, *Bull. Earthq. Eng.*, **8**(3), 519-534.

Albarello, D. & Lunedei, E., 2011. Structure of an ambient vibration wave-field in the frequency range of engineering interest ([0.5, 20] Hz): insights from numerical modelling, *Near Surf. Geophys.*, **9**(6), 543-559.

Arai, H. & Tokimatsu, K., 2005. S-wave velocity profiling by joint inversion of microtremor dispersion curve and horizontal-to-vertical (H/V) spectrum, *Bull. seism. Soc. Am.*, **95**, 1766-1778.

Assimaki, D., Pecker, A., Popescu, R. & Prevost, J., 2003. Effects of spatial variability of soil properties on surface ground motion, *J. Earthq. Eng.*, **7**(S1), 1-44.

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

Bard, P.Y., Campillo, M., Chavez-Garcia, F.J. & Sanchez-Sesma, F., 1988. The Mexico earthquake of September 19, 1985—a theoretical investigation of large-and small-scale amplification effects in the Mexico City Valley, *Earthq. Spectra*, **4**(3), 609-633.

Bard, P.Y. et al., 2010. From non-invasive site characterization to site amplification: recent advances in the use of ambient vibration measurements, in *Earthquake Engineering in Europe, Geotechnical, Geological, and Earthquake Engineering*, Vol. 17, pp. 105-123, Springer.

Bergamo, P., Comina, C., Foti, S. & Maraschini, M., 2011. Seismic characterization of shallow bedrock sites with multimodal Monte Carlo inversion of surface wave data, *Soil Dyn. Earthq. Eng.*, **31**(3), 530-534.

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 (SPAC): application to the Grenoble area, *Bollettino di Geofisica Teorica ed Applicata*, **42**(3-4), 281-304.

Bièvre, G., Kniess, U., Jongmans, D., Pathier, E., Schwartz, S., van Westen, C., Villemin, T. & Zumbo, V., 2011. Paleotopographic control of landslides in lacustrine deposits (Trièves plateau, French western Alps), *Geomorphology*, **125**, 214-224.

Bonnefoy-Claudet, S., Cornou, C., Bard, P.-Y., Cotton, F., Moczo, P., Kristek, J. & Fäh, D., 2006. H/V ratio: a tool for site effects evaluation. Results from 1-D noise simulations, *Geophys. J. Int.*, **167**(2), 827-837.

Boore, D., 2006. Determining subsurface shear-wave velocities: a review, in *Proceeding of the 3rd Int. Symp. on the Effects of Surface Geology on Seismic Motion*, Grenoble, 30 August-01 September, eds Bard, P.Y., Chaljub, E., Cornou, C., Cotton, F. & Guéguen, P., LCPC Editions, paper# 103.

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

Chaljub, E., Moczo, P., Tsuno, S., Bard, P. Y., Kristek, J., Käser, M. & Kristekova, M., 2010. Quantitative comparison of four numerical predictions of 3D ground motion in the Grenoble Valley, France, *Bull. seism. Soc. Am.*, **100**(4), 1427-1455.

Chatelain, J.L., Gueguen, P., Guillier, B., Frechet, J., Bondoux, F., Sarrault, J., Sulpice, P. & Neuville, J.M., 2000. CityShark: a user-friendly instrument dedicated to ambient noise (microtremor) recording for site and building response studies, *Seismol. Res. Lett.*, **71**(6), 698-703.

Chavez-Garcia, F.J. & Bard, P.-Y., 1994. Site effects in Mexico City eight years after the September 1985 Michoacan earthquakes, *Soil Dyn. Earthq. Eng.*, **13**, 229-247.

Chávez-García, F.J. & Faccioli, E., 2000. Complex site effects and building codes: making the leap, *J. Seismol.*, **4**(1), 23-40.
Chin, B.-H. & Aki, K., 1991. Simultaneous study of the source, path, and site effects on strong ground motion during the 1989 Loma Prieta Earthquake: a preliminary result on pervasive nonlinear effects, *Bull. seism. Soc. Am.*, **81**(5), 1859–1884.

Claprood, M., Asten, M.W. & Kristek, J., 2011. Using the SPAC microtremor method to identify 2D effects and evaluate 1D shear-wave velocity profile in valleys, *Bull. seism. Soc. Am.*, **101**(2), 826–847.

Cornforth, D.H., 2005. *Landslides in Practice: Investigation, Analysis and Remedial/Preventative Options in Soils*, Wiley, 624 pp.

Daëron, M., Klinger, Y., Tapponnier, P., Elias, A., Jacques, E. & Sursock, A., 2007. 12 000-year-long record of 10 to 13 paleoearthquakes on the Yammoune fault, Levant fault system, Lebanon, *Bull. seism. Soc. Am.*, **97**(3), 749–771.

Di Giulio, G. et al., 2012. Exploring the model space and ranking a best class of models in surface-wave dispersion inversion: application at European strong-motion sites, *Geophysics*, **77**(3), B147–B166.

Dubertret, L., *Map of Beirut and Surroundings*, Edition June 1944.

Finn, W., 2000. State-of-the-art of geotechnical earthquake engineering practice, *Soil Dyn. Earthq. Eng.*, **20**, 1–15.

Flageollet, J.C., Malet, J.P. & Maqaire, O., 2000. The 3D structure of the Super-Sauze earthflow: a first stage towards modelling its behaviour, *Phys. Chem. Earth, Part B: Hydrol., Oceans Atmos.*, **25**(9), 785–791.

Forbriger, T., 2003. Inversion of shallow seismic wavefields: I. Wavefield transformation, *Geophys. J. Int.*, **153**(3), 719–734.

Foti, S., Parolai, S., Albarello, D. & Picozzi, M., 2011. Application of surface-wave methods for seismic site characterization, *Surv. Geophys.*, **32**(6), 777–825.

Godio, A. & Naldi, M., 2003. Two-dimensional electrical imaging for detection of hydrocarbon contaminants, *Near Surf. Geophys.*, **1**(3), 131–137.

Haghshenas, E., Bard, P.Y. & Theodulidis, N., 2008. Empirical evaluation of microtremor H/V spectral ratio, *Bull. Earthq. Eng.*, **6**(1), 75–108.

Hanks, T.C. & Brady, A.G., 1991. The Loma Prieta earthquake, ground motion, and damage in Oakland, Treasure Island, and San Francisco, *Bull. seism. Soc. Am.*, **81**(5), 2019–2047.

Hizen, K.G., Weber, B. & Scherbaum, F., 2004. On the resolution of H/V measurements to determine sediment thickness, a case study across a normal fault in the Lower Rhine Embayment, Germany, *J. Earthq. Eng.*, **8**(6), 909–926.

Hobiger, M. et al., 2013. Ground structure imaging by inversions of Rayleigh wave ellipticity: sensitivity analysis and application to European strong-motion sites, *Geophys. J. Int.*, **192**(1), 207–229.

Hunter, J.A., Benjumea, B., Harris, J.B., Miller, R. D., Pullan, S.E., Burns, R.A. & Good, R.L., 2002. Surface and downhole shear wave seismic methods for thick soil site investigations, *Soil Dyn. Earthq. Eng.*, **22**(9), 931–941.

Ikeda, T., Matsuoka, T., Tsuji, T. & Hayashi, K., 2012. Multimode inversion with amplitude response of surface waves in the spatial autocorrelation method, *Geophys. J. Int.*, **190**(1), 541–552.

Jongmans, D., 1992. The application of seismic methods for dynamic characterization of soils in earthquake engineering, *Bull. Int. Assoc. Eng. Geol.*, **46**, 63–69.

Jongmans, D. & Campillo, M., 1990. The Liege earthquake of November 8, 1983, Damage distribution and site effects, *Earthq. Spectra*, **6**(4), 713–737.

Jongmans, D., Bièvre, G., Renalier, F., Schwartz, S., Beaurez, N. & Orengo, Y., 2009. Geophysical investigation of a large landslide in glaciolacustrine clays in the Trièves area (French Alps), *Eng. Geol.*, **109**(1), 45–56.

Kennett, B.L.N., 1983, *Seismic Wave Propagation in Stratified Media*, Cambridge Univ. Press.

Komatitsch, D., Liu, Q., Tromp, J., Süss, P., Stidham, C. & Shaw, J. H., 2004. Simulations of ground motion in the Los Angeles basin based upon the spectral-element method, *Bull. seism. Soc. Am.*, **94**(1), 187–206.

Konno, K. & Ohmachi, T., 1998. Ground-motion characteristics estimated from spectral ratio between horizontal and vertical components of microtremor, *Bull. seism. Soc. Am.*, **88**(1), 228–241.

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

Le Roux, O., Cornou, C., Jongmans, D. & Schwartz, S., 2012. 1-D and 2-D resonances in an Alpine valley identified from ambient noise measurements and 3-D modelling, *Geophys. J. Int.*, **191**(2), 579–590.

Lomax, A. & Snieder, R., 1994. Finding sets of acceptable solutions with a genetic algorithm with application to surface wave group dispersion in Europe, *Geophys. Res. Lett.*, **21**(24), 2617–2620.

Loke, M.H. & Barker, R.D., 1996. Least-squares deconvolution of apparent resistivity pseudosections, *Geophysics*, **60**(6), 1682–1690.

Lu, L. & Zhang, B., 2006. Inversion of Rayleigh waves using a genetic algorithm in the presence of a low-velocity layer, *Acoust. Phys.*, **52**(6), 701–712.

Lu, L., Wang, C. & Zhang, B., 2007. Inversion of multimode Rayleigh waves in the presence of a low-velocity layer: numerical and laboratory study, *Geophys. J. Int.*, **168**(3), 1235–1246.

Makra, K., Chávez-García, F.J., Raptakis, D. & Pitilakis, K., 2005. Parametric analysis of the seismic response of a 2D sedimentary valley: implications for code implementations of complex site effects, *Soil Dyn. Earthq. Eng.*, **25**(4), 303–315.

Maraschini, M. & Foti, S., 2010. A Monte Carlo multimodal inversion of surface waves, *Geophys. J. Int.*, **182**(3), 1557–1566.

Maraschini, M., Ernst, F., Foti, S. & Socco, L. V., 2010. A new misfit function for multimodal inversion of surface waves, *Geophysics*, **75**(4), G31–G43.

Nakamura, Y., 1989, A method for dynamic characteristics estimation of subsurface using microtremors on the ground surface, *Quat. Rep. Railway Technol. Res. Inst.*, **30**, 25–33.

Nguyen, F., Rompaey, G., Teerlynck, H., Van Camp, M., Jongmans, D. & Camelbeeck, T., 2004. Use of microtremor measurement for assessing site effects in Northern Belgium—interpretation of the observed intensity during the Ms = 5.0 June 11 1938 earthquake, *J. Seismol.*, **8**, 41–56.

Nogoshi, M. & Igarashi, T., 1971. On the amplitude characteristics of microtremor (Part 2), *J. Seismol. Soc. Jpn.*, **24**, 26–40.

O'Neill, A., 2004, Shear wave velocity model appraisal in shallow surface wave inversion, in *Proceedings of the ISC-2 on Geotechnical and Geophysical Site Characterization*, eds Viana de Fonseca, A. & Mayne, P.W., pp. 539–546, Millpress, Rotterdam, Portugal.

Özalaybey, S., Zor, E., Ergintav, S. & Tapırdamaz, M.C., 2011. Investigation of 3-D basin structures in the Izmit Bay area (Turkey) by single-station microtremor and gravimetric methods, *Geophys. J. Int.*, **186**(2), 883–894.

Pagliaroli, A., Moscatelli, M., Raspa, G. & Naso, G., 2014. Seismic microzonation of the central archaeological area of Rome: results and uncertainties, *Bull. Earthq. Eng.*, **12**, 1405–1428.

Park, C.B., Miller, R.D. & Xia, J., 1999. Multichannel analysis of surface waves, *Geophysics*, **64**(3), 800–808.

Park, C.B., Miller, R.D., Rydén, N., Xia, J. & Ivanov, J., 2005. Combined use of active and passive surface waves, *J. Environ. Eng. Geophys.*, **10**(3), 323–334.

Park, C.B., Miller, R.D., Xia, J. & Ivanov, J., 2007. Multichannel analysis of surface waves (MASW)—active and passive methods, *Leading Edge*, **26**(1), 60–64.

Picozzi, M., Parolai, S. & Richwalski, S.M., 2005. Joint inversion of H/V ratios and dispersion curves from seismic noise: estimating the S-wave velocity of bedrock, *Geophys. Res. Lett.*, **32**(11), L11308, doi: 10.1029/2005GL022878.

Poggi, V. & Fäh, D., 2010. Estimating Rayleigh wave particle motion from three-component array analysis of ambient vibrations, *Geophys. J. Int.*, **180**(1), 251–267.

Poggi, V., Fäh, D., Burjanek, J. & Giardini, D., 2012. The use of Rayleigh wave ellipticity for site-specific hazard assessment and microzonation: application to the city of Lucerne, Switzerland, *Geophys. J. Int.*, **188**(3), 1154–1172.

Popescu, R., Deodatis, G. & Nobahar, A., 2005. Effects of random heterogeneity of soil properties on bearing capacity, *Probab. Eng. Mech.*, **20**, 324–341.

Renalier, F., Jongmans, D., Savvaidis, A., Wathelet, M., Endrun, B. & Cornou, C., 2010. Influence of parameterization on inversion of surface wave dispersion curves and definition of an inversion strategy for sites with a strong Vs contrast, *Geophysics*, **75**(6), B197–B209.
Scherbaum, F., Hinzen, F.K. & Ohrnberger, M., 2003. Determination of shallow shear wave velocity profiles in the Cologne, Germany area using ambient vibrations, *Geophys. J. Int.*, **152**, 597–612.

Socco, L.V. & Strobbia, C., 2004. Surface-wave method for near-surface characterization: a tutorial, *Near Surf. Geophys.*, **2**(4), 165–185.

Socco, L.V., Foti, S. & Boiero, D., 2010. Surface-wave analysis for building near-surface velocity models—established approaches and new perspectives, *Geophysics*, **75**(5), A83–A102.

Souriau, A. et al., 2011. Multimethod characterization of the French-pyrenean valley of Bagnères-de-Bigorre for seismic-hazard evaluation: observations and models, *Bull. seism. Soc. Am.*, **101**(4), 1912–1937.

Thompson, E.M., Baise, L.G., Kayen, R.E. & Guzina, B.B., 2009. Impediments to predicting site response: Seismic property estimation and modeling simplifications, *Bull. seism. Soc. Am.*, **99**(5), 2927–2949.

Uebayashi, H., 2003. Extrapolation of irregular subsurface structures using the horizontal-to-vertical spectral ratio of long-period microtremors, *Bull. seism. Soc. Am.*, **93**(2), 570–582.

Wathelet, M., 2008. An improved neighborhood algorithm: parameter conditions and dynamic sjkscaling, *Geophys. Res. Lett.*, **35**(9), L09301, doi:10.1029/2008GL033256.

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

Wathelet, M., Jongmans, D., Ohrnberger, M. & Bonnefoy-Claudet, S., 2008. Array performances for ambient vibrations on a shallow structure and consequences over V s inversion, *J. Seismol.*, **12**(1), 1–19.

Zhang, B. & Lu, L., 2003a. Guided waves in a stratified half-space, *Acoust. Phys.*, **49**(4), 420–430.

Zhang, B. & Lu, L., 2003b. Rayleigh wave and detection of low-velocity layers in a stratified half-space, *Acoust. Phys.*, **49**(5), 516–528.

Zor, E., Özalaybey, S., Karaaslan, A., Tapirdamaz, M.C., Özalaybey, S.Ç., Tarancioğlu, A. & Erkan, B., 2010. Shear wave velocity structure of the İzmit Bay area (Turkey) estimated from active-passive array surface wave and single-station microtremor methods, *Geophys. J. Int.*, **182**(3), 1603–1618.

## APPENDIX A: CONTAMINATION EFFECT ON ELECTRICAL IMAGES

The blurring effect of contamination has already been observed by Godio & Naldi (2003) who studied long-term diesel oil pollution due to leakage from buried tanks in an alluvial plain. The biological degradation of hydrocarbons created highly conductive anomalies masking the fluvio-glacial sequence made up of silt, sand and gravel deposits. In order to test this effect in the alluvial plain of Beirut, the presence of a contaminant plume was numerically simulated using the RES2DMOD software (http://www.geotomosoft.com/). The four-layer electrical model (Fig. 8c) was built from the geotechnical data at the bus site and from the resistivity values deduced from electrical profile E1. From top to bottom, the model includes (1) a

<figure><img src="image_19.png" /><figcaption>**Figure A1.** (a) Measured pseudosection (E3). (b) Electrical tomography (E3). (c) Synthetic resistivity model of the pollution (see text for details). (d) Pseudosection computed for the synthetic resistivity model. (e) Electrical resistivity tomography obtained for the synthetic model (Abs error = 3.3 per cent after 4 iterations).</figcaption></figure>
superficial conductive layer (12 Ωm) with a thickness of 1 m, (2) a 4-m-thick resistive layer made up of sand and gravel, (3) a 13-m-thick clay layer (12 Ωm) and (4) the marly bedrock (50 Ωm). A conductive contaminant plume (1 Ωm) with geometry close to the one observed in Fig. 8(b) was considered. Some shallow conductive and resistive anomalies were introduced to simulate the heterogeneity of the near-surface. A 1 per cent noise was introduced in the synthetic data and the inversion was performed using the L1 norm. The measured and simulated pseudo-sections are shown in Figs 8(a) and (d), respectively. Although they exhibit differences, they are characterized by the same range in apparent resistivity values and similar

patterns. The inverted synthetic electrical image (Fig. 8e), the features of which are comparable to those of profile E3, shows that the presence of the conductive contaminant plume causes the electrical bedding to be blurred out. In particular, the contact between the resistive coarse soil layer and the underlying clay layer is untraceable. These results show that this anthropogenic effect prevents electrical tomography to be used for determining the alluvial sequence in the polluted area.

APPENDIX B: PARAMETRIZATION

<table><caption>**Table B1.** Parametrization chosen for the surface wave inversion by using the knowledge inferred from geotechnical and geophysical tests at the different array sites. The constraint imposed on the resonance frequency during inversion is indicated in the second column while the frequency range for the dispersion curve inversion is given in the third column. Bottom depth ranges are indicated in the fourth column for the three upper layers: GL, gravel layer; CL, clay layer; SL, sand layer; or WM, weathered marly limestone layer.</caption><thead><tr><th>Array #</th><th>H/V peak frequency (Hz)</th><th>Frequency range for DC inversion (Hz)</th><th>Range of bottom depth (m)</th><th>Range of Vs (m s<sup>-1</sup>)</th><th>Velocity law</th><th>Number of acceptable Vs profiles</th></tr></thead><tbody><tr><td>A2</td><td>1.40 ± 0.10</td><td>3.2–20.0</td><td>GL: 8–16<br>CL: 8–30<br>SL: 50–200<br>Half-space</td><td>300–1000<br>100–300<br>Top: 200–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>8070</td></tr><tr><td>A3</td><td>1.56 ± 0.10</td><td>2.8–14.7</td><td>GL: 8–12<br>CL: 13–20<br>SL: 15–200<br>Half-space</td><td>300–1000<br>100–200<br>Top: 200–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>606</td></tr><tr><td>A5</td><td>1.42 ± 0.13</td><td>3.3–14.0</td><td>GL: 8–10<br>CL: 18–20<br>WM: 18–200<br>Half-space</td><td>300–1000<br>100–200<br>Top: 400–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>861</td></tr><tr><td>A6</td><td>1.70 ± 0.10</td><td>3.0–9.8</td><td>GL: 5–7<br>CL: 18–22<br>WM: 18–200<br>Half-space</td><td>300–1000<br>100–200<br>Top: 400–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>3420</td></tr><tr><td>A7</td><td>1.80 ± 0.70</td><td>2.9–12.0</td><td>GL: 6–8<br>CL: 18–20<br>WM: 18–200<br>Half-space</td><td>300–1000<br>100–200<br>Top: 400–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>22 246</td></tr><tr><td>A8</td><td>3.30 ± 0.90</td><td>4.8–19.0</td><td>GL: 5–8<br>CL: 12–15<br>WM: 15–200<br>Half-space</td><td>300–1000<br>100–200<br>Top: 400–2000<br>Bottom: 400–2000<br>400–3500</td><td>Uniform<br>Uniform<br>Linear increase<br>with five sublayers</td><td>14 095</td></tr></tbody></table>