We have inverted a very large combination of delay times for P, pP, PcP, PKPab, PKPbc from the ISC-EHB catalog, augmented with hand-picked OBS and MERMAID onset data and cross-correlation delays from broadband stations to derive a new global model, UNICA25, of P-wave velocity variations in the mantle. This paper presents the model, summarizes the data processing and inversion, and shows resolution tests at different scales. After rejecting 109,917 outliers (1% of the total), the data set comprises 10,571,152 arrival times. The misfit of the observed delays, with respect to a background model consisting of the crust of LITHO1.0 on top of model AK135, is on average 2.14 standard errors (σ), which is reduced to 0.99σ by inverting for P-wave velocity variations using both ray theory and finite frequency theory, and 5000 iterations of the iterative solver LSQR. Extensive resolution tests show a resolution at 500 km scale or better in most of the lower mantle. An upper mantle resolution of the order of 300 km is obtained under densely instrumented continents, but upper mantle resolution remains deficient under much of the oceans except where covered by MERMAIDs. Linearized tomography is at least three orders of magnitude faster than recent full waveform inversions (FWI), and reaches much higher frequencies for cross-correlated delays. Model UNICA25 has a P-velocity resolution that rivals that of FWI and shows a superior performance in fitting delay times, notably the ISC-EHB P-delays.

We present a statistically and computationally efficient spectral-domain maximum-likelihood procedure to solve for the structure of Gaussian spatial random fields within the Mat & eacute;rn covariance hyperclass. For univariate, stationary and isotropic fields, the three controlling parameters are the process variance, smoothness and range. The debiased Whittle likelihood maximization explicitly treats discretization and edge effects for finite sampled regions in parameter estimation and uncertainty quantification. As even the 'best' parameter estimate may not be 'good enough', we provide a test for whether the model specification itself warrants rejection. Our results are practical and relevant for the study of a variety of geophysical fields, and for spatial interpolation, out-of-sample extension, kriging, machine learning and feature detection of geological data. We present procedural details and high-level results on real-world examples.

Elastic full-waveform inversion (EFWI) is a state-of-the-art seismic tomographic method. Recent advances in technology and instrumentation, combining crosstalk-free source-encoded FWI (SE-FWI) with multicomponent marine data acquisition using ocean-bottom nodes (OBNs), enable full-physics wave propagation and parameter inversion without the computational burden of traditional FWI. With OBN acquisition, P waves, S waves, and P-to-S conversions are recorded. It is not well understood to what extent adding horizontal components to SE-FWI improves the resolution of subsurface modeling. We assess their potential for the reconstruction of shear and compressional wave speeds (V P and V S ) by using a synthetic data set modeled after a recently acquired OBN survey in the North Sea. We perform synthetic inversion tests to design suitable strategies that leverage the information recorded in the horizontal components of the data to improve the reconstructed model resolution laterally and in depth. We advocate for a hierarchical inversion approach to recover the elastic parameters. We exploit the P and P-to-S converted waves recorded on the horizontal components to robustly reconstruct both V P and V S . Adding horizontal components to the SE-FWI modeling workflow results in improved spatial resolution, enhanced depth coverage, and more accurate elastic wave speed estimates.

We present the analysis of more than 9000 hydro-acoustic earthquake records, recorded by a network of 50 instruments called Mobile Earthquake Recording in Marine Areas by Independent Divers (MERMAIDs), which are freely floating in the South Pacific Ocean. This network is part of the collaborative South Pacific Plume Imaging and Modeling (SPPIM) project. Our analysis focuses on evaluating how the MERMAID stations complement the conventional worldwide network of stations reporting to the ISC, comprising records of globally distributed seismic stations. In the context of routine earthquake (re)location, we evaluate the improvement of results of earthquake location estimates, particularly focusing on the Tonga-Kermadec subduction zone, where current distribution of seismic land stations is extremely sparse. MERMAIDs are often the closest “station” to earthquakes occurring in the Tonga subduction zone, and frequently fill significant station azimuthal gaps. We matched the MERMAID records from June 2018 to December 2023 to more than 3000 earthquakes which are reported to the ISC by many international agencies worldwide, allowing us to construct the MERMAID catalogue, complimented with additional parametric data such as MERMAIDs locations, observed direct P-phase picks and travel time residuals with respect to ak135 model, as well as their corrections for bathymetry. These data were integrated into the existing ISC Bulletin, and relocated using the ISC hypocentre relocation algorithm ISCLoc. The results of this successful first integration of MERMAID data into ISC routines presents an important step towards routinely including MERMAID traveltime picks in ISC Bulletin.

