Locating Earthquakes:At WhatDistance Can the EarthNo LongerBe Treatedas Flat?
J. Arthur Snoke
Virginia Tech
John C. Lahr
U.S. Geological Survey, Denver

INTRODUCTION
Earthquake location programs for events recorded at regional to teleseismic distances have traditionally used travel-time tables such as the J-B tables (Jeffreys and Bullen, 1958), which were based on ray tracing in a laterally homo-
geneous spherical earth (e.g., Boyd et al., 1984). Such travel-
time tables may not give accurate locations for events located using stations at local to near-regional distances because of regional variations in the crustal and upper-mantle velocity structure~the regions of the Earth in which the first-arrival P and S waves spend most of their time for such events. Until recently, it was a nontrivial task to create travel-time tables for different spherical-earth velocity structures. Accordingly, for the past 30+ years, most locations for local and nearregionally-recorded events have used programs which used flat-earth ray tracing through constant-velocity layers. Exam-
ples of such programs are HypoTi (Lee and Lahr, i972), Hypoinverse (Klein, 1985), Fasthypo (Herrmann, 1979), and Hypoellipse (Lahr, 1999).
With the advent of digital recording and advances in telemetry, there are an increasing number of data sets which include arrival times for small to medium-sized events at both local and regional distances. Most researchers are aware of the fact that the assumption of lateral homogeneity may be invalid for event-station paths at larger epicentral distances, but when shown results as given in Figure 1 below, seismologists often express surprise that the spherical-earth
travel-time corrections became significant (e.g., > 0.1 s) at
what seemed relatively small epicentral distances. The purpose of this paper is to document the limita-
tions of using a flat-earth velocity model to approximate travel times within the earth and to suggest strategies for minimizing errors in the analysis of events which are observed at distances for which the spherical-earth corrections become significant.

RAYTRACINGIN A SPHERICALEARTHAND IN A FLAT EARTH
A review of the math and physics involved in flat-earth and spherical-earth ray tracing can be found in most books on geophysics which have sections on earthquake seismology,
for example, Chapter 5 in the third edition of Physics of the Earth (Stacey, 1992).
For regions in which the seismicity occurs primarily in the crust, velocity models are usually developed using local or regional networks and flat-earth ray tracing. In a flat-earth model, a down-going ray can get to the surface only if it is reflected at a velocity discontinuity or if it is critically refracted at a discontinuity and travels along it as a head wave (commonly called P, or S, if the discontinuity is the Moho). If the earthquakes are assumed to be in the crust, the mantle is usually given as a single layer (or half space), and it enters into travel-time calculations only as the velocity for the horizontal part of the path of the head wave.
Typically, only the first P-wave and S-wave arrivals are used when locating earthquakes. For crustal earthquakes at short distances, these will correspond to the direct, up-going
waves (Pg or Sg). At greater distances, because the mantle
velocities are higher than crustal velocities, the P, or S, waves come in as the first arrivals. (We use Kennett's definition for P, or S, [Kennett, 1983, p. 204]: body waves which travel in the uppermost mantle. Note that this is not the same definition as given by Aki and Richards [1980, pp. 212-214], who use the term P, for a head wave refracted along the top of the mantle in a spherical Earth. Our usage is consistent with the
notation in the iasp91 travel-time tables [Kennett and
Engdahl, 1991 ].) The distance at which the arrival times for
the direct (Pg) and refracted (P,) rays are the same is the crossover distance. If the P-to-S velocity ratio is not constant with
depth, the crossover distances will differ for P and S. The crossover distance will also decrease with increasing focal

538 SeismologicalResearchLetters Volume72, Number5 September/October2001

