Imaging, modelling and inversion to explore the Earth’s lithosphere

This session will cover applied and theoretical aspects of geophysical imaging, modelling and inversion using active- and passive-source seismic measurements as well as other geophysical techniques (e.g., gravity, magnetic and electromagnetic) to investigate properties of the Earth’s crust and uppermost mantle, and explore the processes involved. We invite contributions focused on methodological developments, theoretical aspects, and applications. Studies across the scales and disciplines are particularly welcome.

Among others, the session may cover the following topics:
- Active- and passive-source imaging using body- and surface-waves;
- Full waveform inversion developments and applications;
- Advancements and case studies in 2D and 3D imaging;
- Interferometry and Marchenko imaging;
- Seismic attenuation and anisotropy;
- Developments and applications of multi-scale and multi-parameter inversion; and,
- Joint inversion of seismic and complementary geophysical data.

Co-organized by GD3
Convener: Matthew AgiusECSECS | Co-conveners: Giovanni DiaferiaECSECS, Michal Malinowski, Monika Ivandic, Cedric Schmelzbach
vPICO presentations
| Mon, 26 Apr, 09:00–12:30 (CEST)

Session assets

Session materials

vPICO presentations: Mon, 26 Apr

Chairpersons: Monika Ivandic, Matthew Agius, Pascal Edme
Ivan Vasconcelos, Wouter Klessens, Yang Jiao, Andre Niemeijer, and Suzanne Hangx

More often than not, important geologic processes occur at micro-scales, e.g., fluid flow, mineral-phase changes, chemically-induced alteration, rock-frame compaction, or even mechanical ruptures/instabilities leading to large earthquakes. However, reliably imaging material properties at such scales from remote long-wavelength information contained in either seismic or EM fields has long been a challenge to the geophysical, engineering and material science communities. In this talk, we present a general framework for the estimation of sub-wavelength material properties from long-scale waves, building on recent advances on statistical microstructure descriptors (SMDs) within the field of material science.


In geoscience, traditional approaches to describing material microheterogeneity rely on either analytical inclusion-based models, or in sample-based digital rocks: each of these having their pros and cons. Here, we instead rely on SMDs, namely two-point correlation and polytope functions, to describe microheterogeneous geo-materials in a manner that is capable of generalizing complex geometrical information hidden in microstructures, while also retaining realism and sample fidelity. Using SMDs, we rely on wave-equation-based Strong Contrast Expansions (SCEs) to predict frequency/scale-dependent effective wave properties for acoustic, elastic and EM waves.  We briefly discuss how SMD-described microstructures affect long-wave properties – and in particular how they not only predict frequency-dependent attenuation due to sub-wavelength scattering, but that attenuation is particularly sensitive to microstructure when compared to effective wavespeeds.


When it comes to the estimation of microstructure properties from wave observations, the problem becomes substantially more difficult because realistic microscale parameters could in principle have far too many degrees of freedom than what is observable from finite-frequency wave data. As such, it is key that any method that aims at realistically retrieving microstructure information from long-scale wave data accounts for uncertainty, while also handling the highly nonlinear nature of microstructure-dependent effective wave properties. To that end, we combine our SMD and SCE approaches for effective wave properties with the supervised machine-learning method of Random Forests to construct a Bayesian approach to infer microstructure properties from effective wave parameters as observables. This method yields full posterior distributions for microstructure parameters (e.g., property contrast, volume fraction, and geometry information) from frequency-dependent observations of wave velocities and attenuation. We present several examples of inference scenarios, showing, for example, that i) attenuation is key to microstructure imaging, and ii) microgeometry information can only be reliably retrieved if contrast and volume fraction are relatively well known a priori. We illustrate of inference approach with several examples of analytical and real microstructures, including data from a  laboratory compaction experiment controlled by microscale CT imaging.

How to cite: Vasconcelos, I., Klessens, W., Jiao, Y., Niemeijer, A., and Hangx, S.: Imaging the Earth’s small structures: AI-driven, Bayesian inference of microstructure descriptors from finite-frequency waves, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-10058,, 2021.

Active sources
Richard Hobbs and Christine Peirce

The transition zone between the more porous upper extrusive layer (2A) and the less porous lower dyke layer (2B) within the oceanic crust is characterised by a high velocity gradient based on inversion of controlled-source, long-offset refraction data. In these data the phase associated with this high velocity gradient, termed the 2A Event, has an anomalously high amplitude over a limited range of offsets and appears to form a triplication with refractions from layer 2A above and 2B below. These characteristics fit the accepted model that this event is a caustic or retrograde phase, generated by a distinct layer whose thickness and velocity gradient can be determined by ray-trace modelling. Hence, a velocity model for Layer 2 (derived from seismic data acquired near ODP 504B) consists of a ~500 m-thick 2A with a velocity gradient of ~1.0 s-1; a ~200 m-thick transition zone with a high velocity gradient of ~4.0 s-1; and a ~1300 m-thick 2B with a velocity gradient of ~0.3 s-1. However, this model is at odds with observation of Layer 2 lithology obtained from coring and ophiolites where the 2A is composed of a mixture of higher velocity basalt flows and lower velocity pillow lavas and breccia, with the transition zone represented by an increasing number of dykes which eventually make up 100% of the section in layer 2B combined and the effects of high-temperature alteration. Starting with a simplified but plausible geologically-based model, we show that it is possible to synthetically generate the observed 2A Event, and gain insight into what controls its visibility and variability in refraction data. Our primary findings show that the 2A Event will only form and propagate in the base of layer 2A, above the level where the higher velocities dominate. We also show that the amplitude of the 2A Event is sensitive to the local velocity structure of the extrusive layer and is most visible when seismic energy is focused by a low velocity layer. Hence, we conclude that the 2A Event is not a simple caustic, as defined by geometrical optics, but instead caused by the incident seismic energy being briefly concentrated in a leaky waveguide close to, but above, the mean depth of the dykes and the onset of high temperature alteration.

How to cite: Hobbs, R. and Peirce, C.: What does the 2A Event observed in the oceanic crust actually represent?, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-5859,, 2021.

Vera Lay, Stefan Buske, Franz Kleine, John Townend, Richard Kellett, Martha Savage, Douglas R. Schmitt, Alexis Constantinou, Jennifer Eccles, Donald Lawton, Malcolm Bertram, Kevin Hall, Randolph Kofman, and Andrew Gorman

The Alpine Fault at the West Coast of the South Island (New Zealand) is a major plate boundary that is expected to rupture in the next 50 years, likely as a magnitude 8 earthquake. The Deep Fault Drilling Project (DFDP) aimed to deliver insight into the geological structure of this fault zone and its evolution by drilling and sampling the Alpine Fault at depth. Here we present results from a seismic survey around the DFDP-2 drill site in the Whataroa Valley where the drillhole almost reached the fault plane. This unique 3D seismic survey includes several 2D lines and a 3D array at the surface as well as borehole recordings. Within the borehole, the unique option to compare two measurement systems is used: conventional three-component borehole geophones and a fibre optic cable (heterodyne Distributed Vibration Sensing system (hDVS)). Both systems show coherent signals but only the hDVS system allowed a recording along the complete length of the borehole.

Despite the challenging conditions for seismic imaging within a glacial valley filled with sediments and steeply dipping valley flanks, several structures related to the valley itself as well as the tectonic fault system are imaged. The pre-processing of the seismic data also includes wavefield separation for the zero-offset borehole data. Seismic images are obtained by prestack depth migration approaches.

Within the glacial valley, particularly steep valley flanks are imaged directly and correlate well with results from the P-wave velocity model obtained by first arrival travel-time tomography. Additionally, a glacially over-deepened trough with nearly horizontally layered sediments is identified about 0.5 km south of the DFDP-2B borehole.

With regard to the expected Alpine fault zone, a set of several reflectors dipping 40-56° to the southeast are identified in a ~600 m wide zone between depths of 0.2 and 1.2 km that is interpreted to be the minimum extent of the damage zone. Different approaches image one distinct reflector dipping at 40°, which is interpreted to be the main Alpine Fault reflector. This reflector is only ~100 m ahead from the lower end of the borehole. At shallower depths (z<0.5 km), additional reflectors are identified as fault segments and generally have steeper dips up to 56°. About 1 km south of the drill site, a major fault is identified at a depth of 0.1-0.5 km that might be caused by the regional tectonics interacting with local valley structures. A good correlation is observed among the separate seismic data sets and with geological results such as the borehole stratigraphy and the expected surface trace of the fault.

In conclusion, several structural details of the fault zone and its environment are seismically imaged and show the complexity of the Alpine Fault at the Whataroa Valley. Thus, a detailed seismic characterization clarifies the subsurface structures, which is crucial to understand the transpressive fault’s tectonic processes.

How to cite: Lay, V., Buske, S., Kleine, F., Townend, J., Kellett, R., Savage, M., Schmitt, D. R., Constantinou, A., Eccles, J., Lawton, D., Bertram, M., Hall, K., Kofman, R., and Gorman, A.: 3D seismic imaging of the Alpine Fault and the glacial valley at Whataroa, New Zealand, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-9743,, 2021.

Damian Pasiecznik, Andrew Greenwood, Ludovic Baron, Florian Bleibinhaus, and György Hetényi

The Ivrea Verbano Zone (IVZ) is one of the most complete crust–upper mantle geological references in the world, and the Drilling the Ivrea-Verbano zone project (DIVE) aims to unravel the uncertainties below this area. Geophysical anomalies detected across the IVZ indicate that dense, mantle-like rocks are located at depths as shallow as ca. 3km. Thus, several geological, geochemical and geophysical studies are planned, including the drilling of a 4km deep borehole that will penetrate the Balmuccia Peridotite (Val Sesia, Italy) to approach and possibly cross the crust–mantle transition zone, and provide, for the first time, geophysical in-situ measurements of the IVZ.

One of the primary requirements before drilling is a seismic site characterization, to define with precision the correct positioning and orientation of the borehole, to assess potential drilling hazards and to allow for the spatial extrapolation of the borehole logs. For that goal, two joint geophysical surveys were performed in October 2020 in a collaboration of GFZ Potsdam, Université de Lausanne and Montanuniversität Leoben. First: a deep seismic survey, entitled SEIZE (SEismic imaging of the Ivrea ZonE), consisting of two approximately 15km-long seismic lines performed by GFZ Potsdam, that aim to resolve the deeper structure of the IVZ in the area, and second: a static seismic survey at the proposed drill site, entitled HiSEIZE (High-resolution SEismic imaging of the Ivrea ZonE), geared towards providing high-resolution seismic images of the uppermost few km at the proposed drill site.

