This introduction is a brief overview of the moment tensor inversion problem as handled by the MTUQ
package. It is directly copied from the Google Colab notebook that I used for the SCOPED 2024 Workshop session on moment tensor inversion. The notebook can be found here.
Please reach out to me if you have any issues or questions about the content of the Colab notebook.
Introduction to Moment tensor inversion with the Cut-and-Paste method
In this brief introduction, we will see what the moment tensor inversion problem is about. This should cover all the theory you should be familiar with before handling MTUQ, starting with:
The seismic moment tensor
The seismic moment tensor $\mathbf{M}$ is a mathematical representation of the motion applied to a point-sized seismic source. It is a $3 \times 3$ second-order symmetric tensor with 6 independent parameters.
In the ideal case of a motion applied onto a perfect fault-plane, it can be represented by a supperposition of force-couples or dipoles (with zero net-torque):
Therefore, by estimating the moment tensor parameters for seismic sources that fall into the point-source approximation we can learn about their geometrical properties, and the general motion that was applied to the subsurface - generating seismic waves.
A moment tensor can be conveniently represented as the so-called beachball, which is a representation of the P-wave first motions in 3D space (binary representation of the P-waves radiation pattern). The color of the beachball reflects the direction of motion applied in different directions. It is generally white for motion pointing inward and black (or color shaded) for motion pointing outward of the unit sphere.
Most of the time, moment tensor “beachballs” are visualized as 2D diagrams. These diagrams are stereonet projections of the lower hemisphere of the unit sphere, as seen from above in a map view. This representation allows for easier interpretation and analysis, as it conveys the general motion on a fault at a glance, or the source type in the case of non double-couple sources.
Moment tensor parameterization and the Lune diagram
Moment tensors attributes can be broadly defined as source type and orientation, which are both informed by the geometrical properties of the tensor. Given a moment tensor $\mathbf{M}$, whose eigvenvalues $(\lambda_1, \lambda_2, \lambda_3)$ satisfy $\lambda_1 \geq \lambda_2 \geq \lambda_3$, we can represent a moment tensor as its eigendecomposition:
$\mathbf{M} = \mathbf{U \Lambda U}^T$
where $\mathbf{U}$ contains the eigenvectors that satisfies $\mathbf{UU}^T = \mathbf{I}$, and $\mathbf{\Lambda}$ is the diagonal matrix containing the eigenvalue triple in descending order.
With this decomposition in mind, we see that $\mathbf{U}$ encodes the source orientation (the tensor’s eigenframe), and $\mathbf{\Lambda}$ encodes the source type. From this decomposition stems a very natural way of representing source-types, by projecting the tensor position on the eigenvalue lune:
In the following example, we show for 4 fixed eigenframes, how the moment-tensor source type changes as a function of $\mathbf{\Lambda}$ only:
In MTUQ
, we adopt a moment-tensor parameterization that exploits the geometrical properties of the moment tensor so that we explore the moment-tensor space in a uniform manner. Instead of exploring the 6 parameters of the moment tensor directly, we will express moment tensors with the following set of coordinates:
- $γ, δ$ for the eigenvalue triple, $γ$ being the lune longitude and $δ$ the lune latitude,
- $κ, σ, θ$ for the orientaion, where $κ$, $σ$, $θ$ are for strike, dip and, slip angles respectively.
- $ρ$ gives the scalar seismic moment $M_0 = ρ \sqrt{2}^{-1}$
such that we have $\mathbf{M}(ρ, γ, δ, κ, σ, θ)$.
For more details, you can read Tape & Tape (2015).
Moment tensor inversion
In MTUQ
we solve the moment tensor inversion problem by finding the source that minimizes the following cost function:
$C(\mathbf{M}) = \sum_{r=1}^{N} \sum_{i=1}^{3} \int_0^T \left[d_i^r(t,\mathbf{M}) - d_i^{r,obs}(t)\right]^2 dt$
which translates to the sum of least-squares waveform difference between the synthetic waveforms $d(t)$ and the observed data $d^{obs}(t)$ for N stations with 3 components.
In order to speed up the optimization process, we rely on a set of pre-computed Green’s functions $\mathbf{G}^r$ between each source-receiver pair $r$, such that we can recast the problem as
$C(\mathbf{M}) = \sum_{r=1}^{N} \sum_{i=1}^{3} \int_0^T \left[\mathbf{G}_i^r \mathbf{m} - d_i^{r,obs}\right]^2 dt$
where the time dependency has been dropped for convenience, and $\mathbf{m}$ is a vector containing the 6 independent parameters of the Moment tensor $\mathbf{M}$. Using pre-computed Green’s functions allows the rapid evaluation of the residual term, which we simplify as:
$\mathbf{r} = \mathbf{Gm} - \mathbf{d}$
which is less expensive than explicitly modeling the waveform for each source in a given earth model.
The Cut-and-paste method
The algorithm underlying MTUQ
is known as the “cut-and-paste” (CAP)method of Zhu & Helmberger (1996).
The main goal of the CAP method is to mitigate the effects of 3-D velocity structures that are unaccounted for with a 1D reference model (or an inexact 3D model).
It does so by splitting all the data into specific window groups (typically P-waves and Surface waves) and allowing a time-shift window for the synthetics. Within this allowable time-shift window, a cross-correlation is performed between $d_{syn}$ and $d_{obs}$, and the residual is computed only for the maximum cross-correlation value (when signals are supposedly aligned).
We do this for each data group (P, Rayleigh, Love waves), for each station and each component (Z-R-T), and return one misfit value corresponding to the maximum cross-correlation value.