Past Seismic Tomography Sessions

Session Recordings

To watch the session recordings, please make sure you are logged in to your SSA account. Then, click the button below to access the session recordings.

WATCH

Thank you to SmartSolo for sponsoring the Virtual Tomography series!

smartsolo.com
smartsolo.com.cn

 

Cutting-edge Methods for Seismic Imaging I

6 October 2020 at 10 AM Pacific
Featured Discussants: Jeroen Ritsema (University of Michigan) and Carl Tape (University of Alaska Fairbanks)

Abstracts

Markov Chain Monte Carlo Regional Tomography Models of Western North America: Findings from Alaska & Southern California

Presenting Author: Elizabeth Berg, University of Utah, Utah, USA

Authors: Elizabeth Berg, University of Utah, Utah, USA; Fan-Chi Lin, University of Colorado Boulder, Colorado, USA;  Vera Schulte-Pelkum, University of Colorado Boulder, Colorado, USA; Weisen Shen, State University of New York at Stony Brook, New York, USA; Amir Allam, University of Utah, Utah, USA

We create isotropic shear velocity models of Southern California and Alaska with sensitivity from the near-surface into the upper mantle and resolve Moho depths. This is possible through data recorded on hundreds of stations provided by regional networks in Southern California in 2015, and from the USArray deployment in Alaska from 2014 to 2018. We use Rayleigh-wave measurements from earthquakes for long-periods (15s-100s), and from ambient noise cross-correlations for shorter periods (2s-24s periods). These measurements include Rayleigh-wave phase velocity via Helmholtz tomography and Eikonal tomography for earthquake and ambient noise data respectively, as well as Rayleigh-wave ellipticity. By combining surface wave measurements with azimuthally-corrected receiver functions we obtain sensitivity to near-station subsurface structure, including horizontal interfaces and sensitivity to regional 3-D subsurface structure. We leverage complementary sensitivity of these datasets in a Markov Chain Monte Carlo joint inversion and are able to quantify model uncertainty and gain understanding of sensitivity. In Southern California, we observe well-known large-scale features, including prominent sediment basins in the Central Valley, Los-Angeles, San Bernardino and Ventura basins. We also create self-consistent Moho depth maps across both regions. In Southern California, we observe thicker crust in the Peninsular and Sierra Nevada Ranges, and thinner crust in the Salton Trough and the Coast Ranges. In Alaska, we resolve every major basin, both shallow and deep, and we image the interaction of the Pacific slab and continental crust and its impact to create more or less prominent volcanism. In Alaska, we find thicker crust beneath the Brooks Range and thinner crust in the Yukon Composite Terrane and in interior Alaska.

Adjoint Tomography of New Zealand’s North Island Using an Automated, Open-source Workflow

Presenting Author: Bryant Chow, Victoria University of Wellington, Wellington, New Zealand

Authors: Bryant Chow, Victoria University of Wellington, Wellington, New Zealand; Yoshihiro Kaneko, GNS Science, Wellington, New Zealand, Kyoto University, Kyoto, Japan; Carl Tape, University of Alaska Fairbanks, Alaska, USA; Ryan Modrak, Los Alamos National Laboratory, New Mexico, USA; John Townend, Victoria University of Wellington, Wellington, New Zealand

We have developed an automated, open-source workflow for full-waveform inversion using spectral-element and adjoint methods (Chow et al., 2020). As a study area we choose the eastern North Island of New Zealand, which encompasses the Hikurangi subduction zone. The chosen domain offers a unique opportunity for imaging material properties near an active subduction zone, due to the availability of well-recorded earthquakes in close proximity to the plate interface. We have performed realistic synthetic inversions using New Zealand source and receiver distributions to determine a set of parameters usable in real-data inversions. We have then undertaken an iterative inversion using 5,000 measurements from a subset of long-period (10–30 s) earthquake waveform data to resolve broad-scale velocity changes with respect to an initial ray-based 3D velocity model of New Zealand. Velocity changes are resolved at shallow crustal depths in tectonically active areas, such as the Taupō Volcanic Zone, regions of geodetically detected slow slip, and in the vicinity of the locked-to-creeping transition within the subduction zone. Improved fits between data and synthetic waveforms motivate an ongoing large-scale inversion using 200 earthquakes recorded on as many as 100 broadband stations. Initial iterations are used to fit long-period waveforms (15–30 s), while progressively shorter periods are introduced to a target period of 3 s, corresponding to features on the order of several kilometers in size. We present ongoing efforts towards an accurate, high-resolution tomographic model of the North Island of New Zealand, which will be used to further understanding of enigmatic tectonic processes.

Seismic Traveltime Tomography Based on Stochastic Voronoi Cells Parameterization: Application from Local to Global Scales

Presenting Author: Hongjian Fang, Massachusetts Institute of Technology, Massachusetts, USA

Authors: Hongjian Fang, Massachusetts Institute of Technology, Massachusetts, USA; Malcolm White, University of Southern California, California, USA; Yehuda Ben-Zion, University of Southern California, California, USA; Rob van der Hilst, Massachusetts Institute of Technology, Massachusetts, USA

The amount of seismic arrival time data has increased dramatically in the past few decades. With the use of machine learning based automatic picking technique, the number is expected to grow more rapidly. However, current travel time tomography methods have to use a reduced data set to homogenize the data distribution and to decrease the computational costs. In this study, we propose a new tomography method that adopts a Poisson-Voronoi cells parameterization, which avoids the use of explicit regularization terms by decomposing the original high-dimensional problem into a series of low-dimensional sub-problems. Moreover, we use a subset of the whole dataset in each sub-problem, which is similar to the mini-batch technique in machine learning. The sub data set can be optimally selected in a way that can make the data sampling more homogeneous while making full use of the available data. We apply the method to seismic arrival times at different scales, from a local scale data set in the Ridgecrest area to a global scale in the Earth mantle. Detailed model features will be discussed in the talk.

Accelerating Global-Scale Full-Waveform Inversion

Presenting Author: Sölvi Thrastarson, ETH Zürich, Zürich, Switzerland