depth; it is about 82 km for a surface focus using the iasp91
velocity model (Moho at 35 km). For events which have focal depths beneath the Moho,
flat-earth ray-tracing can be used in principle as long as the first arrival is the direct arrival. Although the direct arrival is the first arrival for up to several hundred kilometers, we show below that assuming flat-earth geometry results in significant errors in travel times for epicentral distances of only a few hundred kilometers.
DISCUSSION
For our tests, we use the iasp91 velocity model (Kennett and
Engdahl, 1991), which has a two-layer crust with a Moho at 35 km depth and a small increase in velocity with depth in the uppermost mantle (Table 1). We used the same model for the flat-earth ray tracing, except we divided the mantle into constant-velocity layers, as is required by our flat-earth ray tracer.
For epicentral distances less than the crossover distance, the flat-earth and spherical-earth arrival times agree to within a few hundredths of a second. The crossover distances differ in general by at most a few kilometers between flatearth and spherical-earth ray tracing. Beyond the crossover distance, the flat-earth times become increasingly slower than the spherical-earth times (Figure 1). (For a constantvelocity mantle, the chords traveled by the spherical-earth ray become increasingly shorter than the path along the Moho traversed by the flat-earth ray.) This difference in distance, and hence the differential travel time, increases very slowly; if the epicentral distance is less than a few hundred kilometers, spherical-earth P, does not penetrate more than a few kilometers into the mantle.
For a crustal focus, the differential travel time is less than 0.1 s for distances less than 187 km, which is about twice the crossover distance. The effect on event location of using flatearth ray tracing for larger epicentral distances will depend on the details of the event-station geometry. One scenario which could lead to systematic mislocations is using arrivals from distant stations in a small azimuthal range. Of course, if the crustal thickness is not uniform beneath the epicenter and the locator network, there will also be systematic and potentially significant errors due to assuming a uniform structure. Such errors can be minimized by using station corrections or separate velocity structures for different stations. The main point of this paper is that using flat-earth ray tracing at large epicentral distance produces systematic errors which can be calculated and eliminated.
As an example, for some of the larger events within the regional network of the Alaska Earthquake Information Center (AEIC), arrival times for very distant stations in the
Aleutian arc (e.g., SMY, Shemya) or along the Alaskan panhandle (e.g., SIT, Sitka) can be measured. The use of spheri-
cal-earth travel-time tables beyond 1,000 km provided residuals without the large bias that would have been associated with a flat-earth model.

TABLE 1
iasp91 P-wave Velocity Model Down to 210 km
(Kennett and Engdahl, 1991)

Depth (km)

P Velocity (km/s)

0

5.8

20

5.8

20

6.5

35.0

8.04

77.5

8.045

120,0

8.05

120.0

8,05

165.0

8.175

210.0

8.30

Repeated depths indicate discontinuities which are first-order (20 kin, 35 km) or second-order (120 kin). For flat-earth ray tracing, we used the samemodel exceptthat we divided the upper mantle down to 150 km into 12 constant-velocity layers.

Also shown in Figure 1 are two differential travel-time curves for events with focal depths in the mantle. One such data set was collected for events which occurred in the horizontal section of the subducting Nazca Plate beneath central Peru. The focal depths ranged from 100 km to 150 km, and the locator stations in a temporary array deployed in 1985
had epicentral distances of up to 700 km (Norabuena et al.,
1994; James and Snoke, 1994). Compared to the crustalfocus example, the errors become significant at a shorter epicentral distance for a deeper focus. For the Peruvian studies, the locations did not differ by more than a few kilometers in either epicenters or depths because of the favorable geometry of the locator arrays, but correcting for this known source of error allowed the researchers to find the best velocity structure for the region.
Another product of event-station ray tracing is the calculation of body-wave take-off angles, which are used (along with the event-based azimuths and arrival polarities) to
determine focal mechanisms (e.g., Snoke, 1989). For the
runs done to produce the curves in Figure 1, take-off angles differed by at most a couple of degrees throughout the epicentral distances plotted. Experience has shown that a bigger source of bias or error in constraining focal mechanisms
comes from using constant-velocity crustal layers. Teague et al. (i986) used an option in the program Hypoellipse of a
constant gradient layer over a half space to get around the "quantization" of take-off angles for very shallow earthquakes observed at nearby stations, which occurs if one uses only constant-velocity layers.
In principle, it is possible to circumvent the problem for subcrustal events by using earth-flattening transformations
for the velocity structure (e.g., Buland and Chapman, 1983).
For a constant-velocity spherical-earth layer, the transformed

