# Seismological Research Letters

- © 2012

*Online material:* Software source code; examples; sample plots.

## INTRODUCTION

The graphical representation of a focal mechanism as a plot of the radiation pattern on the focal sphere (sometimes referred to as “beachball diagram”) is well established. It consists of a spherical projection on a unit sphere centered at the source point of a seismic event. The sign of the *P*-wave radiation pattern is shown as colored areas on the sphere (as described, *e.g.*, by Udías 1999; Stein and Wysession 2003; and Crotin 2004).

For pure shear motion on a plane, the focal sphere may be fully described by the angles of strike, dip, and slip-rake of this fault plane (*cf.*, Shearer 1999; Aki and Richards 2002). In this case the construction of the spherical projection into the *e.g.*, Stein and Wysession 2003, chapter 5).

Several tools are available to do this work. They exist either as a more or less configurable, independent software package (Snoke 2003; Srivastava *et al.* 2006; Luis 2007; Haeussler and Labay 2007), as supplements for *GMT, MATLAB* (Creager 1997; Center for Earthquake Research and Information 2008) or *Mathematica* (Scherbaum *et al.* 2010), or as an online applet (Helffrich 2007). Almost all of the existing software consider only these pure shear cracks. The visualisation of a more complex source mechanism, *e.g.*, an arbitrary seismic moment tensor with six independent entries, is often not possible.

Only a few tools exist at all that can handle this general case. The most commonly used one is the *psmeca*-tool in *GMT* (Wessel and Smith 1999). However, projections for arbitrary viewpoints are not possible with *psmeca*.

We present the free platform-independent command line tool *MoPaD*, written in the *Python* programming language that allows one to generate, plot, and save scale- and case-independent focal spheres in a variety of standard graphic formats. The tool allows the use of *GMT* for plotting focal mechanisms even in the above mentioned cases by using the *GMT* program *psxy*.

As an additional feature, the moment tensor can be decomposed and transformed into different representations.

Alternatively, the core functionality of *MoPaD* is provided as a *Python* library.

After providing the theoretical basis in the next section of this paper, we will turn to the technical aspects of implementation. We give a short introduction into the handling of *MoPaD* and its methods by examples in the last section.

## THEORY

Many publications describe in detail the theory of source mechanisms (Udías 1999) but do not fully cover their visualisation. The idea of constructing focal spheres is mainly sketched from a geometrical point of view in the case of double-couple sources (*cf.,* Stein and Wysession 2003) or handled as a contour plot of first arrival polarizations (*e.g.*, Riedesel and Jordan 1989).

We start with the description of the stable construction of a general focal sphere diagram (FSD) on the basis of an analytical approach. This is based on the calculation of nodal lines on the focal sphere.

### Definitions and Terminology

The seismic moment tensor M is defined in the form of a real-valued symmetric 3×3 matrix. It is fully diagonalizable, so its characteristics are completely determined by its eigensystem with the eigenvalues {*σ _{P}* ≤

*σ*≤

_{N}*σ*}. The corresponding eigenvectors are the pressure-, null-, and tension-axis. This orthogonal principal axis system is used to derive an exact description of the FSD.

_{T}With **M¯** := (*M*_{11}, *M*_{22}, *M*_{33}, *M*_{12}, *M*_{13}, *M*_{23}) there is a bijection **M¯** is identified with **M** for the remainder of this publication.

Within the following steps, we take advantage of a planar reflective symmetry. The eigenvector that corresponds to the largest absolute eigenvalue is the symmetry plane’s normal vector for all cases except two in which the source has a significant volumetric component. To handle these exceptions, we introduce the *exception-value χ*, using the reduced eigenvalues of **M**

Accommodating the mentioned symmetry, we introduce the nomenclature
*a,b*)↦(*b,a*), which is applied if the exception-value is non-zero. The eigenvalue *σ*_{Σ} defines the symmetry, where Σ is its respective eigenvector. The second eigenvector *N* remains unchanged and the third is denoted by *H* (*cf*. Figure 1).

### Nodal Lines

The *P*-wave radiation pattern is characterized by nodal lines of the polarity of the first arrivals on a unit sphere around the origin on the focal sphere, and we are interested in a closed formula describing these lines.

In the basis system *HN*Σ, the focal sphere is planar symmetric with respect to the *HN*-plane, so the lines are given by a closed curve around the positive Σ-axis (shown in Figure 1) and its reflection on the symmetry-plane.

