Journal of Natural Disaster Science, Volume 23, Number 2, 2001, pp89-100

Design Spectra and Phase Spectrum Modeling to Simulate Design Earthquake Motions:
A Case study through Design Standerds of Railway Facilities in Japan
Tadanobu SATO 1), Yoshitaka MURONO 2), Hai-Bo WANG 2), Akihiko NISHIMURA 2)
1)

Disaster Prevention Research Institute, Kyoto University, Kyoto, Japan2)
Railway Technical Research Institute, Tokyo, Japan(Received for 14 Mar., 2002)

ABSTRACT
Design earthquake motions often are defined by response spectra. We introduce three types of design spectra
given in the Seismic Design standards for Railway Facilities (1999). Both the amplitude and phase characteristics of earthquake motion provide key information needed to evaluate the aseismic capacity of designed structures. We present a simple method for modeling the phase characteristic of earthquake motion based on groupdelay time and wavelet analysis. First the phase spectrum is modeled near the source region. Second we develop
a method to model the phase spectrum for the inter-plate earthquakes. Existing data sets on observed earthquake
motions and wavelet analysis then are used to calculate the group delay time of each earthquake motion on eachcompact support of Meyer’s analyzing wavelet. Regression equations for the mean group delay time and standard
deviation on each compact support are obtained as functions of the earthquake’s magnitude and epicentral distance. Simulating sample phase spectra earthquake motions compatible with the design response spectra then are
simulated, and their characteristics are investigated using ductility demand spectra.

1. INTRODUCTION
Because of the severe damage done to urban facilities by the
1995 Hyogoken Nambu earthquake, many Japanese earthquake
design standards have been revised. After two years of activity, acommittee brought together to revise the aseismic design standards
for railway facilities issued Design Standards for Railway
Facilities: commentaries on aseismic design . After introducing
three types of design response spectra defined in the DesignStandards for Japanese Railway Facilities we here develop methods to simulate phase spectra to be used for generating design
earthquake motions.
2. EARTHQUAKE MOTION AT THE ENGINEERINGBASE LAYEREarthquakes that occur the near Japanese archipelago are classified as intra- and inter-plate earthquakes. The Pacific andPhilippine plates move, respectively, in the west and north-west
directions, and both are subducted along the Pacific periphery of
the Japanese archipelago. The associated tectonic forces cause
magnitude 8 class earthquakes. The return period of earthquake
occurrence in the same area ranges from several tens to two hundred years. Recent inter-plate earthquakes are shown in Figure 1.
In contrast, the average return period of an intra-plate earthquake having a magnitude larger than 7 at a certain active fault isseveral thousands years. If, however, such an event occurs near adensely populated urban area severe damage is done to structures.
Moreover, the average occurrence rate of intra-plate earthquakes ofa magnitude larger than 6.8 within the Japanese archipelago isabout 12 years (Rika Nenpyo 2001) 1). The intensities of intraplate earthquakes therefore were taken into account when establishing the design standards, as were those of inter-plate earthquakes (JSCE 1996) 2).
When establishing design standards, the design earthquake
intensity must be defined through risk analyses at the design point.
Although an intensive survey of inland active faults was madeafter the 1995 Hyogoken Nambu earthquake, accuracy in estimating the return period of earthquakes caused by those faults is notgood enough to use for design standards. Therefore, the positionsof all the detected inland active faults only have been used todefine design earthquake intensity in the new design standards. A
flow chart in Figure 2 shows how design earthquake motion isdefined. This motion is defined at the engineering base groundwhere the N values of the standard penetration test are more than
50 or the shear wave velocity is more than 350 m/sec. The intensity of the design earthquake motion is defined by the response spectrum.
Spectrum I is obtained by taking into account the inter-plate
fault activities that cause magnitude 8 class earthquakes. Because
the epicentral distances of these earthquakes are more than 40 km,
the defined spectrum level also envelopes spectra of the nearsource earthquake motions caused by M6.5 class intra-plate earthquakes. This is the minimum intensity level required in the seismic design standards.
Spectrum II is obtained from regression analysis based on
observed earthquake motions caused by intra-plate earthquakes.
Methods to simulate strong earthquake motions based on theoretical fault models have become popular in the field of civil engineering. The spectrum obtained by theoretically taking into account

KEY WORDS: design earthquake motions, phase spectrum, group delay time, design spectra, regression analysis, spectrum compatible earthquake
motion, ductility demand spectra

90

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

Fig. 1

Fig. 2

Recent interplate earthquakes in Japan

Flowchart showing design seismic motion determination