The HiSEIZE survey, the subject of this study, was performed with a fixed spread of 200 vertical geophones and 160 3C-sensors, spaced at ca. 11m along three sub-parallel lines spaced 50-80m apart. Vibroseis source points were at 22m stations along a 2.4km line utilizing a high-frequency (12-140 Hz) 10s linear sweep with 3s listening time. In addition to this, the HiSEIZE receiver spread was active during the deep SEIZE survey, information that may be useful in determination of a velocity model through the Balmuccia Peridotite.

This project will not only provide site characterization for the DIVE project, but also contribute to understanding the structure of the Balmuccia Peridotite, its changes in depth and its relationship with the crustal-mantle transition.

Here we present the data and discuss the challenges of 3D pre-stack-imaging in an area of extreme topography.

How to cite: Pasiecznik, D., Greenwood, A., Baron, L., Bleibinhaus, F., and Hetényi, G.: A high-resolution seismic survey across the Balmuccia Peridotite, Ivrea Zone, Italy - Project DIVE phase two, site investigation, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-14388,, 2021.

Assel Akimbekova, Paolo Mancinelli, Massimiliano Rinaldo Barchi, Cristina Pauselli, and Giorgio Minelli



In the present study, starting from original measurement stations, we created the Bouguer anomaly map of Southern Italy with a reduction density of 2670 kg m-3. We perform a regional gravity modelling at crustal scale along the trace of the CROP-04 (on-shore) and MB6 (off-shore) deep seismic reflection profiles crossing the Southern Apennines and the Southern Tyrrhenian Sea. Along the 320 km-long modelled profile, we investigate crustal-scale sources for the observed gravity anomalies. 

After a compelling review of the published Moho geometries in the area, that were retrieved from either active or passive seismic methods, we test them in the observed gravity field through forward modeling of the Bouguer gravity anomalies. The comparison between the different Moho interpretations shows that the steepness of the subducting slab, the position of the step between the western (Tyrrhenian) and the eastern (Adriatic) Moho and Moho depth represent the main features influencing the observed Bouguer anomalies at crustal scale.

Finally, we provide a best-fitting model across both onshore and offshore areas. In the proposed best-fitting model, the wide wavelength and strong regional Bouguer anomalies correlate with the geometry of the Moho discontinuity and deep tectonic structures. On the other hand, the small-amplitude oscillations of the gravity anomalies were attributed to the low-density values of the Pliocene-Quaternary deposits both on- (e.g. the Bradanic trough) and off-shore (e.g. recent deposits in the Tyrrhenian sea bottom). Gravity minima correspond to the crustal doubling underneath the Southern Apennines where the Tyrrhenian Moho (~27 km depth) overlies the deeper Adriatic Moho (~50 km depth). The positive trend of the observed anomaly toward NE is related to the shallowing of the Adriatic Moho to depths of ~28 km in the Adriatic. Similarly, towards SW, the observed anomaly follows a positive trend towards the maxima located in the Central Tyrrhenian Sea. We model this trend as representative of crustal thinning and shallowing to values of ~12 km depth of the Tyrrhenian Moho. We also model a crustal transition from geometries and density values typical of a continental crust in the Adriatic domain towards a more oceanic structure and composition in the Tyrrhenian domain. This crustal model locates the westward flexure of the Adriatic Moho, mimicking the subduction of the Adriatic lithosphere beneath the Peri-Tyrrhenian block and locates step between the western (Tyrrhenian) and the eastern (Adriatic) Moho beneath the Apennines range.

The resulted gravity forward model provide contributions to the tectonic settings understanding of the area by providing a robust crustal model ranging from the Tyrrhenian Sea to the Apulian foreland.

 Finally, we believe that the proposed model can serve as a starting point for future studies investigating the upper crustal geometries in the area and addressing open questions about its relations with seismicity distribution.


How to cite: Akimbekova, A., Mancinelli, P., Rinaldo Barchi, M., Pauselli, C., and Minelli, G.: Forward  Modelling of Bouguer Anomalies along a transect of the Southern Apennines and the Tyrrhenian back-arc basin., EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-12246,, 2021.

Laura Gómez de la Peña, Ingo Grevemeyer, Anke Dannowsky, Adrià Meléndez, Christine Peirce, and Harm van Avendonk

About 25% of the Earth’s mid-ocean ridges spread at ultraslow rates of less than 20 mm/yr. However, most of these ultraslow spreading ridges are located in geographically remote areas, which hamper investigation. Consequently, how the crust forms and ages at such spreading centres, which traditional models predict to be magma-starved and cold, remains poorly understood.

CAYSEIS project was proposed to survey the Cayman Trough area in order to obtain new data that constraints the nature of the crust, tectonic structures, lithologies outcropping and hydrothermal processes taking place in this area, which includes the Mid Cayman ultra-slow spreading centre (MSCS) with spreading rates of ~15-17 mm/yr. Understanding the sub-seabed geophysical structure of the MCSC is key to understanding not only the lithologies and structures exposed at the seabed, but more fundamentally, how they are related at depth and what role hydrothermal fluid flow plays in the geodynamics of ultraslow spreading. CAYSEIS was a joint and multidisciplinary programme of German, British and US American top tier scientists designed for the obtaining of a new high-quality dataset, including 3D Wide-Angle Seismic (WAS), magnetic, gravimetric and seismological data.

We took leverage of the CAYSEIS dataset to invert a 3D velocity model of the Cayman Trough lithosphere using the Tomo3D code (Meléndez et al., 2015; 2019). This is one of the first times that the Tomo3D code is used for 3D inversion of real datasets. Thus, we are checking our results comparing them with travel time tomographic inversions of 2D lines and testing the different parameters to obtain the more accurate and higher resolution model as possible. The results of this experiment show not only the lithospheric structure along and across the MSCS, including the exhumed Ocean Core Complexes in the surrounding areas, but the 3D lithospheric configuration of the region which is important to understand the crustal formation processes and evolution of ultra-slow spreading settings.

How to cite: Gómez de la Peña, L., Grevemeyer, I., Dannowsky, A., Meléndez, A., Peirce, C., and van Avendonk, H.: 3D tomographic imaging of the Cayman Trough lithosphere: challenges, ongoing work and first results, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-1545,, 2021.

Seismic (coda) attenuation
Itahisa González Álvarez, Sebastian Rost, Andy Nowacki, and Neil Selby

Heterogeneities on the scale of the seismic wavelength in the Earth's crust and mantle cause complex wavefield fluctuations in time and amplitude which are known to affect velocity and source inversions, as well as other seismic characterisations. However, many seismic models ignore these heterogeneities for simplicity. As part of our longer-term goal to account for these, we attempt to rigorously and probabilistically characterise these lithospheric small-scale heterogeneities by combining a single-layer and a multi-layer energy flux models with a new Bayesian inference algorithm. The first technique characterizes energy losses to the ballistic arrival as intrinsic, diffusion and scattering quality factors, which allows us to compare the effects of these attenuation mechanisms on our data. With the second method, we can obtain synthetic coda envelopes for 1- and 2- layer models with different values of the correlation length and fractional velocity fluctuations in each layer. We then use the Metropolis-Hastings algorithm to sample the likelihood space and obtain the posterior probability distributions for each parameter and layer in the model. Our thorough testing of these methods reveals complicated trade-offs between the parameters and highly non-unique solutions, thus highlighting the importance of the Bayesian approach for scattering studies. Previous studies applying these methods used a more traditional grid search for their coda inversion, which may have affected their results. We applied this approach to a data set of over 300 events from three seismic arrays in Australia: Alice Springs array (ASAR), Warramunga Array (WRA) and Pilbara Seismic Array (PSA). The results from the single-layer energy flux model show that all quality factors take higher values for PSA than for the other two arrays, indicating that the structure beneath this array is less attenuating and heterogeneous than for the other arrays. Intrinsic and diffusion attenuation are strongest for ASAR, while scattering and total attenuation are similarly strong for ASAR and WRA. Our multi-layer model results show the crust is more heterogeneous than the lithospheric mantle for all arrays, with crustal values of the correlation length and velocity fluctuations being lower for PSA than for the other arrays, indicating the presence of weaker and smaller scale heterogeneity beneath this array. We attribute these differences and similarities in the attenuation and heterogeneity structure beneath the arrays to variations in the tectonic history of the areas they are located on. This new Bayesian approach to the multi-layer energy flux model, in combination with the single-layer model, not only allows us to determine and compare the different quality factors, but also gives us detailed information about the trade-offs and uncertainties in the determination of the scattering parameters, making it a useful tool for future scattering and small-scale structure studies. 

How to cite: González Álvarez, I., Rost, S., Nowacki, A., and Selby, N.: Small-scale lithospheric heterogeneity characterization using Bayesian inference , EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-11115,, 2021.

Panayiota Sketsiou, Luca De Siena, and David G. Cornwell

The North Anatolian Fault (NAF), a right-lateral strike-slip fault spanning 1500 km in length, stretches across northern Turkey and it marks the boundary between the Eurasian and Anatolian plates. Nucleating in the east, at the Karliova triple junction and reaching the Aegean Sea at the west, it is a particularly active fault zone with a series of migrating high-magnitude earthquakes. Using 6445 Z-component waveforms from a temporary seismic network in the area (Dense Array for North Anatolia – DANA), this study aims to investigate the western part of the NAF, which splays into a northern and southern branch. Coda attenuation imaging is utilised for imaging the absorption characteristics of the area, as it can be used as a marker for source and dynamic Earth processes due to its higher sensitivity to small variations of lithospheric properties compared to seismic velocity. The absorption structure is recovered by inverting for the coda attenuation quality factor, Qc, at frequencies between 3-18 Hz, using sensitivity kernels. The extensive seismicity in the area, as well as the density of the seismic stations, provide high-resolution models of 0.04-0.05 degrees in spacing. The scattering structure of the region is imaged using peak-delay time, which is used as a direct measure for multiple forward scattering. Preliminary results show a clear change in scattering between the Istanbul and Sakarya zones, north and south of the fault respectively, with the scattering increasing from north to south at lower frequencies and decreasing at higher frequencies. At a smaller scale, absorption and scattering anomalies appear to outline contrasting geological units beneath the DANA network.

How to cite: Sketsiou, P., De Siena, L., and Cornwell, D. G.: Coda-attenuation imaging of the North Anatolian Fault Zone, northern Turkey, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15749,, 2021.

Simona Gabrielli, Aybige Akinci, Ferdinando Napolitano, Luca De Siena, Edoardo Del Pezzo, and Guido Ventura

