MLOC Crustal Models

MLOC Crustal Models

This section discusses the format used for custom crustal velocity models in mloc, and some aspects of developing a crustal model for a specific cluster. Although we call them “crustal models” the models discussed here actually consist of one or more layers in the crust over a representation of the upper mantle which can be either a single thick layer (pseudo-half-space) or several layers extending to a depth that is adequate to account for refracted Moho arrivals to a reasonable distance.

The code used in mloc for this purpose is borrowed from the single-event location program HYPOSAT (Schweitzer, 2001). HYPOSAT is capable of travel-time calculations for a large number of crustal seismic phases with models constructed of gradient layers as well as constant-velocity layers. For mloc, however, many of HYPOSAT’s capabilities have been circumscribed.

Usage

A custom crustal model is declared with the command lmod. The use of a custom crustal model is optional. If a custom crustal model is not requested mloc will determine travel times and derivatives for all seismic phases using the tau-p representation of the ak135 global model (Engdahl et al., 1998).

If a custom crustal model is specified in the command file or interactively by the user when running mloc, that model will be used to calculate travel times and derivatives for the crustal phases Pg and Sg and the Moho-refracted phases Pn and Sn, but only to an epicentral distance that is specified in the first line of the crustal model file. Typically that distance will be close to the Pn/P cross-over distance, around 17° epicentral distance. Beyond the specified distance, the tau-p implementation of ak135 is used for all other phases besides the ones calculated with the custom model (so, Pn will not exist beyond 17° in this example). Teleseismic phases that are observed at short epicentral distances (e.g., PcP) are still calculated from ak135 with the tau-p formulation.

Phases

In general, a custom crustal model is used in mloc only when doing direct calibration, i.e., estimating the hypocentroid of the cluster using arrival time data mainly from direct crustal phases at epicentral distances where they are the first-arriving phase (some Moho-refracted arrivals sometimes sneak in). Experience has shown that this kind of analysis works best when all reported arrival times for crustal phases are treated as Pg or Sg. Identification and nomenclature of crustal phases in standard catalogs (e.g., the ISC Bulletin) is too haphazard to rely on, so many arrivals must be assigned new phase codes for use in mloc. The travel-time differences between different theoretical branches of crustal phases are insignificant, given that the crustal structure is seldom known in any detail and the locations of the sources are, before a calibrated location has been determined, often uncertain by tens of kilometers or more. This philosophy does result in loss (or corruption) of small amounts of data in some cases, but the benefits of a simplified approach to phase identification for crustal phases far outweighs those losses.

Model Structure

Although the implementation of HYPOSAT’s code for calculating travel times is capable of handling models with many layers, including layers with gradients in velocity, experience has shown that crustal models used to determine calibrated locations in mloc should consist of no more than 3 crustal layers. Gradient layers are supported in mloc but I nearly always use constant-velocity layers in order to keep things conceptually simple; I have never seen a dataset that required gradient layers in order to adequately model the arrival time data. Although there is a widespread belief that the key to accurate locations of earthquakes is a very complex (but presumably accurate) velocity model, the reality is that most arrival time datasets do not have the quality and coverage to reliably resolve much more than the average crustal velocity, i.e., a single layer of constant velocity (or a gradient at best), and the average crustal thickness when the source locations are also uncertain. In many regions, of course, some additional constraints on gross crustal structure can be derived from other kinds of studies, but arrival time datasets are in general insensitive to structures more complex than an upper and lower crustal layer, and perhaps a sedimentary layer on top. The strong variations in velocity at very shallow depths that undoubtedly exist in nature are averaged out over an earthquake cluster extending across 100 km or so, and have negligible influence on travel times at local and regional distances in any case.

The details of upper mantle structure in a custom crustal model for mloc have very little impact on the relocation. The model needs to extend to a depth that is adequate to generate Pn arrivals to the maximum distance specified in the first line of the file, several hundred km in the usual case where the crustal model is used out to the cross-over distance from Pn to P, around 17°. Since only the relative arrival times of Pn phases would normally be used for relocation, the absolute accuracy of the theoretical Pn times is of minor concern and in any case observed Pn times generally vary by many seconds as a function of azimuth and distance (i.e., at different stations) so the best one can expect with a 1-D crustal model is average agreement. In many cases a single thick layer representing the upper mantle is adequate, but there is no harm in using some version of the layered upper mantle from ak135 or another favored global or regional model, as long as it does not contain any low-velocity zones. This will result in a slight downward curvature of the theoretical travel-time curve of Pn at larger distances that may or may not be observed in the actual data.

Low-velocity zones are not supported.

Format

Formatting of the file defining a custom crustal model for mloc is exactly the same as defined for HYPOSAT. The format is very specific and small deviations can lead to unexpected results. The subroutine that reads crustal models is named rd_loc_mod, found in the module hyposat_loc.f90. The first line is read by:

   read (iunit,'(a)') line
   read (line(1:4),'(f4.1)') rmax

The distance limit to which this model will be used is read from the first four characters as a real value in f4.1. The rest of the line can be used for a comment about the model. All following lines in the file are read by:

   read (line,'(3f10.3,a4)') h(i), v0(1,i), v0(2,i), azo(i)

