Procedures to build
a subsurface velocity structure model

English version
Translated in March, 2022

Earthquake Research Committee, The Headquarters for Earthquake Research Promotion
This document is a translation from Japanese to English of
published in April, 2017.

1

1. A subsurface velocity structure model in strong ground motion prediction
To evaluate seismic ground motion hazards as a measure for earthquake disaster prevention, the Headquarters
for Earthquake Research Promotion (HERP) has released National Seismic Hazard Maps, which contain the
results of ground motion estimation obtained with the strong ground motion prediction method for earthquakes
with specified source faults (“Recipe1”).
This “Recipe” indicates that one can define a seismic velocity structure model for earthquake ground motion
evaluation as a three-dimensional velocity structure with the main physical parameters of density, P-wave and Swave velocities, their Q-values, and the shape of a layer boundary for each layer. The seismic velocity structure
model is built independently by dividing the structure into the following three domains bounded by the topsurface of the engineering bedrock2 and the top-surface of the seismic basement3. The three domains arecharacterized as follows.
・Shallow soil layers are the sedimentary layers from the surface to the top-surface of a firm layer having an
S-wave velocity of 300–700 m/s, which is regarded as the engineering bedrock. The thicknesses of theselayers are between zero and several tens of meters. Mainly, these layers influence the amplification ofseismic waves in a short-period range of less than 2 s.
・Deep sedimentary layers are the layers from the top of the engineering bedrock to the top-surface of a layer
with an S-wave velocity of approximately 3 km/s, which is known as the seismic basement. Their
thicknesses are between several tens of meters and approximately 34000 m. These layers influence theamplification of seismic waves in a wide-period range (0.1–10 s) that are of interest in the evaluation ofbroadband ground motion, including long-period components at periods longer than 2 s.
・Crustal structure is the structure deeper than the top-surface of the seismic basement. This structure
influences the characteristics of the seismic wave propagation path effect. Because seismic waves
propagate in the upper mantle in some cases, depending on epicentral distances of earthquakes
considered, it is necessary to include the upper mantle in modeling this structure. Here, this is referred
to as the crustal structure.

Combined shallow and deep
velocity structure model

Figure 1-1. Schematic diagram of seismic velocity structure model

1The published “Recipe” contains an idea for building a seismic velocity structure model, and a subsurface velocity structure

model in whole Japan made with this idea has been released in J-SHIS. For the previous development of seismic velocity structure
models, see Appendix 1 (in the following, a reference to the Appendix is denoted by a superscript).
2This is a firm layer in which design earthquake ground motion is defined for a structural seismic design in constructionengineering fields. Its S-wave velocity varies depending on the kind of designed structures, in a range of 300–700 m/sA2 inmany cases.
3This is a bedrock located at the top of the crust, which has an S-wave velocity of approximately (over) 3 km/s (Architectural
Institute of Japan, 2009).

2

Simulation of earthquake ground motion in a wide frequency range is firstly performed for a calculation ofground motion at the top-surface of the engineering bedrock by using a hybrid method considering anamplification of seismic waves in the propagation path from the source fault to the top-surface of the engineering
bedrock. By adding the amplification effects caused by the shallow soil layers to the simulated results, it ispossible to obtain the earthquake ground motion on the surface of the earth. In the “Japan Integrated Velocity
Structure Model (provisional version),” which is the currently available model, the surface of the engineering
bedrock is set approximately at the layer with an S-wave velocity of 500–600 m/s. To evaluate a spectrum andwaveform of an earthquake ground motion at short periods, the engineering bedrock must be set at shallower
layers with a lower S-wave velocity than those of the present ones. Furthermore, earthquake ground motions
estimated in many hybrid simulations with matching periods of 0.5–2 s between long- and short-period groundmotions are likely to be affected by both the shallow and the deep sedimentary layers. Therefore, it is necessary
to connect the shallow domain smoothly with the deep domain for a combined layer model without a velocity
gap.
The shallow soil layers model was improved by using a large number of soil data, including boring data. The
top layer of the engineering bedrock has an S-wave velocity of approximately 300 m/s in the created model,
which has a horizontal resolution of several hundred meters. The combined model from the surface to the seismic
basement consists of the advanced shallow soil layers model and the deep sedimentary layers model withconsidering their smooth connection. The goal of this combined model is evaluating a broadband ground motion
in a period range of 0.1 to 10 s on the surface of the engineering bedrock (free engineering bedrock) at depths
close to the surface of the earth using the hybrid method.
Figure 1-2 shows the detailed procedures for building the seismic velocity structure models. In Chapter 2, themodeling of the crustal structure, including the upper mantle, is introduced. The modeling of deep sedimentarylayers that are deeper than the engineering bedrock is described in Chapter 3. In Chapter 4, the detailed procedurefor developing a model of the shallow soil layers is explained. Finally, a “combined model of shallow and deep
layers” is built for a usage to evaluate broadband ground motions in Chapter 5.

3

Figure 1-2. Procedure for seismic velocity structure model to evaluate broadband groundmotions4