Between August and October 2016, the Central Apennines in Italy have been struck by a long-lasting seismic sequence, known as the Amatrice (Mw 6.0) - Visso (Mw 5.9) - Norcia (Mw 6.5) sequence. The cascading ruptures occurred in this sequence have been considered connected to the fluid migration in the fault network, as suggested by previous studies. The behaviour of fluids in the crust is crucial to understand earthquakes occurrence and stress changes since fluids reduce fault stability. It has long been understood that the seismic attenuation is strongly controlled by the structural irregularity and heterogeneities; micro-cracks and cavities, either fluid-filled or dry, temperature and pressure variations cause a decrease in seismic wave amplitude and pulse broadening. Hence seismic attenuation imagining is a powerful tool to be a relevant provenance of information about the influence and abundance of fluids in a seismic sequence.

The aim of this work is to separate scattering and absorption contributions to the total attenuation of coda waves and to provide their spatial and temporal variations at different frequency bands of these quantities using two datasets: the first one comprising 592 earthquakes occurred before the sequence (March 2013-August 2016) and the second one comprising 763 events (ML > 2.8) from the Amatrice-Visso-Norcia sequence. Scattering and absorption have been measured through peak-delay and coda-wave attenuation parameters (the latter inverted using frequency-dependent sensitivity kernels).

The preliminary results show a clear difference between the pre-sequence and sequence images, mainly at low frequencies (1.5 Hz), where we can define a spatial increase of scattering with time attributed to rock fracturing and fluid circulation. The coda attenuation tomography also demonstrates a clear variation between the pre-sequence and the sequence over series of time windows being before and after the largest main shocks of the seismic sequence, with an increase of the attenuation in space with decreasing time. The peak delay indicates a high scattering area corresponding to the Gran Sasso massif and L’Aquila zone, where an important seismic sequence (Mw 6.3) occurred in 2009.

How to cite: Gabrielli, S., Akinci, A., Napolitano, F., De Siena, L., Del Pezzo, E., and Ventura, G.: Spatio-temporal variation of Anaelastic and Scattering Seismic Attenuation during the 2016-2017 Central Italy Seismic Sequence, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-8802,, 2021.

Chiara Nardoni, Luca De Siena, Fabio Cammarano, Elisabetta Mattei, and Fabrizio Magrini

Strong lateral variations in medium properties affect the response of seismic wavefields. The Tyrrhenian Sea is ideally suited to explore these effects in a mixed continental-oceanic crust that comprises magmatic systems. The study aims at investigating the effects of crustal thinning and sedimentary layers on wave propagation, especially the reverberating (e.g., Lg) phases, across the oceanic basin. We model regional seismograms (600-800 km) using the software tool OpenSWPC (Maeda et al., 2017, EPS) based on the finite difference simulation of the wave equation. The code simulates the seismic wave propagation in heterogeneous viscoelastic media including the statistical velocity fluctuations as well as heterogeneous topography, typical of mixed settings. This approach allows to evaluate the role of interfaces and layer thicknesses on phase arrivals and direct and coda attenuation measurements. The results are compared with previous simulations of the radiative-transfer equations. They provide an improved understanding of the complex wave attenuation and energy leakage in the mantle characterizing the southern part of the Tyrrhenian Sea and the Italian peninsula. The forward modelling is to be embedded in future applications of attenuation, absorption and scattering tomography performed with MuRAT (the Multi-Resolution Attenuation Tomography code – De Siena et al. 2014, JVGR) available at

How to cite: Nardoni, C., De Siena, L., Cammarano, F., Mattei, E., and Magrini, F.: Finite difference forward modelling across the Tyrrhenian basin, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-1745,, 2021.

Traveltime tomography and interpretation
Gregor Rajh, Josip Stipčević, Mladen Živčić, Marijan Herak, and Andrej Gosar

The investigated area of the NW Dinarides is bordered by the Adriatic foreland, the Southern Alps, and the Pannonian basin at the NE corner of the Adriatic Sea. Its complex crustal structure is the result of interactions among different tectonic units. Despite numerous seismic studies taking place in this region, there still exists a need for a detailed, smaller scale study focusing mainly on the brittle part of the Earth's crust. Therefore, we decided to investigate the velocity structure of the crust using concepts of local earthquake tomography (LET) and minimum 1-D velocity model. Here, we present the results of the 1-D velocity modeling and the catalogue of the relocated seismicity. A minimum 1-D velocity model is computed by simultaneous inversion for hypocentral and velocity parameters together with seismic station corrections and represents the best fit to the observed arrival times.

We used 15,579 routinely picked P wave arrival times from 631 well-located earthquakes that occurred in Slovenia and in its immediate surroundings (mainly NW Croatia). Various initial 1-D velocity models, differing in velocity and layering, were used as input for velocity inversion in the VELEST program. We also varied several inversion parameters during the inversion runs. Most of the computed 1-D velocity models converged to a stable solution in the depth range between 0 and 25 km. We evaluated the inversion results using rigorous testing procedures and selected two best performing velocity models. Each of these models will be used independently as the initial model in the simultaneous hypocenter-velocity inversion for a 3-D velocity structure in LET. Based on the results of the 1-D velocity modeling, seismicity distribution, and tectonics, we divided the study area into three parts, redefined the earthquake-station geometry, and performed the inversion for each part separately. This way, we gained a better insight into the shallow velocity structure of each subregion and were able to demonstrate the differences among them.

Besides general structural implications and a potential to improve the results of LET, the new 1-D velocity models along with station corrections can also be used in fast routine earthquake location and to detect systematic travel time errors in seismological bulletins, as already shown by some studies using similar methods.

How to cite: Rajh, G., Stipčević, J., Živčić, M., Herak, M., and Gosar, A.: Crustal velocity structure beneath the NW Dinarides from 1-D hypocenter-velocity inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-8226,, 2021.

Ines Hamak, Piedade Wachilala, José Borges, Nuno Dias, Inês Rio, and Mourad Bezzeghoud

This work puts in light the several steps followed to obtain a 3D velocity model in Arraiolos, a region located in central Portugal. After the earthquake of January 2018 occurred, a set of stations were deployed around the main shock area and has recorded the aftershock sequence during a period of six months.

The first stage of this study used a set of data recorded along the 1st month by 21 temporary seismological stations. 317 aftershocks were used to invert a 3D P and S  velocity model, using LOTOS program, and showing an agglomeration of events in one local point leading to a poor resolution. Therefore, we added more stations and data to the second stage of study by integrating 437 aftershocks recorded during a period of 6 months by a set of 34 stations. The tomographic inversion of this extended aftershock sequence has shown a significant improvement of the 3D velocity model resolution and suggesting an alignment of the seismic events cluster. However, the imaged crustal volume was still too small and possessing low resolution on the edges of the area. To fix this issue, additional data and seismological stations were integrated to the study in order to increase the area of interest and cover it entirely in terms of ray density.

The step which we are currently conducting concerns the location of new events followed by their integration to the tomographic study using IPMA and DOCTAR station network records. Since the later phases PmP and SmS has the potential to increase the ray coverage as similarly as the resolution of an area, we will hopefully obtain, after their integration, significant improvements in terms of accuracy and reliability of the crustal image. The main purpose of this new stage of study is to finally provide significant interpretations and figure out precisely the tectonic processes having generated the Arraiolos seismicity.

Thanks are due to FCT for the financial support to the ICT project (UID/GEO/04683/2013) with the reference POCI-01- 0145-FEDER-007690, to the IDL project (UIDB/50019/2020 – IDL).

How to cite: Hamak, I., Wachilala, P., Borges, J., Dias, N., Rio, I., and Bezzeghoud, M.: Improved three-dimensional image of the tomographic inversion of the Arraiolos aftershock sequence., EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-8079,, 2021.

Donato Talone, Rita de Nardis, Giusy Lavecchia, and Luca De Siena

Seismic tomography can be applied to different scales. Over the last two decades, monitoring systems, technical innovations and methodologies have substantially improved, resulting in accurate tomographic images at the global and local scales. Nowadays it is easy to perform travel-time tomography with local seismicity thanks to the increasing density of seismic stations. Nevertheless, it is unlikely to have earthquakes that properly cover the whole studied area at the requested depths. For this reason, many tomographic images are obtained with teleseisms and both far and local earthquakes.

Here, we realized a Local Earthquake Tomography (LET) in an area of high seismic hazard in central-southern Italy, extending from L’Aquila to Benevento, to benchmark the iterative non-linear Fast-Marching code FMTOMO (Rawlinson and Sambridge, 2004) at intracontinental scale. The primary aim is to analyse and discuss the influence of both the inversion parameters and the grid sizes on the inversion results. Special attention was devoted to setting damping factors and smoothing parameters and to study how they can affect the tomographic images and their reliability.

We used 5712 local events (0.2<ML<5.1) recorded by 38 stations of the Italian Seismic network; we jointly inverted 71221 P and S arrival times to obtain Vp and Vs model. We selected earthquakes having: (1) a root-mean-square (RMS) residual less than 0.5 s, (2) more than 10 phases (P and S), (3) azimuth gap less than 180, (4) residual of each phase less than 0.5 s, (5) a depth between 0.5 and 30 km. We used a single layer of 35 km in depth and a grid area extending 162 km in latitude and 245 km in longitude with a node spacing of about 5 km in each direction. As a starting velocity model, we chose a mono-dimensional one of Trionfera et al. (2020).

Using these well-localized earthquakes, we observed low residuals variability despite a full investigation of damping and smoothing parameters. Furthermore, the regularization parameters we obtained are one or two order of magnitude lower than those estimated at the wider scales.

Because of the uncertainties in the depth of events, the fast-marching code needs several nodes above and below the grid set for earthquakes to move sources during each hypocentral inversion. As a consequence, when inverting for both velocity and hypocentral location, FMTOMO performs the calculation even for a wide boundary area without earthquakes, which causes a loss of computational speed.

After properly tuning the inversion parameters, FMTOMO gives reliable and high-resolution tomographic images. We found a good agreement with surface geology and regional tectonic structures, demonstrating that the code works well in areas with such complex geology.

How to cite: Talone, D., de Nardis, R., Lavecchia, G., and De Siena, L.: Regional fast-marching tomography in intracontinental settings: preliminary analysis of limitations and potential, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-12696,, 2021.

Ha Vinh Long, Hsin-Hua Huang, Le Minh Nguyen, Van Duong Nguyen, Quang Khoi Le, Thi Giang Ha, Dinh Quoc Van, Bor-Shouh Huang, and Tu Son Le