Seismic tomography is a principal method for studying deep mantle plume structure. Imaging Earth's wavespeed anomalies is conditioned byseismic wave sampling, and the uneven distribution of receiving stations worldwide leaves several candidate plumes beneath varioushotspots across the globe poorly resolved. We regionally evaluate two full-waveform global tomography wavespeed models, GLAD-M25 andSEMUCB-WM1, focusing on the mantle below the Pacific Ocean in the region of the South Pacific Superswell. This area contains multiple hotspots which may be anchored in the Large Low Shear-Velocity Province at the base of the mantle. The two models show similarities and differences in the target region. With a goal of guiding subsequent iterations in the GLAD model suite, we assess the quality of GLAD-M25 in the target region relative to its global performance using a regional partition of the seismic waveform data used in its construction. We evaluate synthetic waveforms calculated using the spectral-element method, based on how well they fit the data according to a variety of criteria measured across multiple seismic phases, wavetypes, and frequency bands. The distributions of travel-time anomalies that remain in GLAD-M25 are wider regionally than globally, suggestingcomparatively insufficiently resolved seismic velocity structure in the region of interest. This will motivate regionally focused inversions based on a subset of the global data set, and the additionof data sampling new corridors, especially using ocean sensors. We compare GLAD-M25 and SEMUCB-WM1 by cross-validation with a new,independent, data set. Our results reveal that short- and long-wavelength structure is captured differently by the two models. Global models use misfit criteria that may strive for balance between portions of the data set, but could leaveimportant regional domains underserved. Our results lead us to recommend focusing future model iteration and data addition on andaround the Pacific Superswell to better constrain seismic velocity structure in this area of significant geodynamic complexity.

Floating seismometers (‘MERMAIDs’) operating in the noisy environment of the world’s oceans pose a challenge for picking the time of earthquake first arrivals. We report on an experiment to estimate the errors in picked arrivals from 49 MERMAIDS operating in the South Pacific, using two independent strategies. For 15 events, the same arrivals were redundandly picked by several analysts, allowing for a direct estimate of error distributions. Standard errors in times from MERMAID seismograms vary from 0.2 s for close events at mantle depths in the Kermadec subduction zone to more than 2 s for crustal events at large epicentral distance. In a second experiment we analysed the a posteriori misfits after tomographically inverting all events. The residual traveltime misfit is consistent with the error estimates from the first experiment, but also shows inconsistencies with arrival times from the ISC-EHB and NEIC catalogues, which we attribute to errors in the published hypocentres and/or origin times.

Crosstalk-free source-encoded elastic full-waveform inversion (FWI) using time-domain solvers demonstrates skill and efficiency at conducting seismic inversions involving multiple sources and receivers with limited computational resources. A drawback of common formulations of the procedure is that, by sweeping through the frequency domain randomly at a rate of one or a few sparsely sampled frequencies per shot, it is difficult to simultaneously incorporate time-selective data windows, as necessary for the targeting of arrivals or wave packets during the various stages of the inversion. Here, we solve this problem by using the Laplace transform of the data. Using complex-valued frequencies allows for damping the records with flexible decay rates and temporal offsets that target specific traveltimes. We present the theory of crosstalk-free source-encoded FWI in the Laplace domain, develop the details of its implementation, and illustrate the procedure with numerical examples relevant to exploration-scale scenarios.