As a start we calculate the net particle displacement *u*(*x, y, z*) on the focal sphere. For a general solution, the coordinates are expressed in terms of the eigenvector basis *HN*Σ. Due to the symmetry, only points in the upper hemisphere are considered; so *x _{H}, y_{N}, z*

_{Σ}) with

*z*

_{Σ}≥ 0. Following this notation, one obtains for the oriented particle displacement

The effective net particle displacement, pointing radially outward, is obtained by projecting the displacement vector *u*^{*} onto the position vector, expressed via the scalar product

The nodal lines are defined by *u*(*N* and the horizontal projection of *P* to the *HN*-plane, one obtains an explicit function for the polar opening angle *α* between Σ and *P*. This yields the desired description of the two lines around the positive (*α*_{+}) and negative (*α _{-}*) part of the Σ-axis in a closed, yet parametrized form:

The most frequently occurring types pure double couple (*σ*_{Σ} = -*σ _{H}, σ_{N} *= 0) and pure compensated linear vector dipole (CLVD) (

*σ*=

_{H}*σ*= -1/2

_{N }*σ*

_{Σ}) are special cases of this general solution. If the isotropic part of the moment tensor outweighs the deviatoric content, the discriminant in the right hand side of Equation 6 becomes negative, so that no nodal lines can be defined. This leads to a uniformly colored focal sphere, consistently reflecting the constant sign of the radial particle displacement.

### Transformation

A solution in an abstract eigenvector system and its uniqueness is quite satisfactory from a theoretical point of view. For plotting, an arbitrary basis system of the

Having a basis-independent solution for the curve, we append a transformation into the desired basis system to plot the projection of the sphere and all nodal lines. The output can be given in all implemented basis systems (see Table 1). In the following, M is given in north, east, down (*NED*), unless otherwise noted.

### Projections R 3 → R 2

Three projection algorithms are implemented; the tangential stereographic projection is set as default. Other possibilities are lambert azimuthal equal-area and orthographic projection (*cf.* Bronstein and Semendjajew 1987). Figure 2 illustrates the differences between the projections.

In *MoPaD* it is possible to define which hemisphere shall be projected. The standard back-hemisphere projection can be changed into a front-hemisphere projection; thereby the roles of upper and lower hemisphere (as well as the respective projection poles) are interchanged.

### Viewpoint

Figure 3 illustrates the option to transform the viewpoint. The view on the focal sphere solution is defined by three angles. The two angles θ and ϕ define the observer’s position in relation to the sphere, and Ω sets the viewing angle (azimuth) with respect to a local *North N’*. A common FSD is set up so that one looks from above the center (θ = 0, ϕ = 0), while local *North* is on the upper edge of the plot (Ω = 0).

For changing the viewpoint to *V’*, a local coordinate system is set up at this new position. The negative position vector defines the local *Down* = *D’*, while local *North* = *N’* is obtained by an orthonormalization of *D’* and *N* (Gram-Schmidt method). Finally, an outer product of the first basis elements yields *E’*. The set *N’E’D’* provides a right-handed local basis of the

#### Example

We take a closer look at an arbitrarily chosen source mechanism, defined by **M** = (1, 2, 3, -4, -5, -10).

In Figure 4A, the upper hemisphere of the FSD is shown as seen from above (viewpoint (0,0,0)); four additional viewing positions (a-d) for vertical sections are indicated. The center plot of Figure 4B illustrates the standard vertical view on the lower (back) hemisphere from above the source in stereographic projection. It was generated with the *MoPaD* command

Apart from that, the FSDs a-d visualize the respective vertical back hemispheres that have been generated with the *viewpoint*-option arguments (a:{90,0,180}, b:{0,-90,-90}, c:{-90,0,0}, d:{0,90,90}).

### Decompositions

Within *MoPaD*, four possibilities to decompose a moment tensor into different components are available:

Isotropic + double couple + CLVD,

Isotropic + major double couple + minor double couple,

Isotropic + three double couples, and

Isotropic + CLVD + strike-slip movement + dip-slip movement.

The definition and nomenclature of the first three decompositions follow directly (Jost and Herrmann 1989). The fourth alternative is described in Dahm (1993). The definition of the scalar seismic moments is given in Bowers and Hudson (1999).

## IMPLEMENTATION

*MoPaD* is written in *Python* and can be used as a *Python* module or from the command line. In addition to the *Python* standard installation, the modules *numpy* (Oliphant 2006) and *matplotlib* (Hunter 2007) are required.

Developed under the GNU Lesser General Public License (LGPL), *MoPaD* is designed as an open-source tool.

For use as a command line tool, *MoPaD* is subdivided into four methods:

### plot

This method is the main feature of *MoPaD*; it returns a plot of an FSD in several graphic formats. The theoretical shape of the diagram, obtained as described above, is transformed according to given options. After the consecutive projection, a visualisation is provided using the module *matplotlib* (Hunter 2007).

### describe

A brief overview of the provided source mechanism can be obtained using this method. Both the moment tensor representation as well as the orientation of the main shear component are shown.

### decompose

Returns strings containing information about the components of the original source mechanism. Either a full or a partial decomposition of the internally calculated seismic moment tensor can be returned in different basis systems.

### gmt

A string is returned containing *x, y*)-form, which can be fed into the *psxy* function of *GMT*.