The collision between the Eurasian and Indian plates since about 65 Ma caused extensive deformation in the Tibetan Plateau and surrounding areas. Northern Vietnam located in the southeastern Himalayan syntaxis is directly influenced by the collision and extrusion through the Red River shear zone that runs from southeastern Tibet to the South China Sea. Knowledge of crustal structure characteristics in northern Vietnam and across the Red River shear zone is crucial to improve our understanding not only of seismic hazards in the region but also of the regional Himalayan tectonic evolution as a whole. Seismic tomography is one of few methods that allows to study the subsurface structures effectively. In this study, we perform a joint tomographic inversion for northern Vietnam integrating the P- and S-direct waves traveling in the crust and the head waves along the Moho waves arrival time from more than 1000 earthquakes observed by Vietnamese networks. The obtained velocity model shows a good correlation with shallow geological features but also some complexity at crustal-scale. Several velocity anomalies bounded by and across the fault zones are revealed and discussed


Keyword: Traveltime tomography, Northern Vietnam, VpVs ratio, crustal structure, Vietnam seismic network

How to cite: Long, H. V., Huang, H.-H., Nguyen, L. M., Nguyen, V. D., Le, Q. K., Ha, T. G., Van, D. Q., Huang, B.-S., and Le, T. S.: Crustal seismic structure of Red-River shear zone and surrounding area, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-2001,, 2021.

Receiver functions
Matteo Scarponi, György Hetényi, Jaroslava Plomerová, and Stefano Solarino

We present results from a joint inversion study of new seismic and gravity data to constrain a 2D high-resolution image of one of the most prominent geophysical anomalies of the European Alps: the Ivrea geophysical body (IGB). Our work exploits both new data and multidisciplinary a priori constraints, to better resolve the shallow crustal structure in the Ivrea-Verbano zone (IVZ), where the IGB is known to reach anomalously shallow depths and partially outcrop at the surface.

A variety of previous studies, ranging from gravity surveys to vintage refraction seismics and recent local earthquake tomographies (Solarino et al. 2018, Diehl et al. 2009), provide comprehensive but spatially sparse information on the IGB structure, which we aim at investigating at higher resolution, along a linear profile crossing the IVZ. To this purpose, we deployed 10 broadband seismic stations (MOBNET pool, IG CAS Prague), 5 km spaced along a linear West-East profile, along Val Sesia and crossing Lago Maggiore. This network operated for 27 months and allowed us to produce a new database of ca. 1000 seismic high-quality receiver functions (RFs). In addition, we collected new gravity data in the IVZ, with a data coverage of 1 gravity point every 1-2 km along the seismic profile. The newly collected data was used to set up an inversion scheme, in which RFs and gravity anomalies are jointly used to constrain the shape and the physical property contrasts across the IGB interface.

We model the IGB as a single interface between far-field constraints, whose geometry is defined by the coordinates of four nodes which may vary in space, and  density and VS shear-wave velocity contrasts associated with the interface itself, varying independently. A Markov chain Monte Carlo (MCMC) sampling method with Metropolis-Hastings selection rule was implemented to efficiently explore the model space, directing the search towards better fitting areas.

For each model, we perform ray-tracing and RFs migration using the actual velocity structure both for migration and computation of synthetic RFs, to be compared with the observations via cross-correlation of the migration images. Similarly, forward gravity modelling for a 2D density distribution is implemented and the synthetic gravity anomaly is compared with the observations along the profile. The joint inversion performance is the product of these two misfits.

The inversion results show that the IGB reaches the shallowest depths in the western part of the profile, preferentially locating the IGB interface between 3 and 7 km depth over a horizontal distance of ca. 20 km (between Boccioleto and Civiasco, longitudes 8.1 and 8.3). Within this segment, the shallowest point reaches up to 1 km below sea level. The found density and velocity contrasts are in agreement with rock physics properties of various units observed in the field and characterized in earlier studies.

How to cite: Scarponi, M., Hetényi, G., Plomerová, J., and Solarino, S.: Joint receiver function and gravity inversion for new constraints on the Ivrea body structure: a 2D high-resolution view along the Val Sesia profile (N. Italy)., EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-9882,, 2021.

Harry Telajan Linang, Amy Gilligan, Jennifer Jenkins, Tim Greenfield, Felix Tongkul, Sri Widiyantoro, and Nick Rawlinson

Borneo is located at the centre of Southeast Asia, which is one of the most active tectonic regions on Earth due to the subduction of the Indo-Australian plate in the south and the Philippines Sea plate in the east. Borneo resides on the leading edge of the Sundaland block of the Eurasian plate and exhibits lower rates of seismicity when compared to the surrounding regions due to its intraplate setting. Sulawesi, an island which lies just southeast of Borneo, is characterised by intense seismicity due to multiple subduction zones in its vicinity. The tectonic relationship between the two islands is poorly understood, including the provenance of their respective lithospheres, which may have Eurasian and/or East Gondwana origin.

Here, we present recent receiver function (RF) results from temporary and permanent broadband seismic stations in the region, which can be used to help improve our understanding of the crust and mantle lithosphere beneath Borneo and Sulawesi. We applied H-K stacking, receiver function migration and inversion to obtain reliable estimates of the crustal thickness beneath the seismic stations. Our preliminary results indicate that the crust beneath Sabah (in northern Borneo), which is a post-subduction setting, appears to be much more complex and is overall thicker (more than 35 km) than the rest of the island. In addition, we find that crustal thickness varies between different tectonic blocks defined from previous surface mapping, with the thinnest crust (23 to 25 km) occurring beneath Sarawak in the west-northwest as well as in the east of Kalimantan.

We also present preliminary results from Virtual Deep Seismic Sounding (VDSS) in northern Borneo, where from the RF results we know that there is thick and complex crust. VDSS is able to produce well constrained crustal thickness results in regions where the RF analysis has difficulty recovering the Moho, likely due to complexities such as thick sedimentary basins and obducted ophiolite sequences.

How to cite: Linang, H. T., Gilligan, A., Jenkins, J., Greenfield, T., Tongkul, F., Widiyantoro, S., and Rawlinson, N.: Seismic Structure and Tectonic Evolution of Borneo and Sulawesi, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-16456,, 2021.

Ambient seismic noise
Roberto Manzo, Lucia Nardone, Guido Gaudiosi, Claudio Martino, Danilo Galluzzo, Francesca Bianco, and Rosa Di Maio

Following the MD4.0 (Mw3.9) earthquake of August 21 2017 which occurred on the Ischia island (Naples, southern Italy), the local monitoring seismic network was significantly improved in terms of both number of stations and instrumentation performance. Due to the considerable amount of collected data, in particular of seismic noise recorded at broadband stations, some efforts have been addressed in particular to the definition of a 1D average velocity model effective for the whole island. This is an important scientific step because, in complex volcanic areas, the use of reliable velocity models is essential for an accurate localization of local earthquakes. In this work, the main target is to retrieve a pseudo-3D velocity model of the Ischia island. Specifically, we inverted H/V curves and frequency peaks evaluated at about twenty sites to obtain a velocity profile for each of the investigated measurement points. Taking into account that the H/V frequency peak depends on both velocity and thickness of layers, for each site we applied an inversion process fixing the velocities and modifying the thicknesses in order to obtain the corresponding 1D velocity models. We are quite enough confident about the robustness of models, since during the inversion process, we achieved a good convergence towards the best-fit solutions. Then, a pseudo-3D velocity model was obtained by contouring the 1D models of each station site to highlight possible lateral variations of the layer thicknesses and to reconstruct the morphology of the deeper interface characterized by a high impedance contrast. A good correspondence between the pseudo-3D model and the geological features of the island was observed, especially in the northern sector where most of the stations is installed. In particular, the top of the high-impedance contrast interface appears deeper in the northern coastal areas and shallower in the central sector. This is in agreement with the structural setting of the island likely due to the resurgence of Mount Epomeo.


How to cite: Manzo, R., Nardone, L., Gaudiosi, G., Martino, C., Galluzzo, D., Bianco, F., and Di Maio, R.: A Pseudo-3D Vs Velocity Model of Ischia Island (Italy) by HVSR spectral ratio inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-13610,, 2021.

Juan Pinzon, Susana Custódio, Graça Silveira, Luis Matias, and Frank Krüger

The Gloria fault is a strike-slip oceanic plate boundary fault, which has remained poorly studied due mostly to its remote location in the north Atlantic Ocean. The fault has hosted some of the largest strike-slip earthquakes in the oceanic domain, notably the 1941 M8.3 and the 1975 M8.1 earthquakes, and generating tsunamis is the surrounding areas.

The seismic data used for this study was recorded by 12 broadband ocean bottom seismometers (OBSs) located about 100 km north of the Gloria fault during a 10-month experiment. The dataset has been used before to image crustal and mantle discontinuities using receiver function analysis and to infer the S-wave velocity structure of the oceanic lithosphere north of the Gloria fault from P-wave polarization. These studies indicate a slight crustal thickening towards the Gloria fault, as well as an increase in uppermost mantle S-wave velocities towards the fault.

In this study, we use ambient noise surface wave tomography to find the velocity structure beneath the OBS deployment. First, we present a 1D shear-velocity model obtained from inversion of the average fundamental mode Rayleigh and Love wave group and phase velocities. In addition, the hydrophone is also used to better constrain the inversion at shallow depths, because the hydrophone shows a clear fundamental mode without interference of the first higher mode. Because of the short interstation distances of the array, it is not possible to extract the dispersion curves at periods longer than ~16 s. To compute the Vs inversions, we used the code SURF96 (Herrmann and Ammon, 2004) and consider a water layer in the initial model for Rayleigh waves, because these waves are affected by the water layer. Our results show an upper mantle low-velocity zone, which may be related to serpentinization due to the proximity of the Gloria Fault. Finally, we present the lateral variations of group and phase velocities, as a function of period obtained using FMST (Rawlinson and Sambridge, 2005), which show strong contrast velocity anomalies at the center of the array at short periods (shallow depths).

The authors acknowledge support from the Portuguese FCT – Fundação para a Ciência e a Tecnologia, I.P., within the scope of project UTAP-EXPL/EAC/0056/2017 and with the FCT grant PD/BD/135069/2017 - IDL.

How to cite: Pinzon, J., Custódio, S., Silveira, G., Matias, L., and Krüger, F.: Ambient noise surface wave inversion using OBS data from the north Atlantic, north of the Gloria fault., EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-7188,, 2021.

Mohamadhasan Mohamadian Sarvandani, Emanuel Kästle, Lapo Boschi, Sylvie Leroy, and Mathilde Cannat