Authors: Sölvi Thrastarson, ETH Zürich, Zürich, Switzerland; Dirk-Philip van Herwaarden, ETH Zürich, Zürich, Switzerland; Lion Krischer, ETH Zürich, Zürich, Switzerland; Christian Boehm, ETH Zürich, Zürich, Switzerland; Martin van Driel, ETH Zürich, Zürich, Switzerland; Michael Afanasiev, ETH Zürich, Zürich, Switzerland; Andreas Fichtner, ETH Zürich, Zürich, Switzerland

Full-waveform inversion (FWI) has proven to produce subsurface images of previously unprecedented accuracies. However, with the requirement of multiple accurate seismic wavefield simulations, FWI is a computationally demanding procedure. Incorporating all available data into an FWI study is already a challenge. With the wealth of seismic data increasing exponentially, one has to think of ways to exploit the available data without hitting computational barriers. In this study, we present a novel approach to both numerical simulations and to the adjoint method which reduces the computational cost of FWI by an order of magnitude.

In smooth media it has been observed that wavefield complexity is highly anisotropic. The wavefield is relatively smooth in the lateral direction with respect to the source location, a feature which can be exploited in the meshing process. By adapting the meshing of the wavefield to the expected complexity, the number of required elements to fully resolve the wavefield can be reduced by an order of magnitude. This results in a mesh optimized for a single source location and medium. Accurate gradients can be computed via the discrete adjoint approach.

We will present a real data example of a global FWI where we were able to recover the large scale features of the Earth at the cost of less than 20 regular wavefield simulations. In the example, the wavefield adapted meshes are used in combination with a stochastic mini-batch optimization scheme where only a fraction of the dataset is used in each iteration. The example shows how global FWI can be done with a large dataset, in a computationally tractable manner.

Uncertainty Quantification in Full Waveform Tomography with Ensemble Data Assimilation – Methods and Anisotropy

Presenting Author: Julien Thurin, ISTerre, Université Grenoble Alpes, France

Authors: Julien Thurin, ISTerre, Université Grenoble Alpes, France; Romain Brossier, ISTerre, Université Grenoble Alpes, France; and Ludovic Métivier, LJK, CNRS, Université Grenoble Alpes, France

Uncertainty estimation and quality control are critically missing in most geophysical tomographic applications. The few solutions to cope with that issue are often left out in practical applications when these grow in scale and involve complex modeling. We present a joint full-waveform inversion and ensemble data assimilation scheme, allowing local Bayesian estimation of the solution that brings uncertainty estimation to the tomographic problem. This original methodology relies on a deterministic square root ensemble Kalman filter commonly used in the data assimilation community: the ensemble transform Kalman filter. Combined with a 2D visco-acoustic frequency-domain full-waveform inversion scheme, the resulting method allows access to a low-rank approximation of the posterior covariance matrix of the solution. It yields uncertainty information through an ensemble-representation, that can conveniently be mapped across the physical domain for visualization and interpretation. We present the combination of ensemble transform Kalman filter with full-waveform inversion along with the scheme design and algorithmic details that lead to our mixed application. Both synthetic and field-data results are considered, along with the biases that are associated with the limited rank ensemble representation. Finally, we review the open questions and development perspectives linked with data assimilation applications to the tomographic problem.

Cutting-edge Methods for Seismic Imaging II

10 November 2020 at 10 AM Pacific

Featured Discussants: Andrew Curtis, University of Edinburgh and Kees Wapenaar, Delft University of Technology

Abstracts

Characterization of a Magma Sill in the Mayotte Volcanic Zone Using Converted Waves

Presenting Author: Lise Retailleau, Université de Paris, Institut de Physique du Globe de Paris and Observatoire Volcanologique du Piton de la Fournaise, Institut de Physique du Globe de Paris

Authors: Lise Retailleau, Université de Paris, Institut de Physique du Globe de Paris, Observatoire Volcanologique du Piton de la Fournaise, Institut de Physique du Globe de Paris; Jean-Marie Saurel, Université de Paris, Institut de Physique du Globe de Paris; Jean Letort, Observatoire Midi Pyrénées, Université Paul Sabatier; Aude Lavayssière, Université de Paris, Institut de Physique du Globe de Paris; Claudio Satriano, Université de Paris, Institut de Physique du Globe de Paris; Wayne Crawford, Université de Paris, Institut de Physique du Globe de Paris; Jérome Vergne,Institut de Physique du Globe de Strasbourg, Université de Strasbourg;  Arnaud Lemarchand; Université de Paris, Institut de Physique du Globe de Paris

The magmatic event and subsequent creation of a new volcanic edifice that occurred offshore from the Mayotte island in the North Mozambique channel (Indian Ocean) is unique in its emitted lava volume and duration since the Laki eruption in 1784.

The event was first noticed through the strong seismicity it generated. Indeed, multiple Volcano Tectonic (VT) earthquakes of magnitude superior to 4 were felt on the island during the first month of the crisis (May 2018). The seismic activity is still ongoing and regularly felt by the population of the Mayotte island.

Very long period events (VLP) have also been recorded. The 11 November 2018 VLP was recorded worldwide and highlighted the volcanic origin of the activity before a scientific cruise discovered the volcano in May 2019. On active volcanoes, VLP signals are often associated with fluid movement or fluid resonance in a structure of the volcanic system. Evidence of a reservoir above the main seismicity swarm has been proposed by recent Magneto-telluric and petrologic studies.

The signature of this structure also appears as multiple arrivals on recordings of VT earthquakes located beneath it. The P and S-waves propagating upward and crossing the layer undergo conversions and speed deceleration in the sill, generating the supplementary arrivals.

Combining arrival time analysis and wave modeling in teleseismic and regional settings, we characterize a possible sill and identify its geometry and main properties. The layered structure has a thickness of about 2km, topping at 20km depth.

Imaging Light Structures at Earth’s Core-mantle Boundary – In Heavy Detail

Presenting Author: Paula Koelemeijer, Royal Holloway University of London

Authors: Paula Koelemeijer, Royal Holloway University of London; Alex Robson, University of California, Berkeley; Harriet Lau, University of California, Berkeley; Barbara Romanowicz, University of California, Berkeley & College de France; Arwen Deuss, Utrecht University; Jeroen Ritsema, University of Michigan