2. “Crustal structure”
Seismic velocity models beneath the top surface of the seismic bedrock , including the upper mantle, can beconstructed using the velocity structure models estimated with P-wave and S-wave travel times of earthquake
records, results of deep boring surveysA3 such as the MITI (current Ministry of Economy, Trade and Industry)
Exploratory Test Wells, and physical property values obtained with seismic reflection and refraction
surveysA4 .Considering the characteristics of the seismic basement from the viewpoint of its physical properties,
it is necessary to define the top-surface of the basement with a layer that has an S-wave velocity higher than 3
km/s, which can be found in a wide area. Firstly, a relationship between the basement depth and gravity anomaly
values in the area of interest is determined. Then, a distribution of the basement depth is estimated using therelationship and the gravity anomaly map. The Conrad discontinuity (the boundary between the upper and lower
parts of the crust) and the Mohorovičić discontinuity are modeled using the results of the seismic surveys. The
shapes of these boundaries can be modeled with additional data from a regional gravity anomaly distribution.
As indirect information, one can use crust models used in analyses of hypocenter determinations and source
inversions, three-dimensional seismic velocity structures (seismic tomography models), and so on. The S-wave
velocities of these models are considered for interpolating the depths of the seismic basement and the velocity
structures obtained from the original data, and a three-dimensional velocity structure model beneath the seismic
basement are constructed based on the S-wave velocities. When a seismic velocity structure model in an oceanarea is needed, the shape of an oceanic plate estimated from seismic surveys and seismic activity studies in thearea is also considered. The use of a seismic tomography model can be effective for introducing complex velocity
information, such as mantle wedges and subducting slabs, into the crustal and mantle structure model.
In the development of the “Japan Integrated Velocity Structure Model (provisional version),” the surface ofthe seismic basement was constrained by a distribution pattern of gravity anomaly values with data from boring
and geophysical surveys. Furthermore, in modeling the depth to the plate surface in the building of a crustal
model deeper than the seismic basement, the results of studies by Sato et al. (2005) and Baba et al. (2006) wereincorporated, together with referencing other studies (e.g., Ryoki, 1999; Matsubara et al., 2008). In the crustal
model published in J-SHIS, a seismic tomography model described by Matsubara et al. (2008; 2011) (releasedon the website of the National Research Institute for Earth Science and Disaster Resilience) was used.
Furthermore, the crustal structure in the area of the Nansei Islands was modeled with the results obtained by Zhao
et al. (1994) on the Conrad and the Mohorovičić discontinuities. Examples of seismic tomography models of thecrust and the mantle below the Japanese Archipelago are summarized in the “Publications and Databases” on thewebsite of the Earthquake Research Institute at The University of Tokyo, although the information is slightly old.
Measured gravity (anomaly) data for various regions of Japan were published by the Geological Survey ofJapan, which is a part of the National Institute of Advanced Industrial Science and Technology (a DVD versionof the Japan Gravity Database published in 2013 can be downloaded), and numerical data were released in a
database by Kanazawa University (Honda et al., 2012). The Gravity Database (GALILEO) of the National
Institute of Advanced Industrial Science and Technology enables users to specify and view generated plots.
3. “Deep sedimentary layers”
Procedure (1): Collection of geological information and various geophysical data
Existing information sources are collected. These include regional geological maps covering wide areas andvarious geological cross sections, information about soil profiles obtained by deep boring surveys (boreholeobservatory for crustal activity), including PS loggingA5, results from various geophysical surveys (seismicreflection and refraction surveys, microtremor surveys, gravity surveys, and electromagnetic surveys), and studies
using earthquake ground motion records.
In the geophysical and deep boring surveys, information about the deep sedimentary layers in a large-scale
area was also obtained in basic academic researches, such as for earthquake disaster prevention. These data as
the fundamental materials for modeling deep sedimentary layers are available in databases containing subsurface
5

data maintained by the National Research Institute for Earth Science and Disaster Resilience and the National
Institute of Advanced Industrial Science and Technology. The results of basic geophysical surveys andinformation about the MITI Exploratory Test Wells of the Japan Oil, Gas, and Metals National Corporation are
also available for public purposes, such as natural disaster prevention.
Procedure (2): Estimating one-dimensional structure models4 by comparing geological profiles with velocity
structures
Existing geological maps, geological columnar sections, boring data (particularly data from boreholeobservatory for crustal activity), and regional geological stratigraphic data are used for one-dimensional models(geological layer model) composed of multiple strata with individual geological classifications. A onedimensional velocity structure model that corresponds to the above multiple strata is estimated by comparing
physical property values, such as P-wave velocities, obtained in various surveys with the geologicalclassifications. When the data are insufficient, or when only geological information is available, a geologicallayer model is estimated by the geological classifications from available geological data. Then, it is possible toestimate the velocity values from geological classifications. This is done by considering the relationships between
the geological classifications and the physical property values established in previous studies, or with existing
boring data in the other areas (for example, Suzuki, 1996).
To build a one-dimensional structure model that corresponds to the geological layer model, results such as Swave velocities from PS loggings and microtremor surveys, densities from density loggings, and results of model
analysis using earthquake ground motion records are used. They make it possible to determine the S-wave
velocities, densities, and Q-values for individual layers considering their P-wave velocity structures. If only Pwave velocities are available, S-wave velocities, densities, and Q-values can be derived from empirical and/or
theoretical equations with P-wave velocities.
The equation by Gassmann (1951) is an example of the equations available for addressing the relations between
physical values. One can estimate an S-wave velocity and a density when a P-wave velocity and a void ratio are
specified. Furthermore, empirical equations to convert a P-wave velocity (Vp) to a density (ρ) are also available,
such as an equation for a velocity and a density by Ludwig et al. (1970) and an empirical equation of ρ = 0.31Vp0.25
by Gardner et al. (1974). When these theoretical equations and/or empirical equations are used, it is important toconfirm their applicable ranges by referring to the original data used to generate the equations.
A simple assumption for the Q-values can be used at this stage, such as a value proportional to an S-wave
velocity because there is an adjustment of the values as discussed later.

4This corresponds to the zero-th degree model (Fujiwara et al., 2006) in the strong ground motion prediction method forearthquakes with specified source faults (“Recipe”).

6Figure 3-1. Building of one-dimensional structure model by comparing geological
profile with velocity structure
Procedure (3): Development of three-dimensional model using one-dimensional profiles
To estimate a three-dimensional velocity structure model using the one-dimensional velocity structure models
made in Procedure (2) for individual locations or areas, one categorizes layers with similar geologicalclassifications into the same velocity layer and determine their velocity and thickness. For determining the spatial
extent of each layer, spatial interpolation is performed. In this procedure, boundaries of layers with different
velocities detected by seismic reflection and refraction surveys, density structures estimated from the gravity
anomaly distribution filtered at wavelength components corresponding to the target depth of the structures,
existing geological cross sections, geological contour maps, geological ages, and shapes of faults and folding are
considered. Considering the three-dimensional distributions of individual geological boundaries, an initial modelof the deep sedimentary layers from the top of the engineering bedrock to the top of the seismic basement is built.
As the top of the engineering bedrock, a stratum is chosen to have an S-wave velocity used in the definition of
the engineering bedrock. In an area where the top of the seismic basement is very shallow, weathering layers near
the surface over the basement are also considered.
The final shape of each boundary can be refined through a forward analysis and/or an inversion of observed
gravity anomaly data.
The relationship between the gravity anomaly and geologic structures was discussed in detail by Kimura(2002), and their usefulness as geological structure surveys has been mentioned in other works, such as those by
Kono and Furuse (1989) and Yamamoto and Shichi (2004).

7

