Sources and Receivers

Instaseis uses objects to represent seismic receivers and different types of sources. Source objects are standard moment tensor or double couple sources, ForceSource object are (as the name implies) force sources, and finite source are represented by the FiniteSource class.

Note

Instaseis works with geocentric coordinates. Coordinates in for example QuakeML files are defined on the WGS84 ellipsoid. Thus they need to be converted. Instaseis performs this conversion if it is used to read from any file format - if you create the source objects yourself you have to take care to pass geocentric coordinates. The difference in both can be up to 20 kilometers. Depending on the application, this is not an effect you can savely ignore.

The ruleset within Instaseis is simple:

  1. Source and receiver objects are always in geocentric coordinates.

  2. If parsed from any outside file (QuakeML/SEED/SAC/StationXML/…) the coordinates are assumed to be WGS84 and will be converted so the source/receiver objects are again in geocentric coordinates.

The directions of r, theta, and phi are defined according to the standard spherical coordinate system definition used in physics, which is: r positive outward, theta positive downward, and phi positive counter-clockwise.

Receiver

class instaseis.source.Receiver(latitude, longitude, network=None, station=None, location=None, depth_in_m=None)[source]

Class dealing with seismic receivers.

Parameters
  • latitude (float) – The geocentric latitude of the receiver in degree.

  • longitude (float) – The longitude of the receiver in degree.

  • depth_in_m (float) – The depth of the receiver in meters. Only

  • network (str, optional) – The network id of the receiver.

  • station (str, optional) – The station id of the receiver.

  • location (str) – The location code of the receiver.

>>> from instaseis import Receiver
>>> rec = Receiver(latitude=12.34, longitude=56.78, network="AB",
...                station="CDE", location="SY")
>>> print(rec)  
Instaseis Receiver:
    Longitude :   56.8 deg
    Latitude  :   12.3 deg
    Network   : AB
    Station   : CDE
    Location  : SY
static parse(filename_or_obj, network_code=None)[source]

Attempts to parse anything to a list of Receiver objects. Always returns a list, even if it only contains a single element. It is meant as a single entry point for receiver information from any source.

Supports StationXML, the custom STATIONS fileformat, SAC files, SEED files, and a number of ObsPy objects. This method can furthermore work with anything ObsPy can deal with (filename, URL, memory files, …).

Coordinates are assumed to be defined on the WGS84 ellipsoid and will be converted to geocentric coordinates.

Parameters
  • filename_or_obj – Filename/URL/Python object

  • network_code – Network code needed to parse ObsPy station objects. Likely only needed for the recursive part of this method.

Returns

List of Receiver objects.

The following example parses a StationXML file to a list of Receiver objects.

>>> import instaseis
>>> print(instaseis.Receiver.parse(stationxml_file))
[<instaseis.source.Receiver object at 0x...>]

Source

class instaseis.source.Source(latitude, longitude, depth_in_m=None, m_rr=0.0, m_tt=0.0, m_pp=0.0, m_rt=0.0, m_rp=0.0, m_tp=0.0, time_shift=None, sliprate=None, dt=None, origin_time=UTCDateTime(1970, 1, 1, 0, 0))[source]

Class to handle a seismic moment tensor source including a source time function.

Parameters
  • latitude – geocentric latitude of the source in degree

  • longitude – longitude of the source in degree

  • depth_in_m – source depth in m

  • m_rr – moment tensor components in r, theta, phi in Nm

  • m_tt – moment tensor components in r, theta, phi in Nm

  • m_pp – moment tensor components in r, theta, phi in Nm

  • m_rt – moment tensor components in r, theta, phi in Nm

  • m_rp – moment tensor components in r, theta, phi in Nm

  • m_tp – moment tensor components in r, theta, phi in Nm

  • time_shift – correction of the origin time in seconds. Useful in the context of finite source or user defined source time functions.

  • sliprate – normalized source time function (sliprate)

  • dt – sampling of the source time function

  • origin_time – The origin time of the source. If you don’t reconvolve with another source time function this time is the peak of the source time function used to generate the database. If you reconvolve with another source time function this time is the time of the first sample of the final seismogram.

