3.1. Example usage of SeisSrcMoment to calculate \(M_W\) in the time domain

This example is for an icequake at Vatnajokull, Iceland. The moment magnitude, \(M_W\), is calculated in the time domain, i.e. the long period spectral level is calculated from the integral of the displacement through time.

3.1.1. 1. Specify parameters to use:

import numpy as np
from SeisSrcMoment import moment

# Specify variables:
stations_to_calculate_moment_for = ["SKR01", "SKR02", "SKR03", "SKR04", "SKR05", "SKR06", "SKR07"]
stations_not_to_process = ["SKG08", "SKG09", "SKG10", "SKG11", "SKG12", "SKG13", "GR01", "GR02", "GR03","GR04","BARD"]
mseed_filename = "data/mseed_data/20140629184210331.m"
inventory_fname = None
instruments_gain_filename = "data/instrument_gain_data.txt" # File with instrument name, instrument gains (Z,N,E) and digitaliser gains (Z,N,E)
NLLoc_event_hyp_filename = "data/NLLoc_data/loc.Tom__RunNLLoc000.20140629.184210.grid0.loc.hyp"
window_before_after = [0.004, 0.196] # The time before and after the phase pick to use for calculating the magnitude within
filt_freqs = []
# MT_data_filename = "data/MT_data/20140629184210363MT.mat"
MT_six_tensor = np.array([1.,1.,1.,0.,0.,0.]) # Explosion in this example.
density = 917. # Density of medium, in kg/m3
Vp = 3630. # P-wave velocity in m/s
Q = 150. # Quality factor for the medium
use_full_spectral_method = False # If False, will use the time domain displacement integration method
verbosity_level = 1 # Verbosity level (1 for moment only) (2 for major parameters) (3 for plotting of traces)

3.1.2. Run moment calculation:

# Find seismic moment release:
av_M_0, std_err_seis_M_0, n_obs, event_obs_dict = moment.calc_moment(mseed_filename, NLLoc_event_hyp_filename, stations_to_calculate_moment_for, density, Vp, inventory_fname=inventory_fname, instruments_gain_filename=instruments_gain_filename, Q=Q, window_before_after=window_before_after, filt_freqs=filt_freqs, use_full_spectral_method=use_full_spectral_method, stations_not_to_process=stations_not_to_process, MT_six_tensor=MT_six_tensor, verbosity_level=verbosity_level)
print("Seismic moment release (Nm):", av_M_0)

Processing data for station: SKR01
Overall seismic moment (Nm): 7860258.2629
Processing data for station: SKR02
Overall seismic moment (Nm): 7596240.15988
Processing data for station: SKR03
Overall seismic moment (Nm): 5185040.49289
Processing data for station: SKR04
Overall seismic moment (Nm): 3427612.88854
Processing data for station: SKR05
Overall seismic moment (Nm): 4675845.68622
Processing data for station: SKR06
Overall seismic moment (Nm): 3582958.1162
Processing data for station: SKR07
Overall seismic moment (Nm): 7224947.64632
Average seismic moment for event: 5650414.75042 +/- 662869.178247
Seismic moment release (Nm): 5650414.75042

# And find corresponding moment magnitude, M_w (Hanks and Kanamori 1979):
M_w = (2./3.)*np.log10(av_M_0) - 6.0
print("Local moment magnitude, M:", M_w)

Local moment magnitude, M: -1.49861311535