Accurate estimates of density variations are essential for modeling mantle flow and to distinguish between a thermal or compositional origin of mantle heterogeneity. Particularly, constraining the density of the Large-Low-Velocity-Provinces (LLVPs) in the deep mantle is important as the properties and stability of these structures have major implications for the long term evolution of the mantle. However, resolving LLVP density via seismological observations is not straightforward. Body wave data are only indirectly sensitive to density, while the sensitivity of Earth’s normal modes to density is relatively small and also shows trade-offs with topography on the core-mantle boundary. Consequently, few models of lower mantle density exist, with the most recent studies based on specific normal modes (Koelemeijer et al., 2017) and tidal data (Lau et al., 2017) drawing opposite conclusions.

Here, I will review efforts from the last decades to constrain the nature of the LLVPs, focusing on density and CMB topography. Specifically, I will discuss how the spatial extent of any dense material may be restricted to smaller regions within the LLVPs. In addition, I will discuss our current knowledge of the large-scale CMB topography, which also holds clues as to what the LLVPs may be. Furthermore, I will demonstrate that including CMB topography reconciles recent density studies based on normal modes and Earth’s tides, while the presence of a dense layer is also able to explain both data types. Together, these results point towards primarily thermal LLVPs with deep, spatially confined dense material present at their base.

Imaging Stratified Structures Using Spectral and Temporal Autocorrelations of Earthquake-based Data: Applications to Antarctic Grounded and Floating Ice Caps

Presenting Author: Thanh-Son Phạm, The Australian National University

Authors: Thanh-Son Phạm, The Australian National University; Hrvoje Tkalčić, The Australian National University

The use of passive seismic data to study structures of the ice covers in some parts of the Earth has been limited. This is either due to the lack of over-ice seismic deployment or the inefficiency of existing seismic methodology in the presence of shallow but sharp stratifications often encountered in the icy environment or sedimentary basin. In this talk, we present a new class of single-component only analysis, the autocorrelation method, to image stratified ice structures beneath individual receivers. Instead of using continuous ambient noise records, we advocate the use of P-wave coda – the waveform records immediately following the P arrival at steep incidence – from distant earthquakes to obtain high-resolution images of the ice structures. The method is then demonstrated using seismic data from four major over-ice arrays in Antarctica in the last two decades. The stacked autocorrelogram in time domain over earthquakes is deployed to obtain distinct reflection signals from the ice-rock interface beneath grounded ice stations. Meanwhile, at stations deployed over a floating ice shelf, where a significant water layer must be included into the problem, the stacked autocorrelogram in the spectral domain, which is equivalent to the temporal autocorrelation via a Fourier transform, proves to be effective to obtain thickness of both ice and water layers. Because the autocorrelation method presented here requires single-component recordings only, it has great potential in cryoseismological studies, or imaging applications at places with limited logistical access, such as in space missions to icy planets, polar regions or in ocean-bottom deployments.

Imaging Yellowstone’s Crustal Magmatic System Using Ambient Noise Correlations

Presenting Author: Ross Maguire, Michigan State University, University of New Mexico

Authors: Ross Maguire, Michigan State University, University of New Mexico; Brandon Schmandt, University of New Mexico; Min Chen, Michigan State University; Chengxin Jiang, Australian National University; Justin Wilgus, University of New Mexico

Seismic tomography has provided key constraints on the geometry and melt fraction of Yellowstone’s magmatic system and can provide insight into where Yellowstone is in its volcanic life cycle. Recent tomography studies depict seismic wavespeed reductions consistent with a melt-poor, non-mobile, crustal magma reservoir with melt fractions of roughly 10% or less. However, the seismically determined melt fraction may be underestimated due to imperfect data coverage, inversion regularization and inaccurate forward modeling. Additionally, evidence from teleseismic body wave modeling hints that the melt fraction could be much higher than 10%. Here, we use 3D spectral element waveform modeling to illustrate that traditional surface wave tomography can dramatically underestimate the amplitude of seismic anomalies associated with a melt-rich reservoir, which could lead to misinterpretation of the volume and mobility of melt present in Yellowstone’s crustal reservoir. To overcome these limitations, we are developing a new model based on adjoint ambient noise tomography, which incorporates more realistic wave propagation physics compared with traditional inversion techniques. We present a new Vs model of Yellowstone’s crust and uppermost mantle which will serve as our starting point for an adjoint inversion. The model uses the most complete set of ambient noise cross-correlations assembled for the Yellowstone region thus far, including over 2000 distinct source-receiver pairs with exceptional signal-to-noise ratio. Key features of the model, such as the amplitude and geometry of the low wavespeed anomaly below the Yellowstone caldera, agree well with past studies based on body and surface wave tomography. Adjoint inversion of our ambient noise cross-correlation dataset will likely sharpen these features which should help improve our understanding of the melt fraction and distribution in Yellowstone’s magmatic system.

Matrix Approach Applied to Passive Seismic Data

Presenting Author: Rita Touma, ISTerre, Université Grenoble Alpes, Institut Langevin, ESPCI Paris, PSL University, CNRS

Authors: Rita Touma, ISTerre, Université Grenoble Alpes, Institut Langevin, ESPCI Paris, PSL University, CNRS; Thibaud Blondel, Institut Langevin, ESPCI Paris, PSL University, CNRS; Arnaud Derode, Institut Langevin, ESPCI Paris, PSL University, CNRS; Michel Campillo, ISTerre, Université Grenoble Alpes; Alexandre Aubry, Institut Langevin, ESPCI Paris, PSL University, CNRS

To understand fault systems, it is required to identify the structure of the crust and upper mantle. Seismic investigations have long been relying on active sources generating an incident wave-field from the Earth surface. The reflected wave-field is then recorded by sensors deployed at the surface. Nowadays, passive imaging has been adopted as an alternative of this source-receiver configuration by computing the correlations of ambient noise. This process allows to estimate the Green’s function between two receivers. We here present a passive imaging technique applied to data recorded along the Clark Branch of the San Jacinto Fault Zone (Ben-Zion et al., 2015). Inspired by previous works in optics and acoustics, we introduce a matrix approach of seismic imaging based on seismic noise cross correlations. Our method applies focusing operations at emission and reception (Blondel et al., 2018) allowing to project the reflection matrix recorded at the surface to depth (redatuming). Although seismic noise is dominated by surface waves, focusing operations allow to extract the body wave components that carry information about the reflectivity of in-depth structures. However, complex velocity distribution of the Earth’s crust results in phase distortions, referred to as aberrations in the imaging process. Phase distortions prevent the imaging of the true reflectivity of the subsurface leading to unphysical features and blurry images. To overcome these issues, we introduce a new operator: the so-called distortion matrix. It connects any virtual source induced by focusing at emission with the distorted part of the reflected wave-front in the spatial Fourier domain. A time-reversal analysis of the distortion matrix is applied to correct for high-order aberrations. A 3D image of the fault is retrieved with optimal resolution and contrast. The damage volume extends from the surface to the first few hundred meters. Strata layers located at depths on both sides of the fault were also identified.

 