the fault mechanism is called spectrum III.
Areas in which there is high probability of existing active
faults whose detection is not easy due to deep alluvial and diluvial
deposits and metropolitan areas with complicated crustal structures
are classified as areas of high uncertainty of fault existence. In
these areas the design earthquake intensity is defined by spectrum
II modified by a regional correction factor.
3. DESIGN RESPONSE SPECTRA3)
3.1 Spectrum II (for intra-plate earthquakes)(1) Near Source Earthquake Motion Records
Spectrum II is defined stochastically by the use of response
spectra calculated from observed near source earthquake motions.
We collected observed earthquake motions recorded on hardground to define design response spectrum free from the effect of
topological condition and ground irregularity as well as amplification of soft surface layers. Table 1 gives near source earthquake
motions with maximum accelerations of more than 100gal
observed on hard ground during the 1995 Hyogoken Nambu earth-quake (M7.2), Coyote Lake earthquake (1979, M5.9), Loma Prieta
earthquake (1989, M7.1), Landers earthquake (1992, M7.5) andNorthridge earthquake (1994, M6.7). Although collected earthquake motions are mainly from the earthquakes occurred in the
United States except the 1995 Hyogoken Nambu earthquake weassumed that there is no distinguishable difference between the
U.S. and Japanese earthquake motions. Their acceleration
response spectra are shown in Figure 3.
The observed records at Tarzana during the Northridge earthquake and those at the Kobe Power Transfer Station, as well as at
the Takarazuka Station of the Japan Railway Company, during the
1995 Hyogoken Nambu earthquake were not used because effects
of topological conditions were reported to be large. Acceleration
response spectra calculated from these observed earthquake
motions are shown by dotted lines in Figure 3.
Spectra scattering is very extensive. It originates from uncertainty about the 1) source mechanism, 2) transmitting path, 3) geographical structure and 4) local soil condition. To reduce scattering, the calculated response spectra levels are corrected by the proposed attenuation relationship.

DESIGN SPECTRA AND PHASE SPECTRUM MODELING TO SIMULATE DESIGN EARTHQUAKE MOTIONS:A CASE STUDY THROUGH DESIGN STANDERDS OF RAILWAY FACILITIES IN JAPAN
Table 1.

Fig. 3Near-source seismic records from recent earthquakes

Acceleration response spectra of observed records at near source regions ofinland earthquakes(2)

Response Spectra Adjustment Using the Attenuation
Relation
There are several attenuation relationships of the acceleration
response spectrum that take into account the earthquake source
extent. An attenuation relationship that uses the shortest distanceto the earthquake fault, proposed by Fukushima (1994) 4), is givenbylogS(T) = a1(T)M 2w - a2(T)Mw + b(T)・R-log(R+0.025×100.42Mw) +Σcj(T)Ij(1)
An attenuation relationship that uses the equivalent hypocentral
distance was proposed by Ohno(1994) 5);logS(T) = a(T)M w - logXeq -b(T)Xeq + c(T) +Δs(T)N

Xeq-2 = Σ di2Xi-2
i=1

91

(2)

N

di2
Σ
i=1

in which T is the period, Mw the moment magnitude, R the shortest
distance to the earthquake fault, Ij the correction factor, Xeq theequivalent hypocentral distance, N the number of subdivided elements of the earthquake fault, and di and Xi respectively are the dislocation on a subdivided element and the hypocentral distance tothat element.
The shortest distance from an observation point to the earthquake fault and the equivalent hypocentral distance are shown in
Table 1. We used fault models obtained from the publications ofWang et al. (1997) 6) to calculate these values.
Inland earthquake fault planes usually are vertical, the fault
bottoms being 15-20km from the ground surface, depending on the
thickness of the earth’s crust. An increase in earthquake magnitude increases the length of the earthquake fault. The intensity ofearthquake ground motion near a fault region, however, is not proportional to earthquake magnitude because the effect of rupture
energy released on a fault plane distant from the observation point
decreases with an increase in the distance to the rupture point. The
affected area becomes large as earthquake magnitude increases,
but ground motion intensity asymptotes to a constant value forpoints with the same shortest distance to the fault in the nearsource region. To clarify the characteristics of earthquake motion
just above a fault, the calculated response spectra from observed
earthquake motions were adjusted to those at an equivalent
hypocentral distance, 12km, from the fault (Figure 4). The averageadjusted response spectra level becomes high, but almost all these
spectra are bounded by 2000 gal, except those calculated using

92

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

Table 2.

Fig. 4

Seismic records from recent interplate earthquakes in Japan

Adjusted response spectra of observed records to the equivalent hypocentral distance of 12km