Seismological ResearchLetters September/October2001 September/October2001 539

' i .....~J 3.0

2.5
co I
I, 2.0
0
E
i=
1.5 t21
~m
c" 9 1.0 (1.}
6~
ID
0.5
,
0

.~ .~ o~

oO"
~176176
o.~~ ~176 ~

~ o.,
0

.o o,-

. ..*~176 ~

..... ..... ...., ...,.. ..... .....,

.,... ..,.....-

..... ....., .....
i s

---r ~, ~ 1 , , , 1 . . , 1 , , ....j I J ,

400

800

Epicentral Distance (km)

,,1
1200

,~, Figure 1. Differential travel time versus epicentral distance for three focal depths starting at an epicentrat distance for which the time differential is 0.1
s. The differential travel time is defined as the flat-earth travel time minus the spherical-earth travel time for the Pwave based on the iasp91velocity model.
The curves are labeled by the focal depths in kilometers. For other crustal focal depths, the curves are almost identical to that shown for a surface focus, as it is only in the head-wave (flat earth) and ,on(spherical earth) parts of the ray path that the differences become nonnegligible. In the mantle, the differential travel times increase with focal depth for fixed epicentral distance.

flat-earth velocities are not constant. Most flat-earth location programs require constant-velocity layers, so to simulate this in such programs, many additional constant-velocity layers would be required.
An alternative to using flat-earth ray tracers is to use a modern version of the j-B tables~a travel-time table computed for a spherical earth using an appropriate velocity model. This approach has been adopted for the program
Hypoellipse,which now includes the option to use spherical-
earth travel-time tables either for all distances or for distances greater than a chosen threshold (Lahr and Snoke,
2001). A limitation of this solution is that the iasp91 soft-
ware package (Kennett and Engdahl, 1991) was fine-tuned
for the iasp91 earth model, and there are some models for
which the code will not work, including the J-B S-wave velocity model. An alternative package for calculating spherical-earth travel times is The TauP Toolkit (Crotwei1 et. al., 1999). This package is stable for a wider range of velocity structures. However, since the package is written in Java, it would be more difficult to incorporate it in programs written in Fortran or C.

The Hypoellipse package (Lahr and Snoke, 2001) includes the program Hypotable for generating sphericalearth travel-time tables for use with Hypoellipse. Hypotable uses the iasp91 software package to generate tables of travel
times in seconds versus epicentral distance in kilometers for
a spherical earth with a radius of 6,371 km. To compensate
for variations in the earth's radius from 6,378 km at the
equator to 6,357 km at the poles, Hypoellipsecalculates the
angular epicentral distance using geocentric latitudes and
then uses the local earth radius of the epicenter to convert
the distance to kilometers.
Generally, it is only at local and near-regional distances
that a detailed crustal velocity structure different from the
iasp91 model is required to produce reliable travel times and
take-off angles. Hence, using an appropriate model with flat-
earth ray tracing for small epicentral distances and the
iaspgl-package-generated spherical-earth travel-time tables
beyond a chosen threshold distance (such as the crossover
distance) is a reasonable compromise that provides suffi-
ciently accurate results when locating earthquakes with arriv-
als ranging from local to regional distances. Ell

540 Seismological Research Letters Volume72, Number5 September/October2001