Cutting-edge Methods and Applications in Seismic Tomography

2 February 2021 at 10 a.m. –12 Noon Pacific

Abstracts

“From Raw Seismic Waveforms to Detailed Seismicity and Traveltime Tomography around the 2019 M7.1 Ridgecrest, CA Earthquake.”

Presenting Author: Malcolm White, University of Southern California.

Authors: Malcolm C. A. White, University of Southern California (malcolm.white@usc.edu); Hongjian Fang, Massachusetts Institute of Technology (hfang@mit.edu); Rufus D. Catchings, United States Geological Survey (catchings@usgs.gov); Yehuda Ben-Zion, University of Southern California (benzion@usc.edu)

We use an automated procedure to derive a catalog of nearly 95,000 event locations for the aftershock sequence following the 2019 M6.4 and M7.1 Ridgecrest, CA earthquake pair from raw seismic data recorded by rapid-response dense-deployment sensors following the mainshocks. We then iteratively relocate events and update our P- and S-wave velocity models (Vp and Vs, respectively) using fully 3D ray-based tomography with stochastic Voronoi-cell model parameterizations. The final Vp and Vs models correlate well with surface geology in the top 4 km of the crust and spatial seismicity patterns at depth. Joint interpretation of the derived catalog, velocity models, and surface geology suggests that (a) a compliant low-velocity zone near the Garlock Fault arrested the M7.1 rupture at the SE; (b) a stiff high-velocity zone beneath the Coso Mountains acted as a strong barrier that arrested the rupture at the NW; and (c) isolated seismicity on the Garlock Fault accommodated transtensional stepover strain triggered by the main events.

“Imaging 3D Anisotropic Subduction Zone Structure with Teleseismic P-waves: Method and Application to the Western United States.”

Presenting Author: Brandon VanderBeek, University of Padua

Authors: Brandon P VanderBeek, University of Padua (brandon.p.vanderbeek@gmail.com); Miles Bodmer, Sandia National Laboratories (mbodmer@gmail.com); Manuele Faccenda, University of Padua (manuele.faccenda@gmail.com)

Despite the well-established anisotropic nature of Earth’s upper mantle, the inuence of elastic anisotropy on teleseismic P-wave imaging remains largely ignored. Unmodeled anisotropic heterogeneity can lead to substantial isotropic velocity artefacts that may be misinterpreted as compositional and thermal heterogeneities. Here, we present a new parameterization for imaging arbitrarily oriented hexagonal anisotropy using teleseismic P-wave delays. We evaluate our tomography algorithm by reconstructing geodynamic simulations of subduction that include predictions for mantle mineral fabrics. Our synthetic tests demonstrate that accounting for both the dip and azimuth of anisotropy in the inversion is critical to the accurate recovery of both isotropic and anisotropic structure. We then perform anisotropic inversions using data collected across the western United States and offshore Cascadia. Our preliminary models show a clear circular pattern in the azimuth of anisotropy around the southern edge of the Juan de Fuca slab that is remarkably similar to the toroidal ow pattern inferred from SKS splits. We also image dipping anisotropic domains coincident with the descending Juan de Fuca slab. In contrast to prior isotropic tomographic results, the Juan de Fuca slab in our anisotropic model is characterized by more uniform P-wave speeds and is without an obvious slab hole below ~150 km depth. We also nd a general decrease in the magnitude of mantle low-velocity zones throughout the model relative to prior studies. These results highlight the sensitivity of teleseismic P-waves to anisotropic structure and the importance of accounting for anisotropic heterogeneity in the imaging of subduction zones.

“Bayesian Mapping of the Hawaiian ULVZ with Sdiff Postcursors.” 

Presenting Author: Carl Martin, University of Cambridge

Authors: Carl Martin, University of Cambridge (cm828@cam.ac.uk); Sanne Cottaar, University of Cambridge (sc845@cam.ac.uk); Thomas Bodin, ENS de Lyon (thomas.bodin@ens-lyon.fr)

We present a new inverse Bayesian approach to study extreme non-linear phenomena on the core-mantle boundary (CMB). Such features, known as ultra-low velocity zones (ULVZs), produce strong anomalies in S core-diffracted phases (Sdiff) which result in postcursor phases. As the ULVZs represent extremely reduced velocities (-10 to -30%) and the Sdiff postcursors travel times depend non-linearly on the ULVZ properties, conventional linearised tomography cannot be applied.

3D full waveform synthetic methods up to the frequencies required for the postcursors (0.1 Hz) are computationally expensive, making full-waveform inverse modelling intractable. So far, studies have forward modelled Sdiff postcursors assuming a simplified cylindrical anomaly (Cottaar & Romanowicz, 2012; Yuan & Romanowicz, 2017; Kim et al, 2020).

As a forward model, we employ 2D multi-arrival wavefront tracking simulations (Hauser et al, 2008) which can be computed on the order of seconds, making an inverse ULVZ problem feasible. This tool has so far been used for surface waves. To apply it to Sdiff waves on the core-mantle boundary, we project the events and stations onto the CMB and adjust the travel times for the up- and down-going legs.

Using reference arrivals times extracted from 3D full waveform synthetics, we demonstrate the viability of this technique with a level-set parameterisation to characterise the trade-offs between between size and shear velocity reduction, and the resolvability of the ULVZ morphology.

We then perform a Bayesian inversion on a real data set of Sdiff measurements sampling CMB structure beneath Hawaii. Two sets of crossing Sdiff paths are used, both using the dense station coverage provided from the Transportable Array: one geometry from the New Ireland Region to central US (2010), and the other from the Fiji Islands to Alaska (2018). The results show a quasi-cylindrical anomaly offset to the SW of the Hawaiian hotspot.