earthquake record numbers 4 and 11. Scattering of the response
spectra in the short period range is decreased as compared with theoriginal spectra (Figure 3), but in the long period range there is no
improvement. Earthquake motion characteristics in the long period
range are controlled by the source mechanism and crustal profile
surrounding the earthquake fault.
The acceleration response spectra adjusted to the shortest fault
distance of 2km also were examined. In intra-plate earthquakes inJapan, the main earthquake energy-releasing portion of the fault isdeeper than 2km. This shortest distance therefore means that theobservation point is located just above the fault.
(3) Design Response Spectra
A design response spectrum should be defined by a probabilistic base if there is enough data for stochastic analysis. To
define non-exceedance probability the following items must be
taken into account:1) The safety of railway facilities must be evaluated from thestandpoint of protecting the lives of passengers.
2) Damage to a structure causes malfunctioning of the entire railsystem.

In the new design standard, 90% of the non-exceedance probability was used to define the design acceleration spectrum, on theassumption that there was a lognormal distribution of the calculated response spectra for each period. The dotted line in Figure 4
shows the acceleration response spectrum with a non-exceedance
probability of 90%. On the basis of this line, the design acceleration spectrum is defined by the thick polygonal line shown in
Figure 5. We treat this as the mother spectrum defined just above
the inland earthquake fault and modify its level based on distancefrom the fault.
3.2 Spectrum I (inter-plate earthquakes)
The response spectrum level that takes into account inter-plate
earthquake activities was 1G for the elastic response of structures.
This spectrum was revised in the same manner used to define
Spectrum II.
(1) Earthquake motions analyzed
Recent inter-plate earthquakes of magnitudes larger than 7.4
that occurred near the Japanese archipelago are shown in Figure 1.
The records collected consist of 100 components with maximum

DESIGN SPECTRA AND PHASE SPECTRUM MODELING TO SIMULATE DESIGN EARTHQUAKE MOTIONS:A CASE STUDY THROUGH DESIGN STANDERDS OF RAILWAY FACILITIES IN JAPANFig. 5Comparison of acceleration response spectra, design spectrum, average spectrum
and spectrum with 90% non-exceedance probability.

Fig. 6

Distance distribution between the observation point and seismic source

Fig. 7

Acceleration response spectra of observed records of Japanese interplate earthquakes

93

94

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

accelerations of more than 100 gal. Distributions of hypocenters in
the observed records are shown in Figure 6. Most data are clustered at hypocentral distances of 100-200 km.
We selected earthquake motions that satisfy the following
conditions for statistical analysis;(a) An equivalent hypocentral or shortest distance to the earthquake fault of less than 200km(b) No conspicuous peaks in the Fourier amplitude spectra.
(c) Ground condition at the observation points classified as hardground (of which shear wave velocity is from 350m/sec to
400m/sec).
Twenty-seven records satisfied these conditions. Acceleration
response spectra calculated from these records are shown in
Figure7.
(2) Adjustment of Calculated Response Spectra
The distance to the earthquake fault was used to adjust the calculated response spectra. The concepts of the equivalent hypocentral and shortest distances were used to modify the response spectra. To define Spectrum I, inter-plate earthquakes of magnitudeclass 8 were targeted. All records of magnitudes less than 8 therefore were adjusted to magnitude 8. The distances to the earthquake fault in Table 3 were used to adjust the response spectra.
The average spectra calculated from the adjusted ones usingeqs.(1) and (2) are shown in Figure 8 for cases in which the shortest distance to the earthquake fault was 30km and the equivalent
hypocentral distance 80km. The average spectrum for an equivalent hypocentral distance of 80km is larger than that for the shortest distance of 30km. The 80 km equivalent hypocentral distancetherefore is comparable to the shortest distance of 20km. The zone

Table 3.

with this shortest distance to a magnitude 8 class inter-plate earthquake fault is located just above the periphery of the fault. In the
Kanto and Tokai regions inter-plate earthquake sources are locatedjust below these regions, but the area within the shortest distanceof 20km is not wide. In the other regions this distance is classified
as being very close area to the earthquake fault.
From this information the design response spectrum that takes
into account inter-plate earthquake activities (Spectrum I) is
defined by the broken line in Figure 8. In the period range longer
than 1.0sec, the spectrum level coincides with Spectrum II for theconvenience of design. This spectrum envelopes the adjusted
spectra in the period range longer than 0.65 sec but not in the
shorter period range. Spectra calculated from the observed motion(Figure 7), however, are enveloped by this spectrum.
This Spectrum I multiplied by the regional modification factorof 0.7 also envelopes the spectra calculated by eqs.(1) and (2) atjust above an earthquake fault for an earthquake of magnitude 6.5
and source depth of 10km.
4. MODELING OF THE PHASE SPECTRA
Evaluation of the seismic capacities of structures designed
using the defined design spectra, requires compatible earthquake
motions, and simulation of earthquake motions necessitates bothmodeling of amplitude and phase characteristics. We propose two
methods to define the phase characteristics of design earthquake
motion.
The concept of group delay time tgr (ω) was used to model thephase spectra of earthquake motion. Group delay time is a derivative of the phase spectrum with respect to circular frequency(Papoulis 1962) 7).