Passive seismic interferometry (ambient-noise seismology) is an increasingly popular, eco-friendly, relatively inexpensive exploration geophysics tool, to map S-wave velocity in the Earth’s crust. This method has not yet been applied widely to marine exploration. The purpose of this study is to investigate the crustal structure of a quasi-amagmatic portion of the Southwest Indian Ridge by interferometry, and to examine the performance and reliability of interferometry in marine exploration. To achieve this goal, continuous vertical-component recordings from 43 ocean bottom seismometers (OBS) deployed during the SISMO-SMOOTH cruise (2014) were utilized. Recorded signals span frequencies between 0.1Hz and 3Hz. We show that reliable estimates of the Green’s function are obtained for many station pairs, by cross-correlation in the frequency domain. The comparison of the cross-correlations with the theoretical Green’s (Bessel) function provides one Rayleigh-wave dispersion curve per station pair; dispersion curves are then averaged, and inverted through a conditional neighborhood algorithm to determine a 1D S-wave velocity model, that we estimate to be well constrained within the crust. Our S-wave velocity model is analyzed and interpreted with geological information, and independent geophysical studies in the region of interest, as well as other areas characterized by similar tectonically-dominated, quasi amagmatic spreadings.

How to cite: Mohamadian Sarvandani, M., Kästle, E., Boschi, L., Leroy, S., and Cannat, M.: Passive seismic interferometry of the ultraslow-spreading Southwest Indian Ridge , EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-6132,, 2021.

Zhengyang Zhou, Douglas Wiens, and Andrew Lloyd

The Antarctic continent with its large ice sheets provides a unique environment to investigate the response of the solid Earth to ice mass change. A key requirement of such studies is high-resolution seismic images of the crust and upper mantle, which can be used to estimate the region’s viscous structure. Likewise, these images are key to understanding the region’s geologic history and underlying geodynamic processes. Although the existing transverse isotropic seismic model ANT-20(Lloyd et al., 2020) and azimuthally anisotropic seismic model ANT-30 (Lloyd et al., in prep) have regional-scale resolution from the upper mantle to the transition zone, there is a need for higher resolution within the uppermost mantle (< 75 km) and crust of Antarctica. In this study, we use the ANT-30 model (Lloyd et al., in prep), a 3D seismic model from earthquake data, as a starting model. We seek to improve its resolution within the upper ~100 km of the Antarctic mantle by fitting three-component ambient noise correlograms computed from broadband records collected in Antarctica over the past 20 years. This includes data from recent temporary arrays such as TAMSEIS, AGAP, TAMNNET, RIS, POLENET/ANET, and UKANET. The three-component cross-correlations of station pairs are calculated and properly rotated to extract ambient noise surface waves that include both Rayleigh and Love waves, which show excellent signal-to-noise ratio between 15 to 70 seconds. The benefit of including this data is twofold: (1) it provides surface wave observations down to 15 s, as opposed to 25 s and (2) it provides shorter intercontinental paths, which were absent due to the region’s earthquake distribution. We then use the software package SPECFEM3D_GLOBE to iteratively improve the 3-D earth model, minimizing the nondimensionalized traveltime phase misfit between the observed and synthetic waveforms. The preliminary results indicate a stronger positive radial anisotropy (VSH > VSV) in the lower crust and uppermost mantle for West Antarctica and part of East Antarctica.  With more iterations, smaller-scale detail can be revealed by the new ambient noise data, resulting in a more reliable uppermost mantle and crustal structure.

How to cite: Zhou, Z., Wiens, D., and Lloyd, A.: Adjoint Seismic Tomography of the Antarctic Continent Incorporating Green's Functions from Ambient Noise Correlation, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-9357,, 2021.

omid Bagherpur Mojaver and Fiona Darbyshire

Long-duration stacks of ambient seismic noise cross-correlation can be used to generate high-resolution images of the lithosphere. In this study, we investigate the crustal structure beneath southeastern Canada and the northeastern USA, using an ambient noise tomography technique. Our study area covers the Phanerozoic northern Appalachians and the Proterozoic eastern Grenville Province, recording a complex tectonic history since ~1 Ga. Our datasets include continuous records of vertical component time series, recorded by 69 stations belonging to 7 seismograph networks over a more than two-year period. The ambient seismic noise directionality and seasonality variations of our datasets are analyzed in detail, and possible noise source locations are proposed in the Atlantic and Pacific oceans. Our analysis suggests strong variations of dominant seismic noise sources at both Primary (11-20 s) and Secondary (5-10 s) bands in various months, with different observed patterns at these passband periods. Our tomographic models indicate complex and strong variations of Rayleigh wave phase velocities across the study area, providing us evidence to discuss tectonic implications. The resulting Rayleigh wave phase velocity maps suggest generally slower velocities beneath the Appalachians than the Grenville province. A sharp velocity contrast is observed across the Grenville Province-Appalachian domain boundary at periods sensitive to the lower crust, suggesting a step-like geometry of the Moho interface beneath this area.

How to cite: Bagherpur Mojaver, O. and Darbyshire, F.: Crustal structure, seasonal and directional variations of ambient seismic noise in SE Canada and the NE USA, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-6124,, 2021.

Giovanni Diaferia, Fabrizio Magrini, Matthew Agius, and Fabio Cammarano

The dynamics of crustal extension and the crust-mantle interaction in the Central-Western Mediterranean and Italian peninsula (i.e. Liguro-Provençal and Tyrrhenian Basin), and plate convergence (i.e. Alpine and Apennines chains) are key for the understating of the current geodynamics setting and its evolution in the region. However, open questions such as the style, depth and extent of the deformation still exist despite the wealth of seismological and non-seismological data acquired in the past decades. In this context, it is necessary to provide improved subsurface models in terms of seismic velocities, from which better constraints on the geodynamic models can be derived.

We use seismic ambient noise for retrieving phase velocities of Rayleigh and Love waves in the 4-35 s period range, using private (LiSard network in Sardinia island) and publicly available continuous recordings from more than 500 seismic stations. Considering the excellent coverage and the short period of recovered phase velocities, our study aims to provide an unprecedented, high-resolution image of the shallow crust and uppermost mantle.

We employ a Bayesian trans-dimensional, Monte Carlo Markov chain inversion approach that requires no a-priori model nor a fixed parametrization. In addition to the (isotropic) shear wave velocity structure, we also recover the values of radial anisotropy (ξ=(VSH/VSV)2) as a function of depth, thanks to the joint inversion of both Rayleigh and Love phase velocities.

Focusing on radial anisotropy, this appears clearly uncoupled with respect to the shear wave velocity structure. The largest negative anisotropy anomalies (VSH<VSV, ξ<0.9) are found in the Liguro-Provençal and western Tyrrhenian basins in the top 10-15 km, suggesting a common structural imprint inherited during the extensional phases of such basins. Conversely, the eastern Tyrrhenian basin shows positive radial anisotropy (VSH>VSV, ξ>1.1) within the same depth range. This evidence, combined with the observed shear wave velocities typical of the uppermost mantle, corroborates the presence of a sub-horizontal asthenospheric flow driving the current extension and oceanization of the eastern Tyrrhenian basins.

Moving towards the Italian mainland, a strong anomaly of negative anisotropy appears in the eastern portion of the Apennines chain. We relate such an anisotropic signal with the ongoing compressive regime affecting the area. Here, the high-angle thrust faults and folds, that accommodates the horizontal shortening, obliterate the horizontal layering of the sedimentary deposits, currently constituting the flanks of the fold system.

Our results suggest that the combination of radial anisotropy and shear wave velocities can unravel key characteristics of the crust and uppermost mantle, such as inherited or currently active structures resulting from past or ongoing geodynamic processes.

How to cite: Diaferia, G., Magrini, F., Agius, M., and Cammarano, F.: Seismic radial anisotropy in Central-Western Mediterranean and Italian peninsula from ambient noise recordings, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15238,, 2021.

Matthew Agius, Fabrizio Magrini, Giovanni Diaferia, Fabio Cammarano, Claudio Faccenna, Francesca Funiciello, Emanuel Kästle, and Mark van der Meijde

The Central Mediterranean, the area encompassing Italy, Sardinia, Tunisia and Libya, is characterised by multiple tectonic processes (plate convergence, subduction, and backarc extension). The evolution and interaction of the plate margins within this relatively small area are still being unravelled particularly at the adjacent region known as the Sicily Channel located between Sicily, Tunisia, Libya and Malta. This Channel is characterised by a seismically and volcanically active rift zone. Much of the observations we have today for the southern parts of the Calabrian arc are either limited to the surface and the upper crust, or are broader and deeper from regional seismic tomography, missing important details about the lithospheric structure and dynamics. The project GEOMED ( addresses this issue by processing all the seismic data available in the region in order to understand better the geodynamics of the Central Mediterranean.

We measure Rayleigh- and Love-wave phase velocities from ambient seismic noise recordings to infer the structures of the Central Mediterranean, from the Central Apennines to the African foreland, with a special focus on the Sicily Channel Rift Zone (SCRZ). The phase-velocity dispersion curves have periods ranging from 5 to 100 seconds and sample through the entire lithosphere. We invert the dispersion data for isotropic and polarised shear velocities with depth and infer crustal thickness and patterns of radial anisotropy. We find that continental blocks have thick crust (30-50 km), whereas beneath the SCRZ the crust is thin (<25 km), and thinner beneath the Tyrrhenian Sea. Beneath the SCRZ and the Tyrrhenian Sea, the crustal shear velocities are characterised by positive radial anisotropy (VSH>VSV) indicative of horizontal flow or extension, whereas the uppermost mantle is characterised by slow shear velocities indicative of warmer temperatures and strong negative radial anisotropy (VSH>VSV) indicative of vertical flow. We discuss the relevance of these findings together with other geophysical studies such as the regional seismicity and GPS velocity vectors to identify the rifting process type of the SCRZ.

This project has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 843696.

How to cite: Agius, M., Magrini, F., Diaferia, G., Cammarano, F., Faccenna, C., Funiciello, F., Kästle, E., and van der Meijde, M.: Geodynamics of the Central Mediterranean (GEOMED) inferred from ambient seismic noise, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-6028,, 2021.

Chairpersons: Giovanni Diaferia, Michal Malinowski, Matthew Agius
Interferometry and full-waveform inversion
Sneha Singh, Yann Capdeville, Heiner Igel, Navid Hedjazian, and Thomas Bodin

Wavefield gradient instruments, such as rotational sensors and DAS systems, are becoming more and more accessible in seismology. Their usage for Full Waveform Inversion (FWI) is in sight. Nevertheless, local small-scale heterogeneities, like geological inhomogeneities, surface topographies, and cavities are known to affect wavefield gradients. This effect is in fact measurable with current instruments. For example, the agreement between data and synthetics computed in a tomographic model is often not as good for rotation as it is for displacement.

The theory of homogenization can help us understand why small-scale heterogeneities strongly affect wavefield gradients, but not the wavefield itself. It tells us that at any receiver measuring wavefield gradient, small-scale heterogeneities cause the wavefield gradient to couple with strain through a coupling tensor J. Furthermore, this J is 1) independent of source, 2) independent of time, but 3) only dependent on the receiver location. Consequently, we can invert for J based on an effective model for which synthetics fit displacement data reasonably well. Once inverted, J can be used to correct all other wavefield gradients at that receiver.