“Ps-P Crustal Tomography Enhances Imaging of Lower Crustal Magmatic Systems.” 

Presenting Author: Daniel Portner, Carnegie Institution for Science

Authors: Daniel E Portner, Carnegie Institution for Science (dportner@carnegiescience.edu); Lara S Wagner, Carnegie Institution for Science (lwagner@carnegiescience.edu); Helen A Janiszewski, University of Hawaii at Manoa (hajanisz@hawaii.edu); Diana C Roman, Carnegie Institution for Science (droman@carnegiescience.edu); John A Power, United States Geological Survey (jpower@usgs.gov)

Seismic tomography of the crust is an essential tool for characterizing the interconnected networks of magma storage regions and transport pathways that make up the magmatic plumbing systems feeding active volcanoes. Unfortunately, the resolvable depth range of crustal tomography is often limited by the distributions of local seismicity and/or the local seismic monitoring instruments. In particular, the deep crustal earthquakes required to image the lower crust with local body wave tomography are rare, while the broad aperture seismic networks required to resolve lower crustal velocities with surface wave inversions are commonly unavailable. Alternatively, teleseismic receiver functions can be used to illuminate local structural variations, but they do not typically account for the effects of three-dimensional velocity heterogeneities. Here we present Ps-P crustal tomography which harnesses the complementary strengths of both receiver function analysis and seismic tomography by processing Moho-generated Ps-P delay times derived from receiver functions in a tomographic S wave inversion. Using Ps-P tomography, we reveal the relative velocity structure of the entire crustal column beneath Mount Cleveland Volcano in the central Aleutian arc, Alaska, identifying a vertically extensive high Vp/Vs anomaly beneath the volcano that we attribute to a mid-to-lower crustal magma reservoir. This result illustrates the potential of Ps-P crustal tomography to advance imaging of the lower crust in magmatic systems without broad seismic networks or distributed local seismicity and lays the groundwork for more detailed crustal imaging at remote volcanoes.

“Dynamic Upwelling Beneath the Salton Trough Imaged with Teleseismic Attenuation Tomography.” 

Presenting Author: Joseph Byrnes, University of Minnesota/Northern Arizona University

Authors: Joseph S Byrnes, University of Minnesota, Northern Arizona University (byrnes.joseph@gmail.com); Max Bezada, University of Minnesota (mbezada@umn.edu)

The Salton Trough is one of the few regions on Earth where rifting is sub-aerial instead of sub-marine. We use the relative attenuation of teleseismic P phases recorded by the Salton Trough Seismic Imaging Project to investigate lithospheric and asthenospheric structures that form during extension. Map-view analysis reveals stronger attenuation within the Salton Trough than in the adjacent provinces. We then construct tomographic models for variations in seismic attenuation with depth to discriminate between crustal and mantle signals with a damped least-squares approach and a Bayesian approach. Synthetic tests show that models from damped least-squares significantly under-estimate the strength of attenuation and cannot separate crustal and mantle signals even when the tomographic models are allowed to be discontinuous at the lithosphere-asthenosphere boundary. We show that a Bayesian approach overcomes these problems when inverting the same synthetic datasets, and that shallow and deep signals are more clearly separated when imposing a discontinuity. With greater than 95% confidence, the results reveal first, that attenuation occurs primarily beneath the LAB; second, that the width of the attenuative region is narrower than the rift at 120 km depth; and third, that the strength of attenuation requires that the attenuative feature represents a melting-column similar to those beneath mid-ocean ridges. The narrow width of the melting-column below the volatile-free solidus is inconsistent with models for passive upwelling, where flow is driven only by rifting. Instead, we attribute the generation of incipient oceanic crust to mantle upwelling focused by buoyancy into a narrow diapir.

“Proximal Markov Chain Monte Carlo: Towards Building a Sparse Earth Model.” 

Presenting Author: Augustin Marignier, University College London

Authors: Augustin Marignier, University College London (augustin.marignier.14@ucl.ac.uk); Jason McEwen, University College London; Ana M. G. Ferreira, University College London; Thomas Kitching, University College London

Probabilistic sampling methods such as Markov chain Monte Carlo (MCMC) are a popular way to sample the posterior probability density function of an inverse problem. They are common in seismic tomography on local and regional scales, however they struggle on global scales where the dimensionality of the problem is typically much higher. The appeal of these methods for tomography is twofold: they allow for full uncertainty quantification and can solve non-linear inverse problems. Solutions to high-dimensional problems include the Hamiltonian Monte Carlo (HMC) method, which uses the gradient of the posterior to guide the sampling search. Unfortunately, this prohibits using non-smooth priors such as the Laplace distribution which promotes sparsity. Sparse image reconstructions have been shown to be able to recover both sharp and smooth features, even in underdetermined systems or when data are poorly distributed, the latter of which is a common scenario in seismic tomography.

In this work we use a recent proximal MCMC algorithm to build global fundamental mode Rayleigh wave phase velocity maps. Proximal MCMC leverages proximal calculus to allow non-differentiable priors in a high-dimensional inversion. As such, we adopt a sparsity-promoting prior, promoting sparsity in a spherical wavelet basis. We perform synthetic experiments and compare our results to those of a typical linear least-squares inversion. Finally, we discuss the future possibilities of building a 3D Earth model with uncertainties using a similar wavelet parameterised proximal MCMC approach.

 

 

Cutting-edge Methods and Applications in Seismic Tomography

2 March 2021 at 10 a.m. –12 Noon Pacific

Abstracts

“Getting to the Target with Marchenko-based Virtual Sources and Receivers.” Joeri Brackenhoff, Delft University of Technology (j.a.brackenhoff@tudelft.nl); C.P.A. Wapenaar, Delft University of Technology (c.p.a.wapenaar@tudelft.nl)

In many applications of seismic methods, a frequent issue is the presence of events that have scattered more than once in the subsurface, commonly referred to as internal multiples. The problem is that many methods cannot distinguish the internal multiples from the events that have scattered only once. This can cause artifacts in the final result of a method, especially for deep targets.