Target distance to the earthquake fault and the earthquake scale used for adjustment

tgr(ω)=

dφ
dω

(3)

Its mean value is the central arrival time of ground motion with thefrequency component ω, and its standard deviation corresponds tothe duration of the ground motion concerned. Modeling the groupdelay time is much easier than direct modeling of the phase.

Fig. 8Comparison of adjusted acceleration response spectra

DESIGN SPECTRA AND PHASE SPECTRUM MODELING TO SIMULATE DESIGN EARTHQUAKE MOTIONS:A CASE STUDY THROUGH DESIGN STANDERDS OF RAILWAY FACILITIES IN JAPAN

4.1 Modeling of Group Delay Time for Design Spectrum II(1) Modeling of the Rupture Process of an Earthquake Fault
A fault model of the Haskel type is assumed. The fault plane
is divided into n×n elements, each of which corresponds to thearea of a small event. The motion of the large event g0L (t) at theobservation point is expressed by taking into account the timedelay of each small event, g0S (t) (Irikura 1983) 8);nL nW (nD-1)n' 1

nL nW

gOL(t) = ΣΣ gos (t-tij) + ΣΣ Σ
i=1 j=1i=1 j=1 k=1n'(Rij -R0) ξij
+
,
Vr
Vs

tijk = tij +kτ(nD -1)n'(5)where τ is the rise time; Vs the shear wave velocity of the media;
Vr the rupture velocity;ξij the distance between the starting pointof rupture and the small event point; nL, nW, nD the respective similarity ratios for fault length, width and rise time between the large
and small events; Rij the distance between the representative pointof the small event and the observation point; and n’ an arbitrary
integer to eliminate any artificial period caused by subdivision of
the rise time.
Deconvoluting eq.(4) with respect to the source time function
and train of impulses, p(t), expressing the rupture process of earthquake fault gives
nL nW

nL nW (nD -1)n'i=1 j=1i=1 j=1

Σ n'1 aij δ(t-tijk)p(t) = ΣΣ aij δ(t-tij)+ΣΣ
N

(6)

i=1

N =nL×nW×(nD -1) n'
Because the Fourier amplitude of eq.(6) at zero frequency is proportional to the seismic moment of the earthquake fault, for thecase of multiple events the value a ij is assigned as proportional to
the seismic moment of each fault. The effect of geometrical damping factor 1/Rij is considered to evaluate a ij. The phase spectrum
φp(ω) of p(t) is

(

)

(7)

The group delay time tgrp (ω) of φ(ω) is obtained by differentiating
it with respect to the circular frequency, ω;Ntgrp (ω) =1 ∞ 1n(Ah(y))
dyπ ∫-∞ ω-y(9)

Differentiating eq.(9) with respect to ω, the group delay time
tgrh (ω) is1
1 ∞ 1n(Ah(y))
dy = 2 *{1n(Ah(ω))}π ∫-∞ (ω-y)2ω(10)

in which ＊ is the symhol expressing convolution.
The phase spectrum of the frequency transfer function fromthe base rock to ground surface is expressed by the sum of the minimum phase function and the all pass function (Izumi 1990) 9). The
group delay time of the all pass function nearly coincides with thetravel time of the direct wave from the base rock to ground surface.
The effect of the group delay time of the all path function therefore
is just to shift the arrival time of earthquake motion, and this effectwas neglected in the simulation of earthquake motions compatible
with design response spectra.
(3) Linearity of Group Delay TimeIf an earthquake motion is expressed by the convolution of
two functions (a function of source and transmitting path, p(t), andone of local amplification, h(t), the phase of earthquake motion,φ(ω), is the linear sum of these functions’ phases, φp (ω) andφh (ω). Because the linearity of the group delay time also is preserved, the group delay time of earthquake motion t gr (ω) isexpressed as the linear sum of each function’s group delay timek=1

= Σ aiδ(t-ti)

a sin(ωti )φp (ω) = tan-1 Σ i
Σai cos(ωti )φh (ω)=-(4)in whichtij =spectra;tgrh (ω)=
gos (t-tijk)

95

tgr (ω) = tgrp(ω) + tgrh (ω)(11)(4) Simulation of Earthquake Motion Compatible with Design
Spectrum II
The Fourier amplitude spectra calculated from the observed
earthquake motions and the phase spectra obtained by integratinggroup delay times simulated by using eq.(11) taking into accountthe fault rupture process and local site condition were used to resimulate the earthquake motions of the 1995 Hyogoken Nambu
earthquake. The validity of eq.(11) for estimating the phase spectra was verified by comparing the resimulated earthquake motions
with the observed ones (Sato et al. 1999) 10).
Design earthquake motion compatible with Spectrum II therefore is simulated by using eq.(11) to define phase spectra. The
assumed earthquake fault is a vertical one with the size of 40 ×