Figure 3-2. Example of an expansion of one-dimensional velocity structure models in horizontal directions using theseismic surveys and the gravity anomaly distribution. The red lines indicate surveying lines of seismic surveys for Pwave velocity images obtained by reflection surveys and boundary shapes obtained by refraction surveys. The color
contours show the Bouguer anomaly from the Geological Survey of Japan (2013) with an assumed density of 2.0
g/cm3. Layer boundaries are developed in the horizontal directions using this information as interpolation data.

Procedure (4): Confirmation of three-dimensional model5 using earthquake ground motion records
One can estimate S-wave velocities and the depths to boundaries of the deep sedimentary layers using
earthquake ground motion records, such as surface wave dispersion, horizontal-to-vertical spectral ratios, andreceiver functions. One-dimensional profiles from these results are compared with the three-dimensional
sedimentary layer structure model made in Procedure (3), and/or theoretical values calculated from the threedimensional model are compared with these observed values to validate the three-dimensional model. If necessary,the model is adjusted.
The available methods for estimations of seismic velocity structures using earthquake records are the receiver
function method (e.g., Langston, 1979) to detect velocity boundaries in a one-dimensional S-wave velocity
structure below an observation station and a spectral ratio (H/V spectrum) of a horizontal component to a verticalcomponent of Rayleigh wave (as the fundamental mode) to estimate a thickness and an S-wave velocity of eachlayer (e.g., Arai and Tokimatsu, 2004).
Earthquake ground motion data from strong earthquake observation networks densely placed in Japan, such as
K-NET and KiK-net by the National Research Institute for Earth Science and Disaster Resilience, are easily
obtained from the websites. Furthermore, data from the seismic intensity meters obtained by the Japan
Meteorological Agency (JMA) are available to the public through its website and CD-ROMs.
5The adjusted model corresponds to the 0.5-th degree model (Fujiwara et al., 2009) in the strong ground motion prediction

method for earthquakes with specified source faults (“Recipe”).

8Procedure (5): Validation of velocity structure model with simulation of observed earthquake ground motions6
A simulation of earthquake ground motions observed at exiting strong motion stations is performed using theestimated three-dimensional structure model of the deep sedimentary layers. Synthetic seismograms calculatedwith analytical methods based on the elastic wave theory, such as the finite difference method, are compared withthe observed ones to confirm whether the simulation reproduces observed amplitudes, appearance of distinct
phases, amplification characteristics, and predominant periods. The model can be adjusted (model tuning) toenhance its ability to reproduce the ground motion characteristics. Adjusting the velocities and Q-values throughan inversion scheme is recommended.
It is necessary for ground motion simulation to evaluate various characteristics in observed waveforms, suchas amplitudes, predominant periods, arrival times of distinct phases, and time-variant characteristics, with thethree-dimensional model. In this procedure, it is also important to set a proper Q-value of each layer, which
controls the attenuation features of seismic waves. To establish a high-accuracy Q-value model, one must referinformation related with the regionality (i.e., geological classification) of the physical properties for a target
areaA7. Although conducting of inversion analyses, or parameter tuning, by modeling of observed waveforms to
obtain the proper values of those physical properties would be effective, it is not always possible because of thelimitation of data and computing resources. Without tuning the physical properties such as velocities and Qvalues, the ability of the velocity model to reproduce the observed waveforms is limited.
It is necessary to complete the modeling of the velocity structure using Procedure (5) for a ground motion
prediction considering only the effects of the deep sedimentary layers. However, when one combines the modeldsof the shallow soil layers and deep sedimentary layers, it is desirable to conduct the model tuning and theverification with earthquake records after the combination has been done.
4. “Shallow soil layers”
An amplification of earthquake ground motion at short periods (mainly periods less than 2 s) is significantly
affected by S-wave velocity structures in the shallow soil layers from the top of the engineering bedrock to thesurface of the earth. The principal procedure for modeling shallow soil layers is to build one-dimensional multilayered velocity structural models to calculate a ground response by collecting surface geological data and boringdata from soil surveys. However, a large number of the boring data are required to develop a model with asignificant local lateral variation of the shallow soil, which is not rare in actual conditions. Furthermore, it isrecommended, in the building of a three-dimensional model by connecting one-dimensional structures from theboring data, not to apply a simple mathematical interpolation. Instead, the continuity of the surface geologicalcondition should be considered, and tuning should be performed by considering S-wave velocity profiles fromgeophysical surveys, such as microtremor survey. However, the boring data are insufficient to build three-dimensional model in many areas. In this case, one needs to extrapolate information and could get a spatiallysmooth three-dimensional model. One needs to judge the model appropriateness for applying his frequency
range of the seismic waves to be treated and computer resources．
Here, we explain a method using geomorphologic classification model and a method using the boring data formodeling shallow soil layers. Two procedures are introduced for the latter method. One is an interpolation method
using geomorphologic classification for an area without the sufficient data, and the other is a detailed method thatemphasizes the three-dimensional continuity of the surface geological conditions which is used for the combinedmodel of the shallow and deep layer (hereinafter referred to as “SD model”).
6The adjusted model corresponds to the first-degree modelA6 in the strong ground motion prediction method for earthquakes

with specified source faults (“Recipe”).

84.1 Method for building of the geomorphologic classification model
Geomorphologic classification model
The region of interest is discretized into appropriately sized meshes, considering the scales of existing maps,
such as terrain maps and soil classification maps. When a mesh contains multiple kinds of the geomorphologic
units, we choose a geomorphologic unit that occupies the largest area within the mesh or a geomorphologic unit
at the center of the mesh to represent its classification. The representative topographical and geomorphologic
units in each mesh are recategorized with the given criteria. For the categorized meshes, the site amplification
characteristics are empirically estimated with explanatory variables of an altitude, a slope angle, and a distancefrom mountains or hills formed in old geological ages.
An AVS30 (the mean S-wave velocity in the top 30 m from the surface), which has a high correlation with thesite amplification factors of the shallow soil layers, is estimated using an empirical equation based on PS logging
dataA8. Moreover, the amplification factor is calculated for each geomorphologic classification from therelationships of the AVS30 with the site amplification factor of the peak ground velocity, and so on. Thisapproach can be used to evaluate the site amplification characteristics throughout Japan, including areas withinsufficient soil data.
The “site amplification capability map based on the 7.5-arc-second Japan Engineering Geomorphologic
Classification Map” (Matsuoka and Wakamatsu, 2008) is a geomorphologic classification model. In this model,
an AVS30 is estimated using an empirical equation (Matsuoka et al., 2005) with the geomorphologic
classification in the “7.5-arc-second Japan Engineering Geomorphologic Classification Map” (Wakamatsu andMatsuoka, 2003). Then, a site amplification factor is calculated from the AVS30 using an empirical equation of
the site amplification factors with the mean S-wave velocities of the shallow soil layers (Fujimoto andMidorikawa, 2006). The “7.5-arc-second Japan Engineering Geomorphologic Classification Map based on theuniform national standards” (Figure 4-1) (Wakamatsu and Matsuoka, 2013), which was made by reconsidering
the geomorphologic classifications based on new geomorphological data, is also used in the modeling.