In recent years, the Marchenko method has been employed to deal with this kind of internal multiples. The method employs reflection data, which are data that are measured using active sources and receivers located at the surface of the Earth, to simulate sources and receivers in the subsurface. These simulated sources and receivers are referred to as virtual sources and receivers. In order to achieve this, the only information that is required is the first arriving event from the virtual source or receiver position in the subsurface to the Earth’s surface. Such information can be obtained from a background velocity model. The main advantages of the method are that a virtual source or receiver can be created at any point in the subsurface without the need to resolve the overburden above the target. Other methods can achieve similar results for primaries, however the Marchenko method can handle the internal multiples properly.

We will present applications of seismic methods using the Marchenko method, namely the creation of virtual sources and receivers for the purpose of imaging the subsurface and monitoring wavefields in the subsurface. We will first demonstrate this on numerical data, using a model that contains an overburden with strong scatterers above a target with weak scattering energy. We will show the difference between the results achieved with the Marchenko method and achieved by a method that does not handle the internal multiples. We will show the application of these methods to field data, which have been pre-processed in order to apply the Marchenko method.

“Constraints on the Geometry of the Subducted Gorda Plate From Converted Phases Generated by Local Earthquakes.” Jianhua Gong, MIT/WHOI Joint Program

The largest slip in great megathrust earthquakes often occurs in the 10–30 km depth range, yet seismic imaging of the material properties in this region has proven difficult. We utilize a dense onshore-offshore passive seismic dataset from the southernmost Cascadia subduction zone where seismicity in the mantle of the subducted Gorda Plate produces S-to-P and P-to-S conversions generated within a few km of the plate interface. These conversions typically occur in the 10–20 km depth range at either the top or bottom of a ~5 km thick layer with a high Vp/Vs that we infer to be primarily the subducted crust. We use their arrival times and amplitudes to infer the location of the top and bottom of the subducted crust as well as the velocity contrasts across these discontinuities. Comparing with both the Slab1.0 and the updated Slab2 interface models, the Slab2 model is generally consistent with the converted phases, while the Slab1.0 model is 1–2 km deeper in the 2–20 km depth range and ~6–8 km too deep in the 10–20 km depth range between 40.25°N and 40.4°N. Comparing the amplitudes of the converted phases to synthetics for simplified velocity structures, the amplitude of the converted phases requires models containing a ∼5 km thick zone with at least a ~10%–20% reduction in S wave velocity. Thus, the plate boundary is likely contained within or at the top of this low velocity zone, which potentially indicates a significant porosity and fluid content within the seismogenic zone.

“Fault Damage and Pore Fluid Distribution Regionally around Kaikoura, New Zealand from Body-wave Tomography.” Ben Heath University of Wisconsin – Madison (baheath@wisc.edu); Donna Eberhart-Phillips, University of California, Davis

Individual earthquake ruptures are usually assumed to occur on individual faults and are often associated with narrow regions (< 5 km) of altered physical properties, such as areas of increased fracturing and/or increased pore fluids. Recently, earthquakes such as the 2016 Kaikoura, New Zealand earthquake have ruptured multiple faults with different orientations over regions with widths spanning > 25 km. We test whether such regions hosting these earthquakes are associated with anomalous physical properties. We use seismic arrival-time tomography in the Kaikoura region to investigate lateral variations in Vp and Vp/Vs, using these parameters to infer variation in crustal faulting/fracturing. By modeling the effect of fluid-filled fractures on lateral variations in Vp and Vp/Vs, we are able to attribute the lateral variation in seismic velocities (over scales of > 50 km) to fault damage and pore fluid distribution. We find that the immature fault zones ruptured during the Kaikoura earthquake are on average characterized by decreased Vp and elevated Vp/Vs, features that decay (over distances of 50 km) towards background levels with increased distance from Kaikoura earthquake faults (and increased proximity to more mature fault zones). Drops in Vp in the Kaikoura rupture region are found to linearly relate to increases in Vp/Vs at a rate that is consistent with elevated 0.01 aspect ratio fractures, with highest fracturing within 10 km of the ruptured faults. The broad regional fracture distribution is likely the result of distributed long-term deformation, with increased deformation in the Kaikoura region. In contrast to more mature fault zones, which have localized strain accommodation and limited regional fracture distribution, immature fault zones are characterized by broadened, extensive fractures which contribute to complicated rupture dynamics.

“Level-set Imaging of the Los Angeles Basin using the Hierarchical Ensemble Kalman Sampling.” Jack Muir, California Institute of Technology (jmuir@caltech.edu); Robert W Clayton, Caltech (clay@gps.caltech.edu); Victor C Tsai, Brown University (victor_tsai@brown.edu)

A key challenge facing modern tomographic analyses is the interpretability of the resulting images, both in terms of the geological features obtained from imaging, and in how the results of uncertainty quantification modulate our confidence in the results. An innovative choice of model parametrization that improves both the interpretability and robustness of seismic tomography is therefore highly desirable. We have developed a framework employing level-set parametrizations that seeks to define geological units, which we invert using a recently developed Hierarchical Tikhonov Ensemble Kalman Sampling scheme, a highly efficient derivative-free optimizer which includes uncertainty estimates and optimization of regularization parameters. We have applied this scheme in an effort to better map earthquake hazard and improve velocity models within the Los Angeles Basin by generating a local update to the SCEC CVMS-4.26 model. We utilize high resolution data from the Community Seismic Network, a 400-station permanent urban deployment, to invert Love-wave dispersion, derived from eikonal tomography of two-station cross-correlation travel-time delays, and relative amplification data from the Mw 7.1 July 5 2019 Ridgecrest Earthquake. We find that the data is best explained by a deepening of the LA Basin (compared to the CVMS-4.26 reference model) along its Northwest-Southeast axis relative to its deepest point, just south of downtown LA. Additionally, the deeper basin edge extends further to the East of downtown LA towards East Los Angeles. This result offers new progress towards the parsimonious incorporation of detailed local basin models within regional reference models utilizing an objective inverse-problem framework with uncertainty quantification, and highlights the importance of accurate basin geometry models when accounting for the potentially significant amplification of surface waves from regional earthquakes in the high-rise building frequency band.