A combination of methods decompose and plot allows, for instance, for the graphical decomposition of a moment tensor, as shown in Figure 5.

## EXAMPLES

In order to demonstrate the main features of the methods of *MoPaD*, we present here some small example applications. (All input has to be provided as one line, all linebreaks in the following examples are due to the proper typesetting of this article.)

### plot

In addition to the case presented previously for a varying viewpoint, where we got a very simple plot, we give here an input line for calling *MoPaD* that contains a combination of options. In contrast to the aforementioned case, this plot is not shown but stored directly into a file:

With this input, we have chosen to plot the deviatoric part of the moment tensor with an output size (-s) of 10 cm, a view (-v) from south-east showing the front hemisphere (-U), the orthographic projection (-p), an .svg formatted output file (-f), the tensional area in yellow (-r), and the pressure domain in purple (-w). Additionally, one fault plane of the governing double-couple part is superimposed (-d), and the positions of eigenvectors (-e), as well as the axes of the basis system (-b), are indicated. (The generated file is provided in the electronic supplement of this publication.)

This setup is a highly unusual combination for an everyday focal sphere visualisation, but it serves the purpose to show the action of an arbitrary set of options.

### describe

The method describe provides a fast and simple way to obtain an overview of all parameters of M. In addition, it enables an easy conversion between basis systems (options -i and -o define input and output basis):

### decompose

The moment tensor can be decomposed into different parts. Here we use a partial (-p) decomposition of **M** = (1, 2, 3, -4, -5, -10) in order to determine the double-couple part and isotropic component, including the respective percentages:

### gmt

*MoPaD*’s gmt-method provides the coordinates of the nodal lines of the FSD, as well as the eigenvectors. These are to be piped into the psxy-tool of *GMT*. The tension-areas of the FSD hold the color-key “1” for the psxy-*Z*-option; the pressure areas the key “0.” If requesting only the border lines, the key is also set to “1.”

The complete command for generating the FSD in Figure 4B (center) by using a combination of *MoPaD* and *GMT* is:
*GMT* color tables psxy_lines.cpt and psxy_fill.cpt must be generated externally.) That seems to be rather complicated, but the high flexibility of *MoPaD* proves beneficial in the case of handling focal mechanisms, which are not limited to the standard view and projection. There the commands are not more complex than in the standard case:

This yields a rotated FSD (viewpoint 45°, 45°, 45°) in the orthographic projection. (The files BB1.ps and BB2.ps, as well as the color table files, can be found in the electronic supplement.)

## CONCLUSIONS

The tool *MoPaD* provides a variety of plotting options for the widespread focal mechanism’s graphical representation, visualizing seismic source mechanisms.

The concept given in the Theory section allows for a fast and easy plotting of the FSDs, independent of the source type, orientation, and viewpoint of the observer. In this way, the tool represents a significant improvement for everyday work in seismology. Yielding a simplification for creating arbitrary FSDs in combination with *GMT*’s *psxy, MoPaD* can be a useful supplement.

The software is open source, and a first alpha version is freely available from the electronic supplement for this publication. More recent versions will be available on the homepage of *MoPaD*: http://www.mopad.org. An instruction manual with descriptive examples is provided there together with the source code. *MoPaD* is additionally integrated into *ObsPy,* a *Python* toolbox for seismology (Beyreuther *et al.* 2010).

## Acknowledgments

The development of this new tool was initiated by the lack of easy to use stand-alone tools for plotting general focal mechanisms, not restricted to pure double couples. More intuitive manuals within the seismological community are highly desirable. Moreover, the need for an appropriate *Python*-module for this purpose was a further motivation.

The authors thank their supervisor, Professor Dahm at the Institute of Geophysics at the University of Hamburg, for his idea to implement *MoPaD*’s gmt-method, providing the connection to *psxy*, as well as the fourth non-standard decomposition of the moment tensor. Thanks go also to Stefan Trabs for testing the software.