N-1 N

-Σai2 ti -Σ Σ ai aj (tj+ti)cos{ω(tj-ti)}
i=1

N

Σ
i=1

i=1 j=i+1
N-1 N

ai2 +2ΣΣ ai aj cos{ω(ti-ti)}
i=1 j=i+1(8)(2) Modeling of the phase shift caused by local soil conditions
The amplitude characteristic of the frequency transfer functionthat expresses the effect of local soil conditions is assumed apriori
to be Ah (ω). We applied the minimum phase criterion to define its
phase spectrum,φh (ω). This assumption is not valid when theseismic wave takes multiple paths between the observation point
and base rock. The phase spectrum that expresses the local soil
effect is derived from the Hilbert transformation (Papoulis 1962)5)
which defines the relationship between amplitude and the phaseFig. 9a

Example of the rupture fault model used in the analyses

96

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

Fig. 9b

Relative position between the observation points and fault(a)Fig. 11

Fig. 10Group delay time at point A(b) Group delay time at point B
Simulated group delay times at points A and B

Ground profile used to evaluate local soil effect20(km). In the parameter studies we selected several combinationsof the starting point of fault rupture and distribution of asperity
because these values strongly affect the phase characteristics ofearthquake motion. We considered three types of asperity distributions and three cases of rupture starting points. The size and distribution of asperity were defined from seismological research results(Somerville 1993) 11). From the parameter studies, we selected the
combination of asperity distributions and starting points of fault
rupture, shown in Figure 9a, as the most severe fault rupture model
to generate design phase spectrum. Fourteen observation points
around the fault (Figure 9b) are considered to define the design
earthquake motion compatible with design Spectrum II. The procedure used to simulate earthquake motion compatible with thedesign response spectrum is given in the Appendix. The depth of
the base rock is assumed to be 1500m and the shear wave velocity
Vs=3400(m/s). The amplitude spectrum expressing local soil condition Ah(ω) is calculated by the multi reflection theory of one
dimensional SH wave. Ground profile used for the calculation isshown in Figure 10.
Representative results of the parameter analyses are shown in
Figures 11 and 12, in which the group delay times and simulated
earthquake motions compatible with Spectrum II respectively are
given. Point A is located in the direction of rupture propagation,

Fig. 12

Earthquake motions compatible with design spectrum IIand point B in the opposite direction. Scattering of the group delay
time of the train of impulses is small for the observation point
located in the direction of rupture propagation (point A), increasing in the opposite direction (point B). The directivity of rupture
propagation shows clearly in the group delay time. This phenomenon was confirmed in the simulated earthquake motions (Figure12). Even when the same response spectrum is used, the waveform is strongly affected by the relative position of the observation
point and the rupture configuration. If the observation point islocated in the direction of rupture propagation, wave energy is concentrated over a short period. In the opposite direction, however,the duration of earthquake motion lengthens.
The ductility demand spectrum was used to analyze the effectof these earthquake motion characteristics on the dynamic analysesof structure systems. This spectrum was obtained through nonlin-

DESIGN SPECTRA AND PHASE SPECTRUM MODELING TO SIMULATE DESIGN EARTHQUAKE MOTIONS:A CASE STUDY THROUGH DESIGN STANDERDS OF RAILWAY FACILITIES IN JAPAN97

Fourier transformation ofψj,k (t) has a compact support for eachscale factor j (named the jth compact support) defined by{2 j / 3Td < f < 2 j+2/ 3Td}

Fig. 13Strength demand spectra for earthquake motions compatible with design spectrum II (μ=4)

ear dynamic response analyses of a single degree of freedom system. The yield seismic coefficient of the SDF is controlled until
the SDF response reaches the targeted ductility. The chart that
expresses the relationship between the yield seismic coefficientand elastic natural period is called the ductility demand spectrum,in which the targeted ductility is a parameter. The hysteresis
model used in the analyses is the bi-linear type proposed byClough (1966). The second stiffness was selected as 5% of thefirst, and a 5% damping coefficient assigned. Results for a ductility 4 are shown in Figure 13. The thin lines are ductility demand
spectra calculated from simulated earthquake motions compatible
with Spectrum II using the phase spectra of twenty observed earthquake motions. The ductility demand spectrum for the earthquake
motion at point A envelopes all twenty results. Because this spectrum has a very smooth nature in all the frequency ranges we designated this earthquake motion as the design earthquake motion
compatible with Spectrum II.
4.2 Modeling of Group Delay Time for Design Spectrum I(1) Average Group Delay Time and Its Standard Deviation
Wavelet analysis often is used to determine the nonstational
characteristic of earthquake motion in the time and frequency
domains because resolution in these domains is guaranteed by theuncertainty criterion.
The discrete wavelet transformation of function x(t) and its
inverse transformation have been defined by (Yamaguchi, 1990)12);
aj,k =∫-∞ ψ*j,k・x (t) dt(12)x (t) = Σ x (j) (t) = ΣΣ aj,k・ψj,k (t)(13)ψj,k (t) = 2j/2・ψ(2jt - k) j,k : integer