“Defining a New A-priori Velocity Model for the Seismic Tomography Analysis of South Italy.” Cristina Totaro, University of Messina (ctotaro@unime.it); Giancarlo Neri, University of Messina (geoforum@unime.it); Barbara Orecchio, University of Messina (orecchio@unime.it);  Debora Presti, University of Messina (dpresti@unime.it); Silvia Scolaro, University of Messina, (silscolaro@unime.it)

By integrating data and constraints available in the literature we estimated a new “a-priori” 3D seismic velocity model depicting the lithospheric structure of South Italy, a highly complex area of the Mediterranean region characterized by the coexistence of Africa-Europe NNW-trending plate convergence and SE-ward residual rollback of the Ionian lithospheric slab subducting underneath the Tyrrhenian. The integration of data like velocity values from seismic profiles and/or tomographies, moho depth, and subduction interface geometry, is performed with a procedure derived to that already successfully applied in the area about a decade ago and it is aimed to furnish the simplest 3D velocity structure consistent with all the available data. Taking benefit from studies and analyses of the last decade we enlarged and improved the previous estimated model by adding further data and useful constraints. The so obtained “a-priori” velocity model has then been employed as starting model for the local earthquake tomography of the region. For the tomographic study a set of ca. 10000 earthquakes has been used by selecting from the Italian seismic database (www.ingv.it) all the M≥2 earthquakes that occurred in the time period 2000-2020 at depth less than 60 km and with at least 10 station readings. The obtained 3D velocity structure together with the related hypocenter locations have been interpreted in the frame of the geodynamic models proposed for the region.

“An Indian Ocean Cluster of Mantle Plumes Imaged by Multifrequency P-wave Tomography.” Maria Tsekhmistrenko, DIAS (mariat@cp.dias.ie); Karin Sigloch, University of Oxford (karin.sigloch@earth.ox.ac.uk); Kasra Hosseini, The Alan Turing Institute (khosseini@turing.ac.uk); Guilhem Barruol, IPGP CNRS (barruol@ipgp.fr)

Mantle plumes are commonly envisioned as thin, buoyant conduits rising vertically from the core mantle boundary (CMB) to the earth’s surface, where they produce volcanic hot spots. Most hotspots are located in the sparsely instrumented oceans, creating poor prospects for the seismic resolution of thin conduits in the deep mantle.

The RHUM-RUM experiment remedied this issue around the hotspot island of La Réunion by instrumenting 2000km x 2000km of seaoor for 13 months with 57 broadband ocean bottom seismometers (OBS). We present a 3-D P-wave tomography model that was computed from the RHUM-RUM waveform data, supplemented by a global data set of P-diffracted measurements, and by a selection of ISC picks. Multifrequency traveltimes were measured on the waveforms and inverted in a finite-frequency framework. We achieve high image resolution beneath the Indian Ocean hemisphere, and especially beneath La Réunion, from upper mantle to CMB.
We observe the Large Low Velocity Province (LLVP) rising 800 km above the CMB, forming a cusp beneath South Africa. A low-velocity branch undulates obliquely from this cusp region towards the uppermost mantle beneath La Réunion. Hence La Réunion’s connection to the lower mantle is more complex than previously envisioned, being neither a thin vertical conduit nor projecting down to an edge of the LLVP. The deep-mantle connections of the Afar and Kerguelen hotspots emerge from the same LLVP cusp beneath South Africa and extend towards the surface through tilted low- velocity branches.

Our results provide the first high-resolution image of a western Indian Ocean plume cluster, from the surface to the CMB. This represents a key advance for linking geophysical, geodynamic and geochemical observations.”

 

Full-Waveform Methods and Applications

6 April, 10 a.m.–12:40 Pacific

Abstracts

“Adjoint Tomography of South America Based on 3D Spectral-Element Seismic Wave Simulations.” Caio Ciardelli, University of São Paulo

Advances in computational power during the last decade combined with the ever- growing seismographic coverage worldwide lead to better illumination of Earth’s interior and its dynamics. Full-waveform adjoint tomography has become a standard method in many seismic imaging studies, with the availability of computational facilities for large-scale simulations. We use 3D spectral-element continental-scale seismic wave simulations (Komatitsch & Tromp, 2002) and 112 earthquakes from the Harvard CMT catalog recorded by 1311 seismic stations to construct an adjoint tomography model of South America. We detect and remove noisy & problematic data using our multi-stage algorithm before the time-window selection, reducing the likelihood of discarding useful data or assimilating bad-quality waveforms in inversions. Our misfit function is a complex-exponentiated instantaneous phase, which optimizes the information extracted from each time series without the need for small-time windows, which is a requirement for cross-correlation measurements. Therefore we fine-tuned our window-selection algorithm to select long-time windows as much as the data quality permits. We use a preconditioner based on the pseudo-Hessian kernel and weight our misfit function to balance the source-receiver distribution for faster convergence. We have performed 23 iterations, gradually increasing the frequency content of the data to avoid local minima. Our final model (SAAM23) shows a significant decrease in the misfit. We further assessed the improvement by using cross-correlation measurements using 53 earthquakes that were not included in the adjoint inversion. In the long wavelengths, the model is compatible with previous studies, such as Feng et al. (2007), Celli et al. (2020), and Lei et al. (2020). SAAM23 agrees with many geological structures in the lithosphere, including the subduction along the Andes, the Amazonian craton, the São Francisco craton, and the Paranapanema block.

“A 3D Shear Wave Velocity Model of the Entire South American Subduction Margin.” Meng Liu, University of Massachusetts Amherst

