### Microlocal analysis of seismic body waves and linearized inverse problems – Imaging

A separated scale representation of variations in medium properties leads to the introduction of a background model and a contrast or reflectivity. This in turn leads to the introduction and formulation of a separated inverse problem for estimating these, and is tuned to particular (relative high) frequency windows in the data. This formulation encompasses inverse scattering and ‘wave-equation’ tomography and is focussed on the singularities in the data. In this framework, the phases or constituents may be multiply scattered but can interact with the unknown contrast only once.

One key complication in developing proper imaging techniques for transmission and reflection tomography and inverse scattering are the generation of caustics. This is typically associated with regions of low wavespeeds, for example, hot spots and mantle upwellings. A related complication is associated with anisotropy. To avoid misinterpretation of images, it is critical to understand uniqueness of reconstruction (even though full data coverage will not be available), the possible generation of (singular) artifacts (beyond standard resolution analysis), and partial reconstruction (including illumination). One may think of these techniques, which have grown out of traditional tomography and (Kirchhoff) migration, as providing a ‘leading-order’ identification of variations in material properties and structures in Earth’s interior.

**Inverse scattering**

Even though maximal (5-dimensional) surface data seem to define an overdetermined inverse problem, this is, in fact, not the case in the presence of caustics. In reality, this inverse problem is further complicated by the availability of data (for example, linear arrays) often constituting a lower-dimensional set.

*Characterization of the normal operator, extensions*. The general characterization of the normal operator associated with inverse scattering is still incomplete with so-called*I*^{p,l }classes appearing in the analysis. To ensure that false discontinuities cannot be generated, we invoked the*Bolker condition*[29], which may be violated in particular if the acquisition loses a dimension. Possible mitigation of artifacts caused by such a violation is a current subject of research. We pioneered the introduction of extensions leading to (microlocally) invertible single scattering operators. The propagation of singularities by these extensions are described by canonical transformations with associated evolution equations. With these extensions the imaging of reflectivity in inaccurate background models can be undone; beyond that, we obtained an evolution equation (and an associated Hamiltonian) for extended images under changing background velocities [72]. The extensions play a key role in obtaining estimates for the background model. We constructed explicit extensions via local directional decomposition subject to the introduction of a so-called*curvilinear double-square-root (DSR) condition*[46,74]. Our standard extension involves interior curvilinear offset; however, as an alternative we constructed a (wave-equation) angle transform, which can be related (microlocally) to reflection coefficients [41].

*Annihilators, wave-equation reflection tomography*. We introduced annihilators (extending the notion of differential semblance to allow for the formation of caustics in the unknown background) of the data [55]. Annihilators characterize the range of the single scattering operator [52] and provide a procedure for conditional source-receiver continuation. The range can be used as a criterion for background model recovery, which is under certain conditions unique. We designed an iterative method for this recovery and analyzed the associated sensitivity kernels [57].We adapted and applied our methodology using the generalized Radon transform jointly for PP and PS reflections to constrain TI parameters in a sedimentary environment in the Norwegian sector of the North Sea [48].

