Author Archives: ichang

Directed Spanning Set (DSS), Global Injectivity of ReLU Networks

Multi-scale physics of seismic waves

In building a detailed understanding of the propagation and scattering of waves, that is, their interaction with complex media and boundaries, I have followed a hierarchical description of medium variations or coefficients, starting with a smooth component, a C^1,1 component, general interfaces with regularity constraints, and a class of random fluctuations.

• Multi-scale ray-wave duality. The description of waves up to C^1,1 coefficients is naturally obtained introducing the notion of wavepackets. Wavepackets were introduced in harmonic analysis and later identified as curvelets. A wavepacket can be viewed as a localized ‘fat’ plane wave, or a ‘quantum’. Wavepackets follow a scale description (dyadic parabolic decomposition of phase space) and can be used to form a transform pair (generalizing in some sense the notion of beamforming). We obtained a construction of the solution to the wave equation in two steps: in step (i) we smooth the medium on the scales corresponding with the scales of the wavepackets and obtain an approximate construction using geometry (rays on each scale); in step (ii) we incorporate the scattering via a Volterra integral equation. We obtain a concentration of wave packets result; the degree of concentration reveals the smoothness of the medium [65]. The construction also answers precisely questions about how asymptotic representations are related to full wave solutions.

  We generalized the construction and description from scalar waves to elastic waves, and addressed the question of the minimum regularity in Lame´ parameters for which polarized waves can be distinguished, which appears to be C^1,1[86].

• Multi-scale wave-equation transmission tomography. The above construction can be utilized to analyze the physics behind cross-correlation-based wave-equation transmission tomography, that is, the ‘finite-frequency’ analogue of the geodesic ray transform, by taking Fréchet-type derivatives of our solution construction [78], and resolves the seeming insensitivity of the transform along the unperturbed rays via a multi-scale analysis. Indeed, in the case of smooth wavespeeds, and the limit of infinite bandwdith, the geodesic ray transform is recovered [43].
• Scattering series, interfaces. Reducing the regularity of the medium further, leads us to the introduction of (general) interfaces (conormal singularities). We succeeded to track the order of scattering through the regularity of the wave constituents in the case of C^1,ε (ε > 0) coefficients [108]. In this case, the reflected wave is more regular than the incident wave in a Sobolev sense. Via this analysis, we also obtained the first meaningful introduction of a scattering series for waves (beyond the case of planarly layered media).
• Random fluctuations, moments, cross correlations. Via random fluctuations we reduce the regularity of the medium further. Assuming a wavepacket decomposition of the source, we studied waves in random media essentially via a description in terms of moments. Through the parabolic scaling underlying the construction of wavepackets, one obtains Itô-Schrödinger diffusion processes which reproduce the moments, in the limit of fine scales. We incorporated deterministic interfaces as well, moderately tilted relative to the orientation of the wavepacket source, and obtained theoretically an extended notion of enhanced backscattering [101]. We found that the backscattering cone does not depend on the tilt.

  We utilized the results pertaining to the Itô-Schrödinger diffusion processes and partly coherent waves to study the retrieval of Green’s functions from ‘field-field’ cross correlations in the strongly cluttered regime. We obtain filters via Wigner transforms which satisfy a set of transport equations, and determine their asymptotic behaviors which contain information about the statistics of the random fluctuations. We found that the retrieval is partial only, contains these (blurring) filters, and depends on the directionality of the source [69,102].

  • We studied the ‘field-field’ correlations from data acquired on the Tibetan plateau when focussing on surface waves. We assessed the retrieval from coda and from ambient noise by comparing the results with direct surface wave observations from earthquakes [54] (including directionality [73]). We carried the comparison over to the estimation of phase velocities and enabled a fusion of the results while enlarging the bandwidth. Here, we employed a semi-classical analysis point of view inserting a small parameter in the variations in medium properties in the surface normal coordinate.
  • In highly discontinuous and random planarly layered media, we developed an understanding of tunneling, hidden equivalent medium averaging, and microlocal attenuation which cannot originate from standard visco-elastic models [1,2,5,8].

On planetary scale, that is, a ball, we studied a random radially layered medium including random angular fluctuations on distinct scales as well as deterministic discontinuities. We obtained a stochastic analysis in terms of modes and partly coherent waves excited by an interior source. Through cross correlations of the trace (field at the surface) we obtained a (virtual) data set which can be utilized for the imaging of the deterministic discontinuities.