The goal of this study is to image the seismic characteristics of the slab geometry along the entire South American subduction system, and to explore its spatial correlation with the distributions of volcanism and megathrust earthquakes. The collision between the Nazca plate and the South American continent has resulted in megathrust earthquakes and arc volcanoes, whose distribution patterns are significantly controlled by the geometry of the subducting plate. Previous studies revealed that subducting of the Nazca Ridge and the Juan Fernandez Ridge has contributed to the flat-slab segments in Peru and central Chile, respectively, where seismicity is absent. However, the lack of the lithospheric structures from trench to arc along the entire South American margin limits our understanding of the role of subduction dynamics on seismicity and volcanism. In this study, we construct a high-resolution shear velocity model of the crust and upper mantle within the entire South American subduction system by applying full-wave simulation and inversion method. We use 718 broadband seismic stations along the entire South American margin during 1994-2019. The vertical components of continuous seismic ambient noises between each station pair are cross-correlated and stacked in order to extract Rayleigh-wave empirical Green’s functions. The seismic velocity model is progressively improved by iteratively minimizing the phase delay measurements. Our tomographic results reveal four distinctive features, including (1) strong low-velocity anomalies (< 3.2 km/s) distributed extensively within the Andean continental crust, (2) two high-velocity flat slab segments located beneath southern Peru and central Chile, (3) complex slab geometry at flat-to-normal subduction transitions; and (4) low-velocity mantle wedge in correlation with surface volcanism. The scientific implications of our seismic observations in the South American subduction system will be further explored and discussed.

“Evolutionary Full-Waveform Inversion Applied to the African Plate.” Dirk-Philip van Herwaarden, ETH Zürich

We present a novel approach for the assimilation of new data in full-waveform inversion models. This evolutionary inversion is built on a form of stochastic L-BFGS that we refer to as dynamic mini-batch optimization. The motivation for this work is driven by the observation that seismic data volumes have been growing roughly exponentially. In theory, all this new data should enable us to improve our models. In practice, however, most Earth models are currently built from a static dataset, seeing little updates afterward.

In contrast to conventional waveform inversion, our approach uses dynamically changing mini-batches (subsets of data) that approximate the gradient of the larger dataset at each iteration. This has three major advantages, (1) We exploit redundancies within the dataset, which results in a reduced computational cost for model updates, (2) The size of the complete dataset does not directly impact the computational cost of an iteration, thereby enabling us to work with larger datasets, and (3) The nature of the algorithm makes it trivial to assimilate new data, as the new data can simply be added to the complete dataset from which the mini-batches are sampled.

To show the utility of this approach, we first apply it to an example where we invert for upper-mantle structure beneath the African Plate using data filtered down to 65 seconds. Starting from a 1-D model and with data recorded until 1995, we sequentially grow the dataset by incorporating more and more recent data into this ongoing inversion. Finally, we show our latest application, which started from the Collaborative Seismic Earth Model, where we incorporate periods down to 35 seconds.

“Extended Imaging: Convexification of Full-Waveform Inversion Problems and Target-Orient Elastic Inversion” Ettore Biondi, Stanford University

Full-waveform inversion (FWI) has proven its ability to invert high-resolution model parameters. However, its successful application is bounded by the knowledge of accurate initial model guesses. To avoid local minima and allow for inaccurate initial model estimates, we present a waveform inversion method that uses extended images to avoid the well-known issue of cycle skipping.

Moreover, in the context of seismic exploration, we demonstrate the ability of extended images to preserve the elastic amplitude response of pressure data and how they can be employed to perform a target-oriented elastic FWI method to obtain high-resolution elastic properties.

Finally, we apply the target-oriented methodology on an ocean bottom node field data 3D to characterize a potential target. We present the extended image space in the context of seismic exploration examples, but the theory can be easily applied to global seismology scenarios.

“Large-Scale Bayesian Full-Waveform Inversion using Hamiltonian Monte Carlo.” Lars Gebraad, ETH Zürich

Continental and global scale studies investigating mantle structure have in recent years been enhanced by advances in tomographic methods such as Full-Waveform Inversion (FWI), an approach that relies on the accurate simulation of wavefields using partial differential equations to reconstruct observed data. Methods such as these improve our ability to extract knowledge from data and thereby our capability of studying the inner Earth, at the cost of increased computational expense. Recent studies employing FWI acceleration techniques demonstrate that global scale deterministic FWI is within reach.

FWI problems can be recast as a Bayesian inference problems, using the data and prior knowledge to create a distribution over all possible models. Although this makes evaluating the solution of these problems computationally even more expensive, it also allows for powerful uncertainty quantification. Recent studies on the applicability of gradient-based Monte Carlo algorithms, such as Hamiltonian Monte Carlo (HMC), to seismological problems highlight their potential in this field.
The work presented shows how combining novel optimizations for both FWI and HMC create the powerful machinery needed to give Bayesian answers to tomographic questions on the scale of continents to the globe.

The product of combining these approaches yields efficient non-linear uncertainty quantification for FWI problems with favourable scaling properties w.r.t. many other available probabilistic algorithms. We’ll discuss the major ingedrients of such an approach, as well as preliminary considerations related to real data applications.

One major advantage of treating FWI problems probabilistically is that this allows one to distinguish between the effect of regularization (prior) and resolvability likelihood). This is important for studies involving parameters which are typically ill-recovered, such as density, attenuation and anisotropy.

“Downscaling Seismic Tomography.” Navid Hedjazian, Université de Lyon

In elastic full waveform inversion (FWI), scales smaller than a fraction of the minimum wavelength cannot be resolved, only a smoothed version of the true underlying medium can be recovered. Application of FWI to media containing small and strong heterogeneities therefore remains problematic. This smooth tomographic image is related to the effective elastic properties, which can be exposed with the homogenization theory of wave propagation. We study how this theory can be used in the FWI context. The seismic imaging problem is broken down in a two-stage multiscale approach. In the first step, called homogenized full waveform inversion (HFWI), observed waveforms are inverted for a macro-scale, fully anisotropic effective medium, smooth at the scale of the shortest wavelength present in the wavefield. The solution being an effective medium, it is difficult to directly interpret it. It requires a second step, called downscaling, where the macro-scale image is used as data, and the goal is to recover micro-scale parameters. All the information contained in the waveforms is extracted in the HFWI step. The solution of the downscaling step is highly non-unique as many fine-scale models may share the same long wavelength effective properties. We therefore rely on the introduction of external a priori information. In this step, the forward theory is the homogenization itself. It is computationally cheap, allowing to consider geological models with more complexity.

In a first approach to downscaling, the ensemble of potential fine-scale models is described with an object-based parametrization, and explored with a MCMC algorithm. We illustrate the method with a synthetic cavity detection problem. In a second approach, the prior information is introduced by the means of a training image, and the fine-scale model is recovered with a multi-point statistics algorithm. We apply this method on a subsurface synthetic problem, where the goal is to recover geological facies.