>>> import instaseis
>>> source = instaseis.Source(
...     latitude=89.91, longitude=0.0, depth_in_m=12000,
...     m_rr = 4.71e+17, m_tt = 3.81e+15, m_pp =-4.74e+17,
...     m_rt = 3.99e+16, m_rp =-8.05e+16, m_tp =-1.23e+17)
>>> print(source)  
Instaseis Source:
    Origin Time      : 1970-01-01T00:00:00.000000Z
    Longitude        :    0.0 deg
    Latitude         :   89.9 deg
    Depth            : 1.2e+01 km km
    Moment Magnitude :   5.73
    Scalar Moment    :   4.96e+17 Nm
    Mrr              :   4.71e+17 Nm
    Mtt              :   3.81e+15 Nm
    Mpp              :  -4.74e+17 Nm
    Mrt              :   3.99e+16 Nm
    Mrp              :  -8.05e+16 Nm
    Mtp              :  -1.23e+17 Nm
property M0

Scalar Moment M0 in Nm

classmethod from_strike_dip_rake(latitude, longitude, depth_in_m, strike, dip, rake, M0, time_shift=None, sliprate=None, dt=None, origin_time=UTCDateTime(1970, 1, 1, 0, 0))[source]

Initialize a source object from a shear source parameterized by strike, dip and rake.

Parameters
  • latitude – geocentric latitude of the source in degree

  • longitude – longitude of the source in degree

  • depth_in_m – source depth in m

  • strike – strike of the fault in degree

  • dip – dip of the fault in degree

  • rake – rake of the fault in degree

  • M0 – scalar moment

  • time_shift – correction of the origin time in seconds. only useful in the context of finite sources

  • sliprate – normalized source time function (sliprate)

  • dt – sampling of the source time function

  • origin_time – The origin time of the source. If you don’t reconvolve with another source time function this time is the peak of the source time function used to generate the database. If you reconvolve with another source time function this time is the time of the first sample of the final seismogram.

>>> import instaseis
>>> source = instaseis.Source.from_strike_dip_rake(
...     latitude=10.0, longitude=12.0, depth_in_m=1000, strike=79,
...     dip=10, rake=20, M0=1E17)
>>> print(source)  
Instaseis Source:
    Origin Time      : 1970-01-01T00:00:00.000000Z
    Longitude        :   12.0 deg
    Latitude         :   10.0 deg
    Depth            : 1.0e+00 km km
    Moment Magnitude :   5.27
    Scalar Moment    :   1.00e+17 Nm
    Mrr              :   1.17e+16 Nm
    Mtt              :  -1.74e+16 Nm
    Mpp              :   5.69e+15 Nm
    Mrt              :  -4.92e+16 Nm
    Mrp              :   8.47e+16 Nm
    Mtp              :   1.29e+16 Nm
property moment_magnitude

Moment magnitude M_w

static parse(filename_or_obj)[source]

Attempts to parse anything to a Source object. Can currently read anything ObsPy can read, ObsPy event related objects.

For anything ObsPy related, it must contain a full moment tensor, otherwise it will raise an error.

Coordinates are assumed to be defined on the WGS84 ellipsoid and will be converted to geocentric coordinates.

Parameters

filename_or_obj – The object or filename to parse.

The following example will read a local QuakeML file and return a Source object.

>>> import instaseis
>>> source = instaseis.Source.parse(quakeml_file)
>>> print(source)  
Instaseis Source:
    Origin Time      : 2010-03-24T14:11:31.000000Z
    Longitude        :   40.1 deg
    Latitude         :   38.6 deg
    Depth            : 4.5e+00 km km
    Moment Magnitude :   5.09
    Scalar Moment    :   5.37e+16 Nm
    Mrr              :   5.47e+15 Nm
    Mtt              :  -4.11e+16 Nm
    Mpp              :   3.56e+16 Nm
    Mrt              :   2.26e+16 Nm
    Mrp              :  -2.25e+16 Nm
    Mtp              :   1.92e+16 Nm