Mountain (Pre-Tertiary)
Mountain (Tertiary)
Mountain footslope
Hill
VolcanoVolcanic footslope
Volcanic hill
Rocky strath terrace
Gravelly terrace
Loam terrace
Valley bottom lowland
Alluvial fan
Natural levee

Back marsh
Abandoned river channels /Former pond and swamp
Delta and coastal lowland
Marine sand and gravel bars
Sand dune
Lowland between coastal bars and/or dunes
Reclaimed land
Filled land
Rock shore, Rock reef
Dry riverbed
River bed
Lake
OthersFigure 4-1. Geomorphologic classification map in World Geodetic System with reconsidered geomorphologic
classification (Wakamatsu and Matsuoka, 2013)9

4.2 Method for building of soil model using boring data
Procedure (1): Assignment of various collected data to each mesh and classification of “mountain/hill”
and “plateau/lowland” based on geomorphologic classification data
After a target area is divided into meshes, one marks the positions of available boring data, and assigns
attributes, such as the availability of PS logging data, the number of boreholes with their drilling depths, and thegeomorphologic classifications to each mesh. One also notes availability of geophysical surveys, such as seismic
and microtremor surveys without any boreholes in meshes along surveying lines.
A classification criterion (Wakamatsu et al., 2004) was made by considering geomorphologic terrain
classification standards for engineering uses to evaluate shallow soil layer characteristics. The criterion is used
to classify the meshes into the classification of “mountain/hill” corresponding to “mountain, a mountain footslope,a hill, a volcano, a volcanic foots-lope, or a volcanic hill”, and the classification of “plateau/lowland” defined as
the geomorphologic units other than “mountain/hill.”
Some local governments and some institutions have attempted to release borehole information to the public.
Examples include Kunijiban by the Public Works Research Institute, Geo Station from the National ResearchInstitute for Earth Science and Disaster Resilience, the Research Information Database from the National Instituteof Advanced Industrial Science and Technology, and the Nation-wide Electronic Geotechnical Database Systems
by the Japanese Geotechnical Society.
To build the shallow soil layers model, the following data are often used: the Digital Japan Basic Map andFundamental Geospatial Data (provided online by the Geospatial Information Authority of Japan), geology
information that can be confirmed using the 1:200,000 Seamless Digital Geological Map of Japan (multiple piecesof geological information are provided online at Geomap Navi from the National Institute of Advanced Industrial
Science and Technology), and the “7.5-arc-second Japan Engineering Geomorphologic Classification Map basedon the uniform national standards” (Wakamatsu and Matsuoka, 2013), in addition to the information used in the
deep sedimentary layers model.
Procedure (2): Estimation of representative columnar section of each mesh and layer classification
A borehole columnar section is estimated to represent each mesh. When a mesh contains multiple boring data,
the data with the deepest borehole depth or boring data with PS logging are chosen for “representative columnar
section” of the mesh. Alternatively, the representative columnar section can be determined from soil types andthe mean of N-values7 at each depth of the boring data in a mesh, after removing low-quality boring data. For thisrepresentative columnar section, layers are classified considering the geologic stratigraphy, soil classification,
and N-values available in the literature.
In a mesh classified as “mountain/hill,” it is necessary to set weathering layers and covering layers at a depthof several meters from the surface, because of the low S-wave velocities near the surface resulting from theweathering and stress release in the shallow part of the bedrock and an existence of loam layers. When it ispossible to estimate the representative columnar section, the classifications for the weathering layers and thecovering layers are determined considering the depth variation of the N-value distribution.
The soil type, thickness, and N-value vary in horizontal directions, even within a mesh, because of the positionsof the available boring data. In some cases, the geomorphologic unit and geologic stratigraphy may be also
different in a mesh. Such variation in the soil types and the N-values is caused by their spatial extent within amesh, and one can develop the representative columnar section from the averages of these values (Yasuda et al.,2009).
7This value is used as an index of ground strength, and it is defined as the number of drops of a hammer required to make aStandard Penetration Test sampler penetrate to a distance of 300 mm. It is obtained with in-situ test to determine
engineering characteristics using a soil sample (Standard Penetration Test).

10

