Inversion sismique 3D des formes d'onde complètes pour des cibles terrestres complexes

Seismic imaging of onshore targets is very challenging due to the 3D complex near-surface-related effects. In such areas, the seismic wavefield is dominated by elastic and visco-elastic effects such as highly energetic and dispersive surface waves. The interaction of elastic waves with the rough topography and shallow heterogeneities leads to significant converted and scattering energies, implying that both accurate 3D geometry representation and correct physics of the wave propagation are required for a reliable structured imaging. In this manuscript, we present an efficient and flexible full waveform inversion (FWI) strategy for velocity model building in land, specifically in foothill areas.Viscoelastic FWI is a challenging task for current acquisition deployment at the crustal scale. We propose an efficient formulation based on a time-domain spectral element method (SEM) on a flexible Cartesian-based mesh, in which the topography variation is represented by an accurate high-order geometry interpolation. The wave propagation is described by the anisotropic elasticity and isotropic attenuation physics. The numerical implementation of the forward problem includes efficient matrix-vector products for solving second-order elastodynamic equations, even for completely deformed 3D geometries. Complete misfit gradient expressions including attenuation contribution spread into density, elastic parameters and attenuation factors are given in a consistent way. Combined adjoint and forward fields recomputation from final state and previously saved boundary values allows the estimation of gradients with no I/O efforts. Two-levels parallelism is implemented over sources and domain decomposition, which is necessary for 3D realistic configuration. The gradient preconditioning is performed by a so-called Bessel filter using an efficient differential implementation based on the SEM discretization on the forward mesh instead of the costly convolution often-used approach. A non-linear model constraint on the ratio of compressional and shear velocities is introduced into the optimization process at no extra cost.The challenges of the elastic multi-parameter FWI in complex land areas are highlighted through synthetic and real data applications. A 3D synthetic inverse-crime illustration is considered on a subset of the SEAM phase II Foothills model with 4 lines of 20 sources, providing a complete 3D illumination. As the data is dominated by surface waves, it is mainly sensitive to the S-wave velocity. We propose a two-steps data-windowing strategy, focusing on early body waves before considering the entire wavefield, including surface waves. The use of this data hierarchy together with the structurally-based Bessel preconditioning make possible to reconstruct accurately both P- and S-wavespeeds. The designed inversion strategy is combined with a low-to-high frequency hierarchy, successfully applied to the pseudo-2D dip-line survey of the SEAM II Foothill dataset. Under the limited illumination of a 2D acquisition, the model constraint on the ratio of P- and S-wavespeeds plays an important role to mitigate the ill-posedness of the multi-parameter inversion process. By also considering surface waves, we manage to exploit the maximum amount of information in the observed data to get a reliable model parameters estimation, both in the near-surface and in deeper part.The developed FWI frame and workflow are finally applied on a real foothill dataset. The application is challenging due to sparse acquisition design, especially noisy recording and complex underneath structures. Additional prior information such as the logs data is considered to assist the FWI design. The preliminary results, only relying on body waves, are shown to improve the kinematic fit and follow the expected geological interpretation. Model quality control through data-fit analysis and uncertainty studies help to identify artifacts in the inverted models.L’imagerie sismique du sous-sol à partir de données terrestres est très difficile à effectuer due à la complexité 3D de la proche surface. Dans cette zone, les ondes sismiques sous forme d’un paquet compact de phases souvent imbriquées sont dominées par des effets élastiques et viscoélastiques, couplés aux effets dus à la surface libre qui génèrent des ondes de surface de grande amplitude et dispersives.L’interaction des ondes sismiques avec une topographie plus ou moins complexe dans un contexte de fortes hétérogénéités de la proche surface induit d’importantes conversions des ondes avec de fortes dispersions d’énergie. Il est donc nécessaire de prendre en compte à la fois une représentation tridimensionnelle précise de la topographie et une physique correcte qui rend compte de la propagation du champ d’onde dans le sous-sol au niveau de précision réclamé par l’imagerie sismique. Dans ce manuscrit, nous présentons une stratégie d’inversion des formes d’onde complètes (FWI en anglais) efficace, autonome et donc flexible, pour la construction de modèles de vitesse à partir de données sismiques terrestres, plus particulièrement dans les environnements dits de chevauchements d’arrière pays(foothills en anglais) aux variations de vitesse importantes.Nous proposons une formulation efficace de cette problématique basée sur une méthode d’éléments spectraux en domaine temporel sur une grille cartésienne déformée, dans laquelle les variations de topographie sont représentées par une description détaillée de sa géométrie via une interpolation d’ordre élevé. La propagation du champ d’onde est caractérisée par une élasticité linéaire anisotrope et par une atténuation isotrope du milieu: cette deuxième approximation semble suffisante pour l’imagerie crustale considérée dans ce travail. L’implémentation numérique du problème direct inclut des produits matricevecteurefficaces pour résoudre des équations élastodynamiques composant un système différentielhyperbolique du second ordre, pour les géométries tridimensionnelles rencontrées dans l’exploration sismique. Les expressions explicites des gradients de la fonction écart entre les données et les prédictions sont fournies et inclut les contributions de la densité, des paramètres élastiques et des coefficients d’atténuation. Ces expressions réclament le champ incident venant de la source au même temps de propagation que le champ adjoint. Pour ce faire, lors du calcul du champ adjoint à partir de l’instant final, le champ incident est recalculé au vol à partir de son état final, de conditions aux bords préalablement sauvegardées et de certains états intermédiaires sans stockage sur disques durs. Le gradient est donc estimé à partir de quantités sauvegardées en mémoire vive. Deux niveaux de parallélisme sont implémentés, l’un sur les sources et l’autre sur la décomposition du domaine pour chaque source, cequi est nécessaire pour aborder des configurations tridimensionnelles réalistes. Le préconditionnement de ce gradient est réalisé par un filtre dit de Bessel, utilisant une implémentation différentielle efficace fondée sur la même discrétisation de l’espace du problème direct et formulée par une approche d’éléments spectraux composant un système linéaire symétrique résolu par une technique itérative de gradient conjugué. De plus, une contrainte non-linéaire sur le rapport des vitesses de compression et de cisaillement est introduite dans le processus d’optimisation sans coût supplémentaire: cette introductions’avére nécessaire pour traiter les données en présence de faibles valeurs de vitesse proche de la surface libre.L’inversion élastique multi-paramètres en contexte de chevauchement est illustrée à travers des exemples de données synthétiques dans un premier temps, ce qui met en évidence les difficultés d’une telle reconstruction…