(14)

∞

j

j

k

in which the value with the sign * expresses the complex conjugateof the original value, ψ(t) the wavelet function, j is the scale factor, and k the factor expressing the time shift. If the inverse
wavelet transformation exists,ψ(t) is called the analyzing wavelet.
Although there are several ways to define this function, we used
the method of Meyer (Meyer 1989) 13) to composeψ(t). The(15)

in which Td (=NΔt, N: number of data, Δt: sampling interval) is
the duration of earthquake motion. As seen from eq.(15), if scale
factor j becomes large, the treating frequency range shifts to thehigh frequency region.
For all oberved earthquake earthquake motions with a sampling interval of 0.01(sec) zero data were added until the total
number of sampling data for each earthquake motion, N, became
131071 (=217). The wavelet transformation of each earthquake
motion x(t) was used to decompose each time history of earthquake motion to a component time history of each scale factor j ( j=1〜17. We call this the jth component time history, x (j)(t). Basedon the definition of x (j)(t) given in eq.(13) the frequency band of
this time history is limited, as defined by eq.(15).
The group delay time of this time history, tgr(j)(ω), was calculated and its mean, μ(j)tgr , and standard deviation, σ(j)tgr , on the jthcompact support obtained from eqs.(16) and (17):N (j) t (j) (ω)
gr
i
(j)
i=1 N

(j)
μtgr
=Σ

(16)

(j)

(j)
σtgr
=

1 N (j)(j) 2
tgr (ωi) -μtgr)
N (j) Σ(i=1

(17)

in which j=1,2…,17, N(j) is the number of time history data at thejth scale factor defined by N(j)=2j , and tgr(j)(ωi) the group delay timeof the jth component time history at a circular frequency ofωidefined by equation (18);tgr(j) (ωi)=φ(ωi) -φ(ωi+1)
dφ
=dω
Δω

(18)

(2) Regression Equations of Group Delay Time
Existing data sets for the observed earthquake motions areused to calculate the group delay time of each earthquake motion.
The attenuation relationships of the mean and standard deviation of
the group delay time on each compact support are obtained asfunctions of earthquake magnitude and epicentral distance.
Data sets of the earthquake motions analyzed are from the Off
Southwest Hokkaido earthquake (1993, M7.8), the Off East
Hokkaido earthquake (1994, M8.1), the Off Sanriku earthquake(1994, M7.5), and the Hyogoken Nambu earthquake (1995, M7.2).
The starting time of rupture at the earthquake fault was selected as the origin of the time history of earthquake motion. The theoretically calculated travel time of the S wave from the initial rupture point to the observation site is added to the trigger time. In
standard seismic design procedures the period range concerned is
0.1〜5 sec. Regression equations μ(j)tgr and σ(j)tgr for j=6〜15
therefore were developed as functions of earthquake magnitude M
and epicentral distance Δ (km);μ(j)tgr =α1(j)×10β1M×Δγ1 , σ(j)tgr =α2(j)×10β2M×Δγ2(19)in which α1 (j), α2 (j), β1 (j), and β2 (j) are regression coefficients on thejth support.

98

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

Fig. 14

Table 4.

Relationship between the epicentral distance Δ and the mean and standard deviations of group delay time on scale factor j=10,13

Regression analyses resultscorrelation coeff

α1

α2

β1

β2

γ1

γ2

μ

σ

7

1.011

27.708

0.0

0.0

0.864

0.203

0.86

0.41

8

0.830

14.584

0.040

0.0

0.790

0.337

0.96

0.71

9

0.543

17.9680.086-0.030

0.700

0.344

0.85

0.75

10

0.806

8.451

0.060

-0.005

0.686

0.321

0.97

0.73

11

0.850

2.970

0.026

0.016

0.764

0.366

0.98

0.77

12

0.511

0.392

0.058

0.143

0.744

0.295

0.99

0.83

13

0.367

0.0790

0.077

0.267

0.739

0.201

0.99

0.87

14

0.330

0.0572

0.081

0.287

0.742

0.239

0.99

0.81

15

0.04430.0106

0.256

0.439

0.591

0.204

0.97

0.79

j

Fig. 16

Fig. 15

Earthquake motions compatible with design spectrum I

Ductility demand spectra for earthquake motions shown in Figure 15

DESIGN SPECTRA AND PHASE SPECTRUM MODELING TO SIMULATE DESIGN EARTHQUAKE MOTIONS:A CASE STUDY THROUGH DESIGN STANDERDS OF RAILWAY FACILITIES IN JAPAN99

Fig. 17 The design earthquake motion compatible with design spectrum I

Regression analyses results are given in Table 4, and some areshown in Figure 14. Values of μ(j)tgr and σ(j)tgr increase as the epicentral distance increases, especially in the longer period range(for a small j value) because α1 and α2 values are increases as the jvalue decreases. The effect of epicenter distance on μ(j)tgr andσ(j)tgr is almost same because γ1 and γ2 have constant values forall j values. The effect of earthquake magnitude on μ(j)tgr and σ(j)tgr
is rather small comparing with that of epicenter distance because
values of β1 and β2 are small. The values of β1 and β2 are related
to shorter period motions because those values increase as the jvalue increases but the values of α1 and α2 are related to longer
period motions.
(3) Simulation of Earthquake Motion Compatible with Design
Spectrum I
The values of μtgr(j) and σtgr(j) on the jth compact support wereobtained from eq.(19) for an earthquake of magnitude 8.0 andthree epicentral distances (10, 100, 300km) used to simulate earthquake motions compatible with design Spectrum I (Figure 8).
Because the ferquency range of each compact support j defined byeq.(15) is overlapping with that of neighboring compact supports
j−l and j×l, we assume that the frequency range 2 j−l /T/ < f <2 j/T is effective to generate the group delay time in the compact
support j. A sample group delay time on the jth compact support
for each epicentral distance is simulated by generating a random
value based on a normal distribution N(μtgr(j),σtgr(j)) although
Nigam (1982)14) derived a probability density function for groupdelay time different from the Gauss distribution. The phase spectrum on the jth support is obtained by integrating this sample groupdelay time. To guarantee the continuity of phase spectrum the initial value of phase in the frequency range 2 j−l /T < f < 2 j/T, whichis the representative frequency range of compact support j, isassigned as the final phase of the frequency range 2 j−2 / T < f <2 j−l /T which is the representative frequency range of the compact
support j−l. Simulated earthquake motions compatible withdesign Spectrum I (see Appendix) are shown in Figure 15. This
design spectrum is defined by taking into account only the activityof inter-plate earthquakes. The epicentral distance therefore ismore than 30-40km, but we added a case of a 10km epicentral distance. As the epicentral distance increases the duration lengthens,
and longer period motion predominates.
Ductility demand response spectra (Figure 16) were calculated
from the simulated earthquake motions. Although the same elastic
response spectrum was used to simulate the earthquake motions,the strong effect of epicentral distance can be seen on the elastoplastic structural response characteristics. A short epicentral distance requires large ductility, whereas the demand decreases as theepicentral distance increases.