*Generalized Radon transform (GRT)*. We have constructed a GRT for inverse scattering incorporating the asymptotic inverse (parametrix) of the normal operator, in anisotropic media, assuming that the equations for elastic waves are of principal type [22] and carried out an associated resolution analysis [17]. Indeed, for the purpose of imaging Earth’s deep interior, the GRT can be applied to global network data, in principle; one does not need uniform spatial sampling as we demonstrated by invoking quasi-Monte Carlo techniques [10]. Moreover, upon adapting the order of the Fourier integral operator representing the GRT, we obtained stability and a method for statistical inference of singularities [62]. We developed an approach to precise partial reconstruction from partial data with the GRT and illumination correction using a frame of wavepackets [70].- We applied the GRT to image the
*D*^{ll }layer beneath Central America using ScS (precursors), and provided a new view of the post-perovskite phase transition near the CMB. With thermodynamic properties of phase transitions in mantle silicates, we interpreted the images and estimated in situ temperatures. Exploiting the presence of multiple phase transitions, we estimated the core heat flux in the coldest region (a site of deep subduction) and in a neighborhood away from it [59]. - We also applied the GRT to image the hot spot underneath Hawaii using underside reflections, that is, SS precursors [90]. The detection of a thermal plume has been seismically inconclusive. We explained the depths of the discontinuities in the transition zone below the Central Pacific with olivine and garnet transitions in a pyrolitic mantle. The presence of a thermal anomaly west of Hawaii suggest that hot material accumulates near the base of the transition zone before being entrained in flow toward Hawaii.
- The modern view of Earths lowermost mantle considers a
*D*^{ll}*D*^{ll}region [107]. The occurrence of multiple interfaces is inconsistent with expectations from a thermal response of a single isochemical post-perovskite transition but can be explained with post-perovskite transitions in differentiated slab materials. Candidate compositions for a seismically detectable post-perovskite transition at pressures less than the CMB include mid-oceanic ridge basalt (MORB) and harzburgite components of differentiated oceanic lithosphere transported to the lowermost mantle by subduction.

- We applied the GRT to image the

*Reverse-time-continuation (RTM)-based inverse scattering*. We constructed a new inverse scattering transform based on RTM [93], also in the anisotropic case [chapter in book 5], which automatically removes the so-called low-frequency artifacts. This transform applies to data generated by a finite set of (point) sources (events). The formulation in curvilinear coordinates is straightforward and an elastic wave-equation driven imaging procedure is obtained. However, we also designed a fast, approximate algorithm with associated estimates using the dyadic parabolic decomposition of phase space and bypassing time stepping [preprint 1]. In this algorithm it is straightforward to extract information about the local incidence angle and thus generate angle-dependent images.- We can remove the knowledge of the source in the case of converted waves. In [99] we obtained an inverse scattering procedure which is effectively bilinear in the data and replaces traditional receiver functions, with the resolving power of RTM-based inverse scattering.
- We extended our inverse scattering procedure to data from (known) deterministic to (known or unknown, ambient) noise sources [98] ensuring statistical stability. Exploiting mode conversions, we obtained a related result using the coda (only) of an earthquake.

These procedures, in combination with advanced wavefield retrieval techniques, provide a new platform honoring wave dynamics, for example, for high-resolution imaging and improved detection and characterization of interfaces in the upper mantle, including the Moho, the lithosphere-asthenosphere boundary and the 660km discontinuity, as well as low velocity features underneath North America using USArray (including temporary and ‘flexible’ array) data.

*Downward-continuation-based inverse scattering*. The downward continuation concept is based on directional wavefield decomposition and its implementation has its origin in the so-called double-square-root (DSR) equation. We introduced curvilinear coordinates and a pseudodepth via a Riemannian metric. We constructed the normal operators for the above mentioned extended imaging procedures in this downward continuation approach [55]. We also constructed explicitly the corresponding evolution equation (and thus an extended exploding reflector model) which provided a new insight in the geometry of extended imaging in the presence of caustics via the introduction of extended isochron rays [82].

** (Wave-equation) tomography**

*Condition for uniqueness, multiple scattering*. In the case of transmission tomography it is well known that the presence of caustics obstructs the unique recovery of velocities. With multiply broken rays and using the scattering relation as the data – with a finite set of favorably located (unknown) reflectors – the linearized inverse problem and (partial) reconstruction allows for more general background models in which certain caustics may form.*Anisotropy, shear-wave splitting*. We extended our analysis of wave-equation transmisison tomography to shear-wave splitting intensity tomography. We designed, again, an iterative method in a heterogeneous environment [64], and carried out computational experiments for SKS and SKKS phases to demonstrate the possibility of constraining the deformation and flow beneath the Ryukyu arc, Japan. Flow modelling is an essential part of the regularization in estimating local anisotropy.