Furthermore, when the geomorphologic
N-valueClassification of S-wave velocity layers
unit, geologic stratigraphy, and layer
Covering layers/strong weathering layers :<10classification in boreholes within a mesh are
150-200m/s
different or when layer thicknesses vary
10- 50
Weak weathering layers : 200-350m/s
significantly in boreholes within a mesh,one may select only highly reliable boring
Fresh layers :data as representative, or boring data whose>50
S-wave velocity of top of engineering bedrock
geomorphologic
unit
orgeologic
Figure 4-2. Examples of shallow soil classifications
stratigraphy occupies the largest area in theof mountains/hills.
mesh. The layers belonging to theweathering classification, including the weathering layers, are classified in the SD model with N-values fromboring data obtained in “mountain/hill” regions of the Kanto district and nationwide PS logging data (Figure 4-2).
In a mesh having the representative columnar section that belongs to the “mountain/hill” classification, three layers
are prepared: “covering layer/strong weathering part (N-values less than 10),” “weak weathering part (N-values
between 10 and 50),” and “fresh part (N-values of 50 or higher).” If there are no boring data, the layers areclassified considering their spatial continuity.
Procedure (3): Building of one-dimensional velocity structure model using correlation of physical property values
with N-values and soil classification
P-wave and S-wave velocities obtained by PS logging are assigned to the layer classifications made inProcedure (2). When no PS logging data are available in a mesh, A velocity value is estimated to each layer using
an empirical equation for N-value and physical values (mainly S-wave velocity) for each soil type and eachgeological age, and thereby a one-dimensional velocity structure model is built for each mesh.
Many regression formulas for S-wave velocities from various types of soil data have been proposed in previous
investigations. These include regression formulas with explanatory variables of N-values and soil types (Imai andTonouchi, 1982) and N-value, soil types, geological ages, and depths (Ohta and Goto, 1976; Masaki, 1984;
Fukuwa et al., 1999).
One can build a velocity structure model by using an averaged velocity for each classification based on the
representative columnar section without estimating different velocity structures from N-values and soil types foreach mesh. In the SD model, an S-wave velocity for each layer classification is assigned referring a regression
formula of N-values, soil types, and geological ages by the Central Disaster Prevention Council (2001) when the
representative columnar section is not based on PS logging data. For the weathering classification of the“mountain/hill” part, S-wave velocities of 150–200 and 200–350 m/s are used for the “covering layer/strong
weathering part” and “weak weathering part,” respectively, as shown in Figure 4-2. The S-wave velocity of the"fresh part" is set to the same value as the surface of the engineering bedrock. In the method by the Central
Disaster Prevention Council (2001), the weathering layer with an S-velocity of 300 m/s is set with a thickness of
10 m for Quaternary volcanoes and Tertiary sedimentary layers and 5 m for the other layers.
Procedure (4): Modeling for meshes without boring data (with weak emphasis of layer continuity on building athree-dimensional model)
For a mesh without boring data, one assigns a one-dimensional velocity structure model of the representative
columnar section of one of the neighboring meshes that has the same geomorphologic classification as thegeomorphologic features of the largest area in the mesh. When the mesh is surrounded by multiple meshes
belonging to the same geomorphologic classification in “plateau/lowland” part, a velocity model of the closest
mesh is used for that mesh.
For the mesh belonging to the “mountain/hill” part, one also uses a velocity model of a mesh with the same
geomorphologic classification. However, the surface geological condition can be also considered, because the11thickness of the weathering layers varies depending on the surface geology.
In a region without any boring data, the same one-dimensional model may spread widely. However, theresultant three-dimensional model becomes close to the one-dimensional multi-layered model because the layer
structure models are different for the geomorphologic classifications of the meshes. Although the accuracy of thethree-dimensional model depends on the number of the boring data, this approach is widely used because of theeasy modeling work (the left panel of Figure 4-3).
Procedure (5): Modeling for meshes without boring data (with strong emphasis on layer continuity for buildinga three-dimensional model)
The depths of each layer boundary of the one-dimensional velocity structure models made in Procedure (3) are
interpolated for meshes without any boring data using the information of neighboring meshes with boring data.
In this step, the spatial continuity of the layers is considered referring to the topographical and geologicalcharacteristics.
Because an important boundary of the physical properties can be identified using data from geophysicalsurveys for the shallow soil layers, it is important to perform additional microtremor observations and geophysicalsurveys at locations without boring data. Microtremor observations can be used to determine an S-wave velocity
structure to tune the layer thicknesses. The data from the geophysical surveys can provide spatial information tointerpolate the shape of the layer boundaries of the one-dimensional velocity structure models into the horizontal
directions (to develop a three-dimensional velocity model). Moreover, for an interpolation of the one-dimensional
velocity structures, one can calculate interpolated values in the entire target area with a kriging method, an inverse
distance weighting interpolation method, and a spline function method as a surface model, using altitude valuesof the layer boundaries that have been set in the representative columnar section for each mesh. It is necessary toconsider sufficiently the spatial distribution ranges of the layer interfaces and their regional boundaries. For
example, the distribution of the alluvium layers is limited within the lowland regions, and the distribution of theloam layers is limited within the regions belonging to the plateaus. Therefore, a significant irregularity in the topsurface of the engineering bedrock may be generated in the regional boundary of the areas between “mountain/hill”
and “plateau/lowland” because of the geological differences. It is necessary to express the shape of the connecting
part of these areas with a consideration of an accuracy corresponding to the adopted mesh size.
Two approaches are available to interpolate a velocity value for each layer classification in meshes withoutboring data. One is a method in which the velocities are obtained directly by an interpolation of the velocity
values of the corresponding layer classifications in the representative columnar sections of the neighboringmeshes. In the other method, an S-wave velocity is estimated from an empirical formula between a soil
classification and a N-value which is interpolated from a soil types and N-values at the neighboring meshes(Kimura et al., 2014).
Because the SD model is established to reproduce an observed broadband ground motion, the shallow soil
layers model is built using this advanced three-dimensional approach (the right panel of Figure 4-3). Therefore,a large number of the boring data were collected, and many observations of microtremors were also conductedto obtain the local variation of the layer boundary shapes accurately.

12

Figure 4-3. Differences in approaches for shallow modeling in a mesh with and without boring data: approach with weakemphasis of layer continuity on building a three-dimensional model (left), which generates a block-shaped model having
constant depths to layer boundaries in meshes belonging to same geomorphologic classification, and approach with strongemphasis of layer continuity on building a three-dimensional model (right), in which the depths to boundaries are modeledwith consideration of the spatial information for smooth connection of the layer classifications of neighboring meshes
without boring data.
5. Building of “SD model”
A SD model from a combination of the shallow and the deep layers is built for evaluating broadband groundmotion using the shallow soil layers model and the deep sedimentary layers model. The effects of the nonlinearity
amplification in the shallow soil layers are not considered in the present ground motion prediction using the SDmodel. The shallow part of the SD model is used to calculate the linear amplification. The evaluation of the site
amplification characteristics with methods including the nonlinear soil effects is necessary to consider thevariation of earthquake ground motion features caused by the nonlinear behavior of the soil.
When the shallow soil layers model is built using data in boring that penetrates, in general, to a depth of thestratum with an N-value of 50, one can set the top surface of the engineering bedrock with an N-value of 50 andan S-wave velocity of 300 to 400 m/s as the bottom layer of the shallow soil layers model, referring to the Nvalues. However, the top of the deep sedimentary layers model is often set to be near a stratum with an S-wave
velocity of 500 to 700 m/s considering various types of information. Here, the procedure to build the SD modelis explained. By this procedure the two models explained in Chapters 3 and 4 are connected without an artificial
discontinuity between them, and the combined model is improved with the results of earthquake ground motion
observations and microtremor surveys for the SD model.
Procedure (1): Connection of the shallow soil layers model and the deep sedimentary layers model and its adjustment
The shallow soil layers model is connected to the deep sedimentary layers model directly when the differenceof the S-wave velocities (or the impedance contrast) is small between the bottom of the shallow soil layers modeland the top of the deep sedimentary layers model. When the difference is large, intermediate layers having
intermediate S-wave velocities of 400 to 500 m/s are added between the two models to prevent an artificial
discontinuity as a result of the connection. The thicknesses and the velocities of the additional intermediate layers
are set by referring to existing data, such as PS logging data from deep boring surveys. An artificial discontinuity
should not be also generated three-dimensionally with an adjustment of the thicknesses of the additional layers
considering the regional tendency of the depths to the deep sedimentary layers and the top of the seismic basement
in a wide region including a basin and surrounding mountain areas.