We present a computational technique to model hydroacoustic waveforms from teleseismic earthquakes recorded by mid-column MERMAID floats deployed in the Pacific, taking into consideration bathymetric effects that modify seismo-acoustic conversions at the ocean bottom and acoustic wave propagation in the ocean layer, including reverberations. Our approach couples axisymmetric spectral-element simulations performed for moment-tensor earthquakes in a 1-D solid Earth to a 2-D Cartesian fluid-solid coupled spectral-element simulation that captures the conversion from displacement to acoustic pressure at an ocean-bottom interface with accurate bathymetry. We applied our workflow to 1129 seismograms for 682 earthquakes from 16 MERMAIDs (short for Mobile Earthquake Recording in Marine Areas by Independent Divers) owned by Princeton University that were deployed in the Southern Pacific as part of the South Pacific Plume Imaging and Modeling (SPPIM) project. We compare the modelled synthetic waveforms to the observed records in individually selected frequency bands aimed at reducing local noise levels while maximizing earthquake-generated signal content. The modelled waveforms match the observations very well, with a median correlation coefficient of 0.72, and some as high as 0.95. We compare our correlation-based traveltime measurements to measurements made on the same data set determined by automated arrival-time picking and ray- traced traveltime predictions, with the aim of opening up the use of MERMAID records for global seismic tomography via full-waveform inversion.

Seismic discontinuities in the mantle are indicators of its thermo-chemical state and offer clues to its dynamics. Ray-based seismic methods, though limited by the approximations made, have mapped mantle transition zone discontinuities in detail, but have yet to offer definitive conclusions on the presence and nature ofmid-mantle discontinuities. Here, we show how to use a wave-equation-based imaging method, reverse-time migration of precursors to surface-reflected seismic body waves, to uncover both mantle transition zone and mid-mantle discontinuities, and interpret their physical nature. We observe a thinnedmantle transition zone southeast of Hawaii, and a reduction in impedance contrast around 410 km depth in the same area, suggesting a hotter-than-average mantle in the region. Here, we furthermore reveal a 4000-5000 km-wide reflector in new images of the mid mantle below the central Pacific, at 950-1050 km depth. This deep discontinuity exhibits strong topography and generates reflections with polarity opposite to those originating at the 660 km discontinuity, implying an impedance reversal near 1000 km. We link this mid-mantle discontinuity to the upper reaches of deflected mantle plumes upwelling in the region. Reverse-time migration fullwaveform imaging is a powerful approach to imaging Earth's interior, capable of broadening our understanding of its structure and dynamics and shrinking modeling uncertainties.

The highly used Global Seismographic Network (GSN) is a pillar of the seismological research community and contributes to numerous groundbreaking publications. Despite its wide recognition, this survey found that the GSN is not consistently acknowledged in scientific literature and is underrepresented by roughly a factor of 3 in citation searches. Publication tracking is a key metric that factors into operational decisions and funding support for the network; thus, consistent and proper citation of the GSN is important. This study not only serves as a reminder for researchers using GSN observations to cite the network's digital object identifiers (DOIs) but also promotes a community-wide conversation among researchers, journal editors, network operators, and other stakeholders regarding more standardized policies and review processes to ensure seismic networks are properly and consistently recognized for their contributions to research.