For design purpose it is not practical to use a sample earthquake motion compatible with Spectrum I because the time history
characteristic of a simulated earthquake motion is strongly affectedby the random nature of a generated group delay time. We therefore simulated several tens of design earthquake motions compatible with Spectrum I by supposing earthquake ground motions withmagnitude M=8 at the epicentral distance Δ=80(km). From thesean earthquake motion (Fig. 17) was selected as the standard earthquake motion to be used for the design purpose because it has asmooth shape of ductility demand spectra.
5. CONCLUSION
We introduced the concept of design earthquake motions as
defined in the design standards for Japanese railway facilities, inparticular the procedures used to define design acceleration spectra
and to simulate earthquake motions compatible with the design
spectra.
The following findings were obatined:
1) Three types of design acceleration spectra were defined in
the seismic design standards for Japanese Railway Facilities.
2) The concept of group delay time was used to model thephase spectra of earthquake motions. A method to simulate
phase spectra of earthquake motions at near source regions wasdeveloped by taking into account the effect of rupture processof earthquake faults and local soil conditions.
3) From data sets of observed earthquake motions and thewavelet algorithm, the mean group delay time and standard
deviation on the compact support of the Meyer wavelet wereexpressed as functions of earthquake magnitude and epicentral
distance (regression equations). We developed a method tosimulate a sample phase spectrum for an inter-plate earthquake
motion based on the regression equations.
4) Using the simulated sample phase spectrum a method tosimulate earthquake motions compatible with the design
response spectrum was developed. Its suitability for design
purposes was investigated using ductility demand spectra.
APPENDIX
The method to simulate earthquake motions compatible withdesign spectra is summarized. The simulated phase spectrum on
the jth compact support φ(j)(ω) is used to simulate an earthquake
motion compatible with the response spectrum. To do this, an
initial Fourier amplitude spectrum Aint(ω) is needed. This is
assumed to be given by the velocity response spectrum with azero damping constant converted from the design acceleration