property tensor

List of moment tensor components in r, theta, phi coordinates: [m_rr, m_tt, m_pp, m_rt, m_rp, m_tp]

property tensor_voigt

List of moment tensor components in theta, phi, r coordinates in Voigt notation: [m_tt, m_pp, m_rr, m_rp, m_rt, m_tp]


ForceSource

class instaseis.source.ForceSource(latitude, longitude, depth_in_m=None, f_r=0.0, f_t=0.0, f_p=0.0, origin_time=UTCDateTime(1970, 1, 1, 0, 0), sliprate=None, time_shift=None, dt=None)[source]

Class to handle a seismic force source.

Parameters
  • latitude – geocentric latitude of the source in degree

  • longitude – longitude of the source in degree

  • depth_in_m – source depth in m

  • f_r – force components in r, theta, phi in N

  • f_t – force components in r, theta, phi in N

  • f_p – force components in r, theta, phi in N

  • origin_time – The origin time of the source. If you don’t reconvolve with another source time function this time is the peak of the source time function used to generate the database. If you reconvolve with another source time function this time is the time of the first sample of the final seismogram.

  • sliprate – normalized source time function (sliprate)

>>> import instaseis
>>> source = instaseis.ForceSource(latitude=12.34, longitude=12.34,
...                                f_r=1E10)
>>> print(source)  
Instaseis Force Source:
    Origin Time      : 1970-01-01T00:00:00.000000Z
    Longitude :   12.3 deg
    Latitude  :   12.3 deg
    Fr        :   1.00e+10 N
    Ft        :   0.00e+00 N
    Fp        :   0.00e+00 N
Parameters
  • latitude – geocentric latitude of the source in degree

  • longitude – longitude of the source in degree

  • depth_in_m – source depth in m

  • f_r – force components in r, theta, phi in N

  • f_t – force components in r, theta, phi in N

  • f_p – force components in r, theta, phi in N

  • origin_time – The origin time of the source. If you don’t reconvolve with another source time function this time is the peak of the source time function used to generate the database. If you reconvolve with another source time function this time is the time of the first sample of the final seismogram.

  • sliprate – normalized source time function (sliprate)

  • time_shift – correction of the origin time in seconds. Useful in the context of finite source or user defined source time functions.

  • dt – sampling of the source time function

property force_rtp

List of force components in r, theta, phi, coordinates: [f_r, f_t, f_p]

property force_tpr

List of force components in theta, phi, r coordinates: [f_t, f_p, f_r]


FiniteSource

class instaseis.source.FiniteSource(pointsources=None, CMT=None, magnitude=None, event_duration=None, hypocenter_longitude=None, hypocenter_latitude=None, hypocenter_depth_in_m=None)[source]

A class to handle finite sources represented by a number of point sources.

Parameters
  • pointsources (list of Source objects) – The points sources making up the finite source.

  • CMT (Source, optional) – The centroid of the finite source.

  • magnitude (float, optional) – The total moment magnitude of the source.

  • event_duration (float, optional) – The event duration in seconds.

  • hypocenter_longitude (float, optional) – The hypocentral longitude.

  • hypocenter_latitude (float, optional) – The hypocentral latitude.

  • hypocenter_depth_in_m (float, optional) – The hypocentral depth in m.

property M0

Scalar Moment M0 in Nm

compute_centroid(planet_radius=6371000.0, dt=None, nsamp=None)[source]

computes the centroid moment tensor by summing over all pointsource weihted by their scalar moment

find_hypocenter()[source]

Finds the hypo- and epicenter based on the point source that has the smallest timeshift

classmethod from_Haskell(latitude, longitude, depth_in_m, strike, dip, rake, M0, fault_length, fault_width, rupture_velocity, nl=100, nw=1, trise=1.0, tfall=None, dt=0.1, planet_radius=6371000.0, origin_time=UTCDateTime(1970, 1, 1, 0, 0))[source]

Initialize a source object from a shear source parameterized by strike, dip and rake.