Here, we aim to understand the benefits and drawbacks of wavefield gradient sensors in a FWI context. We show that FWIs performed with rotations and strains are equivalent to that performed with displacements provided that 1) the number of data is sufficient, and 2) the receivers are placed far away from heterogeneities. In the case that receivers are placed near heterogeneities, we find that due to the effect of these heterogeneities, an incorrect model is recovered from inversion. In this case therefore, the coupling tensor J needs to be taken into account for each receiver to get rid of the effect.

How to cite: Singh, S., Capdeville, Y., Igel, H., Hedjazian, N., and Bodin, T.: Full-waveform inversion with wavefield gradients, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15143,, 2021.

Varun Kumar Singla and Ivan Lokmer

While there are seismic techniques which make use of surface waves in imaging the subsurface, there are also those where these types of waves are considered coherent noise. Important examples where the surface waves may significantly degrade the obtained images include different types of reflection seismic surveys (e.g., shallow surveys for engineering, environmental and groundwater investigations, and deep surveys for imaging hydrocarbon reservoirs). In a strongly heterogeneous medium (encountered typically in onshore surveys), the conventional methods for attenuating these waves (such as f-k "velocity" filtering) often do not give satisfactory results.

Seismic interferometry is a data-driven approach that offers a viable alternative for removal of surface waves from active seismic records. In this approach, the reflection data of several sources is considered and for each source, the seismic signals at a pair of receivers are cross-correlated to produce the surface wavefield between the receivers. The cross-correlated waveforms are then summed over all the sources to obtain the "interferometric" signal for the considered receiver pair. During this summation, the reflection and non-physical events cancel out due to the variable differences in the travel times to the considered receiver pair from different sources. The "interferometric" signal consequently contains predominantly the surface waves and this makes it conducive for adaptive subtraction (or filtering) from the original records. This study investigates the efficacy of the commonly used filtering techniques in interferometry to remove the surface waves from active seismic records. For this, the reflection data of a complex 2-D elastic medium is simulated and the filtering techniques are applied to this data. The limitations of these techniques inferred from the quality of the filtered data are discussed and possible remedies to overcome them are suggested.

How to cite: Singla, V. K. and Lokmer, I.: Efficacy of Seismic Interferometry in Removing Surface Waves from Active Seismic Records, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-279,, 2020.

Amr El-Sharkawy, Thomas Meier, Christian Hübscher, Sergei Lebedev, Anke Dannowski, Heidrun Kopp, Jan H. Behrmann, Andrew McGrandle, and Mona Hamada

The Earth’s oldest oceanic lithosphere preserved in-situ is in the eastern Mediterranean Sea. It can offer essential information on the oceanic plate evolution. Yet, its thickness and other properties have been difficult to determine by means of seismic imaging due to the high heterogeneity of the region. Here, we combine a large, new surface wave dataset with published wide-angle data in order to map the properties and lateral variability of the oceanic lithosphere, as well as the ocean-continent transition in the easternmost Mediterranean beneath the Levant Basin. We use stochastic joint inversion of broad band, phase-velocity dispersion measurements and seismic refraction P-wave velocity models to obtain 1-D, shear wave velocity models down to 300 km depth and compare the structure beneath the Ionian Sea and the Levant Basin. The thickness of the crust is about 16.4 ± 3 km and 22.3 ± 2 km beneath the chosen locations within the Ionian Sea and the Levant Basin, respectively. The Poisson’s ratio of about 0.32 and Vp/Vs of about 1.93 in the crystalline crust, yielded by the inversion, confirm the presence of oceanic crust beneath the Ionian Sea. The thickness of the Ionian oceanic lithosphere is around 180 km, whereas the continental lithosphere beneath the eastern Levant Basin is ~70 km thick, with low crustal Vp/Vs (~1.7) and Poisson’s (~0.24) ratios. According to 3-D shear wave velocity tomography using the surface wave data, the thickness of the oceanic lithosphere increases from the Triassic Ionian Sea towards the Permian-Carboniferous Libyan Sea and Herodotus Basin. Thicknesses of the Permo-Triassic oceanic lithosphere considerably larger than 100 km indicate that oceanic lithosphere can thicken by cooling substantially beyond the limits suggested by the plate cooling model. The transition from oceanic to continental lithosphere occurs at about 31°E in the crust, as indicated by magnetic and gravity measurements. The continental mantle lithosphere further to the east of this boundary is ~150 km thick beneath the westernmost Levant Basin, as indicated by shear wave velocity tomography and long wavelength gravity anomalies, and strongly thins eastward towards the area of the Levantine Coast and the Dead Sea Fault. The localization of the lithospheric deformation and crustal seismicity along the Dead Sea Fault correlates spatially with the thinning of the underlying continental lithosphere.

Key words: surface wave tomography, wide angle seismic imaging, joint inversion, Vp/Vs and Poisson’s ratios, eastern Mediterranean, Oceanic Lithosphere, Continental Lithosphere, Dead Sea Fault.

How to cite: El-Sharkawy, A., Meier, T., Hübscher, C., Lebedev, S., Dannowski, A., Kopp, H., Behrmann, J. H., McGrandle, A., and Hamada, M.: Permo-Triassic Oceanic and Adjacent Continental Lithosphere in the Eastern Mediterranean, from surface-wave and wide-angle imaging, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-7231,, 2021.

Senlin Yang, Peng Jiang, Yuxiao Ren, and Xinji Xu

The seismic full waveform inversion (FWI), as one of important ways to obtain the seismic wave velocity, has made rapid development in the last decade. In response to problems of cycle-skipping artifacts, dependence on the initial model, and low-frequency information in FWI, researchers have made many improvements, such as multi-scale envelope inversion and low-frequency extension. Recently, deep learning has been also adopted seismic data processing and interpretation, because of its strong nonlinear mapping ability. However, these works depend on labels used for training heavily, especially for the velocity model in the inversion, which prevents them from real application. Referring to these studies, this work combines low-frequency extension commonly as well as multiscale inversion with deep learning, and proposes a multi-scale FWI gradient optimization method based on CNN. CNN we designed is trained to predict the inversion gradient corresponding to the low-frequency band data in FWI, so that multi-scale gradient optimization can be directly used in multi-scale inversion, expanding the low-frequency information in the actual data and reducing the calculation in FWI. With a specially designed dataset, CNN is trained to optimize the gradients computed from the high-frequency band data by predicting the gradients corresponding to the low-frequency band data and the gradients corresponding to the mid-frequency band data, respectively. The predicted gradients are used in different stages of the multi-scale inversion. The low-frequency gradients are used to invert the initial structural construction so as not to rely on a good initial model, and the high-frequency gradients are used to improve the accuracy of the inversion results. In this way, low-frequency expansion and multiscale inversion can be achieved simultaneously. Our method achieves good results on the initial model for a given uniform wave velocity, effectively alleviating the reliance on the initial model in FWI. This study provides a new idea of combining deep learning and full waveform inversion, which will be effectively used in seismic data processing.

How to cite: Yang, S., Jiang, P., Ren, Y., and Xu, X.: CNN-based Multi-Scale Gradient Optimization for Full-Waveform Inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-11916,, 2021.

Raffaele Bonadio, Sergei Lebedev, Thomas Meier, Pierre Arroucau, Andrew Schaeffer, Andrea Licciardi, Matthew Agius, Clare Horan, Louise Collins, Brian O'Reilly, Peter Readman, and the Ireland Array Working Group

Spatial resolution, as the ability to distinguish different features that are close together, is a fundamental concept in seismic tomography and other imaging fields. In contrast with microscopy or telescopy, seismic tomography’s images are computed, and their resolution has a complex, non-linear dependence on the data sampling and errors. Linear inverse theory provides a useful resolution-analysis approach, defining resolution in terms of the closeness of the resolution matrix to the identity matrix. This definition is similar to the universal, multi-disciplinary one in some contexts but diverges from it markedly in others. In this work, we adopt the universal definition of resolution (the minimum distance at which two spike anomalies can be resolved). The highest achievable resolution of a tomographic model then varies spatially and depends on the data sampling and errors in the data. We show that the propagation of systematic errors is resistant to data redundancy and results in models dominated by noise if the target resolution is too high. This forces one to look for smoother models and effectively limits the resolution. Here, we develop a surface-wave tomography method that finds optimal lateral resolution at every point by means of error tracking.
We first measure interstation phase velocities at simultaneously recording station pairs and compute phase-velocity maps at densely, logarithmically spaced periods. Multiple versions of the maps with varying smoothness are computed, ranging from very rough to very smooth. Phase-velocity curves extracted from the maps at every point are then inverted for shear-velocity (VS) profiles. As we show, errors in these phase-velocity curves increase nearly monotonically with the map roughness. Very smooth VS models computed from very smooth phase-velocity maps will be the most accurate, but at a cost of a loss of most structural information. At the other extreme, models that are too rough will be dominated by noise. We define the optimal resolution at a point such that the error of the local phase-velocity curve is below an empirical threshold. The error is estimated by isolating the roughness of the phase-velocity curve that cannot be explained by any Earth structure.
A 3D VS model is then computed by the inversion of the phase-velocity maps with the optimal resolution at every point. The estimated optimal resolution shows smooth lateral variations, confirming the robustness of the procedure. Importantly, optimal resolution does not scale with the density of the data coverage: some of the best-sampled locations require relatively low lateral resolution, probably due to systematic errors in the data.
We apply the method to image the lithosphere and underlying mantle beneath Ireland and Britain, using 11238 newly measured, broadband, inter-station dispersion curves. The lateral resolution of the 3D model is computed explicitly and varies from 39 km in central Ireland to over 800 km at the region boundaries, where the data coverage declines. Our tomography reveals pronounced, previously unknown variations in the lithospheric thickness beneath the region, with implications for the Caledonian assembly of the islands’ landmass and the mechanism of the magmatism of the British Tertiary Igneous Province.

