Theoretical Background¶
Overview¶
source_spec
inverts the P- or S-wave displacement amplitude spectra from
station recordings of a single event.
For each station, the code computes P- or S-wave displacement amplitude spectra for each component (e.g., Z, N, E), within a predefined time window.
The same thing is done for a noise time window: noise spectrum is used to compute spectral signal-to-noise ratio (and possibily reject low SNR spectra) and, optionnaly, to weight the spectral inversion.
The component spectra are combined through the root sum of squares:
where \(f\) is the frequency and \(S_x(f)\) is the P- or S-wave spectrum for component \(x\).
Spectral model¶
The Fourier amplitude spectrum of the S-wave displacement in far field can be modelled as the product of a source term [Brune, 1970] and a propagation term (geometric and anelastic attenuation of body waves):
where:
\(\mathcal{G}(r)\) is the geometrical spreading coefficient (see below) and \(r\) is the hypocentral distance;
\(R_{\Theta\Phi}\) is the radiation pattern coefficient for P- or S-waves (average or computed from focal mechanism, if available);
\(\rho_h\) and \(\rho_r\) are the medium densities at the hypocenter and at the receiver, respectively;
\(c_h\) and \(c_r\) are the P- or S-wave velocities at the hypocenter and at the receiver, respectively;
\(M_O\) is the seismic moment;
\(f\) is the frequency;
\(f_c\) is the corner frequency;
\(t^*\) is an attenuation parameter which includes anelastic path attenuation (quality factor) and station-specific effects.
Geometrical spreading¶
The geometrical spreading coefficient \(\mathcal{G}(r)\) can be defined in one of the following ways:
\(\mathcal{G}(r) = r^n\): \(n\) can be any positive number. \(n=1\) (default value) is the theoretical value for a body wave in a homogeneous full-space; \(n=0.5\) is the theoretical value for a surface wave in a homogeneous half-space.
Follwing Boatwright et al. [2002], eq. 8:
body wave spreading (\(\mathcal{G}(r) = r\)) for hypocentral distances below a cutoff distance \(r_0\);
frequency dependent spreading for hypocentral distances above the cutoff distance \(r_0\).
More precisely, the expression derived from Boatwright et al. [2002] is:
with
Note that here we use the square root of eq. 8 in Boatwright et al. (2002), since we correct the spectral amplitude and not the energy.
Building spectra¶
In source_spec
, the observed spectrum of component \(x\),
\(S_x(f)\) is converted into moment magnitude units \(M_w\).
The first step is to multiply the spectrum for the geometrical spreading coefficient and convert it to seismic moment units:
Then the spectrum is converted in units of magnitude (the \(Y_x (f)\) vector used in the inversion):
The data vector is compared to the teoretical model:
Finally coming to the following model used for the inversion:
where \(M_w \equiv \frac{2}{3} (\log_{10} M_0 - 9.1)\).
Inversion procedure¶
The parameters to determine are \(M_w\), \(f_c\) and \(t^*\).
The inversion is performed in moment magnitude \(M_w\) units (logarithmic amplitude). Different inversion algorithms can be used:
TNC: truncated Newton algorithm (with bounds)
LM: Levenberg-Marquardt algorithm (warning: Trust Region Reflective algorithm will be used instead if bounds are provided)
GS: grid search
IS: importance sampling of misfit grid, using k-d tree
Starting from the inverted parameters \(M_0\) ( \(M_w\) ), \(fc\), \(t^*\) and following the equations in Madariaga [2011], other quantities are computed for each station:
the Brune stress drop
the source radius
the quality factor \(Q_0\) of P- or S-waves
Finally, the radiated energy \(E_r\) can be mesured from the displacement spectra, following the approach described in Lancieri et al. [2012].
As a bonus, local magnitude \(M_l\) can be computed as well.
Event averages are computed from single station estimates. Outliers are rejected based on the interquartile range rule.
Attenuation¶
The retrieved attenuation parameter \(t^*\) is converted to the P- or S-wave quality factor \(Q_0^{[P|S]}\) using the following expression:
where \(tt_{[P|S]}(r)\) is the P- or S-wave travel time from source to station and \(r\) is the hypocentral distance.
Station-specific effects can be determined by running source_spec
on several
events and computing the average of station residuals between observed and
inverted spectra. These averages are obtained through the command
source_residuals
; the resulting residuals file can be used for a second run
of source_spec
(see the residuals_filepath
option in
Configuration File).