13

The intermediate layers are set with a careful consideration of the three-dimensional shapes of the interfaces
and thicknesses of the connecting part between the shallow soil layers model and the deep sedimentary layers
model (Figure 5-1).

Figure 5-1. Image of setting additional layers with intermediate velocities.
The three-dimensional shapes of the top surfaces of the engineering bedrock and the seismic basement
are considered in setting the thicknesses and the velocities of the additional layers with intermediate Swave velocities by referring to PS logging data from deep boring surveys.
Procedure (2): Adjustment of the SD model using earthquake records and microtremor data
Earthquake ground motion records and microtremor data are collected to adjust the SD model made inProcedure (1). The SD model is tuned using information on seismic velocity structures beneath the observation
points extracted from the earthquake records and microtremor data. If there is difficulty in tuning the model, onemust collect additional data and perform microtremor observations.
The thicknesses and the velocities of the layers in the SD model are tuned, including the additional
intermediate layers applied in Procedure (1). One can reconsider the modeling of the connecting part if necessary.
To extract information on seismic velocity structures from earthquake records and microtremor observation
data, various analysis methods are available, such as H/V spectrum, which is the Fourier spectrum ratio ofhorizontal motion to vertical motion of earthquake records (called the “R/V spectrum” when applied to theRayleigh wave part), H/V spectrum of microtremor data, spectral ratio of a seismic record on the surface to aborehole record in a vertical array , and surface wave phase velocity from a microtremor array survey. For
example, one can determine an S-wave velocity and layer thickness from comparison of the peak period andshape of the observed H/V spectrum of microtremor data at a single point with the theoretical ellipticity ofRayleigh waves with an assumption of a homogeneous layered velocity structure. It is also often applied todetermine an S-wave velocity and layer thickness with an inverse analysis of surface-wave phase velocity
obtained from microtremor array records. A joint inversion analysis in which multiple observed data aresimultaneously analyzed is also an effective method for the model tuning.

14

Procedure (3): Validation and adjustment of the SD model with short-period empirical site effect
After estimating site effects at short periods (for example, periods of less than 2 s) by means of an empiricalmethod using ground motion records from moderate to small earthquakes, the estimated site effects are comparedwith the calculated ones using a theoretical method to validate the SD model adjusted in Procedure (2). When
there is difficulty in reproducing the site effects, it is necessary to return to Procedure (2) for readjustment.
Spectral separation inversion (for example, Iwata and Irikura, 1986) can be used as an evaluation method forthe empirical site effects at short periods. Q-values for the soil layers and the propagation path can also bevalidated (for example, Fukushima and Midorikawa, 1994; Yamanaka et al., 2009). The one-dimensional elastic
wave propagation theory is mainly used for the one-dimensional theoretical evaluation of the site effects.
Procedure (4): Comparison with reproduction of observed earthquake ground motion using a three-dimensional
calculation method
The validated or adjusted SD model in Procedure (3) is used to calculate earthquake ground motions formoderate to small earthquakes by means of a theoretical method based on the elastic wave theory, such as thefinite difference or finite element methods, considering the influence of the three-dimensional seismic velocity
structure. Then, the overall characteristics of the observed and calculated waveforms, such as P-S travel times,
amplitudes, peak periods, spectral shapes, time-variant characteristics, and durations, are compared. If there is adifficulty in the reproducibility, one must return to Procedure (2) to adjust the SD model again.
Depending on the purpose of ground motion simulation, a comparison of ground motions on the surface of theearth and the top-surface of the engineering bedrock is performed. For the comparison of the motions on the topsurface of the engineering bedrock, the shallow soil layers are deleted from the SD model established in Procedure(3), and the modified model is used in the theoretical calculations. Besides, earthquake ground motion on thesurface of the engineering bedrock is estimated from the observed waveforms using correction factors calculated
from the one-dimensional soil amplification of the shallow soil layers model that has been removed in the model
used in the three-dimensional simulation.
It may be possible to perform the validations simultaneously in Procedures (3) and (4) to compare thebroadband ground motion calculated by hybrid methods.

