This tutorial is presented as Python code running inside a Jupyter Notebook, the recommended way to use Salvus. To run it yourself you can copy/type each individual cell or directly download the full notebook, including all required files.

Data Adaptive Mesh Masking

In many seismologic applications, the source reveiver distribution does not follow a simple geometry such as a square or a circle. Due to the fact that sources are mostly confined to plate boundaries and stations are almost always on land, very complicated domain shapes are common. See [1] for an application.

In this tutorial, we show two different ways of adapting a mesh to the source receiver distribution my removing elements from a larger mesh that are not passed by waves of interest. The example we use here is a quake in Turkey recorded on the Search Results USArray Reference Network (_US-REF).

%matplotlib inline
%config Completer.use_jedi = False

from salvus.mesh.simple_mesh import SmoothieSEM, Globe3D
from salvus import namespace as sn

from obspy.clients.fdsn import Client

# a quake in Turkey that we will use in this tutorial, original data from IRIS spud:
source = sn.simple_config.source.seismology.SideSetMomentTensorPoint3D(

# get USarray stations from iris
inv = Client("IRIS").get_stations(
    network="_US-REF", level="station", format="text"

receivers = sn.simple_config.receiver.seismology.parse(
    inv, dimensions=3, fields=["displacement"]

# prepare an event collection that will later be used to mask the mesh to a region of interest
event_collection = sn.EventCollection.from_sources(
    sources=[source], receivers=receivers

Method 1: Surface Based

In the first approach, we only mask out element based on there lateral position on the sphere and ignore the depth (this approach is used in [1]). To achieve this, a covex hull is built from the sources and receivers and all elements within that hull as well as those within a specified distance are retained in the mesh.

from salvus.mesh.mask_generators import SurfaceMaskGenerator

sm = Globe3D()
sm.basic.model = "prem_iso_one_crust"
sm.basic.min_period_in_seconds = 100.0
sm.basic.elements_per_wavelength = 2.0
sm.spherical.min_radius = 4000.0

# use event collection to create a surface mask
smg = SurfaceMaskGenerator(

# hand over the mask as a callback funcion
<salvus.mesh.unstructured_mesh.UnstructuredMesh at 0x7f2f25e8f810>

Alternatively, we can use rays computed with obspy.TauPy to mask out all elements beyond a specified distance from the rays. E.g. for teleseismic body wave tomography, this reduces the number of elements further in comparison to the surface based masking. This may be particularly interesting in the context of SmoothieSEM meshes, however, this works identically with standard cubed sphere type meshes.

In this particular example we only use the P-waves, as would be applicable in a teleseismic P-wave tomography:

from salvus.mesh.mask_generators import RayMaskGenerator

sm = SmoothieSEM()
sm.basic.model = "prem_iso_one_crust"
sm.basic.min_period_in_seconds = 50.0
sm.basic.elements_per_wavelength = 2.0
sm.basic.number_of_lateral_elements = 12
sm.spherical.min_radius = 1500.0

# no refinement, rotation to a source in Turkey
sm.basic.number_of_lateral_elements = 10
sm.source.latitude = source._initial_arguments["latitude"]
sm.source.longitude = source._initial_arguments["longitude"]

# use event collection to create a ray mask
rmg = RayMaskGenerator(