ACKNOWLEDGMENTS
This paper has benefited from USGS internal reviews by Jim Dewey, George Choy, and Ray Buland. The iaspgl TauP routines used in program Hypotable are based on the source code made available by Ray Buland (ffp://ghtftp.cr.usgs.gov/
pub/tau/exped mentat/tau.ta r.Z).
REFERENCES
Aki, K. and P. G. Richards (1980). Quantitative Seismology:Theoryand Methods, W. H. Freeman, San Francisco, CA.
Boyd,T. M., J. A. Snoke, I. S. Sacks,and A. RodriguezB. (1984). High resolution determination of the Benioffzone geometry beneath southern Peru, Bull. Seism. Soc.Am. 74, 559-568.
Buland, R. and C. H. Chapman (1983). The computation of seismic travel times, Bull. Seism. Soc.Am. 73, 1,271-1,302.
Crotwell, H. E, T. J. Owens, and J. Ritsema (1999). The TauP Toolkit: Flexible seismictravel-time and ray-path utilities, Seism. Res. Lett. 70, 154-160.
Herrmann, R. B. (1979). FASTHYPO: A hypocenter location program, EarthquakeNotes50(2), 25-37.
James, D. E. and J. A. Snoke (1994). Structure and tectonics in the region of flat subduction beneath central Peru. Part I: Crust and uppermost mantle,J. Geophys.Res.99, 6,899-6,912.
Jeffreys, H. and K. E. Bullen (I 958). SeismologicalTables,British Association for the Advancement of Science, GrayMilneTrust, Office of the British Association, Burlington House, London, England.
Kennett, B. L. N. (1983). SeismicWavePropagationin StratifiedMedia, Cambridge University Press, Cambridge.
Kennett, B. L. N. and E. R. Engdahl (1991). Travel times for global earthquake location and phase identification, Geophys.J. Int. 122, 429--465.
Klein, E W. (1985). User'sGuide to HYPOINVERSE, a Programfor VAXand Professional350 Computersto Solvefor Earthquake Locations, U.S. GeologicalSurveyOpen-File Report 85-515, 53 pp.

Lahr, J. C. (1999). HYPOELLIPSE Y2K: A Computer Programfor Determining LocaIEarthquake HypocentralParameters,Magnitude, imd First-motion Pattern, U.S. Geological Survey Open-file Report 99-023, paper (112 pp.) and online (http://
greenwood,cr.usgs.gov/pub/open-fite-reports/ofr-99-O023).
Lahr, J. C. and J. A. Snoke (2001). The HYPOELLIPSE earthquake location program, in Lee,W. H. K., H. Kanamori, P. C. Jennings, and C. Kisslinger (editors), IASPEI Centennial International Handbook of Earthquake and Engineering Seismolog~ Academic Press, San Diego (in press).
Lee, W. H. K. and J. C. Lahr (1972). HYPO71' A ComputerProgram for DeterminingHypocenter,Magnitude and First-motionPatternof LocalEarthquakes, U.S. GeologicalSurveyOpen-File Report, 100
PP. Norabuena, E. O., j. A. Snoke, and D. E. James (1994). Structure of
the subducting Nazca Plate beneath Peru, J. Geophys. Res. 99, 9,215-9,226. Snoke, J. A. (1989). Earthquake mechanisms, in James, D. E. (editor), Encyclopedia of Geophysics, Van Nostrand Reinhold Company, 239-245. Stacey,E D. (1992). Physicsofthe Earth, 3rd ed., BrookfieldPress, Brisbane. Teague, A. G., G. A. Bollinger, and A. C. Johnston (1986). Focal mechanism analyses for eastern Tennessee earthquakes (19811983), Bull. Seism. Soc.Am. 76, 95-109.
Department of Geological Sciences Virginia Tech
Blacksburg, VA, 24061 USA snoke@vt.edu
(J.A.S.)
U.S. Geological Survey DFC
Denver, CO, 80225, USA Iahr@usgs.gov
(J.C.L.)

Seismological ResearchLetters September/October2001 September/October2001 541