• Directional wavefield decomposition and propagation, sufficiently smooth coefficients. In the case of sufficient medium regularity one can interchange time locally with a (curvilinear) spatial coordinate (‘pseudodepth’) to describe the propagation of waves in terms of an evolution or one-way wave equation. In certain non-trivial wavespeed profiles this can be accomplised exactly, using spectral theory (modes) and special functions [23]. Under certain conditions, we can incorporate and analyze the role of ‘evanescent’ wave constituents [24]. In general, this technique is desgned for relatively high frequency wave propagation. In appropriate coordinates, we derived a novel paraxial equation and obtained estimates for its solution in terms of an energy norm converging to the full wave solution as the frequency becomes large. This is accomplished, again, via parabolic scaling. The ‘up-down’ coupling can be incorporated through a generalized Bremmer series [11]. These consistent paraxial equations enable the development of fast algorithms [15], and are of current interest in, for example, long-range propagation of and imaging with teleseismic body waves.
• Anisotropy. We studied medium anisotropy from the point of view of polarized waves and a corresponding local slowness surface characterization (for example, by introducing anellipticity) and expansions [20], and the point of view of classification, characterization and construction of waveforms via Picard-Fuchs equations [7].

The behaviors and phenomena described above provide fundamental insight in developing multi-scale methods to estimate local continuity, anisotropy, regularity etc. of medium properties and substructures from broadband seismic observations, and can lead to a deeper understanding of the coupling of geological and geodynamical processes, (complex) phase boundaries etc. and the interaction (and, hence, probing) with seismic waves. In this context, we obtained some basic results in the case of viscous flow in the mantle wedge induced by a subducting lithospheric plate; we approximated the subduction of the Philippine Sea plate beneath Eurasia, explored the effects of different rheological and thermal constraints, and modelled the evolution of the finite strain ellipses [61].

Algorithms

Following the multi-scale physics of waves a family of new computational algorithms has emerged:

• We designed and implemented a fast discrete almost symmetric wavepacket transform pair in three dimensions [80]. This transform pair can also be utilized for wave-based signal processing. Waves (as well as coefficient models containing interfaces) can be compressed with this transform.
• Following the dyadic parabolic decomposition of phase space underlying the construction of wave packets, we developed a fast algorithm (with error estimates) for a class of so-called Fourier integral operators comprising propagators, asymptotic solution operators of the wave equations and (inverse) scattering transforms [92]. This algorithm is based on accurate low-rank separated representations, prolate spheroidal wave functions, and ray-geometric computations. (From a physics of waves point of view, in the case of one-way wave propagators, these representations are related to generalizations of phase screens [34,35].) In the application to the half-wave equation, we obtained a theoretical, multi-scale description of the Gouy phase which affects the extraction of travel time in long-range propagation. We designed a special procedure to incorporate the (numerical) generation of general caustics[104].
• We estimated the regularity of the kernel of the Volterra operator in the above mentioned construction of weak, full-wave solutions in media of limited smoothness. This regularity determined the optimal choice of quadrature to discretize the integral. Moreover we estimated the approximation of solutions by finite sets or clusters of wave packets to complete the discretization of the Volterra equation and obtain a finite matrix representation [95].

We constructed frames of wave packets produced by parabolic dilation, rotation and translation of (a finite sum of) Gaussians and obtained asymptotics on the analogue of Daubechies’ frame criterion [110]. We also obtained a frame of wave-atom-like packets generated from pure Gaussians; these can be identified as multi-scale Gaussian beams. We developed an associated parametrix construction and algorithm.

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-calledcurvilinear 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 region of enhanced (seismologically inferred) heterogeneity bounded by the core mantle boundary (CMB) and an interface some 150-300 km above it, with the latter attributed to the post-perovskite phase transition. We carried out a deep-Earth exploration with ScS and SKKS wavefields, which probe this region from above and below, respectively, which suggested that this view be modified. Inverse scattering revealed that, beneath Central America and East Asia, two areas known for subduction of oceanic plates deep into Earths mantle, seismic reflectors exist above the conventional D^llregion [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.

• 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.

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-calledcurvilinear 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 region of enhanced (seismologically inferred) heterogeneity bounded by the core mantle boundary (CMB) and an interface some 150-300 km above it, with the latter attributed to the post-perovskite phase transition. We carried out a deep-Earth exploration with ScS and SKKS wavefields, which probe this region from above and below, respectively, which suggested that this view be modified. Inverse scattering revealed that, beneath Central America and East Asia, two areas known for subduction of oceanic plates deep into Earths mantle, seismic reflectors exist above the conventional D^llregion [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.

• 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.

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-calledcurvilinear 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 region of enhanced (seismologically inferred) heterogeneity bounded by the core mantle boundary (CMB) and an interface some 150-300 km above it, with the latter attributed to the post-perovskite phase transition. We carried out a deep-Earth exploration with ScS and SKKS wavefields, which probe this region from above and below, respectively, which suggested that this view be modified. Inverse scattering revealed that, beneath Central America and East Asia, two areas known for subduction of oceanic plates deep into Earths mantle, seismic reflectors exist above the conventional D^llregion [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.

• 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.