Since the launch of the first third-generation MERMAID in 2018, sixty-seven autonomous freely-drifting mid-column hydrophones have been afloat in the Pacific Ocean, the South China Sea, and the Mediterranean. Over fifty-five instruments remain alive and well, and are continuing to report short, triggered, waveform segments of acoustic pressure variations, with six instruments directly reporting acoustic spectral densities. A new model equipped with a conductivity-temperature-depth (CTD) sensor is due for deployment (four new units). Over these last few years, many thousands of teleseismic arrival-times have been reported, and associated with global earthquake catalogs, so that their travel-time residuals with respect to global reference models can be determined in view of making tomographic models, especially of the area around French Polynesia, the original launch focus of the EarthScope-Oceans (ESO) fleet. The data stream has been flowing into the EarthScope (IRIS) data management center (DMC), and reporting short segments around first-arriving phases and the calculation of travel-times have become routine applications.We briefly review those successes, but we focus on the latest data types (spectral densities computed in-situ rather than time-domain seismograms transmitted via satellite), on the codesign of MERMAID as an acoustic float with an environmental CTD sensor, on the latest modeling efforts at matching waveforms with synthetically computed pressure seismograms, and on the rich set of continuous records that were requested from MERMAID's one-year buffer in the hours and days immediately following the Hunga Tonga Hunga Ha'apai eruption, recorded by over twenty instruments over a wide epicentral distance and backazimuthal range.Waveform modeling was not a design goal for MERMAID, but we discuss an innovative approach to circumvent the computational burden involved in matching the global earth ocean-bottom response to teleseismic earthquakes (computed using an axisymmetric spectral-element code) to the local oceanic mid-column pressure response, including the effects of bathymetry (computed using a local two-dimensional spectral-element modeling step). We applied our method, which is based on precomputed Green's functions via the publicly available code Instaseis and on a custom data base of ocean-floor-to-water-column response functions computed using SPECFEM-2D, to a set of over one thousand waveforms, after dynamically selecting the optimal bandwidth based on adaptive signal-to-noise considerations to steer clear of the noise generated by the ocean wave heave. The correlation between synthetics and observations is as high as 0.98, with a median of 0.72, and very coherent across the array, allowing for the determination of cross-correlation travel times and opening up MERMAID seismograms to conduct full-waveform tomography of Earth's mantle.Similarly, volcanic monitoring was not a design goal for MERMAID, but in recovering, on-demand, the many hours of continuous records of the Hunga Tonga Hunga-Ha`apai eruption, we have obtained a unique data set from which are piecing together a detailed picture of the eruption sequence. We discuss signal correlations and disparities within and across the South Pacific Plume Imaging and Modeling (SPPIM) array, and address, in particular, path-dependent effects due to "bathymetric occlusion", the influence of seafloor topography on the coherent propagation of hydroacoustic energy over large distance ranges.

The eastern continental margin of North America, despite being a passive margin at present, records a comprehensive tectonic history of both mountain building and rifting events. This record is punctuated by several igneous events, including those associated with the Great Meteor and Bermuda hotspots. To gain a better understanding of the state of the mantle beneath this region, we employ the massive quantity of seismic data recorded by the USArray to image the mantle transition zone beneath eastern North America. To construct these images, we first calculate P-to-s receiver functions using an iterative time-domain deconvolution algorithm. These receiver functions are then automatically filtered by their quality, using a set of rigorous criteria, and subsequently summed using common conversion point stacking. We present several cross sections through these stacks, which show remarkable features such as a thinned transition zone beneath the independently observed northern Appalachian and central Appalachian low-wavespeed anomalies, as well as a thickened transition zone beneath western Tennessee associated with the Laramide slab stagnating at depth. In addition to discussing these geologically relevant features, we perform a technical analysis of the effects of using various seismic velocity models for the moveout correction of our receiver functions. We find that the thickness of the mantle transition zone under eastern North America is a robust measurement, while the resolved depths of the 410 and 660 km discontinuities are model dependent.

<p><br>MERMAIDs (Mobile Earthquake Recording in Marine Areas by Independent Divers) are seismic instruments which record local, regional and teleseismic earthquakes, and other signals, in the oceans. In the Southern Pacific Ocean, some fifty MERMAIDs are collecting acoustic pressure time series as part of the South Pacific Plume Imaging and Modeling (SPPIM) project. Deployed in 2018 and 2019 by members of the EarthScope-Oceans consortium, these instruments record continuous time series data on a one-year buffer and autonomously report a wealth of waveforms, selectively triggered mostly by teleseismic events, suitable for mantle tomography. </p><p>Listening for signals while roughly 1.5 km below the ocean surface, MERMAIDs' primary mission is to detect and deliver records of P-waves generated by distant earthquakes, and collectively they have returned many thousands of seismograms corresponding to such signals. We present highlights from our earthquake catalog and discuss the changing character and causes of the background noise. Whilst the South Pacific fleet is programmed to only send short seismic records, corresponding to confident identifications of teleseismic first arrivals, some records contain later-arriving phases. In addition to <em>P</em>-waves, we present a miscellany of observations of other signals, including core phases, converted <em>S</em>- and surface waves, and <em>T</em>-phases. Furthermore, we illustrate that we are able to obtain other recorded data segments through buffer requests via satellite. </p><p>Data from the MERMAIDs owned by Geoazur and Princeton University that we report on here are being archived by IRIS—those from three instruments is available without embargo. We highlight MERMAID waveform availability and its utility to the scientific community via examples of their modeling and preliminary interpretations that can be made regarding wavespeed heterogeneity in the dynamic mantle below the Pacific Ocean. </p>

We establish a theoretical framework, an algorithmic basis, and a computational workflow for the statistical analysis of multi-variate multi-dimensional random fields - sampled (possibly irregularly, with missing data) and finite (possibly bounded irregularly). Our research is practically motivated by geodetic and scientific problems of topography and gravity analysis in geophysics and planetary physics, but our solutions fulfill the more general need for sophisticated methods of inference that can be applied to massive remote-sensing data sets, and as such, our mathematical, statistical, and computational solutions transcend any particular application. The generic problem that we are addressing is: two (or more) spatial fields are observed, e.g., by passive or active sensing, and we desire a parsimonious statistical description of them, individually and in their relation to one another. We consider the fields to be realizations of a random process, parameterized as a Matern covariance structure, a very flexible description that includes, as special cases, many of the known models in popular use (e.g. exponential, autoregressive, von Karman, Gaussian, Whittle, ...) Our fundamental question is how to find estimates of the parameters of a Matern process, and the distribution of those estimates for uncertainty quantification. Our answer is, fundamentally: via maximum-likelihood estimation. We now provide a computationally and statistically efficient method for estimating the parameters of a stochastic covariance model observed on a regular spatial grid in any number of dimensions. Our proposed method, which we call the Debiased Spatial Whittle likelihood, makes important corrections to the well-known Whittle likelihood to account for large sources of bias caused by boundary effects and aliasing. We generalise the approach to flexibly allow for significant volumes of missing data including those with lower-dimensional substructure, and for irregular sampling boundaries. We build a theoretical framework under relatively weak assumptions which ensures consistency and asymptotic normality in numerous practical settings including missing data and non-Gaussian processes. We also extend our consistency results to multivariate processes. We provide detailed implementation guidelines which ensure the estimation procedure can still be conducted in O(n log n) operations, where n is the number of points of the encapsulating rectangular grid, thus keeping the computational scalability of Fourier and Whittle-based methods for large data sets. We validate our procedure over a range of simulated and real world settings, and compare with state-of-the-art alternatives, demonstrating the enduring practical appeal of Fourier-based methods, provided they are corrected and augmented by the procedures that we developed.

We provide a computationally and statistically efficient method for estimating the parameters of a stochastic covariance model observed on a regular spatial grid in any number of dimensions. Our proposed method, which we call the Debiased Spatial Whittle likelihood, makes important corrections to the well-known Whittle likelihood to account for large sources of bias caused by boundary effects and aliasing. We generalize the approach to flexibly allow for significant volumes of missing data including those with lower-dimensional substructure, and for irregular sampling boundaries. We build a theoretical framework under relatively weak assumptions which ensures consistency and asymptotic normality in numerous practical settings including missing data and non-Gaussian processes. We also extend our consistency results to multivariate processes. We provide detailed implementation guidelines which ensure the estimation procedure can be conducted in O(nlogn) operations, where n is the number of points of the encapsulating rectangular grid, thus keeping the computational scalability of Fourier and Whittle-based methods for large data sets. We validate our procedure over a range of simulated and realworld settings, and compare with state-of-the-art alternatives, demonstrating the enduring practical appeal of Fourier-based methods, provided they are corrected by the procedures developed in this paper.

Mid-plate volcanoes are well known as hotspots. They represent the surface signature of mantle plumes, nevertheless their origin and their role in geodynamics are still a challenge in the Earth sciences. Even though plate tectonics and mantle plumes were discovered at the same time, the latter cannot be explained by the former. Plumes’ birth, life and death play a fundamental role on the evolution of life on Earth and on plate-tectonic reorganization. La Réunion hotspot is known as one of the largest on the Earth, that created the Deccan volcanic traps in India (almost 2 million km2) and the death of more than 90% of life on the Earth including dinosaurs ~65Ma ago. So far the origin of the mantle plumes and their role in geodynamics are still unclear in Earth sciences. In that respect, we use the dataset from the French-German RHUM-RUM experiment around La Réunion hotspot (2012-2013), from IRIS data center and FDSN to extensively investigate the deep structure of the plume along its complete track from its birth to its present stage, as well as from the upper mantle to the lowermost mantle. Several shear-wave anomalies are resolved underneath Indian Ocean and the upper mantle beneath this region is fed by mantle plume rising from the core-mantle boundary. The lower mantle thermochemical dome associated to the South-African Large Low-Shear Velocity Province (LLSVP) is found to be composed of several conduits. Plume branches are highlighted at ~900 km depth. Thermal instability and thermochemical heterogeneities in the D" layer are likely the principal reasons of the plumes birth at the core-mantle boundary, and therefore an indicator of long-life of the Réunion hotspot.

SUMMARY Passive seismic inversion at the reservoir scale offers the advantages of low cost, negligible environmental impact and the ability to probe a target area with low-frequency energy not afforded by even the most modern active-source seismic technology. In order to build starting models suitable for full-waveform wave speed tomography, characterization of earthquake sources is an indispensable first step. We present a workflow for the centroid moment tensor (CMT) inversion of seismic events identified in a passive seismic data set acquired by a large and dense array of three-component broad-band seismic sensors in a mountainous setting in the Himalayan foothills. The data set comprised 256 instruments operating for 2×4 months over an area of 8000 km2. An initial 3-D wave speed model was determined for the region via the analysis of first-arriving traveltime picks. Of the 2607 identified seismic events that were well recorded at frequencies between 0.2–50 Hz, 86 with magnitudes 1.3 ≤ M ≤ 3.0 initially had their CMT focal mechanisms determined by a waveform fitting procedure built on a Green’s function approach in a 1-D layered average wave speed model, for stations within an offset of 10 km, in the frequency range 0.2–1.4 Hz. Here, we obtain updated CMT mechanisms for the 86 events in that catalogue via multicomponent full-waveform inversion in the 3-D wave speed model. Our workflow includes automated data- and model-driven data selection using a combination of different metrics derived from signal-to-noise considerations and waveform-fitting criteria, and relies upon spectral-element simulations of elastic wave propagation in the 3-D wave speed model, honouring topography. Starting from the initial CMT solutions, we seek improvement to the data fit within the frequency band 0.5–2.5 Hz by minimizing the waveform difference between observed and synthetic data, while accommodating wave speed-model errors by allowing for small time-shifts. We balance uneven data coverage and tune their contributions via data-space weighting functions. We quantify the improvements to the data fit in terms of different metrics. We summarize the changes to the CMT solutions, and present and analyse the resulting catalogue for the region, including their breakdown into double-couple and non-double couple components, and their relation to mapped faults.

We present the first 16 months of data returned from a mobile array of 16 freely floating diving instruments, named MERMAID for Mobile Earthquake Recording in Marine Areas by Independent Divers, launched in French Polynesia in late 2018. Our 16 are a subset of the 50 MERMAID deployed over a number of cruises in this vast and understudied oceanic province as part of the collaborative South Pacific Plume Imaging and Modeling (SPPIM) project, under the aegis of the international EarthScope-Oceans consortium. Our objective is the hydroacoustic recording, from within the oceanic water column, of the seismic wavefield generated by earthquakes worldwide, and the nearly real-time transmission by satellite of these data, collected above and in the periphery of the South Pacific Superswell. This region, characterized by anomalously elevated oceanic crust and myriad seamounts, is believed to be the surface expression of deeply rooted mantle upwellings. Tomographically imaging Earth’s mantle under the South Pacific with data from these novel instruments requires a careful examination of the earthquake-to-MERMAID traveltimes of the high-frequency P-wave detections within the windows selected for reporting by the discrimination algorithms on board. We discuss a workflow suitable for a fast-growing mobile sensor database to pick the relevant arrivals, match them to known earthquakes in global earthquake catalogues, calculate their traveltime residuals with respect to global seismic reference models, characterize their quality and estimate their uncertainty. We detail seismicity rates as recorded by MERMAID over 16 months, quantify the completeness of our catalogue and discuss magnitude–distance relations of detectability for our network. The projected lifespan of an individual MERMAID is 5 yr, allowing us to estimate the final size of the data set that will be available for future study. To prove their utility for seismic tomography we compare MERMAID data quality against ‘traditional’ land seismometers and their low-cost Raspberry Shake counterparts, using waveforms recovered from instrumented island stations in the geographic neighbourhood of our floats. Finally, we provide the first analyses of traveltime anomalies for the new ray paths sampling the mantle under the South Pacific.