How to cite: Bonadio, R., Lebedev, S., Meier, T., Arroucau, P., Schaeffer, A., Licciardi, A., Agius, M., Horan, C., Collins, L., O'Reilly, B., Readman, P., and Ireland Array Working Group, T.: Optimal resolution seismic tomography with error tracking: imaging the upper mantle beneath Ireland and Britain, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-11754,, 2021.

Tong Zhou, Min Chen, Ziyi Xi, and Jiaqi Li

Cratonic lithosphere is believed to be rigid and less deformed during a long period of time. However, the detailed structure of Cratons may bring information of the complex formation and assemblage process of the continental lithosphere. Here, we present the seismic radial anisotropic structure of the North American Craton (NAC) constrained by a regional full-waveform inversion (FWI) with 465,422 high-quality frequency-dependent travel time misfit measurements with the shortest period of 15 s from both the body wave and surface wave recordings of 5,120 stations and 160 earthquakes located in the contiguous U.S and surrounding regions. Started from an initial model constructed by combining US.2016 and Crust1.0 in the crust and S40RTS (isotropic) in the mantle, we are able to have the optimized crustal structure in terms of initial waveform similarity and get rid of existing features from other radially anisotropic mantle models.

Our new model reveals the NAC lithosphere with about +2% voigt shear wave speed anomaly and an average thickness of 200–250 km beneath the Superior Craton, and becomes thinner towards the eastern, the southern, and the southwestern margins with a thickness decreased to 100–150 km. The radial anisotropy manifests a layer of higher horizontal shear wave speed VSH (ξ=VSH2/VSV2>1) beneath the core of Superior Craton down to around 160 km, where the higher vertical shear wave speed VSV (ξ<1) is observed beneath 160 km. Such radial anisotropy layering is also observed in the margin of continental lithosphere but with shallower depth. The radial anisotropic layer matches the receiver function results of mid-lithosphere discontinuities of the Craton cores, and the lithosphere conductivity result. The radial anisotropy layering observation confirms the two-layered lithosphere structure of the NAC, where the upper layer likely represents the original radial anisotropy fabric related to the cooling of the craton core, while the lower layer might be related to the tectonic processes more recently, e.g., accretion . The lithospheric thinning beneath the NAC margins indicates the deformation of the lithosphere and is likely controlled by the large-scale mantle convection, therefore relates to the further modification process of the NAC.

How to cite: Zhou, T., Chen, M., Xi, Z., and Li, J.: Lithospheric structure of the North American Craton constrained by full waveform inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-14319,, 2021.

Hui Dou, Sergei Lebedev, Bruna Chagas de Melo, Baoshan Wang, and Weitao Wang

We present a new shear-wave velocity model of the upper mantle beneath the East Asia, ASIA2021, derived using the Automatic Multimode Inversion technique. We use waveform fits of over 1.3 million seismograms, comprising waveforms of surface waves, S and multiple S waves.  In total, data from 9351 stations and 23344 events constrain ASIA2021, which maps in detail the structure of the lithosphere and underlying mantle beneath the region. Our model reveals deep structure beneath the tectonic units that make up East Asia. It shows agreement with previous models at larger scales and, also, sharper and stronger velocity anomalies at smaller regional scales. High-velocity continent roots are mapped in detail beneath the Sichuan Basin, Tarim Basin, Ordos Block, and Siberian Craton, extending to over 200 km depths. The lack of a high-velocity continental root beneath the Eastern North China Craton (ENCC), underlain, instead, by a low-velocity anomaly, is consistent with the destruction of this Archean nucleus. Strong low-velocity anomalies are mapped within the top 100 km beneath Tibet, Pamir, Altay-Sayan area, and back-arc basins. At greater depths, ASIA2021 shows high-velocity anomalies related to the subducted and underthrusted lithosphere of India beneath Tibet and the subduction of the Pacific and other plates in the upper mantle. In the mantle transition zone (MTZ), we find high-velocity anomalies probably related to deflected subducted slabs or detached portions of ancient continent cratons. In particularly, ASIA2021 reveals separate bodies, probably originating from the Indian Plate lithosphere beneath central Tibet, with one at 100-200 km beneath Songpan-Ganzi Block (SGFB) and the other in the MTZ. A strong low-velocity anomaly extending from the surface to the lower mantle beneath Hainan volcano and South China Sea is consistent with the hypothesis of the Hainan mantle plume. The high-velocity anomaly beneath ENCC in MTZ can be interpreted as a detached Archean continent root. The Pacific Plate subducts beneath the eastern margin of Asia into the MTZ and appears to deflect and extend horizontally as far west as the Songliao Basin. The absence of major gaps in the stagnant slab is consistent with the origin of Changbaishan volcano above being related to the Big Mantle Wedge, proposed previously. The low-velocity anomalies down to ~ 700 km depth beneath the Lake Baikal area suggest a hot upwelling (mantle plume) feeding the widely distributed Cenozoic volcanoes in central and western Mongolia.

How to cite: Dou, H., Lebedev, S., Chagas de Melo, B., Wang, B., and Wang, W.: Imaging the East Asia using waveform tomography with massive datasets , EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-13665,, 2021.

Deborah Wehner, Nienke Blom, Nicholas Rawlinson, Meghan Miller, Sri Widiyantoro, and Mr Daryono

Southeast Asia is one of the most complex tectonic regions on Earth. This is mainly a result of its location within the triple junction of the Australian, Eurasian and Philippine Sea plates which has created a complicated configuration of active plate tectonic boundaries. Adjoint waveform tomography is especially suitable for imaging such complex regions. By simulating the 3D wavefield, it is possible to directly compare observed and simulated seismograms, thereby taking into account both body and surface waves. The method can account for the effects of anisotropy, anelasticity, wavefront healing, interference and (de)focusing that can hamper other seismological methods.

To date, sparse instrument coverage in the region has contributed to a heterogeneous path coverage. In this project, we make use of publicly available data as well as our recently deployed networks of broadband seismometers on Borneo and Sulawesi. This, in addition to access to national permanent networks, provides data from over 300 stations which promises a significant improvement in data coverage around the Banda Arc, Borneo and Sulawesi. We employ a geographical weighting scheme to minimise the effect of dense regional arrays and compile a catalogue of 118 well-constrained earthquakes, optimising for coverage, signal-to-noise ratio and data availability. An optimised window selection algorithm allows us to balance amplitude differences and include as much signal as possible while avoiding noisy data.

Here, we present a seismic waveform tomography for upper mantle structure in Southeast Asia, imaging radially anisotropic S velocity, P velocity and density. We use a gradient-based optimisation scheme (L-BFGS) and adjoint methods to obtain sensitivity kernels as the corresponding gradients. In the first part of the inversion, periods down to 50 s are used to update a 1D initial model, adapting a multi-scale approach in which long periods are inverted for first to avoid cycle skipping. In our long-period results, we observe a strong regional low S-velocity structure with an underlying high-velocity anomaly. The results are consistent with the global S40RTS model. 

How to cite: Wehner, D., Blom, N., Rawlinson, N., Miller, M., Widiyantoro, S., and Daryono, M.: Towards 3D Multiscale Adjoint Waveform Tomography of the Lithosphere and Underlying Mantle beneath Southeast Asia, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-1586,, 2021.

Janneke de Laat, Sergei Lebedev, Bruna Chagas de Melo, Nicolas Celli, and Raffaele Bonadio

We present a new S-wave velocity tomographic model of the Australian Plate, Aus21.  It is constrained by waveforms of 0.9 million seismograms with both the corresponding sources and stations located within the half-hemisphere centred at the Australian continent. Waveform inversion extracts structural information from surface, S- and multiple S-waves on the seismograms in the form of a set of linear equations. These equations are then combined in a large linear system and inverted jointly to obtain a tomographic model of S- and P-wave speeds and S-wave azimuthal anisotropy of the crust and upper mantle. The model has been validated by resolution tests and, for particular locations in Australia with notable differences with previous models, by independent inter-station measurements of surface-wave phase velocities, which we performed using available array data. 


Aus21 offers new insights into the structure and evolution of the Australian Plate and its boundaries. The Australian cratonic lithosphere occupies nearly all of the western and central Australia but shows substantial lateral heterogeneity. It extends up to the northern edge of the plate, where it is colliding with island arcs, without subducting. The rugged eastern boundary of the cratonic lithosphere provides a lithospheric definition of the Tasman Line. The thin, warm lithosphere below the eastern part of the continent, east of the Tasman Line, underlies the Cenozoic volcanism locations in the area. The lithosphere is also thin and warm below much of the Tasman Sea, underlying the Lord Howe hotspot and the submerged part of western Zealandia. A low velocity anomaly that may indicate the single source of the Lord Howe and Tasmanid hotspots is observed in the transition zone offshore the Australian continent, possibly also sourcing the East Australia hotspot. Another potential hotspot source is identified below the Kermadec Trench, causing an apparent slab gap in the overlying slab and possibly related to the Samoa Hotspot to the north. Below a portion of the South East Indian Ridge (the southern boundary of the Australian Plate) a pronounced high velocity anomaly is present in the 200-400 km depth range just east of the Australian-Antarctic Discordance (AAD), probably linked to the evolution of this chaotic ridge system.

How to cite: de Laat, J., Lebedev, S., Chagas de Melo, B., Celli, N., and Bonadio, R.: The Australian Plate and underlying mantle from waveform tomography with massive datasets , EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-5272,, 2021.

Maryam Rezaei and Taghi Shirzad

Classical surface wave tomography based on waveform scattering through seismic data has an important role in studying the structure of the Earth‘s crust and upper mantle on regional and global scale. The shallow crustal velocity structure is studied using earthquake waveforms in Khorasan/E-NE Iran. For this purpose, 522 local recorded waveforms with M≥4, which occurred between 53°-63°E and 30°-42°N, were selected. Therefore, all available vertical components of waveforms recorded at the stations in the Iranian Seismological Center (IrSC), the International Institute of Earthquake Engineering and Seismology (IIEES), and the IRIS global network collected in the period between January 2006 to October 2020. Then, some data selection criteria were applied for each waveform, including (i)SNR>4, (ii) the gap time less than 2 s within the expected signal window (1.5-4.5 km/s), (iii) epicentral distance greater than 20 km. The multiple-filter analysis technique was then applied by the computer program in seismology Hermann and Ammon (2013) to measure Rayleigh wave dispersion curves in the period range of 3-50 s. Finally, Rayleigh wave 2D horizontal group velocity maps are calculated by the fast marching surface wave tomography method. Our tomographic results indicate some local low velocity anomalies appeared in the W-NW of the study area where it connected with the East-Alborz tectonic structure. Also, the Doruneh Fault System clearly separates Kopeh-Dagh tectonic zone and Central Iran micro-plateau. However, a high velocity anomaly appears in Kopeh-Dagh tectonic zone at periods larger than 16 s and 30 s.