15Appendices1. Previous developments in building of seismic velocity structure models
The Headquarters for Earthquake Research Promotion summarized the procedures for building a seismic velocity structure modelby referring to Koketsu et al. (2009), Fujiwara et al. (2009), and so on. These procedures were published in the strong groundmotion prediction method for earthquakes with specified source faults (“Recipe”). The sedimentary layer model for the wide
area made with the “Recipe” was published in J-SHIS by the National Research Institute for Earth Science and Disaster
Resilience in 2010. Furthermore, the “Japan Integrated Velocity Structure Model (provisional version)” and its background idea
for the model building for long-period ground motion predictions were released in “Long-period ground motion hazard map”
(prototype version in 2012).
2. Engineering bedrock
This is defined as a firm soil layer that is used as the position for setting predicted earthquake ground motion in building designin construction engineering fields, such as building engineering and civil engineering. Its S-wave velocity varies depending on
the kind of man-made structures and soil conditions. However, it is within a range of 300 to 700 m/s in many cases(Architectural Institute of Japan, 2009). Furthermore, the stratum with an N-value of more than 50 is often considered as the topof the deep sedimentary layers. Careful attention must be paid to determine how deep the surface of the engineering bedrock should
be set, because a stratum with an N-value of less than 50 may appear below the above-defined engineering bedrock.
3. Seismic reflection and refraction surveys
The seismic survey methods are used to explore seismic velocity structures from surface observation of ground vibration
generated by a controlled seismic source (or an artificial seismic source). Physical phenomena of elastic seismic wave propagating
to the boundary of geological layers with different physical properties, such as refraction and reflection are utilized in the methods.
There are two survey methods. A seismic refraction survey captures the wave propagation velocities and attenuation factors
mainly by using refracted waves, whereas a seismic reflection survey captures the distribution and shapes of reflection interfaces
using reflected waves.
4. Deep boring surveys
MITI exploratory test wells were conducted to discover oil and gas fields on land. There are also examples of deep boring
surveys down to 3000 m in the Kashiwazaki deep ground motion observation project. In sea areas, b o r i n g surveys have been
conducted with the research vessel “Chikyu” (which drilled down to approximately 2200 m in the Nankai Trough) as an example.
5. PS logging
Velocity logging is a survey method in which a velocity structure is deduced by utilizing a borehole to measure the elastic
waves propagating in soil. The elastic wave velocities that can be obtained include those for P-waves and S-waves. When both
velocities are obtained, it is called “PS logging.”
6. First-degree modelThis is a seismic velocity structure model for an assumed Tokai earthquake and Tonankai earthquake (“Long-period groundmotion hazard map,” Headquarters for Earthquake Research Promotion, 2009); “ Japan Integrated Velocity Structure Model(provisional version)”; ( “Long- period ground motion hazard map,” Headquarters for Earthquake Research Promotion, 2012;
Koketsu et al., 2008, 2012), which can be found on the website of the Headquarters for Earthquake Research Promotion(https://www.jishin.go.jp/evaluation/seismic_hazard_map/lpshm/).
7. Important notice about ground motion simulations for moderate to small earthquakes
Even if one uses a velocity model revised by repeated verification through the reproductions of earthquake ground motions ofmoderate to small earthquakes, discrepancies between the observed motions and the calculated ones may be found for anearthquake that was not used for the verification. In such cases, one must be careful with a consideration of non-neglectable effects
due to the hypocenter parameters, the crustal structures and anisotropy of the propagation path from the hypocenter of each
earthquake to an observation point as factors that are not related with the characteristics of the shallow soil layers and the deep
sedimentary layers.
8. Empirical equations for AVS30 (mean S-wave velocity from the surface of the earth to a depth of 30 m)
The followings are examples of the empirical equations for an AVS30 and a site amplification factor.
・ Matsuoka and Midorikawa (1994): A mean S-wave velocity of soil layers was empirically obtained from Digital National
Land Information data which include the topographical classifications and altitudes in the Kanto region. The empiricalequation is mainly applicable to the Kanto region.
・ Midorikawa and Matsuoka (1995): The relationship between an AVS30 and an amplification factor for the peak ground
velocities of strong earthquake records was obtained. Then, an AVS30 was generated from the empirical equation withthe topography, altitude and distance from a river on the basis of the topographic information, soil data and so on.
・ Matsuoka et al. (2005): An empirical method for the AVS30s was proposed with a consideration of altitudes, slope angles,
and distances from mountains and hills formed in the old geological ages, considering the significant differences of the
AVS30s in each geomorphologic classification caused by the soil formation processes and sedimentary environments.
The empirical equation was applied for the 7.5-arc-second Japan Engineering Geomorphologic Classification Map.

16

References (alphabetical order)
Arai, H. and Tokimatsu, K. (2004): S-wave velocity profiling by inversion of microtremor H/V spectra, Bull.
Seism. Soc. Am., 94, 1, 53-63.
Architectural Institute of Japan (2009): Generation guide for seismic input motions based on the recent
advancement of ground motion studies, Maruzen Inc., 163 (in Japanese).
Baba, T., Ito, A., Kaneda, Y., Hayakawa, T. and Furumura, T. (2006): 3-D seismic wave velocity structures in
the Nankai and Japan Trench subduction zones derived from marine seismic surveys, Abstracts for Japan
Geoscience Union Meeting, S111-006.
Central Disaster Management Council (2001): Estimate of the mean S-wave velocity of a ground using Digital
National Land Information, “Expert examination committee for Tokai Earthquake” (the 11th), material 31, related materials (about earthquake ground motion), 1-33 (in Japanese).
Fujimoto, K. and Midorikawa, S. (2006): Relationship between average shear-wave velocity and site
amplification inferred from strong motion records at nearby station pairs, Journal of Japan Association
for Earthquake Engineering, 6, 1, 11-22 (in Japanese).
Fujiwara, H., Kawai, S., Aoi, S., Senna, S., Ooi, M., Matsuyama, H., Iwamoto, K., Suzuki, H. and Hayakawa, Y.
(2006): A subsurface structure modeling in all over Japan for strong-motion evaluation, The 12th Japan
Earthquake Engineering Symposium, 1466-1469 (in Japanese).
Fujiwara, H., Kawai, S., Aoi, S., Morikawa, N., Senna, S., Kudo, N., Ooi, M., Hao, K. X., Wakamatsu, K.,
Ishikawa, Y., Okumura, T., Ishii, T., Matsushima, S., Hayakawa, Y., Toyama, N. and Narita, A. (2009):
Technical Reports on National Seismic Hazard Maps for Japan, Technical Note of the National ResearchInstitute for Earth Science and Disaster Prevention, 336.
Fukushima, Y. and Midorikawa, S. (1994): Evaluation of site amplification factors based on averagecharacteristics of frequency dependent Q-1 of sedimentary strata, Journal of Structural and Construction
Engineering, AIJ, 460, 37-46 (in Japanese).
Fukuwa, N., Arakawa, M., Koide, E. and Ishida, E. (1999): Dynamic soil modeling based on existing soil datain urban area using GIS, AIJ Journal of Technology and Design, 9, 249-254 (in Japanese).
Gardner, G. H. F., Gardner, L. W. and Gregory A. R. (1974): Formation velocity and density - The diagnostic
basics for stratigraphic traps, Geophysics, 39, 6, 770-780.
Gassmann, F. (1951): Über die elastizität poröser medien, Vierteljahrss-chrift der Naturforschenden Gesellschaft
in Zurich, 96, 1-23.
Geological Survey of Japan, edited (2013): Japan Gravity Database, DVD version, Digital geological map, P-2,
Geological Survey of Japan at National Institute of Advanced Industrial Science and Technology.
Honda, R., Sawada, A., Furuse, N., Kudo, T., Tanaka, T. and Hiramatsu, Y. (2012): Release of gravity databaseof the Kanazawa University, Journal of the Geodetic Society of Japan, 58, 4, 153-160 (in Japanese).
Imai, T. and Tonouchi, K. (1982): Relationship between the N-value and the S-wave velocity and its usage
examples, Journal of Foundation Engineering and Equipment, 16, 6, 70-76 (in Japanese).
Iwata, T. and Irikura, K. (1986): Separation of source, propagation and site effects from observed S-waves,
Journal of the Seismological Society of Japan. 2nd ser., 39, 579-593 (in Japanese).
Kimura, T. (2002): Crustal movements in the Japanese Archipelago – from a new viewpoint-, Aichi Shuppan,
470 (in Japanese).
Kimura, K., Hanashima, Y., Nishiyama, S. and Ishihara, Y. (2014): Shallow-level subsurface three-dimensional
geological modeling based on analysis of borehole data: A case of the northern Tokyo and the southern
Nakagawa lowland area, Geological Survey of Japan, Miscellaneous Map Series, No. 40, “Integral
research and application of the subsurface geology of the central Kanto Plain,” 56-113 (in Japanese).
Koketsu, K., Miyake, H., Fujiwara, H. and Hashimoto, T. (2008): Progress towards a Japan integrated velocity
structure model and long-period ground motion hazard map, Proceedings of the 14th World Conference17

on Earthquake Engineering, Paper No. S10-038.
Koketsu, K., Miyake, H., Afnimar and Tanaka Y. (2009): A proposal for a standard procedure of modeling 3-D
velocity structures and its application to the Tokyo metropolitan area, Japan, Tectonophysics, 472, 290300.
Koketsu, K., Miyake, H. and Suzuki, H. (2012): Japan Integrated Velocity Structure Model Version 1,
Proceedings of the 15th World Conference on Earthquake Engineering, Paper No. 1773.
Kono, Y. and Furuse N. (1989): 1:1,000,000 scale gravity anomaly map in the Japanese Archipelago, Universityof Tokyo Press, 88.
Langston, C. A. (1979): Structure under Mount Rainier, Washington, inferred from teleseismic body waves, J.
Geophys. Res., 84(B9), 4749-4762.
Ludwig, W. J., Nafe, J. E. and Drake, C. L. (1970): Seismic refraction, in The Sea, vol. 4, edited by A. E. Maxwell,pp. 53-84, Wiley-Interscience, New York.
Masaki, K. (1984): A study about evaluation of ground motion characteristics for earthquake disaster prevention,
Doctoral dissertation at Tokyo Institute of Technology, 167 (in Japanese).
Matsubara, M., Obara, K. and Kasahara, K. (2008): Three-dimensional P- and S-wave velocity structures beneath
the Japan Islands obtained by high-density seismic stations by seismic tomography, Tectonophysics, 454,86-103.
Matsubara, M. and Obara, K. (2011): The 2011 Off the Pacific Coast of Tohoku earthquake related to a strong
velocity gradient with the Pacific plate, Earth, Planets and Space, 63, 663-667.
Matsuoka, M. and Midorikawa, S. (1994): The Digital National Land Information and seismic microzoning,
Reference materials of the 22nd Ground Motion Symposium, 23-34 (in Japanese).
Matsuoka, M., Wakamatsu, K., Fujimoto, K. and Midorikawa, S. (2005): Average shear-wave velocity mapping
using Japan engineering geomorphologic classification map, Transactions of Japan Society of Civil
Engineers, 794/I-72, 239-251 (in Japanese).
Matsuoka, M. and Wakamatsu, K. (2008): Site amplification capability map based on the 7.5-arc-second Japan
Engineering Geomorphologic Classification Map, National Institute of Advanced Industrial Science andTechnology, Intellectual property management, H20PRO-936.
Midorikawa, S. and Matsuoka, M. (1995): GIS-based integrated seismic hazard evaluation using the Digital
National Land Information, Butsuri-Tansa, 48, 6, 519-529 (in Japanese).
Ohta, Y. and Goto, N. (1976): Estimation of S-wave velocity in terms of characteristic indices of soil, ButsuriTanko (Geophysical Exploration), 29, 4, 251-261 (in Japanese).
Ryoki, K. (1999): Three-dimensional depth structure of the crust and uppermost mantle beneath SouthwesternJapan and its regional gravity anomalies, Journal of the Seismological Society of Japan. 2nd ser., 52, 5163 (in Japanese).
Sato, H., Hirata, H., Koketsu, K., Okaya, D., Abe, S., Kobayashi, R., Matsubara, M., Iwasaki, T., Ito, T., Ikawa,T., Kawanaka, T., Kasahara, K. and Harder, S. (2005): Earthquake source fault beneath Tokyo, Science,
309, 462-464.
Suzuki, H. (1996): Geology of the Koto deep borehole observatory and geological structure beneath themetropolitan area, Japan, Report of the National Research Institute for Earth Science and Disaster
Prevention, 56, 77-123 (in Japanese).
Wakamatsu, K. and Matsuoka, M. (2003): Development of the 7.5-arc-second Japan engineering
geomorphologic classification map for metropolitan areas, JSCE journal of Earthquake Engineering, 27,050.pdf, CD-ROM (in Japanese).
Wakamatsu, K., Matsuoka, M., Kubo, S., Hasegawa, K. and Sugiura, M. (2004): Development of GIS-based
Japan engineering geomorphologic classification map, Transactions of Japan Society of Civil
Engineering, .759/I-67, 213-232 (in Japanese).
18

Wakamatsu, K. and Matsuoka, M. (2013): Development of the 7.5-arc-second Japan engineering
geomorphologic classification map based on uniform national standards and the use of the map, Bulletinof Japan Association for Earthquake Engineering, 18, 33-38 (in Japanese).
Yamamoto, A. and Shichi, R. (2004): Gravity atlas of the Japanese Archipelago – Southwest Japan and Central
Japan -, University of Tokyo Press, 136, with CD-ROM.
Yamanaka, H., Ohori, M. and Midorikawa, S. (2009): Estimation of seamless Vs structure for Kanto basin usingsite amplification from spectral inversion of earthquake motion records, Summaries of technical papersof Annual Meeting Architectural Institute of Japan. B-2, Structures II, Vibration, 47-48 (in Japanese).
Yasuda, S., Todo, H., Mimura, M. and Yamamoto, K. (2009): Developing collaborative linkage scheme among
geotechnical databases throughout Japan, Abstracts for the 3rd Symposium “Development of Integrated
Geophysical and Geological Information Databases,” 49-58 (in Japanese).
Zhao, D., Hasegawa, A. and Kanamori, H. (1994): Deep structure of Japan subduction zone as derived from local,
regional, and teleseismic events, J. Geophys. Res., 99(B11), 22313-22329.

19