where h is a depth in km. v0 is an array with two columns that hold Vp and Vs for the layer in km/s, and azo is a character variable with the value “CONR” or “MOHO” that identifies the major crustal discontinuities. The label is used only in the line that defines the “upper” side of the discontinuity. First-order discontinuities require two lines at the same depth in the input file.

The code recognizes the “CONR” flag but for the reasons discussed above, the declaration of a Conrad Discontinuity in the model is strongly advised against. It is possible (probably, it has not been much tested) to configure mloc to use the associated Pb and Sb phases but no case has ever been encountered where it was helpful and it is usually very unhelpful. On the other hand, the Moho discontinuity should always be declared.

Examples

Here is a crustal model based on the ak135 global model:

17.0 Crust and upper mantle of ak135
     0.000     5.800     3.460    
    20.000     5.800     3.460    
    20.000     6.500     3.850    
    35.000     6.500     3.850MOHO
    35.000     8.040     4.480    
    77.500     8.045     4.490    
   120.000     8.050     4.500    
   165.000     8.175     4.509    
   210.000     8.300     4.518    
   210.000     8.300     4.523    
   260.000     8.483     4.609    
   310.000     8.665     4.696

This model uses gradient layers in the upper mantle, following the original, but a simplified version using a single layer of constant velocity would work as well:

17.0 Crust of ak135 with simplified upper mantle
     0.000     5.800     3.460    
    20.000     5.800     3.460    
    20.000     6.500     3.850    
    35.000     6.500     3.850MOHO
    35.000     8.040     4.480    
   210.000     8.040     4.480

The fit to observed arrival times of Pn at larger distances might be somewhat poorer but that would not compromise a calibrated relocation analysis.

Location

As discussed in the section on data files, custom crustal models are normally stored either in the cluster directory (along with the event data files, command files and output files from individual runs), or in the tables/crust directory. Because they are referenced by the full path name in a command file they can actually reside anywhere in the user’s directory structure.

Forensics

The easiest way to determine the details of the crustal model that was used for a specific run of mloc is to have access to the cluster directory, find the “.cr” model in the directory and check in the header section of the “.summary” file that it is indeed the model that was used; the pathname will be given in parentheses after the global model that was used in the line that begins with “Travel times: “. If for some reason you do not have access to the crustal model file, the full crustal model (if one was used) is also logged in the “.log” file of the run.

Model Refinement

For a calibrated relocation analysis, I normally make the first run with the ak135 model distributed in Data Files. The only difference between this and running without a custom crustal model (using tau-p formulation of ak135 for everything) is that the only crustal phases will be Pg and Sg. Decisions about refinement to the model are based on examination of several output plots of travel-time vs distance (see Summary Plots):

  • Focal depths will need to be considered in parallel with the crustal part of velocity model. They don’t need to be completely finalized but they need to be in the ballpark. For most clusters a well-chosen default (average) depth for the cluster will suffice. Crustal thickness cannot be established until focal depths are stabilized, and systematic changes to focal depths at a later stage of analysis will require adjustment of crustal thickness.
  • Changes to Vp and Vs in the crust are primarily based on observation of distance-dependent residuals in the local distance TT (tt5) plot and to a lesser extent the near-source TT (tt4) plot. It is usually very difficult to constrain the relative thicknesses of crustal layers. Determining whether to make changes in upper or lower layers of the crust can be guided to some extent by attention to the pattern of residuals in the near-source TT plot, but the effects are similar to the signature of errors in assumed focal depth.
  • Changes to the Vp/Vs ratio are based primarily on the offset from zero mean of the P and S residuals of near-source TT (tt4) plot and to a lesser extent on the pattern of residuals in the local distance TT (tt5) plot. They should both have offsets with absolute value less than ~0.1 s in the near-source TT plot.
  • Changes to crustal thickness are based on the fit to Pn arrivals in the local-regional TT (tt6) plot. When the hypocentroid is constrained using crustal phase alone, as it is in direct calibration, the arrival times of Pn are controlled by crustal thickness. Thicker crust = later Pn arrivals. Sn arrivals (local-regional S plot, tt7) can be used for this too, of course, but they are often a bit inconsistent with the Pn data, which should be preferred.
  • Changes to Pn velocity are best made through reference to the local-regional TT (tt6) plot, as upper mantle velocity variations are more clearly reflected at greater distances. However, if there is a discrepancy between closer and farther distances (often around 10°) I would favor the near-regional data for refining the model.

It should be emphasized that in direct calibration the most important aspect of the crustal model is the average crustal velocity, as measured by the fit of the observed arrival time data to the theoretical times over the epicentral distance specified for the location of the hypocentroid, which can range from as close as 0.2° to as much as 1.5° or more. The Pg/Pn cross-over distance (a function of average focal depth and crustal thickness) generally establishes an upper bound on this distance, but if the data permit (sufficient number of readings with good azimuthal coverage), restricting the distance range reduces the importance of the velocity model for the calibration. Beyond the cut-off distance for estimation of the hypocentroid, the fit of the observed data to the theoretical TT curves is of little consequence for the relocation analysis, but the the average crustal thickness obtained by fitting the Pn data closely is certainly a parameter of broader interest, worth a careful investigation. Although a calibrated relocation analysis with mloc can reveal useful information about certain aspects of the velocity structure in the source region, it must be remembered that the model was developed only to predict travel times of certain seismic phases. Those phases sample the crust in only limited ways and do not constitute a credible investigation of crustal structure as a whole.