\makeatletter \@ifundefined{HCode} {\documentclass[screen, CRGEOS, Unicode, Thematic]{cedram} %\usepackage[sort&compress,square,comma,authoryear]{natbib} \usepackage[square,comma,authoryear]{natbib} \newenvironment{Table}{\begin{table}}{\end{table}} \renewcommand\bibsection{\section*{\refname}} \newenvironment{noXML}{}{} \def\thead{\noalign{\relax}\hline} \def\endthead{\noalign{\relax}\hline} \def\tbody{\noalign{\relax}} \def\bcaption#1#2#3{\caption{Caption continued on next page.}\end{figure*}\setcounter{figure}{#2}\begin{figure*}\vspace*{-1pc}\caption{(\textbf{cont}.)\space#1\space #3}} \def\botline{\\\hline} \def\0{\phantom{0}} \def\tsup#1{$^{{#1}}$} \def\tsub#1{$_{{#1}}$} \RequirePackage{etoolbox} \usepackage{upgreek} \def\jobid{crgeos20230720} %\graphicspath{{/tmp/\jobid_figs/web/}} \graphicspath{{./figures/}} \newcounter{runlevel} \usepackage{multirow} \def\nrow#1{\@tempcnta #1\relax% \advance\@tempcnta by 1\relax% \xdef\lenrow{\the\@tempcnta}} \def\morerows#1#2{\nrow{#1}\multirow{\lenrow}{*}{#2}} \let\MakeYrStrItalic\relax \def\refinput#1{} \def\back#1{} \csdef{Seqnsplit}{\\} \DOI{10.5802/crgeos.243} \datereceived{2023-08-10} \daterevised{2024-03-29} \dateaccepted{2024-05-27} \ItHasTeXPublished } { \PassOptionsToPackage{authoryear}{natbib} \documentclass[crgeos]{article} \usepackage[T1]{fontenc} \def\CDRdoi{10.5802/crgeos.243} \makeatletter \def\CDRsupplementaryTwotypes#1#2{} \def\selectlanguage#1{} \def\bcaption#1#2#3{\caption{#1\space#3}} \def\href#1#2{\url[#1]{#2}} \let\newline\break } \makeatother \dateposted{2024-07-22} \begin{document} \begin{noXML} \shortrunauthors \CDRsetmeta{articletype}{research-article} \makeatletter \def\TITREspecial{\relax} \def\cdr@specialtitle@english{New Developments in Passive Seismic Imaging and Monitoring} \def\cdr@specialtitle@french{Nouveaux d\'eveloppements en mati\`ere d'imagerie sismique passive et de surveillance} \makeatother \title{Potential and limitations of noise-based surface-wave tomography for numerical site effect estimation: a case study in the French Rh\^{o}ne valley} \alttitle{Potentiel et limites de l'imagerie sismique passive en ondes de surface pour l'estimation num\'erique des effets de site : un cas d'\'etude dans la vall\'ee du Rh\^one (France)} \author{\firstname{Fran\c{c}ois} \lastname{Lavou\'{e}}\CDRorcid{0000-0002-1723-4730}\IsCorresp} \address{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \address{Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, UGE, ISTerre, Grenoble, France} \email[F. Lavou\'{e}]{francois.lavoue@univ-grenoble-alpes.fr} \author{\firstname{B\'{e}r\'{e}nice} \lastname{Froment}\CDRorcid{0000-0001-5570-9365}} \addressSameAs{1}{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \author{\firstname{C\'{e}line} \lastname{G\'{e}lis}\CDRorcid{0000-0002-0915-3535}} \addressSameAs{1}{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \author{\firstname{Pierre} \lastname{Bou\'{e}}\CDRorcid{0000-0001-9153-4048}} \addressSameAs{2}{Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, UGE, ISTerre, Grenoble, France} \author{\firstname{Emmanuel} \lastname{~Chaljub}\CDRorcid{0000-0001-5498-885X}} \addressSameAs{2}{Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, UGE, ISTerre, Grenoble, France} \author{\firstname{Laurent} \lastname{Stehly}\CDRorcid{0000-0002-1854-7157}} \addressSameAs{2}{Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, UGE, ISTerre, Grenoble, France} \author{\firstname{Sophie} \lastname{Beaupr\^{e}tre}} \address{Sisprobe/EGIS, 3 rue du Dr Schweitzer, 38180 Seyssins, France} \author{\firstname{Florent} \lastname{De Martin}\CDRorcid{0000-0003-3746-9067}} \address{BRGM (French Geological Survey), Risks and Prevention Division, Orl\'{e}ans, France} \author{\firstname{Lo\"{i}c} \lastname{~Gisselbrecht}\CDRorcid{0000-0002-1943-1076}} \addressSameAs{1}{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \addressSameAs{2}{Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, UGE, ISTerre, Grenoble, France} \author{\firstname{Marco} \lastname{Pilz}\CDRorcid{0000-0002-8575-579X}} \address{GFZ, German Research Center for Geosciences, Potsdam, Germany} \author{\firstname{Denis} \lastname{Moiriat}\CDRorcid{0000-0003-1264-122X}} \addressSameAs{1}{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \author{\firstname{Edward Marc} \lastname{Cushing}\CDRorcid{0000-0002-1601-8323}} \addressSameAs{1}{Institut de Radioprotection et S\^{u}ret\'{e} Nucl\'{e}aire (IRSN), PSE-ENV, SCAN, BERSSIN/BEHRIG, Fontenay-aux-Roses, France} \keywords{\kwd{Seismology} \kwd{Seismic hazard} \kwd{Site effects} \kwd{Passive seismic imaging} \kwd{Numerical simulations}} \altkeywords{\kwd{Sismologie} \kwd{Al\'ea sismique} \kwd{Effets de site} \kwd{Imagerie sismique passive} \kwd{Simulations num\'eriques}} \thanks{French National Research Agency (Agence Nationale pour la Recherche, ANR, grant number ANR19-CE31-0029), German Research Foundation (Deutsche Forschungsgemeinschaft, DFG, project number 431362334), HPC resources of TGCC under the allocation 2022-gen13461 made by GENCI, GRICAD infrastructure} \begin{abstract} In certain geological settings such as sedimentary basins, the ground motion induced by an earthquake may be amplified by local site conditions. Estimating these site effects is important for seismic hazard assessment but can be difficult to do empirically due to the scarcity of site-specific field data in time and space, especially in low-to-moderate seismicity regions where the earthquakes needed for measuring the site effects have long return periods. In this study, we try to overcome these limitations and investigate an alternative approach based on ambient seismic noise and numerical simulations. More specifically, we use a 3D numerical model of seismic properties derived from Ambient Noise Surface-Wave Tomography (ANSWT) for 3D numerical simulations of seismic wave propagation, and consequently for a numerical estimation of seismic amplification in the basin. We illustrate the approach on a target site located in the French Rh\^{o}ne valley, where the Messinian salinity crisis has dug a paleo-canyon which is now filled by soft sediments in direct contact with a harder substratum, thereby providing typical conditions for significant site effects, as also observed by previous studies in the area. This work makes use of two dedicated datasets. On one hand, we use earthquake recordings acquired by a network of broadband stations deployed over the target site over 8 months, in order to estimate seismic amplification in the basin with respect to a rock-site reference via Standard Spectral Ratios (SSR), which we consider as our reference for evaluating our numerical results. On the other hand, we exploit one-month-long ambient noise recordings acquired by a dense array of 400 3C sensors. Prior to this work, this noise data was used to build a 3D shear-wave velocity (\textit{V\textsubscript{S}}) model of the target site via ANSWT, and also to estimate seismic amplification via noise-based Standard Spectral Ratios (SSRn). The obtained ANSWT model well reproduces the main geological structures of the basin, with lateral variations of velocities at depth depicting the deeper parts of the basin. However, our simulation results also show that some of its limitations related to surface wave sensitivity and resolution capability have an impact on the numerical amplification predicted in the basin. In particular, this ANSWT model lacks clear basin edges in order to efficiently trap seismic waves in the basin and to generate significant 3D wave propagation effects (diffractions, reflections, and generation of laterally propagating surface waves at the edges of the basin). As a result, the numerical amplification predicted in the ANSWT model remains dominated by a 1D response and does not reproduce the broadband character of the observed amplification at locations affected by significant 3D propagation effects. On the other hand, the numerical amplification predicted in the ANSWT model shows a good agreement with the observations at locations that seem less affected by 3D propagation effects, including in complex regions of the model where lateral variations must be taken into account. Our results therefore contribute to identify and better understand the potential and limitations of using ANSWT models for numerical site effect estimation. This study allows us to propose perspectives for future work to improve the approach, which remains promising for site effect assessment in low- to moderate-seismicity contexts. \end{abstract} \begin{altabstract} Dans les bassins s\'{e}dimentaires, le mouvement du sol suite \`{a} un s\'{e}isme peut \^{e}tre amplifi\'{e} selon les conditions locales de site. L'estimation de ces effets de site est importante pour l'\'{e}valuation correcte et pr\'{e}cise de l'al\'{e}a sismique, mais est souvent rendue difficile \`{a} faire exp\'{e}rimentalement du fait du peu de donn\'{e}es disponibles sur un site sp\'{e}cifique, surtout dans des r\'{e}gions de sismicit\'{e} faible \`{a} mod\'{e}r\'{e}e o\`{u} les s\'{e}ismes n\'{e}cessaires \`{a} la mesure de ces effets de site ont des longues p\'{e}riodes de retour. Dans cette \'{e}tude, nous explorons une approche alternative aux m\'{e}thodes classiques d'estimation des effets de site, en proposant de baser cette estimation sur des donn\'{e}es de bruit sismique ambiant et sur des simulations num\'{e}riques. Plus pr\'{e}cis\'{e}ment, nous proposons d'utiliser des m\'{e}thodes d'imagerie passive, en particulier la tomographie d'ondes de surface extraites du bruit sismique (\textit{Ambient Noise Surface-Wave Tomography}, ANSWT) pour construire des mod\`{e}les num\'{e}riques 3D des propri\'{e}t\'{e}s sismiques d'une zone d'\'{e}tude, mod\`{e}les que nous utilisons ensuite au sein de simulations num\'{e}riques de propagation d'ondes pour calculer une estimation num\'{e}rique d'amplification sismique dans un bassin. Nous illustrons cette approche avec un cas d'\'{e}tude situ\'{e} dans le bassin de Tricastin, en vall\'{e}e du Rh\^{o}ne (France), o\`{u} la crise de salinit\'{e} messinienne a creus\'{e} un pal\'{e}o-canyon maintenant rempli de s\'{e}diments meubles en contact avec un substratum indur\'{e}, cr\'{e}ant ainsi des conditions propices \`{a} l'\'{e}mergence d'effets de site significatifs document\'{e}s par des \'{e}tudes pr\'{e}c\'{e}dentes. Ce travail utilise deux jeux de donn\'{e}es compl\'{e}mentaires. Nous utilisons d'une part des donn\'{e}es de sismicit\'{e} enregistr\'{e}es par un r\'{e}seau de stations large bande d\'{e}ploy\'{e} pendant 8 mois sur la zone d'\'{e}tude. Ces donn\'{e}es de sismicit\'{e} nous permettent de mesurer l'amplification sismique dans le bassin par rapport \`{a} une r\'{e}f\'{e}rence au rocher via le calcul de rapports spectraux standards (\textit{Standard Spectral Ratios}, SSR) qui nous servent de points d'observation de r\'{e}f\'{e}rence pour \'{e}valuer nos r\'{e}sultats num\'{e}riques. D'autre part, nous exploitons des enregistrements continus de bruit sismique acquis par un r\'{e}seau dense de 400 capteurs autonomes 3-composantes. En amont de ce travail, ces donn\'{e}es de bruit ont \'{e}t\'{e} utilis\'{e}es pour construire un mod\`{e}le 3D de vitesse des ondes de cisaillement (\textit{V\textsubscript{S}}) par ANSWT, ainsi que pour estimer l'amplification sismique dans le bassin via des SSR bas\'{e}s sur le bruit sismique (\textit{noise-based} SSR, SSRn). Le mod\`{e}le ANSWT obtenu reproduit bien les structures g\'{e}ologiques attendues \`{a} grande \'{e}chelle, avec des variations lat\'{e}rales de vitesses en profondeur correspondant aux parties profondes du bassin. Cependant, nos r\'{e}sultats de simulation montrent aussi que certaines limites de ce mod\`{e}le li\'{e}es \`{a} la sensibilit\'{e} et \`{a} la capacit\'{e} de r\'{e}solution des ondes de surface utilis\'{e}es pour l'imagerie ont un effet sur l'amplification num\'{e}rique pr\'{e}dite dans le bassin. En particulier, ce mod\`{e}le ANSWT ne comporte pas de limites franches entre le bassin et les sites au rocher affleurant en surface. Par cons\'{e}quent, dans les simulations num\'{e}riques effectu\'{e}es dans ce mod\`{e}le, les ondes sismiques ne sont pas pi\'{e}g\'{e}es de mani\`{e}re efficace dans le bassin, dont les bords --- ou plus pr\'{e}cis\'{e}ment, le manque de bords --- ne g\'{e}n\`{e}rent pas de diffractions significatives d'ondes de surface susceptibles de se propager lat\'{e}ralement dans le bassin. L'amplification sismique pr\'{e}dite dans le bassin reste donc domin\'{e}e par une r\'{e}ponse 1D, et ne reproduit pas le caract\`{e}re large bande de l'amplification observ\'{e}e en certains points du bassin, visiblement affect\'{e}s par des effets de propagation 3D. En revanche, l'amplification num\'{e}rique pr\'{e}dite dans le mod\`{e}le montre un bon accord avec les observations en d'autres points du bassin, moins affect\'{e}s par ces effets de propagation 3D, notamment dans des zones relativement complexes o\`{u} les variations lat\'{e}rales de la structure du sous-sol doivent \^{e}tre prises en compte. Dans ces conditions, l'imagerie passive semble fournir une caract\'{e}risation satisfaisante du sous-sol \`{a} m\^{e}me d'expliquer l'amplification observ\'{e}e. Nos r\'{e}sultats permettent donc d'identifier et de mieux comprendre le potentiel et les limites de l'utilisation de mod\`{e}les sismiques 3D construits par imagerie sismique passive pour l'estimation num\'{e}rique des effets de site. Ces r\'{e}sultats nous am\`{e}nent \`{a} proposer plusieurs perspectives d'am\'{e}lioration de cette approche, qui reste prometteuse pour l'estimation des effets de site en contexte de sismicit\'{e} faible \`{a} mod\'{e}r\'{e}e. \end{altabstract} \maketitle \vspace*{1pc} \twocolumngrid \end{noXML} \section{Introduction}\label{sec1} In sedimentary basins, the impedance contrast at the interface between soft sedimentary layers and the underlying bedrock leads to the trapping of seismic waves within the sedimentary in-filling. This gives rise to complex wave phenomena (body wave resonance, generation and diffraction of surface waves at the edges of the basin, vertical and lateral reverberations, focusing effects) which highly depend on the three-dimensional (3D) geometry of the basin and often result in an increased amplitude (or ``amplification'') and duration of the ground motion (or ``duration lengthening''). This modification of ground motion due to local geology is referred to as \textit{site effects}. Site effects have been the subject of many studies, especially since the devastating 1985 Mexico earthquake that brought to light the influence of local soil conditions on the strong amplification of ground motion observed in Mexico City, despite the long distance to the seismic source \citep[e.g.,][]{Campilloetal1989}. This example, along with many other observations around the world \citep[e.g.,][]{Kawase1996,Gravesetal1998,Lebrunetal2001, Rotenetal2008,Bindietal2011,Ktenidouetal2016}, makes it clear that the quantification of these site effects is essential for seismic hazard assessment (SHA). Because they are related to local soil conditions, site effects can be highly variable from one site to another, and therefore require site-specific studies for a robust estimation that accounts for the whole\unskip\break complexity of wave phenomena, in particular in 3D geological structures. The French-German DARE project (\textit{Dense ARray for site effect Estimation}) has been conceived and designed in this line. The idea is to implement various and complementary approaches to perform a detailed study of site effects at a target site located in the French Rh\^{o}ne valley. The area hosts critical facilities including nuclear installations, thereby motivating the need for robust SHA studies locally. This site is located on the deep and elongated Messinian Rh\^{o}ne Canyon, whose geometry and lithological characteristics make it a good candidate for generating multidimensional site effects. DARE is centred on the exploitation of dense and complementary datasets acquired in the area \citep{Fromentetal2022a}. The project proposes to investigate the interest of using such datasets for a robust estimation of site effects, especially in low-to-moderate seismicity areas such as metropolitan France. One first, standard approach considered in DARE relies on earthquake recordings through the calculation of so-called site/reference Standard Spectral Ratios \citep[SSR,][]{Borcherdt1970}. SSR estimate the local seismic amplification by direct comparison between earthquake seismograms simultaneously recorded at a given \textit{site} station laying on a sedimentary basin (subject to \textit{site effects)} with respect to a nearby \textit{reference} station (typically on a bedrock outcrop, considered free of site effects). This empirical method has proven to be efficient for a robust quantification of site effects in various configurations. The implementation of this method may however show some difficulties in low-to-moderate seismicity areas, such as mainland France \citep{Traversaetal2020}, where moderate to large earthquakes $(M_w>5.0)$ have long return periods, therefore requiring long deployments. Alternative approaches may thus be considered to complement such seismicity-based analysis of site effects. To this end, the DARE project investigates two main research tracks in order to explore the potential of using weak but ubiquitous vibrations---known as ambient seismic noise---as an alternative source of data for estimating site effects due to complex wave propagation in sedimentary basins \citep[e.g.,][]{Boueetal2016}. On one hand, seismic noise can be used directly for empirical estimations of site effects, via H/V analysis \citep[e.g.,][]{Bonnefoy-Claudetetal2006,Spicaetal2017,vanGinkeletal2019} or noise-based SSR estimations \citep[e.g.,][]{Perronetal2018,Gisselbrechtetal2023}. On the other hand, ambient noise could also be used in a more indirect way as the initial ingredient to build seismic models of the target site that could then be used for numerical prediction of the ground motion and thus of seismic amplification. The present paper focuses on this numerical\unskip\break aspect. By giving access to the complete wavefield, simulations help to evaluate the spatial variability of site effects and more importantly to understand the underlying physical parameters to which they are sensitive \citep[e.g.,][]{DeMartinetal2021}. 1D modelling of body-wave resonance phenomena is a first way to calculate the response of a sedimentary layer overlaying a rigid bedrock \citep[e.g.,][]{Thomson1950, Haskell1953}. However, in the presence of complex geological structures, several studies have shown the limitation of 1D numerical simulations, and the necessity of 2D or 3D simulations of seismic wave propagation to reproduce the observed amplifications \citep[e.g.,][]{Kawase1996,Smerzinietal2011,Matsushimaetal2014,Ktenidouetal2016}. \citet{Gelisetal2022} reach the same conclusions for our target site, as they observe that 1D simulations do not reproduce the amplification measured in the Tricastin basin, in particular regarding its maximum amplitude (up to a factor 8) and its broadband spectral character, typical of 3D wave propagation effects \citep{Chavez-Garciaetal2000,CornouBard2003,Bindietal2009,Micheletal2014}. The conclusions of this previous study by \citet{Gelisetal2022} form the motivation for the use of 3D numerical simulations in our work. 3D simulations have become more affordable lately thanks to the rapid increase of computational resources and to the development of dedicated software, using in particular spectral-element methods \citep[SEM, e.g.,][]{KomatitschVilotte1998,DeMartin2011,Trinhetal2019}. These developments led to many applications, notably in sedimentary basins with complex geometries that require detailed simulations \citep[e.g.,][]{Komatitschetal2004, Maufroyetal2015, Maufroyetal2016,Chaljubetal2010,Chaljubetal2015,Paoluccietal2015,Thompsonetal2020,DeMartinetal2021, Panzeraetal2022}. These simulations however require an accurate knowledge of the subsurface, both in terms of geometry and of seismic properties \mbox{($S$-wave} velocity $V_S$, $P$-wave velocity $V_P$, density ${\rho}$, and $S$- and $P$-wave attenuation factors $Q_S$ and $Q_P$). As a consequence, these numerical approaches also rely on field data and geophysical surveys in order to constrain numerical models and define reliable input parameters for the simulations. On the other hand, seismic data (earthquake and noise recordings) are also essential to provide observations which the outputs of the simulations can be compared and calibrated with. Most studies involving 3D numerical simulations for site effect estimation rely on layered models \citep[e.g.,][]{TabordaBielak2013,Maufroyetal2016,DeMartinetal2021, Panzeraetal2022}. In these models, layers are separated by interfaces associated to sharp impedance contrasts, and seismic properties are usually assumed either homogeneous or varying only vertically within each layer. The geometry of the interfaces is often derived from geological knowledge (borehole data, field campaigns, surface mapping or interpreted cross-sections), from the interpretation of active seismic migrated images, or from H/V spectral ratio analysis, while seismic properties are often estimated from a limited number of local measurements such as Ambient Vibration Analysis (AVA) and extrapolated to the entire layers, assuming vertical but no lateral variations within the layers \citep[e.g.,][]{Manakouetal2010, Molinarietal2015, Cushingetal2020, Panzeraetal2022}. In the present paper, we propose an alternative approach that uses ambient-noise surface-wave tomography (ANSWT) for building 3D seismic models of the target area. ANSWT is a tomographic method that has proven very successful for imaging shear-wave velocities $(V_S)$ at the crustal and lithospheric scales \citep[e.g.,][]{Shapiroetal2005} and more recently at smaller, basin scales \citep[e.g.,][]{Boueetal2016,Chmieletal2019}. Taking full advantage of the deployment of dense seismic networks, this passive seismic imaging technique is particularly attractive. ANSWT has the advantage to provide a quantitative model of shear-wave velocity structure in 3D, including both vertical and lateral variations, at relatively low cost. However, ANSWT provides smoother models compared to active seismics and to the geology-based layered models usually considered in numerical simulations. The aim of this paper is to investigate the use of standard ANSWT models for the numerical estimation of seismic amplification and for site effect assessment. More precisely, the question we address here is the following: what are the potential and limitations of 3D numerical simulations based on standard ANSWT models to assess seismic amplification in complex sedimentary basins? The paper is organized as follows. Section~\ref{sec2} presents the target zone, its geological context, and the seismic data acquired in the frame of the DARE project. Section~\ref{sec3} reminds the steps of the standard ANSWT workflow used to build a 3D $V_S$ model of the target site, and explains how we derive other seismic properties ($V_P$, ${\rho}$, $Q_S$, $Q_P$) in order to build a full 3D seismic model. In Section~\ref{sec4}, we present our numerical results and the seismic amplification predicted in the 3D model, which we compare with observed, earthquake-based SSR, as well as with noise-based SSRn and with 1D approximations. Finally, we discuss the results in Section~\ref{sec6}, highlighting the potential of using ANSWT models for the seismic characterization of sedimentary basins, but also underlining some limitations in their use for the numerical estimation of site effects, which leads us to propose perspectives for future work. \section{Target site and data}\label{sec2} The DARE project targets the area of the Tricastin Nuclear Site (TNS) in the French Rh\^{o}ne Valley (Figure~\ref{fig1}). The TNS is located on the Messinian Rh\^{o}ne canyon that was dug about 6~million years ago during the Messinian Salinity Crisis (MSC) by the paleo-Rh\^{o}ne river in older geological formations. In this area, the geological basement is mainly\unskip\break constituted by hard and thick (several hundreds of meters) lower Cretaceous limestones (Barremian), so-called ``Urgonian'' limestones, which we consider to be the reference \textit{bedrock} unit in the region. These Urgonian limestones are overlain by more detrital lower to upper Cretaceous (Aptian to Turonian) formations (sands, sandstones, and marls), and then by Tertiary marine and continental detrital formations. These units are called hereafter ``post-Urgonian Cretaceous (and/or Tertiary) formations'' in order to distinguish them both from the Urgonian bedrock and from post-Messinian sediments. Indeed, after the MSC, the canyon has been filled with Pliocene sediments of marine (sands and clays) and continental (fluviatile conglomerates) origins, nowadays covered by the Rh\^{o}ne lower-to-recent Quaternary terraces. In 2019, when the DARE project was initiated, the local geology of the Messinian canyon remained poorly documented in the region of Tricastin. \citet{Gelisetal2022} provide some first insights about the canyon rims and the local subsurface characteristics based on borehole data, geological study, and local 1D geophysical characterization campaigns (H/V and AVA measurements). This first study provided local knowledge about the characteristics (thickness and $V_S$ velocities) of the sedimentary canyon in-filling and of the underlying bedrock at two sites in the area. A first site is located 2--3~km south of the TNS (seismic station E1/BOLL in Figure~\ref{fig1}), on top of the sedimentary basin. At this site, \citet{Gelisetal2022} show that the base of the canyon reaches a depth of at least 500~m, incising---or, at least, lying directly on top of---high-velocity Urgonian limestones. The 2nd site is station G6/ADHE on a nearby outcrop of Urgonian limestones, which \citet{Gelisetal2022} characterize as a hard rock site with an estimated $V_S30$ (the average shear-wave velocity in the first 30~m) of about 2000~m/s. At larger depth, the 1D $V_S$ profiles show velocity rapidly increasing with depth and reaching about 3000~m/s beyond 50~m depth. \citet{Gelisetal2022} present the various criteria for considering G6/ADHE as a good reference for SSR calculations. It is worth noting that $V_S$ profiles obtained at ADHE and BOLL are consistent (i) with available geological data \citep[e.g.,][]{Bagayoko2021,DoCoutoetal2024} and (ii) with each other in terms of mean velocities at depth, therefore giving a reference velocity of $V_S$ ${\cong}$ 3000~m/s for the deep Urgonian\unskip\break substratum. \begin{figure*} \includegraphics{fig01} \vspace*{-4pt} \caption{\label{fig1}Simplified geological map of the target zone and localisation of the two deployments (red circles: nodal array, blue triangles: broadband network). Background colors correspond to main geological units (see legend). The extent of the map corresponds to the domain of interest for numerical simulations (excluding 15-km-wide margins).} \end{figure*} From this first knowledge of the canyon, a 10~km $\times$ 10~km area surrounding the imprint of Pliocene and Quaternary sediment deposits around the TNS\unskip\break (Figure~\ref{fig1}) was targeted. This extension allows us to embed the edges of nearby outcrops of Cretaceous series incised by the canyon, that constitute the basement of the canyon sedimentary in-filling. Most of this target zone is located in a heavily industrialized area, including the widespread TNS, a hydroelectric dam, towns, several railroads and a highway. The associated anthropogenic activity controls the distribution of high-frequency noise sources locally \citep{Gisselbrechtetal2023}. Two complementary seismic campaigns were carried out in the framework of the DARE project \citep{Fromentetal2022a}. The first campaign consisted of deploying 400 3-component seismic nodes over a 10~km $\times$ 10~km area for one month (red dots in\unskip\break Figure~\ref{fig1}). This campaign targeted the recording of seismic ambient noise. A second campaign consisted of deploying about 50 broadband stations over the same target area for more than eight months and targeted the recording of seismicity (including teleseismic events, regional, and local seismicity). These two datasets \citep{Pilzetal2021,Fromentetal2023b} are presented in detail in a data paper \citep{Fromentetal2022a} and are publicly available (see section \textit{Data and software availability}). The present study exploits the ambient noise data recorded by the 400-node network. In the methodological approach adopted in this work, seismicity data recorded by the 50-station network will only be used here to compare our numerical estimations of seismic amplification with observations. Three-component Geospace GSX nodes have been used for the dense nodal experiment whose design is shown in Figure~\ref{fig1}. The average inter-node distance is about 800~m over the area. About half of the stations spaced 200--250~m apart are used to form a denser grid located south of the TNS. Similarly, two dense east--west lines following two roads are located north of the TNS. Information about data completeness, noise levels, and overall data quality is available in \citet{Fromentetal2022a}. A detailed characterization of seismic noise sources recorded by the nodal network is given in \citet{Gisselbrechtetal2023}. \begin{figure*} \vspace*{-4pt} \includegraphics{fig02} \caption{\label{fig2}Regionalization of the tomographic domain in sub-areas based on the clustering of local dispersion curves.} \end{figure*} \section{Soil model from noise-based medium characterization}\label{sec3} \subsection{3D $V_S$ model from ANSWT}\label{sec3.1} As mentioned in the introduction, the dataset of continuous seismic noise recorded by the nodal array was processed using a standard ANSWT workflow to obtain a 3D model of shear-wave velocity, as presented in \citet{Fromentetal2022b}. In order to bring elements for the discussion of the results obtained by using this model, we detail the workflow used in this previous work hereafter: \begin{enumerate}[(1)] \leftskip=-1.2pc \item Cross-correlations of seismic signals between all pairs of stations of the nodal array (except 6 nodes located outside the dense $10\times 10$~km grid, and 24 nodes with unusable data, resulting in 69,192 valid cross-correlations). Continuous signals were cross-correlated over 30-min-long time windows after spectral whitening, and stacked over the one-month duration of the acquisition. The cross-correlations of all three components (N,E,Z) provided the full cross-correlation tensor. Inter-station cross-correlations were then rotated in terms of radial (RR) and transverse (TT) components assuming straight inter-station paths. The vertical components (ZZ) of the cross-correlations were used for exploiting Rayleigh waves and the transverse (TT) for exploiting Love waves. \item Semi-automatic picking of fundamental-mode group-velocity dispersion curves using frequency-time analysis \citep[FTAN,][]{Dziewonskietal1969,Levshinetal1989}. A statistical quality control of the picked dispersion curves was used to reject outliers falling outside two standard deviations of the distribution of picked values. After this quality\unskip\break control, a total of 17,031 Love dispersion curves between 0.4 and 3~Hz and 29,719 Rayleigh dispersion curves between 0.35 and 6.5~Hz were kept for tomography (i.e., 25\% and 43\% of the full dataset, respectively). The tomography is therefore considered to be well constrained up to 3~Hz (by both Love and Rayleigh data) and partially constrained up to 5--6~Hz (only by Rayleigh data). \item Frequency-dependent 2D traveltime tomography \citep{Barminetal2001} in order to convert inter-station dispersion curves into local dispersion curves, i.e.\ build group velocity maps. This step was performed via a linearized inversion involving regularization in the form of norm damping and lateral smoothing. The choice of these regularization parameters plays a role in the resolution of the final model. \item Inversion of local group-velocity dispersion curves into local 1D $V_S$ profiles using a Neighbourhood Algorithm \citep{Sambridge1999,Mordretetal2014}. Here the use of a global optimization scheme allowed for a statistical exploration of the model space and provided an average of best-fitting models for each 1D $V_S$ profiles, which were then linearly interpolated into a 3D $V_S$ model. \end{enumerate} The tomographic process was guided by assumptions derived from the geological knowledge at the time of this first imaging. In particular, for the 1D depth inversion (step 4), the expectation of a strong velocity contrast between sediments and bedrock led to a parameterization of the 1D profiles consisting in two smooth layers (represented by splines functions) potentially separated by a velocity discontinuity (if required by the data). This parameterization was adapted locally within sub-areas defined by clustering the local dispersion curves after the 2D tomography stage (step 3), using a data-driven \textit{K-means} algorithm \citep{MacQueen1967}. The resulting four sub-areas are shown in Figure~\ref{fig2} and turn out to well coincide with previous geological knowledge. Area \#1 (in blue in Figure~\ref{fig2}) corresponds to the deepest parts of the basin. Area \#2 (in green) corresponds to the shallower northwestern edge of the basin, including the northwesternmost corner where Urgonian limestones are outcropping. Area \#3 (in yellow) corresponds to the eastern edge and its outcrops of post-Urgonian/pre-Messinian formations (limestones, sands, sandstones, and marls of Aptian to Miocene ages). Finally, Area \#4 corresponds to a complex zone, including the so-called Lapalud island, formed by post-Urgonian Cretaceous units (limestones, sandstones and marls of Aptian to Turonian ages) expected to have lower seismic velocities than the hard Urgonian limestones \citep[e.g.,][]{Bagayoko2021,DoCoutoetal2024}. In addition to adapting the parameterization of the 1D inversion in these sub-areas, the sub-arrays were also used to perform FK analysis and estimate an average phase velocity curve that was used to constrain the 1D depth inversion of local group velocity curves within each sub-area. Besides the assumption of a two-layer model with smooth velocity variations within each layer, we did not impose any strong constraint to the inversion. This preliminary 3D $V_S$ model is therefore mainly data driven, its purpose being precisely to investigate how such a model, based on seismic noise only via a \textit{blind} ANSWT workflow, can be used to predict seismic amplification in the basin. This remains a fairly open question, considering that the ANSWT procedure uses phase information only (amplitudes are discarded by spectral whitening in the cross-correlations). In order to use the ANSWT model for numerical simulations, the 3D $V_S$ volume is extrapolated outside the tomographic domain, both laterally and vertically, to define seismic properties in the full simulation domain, which extends laterally by 15~km from the limits of the domain of interest (\textit{SEM domain} in Figure~\ref{fig2}), and vertically down to 30~km depth.\unskip\break Because, in a first time, we want to rely on our data-driven ANSWT model and avoid making strong \mbox{assumptions} about the subsurface, lateral extrapolation is performed invariantly at a given depth. Downward vertical extrapolation consists in a smooth transition from the bottom of the tomographic model (1000~m below ground surface) to a constant velocity of 3000~m/s at 1500~m depth. This value of $V_S = 3000$~m/s is chosen for consistency with prior measurements in the area \citep{Gelisetal2022}. It also roughly corresponds to the maximum velocity of the tomographic model at 1~km depth, and therefore ensures a positive velocity gradient with respect to depth. Outside the basin, we also perform a vertically invariant extrapolation upwards in order to assign seismic properties to topographical heights above basin level. Figure~\ref{fig3} gives several views of the obtained 3D $V_S$ model. Figure~\ref{fig3}a gives a 3D view that includes N--S and E--W vertical cross-sections through the $V_S$ model (left and back panels, first colorbar), as well as the 1200-m/s iso-velocity surface that can be associated to the interface between sediments and bedrock (second colorbar). Figure~\ref{fig3}b gives a 3D view of the numerical model, after extrapolation of the ANSWT model outside the tomographic domain over the full domain of interest (but excluding 15-km-wide margins, which just consist of further extrapolation).\unskip\break Figure~\ref{fig3}c shows a vertical N--S cross-section along the expected axis of the Rh\^{o}ne paleo-canyon \citep{Gelisetal2022,Fromentetal2022b}, passing through the TNS and through station E1/BOLL [on which we will focus later on to illustrate our results, and compare them with those of \citealp{Gelisetal2022}]. Figure~\ref{fig3}d shows a vertical E--W cross-section in the south of the tomographic domain, through the Lapalud island and station D0 (on which we will also focus later on to illustrate our results). \begin{figure*} \includegraphics{fig03} \bcaption{\label{fig3}}{2}{3D Vs model. (a)~3D view showing the interface between sediments and bedrock \citep[adapted from][]{Fromentetal2022b}. (b)~3D view of the numerical model after extrapolation outside the tomography domain (excluding 15-km-wide margins). Green dots represent nodes and broadband stations. (c)~South--north vertical cross-section along the axis of the paleo-Rh\^{o}ne canyon. (d)~West--east vertical cross-section through the northern part of the array and station A4. (e)~West--east vertical cross-section through the southern part of the array and the Lapalud island. In (c--e), the dashed box corresponds to the tomography domain. Red dots and blue triangles represent nodes and broadband stations located within 500~m of the sections (see Figure~\ref{fig2} for the location of the sections).} \end{figure*} As already described by \citet{Fromentetal2022b}, the ANSWT $V_S$ model agrees well with the main geological structures expected in the area. In particular, the range of estimated shear-wave velocity values are consistent with previous studies \citep[e.g.,][]{Gelisetal2022}, with values ranging from 500 to about 1200--1400~m/s in the sediments, and from about 1700--2000 to 3000~m/s in the underlying bedrock. Moreover, the velocity discontinuity between the two layers (which roughly corresponds to the 1200-m/s iso-velocity surface shown in Figure~\ref{fig3}a) coincides well with the expected depth of the paleo-canyon, at least for its deeper parts, with a north--south axis corresponding to the paleo-Rh\^{o}ne and a southwestern branch corresponding to the paleo-Ard\`{e}che (Figure~\ref{fig3}a). Between these two branches, the model also depicts the so-called Lapalud island with higher velocities reaching shallower depths, representing post-Urgonian Cretaceous units\unskip\break (Figures~\ref{fig3}a,c). However, when looking in more detail, we notice that some features of the model are questionable. In particular, the model does not display high velocities reaching the surface, even in regions where we expect surface outcrops of Cretaceous formations (especially in the northwestern corner and on the eastern edge of the basin, Figures~\ref{fig3}a, b). Instead, the model exhibits a lower-velocity layer ($500 < V_s < 1000$~m/s) with a thickness of at least 200~m over the entire domain (Figures~\ref{fig3}b--d). This low-velocity layer is not restricted to the basin, which is therefore not well delimited laterally. We will see in the following that this lack of basin edges will have consequences in terms of seismic amplification. In spite of these limitations, the fact that the 3D $V_S$ model derived from ANSWT well depicts the expected large-scale geometry of the Messinian paleo-valley at depth, including some complex structures such as the Lapalud island, motivates us to use this model for numerical simulations, in order to look at the seismic amplification that it may generate. This, however, first requires the definition of other seismic properties in the considered 3D volume. We will now present how we define these properties, and we will distinguish between properties that are constrained by the same seismic noise data as the $V_S$ model (namely $Q_S$) and other properties that are estimated by other means ($V_P$, density, $Q_P$). \subsection{Estimation of shear-wave quality factors $Q_S$}\label{sec3.2} Besides ANSWT, the seismic noise recorded by the dense array is processed via Q-SPAC analysis \citep{Prietoetal2009} to estimate seismic attenuation parameters $(Q_S)$, following the methodology of \citet{Boxbergeretal2017}. To this end, the study area is subdivided into 15 sub-arrays, each of them containing 15 to 25 seismic nodes. \begin{figure*} \includegraphics{fig04} \caption{\label{fig4}Calibration of a site-specific $V_S/Q_S$ relationship using Q-SPAC inversion results.} \end{figure*} In a first step, the Extended Spatial AutoCorrelation (ESAC) method is adapted to first obtain mean 1D $V_S$ profiles below each sub-array by joint inversion of Rayleigh wave dispersion and H/V \mbox{spectral} ratios, and then estimate frequency-dependent Rayleigh-wave attenuation factors from the mean 1D $V_S$ profiles \citep{Ohorietal2002,Boxbergeretal2011}. In a second step, individual 1D $Q_S$ profiles are obtained from the inversion of the Rayleigh-wave attenuation coefficients by constraining the $V_{S}$ profiles to the values obtained in the first step, following \citet{Xia2014}. It is worth noting that the seismic attenuation discussed here does not distinguish between intrinsic and scattering attenuation. Furthermore, while the Rayleigh-wave attenuation factors of the input data depend on frequency, the $Q_S$ parameter of the obtained layered model is assumed to be frequency-independent. Finally, following \citet{Xiaetal2002}, we disregard the contributions of P waves on the Rayleigh-wave attenuation\break factors. In the end, 15 layered 1D profiles representative of average $V_S$ and $Q_S$ values as a function of depth are derived, one for each of the 15 sub-arrays of the nodal network. $V_S$ profiles are found to be in general agreement with the ANSWT model (with\unskip\break deviations smaller than 15\%). These profiles are used to calibrate a relationship between $V_S$ and $Q_S$ values in the form of a 6-order polynom [after \citealp{TabordaBielak2013}, see Figure~\ref{fig4}] which is then used to derive a 3D $Q_S$ model from the 3D ANSWT $V_S$ model. This is expected to provide more realistic and site-specific $Q_S$ values than assuming generic relationships (e.g., $V_S/Q_S=10$), or other relationships from the literature calibrated for other sites [e.g., \citealp{TabordaBielak2013}, see Figure~\ref{fig4}]. The $Q_S$ values obtained from our site-specific relationship range from 21 for small $V_S$ values in shallow sediments to 360 for high $V_S$ values in the deep substratum\break (Table~\ref{tab1}). \begin{table*} \caption{\label{tab1}Ranges of seismic properties in the two layers of the 3D $V_S$ model} \begin{tabular}{ccccccccccc} \thead & \multicolumn{2}{c}{$V_S$~(m/s)} & \multicolumn{2}{c}{$V_P$~(m/s)} & \multicolumn{2}{c}{Density~(kg/m$^3$)} & \multicolumn{2}{c}{$Q_S$} & \multicolumn{2}{c}{$Q_P$}\\ \cline{2-3} \cline{4-5} \cline{6-7} \cline{8-9} \cline{10-11} & min & max & min & max & min & max & min & max & min & max \\ \endthead 1st layer (sediments) & \0500 & 1400 & 2000 & 2675 & 2075 & 2225 & \021 & 130 & 260 & 350 \\ 2nd layer (bedrock) & 1700 & 3070 & 3075 & 5300 & 2300 & 2640 & 170 & 340 & 420 & 760 \botline \end{tabular} \end{table*} \subsection{Other seismic parameters ($V_P$, $Q_P$, density)}\label{sec3.3} Unlike $V_S$ and $Q_S$ properties, $P$-wave velocity $V_P$ and quality factor $Q_P$, and density parameters are not (or poorly) constrained by our ambient noise data, which are typically dominated by surface waves. \textit{In situ} geotechnical measurements performed in the vicinity of the TNS down to approximately 30~m depth provide us with local estimations of $V_P$, $V_S$, and density \citep{Moiriat2019}. These measurements show $V_P/V_S$ ratios of up to 8 in very shallow sediments (Quaternary limons and alluvions forming thin (${<}$20~m) layers and lenses, which are not included in the ANSWT due to its lack of sensitivity to these shallow layers), and of about 4 for values of $V_S\approx500$~m/s in the blue marls encountered below 15 to 20~m depth. It is worth noting that these values of $V_S\approx500$~m/s are in very good agreement with the $V_S$ values of the shallow part of the basin in the ANSWT model which, according to the smallest wavelength considered in the tomography (${\approx}$100~m), we can regard as effective $V_S$ values for the basin sediments in the first 30 to 50~m. Based on these geotechnical measurements, we calibrate a site-specific $V_P/V_S$ relationship of the form {\begin{equation}\label{eq1} \frac{V_{P}}{V_{S}}= 1.73+a\,\mathrm{e}^{-bV_{S}}, \end{equation}}\unskip with $a = 9.46$ and $b = -2.82\times 10^{-3}$. This \textit{ad hoc} relationship is designed such as to yield a $V_P/V_S$ ratio of 1.73 for large $V_S$ values (deep bedrock), corresponding to the usual assumption of a Poisson's solid, while the calibration yields $V_P/V_S$ ratios of about 4 for $V_S\approx500$~m/s, as indicated by \textit{in situ} geotechnical data. In lack of sufficient density measurements, we use Gardner's empirical law [\citealp{Gardneretal1974}, equation~(7); \citealp{Brocher2005}, equation~(2)] in order to relate density to $P$-wave velocity $V_P$. This law is supposed to be valid for sedimentary rocks such as limestones. We verify that the density values obtained for the range of seismic velocities encountered in our ANSWT model roughly coincide with expected values for the known lithologies in the area (limestones, sandstones, and marls), as well as with the abovementioned geotechnical measurements. Finally, we must define values for the $P$-wave quality factor $Q_P$ in our 3D model, although this parameter is expected to have very little impact on the amplification of the SH excitation that we will simulate. In lack of any constraint on these $Q_P$ values, we derive them from $Q_S$, $V_S$, and $V_P$ values, under the assumption that the compressibility quality factor is much larger than the shear quality factor [$Q_\kappa\gg Q_\mu$, \citealp{DahlenTromp1998}, equation~(9.59), p.~350]. Table~\ref{tab1} summarizes the ranges of values of the different seismic properties in the final 3D seismic model. \section{Numerical estimation of seismic\newline amplification}\label{sec4} \subsection{Simulation of seismic wave propagation}\label{sec4.1} Simulations of seismic wave propagation in the constructed 3D model are performed using the EFISPEC software \citep{DeMartin2011}, which makes use of a time-domain spectral-element method (SEM) solving the 3D equation of motion in visco-elastic media. A hexahedral mesh is designed for simulations valid up to 5~Hz, which is about the maximum frequency used to constrain the ANSWT model (at least with Rayleigh data, see Section~\ref{sec3.1}). The simulation domain extends laterally 15~km further than our domain of interest (Figures~\ref{fig1} and~\ref{fig2}) and vertically down to 30~km, in order to mitigate parasite reflections on domain boundaries, where the absorbing condition is based on a paraxial approximation \citep{Stacey1988}. This leads to a total of 6.5~M elements, with an element size that varies between 100~m and 300~m from the shallow to the deeper parts of the model. In order to reproduce the assumptions underlying SSR calculations, the simulated source is a vertically-incident plane wave with SH polarization in the $X$- (east--west) or $Y$- (north--south) direction, injected at 5~km depth. We perform two simulations, one for each polarization of the plane wave, for a duration of 60~s. Each simulation costs 3500 CPUh and is parallelized on 480 cores. The source time function is a low-pass-filtered Dirac delta function, filtered\unskip\break below 5~Hz, such as to remain in the domain of validity of the simulation (and of the ANSWT model) and to be able to visualize and exploit simulated outputs directly, without filtering them at the post-processing stage. Figure~\ref{fig5} shows snapshots of the simulated wavefield on the free surface at different times. The vertically-incident plane wave first arrives after approximately 2.2~s of simulation in the northwestern corner, which has slightly higher velocities than the rest of the domain (Figure~\ref{fig5}a). At 2.5~s, amplitudes are saturated over the entire domain while the plane wave is reflecting on the free surface (Figure~\ref{fig5}b). At 3.5 and 4~s (Figures~\ref{fig5}c, d), high amplitudes are visible in the deepest parts of the basin, reflecting 1D basin amplification due to vertical body-wave resonance, but also outside the tomographic domain, due to topographic effects, in particular on the margins that are not constrained by the tomography. After 4.5~s (Figures~\ref{fig5}e, f), late reverberations are slightly visible in the basin (likely due to reflections on the edges of the basin), but also in the margins of the model outside the tomographic domain (due to topographic effects). These snapshots show that waves are not efficiently trapped in the basin, as we would expect in a well-delimited basin \citep[e.g.,][]{Chaljubetal2015,DeMartinetal2021}. Instead, waves escape in the low-velocity shallow layer, resulting in an amplification over the entire domain, in particular in the margins of the simulation domain that are not constrained by the tomography. \begin{figure*} \includegraphics{fig05} \vspace*{-4pt} \caption{\label{fig5}Snapshots of the simulated wavefield ($v_y$ component excited by a vertically-incident SH plane wave polarized in the $Y$ direction) on the free surface at different times. The black contour corresponds to the iso-depth 375~m of the velocity discontinuity between sediments and bedrock, delineating the deepest part of the basin. The dashed line corresponds to the tomographic domain. The full movie of wave propagation is provided in Supplementary Material.} \vspace*{-4pt} \end{figure*} Based on these first observations and despite some limitations, we will now continue to analyse our simulation results, especially in terms of predicted amplification, with the objective to better understand the limitations of our ANSWT model, which we just started to point out, and thus to draw perspectives about how such ANSWT models could be improved to better reproduce basin amplification. We propose to start by looking into the seismograms simulated at the locations of known stations, including E1/BOLL that has been investigated by \citet{Gelisetal2022}, as well as A4 that serves as a reference station for the calculation of our empirical, earthquake-based SSR (see Supplementary Material for more details). We shall specify here that station A4 is located on the same kind of Urgonian outcrop and presents similar H/V characteristics as station G6/ADHE that was used as a reference station by \citet{Gelisetal2022}. Figure~\ref{fig6} shows seismograms (Figure~\ref{fig6}a) simulated at stations E1/BOLL (in blue) and A4 (in orange) and their amplitude spectra (Figure~\ref{fig6}b). In Figure~\ref{fig6}b, thin lines represent raw spectra and thick lines represent spectra smoothed using Konno-Omachi smoothing with a bandwidth $b= 40$ \citep{KonnoOhmachi1998} which are then used to compute spectral ratios\unskip\break (Figure~\ref{fig6}c). As expected, the signal simulated at station E1/BOLL has a long duration because of wave reverberations in the basin (only the first 11~s of the signals are shown in Figure~\ref{fig6}, but reverberations generate non-negligible amplitudes---of the order of about 10\% of the maximum peak amplitude---over the entire 60-s duration of the simulation). The signal simulated at station A4, on the other hand, does not correspond to our expectations for a reference rock site precisely because it also contains late reverberations of similar, non-negligible amplitudes as station E1. Apart from a high-amplitude arrival at about 4~s, probably due to a vertical body-wave reflection in the basin, the ground motion in E1 does not seem much amplified compared to the one in A4. In the Fourier domain (Figure~\ref{fig6}b), the amplitude spectra at stations E1 and A4 differ mostly by the frequencies of their respective peaks at low frequencies: around 0.5~Hz for E1 (blue line in Figure~\ref{fig6}b), which corresponds well to the 1D resonance of the basin \citep{Gelisetal2022}, and around 1~Hz for A4 (orange line in Figure~\ref{fig6}b), which rather corresponds to the 1D resonance of the shallow low-velocity layer extrapolated outside the actual extent of the basin (if we consider a relation between resonance frequency $f_{\mathrm{res}}$ and layer thickness $h$ of the form $f_{\mathrm{res}} = V_S^{\mathrm{avg}}/(4h)$, with $V_S^{\mathrm{avg}}$ the average $V_S$ velocity in the layer). \begin{figure*} \includegraphics{fig06} \vspace*{4pt} \caption{\label{fig6}(a)~Seismograms simulated at stations E1/BOLL (in blue) and A4 (in orange) and their amplitude spectra (b), for a vertically-incident plane wave source polarized in the $Y$ direction ($v_y$ components are shown). Green curves show the reference signal considered for computing amplification, which corresponds to the recording of the incident plane wave at 10~km depth. (c)~Spectral ratios.} \vspace*{2pt} \end{figure*} As a consequence, we cannot consider the signal simulated at A4 as a reference for estimating amplification via spectral ratios: the resulting amplification would be dramatically underestimated, especially at frequencies corresponding to the 1D resonance of the shallow low-velocity layer (e.g., 1~Hz, see red curve in Figure~\ref{fig6}c). Instead, we consider a deep reference point, located at 10.050~km depth, i.e., at the same distance from the depth at which the plane wave source is injected (5~km) as the free surface in the basin (which has an average elevation of 50~m asl). The seismogram extracted from this deep reference point is time-windowed, such as to retain only the incident plane wave, which has travelled in a homogeneous medium between its injection\unskip\break point at 5~km depth and this reference point at 10.050~km depth, therefore yielding a clean signal (green curves in Figure~\ref{fig6}). The amplitude spectrum of this windowed signal is then multiplied by 2, which corresponds to the amplitude spectrum of a point that would be located on the free surface. In doing so, we therefore define a reference amplitude spectrum that is equivalent to the theoretical response of a point on the free surface above a homogeneous halfspace that would have the same seismic properties as the deep part of our 3D numerical model (including visco-elastic attenuation, which makes this deep reference spectrum different from the spectrum of the injected source time function) and that is consistent with the wavefield simulated at the free surface (since it is extracted from the simulated wavefields, and not computed separately). We refer to the spectral ratio computed with this deep reference as the amplification function (AF), as opposed to the standard spectral ratio (SSR) computed with a\unskip\break reference point located on the free surface. We have verified that AF and SSR are identical in the case of a surface reference point located above a homogeneous subsurface, and far away from basin edges (clean numerical reference). As a consequence, our numerical amplification functions remain comparable to empirical SSRs computed with respect to a clean empirical reference on the field, which is the case of station A4. Figure~\ref{fig6}c shows the resulting amplification functions for station E1 (in blue), with a clear peak around 0.55~Hz corresponding to basin resonance, and for station A4 (in orange), with a peak at 1~Hz, related to the resonance of the artificially extrapolated low-velocity layer. In the following, we will now analyse in more detail these numerical amplification functions, with a particular focus on points located in the basin (since we know that rock sites outside the basin are misrepresented in our model, and therefore may provide less relevant amplification\unskip\break results). \subsection{Numerical amplification}\label{sec4.2.} In this section, we will analyse in three different ways the numerical amplification simulated in our ANSWT model. First, we will compare our numerical amplification functions and their spectral characteristics to empirical SSR measurements based on earthquake data recorded during the deployment of the DARE broadband network. We will exemplify this comparison at two specific locations in the basin, which we have identified as representative of two different signatures that we find of particular interest for our scope, especially when compared to a 1D SH response approximation, which we will also provide as another point of comparison. Second, we will comment on the spatial variability of the numerical amplification predicted in our model, for a few example frequencies. Third, we will have a closer look at this spatial variability at low frequency, and more precisely at the spatial pattern of the (mis)match between our numerical amplification and two empirical amplification estimations based on earthquake and noise data. \subsubsection{Local comparison of numerical amplification with earthquake-based SSR and 1D\newline amplification} Figure~\ref{fig7} shows the amplification at stations E1 (a) and D0 (b) as a function of frequency (see Figure~\ref{fig1}\unskip\break for the location of these two stations). Purple lines represent the empirical, earthquake-based SSR computed with respect to reference station A4 and their uncertainties (purple dashed lines). A catalog of these empirical SSR for all stations of the DARE broadband network is provided in Supplementary Material, together with details on the computation of these SSR. The green line corresponds to the SH amplification computed in a 1D model extracted from the 3D model below the station of interest (right panels), using the \textit{gpsh} function of the Geopsy software \citep{Watheletetal2020}. The numerical amplification of the horizontal component $(\mathrm{AF}_H)$ computed in the 3D tomographic model is shown as dark blue and red dashed lines for vertically-incident plane wave sources polarized in the E--W $(\mathrm{PW}_X)$ and N--S $(\mathrm{PW}_Y)$ directions, respectively. The set of grey curves corresponds to the amplification of the horizontal component for plane wave sources polarized in different directions ($\mathrm{PW}_\theta$, with $\theta$ the polarization azimuth), showing the variability of the amplification with respect to source polarization. These curves are obtained by a linear combination of the two wavefields simulated for plane wave sources polarized in the $X$ and $Y$ directions, in order to reproduce, by linearity of the wave equation with respect to the source, the wavefield generated by a plane wave source polarized in any given direction. \begin{figure*} \includegraphics{fig07} \caption{\label{fig7}Amplification spectra at stations E1 (a) and D0 (b). Dark blue and red dashed lines: numerical amplification of the horizontal component $(\mathrm{AF}_H)$ in the 3D tomographic model, for vertically-incident plane wave sources polarized in the E--W $(\mathrm{PW}_X)$ and N--S $(\mathrm{PW}_Y)$ directions, respectively. Grey curves: amplification of the horizontal component for plane wave sources polarized in different directions $(\mathrm{PW}_{\theta})$. Green line: 1D SH amplification $(\mathrm{AF}_{\mathrm{SH}}^{1\mathrm{D}})$ computed in a 1D model extracted from the 3D model below the station of interest (right panels, where the dashed line below 1~km depth corresponds to the part of the model that is not constrained by ANSWT and is extrapolated to a constant $V_S=3000$~m/s). Purple lines: empirical, earthquake-based SSR$_EQ$ computed with respect to reference station A4, and their uncertainties (purple dashed lines). Note that uncertainties do not appear for high frequencies where the SSR is often constrained by a single event, which does not enable to compute a standard deviation (see Supplementary Material for details).} \end{figure*} The two stations considered in Figure~\ref{fig7}, E1 and D0, exhibit two different amplification patterns that can be linked to their respective locations in the basin, and eventually to the corresponding subsurface structures. Below station E1, we know from previous studies that the paleo-canyon is quite deep [about 560~m according to \citealp{Gelisetal2022}], and that a thick pile of Pliocene sediments probably lies directly on top of the hard substratum of Urgonian limestones. This results in a significant observed\unskip\break amplification (up to $7\pm2$) above the main resonance frequency of the canyon (about 0.5~Hz for this location/depth, according to \citet{Gelisetal2022}. Here it is worth noting that the numerical results do not retrieve this observed level of amplification (with a maximum of about $4\pm 1$), nor its broadband character: while the observed amplification exhibits a plateau above 0.55~Hz, which is commonly interpreted as the signature of 3D wave propagation effects [due, in particular, to laterally propagating surface waves generated at the edges of the basin, see e.g.\ \citealp{CornouBard2003,Bindietal2009}], the numerical amplification computed from 3D numerical simulations (envelope of gray curves) remains quite close from the 1D resonance peaks predicted in a local 1D model extracted below the considered station (green curve). This suggests that the ANSWT model does not generate significant amounts of laterally propagating surface waves within the basin, which can be easily understood from our observations of the model (Figure~\ref{fig3}) and of the wavefield (Figure~\ref{fig4}), where we already noted the lack of clear lateral basin edges, resulting in a lack of wave trapping within the basin. Nevertheless, the numerical amplification computed in the 3D ANSWT model is not strictly identical to a 1D amplification. In particular, it exhibits non-negligible variations with respect to source polarization for some frequencies, including in the vicinity of the main resonance frequency (0.45--0.55~Hz). While this source-related variability would probably deserve a more detailed investigation \citep[e.g.,][]{Maufroyetal2017}, it is a clear indication that the ANSWT model well includes some 3D characteristics of the basin that cannot be captured by purely 1D approximations. Moreover, it should also be noted that the local 1D profile extracted from the 3D ANSWT model below station E1 and the associated 1D amplification (green curves in Figure~\ref{fig7}a) are very similar to the results obtained by \citet[][their figure~8d]{Gelisetal2022} for this location. This suggests that the ANSWT well played its role in terms of local $V_S$ estimation, just as well as the local Ambient Vibration Analysis performed by \citet{Gelisetal2022}, but with the added value of imaging the lateral variations of the 1D $V_S$ structure in the basin. \begin{figure*} \includegraphics{fig08} \vspace*{-4pt} \caption{\label{fig8}Maps of the amplification of the horizontal component at 0.5, 1, 2, and 4~Hz, for a vertically-incident SH plane wave polarized in the $Y$ direction (N--S). The black line represents the iso-375-m depth contour of the interface between the two layers considered in the 1D depth inversion. Margins outside the tomographic domain have been masked.} \vspace*{-4pt} \end{figure*} Regarding station D0, it is located at the border of the Lapalud island, where the subsurface has a complex structure due to the presence of the late lower Cretaceous units in between the Pliocene\unskip\break sediments and the Urgonian bedrock \citep{Bagayoko2021,DoCoutoetal2024}. These units are made of sandstones and marls, and have intermediary seismic velocities compared to the Pliocene and Urgonian lithologies, which in the $V_S$ profiles is represented by a smooth gradient of velocity progressively increasing with depth. This results in a lower level of observed amplification (up to $3 \pm 0.5$). This amplification also displays a broadband character as a function of frequency, but for a different reason as previously: in this case, the amplification plateau above 0.5~Hz is due to the smooth gradient of increasing\unskip\break velocity with depth, and not to 3D wave propagation effects (since 1D amplification also shows this broadband character). Interestingly, both 3D and 1D numerical results show a good agreement with the observations in this case, which again suggests that the amplification is not dominated by 3D propagation effects at this location, but by 1D vertical resonance, so that the 1D approximation is sufficient to explain the observations. It should be noted, however, that this agreement is only possible because the tomography well played its role in terms of (3D) characterization of the local (1D) velocity structure. \subsubsection{Spatial variability of the numerical\newline amplification} One of the interests of numerical simulations is to provide access to the spatial variability of the amplification in the basin. Figure~\ref{fig8} shows amplification maps at different frequencies. At 0.5~Hz, most of the amplification corresponds to the deepest part of the basin (delineated by the iso-375-m depth contour of the interface between sediments and bedrock), as expected. At higher frequencies, however, amplification does not behave as expected. At 1~Hz, amplification is larger outside the basin than inside, which is likely due to the resonance of the shallow 200-m-thick low-velocity layer at this frequency. At 2~Hz, the amplification pattern is less clear but seems to concern mainly the slopes of the basin. At 4~Hz, results simply do not seem to be relevant, and are probably dominated by the topographic effects identified in the wavefield snapshots (Figure~\ref{fig5}). In the following, we focus our analysis to the low-frequency part of the amplification, in the frequency range 0.3--0.7~Hz, related to the deepest parts of the basin (note however that this frequency range is wide enough to cover a relatively wide range of canyon depths, from about 300~m to 700~m). More specifically, we propose to look at the ratio between our numerical amplification and the observed amplification in this frequency range, for all broadband stations (earthquake-based SSR$_{\mathrm{EQ}}$, Figure~\ref{fig9}) and all nodes (noise-based SSRn, Figure~\ref{fig10}). \begin{figure*} \includegraphics{fig09} \vspace*{-5pt} \caption{\label{fig9}Mean ratio between numerical amplification and earthquake-based SSR$_{\mathrm{EQ}}$ between 0.3 and 0.7~Hz (note the logarithmic color scale).} \vspace*{-7pt} \end{figure*} \begin{figure*} \vspace*{3pt} \includegraphics{fig10} \vspace*{5pt} \caption{\label{fig10}Mean ratio between numerical amplification and noise-based SSRn \citep{Gisselbrechtetal2023} between 0.3 and 0.7~Hz (note the logarithmic color scale).} \vspace*{4pt} \end{figure*} \subsubsection{Spatial variability of the misfit between\newline numerical amplification and observed SSR} \paragraph{Comparison to earthquake-based SSR$_{\mathrm{EQ}}$} Figure~\ref{fig9} shows the spatial distribution, for all broadband stations, of the mean ratio between our numerical amplification $\mathrm{AF}_{\mathrm{num}}$ (log-average of the gray curves in Figure~\ref{fig7}) and the observed earthquake-based amplification SSR$_{\mathrm{EQ}}$ with respect to reference station A4 in the frequency range 0.3--0.7~Hz. According to the logarithmic color scale, dark red and dark blue colors represent stations where the numerical amplification underestimate and overestimate the observations by a factor of 4, respectively. Empty triangles correspond to stations for which we could not estimate empirical SSR$_{\mathrm{EQ}}$, due to an insufficient signal-to-noise ratio of the earthquake recordings or to a lack of data (see Supplementary Material for details). The dashed-dotted line represents the surface imprint of the canyon rim, as interpreted from geological mapping (Figure~\ref{fig1}), which provides a visual guide to locate the results with respect to the basin. The dotted line delimits the domain constrained by the tomography, outside which laterally-invariant extrapolation is performed (Figure~\ref{fig3}) and makes the results hardly interpretable, in particular at rock sites (e.g., A4, G6/ADHE). Within the tomography domain, however, the ratio between numerical and observed amplification is always comprised between 0.5 and 2, and is close to 1 for some stations. This suggests that the ANSWT model, although not perfect and still improvable, is able to explain, in some extent, the amplification observed in this low-frequency range, and in various parts of the basin (where, incidentally, the frequency range 0.3--0.7~Hz may include different amplification phenomena depending on the location in the basin and on the underlying/surrounding geological structures). In order to look in more detail at the spatial pattern of the (mis)match between\unskip\break numerical and observed amplification, we now propose to look at the ratio between numerical and noise-based amplification. \paragraph{Comparison to noise-based SSRn} Similarly to\unskip\break Figure~\ref{fig9}, Figure~\ref{fig10} shows the spatial distribution, for all nodes, of the mean ratio between our numerical amplification $\mathrm{AF}_{\mathrm{num}}$ and the observed noised-based amplification SSRn in the frequency range 0.3--0.7~Hz, in which we consider SSRn estimations to be reliable \citep{Gisselbrechtetal2023}. Indeed, several studies have reported a good agreement between noise-based SSRn and earthquake-based SSR in this frequency band \citep[e.g.,][]{LermoChavez-Garcia1994, Bard1999}, mostly because distant sources dominate the ambient noise wavefield at these frequencies \citep[e.g.,][]{Bonnefoy-Claudetetal2006}. In order to further mitigate undesired effects of local and transient low-frequency sources (such as wind), \citet{Gisselbrechtetal2023} computed these SSRn using short time windows (1~min) restricted to quieter nighttime recordings. Note also that these empirical SSRn have been estimated with respect to a node co-located with reference station G6/ADHE. The similarities between Figures~\ref{fig9} and~\ref{fig10} in the vicinity of broadband stations suggests that the two empirical estimations, based on earthquake and ambient noise, are overall consistent, and that the choice of the reference station, A4 or G6/ADHE, has very little effect on the results in this low-frequency band, as also suggested by the SSR ${\cong}$ 1 between the two stations (see Supplementary Material). Thanks to the density of the nodal array, Figure~\ref{fig10} gives a nice spatial view of the (mis)match between\unskip\break empirical and numerical amplification, and enables to look at its spatial pattern in more detail. In particular, we see again that our numerical amplification over-estimates the observations in the northwesternmost corner of the domain, including at rock sites where we know that Urgonian limestones are outcropping, but also in the north-west of the basin (south-west of Pierrelatte), suggesting that the ANSWT model might over-estimate sediment thickness in this area. In contrast, amplification seems under-estimated in the deepest part of the basin, along the north--south axis of the paleo-Rh\^{o}ne canyon, as well as along the paleo-Ard\`{e}che in the south-west of the domain. Knowing that ANSWT provides a similar estimate of the canyon depth at station E1 as AVA \citep{Gelisetal2022}, this\unskip\break underestimation could suggest a lack of 3D wave propagation effects in these parts of the model, all the more that they also correspond to regions susceptible to be affected by laterally propagating surface waves generated at basin edges. Interestingly, one of the areas where the match between numerical and empirical amplification seems the most satisfying corresponds to the Lapalud island in the South of the domain. Again, this suggests that the ANSWT model is successful in estimating the subsurface velocity structure in this region, despite its relative complexity due to the presence of heterogeneous post-Urgonian Cretaceous units. It also suggests that local 1D effects may be dominant in the amplification at this location, which may be intuitively understood as the fact that the top of the island itself is not affected by waves that are diffracted by (/away from) the island. The ratio between numerical and noise-based amplification is also quite satisfying north-west from Lapalud, where the basin is shallower \citep{DoCoutoetal2024}, as well as on its eastern margin, which suggests that the post-Urgonian Cretaceous and Tertiary units that form this margin do not have the same properties as the rock sites located on Urgonian outcrops (A4, G6/ADHE), and experience some amplification compared to these hard rock sites. These lateral variations of bedrock properties, rarely taken into account in numerical studies, are important to consider, as they affect not only the local response, and thus the choice of a reference for SSR calculation \citep{Steidletal1996}, but also the impedance contrast at the sediment/bedrock interface. \section{Discussion}\label{sec5} The aim of this paper was to investigate if and how a standard ANSWT model could be used for the numerical estimation of seismic amplification in sedimentary basins. In light of our results, we can state that ANSWT is surely a method of choice for depicting the large-scale velocity structure of the subsurface, but that a standard, purely data-driven, ANSWT workflow alone is visibly not sufficient to generate seismic models that well reproduce the observed amplification in complex sedimentary basins where 3D wave propagation effects are significant. In particular, we clearly saw the consequences of the lack of well-marked lateral basin edges, and thus of laterally propagating surface waves in the basin. There are several reasons for this lack of shallow lateral contrasts in the model, in relation to the sensitivity and resolution capacity of the surface wave data used for tomography. First, the upper layers of the subsurface (first 50~m in the basin, where the minimum wavelength is about 100~m; first 200~m outside the basin, where the minimum wavelength is about 400~m) are poorly constrained by our measurements of the dispersion of fundamental Love and Rayleigh modes in the considered frequency range (0.3--5~Hz). Second, due to this weak data sensitivity and to the relative weight of regularization (smoothing), there is a smearing of low-velocity anomalies from the center to the edges of the basin at the 2D tomography stage, all the more that the edges of the model are less well constrained by data coverage, given the imprint of the array. This understanding is important because it gives precise clues about how to improve future ANSWT applications for site effect assessment in sedimentary basins. In terms of acquisition design, it suggests that the seismic array should probably extend more outside the basin, in order to better constrain basin edges and surrounding bedrock properties. More generally, it suggests that we should seek more information about basin edges, either in the seismic data themselves or elsewhere. In this preliminary study, indeed, we deliberately adopted a purely data-driven ANSWT approach, in order to assess its potential and limitations. But if surface wave data alone are not sufficient to generate sharp lateral variations, then we should turn to other sources of information which we first need to collect and then to integrate in the tomographic process in order to constrain our models. In our case, the surface imprint of the Messinian paleo-canyon is known from geological mapping (Figure~\ref{fig1}), field campaigns, and borehole data. Besides, active seismic profiles are also available in the area. The interpretation of these profiles (time-domain migrated sections), together with geological field campaigns and borehole data, allows for a high-resolution identification of geological horizons \citep{Bagayoko2021,DoCoutoetal2024}. These horizons, however, cannot be considered as a strict ground truth, because the inference of their depth depends on $P$-wave velocity values assumed for time-to-depth conversion. Moreover, these geological horizons do not necessarily correspond to seismic impedance contrasts, as discussed in \citet{Fromentetal2022b}, and need to be interpolated over the entire domain of interest. Nevertheless, these horizons could be used as prior information in order to constrain the geometry of subsurface layers in the tomographic process. This will be the subject of further work following directly from this study, which enables us to pinpoint more precisely where this prior information should be introduced: namely both at the 2D tomography stage, where lateral coherency could be enforced using these geological constraints rather than blind smoothing, and at the 1D inversion stage, where interface depths could be constrained directly. Finally, the computation of H/V ratios from the ambient noise data recorded by the dense nodal array also results in the estimation of an interface related to an impedance contrast that controls seismic amplification, and this information could also be worth taken into account in the tomography and in our simulations. In practice, however, how to take this interface into account is not completely obvious, given that (i)~the interpretation of this interface in terms of depth also implies some assumption about $V_S$ values in the sediments, (ii)~this interface may not correspond to geological stratigraphic horizons, or not always to the same horizons, and may even not be continuous \citep{Fromentetal2022b}. In all these cases, the consideration of such prior information in the tomographic process would be best handled in the frame of probabilistic methods, such as trans-dimensional tomography \citep[e.g.,][]{Bodinetal2012,Galettietal2017} or Bayesian approaches \citep[e.g.,][]{Luetal2020,Nouibatetal2022}. This would at least require the use of state-of-the-art imaging techniques, and probably some dedicated methodological developments as well, in particular to properly evaluate and propagate uncertainties throughout the whole ANSWT workflow, starting from reliable prior data uncertainty related to dispersion curve picking, and ending with robust posterior uncertainties on model parameters, which could then be used in our numerical simulations to explore the variability of the predicted amplification within the range of model uncertainties. As a direct perspective of this work, we are currently investigating such probabilistic methods in order to improve our tomographic model of the Tricastin basin, and explore the sensitivity of our numerical amplification to input data and model features, within realistic uncertainty ranges. Another important---and yet open---question related to data and model uncertainties concerns the frequency range in which we can reliably predict seismic amplification, based on a certain set of informations to design our numerical models. Finally, the lack of 3D characteristics in the ANSWT model likely comes in part from the 1D approximation which the standard 2-step ANSWT procedure relies on for interpreting surface-wave dispersion in terms of velocity structure as a function of depth. This assumption could be overcome by using other imaging techniques, such as 1-step 3D ANSWT \citep[e.g.,][]{Zhangetal2018} or full waveform inversion \citep[FWI, e.g.,][]{VirieuxOperto2009,Luetal2020}. FWI, in particular, would enable one to fully consider 3D wave propagation, as well as higher-mode surface waves, and eventually body waves (if present in the cross-correlations), that would bring additional information on the impedance contrast between sediments and bedrock, and potentially provide higher-resolution models with sharper lateral variations associated to basin edges. We would recommend the investigation of such waveform tomography techniques for building seismic models for site effect assessment as a long-term research track, keeping in mind that FWI applications in complex subsurface settings are always challenging, and that noise-based FWI is still the subject of current developments and raises challenges of its own, in particular regarding the effect of noise sources on parameter estimation \citep[e.g.,][]{Trompetal2010,Sageretal2018}. Standard ANSWT models such as the one considered in this study could be used as initial models for further waveform inversion approaches. \section{Conclusion}\label{sec6} In this paper, we have used a 3D $V_S$ model of a complex sedimentary basin built via a standard ANSWT workflow to perform numerical simulations of seismic wave propagation, and compute ground motion amplification within the basin. The tomographic model well depicts the main geological structures of the basin via its lateral variations of seismic velocities at depth, but it also suffers from limitations related to surface wave sensitivity and resolution capabilities, especially in its shallow part. In particular, this ANSWT model lacks of clear basin edges in order to efficiently trap seismic waves and generate significant amounts of surface waves susceptible to propagate laterally in the basin. As a result, the numerical amplification predicted in this ANSWT model presents some 3D characteristics, such as non-negligible variations with respect to the polarization of the incident seismic wave, but also lacks of other expected 3D features, such as the broadband signature of observed amplification at locations affected by laterally propagating surface waves coming from basin edges. On the other hand, in basin sites where a 1D response is sufficient to explain the local amplification, the numerical amplification predicted in the ANSWT model show a good agreement with the observations, including in complex regions of the basin where lateral variations of the deeper subsurface velocity structure must be taken into account in order to reproduce the local 1D response. The tomography also enables one to image lateral variations of bedrock properties which are rarely taken into account in numerical site effect estimations. These observations still concur in considering noise-based tomography as a promising method for building seismic models for site effect assessment in sedimentary basins, provided we can compensate the lack of data sensitivity to shallow and sharp lateral variations, and improve the reconstruction of basin edges. Our results therefore show the potential of the approach, while identifying its limitations and giving perspectives for future work, which should aim at considering other imaging paradigms (e.g., noise-based waveform inversion) and at including prior information from other sources of geological and geophysical data into the tomographic process, in order to better constrain the geometry of the basin, in particular its surface boundaries in the shallow subsurface which is poorly constrained by surface-wave data alone. \section*{Declaration of interests} The authors do not work for, advise, own shares in, or receive funds from any organization that could benefit from this article, and have declared no affiliations other than their research organizations. \section*{Data and software availability} The datasets acquired in the frame of the DARE project are freely available via their doi: \href{https://doi.org/10.14470/L27575187372}{10.14470/L27575187372} \citep{Pilzetal2021}, \href{https://doi.org/10.15778/RESIF.3T2019}{10.15778/RESIF.3T2019} \citep{Fromentetal2023a}, and \href{https://doi.org/10.15778/RESIF.XG2020}{10.15778/RESIF.XG2020} [\citealp{Fromentetal2023b}, see also \citealp{Fromentetal2022a}, for details]. The EFISPEC software \citep{DeMartin2011} is freely available at \url{https://gitlab.brgm.fr/brgm/efispec3d} (last accessed 25 January 2024). \section*{Funding} This work is part of the DARE project (\textit{Dense ARray for site effect Estimation}, \url{https://anr.fr/Projet-ANR-19-CE31-0029}, last accessed 25 January 2024), funded by the French National Research Agency (\textit{Agence Nationale pour la Recherche}, ANR, grant number ANR19-CE31-0029) and by the German Research Foundation (\textit{Deutsche Forschungsgemeinschaft}, DFG, project number 431362334). ISTerre is part of Labex OSUG@2020 (ANR10 LABX56). Numerical simulations have been performed with the EFISPEC software developed at BRGM [\citealp{DeMartin2011}, \url{https://gitlab.brgm.fr/brgm/efispec3d}, last accessed 25 January 2024]. This work was granted access to the HPC resources of TGCC under the allocation 2022-gen13461 made by GENCI, and of the GRICAD infrastructure (\url{https://gricad.univ-grenoble-alpes.fr}, last accessed 25 January 2024), which is supported by Grenoble research communities. \CDRGrant[ANR]{ANR19-CE31-0029} \CDRGrant[DFG]{431362334} \CDRGrant[ANR]{ANR10 LABX56} \CDRGrant[TGCC]{2022-gen13461} \section*{Acknowledgements} The geological context of the studied area was discussed and completed with the help of Damien Do Couto, Nanaba Bagayoko, Ludovic Mocochain and Jean-Loup Rubino. We warmly acknowledge them. \back{} \section*{Supplementary data} Supporting information for this article is available on the journals website under \printDOI\ or from the author. \CDRsupplementaryTwotypes{supplementary-material}{\cdrattach{crgeos-243-suppl.pdf}} \CDRsupplementaryTwotypes{supplementary-material}{\cdrattach{File-S1.mp4}} \CDRsupplementaryTwotypes{supplementary-material}{\cdrattach{File-S2.pdf}} \bibliographystyle{crgeos} \bibliography{crgeos20230720} \refinput{crgeos20230720-reference.tex} \end{document}