How to cite: Rezaei, M. and Shirzad, T.: Crustal structure of Khorasan (E-NE Iran) using the Rayleigh wave tomography, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-1186,, 2021.

Geothermal studies
Emma L. Chambers, Raffaele Bonadio, Sergei Lebedev, Javier Fullea, Duygu Kiyan, Christopher Bean, Brian O'Reilly, and Patrick Meere

Deep geothermal resources in low- to medium-temperature settings remain poorly understood and untapped in Ireland and much of Europe. Our new project DIG (De-risking Ireland’s Geothermal Potential) integrates multi-disciplinary, multi-scale datasets in order to investigate Ireland’s low-enthalpy geothermal energy potential. Seismic measurements constrain the distributions of seismic velocities and, through them, the composition and temperature within the lithosphere and underlying mantle. Recent deployments of broadband seismic stations and the surface-wave measurements using the new data yield an unprecedentedly dense data sampling of the crust and upper mantle beneath Ireland and neighbouring Britain. These data form a foundation for the region-scale, multi-parameter modelling of the thermal state of the lithosphere.

We use the recently assembled dataset of over 11,000 Rayleigh-wave, phase-velocity curves, measured for pairs of stations across Ireland and Britain (Bonadio et al. 2021) and complement it with new interstation measurements of Love-wave phase velocities. The measurements were performed using two methods with complementary period ranges, the teleseismic cross-correlation method and waveform inversion. Spanning a very broad period range (from as short as 4 s to as long as 500 s), the phase velocities provide resolution from the upper-middle crust to the asthenosphere. The joint analysis of Rayleigh and Love measurements constrains the isotropic-average shear-wave velocity, relatable to temperature and composition. The optimal-resolution, phase-velocity maps of Bonadio et al. (2021) for Rayleigh waves and the new maps for Love waves computed in this study provide essential constraints on the thermal structure of the region’s lithosphere. We demonstrate this by inverting the data using an integrated geophysical-petrological thermodynamically self-consistent approach (Fullea et al., 2021). The multi-parameter models produced by the integrated inversions fit the surface-wave and surface-elevation data and reveal the temperatures and geothermal gradients within the crust. 

Bonadio, R., Lebedev, S., Meier, T., Arroucau, P., Schaeffer, A. J., Licciardi, A., Agius, M. R., Horan, C., Collins, L., O'Reilly, B. M., Readman, P. and the Ireland Array Working Group (2021). Optimal resolution tomography with error tracking: imaging the upper mantle beneath Ireland and Britain. Geophys. J. Int., in revision.

Fullea, J., Lebedev, S., Martinec, Z., & Celli, N. L. (2021). WINTERC-G: mapping the upper mantle thermochemical heterogeneity from coupled geophysical-petrological inversion of seismic waveforms, heat flow, surface elevation and gravity satellite data. Geophys. J. Int., revised version under evaluation.

How to cite: Chambers, E. L., Bonadio, R., Lebedev, S., Fullea, J., Kiyan, D., Bean, C., O'Reilly, B., and Meere, P.: Imaging the temperature beneath Ireland and Britain using massive, broadband surface-wave datasets and petrological inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-5381,, 2021.

Riddhi Dave, Fiona Ann Darbyshire, Juan Carlos Afonso, and Khaled Ali

We present new thermochemical models of the lithosphere and upper mantle beneath the Superior craton and surrounding regions. The study area is dominated by the Archean Superior Province, surrounded by Proterozoic orogenic belts such as the Trans-Hudson Orogen (THO) to the north and the Grenville Orogen to the southeast. Portions of the Rae and Hearne cratons north of the THO are also studied, as is the Mid-continent Rift to the south. Over a period of ∼3 Ga, the region has seen assembly and modification by accretionary and orogenic events, periods of rifting, and the influence of a number of mantle hotspots. Here, we use a probabilistic inverse method to jointly invert Rayleigh wave dispersion data, Vp data, geoid anomalies, surface heat flow, and absolute elevation. The output is a 3D model of the seismic, temperature, bulk density, and compositional structure of the whole lithosphere beneath the Superior craton.

The resulting model will provide new opportunities for joint studies of the structure of the upper mantle and will shed light on the thermal and compositional variations beneath the region. In this presentation, we will discuss the results from our model and several robust features that carry important geological and geodynamical implications for this region.

How to cite: Dave, R., Darbyshire, F. A., Afonso, J. C., and Ali, K.: Thermochemical structure of the lithosphere and upper mantle beneath Superior craton: Results from multi-observable probabilistic inversion, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-6456,, 2021.

Cécile Ducrocq, Halldór Geirsson, Alex Hobé, Gylfi Páll Hersir, Thóra Árnadóttir, Freysteinn Sigmundsson, Ari Tryggvason, and Vincent Drouin

Crustal deformation in volcanic areas relates ground motions, measured by geodetic techniques, to physical (e.g. pressure or volumetric) changes of magmatic sources below the surface. These measurements contribute to studies of ongoing processes at the source of possible unrest, and are thus key for hazard assessment in active volcanic areas around the globe. However, such assessments often rely on geodetic-based models with quite simplistic assumptions of the physical structure of the volcanic complex. Particularly, constant values of elastic parameters (e.g. Poisson’s ratio and shear moduli) are commonly used for entire active volcanic areas, thus overlooking the spatial effects of lithological properties, depth-dependant compression and temperature variations on those parameters. These simplifications may lead to inaccurate interpretation of the location, shape, and volume change of deformation sources.


In this study we ask how the 3-D heterogeneities of the elastic crustal structure beneath the Hengill volcanic system, SW Iceland, affects models of deformation sources in the area. The Hengill area hosts two active volcanic systems (Hengill and Hrómundartindur), and two high-enthalpy geothermal power plants (Nesjavellir and Hellisheiði), which provide thermal and electrical power to Reykjavík, the capital of Iceland, only 30 km away. To retrieve information on the spatial heterogeneities in the shear moduli and Poisson’s ratio beneath the Hengill area, we first estimate the 3-D shallow density structure of the area using results from regional and local gravimetric surveys. We implement this structure, along with seismic tomographic studies of the SW Iceland, in a Finite Element Model to solve, using forward models, for the 3-D heterogeneities in the shear moduli and Poisson’s ratio beneath the Hengill area. Furthermore, we discuss the difference between static and kinematic elastic moduli, which is important when building deformation models from seismic tomography. The new 3-D inferred elastic model is then used to re-estimate parameters for different sources of deformation causing uplift and subsidence in the area in the past decades. This study shows the importance of accounting for heterogeneities in the crustal elastic structure to estimate with higher accuracy the sources of deformation in volcanic areas around the world.

How to cite: Ducrocq, C., Geirsson, H., Hobé, A., Hersir, G. P., Árnadóttir, T., Sigmundsson, F., Tryggvason, A., and Drouin, V.: 3-D Heterogeneous Elastic Crustal Structure for Deformation Models in the Hengill Area, SW Iceland, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-7553,, 2021.

Bing Xia, Irina Artemieva, and Hans Thybo

We present a thermal model for the lithosphere in Tibet and adjacent regions based on the new thermal isostasy method and our compilation of the Moho depth based on published seismic models. The predicted surface heat flow is in agreement with the few available, reliable borehole measurements. Cratonic-type cold and thick lithosphere (200-240 km) with a surface heat flow of 40-50 mW/m2 typifies the Tarim craton, the north-western Yangtze craton, and most of the Lhasa Block that is possibly refrigerated by underthrusting Indian lithosphere. The thick lithosphere of the Lhasa block extends further north in its western and eastern segments than in its central section. We identify a North Tibet anomaly with a thin (<80 km) lithosphere and high surface heat flow (>80-100 mW/m2), possibly associated with the removal of lithospheric mantle and asthenospheric upwelling. Other parts of Tibet have an intermediate lithosphere thickness of 120-160 km and a surface heat flow of 45-60 mW/m2, with a patchy style in eastern Tibet. In the Qaidam deep sedimentary basin the lithosphere is about 100-120 km thick. The heterogeneous thermal lithosphere beneath Tibet suggests an interplay of several mechanisms as the driver of the observed uplift.

How to cite: Xia, B., Artemieva, I., and Thybo, H.: Extreme variability of Tibetan thermal lithosphere , EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-14135,, 2021.

Azimuthal anisotropy
Hamzeh Mohammadigheymasi, Mohammad Reza Ebrahimi, Graça Silveira, and David schlaphorst

Shear wave splitting analysis is a frequently used tool to study elastic anisotropy from the lower mantle to the crust. Several methods have been developed to evaluate the splitting parameters, Φ (fast axis) and δt (delay time), including the correlation of wave components, minimization of covariance matrix eigenvalues, and minimizing energy on the transverse component. Despite massive progress in introducing sophisticated methods, still fundamental problems, related mainly to noisy data, interfering phases, length of the analyzed waveform, and stability and reliability of results, remain. This study presents a sparsity-based adaptive filtering method to magnify the SKS waveforms and suppress the unwanted noise and interfering phases. The study is an extension of Jurkevics (1988), computing the semi-minor and semi-minor axis of the polarized motion in the time-frequency domain using a regularized inversion-based approach imposing a sparsity constraint. Afterward, the elliptical particle motion caused by the split shear waves and correspond to high semi-minor amplitude is derived in the time-frequency domain. The information is used to design an adaptive filter in the time domain to amplify the SKS phase and suppress the noise and other phases having non-elliptical polarization. The regularized inversion-based approach enables obtaining a sparse time-frequency semi-minor map while handling noise problems in the time-frequency decomposition. Conducting synthetic simulations, we show that the proposed method increases the signal-to-noise ratio of the SKS phase in radial and transverse components, giving a better estimation of anisotropy parameters in the presence of noise and other interfering phases. Future work involves implementing the processing algorithm on real data recorded in São Tomé and Prı́ncipe, Madeira, and Canary islands. This research contributes to the FCT-funded SHAZAM (Ref. PTDC/CTA-GEO/31475/2017) and SIGHT (Ref. PTDC/CTA-GEF/30264/2017) projects.

How to cite: Mohammadigheymasi, H., Ebrahimi, M. R., Silveira, G., and schlaphorst, D.: A sparsity-based adaptive filtering approach to Shear Wave Splitting, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-13988,, 2021.