Parameters
  • latitude – geocentric latitude of the source centroid in degree

  • longitude – longitude of the source centroid in degree

  • depth_in_m – source centroid depth in m

  • strike – strike of the fault in degree

  • dip – dip of the fault in degree

  • rake – rake of the fault in degree

  • M0 – scalar moment

  • fault_length – fault length in m

  • fault_width – fault width in m

  • rupture_velocity – rupture velocity in m / s. Use negative value to propagate the rupture in negative rake direction.

  • nl – number of point sources along strike direction

  • nw – number of point sources perpendicular strike direction

  • trise – rise time

  • tfall – fall time

  • dt – sampling of the source time function

  • planet_radius – radius of the planet, default to Earth.

  • origin_time – The origin time of the first patch breaking.

classmethod from_srf_file(filename, normalize=False)[source]

Initialize a finite source object from a ‘standard rupture format’ (.srf) file

Coordinates are assumed to be defined on the WGS84 ellipsoid and will be converted to geocentric coordinates.

Parameters
  • filename (str) – path to the .srf file

  • normalize (bool, optional) – normalize the sliprate to 1

>>> import instaseis
>>> source = instaseis.FiniteSource.from_srf_file(srf_file)
>>> print(source)  
Instaseis Finite Source:
    Moment Magnitude     : 7.60
    Scalar Moment        :   3.20e+20 Nm
    #Point Sources       : 10
    Rupture Duration     :  222.2 s
    Time Shift           :    0.0 s
    Min Depth            : 50000.0 m
    Max Depth            : 50000.0 m
    Hypocenter Depth     : 50000.0 m
    Min Latitude         :    0.0 deg
    Max Latitude         :    0.0 deg
    Hypocenter Latitude  :    0.0 deg
    Min Longitude        :    0.0 deg
    Max Longitude        :    9.0 deg
    Hypocenter Longitude :    0.0 deg
classmethod from_usgs_param_file(filename_or_obj, npts=10000, dt=0.1, trise_min=1.0)[source]

Initialize a finite source object from a (.param) file available from the USGS website

Coordinates are assumed to be defined on the WGS84 ellipsoid and will be converted to geocentric coordinates.

Parameters
  • filename_or_obj (str) – path to the .param file

  • npts (int) – number of samples for the source time functions. Should be enough to accommodate long rise times plus optional filter responses.

  • dt (float) – sample interval for the source time functions. Needs to be small enough to sample short rise times and can then be low pass filtered and downsampled before extracting seismograms.

  • trise_min (float) – minimum rise time. If rise time or fall time is smaller than this value, it is replaced by the minimum value. Mostly meant to handle zero rise times in some USGS files.

>>> import instaseis
>>> source = instaseis.FiniteSource.from_usgs_param_file(param_file)
>>> print(source)  
Instaseis Finite Source:
    Moment Magnitude     : 7.87
    Scalar Moment        :   8.06e+20 Nm
    #Point Sources       : 121
    Rupture Duration     :  107.6 s
    Time Shift           :    0.0 s
    Min Depth            :  860.1 m
    Max Depth            : 26907.3 m
    Hypocenter Depth     : 26907.3 m
    Min Latitude         :   26.7 deg
    Max Latitude         :   28.7 deg
    Hypocenter Latitude  :   27.9 deg
    Min Longitude        :   84.1 deg
    Max Longitude        :   86.6 deg
    Hypocenter Longitude :   84.8 deg
property moment_magnitude

Moment magnitude M_w

next()[source]

For Py2K compat.

normalize_sliprate()[source]

normalize the sliprate using trapezoidal rule

resample_sliprate(dt, nsamp)[source]

For convolution, the sliprate is needed at the sampling of the fields in the database. This function resamples the sliprate using linear interpolation for all pointsources in the finite source.

Parameters
  • dt – desired sampling

  • nsamp – desired number of samples

set_sliprate_dirac(dt, nsamp)[source]
Parameters
  • dt – desired sampling

  • nsamp – desired number of samples

set_sliprate_lp(dt, nsamp, freq, corners=4, zerophase=False)[source]
Parameters
  • dt – desired sampling

  • nsamp – desired number of samples