100

T. SATO, Y. MURONO, H. WANG AND A. NISHIMURA

response spectrum given by
Aint (ω) =β+γ)ωS (ω,0.05),
( 1+α
h

(20)

in which S(ω,0.05) is the targeted design spectrum with 5%
damping, h the damping constant, and α,β and γ are given parameters (Kawashima and Aizawa, 1986) 15)
The initial Fourier amplitude of the jth compact support
A(j)(ω) is assumed to coincide with the initial Fourier amplitude
Aint(ω) in the frequency range {2j-1/Td ≦f≦2j /Td} and to be zero
outside of this region; i.e., in the frequency regions {2j/3Td≦
f<2j-1/Td} and {2j /Td <f≦2j+2/3Td}. A jth component time history
can be simulated using φ(j)(ω) and A(j)(ω). Earthquake motion of
the first step is obtained by summing up all the component time
histories. Because the calculated response spectrum based on firststep earthquake motion does not coincide with the design response
spectrum, we modified the amplitude of the initial Fourier amplitude, Aint(ω), at the circular frequency, ω, by multiplying the ratior(ω) of the design response spectrum and calculated response spectrum. Assuming the calculated r(ω) Aint(ω) as the new initial
Fourier amplitude of earthquake motion and using the original
phase spectrum of the jth component φ (j)(ω), the earthquake
motion of the second step can be simulated. Repetition of thisprocess until the ratio r(ω) nearly equals 1.0, gives simulated
earthquake motion compatible with the design response spectrum.
REFERENCES
1) Rika Nenpyo (Chronological Scientific Tables 2001),pp.800-831,
National Astronomical Observatory, Japan. Maruzen Co., Ltd., 2000
2) JSCE, ‘Proposal on earthquake resistance for civil engineering structures’, Special task committee of earthquake resistance of civil engineering structures, Japan, (1996)
3) Ministry of Transportation and Railway Research Tustitute, ‘Seismic
Design Standards for Railway Facitities’, Marzen, (1999), (inJapanese).

4) Y. Fukushima, ‘Empirical prediction for strong ground motion reflectedon theoretical backgrounds of source and propagation of seismic
wave’, ORI Report 93-07, Ohsaki Research Institute, (1994), (inJapanese).
5) S. Ohno and K. Takahashi, ‘Evaluation of strong-motion attenuation
relation ssing near-source data in California’, Proc. of the 9th Japan
Earthquake Engineering Symposium, (1994), (in Japanese).
6) H. Wang, Y. Murono and A. Nisimura, ‘On evaluation of ground motion
near seismic faults based on strong seismic records’, Proc. of the 2nd
Urban Earthquake Disaster Symposium, (1997), (in Japanese).
7) A. Papoulis, ‘ The Fourier integral and its application’, McGRAW-Hill,(1962).
8) K. Irikura, ‘Semi-empirical estimation of strong ground motions during
a large earthquake’, Bulletin of Disaster Prevention ResearchInstitute, Kyoto University, Vol.33, Part2, No.297, 63-104, (1983).
9) M. Izum, S. Kurita, Y. Endo, J. Tobita and T. Hanzawa, ‘ Study on
causality of transfer functions and components of causal transfer functions in systems of seismic wave propagation in soil’, J. of Structure
Construction Engineering, AIJ., No.412, 31-41, (1990), (in Japanese)10) T. Sato, Y. Murono and A. Nishimura, ‘Modeling of phase characteristics of strong earthquake motion’, J. Struct. Mech. Earthquake Eng.,
JSCE, No.612/I-46, 201-213, (1999), (in Japanese)
11) P. G. Somerville, K. Irikura, S. Sawada, Y. Iwasaki, M. Tai and M.
Fushimi, ‘A study on empirical modeling of slip distribution on
faults’, Proc. of the 22th JSCE Earhtquake Engineering Symposium,
291-294, (1993).
12) M. Yamaguchi and M. Yamada, ‘Wavelet analysis’, Science, Vol.60
No.6, 398-405, (1990) (in Japanese).
13) Y. Meyer, ‘Wavelets and operator’, Cambridge University Press,(1992).
14) N. C. Nigam, ‘Phase properties of a class of random process’
Earthguake Engineering and Structural Dynamics, Vol. 10, 711-714,(1982).
15) K.Kawashima, K. Aizawa and K. Takahashi, ‘Attenuation of peakground motion and absolute acceleration response spectra’ Proc.of 8thWorld Conference on Earthquake Engineering, Vol.II, 257-264,(1984).