pysep.utils.mt
#
Moment Tensor related functions for grabbing moment tensors, appending them to existing event objects, and writing them out in specific formats such as the CMTSOLUTION format required by SPECFEM
Module Contents#
Functions#
|
Return the seismic moment based on a moment tensor. |
|
Return the moment magitude, M_w, based on a seismic moment. Equation from |
|
Empirical formula for half duration used by Harvard CMT, stated in |
|
Transform moment tensor between XYZ and RTP coordinates. Used primarily |
|
Query online GCMT moment tensor catalog via URL access for moment tensor |
|
Query FDSN webservices USGS client for moment tensors using the current |
|
Focal mechanisms created by John Ristau are written to a .csv file |
|
Get moment tensor information from a internal csv file, |
|
Attempt to find focal mechanism information with a given ObsPy Event object. |
- pysep.utils.mt.seismic_moment(mt)[source]#
Return the seismic moment based on a moment tensor. Can take a list of tensor components, or a Tensor object from ObsPy.
- Parameters:
mt (list of floats or obspy.core.event.source.Tensor) – the components of the moment tensor M_ij
- Return type:
float
- Returns:
the seismic moment, in units of N*m
- pysep.utils.mt.moment_magnitude(moment, c=10.7)[source]#
Return the moment magitude, M_w, based on a seismic moment. Equation from Hanks & Kanamori (1979)
- Parameters:
c (float) – correction factor for conversion, 10.7 for units of N*m, 16.1 for units of dyne*cm
moment (float) – the seismic moment, in units of N*m
- Return type:
float
- Returns:
moment magnitude, M_w
- pysep.utils.mt.half_duration_from_m0(moment)[source]#
Empirical formula for half duration used by Harvard CMT, stated in Dahlen and Tromp (1998, p.178).
- Parameters:
moment (float) – seismic moment in N*m
- Return type:
float
- Returns:
empirically scaled half duration in unit seconds
- pysep.utils.mt.mt_transform(mt, method)[source]#
Transform moment tensor between XYZ and RTP coordinates. Used primarily to transform GeoNet (John Ristau) moment tensors into the correct coordinate system for use in SPECFEM
Based on Equation from ‘Aki and Richards: Quantitative Seismology’ book
Note
Acceptable formats for the parameter mt:
[m11, m22, m33, m12, m13, m23]
[mxx, myy, mzz, mxy, mxz, myz]
[mrr, mtt, mpp, mrt, mrp, mtp]
- Parameters:
mt (dict) – moment tensor in format above
method (str) – type of conversion, “rtp2xyz” or “xyz2rtp”
- Return type:
dict
- Returns:
converted moment tensor dictionary
- pysep.utils.mt.get_gcmt_moment_tensor(event=None, origintime=None, magnitude=None, time_wiggle_sec=120, magnitude_wiggle=0.5)[source]#
Query online GCMT moment tensor catalog via URL access for moment tensor components of a given event. Searches based on origin time and magnitude of an event with a given amount of wiggle room for catalog mismatch of origin time and magnitude.
Note
input is either event OR origintime AND magnitude
- Parameters:
event (obspy.core.event.Event) – Event to use to query for moment tensor
origintime (UTCDateTime or str) – event origin time
magnitude (float) – centroid moment magnitude for event lookup
time_wiggle_sec (int) – padding on catalog filtering criteria realted to event origin time
magnitude_wiggle (float) – padding on catalog filter for magnitude
- Return type:
obspy.core.event.Event
- Returns:
event object for given earthquake
- pysep.utils.mt.get_usgs_moment_tensor(event, time_wiggle_sec=120.0, magnitude_wiggle=0.5, latitude_wiggle_deg=1.0, longitude_wiggle_deg=1.0, depth_wiggle_km=5.0, **kwargs)[source]#
Query FDSN webservices USGS client for moment tensors using the current event definition, which may or may not have been collected via USGS.
Kwargs passed to Client.get_events() for additional event constraint pars
- Parameters:
event (obspy.core.event.Event) – Event to use to query for moment tensor
time_wiggle_sec (float) – padding on catalog filtering criteria realted to event origin time
magnitude_wiggle (float) – +/- padding on magnitude search
latitude_wiggle_deg (float) – +/- padding on latitude search
longitude_wiggle_deg (float) – +/- padding on longitude search
depth_wiggle_km (float) – +/- padding on depth search
- Return type:
obspy.core.event.Event
- Returns:
event object for given earthquake
- pysep.utils.mt.get_geonet_mt(event_id, units, csv_fid=None)[source]#
Focal mechanisms created by John Ristau are written to a .csv file located on Github. This function will append information from the .csv file onto the Obspy event object so that all the information can be located in a single object
- Parameters:
event_id (str) – unique event identifier
units (str) – output units of the focal mechanism, either: ‘dynecm’: for dyne*cm or ‘nm’: for Newton*meter
csv_fid (str) – optional local path to .csv file containing Ristau catalog of moment tensors. If not given, will search a default online URL where the catalog is assumed to be stored
- Rtype focal_mechanism:
obspy.core.event.FocalMechanism
- Return focal_mechanism:
generated focal mechanism
- pysep.utils.mt.query_geonet_mt_catalog(event_id, csv_fid=None)[source]#
Get moment tensor information from a internal csv file, or from an external github repository query. Only relevant to the new zealand tomography problem. Geonet moment tensors stored with a specific column format.
- Parameters:
event_id (str) – unique event identifier
csv_fid (str) – optional path to GeoNet CMT solution file that is stored locally on disk, will be accessed before querying web service
- Rtype moment_tensor:
dict
- Return moment_tensor:
dictionary created from rows of csv file
- pysep.utils.mt.append_focal_mechanism_to_event(event, method='all', overwrite_focmec=False, overwrite_event=False, client=None)[source]#
Attempt to find focal mechanism information with a given ObsPy Event object.
Note
FDSN fetched events are devoid of a few bits of information that are useful for our applications, e.g. moment tensor, focal mechanisms. This function will perform the conversions and append the necessary information to the event located in the dataset.
- Parameters:
event (obspy.core.event.Event) – Event object to append a focal mechanism to.
method (bool) – try to find correspondig focal mechanism using various public catalogs. Currently available: ‘all’: Try all available options in order until MT is found ‘USGS’: Search the USGS moment tensor catalog ‘GCMT’: Search the GCMT moment tensor catalog False: Don’t attempt to search for moment tensors
client (str) – Specific `client`s come built-in with specific MT catalogs If matching client, will ignore other MT choices: ‘GEONET’: will search John Ristau catalog for moment tensors,
overwrite_focmec (bool) – If the event already has a focal mechanism, overwrite the existing focal mechanism.
overwrite_event (bool) – A new event object is usually retrieved when gathering MT from USGS or GCMT. Often the locations/timing of this event are less accurate than the input event (which is usually sourced from a regional catalog). This parameter controls which event object is taken. If True, takes the USGS or GCMT catalog information, if False only takes the focal mechanism attribute.
- Rtype event:
obspy.core.event.Event
- Return event:
Event with a new focal mechanism if one was found
- Raises:
TypeError – if event is not provided as an obspy.core.event.Event