University of Southern California Dissertations and Theses

Determining fault zone structure and examining earthquake early warning signals using large datasets of seismograms

(USC Thesis Other)

Determining fault zone structure and examining earthquake early warning signals using large datasets of seismograms

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content DETERMINING FAULT ZONE STRUCTURE AND EXAMINING EARTHQUAKE
EARLY WARNING SIGNALS USING LARGE DATASETS OF SEISMOGRAMS

by

Michael Antony Lewis

A Dissertation Presented to the
FACULTY OF THE GRADUATE SCHOOL
UNIVERSTIY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(EARTH SCIENCES)

DECEMBER 2008

Copyright 2008 Michael Antony Lewis

ii
Acknowledgments
Many figures in the paper are generated using GMT (Wessel & Smith 1995).
For Chapter 1 I would like to thank Jennifer Eakins and other staff members of the
Anza seismic network for their hard work to make the data available. The manuscript
benefited from useful comments by two anonymous referees and editor Rob van der
Hilst (GJI). The work was funded by the National Science Foundation (grant EAR-
0409048).
In Chapter 2 I thank the IRIS Data Management Center and the Northern
California Earthquake Data Center for distributing the employed dataset. The
manuscript benefited from useful comments by the GJI editor J. Trampert and three
anonymous referees.
Chapter 3 benefitted from discussions with Frank Scherbaum who also provided
the parameters for removing the FIR filter artifacts from the data. The study was funded
by the National Science Foundation (grant EAR-0409605). The paper benefited from
constructive comments by the GJI editor J. Trampert and two anonymous referees.
For Chapter 4 I need to thank ISS international and Margaret Boettcher for
making the waveform data available. The paper benefited from useful discussions with
Margaret Boettcher and comments from Richard Allen, Michel Bouchon and an
anonymous associated editor (JGR). The study was funded by the National Science
Foundation (grant EAR-0409605).

iii
Table of contents
Acknowledgments

ii
List of figures

v
Abstract

viii
Introduction

1
Chapter 1: Shallow Seismic Trapping Structure in the San Jacinto Fault Zone
Near Anza, California
4
Summary 4
1-1: Introduction 5
1-2: Analysis 8
1-3: Discussion 31
1-4: Chapter 1 References

36
Chapter 2: Imaging the Deep Structure of the San Andreas Fault South of
Hollister With Joint Analysis of Fault-Zone Head and Direct P
Arrivals
42
Summary 42
2-1: Introduction 43
2-2: Method 54
2-3: Results 64
2-4: Discussion 73
2-5: Chapter 2 References

81
Chapter 3: Examination of Scaling Between Proposed Early Signals in P
Waveforms and Earthquake Magnitudes
85
Summary 85
3-1: Introduction 86
3-2: Data 90
3-3: Methods 96

3-4 Results 102
3-5: Discussion 107
3-6: Chapter 3 References

112
Chapter 4: Examination of Scaling Between Earthquake Magnitude and
Proposed Early Signals in P Waveforms From Very Near-Source
Stations in South African Gold Mines
115
Summary 115
4-1: Introduction 116
4-2: Data and Methods 120
4-3: Results 131

iv
4-4: Discussion 139
4-5: Chapter 4 References

145
Recapitulation

149
References 150

v
List of Figures
Figure 1-1: Location map
6
Figure 1-2: Spectral ratio
10
Figure 1-3: Example waveforms array CC
11
Figure 1-4: Example waveforms array CL
12
Figure 1-5: Example waveforms array BR
13
Figure 1-6: Event Quality Maps
15
Figure 1-7: Ratio horizontal components
17
Figure 1-8: Moveout with depth
19
Figure 1-9: Moveout with depth for array CL and BR.
20
Figure 1-10: Trapped wave quality array CC
21
Figure 1-11: Trapped wave quality array CL
22
Figure 1-12: Travel time moveout
24
Figure 1-13: Waveform inversion CC
27
Figure 1-14: Waveform inversion CL
28
Figure 1-15: Waveform inversion BR
29
Figure 2-1: Fault model
46
Figure 2-2: Synthetic head waves
47
Figure 2-3: Location map
49
Figure 2-4: Seismicity in a around a surface trace complexity.
52
Figure 2-5: Vertical component velocity seismograms recorded at station SUM
55
Figure 2-6: Differential arrival time analysis
56

vi
Figure 2-7: Optimizing the inversion
60
Figure 2-8: Inversion of synthetic data
62
Figure 2-9: Depth averaged velocity contrasts
63
Figure 2-10: Inversion Results
65
Figure 2-11: Variations in contrast with distance from the fault.
67
Figure 2-12: Inversion using events from northwest
69
Figure 2-13: Inversion using events from southeast
70
Figure 2-14: Waveform fitting
72
Figure 2-15: A schematic diagram of the inferred velocity structure in the study
area
74
Figure 3-1: A location map of the study area
88
Figure 3-2: Similar waveforms
91
Figure 3-3: FIR artifacts
94
Figure 3-4: Example seismograms before the removal of FIR filters
95
Figure 3-5: Same as 3-4 but with the removal of acausal response of FIR filters
96
Figure 3-6: A seismogram and the calculated values of τ
p
98
Figure 3-7: An Example seismogram with an initial low amplitude arrival
101
Figure 3-8: calculation of t
sip
102
Figure 3-9: Calculated signals from early in the P waveforms
104
Figure 3-10: similar to figure 3-9 for data that has been processed to remove FIR
filter artifacts
105
Figure 3-11: Calculating corrected Pd
106
Figure 3-12: Normalized velocity seismograms of the 5 similar earthquakes of
cluster 30 109

vii

Figure 4-1: Radial component normalized velocity seismograms at station K1
118
Figure 4-2: Radial component normalized velocity seismograms at 5 station for a
magnitude 1.5 event
121
Figure 4-3: Obtaining signals
122
Figure 4-4: Calculating corrected Pd
126
Figure 4-5: Dependence of τ
pmax
with time window length
128
Figure 4-6: Similar to figure 4-5 for a smaller event
129
Figure 4-7: Plots of the 6 investigated parameters as a function of magnitude
132
Figure 4-8: Variety of P wave onset style
133
Figure 4-9: The effects of time window change on τ
pmax
136
Figure 4-10: The effects of time window change on corrected Pd
137
Figure 4-11: Comparison between estimated rupture duration and time interval
139
Figure 4-12: A compilation of estimated rupture duration divided by time window

143

viii
Abstract
Seismic signals associated with near-fault waveforms are examined to determine
fault zone structure and scaling of earthquake properties with event magnitude. The
subsurface structure of faults is explored using fault zone head and/or trapped waves,
while various signals from the early parts of seismograms are investigated to find out
the extent to which they scale with magnitude.
Fault zone trapped waves are observed in three arrays of instruments across
segments of the San Jacinto fault. Similarly to previous fault zone trapped wave studies,
the low velocity damage zones are found to be 100-200m wide and extend to a depth of
~3-5km. Observation and modeling indicate that the damage zone was asymmetric
around the fault trace. A similar sense of damage asymmetry was observed using detailed
geological mapping by Dor et al. (2006) nearby on the San Jacinto fault at Anza. Travel time
analysis and arrival time inversions of fault zone head waves were used to produce high
resolution images of the fault structure of the San Andreas fault south of Hollister. The contrast
of P wave velocities across the fault was found to be ~50% in the shallow section,
lowering to 10-20% below 3 km, with the southwest side having faster velocities.
Inversions making use of different subsets of stations suggest that a low velocity damage zone
also exists in this area and that it is more prominent on the faster velocity side of the fault. The
patterns of damage from these studies of fault zone head waves and trapped waves are
consistent (Ben-Zion and Shi, 2005) with the theoretical prediction that earthquake ruptures on
these fault sections have statistically-preferred propagation directions.
The early parts of P waveforms are examined for signals that have previously
been proposed to scale with the final event magnitude. Data from Turkey and a deep

ix
South African gold mine show that scaling is present in signals related to the maximum
displacement amplitude and frequency content. The high sampling rate of the
instruments in the gold mine enables the reduction of the time window in which
measurements are made to below the estimated rupture duration of the largest events.
Using increasingly small time windows has only a minimal effect on the scaling of the
signals with event magnitude, implying that the size of earthquakes is affected
statistically by some property of the early part of the rupture.

1
Introduction
As the number and quality of seismometers increases and generates more near
fault and near-source data, this facilitates the improved investigation of properties of the
structure of faults and properties of earthquake rupture. The purpose of this thesis is to
employ large seismic datasets in near-fault regions to 1) image fault zone structures
using fault zone head and trapped waves and 2) examine earthquake early warning
signals to gain insight into the properties of rupture initiation.
This work looks at seismic signals generated by reflections and refractions along
and within the internal structure of large faults. The two fault zone phases discussed in
this work are fault zone trapped waves (FZTW) and fault zone head waves (FZHW).
Unlike the P or S waves which are traditionally used to image the earth, these fault zone
phases are generated by interaction with the fault structure and as such can be used to
glean high resolution details of that structure. The properties of fault zones may have
ramifications for earthquake hazard during large earthquakes, act as a record of the long
term rupture properties of large earthquakes in the area and have consequences for
many aspects of earthquake physics. For example, understanding where differing
material contrasts exist may predict a preferred rupture direction for future large
earthquakes (e.g. Weertman 1980 Ben-Zion 2001).
This thesis also investigates a number of signals that were previously proposed
to scale with the final event size and what these can tell us about the nature of
earthquake rupture. These signals include some that were purported to relate to
nucleation process and others that were suggested for possible use in earthquake early

2
warning systems. The nature of the onset of earthquake rupture and whether the rupture
process is intrinsically predictable is a fundamental question within seismology.
The four chapters of this thesis are each a published paper. Chapter 1 (Lewis et
al., 2005) images the structure of three branches of the San Jacinto fault using FZTW.
These phases are generated by the total internal reflection of seismic waves within a low
velocity damage zone associated with the fault. In this chapter observation and
waveform inversion are used to determine the physical properties of the damage zones
associated with the branches of the faults across which the arrays of instruments are
deployed.
In Chapter 2 (Lewis et al., 2007) the 1-D velocity structure on each side of the
San Andreas Fault is obtained using joint arrival time inversion of direct P waves and
FZHW. Further variations in the velocities along strike and away from the fault are also
obtained by using subsets of stations and events. This illustrates that when the fault is
also the interface between two materials, head waves that refract along the fault can
detail the near fault velocity structure better than direct arrivals alone
Chapter 3 (Lewis and Ben-Zion, 2007) examines seven different signals from
early parts the seismograms that were proposed to scale with the final event magnitude
and tests whether artifacts introduced in the data processing and recording of
seismograms can impact these signals. The signals can be divided into two types: first,
those which were reported to reflect signatures of seismic nucleation phases and second,
those suggested for use in earthquake early warning systems. The data used in this
chapter consists primarily of clusters of repeating earthquakes that allow for the

3
performance of analysis without the mixing of receiver site and path effects that is
present in all the previous work.
Finally, Chapter 4 (Lewis et al., 2004) examines whether the final event
magnitude can be calculated from a time window in the seismogram that ends before
the cessation of rupture. The same methods and techniques are employed as in Chapter
3 but using high sampling rate data from a deep South African gold mine. This provides
minimal distance between the source and receiver and the ability to reduce the time
window in which observations are made to below the estimated rupture duration.

4
Chapter 1: Shallow Seismic Trapping Structure in the San Jacinto
Fault Zone Near Anza, California
Summary
We analyze fault zone trapped waves, generated by ~500 small earthquakes, for
high-resolution imaging of the subsurface structure of the Coyote Creek, Clark Valley,
and Buck Ridge branches of the San Jacinto fault zone near Anza, California. Based on
a small number of selected trapped waves within this dataset, a previous study
concluded on the existence of a low-velocity waveguide that is continuous to a depth of
15-20 km. In contrast, our systematic analysis of the larger dataset indicates a shallow
trapping structure that extends only to a depth of 3 to 5 km. This is based on the
following lines of evidence. (1) Earthquakes clearly outside these fault branches
generate fault zone trapped waves that are recorded by stations within the fault zones.
(2) A travel time analysis of the difference between the direct S arrivals and trapped
wave groups shows no systematic increase (moveout) with increasing hypocentral
distance or event depth. Estimates based on the observed average moveout values
indicate that the propagation distances within the low velocity fault zone layers are 3-5
km. (3) Quantitative waveform inversions of trapped wave data indicate similar short
propagation distances within the low velocity fault zone layers. The results are
compatible with recent inferences on shallow trapping structures along several other
faults and rupture zones. The waveform inversions also indicate that the shallow
trapping structures are offset to the northeast from the surface trace of each fault branch.
This may result from a preferred propagation direction of large earthquake ruptures on
the San Jacinto fault.

5
1-1. Introduction
The San Jacinto Fault Zone (SJFZ) is one of the most active fault zones in
southern California (e.g. Sanders & Kanamori 1984). It extends ~230 km southeastward
from the San Gabriel Mountains to the Salton Trough at its southern extent (inset of
Figure 1-1). The cumulative slip across the SJFZ is estimated at between 24 km (Sharp
1967) and 29 km (Scott et al. 1994) over about 2 million years. Large faults that have
accommodated such significant slip typically generate zones of intensely damaged
material with lower seismic velocity than the surrounding rocks (Ben-Zion & Sammis
2003, and references therein). If spatially persistent and sufficiently uniform, the low
velocity fault zone (FZ) rock can act as a waveguide for seismic energy. The
constructive interference of critically reflected phases, which are generated within a low
velocity layer and arrive after the S waves, are generally termed FZ trapped waves (e.g.
Ben-Zion 1998; Li & Leary 1990; Ben-Zion & Aki 1990). Since the amplitude and
frequency content of the trapped waves depend strongly on the properties of the FZ
waveguide within which they propagate, FZ trapped waves can provide high resolution
information on geometrical and seismic properties of the FZ waveguide at depth.
Various previous FZ trapped wave studies relied on a small number of selected
waveforms in drawing conclusions on approximately 100 m wide low-velocity FZ
layers that extend to the base of the seismogenic zone (e.g. Li & Leary 1990; Li et al.
1990; Li et al. 1997). Li & Vernon (2001) recorded about 1500 small earthquakes on
FZ arrays located across the Buck Ridge (BR), Clark Valley (CL) and the Coyote Creek
(CC) branches of the SJFZ near Anza, California.

6

Figure 1-1: Location map. A map of the San Jacinto fault in southern California. Black triangles
represent stations in the Anza seismic network. The enlarged blue triangles show the locations of the
three dense FZ arrays. Circles denote events that are located by the Anza network during the operation
period of the dense arrays. The size of each circle is scaled with the event magnitude and is color-coded
with depth. The red lines denote the surface traces of the SJFZ and other major faults in Southern
California. Events located by arrows are plotted in Figures 1-3 to 1-5. The inset shows a broader view of
the study area. The cross section A-B marked by a box is used in Figures 1-10 and 1-11.

7
They identified FZ trapped waves on the three fault branches from a small
waveform subset of 25 on-fault events and 1 off-fault event. The FZ structures were
modeled as waveguides that are continuous to the base of the seismicity (~18 km) and
have depth varying properties (width of 75-100 m, shear velocity reduction of 20-30%
from that of the surrounding rock, and Q values of 40-90). They further inferred that the
BR fault zone dips to the southwest to join with the northeastward dipping CL fault
zone at depth. The CC fault was inferred as near vertical and unconnected to the others.
Igel et al. (1997) and Jahnke et al. (2002) showed with extensive 3D numerical
simulations that smoothly varying width and velocity with depth, as well as other
internal variations on a scale smaller than the FZ width, are not resolvable parameters
from trapped waves. In a continuous waveguide structure, FZ trapped waves are
sensitive only to the average properties of the waveguide. Thus the models derived by
Li & Vernon (2001) include many details that can not be resolved reliably from analysis
of trapped waves. In this work we perform systematic analysis of a much larger dataset,
modeling quantitatively trapped waves data with a considerably simpler model having
only parameters that are sensitive degrees of freedom of trapped waves.
Igel et al. (2002), Jahnke et al. (2002) and Fohrmann et al. (2004) showed that a
shallow FZ layer can trap seismic energy from events which are deeper than and well
outside the FZ. In contrast, a FZ waveguide that is continuous with depth can only trap
seismic waves of events very close to or within the FZ layer. Constraining the FZ
structures thus requires a systematic analysis of a large number of broadly distributed
events. The use of such large datasets on the San Andreas fault (SAF) at Parkfield

8
(Michael & Ben-Zion 1998; Korneev et al., 2003), the rupture zone of the 1992
Landers, California, earthquake (Peng et al. 2003), an inactive fault in central Italy
(Marra et al. 2000; Rovelli et al. 2002) and the Karadere-Duzce branch of the North
Anatolian Fault (NAF) in Turkey (Ben-Zion et al. 2003) resulted in conclusions of
shallow trapping structures.
In the following sections, we first describe briefly general aspects of the
earthquake data used in this study. Then we assign a quality of trapped wave generation
to each earthquake using a spectral ratio method and examine the spatial distribution of
the events with high quality. We further analyze the travel time moveout between body
S and trapped waves in various cross sections, and finally perform synthetic waveform
modeling of several sets of FZ waves. In contrast to the previous conclusion of
waveguides continuous to the bottom of the seismogenic zone in the SJFZ (Li &
Vernon 2001), our analysis based on a larger dataset indicates that the trapping
structures extend only over the (largely aseismic) top 3-5 km of the crust.

1-2. Analysis
1-2.1 Data collection
Three dense FZ arrays were deployed across the CC, CL and BR fault branches
of the SJFZ south of Anza to record FZ trapped waves from small earthquakes (Li &
Vernon 2001). Each array consists of 12 three-component short-period instruments that
are aligned vertical, parallel and perpendicular to the fault. Station separation is on the
order of 25 m in the center and ~50 m near the end, leading to a total array length of

9
~350 m. The arrays perpendicular to the BR and CL fault branches were deployed for
only 2 months, while the array across the CC fault recorded 6 months of data. During
the operational period, waveforms of ~1500 events within 100 km were recorded by the
arrays. In this study, we analyze a subset of ~500 events that were recorded by at least 7
instruments (including the central FZ station) within one of the arrays. The event
locations are obtained from the Anza seismic network (e.g. Berger et al. 1984). The
estimated location errors are less than 1.5 km horizontally and less than 3 km vertically
(Fletcher et al. 1987). These errors are sufficiently small for the analysis performed in
this work. Additional details on the experiment and dataset are given by Li & Vernon
(2001).
1-2.2 Spatial distribution of events generating trapped waves
The spatial distribution of earthquakes generating FZ trapped waves provides
important constraints on the overall properties of the trapping structure. In previous
studies (e.g. Li & Vernon 2001; Ben-Zion et al. 2003), the quality of the trapped waves
has generally been assigned by visual inspection. However, an automatic method is
needed in the analysis of a large number of events. Here we use a spectral ratio
technique to automatically assign quality to the FZ trapped waves generation.
The procedure employed is as follows. We first compute the amplitude spectrum
of each fault-parallel seismogram over a 2.5 s time window, starting 0.5 s before the S
arrival (Figure 1-2a). The spectra for the seismograms recorded by the stations near the
FZ (here station NE1) are divided by the mean spectra of the 4 stations (2 on each end)
that are furthest from the center of the array (Figure 1-2b). The area under the resulting

10

Figure 1-2: Spectral ratio. (a) Example fault-parallel-component seismogram recorded at array CC. The
solid vertical lines mark the portion used in the spectral analysis, while the dashed boxes denote the
stations used in the spectral ratio calculation. (b) The average velocity spectra recorded at array stations
relatively off the fault (dashed) and stations close to the fault (solid). (c) The spectral ratio obtained from
the average velocity spectra in (b). The shaded area after normalization is used as a measure of the
trapped wave quality.
spectral ratio in the range 2-12 Hz (Figure 1-2c) is normalized by dividing by the
frequency range (i.e. 10), and this value is used as a measure of the quality of the
trapped waves generated for the event.
Figures 1-3 to 1-5 show examples of seismograms and corresponding spectra
generated by events that have high quality of trapped wave generation and are located
relatively close to and off the fault. In general, stations close to the surface trace of the
faults (marked in bold) record large-amplitude low-frequency phases after the S wave

11
arrival, while these characteristics are much weaker or missing for stations further away
from the faults.

Figure 1-3: Example waveforms array CC. Examples of fault-parallel-component seismograms recorded
at array CC and generated by events located relatively close to (a) and off (b) the fault. The event ID,
depth, epicentral distance (Dist), and back azimuth (BAZ) are marked on top of each panel. The event
locations are shown in Figure 1. Red vertical bars mark the catalogue picks of the P and S wave arrivals.
The station name is marked on the left hand side of each trace, and stations close to the surface trace are
highlighted in bold. (c-d) Normalized amplitude spectra versus position of stations across array CC from
waveforms shown in (a-b).

12

Figure 1-4: Example waveforms array CL. Examples of fault-parallel-component seismograms recorded
at array CL and generated by events located relatively close to (a) and off (b) the fault. (c-d) Normalized
amplitude spectra versus position of stations across array CL from waveforms shown in (a-b). Other
notations are the same as in Figure 1-3. The less obvious phase following the S wave arrivals in (a-b) and
more broadly distributed energy in (c-d) than corresponding results in Figure 1-3 suggest a less developed
waveguide on this fault branch.

13

Figure 1-5: Example waveforms array BR. Examples of fault-parallel-component seismograms recorded
at array BR and generated by events located relatively close to (a) and off (b) the fault. (c-d) Normalized
amplitude spectra versus position of stations across array CL from waveforms shown in (a-b). Other
notations are the same as in Figure 1-3.
Figures 1-3 to 1-5(c to d) give contour maps of normalized amplitude spectra versus
position of the stations. The concentration of spectral energy between about 4-10 Hz
within ~100 m of the surface trace of these faults is associated with FZ trapped waves.
The differences of the spectral distributions for events close to and off the faults are

14
smaller than the differences between data recorded by the different arrays, suggesting
that the near-station fault structures play an important role in generating the observed
FZ trapped waves. We note that the maximum amplitude of trapped waves and
corresponding peaks of spectral energy are offset to the NE from the stations located
closest to the surface traces of each fault. We return to this issue in later sections.
Figure 1-6 shows the locations, coded with the value of spectral ratios, for all events
recorded by each of the three arrays. The top 25%, middle 50% and bottom 25% of the
spectral ratios for all of the arrays are assigned, respectively, quality A, B and C for FZ
trapped waves generation. In general, events near the faults have higher spectral ratios
than those that are further away. However, some events clearly off the faults generate
FZ trapped waves with quality A or B. In section 2.3 we show clear examples of events
outside the faults, along cross-section A-B of Figure 1-1, with quality A or B of trapped
waves generation (Figures 1-10a and 1-11a). The observation of events outside the FZ
generating considerable trapped wave energy indicates the existence of a shallow
trapping structure at these faults. In addition, the overall spectral ratios are higher for
the CC array (with a mean value of 1.68) than those for the CL and BR arrays (with
mean values of 1.53 and 1.50, respectively), suggesting that the trapping efficiencies for
the three faults are different.

15

Figure 1-6: Event Quality maps. Maps of the
events showing the calculated quality of the
trapped waves recorded at array CC (a), CL (b)
and BR (c). The shading represents the value of
the assigned quality with darker representing
higher values, and thus larger trapped waves.
The top 25 per cent, middle 50 per cent and
bottom 25 per cent of the spectral ratios for all of
the data are assigned quality A, B and C for FZ
trapped wave generation, and are marked as
stars, triangles and circles, respectively. The
fault surface traces are marked as thin lines and
the location of the array as a large red triangle.
Waveforms generated by events that are pointed
to by arrows and highlighted are used in
waveform modeling, Figs 1-13 to 1-15.

16
One theoretical characteristic of trapped waves is a predominance of motion
parallel to the low-velocity FZ layer (Ben-Zion and Aki, 1990). Our method of
assigning quality of trapped wave generation does not differentiate between trapped
waves and other FZ related site effects. To verify the validity of the method, the
obtained quality values are compared with the spectral ratios of the fault-parallel
divided by fault-perpendicular components of the near-fault records. The results (Figure
1-7) show that as the events quality increase, so do the ratios of amplification in the
fault-parallel motion relative to the fault-perpendicular component. This lends support
to our automatic procedure of assigning quality of trapped wave generation.

17

Figure 1-7: Ratio of horizontal components. The ratio between fault-parallel and fault-perpendicular
components plotted against the quality of the trapped wave generation for the events recorded at array CC
(a), CL (b) and BR (c).

18
1-2.3 Travel time moveout analysis
To further constrain the depth extent of the trapping structure, we examine the
dependence of the travel time delays, or moveout, between the FZ trapped waves and S
phases on the propagation distance inside the FZ. For a low-velocity FZ layer in a
homogenous half space, the time difference ∆t between the S phase and trapped waves
group increases with propagation distance r
s
in the FZ layer (Ben-Zion et al. 2003) as
t r
FZ HS
FZ HS
s
∆
−
=
β β
β β 2
, (1)
where β
HS
and

β
FZ
are the S wave velocity of the half space (HS) and FZ, respectively. If
the waveguide is continuous, we would expect to see a moveout between the S phase
and trapped waves for events with increasing propagation distance within the FZ.
Figure 1-8a shows the fault-parallel seismograms at FZ station CCNE1
generated by events within a 5 km epicentral distance from the station. Most of the
events are located inside or close to the CC fault (Figures 1-8b and 1-8c) and are
assigned quality A of trapped waves generation. Figure 1-8d shows the time delay
between the S arrival and the center of the trapped wave group as a function of depth.
The center of the trapped wave group is calculated as the mid point between the average
time of one oscillation after the S arrival and the end of the trapped wave group. The
latter are picked by hand on the fault-parallel seismograms as the point when the
amplitude reduces to below that of the S wave. Clearly, there is no persistent moveout
between the S and FZ trapped waves with increasing depth.

19

Figure 1-8: Moveout with depth. (a) Fault-parallel seismograms at station CCNE1 for 20 events with
trapped wave quality A. The waveforms are aligned by their S waves at 4 s and marked by a red line. The
horizontal bars above and below the seismograms bracket the trapped wave energy and the plus signs
show the estimated trapped wave group centre. The vertical dashed lines mark the picked end of the
trapped wave group. The events that produce these seismograms are a subset of the highest quality
events, located within 5 km of the station. (b) The location of all the events within 5 km of array CC.The
waveforms generated by the events marked by stars (quality A) are shown in (a). The shading and
symbols represent the quality as described for Fig. 6. (c) Hypocentres projected on a cross-section
perpendicular to the CC fault. (d) Time delay between the S arrivals and the mid point of the trapped
waves group plotted against their depth for the events in (b).

20

Figure 1-9 shows similar results for FZ stations within arrays CL and BR and events
located nearby. Again there is no clear increasing trend between the time delays and
hypocentral depths. The results of Figures 1-8a to d and 1-9 indicate similar propagation
distances inside a low-velocity FZ waveguide above the shallowest events at each array.

Figure 1-9: Moveout with depth for array CL and BR. Time delay between the S arrival and the mid
point of the trapped wave group plotted against the depths for a subset of events located directly under
array CL (a), BR (b). The shading and symbols represent the quality as described for Figure 1-6.
Figures 1-10 and 1-11 give results analogous to those of Figures 1-8 and 1-9, but
for events generated along the cross-section A-B of Figure 1-1. The locations of these
events may delineate a cross-fault, which have been suggested to exist in the SJFZ
(Mori 1993; Sharp 1967). If we assume a deep and continuous waveguide, then only
events inside or close to the damage zone will produce clear FZ trapped waves.

21

Figure 1-10: Trapped wave quality array CC. (a) Quality of trapped waves generation for a subset of
events across the FZ in the box A-B of Figure 1-1 and recorded at array CC (large triangle). Other
symbols are the same as in Figure 6. (b) Hypocenters projected from the box in Figure 1-1 along the
cross-section A-B. (c) Time delay between the S arrivals and the mid point of the trapped waves group
plotted against the distances from the fault.

22

Figure 1-11: Trapped wave quality array CL (a) Quality of trapped waves generation for a subset of
events across the FZ in the box A-B of Figure 1-1 and recorded at array CL (large triangle). Other
symbols are the same as in Figure 6. (b) Hypocenters projected from the box in Figure 1-1 along the
cross-section AB. (c) Time delay between the S arrivals and the mid point of the trapped waves group
plotted against the distances from the fault. (Same as in Figure 1-10).
Such a continuous waveguide would result in a peak in the quality of trapped waves and
time delay for events with zero or small horizontal offset from the fault. However, if a

23
shallow discontinuous waveguide is responsible for generating the trapped waves, then
the quality of the trapped waves and the time delay between S and trapped waves will
not have a strong dependence on the horizontal offset for events below the trapping
structure and at different distances from the fault. The data in Figures 1-10(c) and 1-
11(c) support the latter interpretation, and again suggest shallow waveguide structures
in our study area.
Figure 1-12 shows the time delay between the S waves and trapped wave groups
for all events recorded at each array. For every event where a time delay can be
measured, the value is approximately constant and no persistent moveout with distance
exists. The variance of the data points does illustrate how using only a selected number
of events could lead to the conclusion that a trend with distance exists.
The mean time differences between the S phases and trapped waves for all
events recorded at the arrays CC, CL and BR are 0.58, 0.52 and 0.58 seconds,
respectively. The synthetic waveform modeling in the next section gives an average S
wave velocity for the HS of 3 km/s, and FZ velocities of 1.8 km/s for the three different
branches. Using these values in equation (1), the propagation distances inside the FZ
layer of the CC, CL and BR branches are estimated to be 5.2 km, 4.7 km and 5.2 km,
respectively. Since the propagation path of the FZ trapped waves includes both vertical
and along-strike components, these values are upper bounds of the depth extent of the
waveguides at these faults.

24

Figure 1-12: Travel time moveout. Time delay between the S arrival and the mid point of the trapped
wave group plotted against distance for events recorded at array CC (a), CL (b), and BR (c). The shading
denotes the spectral ratio between FZ trapped waves and S waves. Stars and triangles denote quality A
and B respectively.

25
1-2.4 Synthetic waveform modeling of the FZ waves
The observation of a broad distribution of events generating FZ waves in
Section 2.2 suggests that the trapping structure is shallow. This interpretation is
supported in Section 2.3 by a lack of travel time moveout between the FZ trapped
waves and S arrivals with respect to increasing propagation distances inside the FZ. To
further quantify the seismic properties of the FZ waveguide, we perform synthetic
waveform modeling of observed FZ waves using the 2D analytical solution of Ben-Zion
& Aki (1990) for an SH line dislocation with a unit step source-time function in a plane-
parallel layered FZ structure between two quarter-spaces. The analysis done here
employs a simplified model configuration consisting of a single vertical FZ layer in a
surrounding half-space (e.g. Ben-Zion et al. 2003). The 2D analytical solution provides
an effective modeling tool for trapped waves in a FZ structure with width much smaller
than the length and depth dimensions and larger than the length scale of internal FZ
heterogeneities (Igel et al. 1997; Jahnke et al. 2002; Ben-Zion et al. 2003). Previously,
using the same model, good waveform fits were obtained for trapped waves observed
along the Parkfield section of the SAF (Michael & Ben-Zion 1998), the Karadere-
Duzce branch of the NAF (Ben-Zion et al., 2003), the rupture zone of the Landers
earthquake (Peng et al. 2003) and other locations (Haberland et al. 2003; Mizuno et al.
2004). While trapped waves are not sensitive to gradual internal variations in the FZ
structure, they are highly sensitive (Ben-Zion 1998) to changes in the overall average
properties of the waveguide. The employed model accounts for these sensitive
parameters.

26
Ben-Zion (1998) showed that there are strong non-orthogonal trade-offs
between sets of average 2D fault zone parameters. To quantitatively account for the
trade-offs, we model the observed FZ trapped waves recorded at multiple stations using
a genetic inversion algorithm (GIA) that employs the 2D analytical solution as a
forward kernel (Michael & Ben-Zion 1998). The free parameters in the inversion are the
S-wave velocities of the FZ and host rock, the S-wave attenuation coefficient of the FZ
material, the width and propagation distance inside the FZ layer, the source position,
and the center of the FZ layer with respect to the instrument located at the surface trace
of the fault. To reduce the number of parameters, we fix the attenuation coefficient the
host rock to be 1000 (e.g., Ben-Zion et al. 2003). Additional details on the method can
be found in Ben-Zion et al. (2003).
Prior to the inversion, we remove the instrument response and convert the fault-
parallel seismograms into displacement. We also convolve the seismograms with 1/t
1/2

to obtain the equivalent 2D line-source seismograms (e.g. Igel et al. 2002; Ben-Zion et
al. 2003). Figure 1-13(a) shows examples of synthetic waveforms fits for the CC array
along with the fitness values calculated by the GIA for different FZ parameters. The
fitness for a given set of parameters is defined as (1+C)/2, where C is the cross-
correlation coefficient between the employed sets of observed and synthetic waveforms
at multiple stations. The thin curves in Figure 1-13(b) give probability density functions
for the various model parameters, calculated by summing the fitness values of the final
2000 inversion iterations and normalizing the results to have unit sums. The synthetic
waveforms fits were generated using the best fitting parameters associated with the

27
highest fitness values during 10,000 iterations. The ranges of parameters with relatively
high fitness values provide estimates for the errors associated with the best fitting
parameters. Figures 1-14 and 1-15 show corresponding waveform modeling results for
data recorded by arrays CL and BR.

Figure 1-13: Waveform inversion CC (a) Synthetic (red lines) waveform fit of fault-parallel displacement
seismograms (dark lines) at array CC for the event located in Figure 1-6(a). (b) Fitness values (dots)
associated with different FZ parameters tested by the genetic inversion algorithm. The model parameters
associated with the highest fitness values (solid circles) were used to generate the synthetic waveforms in
(a). The curves give probability density functions for the various model parameters.

28

Figure 1-14: Waveform inversion CL (a) Synthetic (red lines) waveform fits of fault-parallel
displacement seismograms (dark lines) at array CL for the event located in Figure 1-6(b). (b) Fitness
values (dots) associated with different FZ parameters tested by the genetic inversion algorithm for array.
The model parameters associated with the highest fitness values (solid circles) were used to generate the
synthetic waveforms in (a). The curves give probability density functions for the various model
parameters.

29

Figure 1-15: Waveform inversion array BR. (a) Synthetic (red lines) waveform fits of fault-parallel
displacement seismograms (dark lines) at array BR. For the event located in Figure 1-6(c). (b) Fitness
values (dots) associated with different FZ parameters tested by the genetic inversion algorithm. The
model parameters associated with the highest fitness values (solid circles) were used to generate the
synthetic waveforms in (a). The curves give probability density functions for the various model
parameters.

Our waveform modeling suggests a waveguide propagation distance of ~5 km
for all three faults. This value is compatible with the results obtained in the previous
sections and is associated, as mentioned earlier, with a combination of along-strike and
vertical components of the propagation path. If we assume for simplicity that the

30
average along-strike and vertical components are the same, we get waveguides depth of
~3.5 km. The other average estimated FZ properties for the CC fault are FZ width of
~125 m, a reduced S-wave velocity relative to the host rock of approximately 35-45%,
and S-wave quality factor of around 30-40. For fault CL, we obtain a width of ~180 m,
a reduced S-wave velocity relative to the host rock of 40%, and S-wave quality factor of
around 20-30. The corresponding results for fault BR are a width of ~140 m, a reduced
S-wave velocity relative to the host rock of approximately 30%, and S-wave quality
factor of around less than 20. The width of the damage zone at the CC branch, with the
highest quality of trapped wave recordings, is the narrowest among the three branches
in the study area. The inversion results also indicate that the waveguides are not
symmetric about the surface traces of the faults but are offset to the northeast by 50-100
m. This can be seen directly from the waveform records of Figures 1-13 to 1-15 (and
the results of Figure 1-3 to 1-5). Waveform inversions done for several other sets of
seismograms recorded at each FZ array lead to similar overall results to those associated
with Figures 1-13 to 1-15.
The flat regions with relatively high fitness values in some of the panels of
Figures 1-13(b) to 1-15(b) reflect the strong trade-offs that exist between the parameters
governing the trapped waves (Ben-Zion, 1998; Ben-Zion et al., 2003). We can, for
example, fit the data well with larger FZ propagation distances if we adjust properly the
values of FZ width, velocity and attenuation coefficient. This may explain why Li and
Vernon (2001) were able to fit observed FZ trapped waves with a FZ layer that spans
the entire seismogenic zone (although we note that they did not give quantitative

31
measures for the fitting quality and parameter-space results as done in Figures 1-13 to
1-15). As mentioned above, waveform fits of FZ trapped can not provide a unique
image of the FZ structure because of the strong trade-offs between model parameters.
For this reason, it is essential to obtain additional constraints on the likely ranges of
parameters from analysis of the spatial distribution of events generating trapped waves,
and travel time moveout of phases, as done here and our previous related studies (Ben-
Zion et al., 2003; Peng et al., 2003).

1-3. Discussion
We analyzed waveform and travel time data from a set of ~500 events recorded
on three linear arrays across the CC, CL and BR branches of the SJFZ. We first assign a
quality of the FZ trapped waves generation based on a spectral ratio method. Many
events, including some clearly off the fault, have a high value of trapped wave quality,
indicating that the trapping structures are shallow (Igel et al. 2002; Ben-Zion et al.
2003; Fohrmann et al. 2004). This interpretation is further supported by a lack of
moveout for the time delays between the S wave and trapped waves group with
increasing propagation distance. The average time delay between the S and trapped
waves gives a total propagation distance in the FZ waveguide of less than 5 km,
suggesting a waveguides depth of about 3.5 km. The synthetic waveform modeling
provides additional information on the depth and other properties of the trapping
structure. The total propagation distance calculated by the GIA for each fault is ~5 km,
compatible with the travel time analysis and waveguide depth of 3 to 5 km.

32
The inversion results indicate that the waveguide at each fault is not centered on
the surface trace but is offset by 50-100 m to the northeast. Interestingly, detailed
geological mapping of the SJFZ at Anza show that the surface structure of the FZ has a
similar asymmetry at a scale of meters, with the gouge and host rock to the northeast of
the principal slip zone having considerably more damage than the corresponding rock
units to the southwest (Dor et al. 2006). These asymmetric patterns of damaged rocks
across the main traces of the faults may result from a preferred northwest propagation
direction of earthquake ruptures on the SJFZ (Ben-Zion & Shi 2005, and references
therein). Such a preferred propagation direction may, in turn, result from a velocity
contrast across the SJFZ (e.g., Weertman, 1980; Andrews and Ben-Zion 1997), with the
northeast side having higher seismic velocities. Local travel time tomography in the
Anza area (Scott et al., 1994) and regional imaging studies (Magistrale & Sanders 1995;
Shapiro et al., 2005) show that this is indeed the case. The analysis of trapped waves is
not sensitive to a small velocity contrast between the crustal blocks on the opposite
sides of the SJFZ, but such a contrast can have significant implications for many aspects
of earthquake dynamics (e.g., Ben-Zion, 2001, and references therein). It is therefore
important to obtain in future studies stronger constraints on the velocity contrast across
the SJFZ from detailed tomography (e.g., Thurber et al. 2004) and analysis of FZ head
waves that refract along material interfaces in the FZ structure (e.g., Ben-Zion et al.
1992; McGuire & Ben-Zion 2005).
In various recent studies, similar shallow trapping structure have been found in
the Karadere-Düzce branch of the NAF (Ben-Zion et al. 2003), the Parkfield segment of

33
the SAF (Michael & Ben-Zion 1998; Korneev et al. 2003), the rupture zone of the
Landers earthquake (Peng et al. 2003), and an inactive FZ near Norcera Umbra, central
Italy (Rovelli et al. 2002). This wide global distribution of faults with shallow trapping
waveguides indicates that it is a common element of FZ structures. The fact that such a
structure exists in inactive FZ (Rovelli et al. 2002) implies that the shallow FZ
waveguide is a long-lasting structure. Ben-Zion et al. (2003) suggested that the shallow
trapping waveguide corresponds to the top part of a flower-type structure that is
mechanically stable and is associated with relatively broad zone of deformation. This
interpretation is supported by recent observations of shallow crustal anisotropy (3-4 km)
around the rupture zone of 1999 Hector Mine, California, earthquake (Cochran et al.
2003), near the 1999 Chi-Chi, Taiwan, earthquake (Liu et al. 2004), and along the
Karadere-Düzce branches of the North Anatolia fault (Peng & Ben-Zion 2004). Since
the observed FZ trapped waves (and anisotropy) are dominated by the upper 3-4 km of
the shallow crust, the results do not provide in general information on properties of the
FZ structure at seismogenic depth, where earthquakes nucleate and the bulk of seismic
slip occurs.
Ben-Zion and Aki (1990) showed with analytical model calculations that a low
velocity FZ layer with realistic parameters can produce motion amplification of factor
30 and more in the vicinity of the fault. Cormier & Spudich (1984) and Spudich &
Olson (2001) found a large amplification for 0.6-1.0 Hz waves within ~1-2 km wide
low-velocity zone around the rupture of the 1984 Morgan Hill earthquake. Seeber et al.
(2000) observed a factor 5 amplification of acceleration in a station located in the

34
rupture zone of the 1999 Izmit earthquake on the Karadere branch of the NAF with
respect to nearby off-fault station. Rovelli et al. (2002) reported large amplification of
ground motion at stations inside an inactive fault for an M
W
5.3 subcrustal earthquake
occurred at a depth of 48 km and an epicentral distance of 10 km from Nocera Umbra,
Italy. Davis et al. (2000) reported an anomalously concentrated damage around the city
of Santa Monica during the 1994 Northridge, California earthquake, and interpreted it
as a focusing of seismic waves by subsurface faults. In this study, we found shallow (~
3 to 5 km) and ~100 m wide waveguides that can trap seismic energy generated by
earthquakes clearly outside the FZ. Because of the large volume of potential sources
that can produce motion amplifications around shallow FZ waveguides, such structures
can have important implications for seismic shaking hazard (e.g. Spudich & Olson
2001; Ben-Zion et al. 2003; Fohrmann et al. 2004).
The numerical simulations of Fohrmann et al. (2004) indicate that the volume of
sources capable of generating trapped waves at shallow waveguide structures increases
with the hypocentral depth, and that the amount of observed trapped energy depends on
the source types (e.g., orientation of the rupture plane and slip direction) and receiver
position within the radiation pattern of the sources. The SJFZ is a complex structure
with many discontinuities and subsidiary faults (e.g. Sanders & Kanamori 1984;
Sanders & Magistrale 1997). The focal mechanisms of earthquakes around the SJFZ are
highly diverse (e.g. Hauksson 2000) and the event depths are in the range of ~5-20 km.
The large diversity of focal mechanisms and range of event depths explain why many
earthquakes in our study area are capable of producing trapped waves in the FZ stations.

35
A similar situation was found to exist along the Karadere-Duzce branch of the North
Anatolian fault (Ben-Zion et al., 2003), and to a somewhat lesser extent also along the
rupture zone of the 1992 Landers earthquake (Peng et al., 2003).
The obtained FZ width of ~100 m is associated with observations of FZ trapped
waves that are in the frequency range of ~2-10 Hz, and produced by microearthquakes
typically less than magnitude 3. The width associated with FZ amplification during a
major earthquake can be much larger than ~100 m, because of the relative low-
frequency seismic waves radiated by large earthquakes (e.g. Spudich & Olson 2001).
On the other hand, non-linear amplification effects (e.g. Field et al. 1997) may reduce
the level of shaking inside the FZ during a major earthquake. The above issues should
be addressed when extrapolating the results of FZ trapped waves generated by
microearthquakes to evaluate the likely amplification that can be generated by a FZ
during a major earthquake (e.g. Michael et al. 2002).

36
Chapter 1 References
Andrews, D.J. and Y. Ben-Zion, (1997). Wrinkle-like Slip Pulse on a Fault between
Different Materials, J.Geophys. Res., 102, 553-571.
Ben-Zion, Y., (1998). Properties of seismic fault zone waves and their utility for
imaging low velocity structure, J. Geophys. Res., 103, 12,567-12,585.
Ben-Zion, Y., (2001). Dynamic Rupture in Recent Models of Earthquake Faults, J.
Mech. Phys. Solids, 49, 2209-2244.
Ben-Zion, Y. And K. Aki, (1990). Seismic radiation from an SH line source in a
laterally heterogeneous planar fault zone, Bull. Seism. Soc. Am., 80, 971-994.
Ben-Zion, Y., S. Katz and P. Leary, (1992). Joint inversion of fault zone head waves
and direct P arrivals for crustal structure near major faults, J. Geophys. Res., 97, 2,
1943-1951.
Ben-Zion, Y., Z. Peng, D. Okaya, L. Seeber, J.G. Armbruster, N. Ozer, A.J. Michael, S.
Baris and M. Aktar, (2003). A shallow fault zone structure illuminated by
trapped waves in the Karadere-Dusce branch of the north Anatolian Fault,
Western Turkey, Geophys. J. Int., 152, 699-717.
Ben-Zion, Y. and C.G. Sammis, (2003). Characterization of fault zones, Pure Appl.
Geophys., 160, 677-715.
Ben-Zion, Y. and Z. Shi, (2005). Dynamic rupture on a material interface with
spontaneous generation of plastic strain in the bulk, Earth. Plant. Sci. Lett., in
press.
Berger, J., L.N. Baker, J.N. Brune, J.B. Fletcher, T.C. Hanks and F.L.Vernon, (1984).
The Anza array: a high dynamic-range, broad-band, digitally radio-telemetered
seismic array, Bull. Seism. Soc. Am., 74, 1469-1682.
Crase, E., A. Pica, M. Noble, J. McDonald and A. Tarantola, (1990). Robust elastic
non-linear inversions: applications to real data, Geophysics, 55, 527-538.
Cormier, V. F. and P. Spudich, (1984). Amplification of ground motion and waveform
complexities in fault zones: examples from the San Andreas and the Calaveras
faults, Geophys. J. R. Astr. Soc., 79, 135-152.
Cochran, E.S., J.E. Vidale and Y-G. Li, (2003). Near-fault anisotropy following the
Hector Mine earthquake, J. Geophys. Res., 108, no. B9, 2436,
doi:10.1029/2002JB002352.

37
Davis, P.M., J.L. Rubinstein, K.H. Liu, S.S. Gao and L. Knopoff, (2000). Northridge
earthquake damage caused by geologic focusing of seismic waves, Science, 289,
1746-1750.
Dor O., T. Rockwell and Y. Ben-Zion, (2006). Geologic observations of damage
asymmetry in the structure of the San Jacinto, San Andreas and Punchbowl faults
in southern California: A possible indicator for preferred rupture propagation
direction, Pure Appl. Geophys, 163, 301-349.
Fletcher, J., L. Haar, T. Hanks, L. Baker, F. Vernon, J. Berger and J. Brune, (1987). The
digital array at Anza, California: Processing and initial interpretation of source
parameters, J. Geophys. Res., 92, 369-382.
Field, E.H., P.A. Johnson, I.A. Beresnev and Y. Zeng, (1997). Nonlinear ground-motion
amplification by sediments during the 1994 Northridge earthquake, Nature, 390,
599–602.
Fohrmann, M., G. Jahnke, H. Igel and Y. Ben-Zion, (2004). Guided waves from sources
out side faults: an indication of shallow trapping structure?, Pure appl, Geophys.,
161, 2125-2137.
Haberland, C. et al., (2003). Modeling of seismic guided waves at the Dead Sea
Transform, J. Geophys. Res., 108, no. B7, 2342, doi:10.1029/2002JB002309.
Hauksson, E., (2000). Crustal structure and seismicity distribution adjacent to the
Pacific and North America plate boundary in southern California, J. Geophys.
Res., 105, 13,875-13,903.
Igel, H., G. Jahnke and Y. Ben-Zion, (1997). Simulation of SH and P-SV wave
propagation in fault zones, Geophys. J. int., 128, 533-546.
Igel, H., G. Jahnke and Y. Ben-Zion, (2002). Numerical simulation of fault zone guided
waves: accuracy and 3-D effects, Pure appl. Geophys., 159, 2067-2083.
Jahnke, G., H. Igel and Y. Ben-Zion, (2002). Three-dimensional calculations of fault
zone guided waves in various irregular structures, Geophys. J. Int., 151, 416-426.
Korneev, V.A., R.M. Nadeau and T.V. McEvilly, (2003). Seismological Studies at
Parkfield IX: Fault-one imaging using Guided Wave Attenuation, Bull. Seism.
Soc. Am., 93, 1415-1426.
Li, Y-G. and P.C. Leary, (1990). Fault zone trapped waves, Bull. Seis. Soc. Am., 80, 5
1245-1271.
Li, Y-G., P.C. Leary, K. Aki and P. Malin, (1990). Seismic trapped modes in the Orville
and San Andreas Fault zones, Science, 249, 763-766.

38
Li, Y-G., J.E. Vidale and K. Aki, (1997). San Jacinto fault zone guided waves: a
discrimination for recently active fault strands near Anza, California, J. Geophys.
Res., 102, 11,689-11,701.
Li, Y-G. and F.L. Vernon, (2001). Characterization of the San Jacinto fault zone near
Anza, California, by fault zone trapped waves, J. Geophys. Res., 106, 30,671-
30,688.
Liu, Y., T-L. Teng and Y. Ben-Zion, (2004). Systematic analysis of shear-wave
splitting in the aftershock zone of the 1999 Chi-Chi earthquake: shallow crustal
anisotropy and lack of precursory variations, Bull. Seismol. Soc. Am., 94, 2330–
2347.
Magistrale, H. and C. Sanders, (1995), P wave image of the Peninsular Ranges
batholith, southern California, Geophys. Res. Lett., 22, 2549-2552.
Marra, F., R. Azzara, F. Bellucci, A. Caserta, G. Cultrera, G. Mele, B. Palombo, A.
Rovelli and E. Boschi, (2000). Large amplification of ground motion at rock sites
within a fault zone in Nocera Umbra (Central Italy), J. of Seis., 4, 534-554.
McGuire, J. and Y. Ben-Zion, (2005). High-resolution imaging of the Bear Valley
section of the San Andreas Fault at seismogenic depths with fault-zone head
waves and relocated seismicity, Geophys. J. Int., 163, 152-164.
Michael, A.J. and Y. Ben-Zion, (1998). Inverting Fault Zone Trapped Waves with
Genetic Algorithm, EOS, Trans. Am. geophys. Un., 79, F584.
Michael A.J., S.L. Ross and H.D. Stenner, (2002). Displaced rocks, strong motion, and
the mechanics of shallow faulting associated with the 1999 Hector Mine,
California, earthquake, Bull. Seis. Soc. Am., 92 (4): 1561-1569.
Mizuno, T., K. Nishigami, H. Ito and Y. Kuwahara, (2004). Deep structure of the
Mozumi-Sukenbu fault, central Japan, estimated from the subsurface array
observation of the fault zone trapped waves, Geophys. J. Int., 159, 622-642.
Mori, J., (1993). fault plane determinations for three small earthquakes along the San
Jacinto fault, California: search for cross faults, J. Geophys. Res., 98, 17, 711-17,
722.
Peng, Z. And Y. Ben-Zion, (2004). Systematic analysis of crustal anisotropy along the
Karadere-Düzce branch of the north Anatolian fault, Geophys. J. Int., 159, 253–
274.

39
Peng, Z., Y. Ben-Zion, L. Zhu and A.J. Michael, (2003). Inference of a shallow fault
zone layer in the rupture zone of the 1992 Landers, California earthquake fro
locations of events generating trapped waves and travel time analysis, Geophys. J.
Int., 155, 1021-1041.
Rovelli, A., A. Caserta, F. Marra and V. Ruggiero, (2002). Can seismic waves be
trapped inside an inactive fault zone? The case study of Nocera Umbra, central
Italy, Bull. Seis. Soc. Am,. 92, 2217-2232.
Sanders, C.O. and H. Kanamori, (1984). A seismo-tectonic analysis of the Anza seismic
gap, San Jacinto fault zone, southern California, J. Geophys. Res., 89, 5,873–
5,890.
Sanders C.O. and H. Magistrale, (1997). Segmentation of the northern San Jacinto fault
zone, southern California, J. Geophys. Res., 102(B12), 27453-27467.
Seeber, L., J.G. Armbruster, N. Ozer, M. Aktar, S. Baris, D. Okaya, Y. Ben-Zion and E.
Field, (2000). The 1999 Earthquake Sequence along the North Anatolia Transform
at the Juncture between the Two Main Ruptures, in Barka et al. edit. The 1999
Izmit and Duzce Earthquakes: preliminary results, Istanbul technical universit,
209-223.
Scott, J. S., T.G. Masters and F.L. Vernon, (1994). 3-D velocity structure of the San
Jacinto Fault zone near Anza, California, Geophys. J. Int., 119, 611-626.
Sharp, R.J., (1967). San Jacinto fault zone in the peninsular ranges of southern
California, Geol. Soc. Am. Bull., 78, 705-730.
Shapiro, N.M., M. Campillo, L. Stehly and M.H. Ritzwoller, (2005). High resolution
surface wave tomography from ambient seismic noise, Science, 307, 1615-1618.
Spudich, P. and K.B. Olson, (2001). Fault zone amplified waves as a possible seismic
hazard along the Calaveras fault in central California, Geophys. Res. Lett., 28,
2533-2536.
Thurber, C., S. Roecker, H. Zhang, S. Baher and W. Ellsworth, (2004). Fine-scale
structure of the San Andreas fault and location of the SAFOD target earthquakes,
Geophys. Res. Lett., 31, L12S02, doi 10.1029/2003GL019398.
Vidale, J.E., D.V. Helmberger and R.W. Clayton, (1985). Finite-difference
seismograms for SH waves, Bull. Seis. Soc. Am., 75, 1765-1782.
Wagner, G. S., (1998). Local wave propagation near the San Jacinto fault zone,
Southern California: observations from a three component seismic array, J.
Geophys. Res., 103, 7231-7246.
42
Chapter 2: Imaging the Deep Structure of the San Andreas Fault
South of Hollister With Joint Analysis of Fault-Zone Head and Direct
P Arrivals
Summary
We perform a joint inversion of arrival time data generated by direct P and fault
zone head waves in the San Andreas Fault south of Hollister, CA, to obtain a high-
resolution local velocity structure. The incorporation of head waves allows us to obtain
a sharp image of the overall velocity contrast across the fault as a function of depth,
while the use of near-fault data allows us to resolve internal variations in the fault zone
structure. The data consist of over 9800 direct P and over 2700 head wave arrival times
from 450 events at up to 54 stations of a dense temporary seismic array and the
permanent northern California seismic network in the area. One set of inversions is
performed upon the whole dataset, and 5 inversion sets are performed on various data
subsets in an effort to resolve details of the fault zone structure. The results imply a
strong contrast of P wave velocities across the fault of ~50% in the shallow section, and
lower contrasts of 10-20% below 3 km, with the southwest being the side with faster
velocities. The presence of a shallow low velocity zone around the fault, which could
corresponds to the damage structures imaged in trapped wave studies, is detected by
inversions using subsets of the data made up of only stations close to the fault. The
faster southwest side of the fault shows the development of a shallow low velocity fault
zone layer in inversions using instruments closer and closer to the fault (<5 km and <2
km). Such a feature is not present in results of inversions using only stations at greater
distances from the fault. On the slower northeast side of the fault, the presence of a low
velocity shallow layer is only detected in the inversions using the stations within 2 km
43
of the fault. We interpret this asymmetry across the fault as a possible indication of a
preferred propagation direction of earthquake ruptures in the region. Using events from
different portions of the fault, the head wave inversions also resolve small-scale features
of the fault visible in the surface geology and relocated seismicity.
2-1. Introduction
2-1.1 Fault zone structure and head waves
Understanding the structure of large faults is an important step towards the
understanding of earthquake processes on those faults. Various studies have made use
of closely spaced networks of seismometers to produce tomographic images of the
crustal seismic velocities around large continental faults (e.g., Scott et al., 1994;
Eberhart-Phillips et al., 1993; Thurber et al., 1997). However, direct P or S wave arrival
time tomography can not resolve well the short length scales that are important for
earthquake physics, and the same holds for reflection/refraction studies in the case of
near-vertical faults. This is because (1) the arrival time picks of P and S body waves are
not sensitive to low velocity zones, and (2) the travel time is a function of the entire
source-receiver path, of which a narrow fault zone (FZ) structure would be a very small
part. These methods are thus inadequate for revealing details of the localized FZ
structure at depth.
A number of studies have suggested that FZ trapped waves can be used to obtain
a high resolution image of the seismogenic structure of specific fault segments (e.g., Li
et al., 1990; 1997b). FZ trapped waves result from the constructive interference of
critically reflected phases within a sufficiently uniform zone of low velocity (damaged
44
rock) that acts as a waveguide (Ben-Zion and Aki, 1990; Igel et al., 1997). As a result,
the characteristics of observed FZ trapped waves depend strongly on the average
geometry and seismic properties of the generating low velocity layer (Ben-Zion, 1998;
Jahnke et al., 2002). In the last 15 years, trapped waves were used to image FZ layers
with width of 10’s to a few 100’s of meters in the structure of several large fault and
rupture zones. These include the San Andreas fault (SAF) near Parkfield (Li et al.,
1991, 1997a; Michael & Ben-Zion, 1998; Korneev et al., 2003), the San Jacinto fault
near Anza (Li et al., 1997b; Lewis et al., 2005), an inactive fault in central Italy (Marra
et al., 2000; Rovelli et al., 2002), the Landers rupture zone (Li et al., 1994; Peng et al.,
2003) and the Karadere-Duzce branch of the North Anatolian fault (Ben-Zion et al.,
2003). The modeling of trapped waves is, however, associated with large trade-offs
between key FZ parameters (Ben-Zion, 1998). For example, considerably different
propagation distances within the FZ can fit the same waveforms equally well if
appropriate changes in the FZ width, velocity and attenuation coefficients are made. As
a result, robust conclusions based on waveform modeling of FZ trapped waves, without
constraining evidence on some of the parameters, are very hard to draw. Furthermore, a
growing body of evidence based on comprehensive multi-signal analysis of large
datasets suggests that the damage zones producing trapped waves are generally limited
to the upper ~3-4 km (Peng et al., 2003; Ben-Zion et al., 2003; Korneev et al., 2003;
Lewis et al., 2005; Rovelli et al., 2002). It thus appears that trapped wave studies cannot
provide detailed information on the fault structure at seimogenic depths, except that it is
not characterized by a uniform damage zone extending continuously from the surface.
45
Plate boundaries and other major faults may juxtapose rock bodies with different
elastic properties. Two regions where juxtaposition of different crustal blocks has been
imaged using seismic tomography are the Cienega Valley area (Thurber et al., 1997)
and Parkfield section (Eberhart-Phillips et al., 1993) of the SAF, with the contrast of
seismic velocities extending throughout the seismogenic zone. A sharp material contrast
separating different crustal blocks will generate (Ben-Zion, 1989, 1990; Ben-Zion and
Aki, 1990) fault zone head waves (FZHW) that spend a large portion of their
propagation paths refracting along the material interface that defines the fault (Figure 2-
1). The head waves propagate along the material interface with the velocity and motion
polarity of the faster block, and are radiated from the fault to the slower velocity block
where they are characterized (Figure 2-2) by an emergent waveform with opposite
motion polarity to that of the direct body waves. The FZHW are the first arriving
seismic energy at locations on the slower side of the fault with normal distance to the
fault (Ben-Zion, 1989) less than a critical distance xc















• =
−
1
2 1
cos tan
α
α
r x
c
, (1)

where r is the along-fault propagation distance and α
2
and α
1
are the average P wave
velocities of the slower and faster media, respectively.
46

Figure 2-1: Fault model. (a) Cross section perpendicular to the fault showing the structure of the
layered quarter spaces fault zone velocity model. The P wave velocities and layer thickness are denoted
by α and h, respectively, with subscripts f
j
and s
j
representing the fast and slow sides of the fault for each
layer. The source (star) is on the fault at depth and the receivers (triangles) are located in the x-y plane on
the surface or at some depth z, less than the thickness h
1
of the shallowest layer. (b) A view in the
direction perpendicular to the fault, showing the ray path (arrows) of the head wave along the fault plane,
which is the coordinate used for the head wave travel time calculation. The vertical components of rays
and P wave velocities in layer j on the faster side of the fault are denoted by ζ
fj
and α
fj
, respectively.
From Ben-Zion et al. (1992).

47
-1
-0.5
0
0.5
1
1.76 1.78 1.8 1.82 1.84 1.86 1.88 1.9 1.92 1.94
-1
-0.5
0
0.5
1
Normalized amplitude
Time from origin (sec)
a)
b)
head
wave
direct P
wave
direct P
wave

Figure 2-2: Synthetic head waves. Synthetic seismograms based on the 2D analytical solution of Ben-
Zion (1989) for line dislocation with a unit step function in time at the interface between two quarter-
spaces with different materials properties. The contrast of P wave velocities between the two quarter-
spaces is set at 25%. (a) Normalized velocity (solid) and displacement (dashed) seismograms for an
instrument 10 km along and 2 km normal to the fault on the quarter-space with slower velocity. The
motion at this position consists of a first arriving emergent head wave followed by a sharp direct P
arrival. (b) Same as (a) but for a receiver 2 km normal to the fault on the quarter-space with faster
velocity. The motion at this position consists only of a direct P arrival.

The existence of material interfaces in fault zone structures can have important
implications for many aspects of earthquake physics, including effective constitutive
laws, suppression of branching, frictional heat, short rise-time of earthquake slip,
earthquake interaction, and seismic shaking hazard (Ben-Zion and Andrews, 1998; Ben-
Zion, 2001). Since the FZHW are generated by and travel along material interfaces,
they provide the best diagnostic and imaging tool of material contrasts at seismogenic
depths of faults. Ben-Zion and Malin (1991) and Ben-Zion et al. (1992) identified
FZHW along the Parkfield section of the SAF and inverted arrival time data of the FZ
head and direct P waves for the local velocity structure. Ben-Zion et al. (1992) also
48
performed numerical experiments showing that the use of head waves can significantly
improve the resolution of the velocity structure across the fault. More recently, McGuire
and Ben-Zion (2005) analyzed arrival times and waveforms characteristics of FZHWs,
direct P waves and additional FZ phases in the early portion of seismograms, using data
collected from a dense temporary array (Thurber et al., 1997) around the SAF south of
Hollister. The results demonstrated clearly the existence of pronounced material
interfaces, sufficiently planer and sharp to generate head waves, in the seismogenic
structure of the SAF in that area.
In the present work we derive additional information on the fine velocity
structure of the SAF south of Hollister, based on joint inversions of FZHW and direct P
arrival times. The inversions employ a simple model geared to focus on the arrival
times of seismic phases generated in a layered structure with a sharp planar material
interface. Using this model with different station and earthquake combinations allows
some exploration of how the structure varies along strike and with distance from the
fault. As shown in the following sections, the employed simple model and a related
complementary framework can explain well the major features (arrival time, amplitude
and polarity) of the seismic phases in the early portions of the observed seismograms
While a fully 3D model would be required to explain all the observed features, our
focus is on imaging the velocity structure in the immediate vicinity of the SAF. This
narrow structural component can be imaged well by analysis of the observed near-fault
P phases and FZ head waves within the simplified framework of our model.
49
121˚ 25' 121˚ 20' 121˚ 15' 121˚ 10'
36˚ 35'
36˚ 40'
36˚ 45'
121˚ 25' 121˚ 20' 121˚ 15' 121˚ 10'
36˚ 35'
36˚ 40'
36˚ 45'
0 5
km
PIT
STS
TTN
SET
YEL
NAK
CAR
SUL
TMP
TRP
LIB
COW
STK
CRH
SLD IND
PAR
FOX
DOZ
PEA
VLT
TOP
RES
PBB
NRA
RNZ
CLF
GIB
AGY
VIN
UUP
SAS SA2
AKC
PIP
TBS
ATK
FWD
H20
PON
FOR
WIN
MUD
BAT
KEY
SUM
HVL
BRG
BCG
BEH
BHR
BJO
BLR
BSL
BSC
Depth Range
0 ~ 3 km
3 ~ 6 km
6 ~ 9 km
9 ~ 12 km
> 12 km
M range
1
2
3
Los Angeles
San Franciso
Pacific Ocean
California

Figure 2-3: Location map. A map of the study area with gray shading used to represent topography and
red lines the surface traces of significant faults. The black triangles represent stations in the temporary
array while the red triangles are permanent station from the NCSN, all of which labeled with their three
letter code names. The circles are positions of 12 years of events relocated by McGuire and Ben-Zion
(2005). The size of each circle scales with the event magnitude and the colors represent their depth. The
inset on the bottom left shows a large scale map, with major faults marked in red, major cities, and the
location of the main map inside a small box.
2-1.2 Data and geological setting
The Cienega Valley to Bear Valley region of the SAF south of Hollister lies at
the northern end of the creeping section (Figure 2-3) and is characterized at the surface
by granites plus Pliocene and Cretaceous sediments (Thurber et al., 1997). The SAF is
50
generally the interface between the higher velocity igneous rocks to the SW and lower
velocity sedimentary rocks to the NE. Thurber et al. (1997) deployed a 48-station
PASSCAL array around the SAF (Figure 2-3), with the aim of performing local P wave
arrival time tomography. The instruments were a mix of standard short period
seismometers, primarily Mark products L22s, L4s and Teledyne S-13s, with 33 located
on the NE side between ~200 m and as much as 10 km away from the surface trace of
the fault. In the 6 months that the temporary array was deployed, some 1200 local
events were recorded and have locations in the Northern California Seismic network
(NCSN) catalogue. With instruments positioned closely around the fault on both the fast
SW and slow NE sides, and with a clear velocity contrast across the fault, the
deployment is well suited for FZHW analysis.
In the present work we analyze data recorded by Thurber et al. (1997) and
additional 7 permanent NCSN stations in the study area. We use accurate earthquake
locations from McGuire and Ben-Zion (2005) based on the double difference technique
(Waldhauser and Ellsworth, 2000) when available, and otherwise the standard NCSN
catalog locations from the Northern California Earthquake Data Center. The relocated
events generally collapse onto an approximately 100-200 m wide near-vertical zone that
is offset to the SW from the surface trace of the fault. There is a long standing debate on
whether this offset is produced by the velocity contrast across the fault (e.g., Bakun et
al., 1980; Ellsworth, 1975), or whether it reflects a fault dip (e.g., Aki and Lee 1976;
Thurber et al., 1997). To reconcile the narrow near-vertical zone containing the relative
locations of the relocated seismicity and the surface trace of the fault, we follow
McGuire and Ben-Zion (2005) and shift the double-difference locations 0.7 km to the
51
NE. This distance is chosen so that the seismicity is 100-200 m to the SW of station
PIT, a location consistent with the surface trace of the fault. This shift is also necessary
to produce first motion polarities and frequency contents consistent with right lateral
focal mechanisms on a near-vertical fault at stations on the southwest side near the
fault.
The offset between the surface trace and the inferred event locations is not
consistently 0.7 km throughout the area but rather increases to the NW. Nevertheless,
the seismicity in our analysis is approximated as a plane with the assumed fault model
consistent with the average locations, after they are shifted 0.7 km to the NE. We note
that the uncertainty in the absolute locations normal to the fault has a small effect on our
results, which are based on arrival time data that are sensitive primarily to the distance
traveled along the fault. As discussed in section 2.2, sets of inversions in which the
assumed fault position is varied within a zone of 1000 m show that the resulting change
to the obtained velocity contrast is only about 1% (which is less than the variations
between inversion runs in the same set).
To the southeast of the network the relocated hypocenters do not collapse onto a
single plane, but to a pair of parallel strands that are active in different depth ranges
(Figure 2-4). The strands are also shown in the geological map of the area (Dibblee,
1974), with the surface trace of the SAF splitting into two branches separated by a sliver
52
-6
-5
-4
-3
-2
-1
0
1
2
Distnace along strike (km)
-2 -1 0 1 2
Distance from PIT (km)
PIT PI T
STS STS
TTN TTN
SET SET
-10
-9
-8
-7
-6
-5
-4
-3
-2
-1
0
Depth (km)
-2 -1 0 1 2
Distance from PIT
a) b)

Figure 2-4: Seismicity in a around a surface trace complexity. Relocated earthquakes in the region of the
structural complexity. (a) A map view of the earthquakes, black dots, and the southeasterly stations of the
network, triangles, rotated to have the fault running straight from top to bottom and on a local coordinate
system centered on station PIT. The box denotes an area, shown in cross section in (b), where the fault
splits into two branches. (b) The events from the box in (a) plotted onto a vertical plane. The two fault
branches can be seen to be active in different depth ranges.
of granitic rocks. The double-strand structure also manifests itself in the waveforms of
events with origins to the south of the temporary network, with at least one extra phase
between the head and direct wave arrivals. When performing waveform modeling,
McGuire and Ben-Zion (2005) found that the energy between the head and direct waves
(Figure 2-5b) requires an intermediate velocity layer, consistent with a sliver of granite,
between two quarter spaces with different elastic properties. This structure impinges our
study in two main ways. First, picking the arrival times of the phases, particularly the
direct wave, for events to the SE of the network is less accurate. Second, the model that
53
we use for the arrival time inversions (Figure 2-1) consists of two layered quarter-
spaces separated by a single interface. Thus data generated by events on both the
parallel strands are mapped by our ray-tracing model to a single effective interface,
rather than being fully accounted for. As discussed above, small changes of locations
normal to the fault have little impact on our results.
Taking into account the structural complexity to the southeast, it would be
preferable to concentrate on performing inversions for the part of the SAF within and to
the NW of the network. However in the period of operation of the dense temporary
array, only about 120 (10%) of the events occur on the simple (single-trace) portion of
the fault. Of the 300 events for which arrival times have been picked on instruments of
the temporary array, more than 50% come from within, or to the SE of, the section of
the fault with the double-strand structure. Results based on such data would therefore be
skewed towards this area where our model simplifies an important structural
component. To work with a more balanced dataset that averages small-scale
complexities not accounted for by the model, 150 events that occurred before the
operation of the array to the north of the double-strand area and were recorded by the
seven NCSN stations are also used. This gives an approximately even distribution of the
number of arrival times from the regions of the fault within and to the NW and SE of
the dense temporary network. As demonstrated in the subsequent sections, the assumed
model of two juxtaposed vertically-layered quarter-spaces can be used to invert the
observed P and head arrival times at the scale of the observations to a stable sharp
image of the velocity contrast across the fault at depth. To explore how additional
complexities affect the results, we perform a suite of inversions using various data
54
subsets, along with some synthetic waveform fits in related model that includes fault
zone layers.
2-2. Method
2-2.1 Travel times
The FZ head and direct P wave arrival times were picked by hand with each
waveform being examined at least twice to increase consistency. For events on the slow
side of the fault within the critical distance (equation 1), the waveforms are
characterized by an emergent first arriving head wave, opposite in polarity to that
expected for a right-lateral strike-slip event located on the SAF, followed by a larger
amplitude direct P wave with the expected polarity. Figure 2-5a shows a set of velocity
seismograms from events at various distances to the north of station SUM, located 1-2
km from the fault on the slow side (NE) of the fault. The waveforms are aligned on their
P wave arrival time and display an increase in the time between the head and direct P
wave arrival with propagation distance along the fault. For an interface between two
different quarter spaces, this differential arrival time ( ∆t) grows with the distance
traveled along the fault (r) as






 ∆








− ∆
2
1 2
~
1 1
~
α
α
α α
r r t
, (2)

with α and ∆α denoting the average and differential P waves velocities, respectively
(Ben-Zion and Malin, 1991). Thus, if an average P wave velocity can be accurately
55
estimated, it can be used along with the gradient of the differential arrival time versus
propagation distance curve to obtain a first order approximation for the strength of the
velocity contrast (∆α).
SUM P Wave Arrivals From the NW and SE
HW
HW
0.6 0.4 0.2 0
8
9
10
11
Time (s)
Distance Along Fault to SE (km)
0.6 0.4 0.2 0
5
6
7
8
9
Distance Along Fault to NW (km)
Time (s)
(a)
(b)

Figure 2-5: Vertical component velocity seismograms recorded at station SUM and aligned on their P
waves arrival plotted at the distance of that event along strike from the station. (a) Seismograms of events
to the northwest of the station with the head wave arrivals delineated by a black line. A clear moveout
between the direct P and head waves as the distance traveled increases can be seen. (b) Seismograms of
events to the southwest of the station showing additional complexities than those in (a). In addition to the
head wave marked by the first black line, there is coherent energy between the head and direct P arrival,
outlined by black lines at distances greater than 9.5 km. From McGuire and Ben-Zion (2005).
56
Figure 2-6 shows the differences between the direct P and head wave arrival
versus hypocenter distances for 32 of the stations on the slow side of the fault.
121°20' 121°16' 121°12'
36°40'
36°44'
36°48'
0 5
km
SU L
PEA
HVL
COW
TOP
SUM
PIT
TTN
STK
CRH
SLD
FO X
DOZ
VLT
BOB
RES
PBB
NRA
RNZ
CLF
GIB
AGY
VIN
SAS
TBS
ATK
H20
PON
FOR
WIN
MJD
KEY
=20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20% =20%
a)
b)
Propagation distance along fault (km)
Differential arrival time (sec)
0
0.5
1
1.5
SUL
14%13%
PEA
8% 3%
HVL
26% 18%
COW
10% 11%
TOP
3% 2%
SUM
17% 12%
PIT
35%14%
TTN
18% 12%
0
0.5
1
1.5
STK
7% 11%
CRH
17% 12%
SLD
16% 18%
FOX
6% 7%
DOZ
10% 3%
VLT
1% 4%
BOB
0% 0%
RES
0% 3%
0
0.5
1
1.5
PBB
0% 0%
NRA
3% 0%
RNZ
4% 1%
CLF
10% 3%
GIB
7% 12%
AGY
6% 9%
VIN
14% 12%
SAS
0% 0%
50 0 50
0
0.5
1
1.5
TBS
0% 0%
50 0 50
ATK
0% 2%
50 0 50
H20
15% 6%
50 0 50
PON
4% 2%
50 0 50
FOR
3% 4%
50 0 50
WIN
1% 2%
50 0 50
MJD
0% 2%
50 0 50
KEY
24%14%

Figure 2-6: Differential arrival time analysis. (a) Differential arrival times between the direct and head P
waves as a function of along-strike propagation distance for 32 stations on the slow side of the fault. The
red circles give results from events located to the northwest and blue circles from events to the southeast
of the station. The black lines are least-squares fits to the points, giving (for those with enough data)
estimates of the average velocity contrasts to the NW and SE. The obtained average velocity contrasts in
percent are marked at the top of each panel. (b) A map view with the stations as triangles and the major
faults as red lines. The velocity contrasts estimated in (a) at each station position, plotted as an arrow
whose length corresponds to the size of the contrast estimated from events in the indicated direction
(using same colors as in (a)).
57
The moveout gradients were estimated with a least-squares fit. The maximum
differences of differential arrival times for events at the same station and similar
distances range from ~-0.5 to +0.5 seconds. Stations closest to the fault, such as KEY
and PIT have the highest residual to the least-squares fits. It is typically quite difficult to
accurately pick arrival times at these stations because the signal/noise ratio is lower and
the waveforms are more emergent and/or complicated. This can perhaps be attributed to
the damage zone around the fault and near-nodal position of the stations in the P wave
radiation pattern. For some stations, usually those further from the fault, there are
insufficient head wave first arrival picks to constrain the moveout velocity. Taking
these available gradient estimates and an average velocity of 5 km/s (McGuire and Ben-
Zion 2005), the average velocity contrast for each station can be estimated using eq. (2).
The velocity contrast estimates (labeled in Figure 2-6a) show that stations generally
measure a higher contrast from events to the north of a given station (shown in red), and
that stations closer to the fault are affected by higher contrast values than those further
away (Figure 2-6b). The largest average contrast is about 35% at station PIT which is
within 200 m of the fault. The largest value might reflect the contrast between the low
velocity damage zone associated with the fault, on which this station likely resides, and
the faster side of the fault. More typical values are in the 10-20% range, which perhaps
better represents the true contrast in velocity between the host rocks on the opposite
sides of the fault.
58
2-2.2 Joint inversion of FZHW and direct P arrivals
To obtain more accurate estimates of the velocity contrast across the SAF,
accounting for depth-variations, we perform a joint inversion of the FZHW and P
arrival time data. The model consists of two vertically-layered quarter-spaces that
represent the two crustal blocks on the opposite sides of the fault (Figure 2-1). For
simplicity, the number of layers in each block is kept constant leaving the thicknesses
and velocities of the layers as free model parameters. We employ the ray tracing
procedure and the adjusting random search inversion algorithm of Ben-Zion et al.
(1992). The adjusting random search is similar to simulated annealing, but uses a fixed
(rather than progressively decreasing) maximum perturbation of parameters. In an
earlier phase of the work we found that the adjusting random search method, with a
maximum perturbation size tailored to properties of the data by trial and error, produces
a smaller misfit than inversions with the same number of iterations using progressively
decreasing maximum perturbation similar to simulated annealing. Additional details on
the inversion algorithm can be found in Ben-Zion et al. (1992).
The travel times for the direct and head waves to the various stations are
calculated for the model shown in Figure 2-1 with the following expressions (Ben-Zion
et al., 1992). For the direct P wave, the travel time td is

∑
− + ⋅ =
− 2 2
p p t
j j d
α ζ ξ , (3a)

where ξ is the horizontal component of the straight source-receiver line, p is the
component of slowness along ξ, and ζ
j
and α
j
are the vertical component of the ray path
59
and the P wave velocity in medium j, respectively. For the head wave, the travel time t
h

is

[ ]
2
1
2
1
2 2 − − −
− + − + ⋅ =
∑ f s fj fj h
x p p y t α α α ζ , (3b)

where y is the component of the source-receiver separation along the fault, p is the y
component of slowness, ζ
fj
and α
fj
are, respectively, the vertical component of the ray
and P wave velocity in medium j on the faster side of the fault, and α
f1
, α
s1
denote the P
wave velocity of the top layers on the faster and slower sides of the fault, respectively.
The calculated travel times are used to back-calculate model origin times of the recoded
earthquakes via subtraction from the measured arrival times. The difference between the
modeled origin times of each earthquake, calculated with each of the available direct
and head wave arrival times, gives a measure of how well the model approximates the
true structure. The model origin times based on all the available data form a model error
function Merr given by:

[]
∑∑∑
=
+
=
−
=
− =
eq
i
i h i d
j
j
k
k i O j i O abs F Merr
1
) ( ) (
2
1
1
) , ( ) , ( )( , (4)

where F is a vector containing the model parameters and eq is the total number of
earthquakes. The indexes d(i) and h(i) denote, respectively, the number of direct and
head wave time picks made for earthquake i, O(i, j) is the origin time of earthquake i
calculated using the station j observed direct wave (j ≤ d(i)) or head waves (j > d(i))
60
arrival times and corresponding model calculated travel times. Since this model error
would increase with an increasing number of data points, it is normalized by the number
of origin times used in its calculation to create NMerr, the mean mismatch per a single
used pair of arrival time data for the same event. The search of best fits to both head and
direct wave arrivals involves minimization of the normalized error NMerr, using a self-
adjusting random search where the layer thicknesses and velocities are perturbed at
every inversion iteration (see additional details in Ben-Zion et al., 1992).
To maximize the efficiency of the inversion and to tailor the values of maximum
perturbation size to properties of parameters in the inversion, we perform tests using a
dataset of synthetic arrival times, calculated using the real event and station locations.
The inversion parameters that can be tuned to lead most efficiently to the minimum
model misfit are the number of iterations, the size of the maximum random
perturbations to the model parameters (layer velocities and thicknesses), and the initial
values of the parameters. While the latter should have little impact in the limit of a large
number of iterations, the effects of the other choices should be explored.
Figure 2-7 shows how the size of the maximum employed random perturbation
changes the size of the final normalized error for different numbers of iterations. The
calculations are done using the two different layered quarter-spaces obtained by Ben-
Zion et al. (1992) for the Parkfield section of the SAF. When the absolute value of the
maximum perturbation is above 0.2, the final error is increasing. The lowest final error
is for a maximum perturbation between 0.1 and 0.15. Increasing the number of
iterations from 1000 to 2000 improves the fit and reduces the error. There is some
61
improvement between 2000 and 3000 iterations, but no further reduction is found in the
region with the best maximum perturbation beyond 3000 iterations.
0.05 0.1 0.15 0.2 0.25 0.3
1000
1500
2000
2500
3000
3500
4000
4500
5000
0.026
0.028
0.03
0.032
0.034
0.036
0.038
0.04
0.042
Maximum perturbation size
Number of iterations
Normalized error

Figure 2-7: Optimising the inversion. A contour plot summarizing the ability of the inversion to
minimize the error in synthetic datasets. The absolute value of the maximum perturbation size (x-axis)
and number of iterations (y-axis) are adjustable parameters of the inversion. The local minimum of the
normalized error for different sets of values is plotted as a red line.
Assuming that the synthetic calculations are representative of the features of the data,
values of 3000 iterations and maximum perturbation of 0.1 appear to offer the best
chance of rapidly finding the global minimum. The assumed velocity model used in the
synthetic calculation and the best fitting model obtained by an inversion using the above
parameters and starting with the standard 5-layers 1-D velocity model used to locate
events in that area is shown in Figure 2-8. Further synthetic tests and explanation of the
improvement in resolution from the introduction of head waves can be found in Ben-
Zion et al. (1992).
62
3.5 4 4.5 5 5.5 6
12
10
8
6
4
2
0
Velocity (km/s)
Depth (km)

Figure 2-8: Inversion of synthetic data generated using the event and station locations of temporary
deployment and the velocity structure obtained by Ben-Zion et al. (1992) for the Parkfield section of the
SAF. The dashed lines give the velocity profiles used in the generation of the synthetic arrival times for
the fast (blue) and slow (red) sides of the fault. The solid lines are the recovered velocity structure using
the inversion method outlined in the text.
As discussed earlier, the relocated seismicity in the study area falls onto a
narrow near-vertical zone to the southwest of the surface trace of the SAF. However,
these locations are not consistent with the first motion polarities at instruments within
100 m of the fault. The first motion polarities are fit properly by assuming that the
events are on a vertical plane under the surface trace of the fault. As mentioned in the
introduction, to resolve the discrepancy between the catalog locations and first motion
polarities, we moved the vertical plane defined by the seismicity 700 m northeast such
that it intersects with the known surface trace location relative to station PIT. To test the
possible effects of this shift on our results, a series of inversions were carried out on the
observed data assuming event locations at various distances from the surface trace. The
63
set up for the inversions is the same as that described in section 3 as set (B3), and it
consists of 24 stations located at over 2 km from the surface trace, selected to have
approximately equal numbers of phase arrivals on the slow and fast sides of the fault.
The event locations are shifted within a zone that is 1000 m wide and 10 inversion runs
are performed for each assumed configuration. Figure 2-9 shows how the depth-average
velocity contrast for each set of inversions changes as a function of the distance of the
assumed event locations from the surface trace of the SAF. Within the range over which
the events are shifted, the average velocity contrasts changes by only 1% out of average
contrasts of about 26%. Hence the 700 m shift required for consistency of motion
polarities has very small effect on the inversion results.
-800 -600 -400 -200 0 200
23
23.5
24
24.5
25
25.5
26
26.5
27
Normal distance of events from fault (m)
Average velocity contrast (%)

Figure 2-9: Depth-average velocity contrasts (circles) obtained by sets of inversion assuming the fault at
depth is located at different normal distances from the surface trace of the San Andreas Fault. Negative
and positive horizontal coordinates are distances perpendicular to the surface trace (black dashed line) in
the southwest and northeast directions, respectively. The events are assumed in this study to be below the
surface trace, some 700 m from their catalog locations.
64
2-3. Results
We perform inversions using data generated by 450 events picked at up to 53
stations, resulting in over 9800 direct and 2700 head P wave arrival time picks. These
are utilized in 6 sets of inversions that are performed on the whole dataset (A), using
stations which are within 5 km (B1), within 2 km (B2) and further than 2 km (B3) from
the fault, and finally using events to the northeast of the complexity (C1) or events
within and to the southeast of the complexity (C2). Initially, a 5 layer model was used
but the inversions always reduced the thickness of one of the layers to near 0. Thus all
results presented are for 4 layers. The starting model in all cases (dashed black line in
Figure 2-10a) was the Northern California velocity structure given out with the hypoDD
program by the USGS.
Figure 2-10a shows the results from 10 inversion runs using all the data (set A),
with each inversion run performed using a different sequence of random numbers for
the self-adjusting random search. The solid colored lines represent the model with the
lowest misfit, while the dashed colored lines and grey shaded areas represent the mean
and standard deviation associated with the10 inversion runs. A general feature of the
models is a subsurface layer (depth < 1 km) on the fast side of the fault. Although this
layer is very thin, it is a common feature of the inversion results (Figure 2-10b).
However, the ranges of depths and velocities covered by the uncertainties in the top two
layers of Figure 2-10 imply that the thin subsurface layer is not well constrained by the
inversions.
65
0.5 1 1.5 2 2.5 3 3.5 4 4.5
10
9
8
7
6
5
4
3
2
1
0
Velocity Difference (km/s)
Depth (km)
(b)
(a)
0
Number of events
0 40 80 120
10
9
8
7
6
5
4
3
2
1
Depth (km)
2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
Velocity (km/s)

Figure 2-10: Inversion results. (a) Inversion results of P wave velocities versus depth on the SAF using
all the stations and data. (a) The solid red and blue lines are the best-fitting velocity profiles for the slow
and fast sides of the fault, respectively, out of 10 inversion runs performed on the data. The dashed red
and blue lines are the average depth and velocity of each layer from all 10 inversions, and the grey shaded
area around the mean represent the standard deviation of the depth and velocity of that layer. The dashed
black line is the initial velocity and depth model used in the inversion. The dashed gray line and
horizontal axis on top give the number of events as a function of depth. (b) The differences between the
velocities on the fast and slow sides of the fault, plotted at the depth where that contrast occur for each of
the 10 inversion runs. The mean and standard deviation of the values at that depth are plotted with the red
symbols.
66
Between 1 and 3 km the velocity contrast between the two sides of the fault increases to
its highest values, with the velocity of the slow side only 50% of the fast side. At
greater depths the velocity difference across the fault remains fairly constant at
approximately 1 km/s or between 10-20%. The average values of the depths and
velocities are consistent overall with the previous tomographic study of Thurber et al.
(1997), as well as with synthetic waveform fits obtained by forward modeling in the end
of the section.
To constrain further the velocity structure in the immediate vicinity of the fault,
additional inversions are performed (Figure 2-11) using subsets of instruments located
closer and closer to the fault (sets B1-B3). More prominent shallow low velocity layers
emerge for both sides of the fault, and become significant features outside the range of
the uncertainties when inversions are performed using only near fault instruments
(Figure 2-11d). In the inversion using all the stations, those further from the fault are
not affected by a low velocity fault zone layer above the seismicity of the type imaged
by recent studies with FZ trapped waves (Ben-Zion et al., 2003; Peng et al., 2003;
Lewis et al., 2005). In contrast, near-fault stations would be affected by such a FZ layer
and will produce a modification in the inversion results.
When subsets of stations within 5 km and 2 km of the fault (Figure 2-11c and 2-
11d) are used, a larger proportion of the ray paths are affected by the internal fault zone
structure, and the results more strongly reflect the structure in the immediate vicinity of
the fault. The results from those inversions imply a shallow low velocity layer on the
faster side of the fault. This feature becomes increasingly more evident as the number of
stations closer to the fault increases (Figure 2-11a to 2-11d), and it appears to be well-
67
2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
Velocity (km/s)
Depth (km)
1.5
10
9
8
7
6
5
4
3
2
1
0
(e)
0
0.2
0.4
0.6
0.8
1
1.2
(B3) > 2 km
(A) all
(B1) < 5 km
(B2) < 2 km
Number of stations used in inversion
Normalized velocity difference x depth of
shallowest layer
2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
Velocity (km/s)
1.5
Depth (km)
10
9
8
7
6
5
4
3
2
1
0
Stations
>2km (B3)
Stations
all (A)
Stations
<5km (B1)
Stations
<2km (B2)
(a) (b)
(d)
(c)

Figure 2-11: Variations in contrast with distance from the fault. (a) Results from inversions on four
different sets of data associated with increasing numbers of near-fault stations. (a) Stations further than 2
km from the fault, (b) all stations, (c) only stations within 5 km and (d) only stations within 2 km. The
solid and dashed lines represent the minimum misfit and mean values from 10 sets of inversions, while
the grey area gives the standard deviation of depth and velocity. (e) The value of the velocity difference
between the top layer and the one below it, normalized to the maximum velocity difference, and
multiplied by the depth of the top layer. This quantity represents the characteristics of the top most layer
of the velocity structure and it changes in a systematic manner in the inversion sets (a-d). The triangles
(red) and circles (blue) are for the slow and fast sides, respectively. The dashed and solid symbols
represent the best fitting and mean model results, respectively.
68
constrained outside the variations of model results obtained by different inversion runs.
For the slower side of the fault, the results show a corresponding shallow low velocity
zone only with inversions using just the stations within 2 km of the fault (Figure 2- 11d
and 2-11e). The results may be affected in part by the different takeoff angles and ray
paths from the earthquakes to the different stations. However, synthetic model
calculations with two layered quarter-spaces and the used earthquake and station
locations indicate that the trends seen in Figure 2-11 are unlikely to be produced as an
artifact of the employed subsets of stations.
The final two sets of inversions involve events to the northwest (C1) and to the
southeast (C2) of the array. While the results based on events to the northwest (Figure
2-12) are very similar to the inversion results of all events (Figure 2-10), the results
based on the subset to the southeast (Figure 2-13) have some distinct features. The
velocity of the faster side of the fault below 4-5 km is between 0.5-1 km/s faster, while
the shallowest layer extends to a greater depth (1-1.5 km) than in the other inversions.
In the depth range 6 to 8 km, the velocities on both the fast and slow sides of the fault
show a reduction with depth that is significant enough to be visible within the variations
of the different inversion runs. These features probably reflect properties of the sliver of
high velocity granite in the southeast area discussed by McGuire and Ben-Zion (2005).
The variations of results from the inversion runs using the dataset C2 are larger (Figure
2-13b) than in the other cases, showing as expected that the data from the southeast are
harder to fit with our model.
69
2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
10
9
8
7
6
5
4
3
2
1
0
Velocity (km/s)
Depth (km)
(a)
Velocity Difference (km/s)
Depth (km)
0 0.5 1 1.5 2 2.5 3 3.5 4
10
9
8
7
6
5
4
3
2
1
0
(b)

Figure 2-12: Inversion using events from the northwest. (a) The 4-layer velocity models for the opposite
sides of the fault produced by inversion for a subset of events in the central and northwestern part of the
study area. The solid lines are the results with the lowest misfit, the dashed lines are the average of the
results from 10 inversion runs using the same setup, and the gray area is the standard deviation. (b) The
velocity difference between the layers in all of the inversion runs used to obtain the average result in (a)
as a function of depth. The mean and standard deviation of the values at a given depth are plotted with the
red symbols.
70
(b)
2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5
10
9
8
7
6
5
4
3
2
1
0
Velocity (km/s)
Depth (km)
(a)
Velocity Difference (km/s)
Depth (km)
0.5 1 1.5 2 2.5 3 3.5 4 4.5
10
9
8
7
6
5
4
3
2
1
0

Figure 2-13: Inversion with events from the southeast (a) The 4-layer velocity models for the opposite
sides of the fault produced by inversion for a subset of events in the southeastern part of the study area.
(b) The velocity difference between the layers in all of the inversion runs used to obtain the average result
in (a). The notations are the same as in Figure 2-12.
71

The inversions of the arrival time data discussed so far give the depth variations
of laterally-averaged velocity contrasts between regions on the opposite sides of the
fault containing the employed events and stations. We now discuss waveform modeling
of the seismic phases in the early portions of the observed waveforms (Figure 2-14),
using a complementary model representing the depth-averaged lateral velocity
variations in the immediate vicinity of the fault. The synthetic calculations are
performed with the 2D analytical solution of Ben-Zion and Aki (1990) and Ben-Zion
(1998) for the scalar wave equation in a model having one or two vertical layers
between two quarter spaces. The source consists of an SH line dislocation with a unit
step in time and the solution assumes particle motion parallel to the structural interfaces,
allowing fault parallel and vertical synthetic to be generated. This results in the solution
being suitable only up to the direct P wave arrival, after which P-SV conversions and
additional phases become important parts of the seismograms.
The stations used in Figure 2-14 were chosen to have similar source-receiver
distances along the fault but different normal distances from the fault. The arrival times
and amplitudes of the direct P waves at the stations are controlled by the average
velocity and Q values of the crustal blocks on which they reside, while the properties of
the head waves at stations on the slow side depends on the attributes of both crustal
blocks. The relative times, amplitudes and motion polarities of the various phases
depend also on the propagation distance along the fault and the receiver offset from the
fault. These parameters are determined from the assumed values of the source and
receiver locations that were used in the previous arrival time analysis.

72
0.5 1 1.5 2 2.5
-3
-2
-1
0
1
2
3
4
5
6
FOR (6.70)
(4.00, 3.00)
TOP (7.42)
(4.00, 3.08)
SUM (6.61)
(4.05, 2.79)
PIT (5.63)
(3.95, 2.97)
AKC (4.30)
(4.69, 3.00)
STS (5.70)
(4.48, 3.00)
Time (sec)
Distance from fault (km)
a)
Direct P
wave
FZHW
0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.2
-1
-0.5
0
0.5
1
1.5
2
2.5
3
HVL (5.06)
(4.20, 2.81)
COW (4.92)
(3.86,2.80)
STK (5.25)
(4.00, 2.98)
KEY (5.68)
(4.21, 2.98)
YEL (4.60)
(4.28, 2.98)
BRG (5.03)
(4.72, 2.90)
Time (sec)
Distance from fault (km)
b)
FZHW
Direct P
wave
FZHW
Direct P
wave

Figure 2-14: Waveform fitting. Observed (blue) and synthetic (red) velocity seismograms for stations on
the fast and slow sides of the fault. The synthetic waveforms are calculated using the 2D analytical
solution of Ben-Zion (1998) for line dislocation with a unit step function in time at the interface between
two different quarter-spaces with possible vertical layer in between. Each pair of seismograms is labeled
with the station name followed by the hypocentral distance in parenthesis. Some of the direct P and
FZHW arrivals are indicated by arrows with corresponding labels. The P wave velocities in km/s of the
quarter-space used to generate the synthetic seismograms are given in parenthesis under the waveforms.
The employed Q values are 60 and 30 for the fast and slow sides of the fault, respectively. (a)
Seismograms associated with an event in the complex area to the SW of the array. (b) Seismograms
associated with an event to in the simpler central area.
73
To increase the spatial resolution of the results, the fit to each seismogram in Figure 2-
14 is produced individually with different structural parameters. In most cases, the
synthetic waveform fits to the early portions of the P waveforms require only two
quarter spaces. As discussed by McGuire and Ben-Zion (2005), for arrivals at stations
near the fault and events located in the region with the fault complexity, an intermediate
velocity layer is needed to account for the energy between the head and direct P waves.
Incorporation of such a layer in the model produces an additional fault zone phase that
is seen clearly in the both the observed and synthetic waveforms at station PIT in Figure
2-14. In the next section we discuss the results obtained by the arrival time inversions
and waveform modeling in the context of a physical model of a fault separating two
different crustal blocks and other available velocity images for the region.

2-4. Discussion
Figure 2-15 provides a schematic representation of the general features of the
near-fault velocity structure in the Cienega Valley to Bear Valley region of the SAF,
based on the arrival time inversions and waveform modeling done in this work. The
abundant observations of head waves imply that a sharp near-planar interface between
materials with differing properties exists on this portion of the fault. We further observe
a high velocity contrast in the shallow structure, a large increase in velocity around 3
km on the slower side of the fault, and a decreasing velocity contrast with depth. All the
sets of inversion consistently show high velocity contrasts (~50%) above 3 km, while
below this depth the variation across the fault reduces to about 10-20%. Between these
two depth ranges, the velocity on the slower side increases in a single large step of 2-2.5
74
km/s. Sets of inversions with 5, 6 and 8 layers were attempted to estimate how sharp or
gradual this transition is. However, the inversion results always produce profiles of the
type shown in Figure 2-10. This feature of the results remains largely unchanged within
the inversion sets (B1)-(B3), implying that it is not a structure associated with the
immediate vicinity of the SAF such as a low velocity trapping structure.

Figure 2-15: A schematic diagram of the inferred velocity structure in the study area, consisting of two
layered quarter-spaces, joined along a sharp material interface, and a shallow asymmetric low velocity
zone around the fault. The colors represent the media velocities, with violet being the fastest and orange
being the slowest.
The arrival time inversion of Thurber et al. (1997) shows a region of 3-4 km/s P
wave velocity extending 6 km to the NE of the fault and down to 3-4 km depth. This
75
region may represent the same structural feature seen in our inversions as a large
velocity step at 3km depth on the slow side. The sharp velocity increases in our
inversions probably reflect geological transitions, perhaps exaggerated somewhat by the
assumed layered structure and limited inversion resolution. Additional data generated
by events around 3 km depth could lead to a more gradual velocity changes with depth.
Other possible structural elements that could explain the top section include a region of
high fluid content, and thus presumably cracks and damage, that was imaged during a
magnetotelluric study down to a depth of 2-4 km in the 6 km interval between the SAF
and the Calaveras fault to the northeast (Bedrosain et al., 2004).
The inversions also show finer results, such as structural complexities in the
fault zone and the possible existence of an asymmetry of rock damage across the fault.
To the southeast of the array the fault splits into two distinct branches, with shallow
seismicity on the east branch and at approximately 6 km depth a switch to seismicity
located on the more westerly branch (Figure 2-4). Based on waveform and travel time
analysis, combined with geological surface mapping, McGuire and Ben-Zion (2005)
conclude that this represents a sliver of granite between the two fault branches
separating the granites and volcanic rock of the southwest (fast side) from the
sedimentary rocks on the northeast (slow) side of the fault. The inversion results based
on the 150 events to the southeast (set C2) that are within or have paths through this
complexity show at 6-8 km depth a reduction in velocity (Figure 2-13). This
corresponds to the depth at which the seismicity changes from one fault branch to the
other. The character of waveforms that pass through this area is also different (Figure 2-
5), with additional phases occurring between the FZHW and direct P arrivals. The
76
synthetic waveform modeling of Figure 2-14 (see also McGuire and Ben-Zion, 2005)
demonstrates that the additional phases can be explained by a vertical layer with
intermediate velocity between those of the bounding quarter-spaces, which represents
the sliver of granite.
While it is interesting that our simple model can distinguish features on the scale
of approximately 200 m, and further confirm the existence of the complexity in the
structure, it can provide little more quantitative results on the nature of the complexity
without becoming more complex itself. The mapping of event locations in the model
onto a single interface affects the accuracy of the results, primarily for the structure on
the fast side of the fault. The propagation distances from the events on the deeper fault
branch to the stations on the fast side are in reality smaller than those used in the model
with dataset C2 (Figure 2-13). This produces somewhat unrealistically high velocities
for the deepest layer of ~7.2-7.4 km/s compared to the ~6.6-6.8 km/s for the same depth
range obtained by Thurber et al. (1997). This would also impact the inversion results
using the set A containing all the data, so the velocities of the deepest layers are
probably slightly increased. The velocity below 4 km of 6.2-6.4 km/s in Figure 2-12
from inversions using the dataset C1 that does not contain the branched section of the
fault might provide the closest representation of the true structure.
When using subsets of the data with increasing proportions of stations close to
the fault (sets B1-B3), the inversions show increasingly clearer near-surface low
velocity (damage) structures, with the low velocity surface layer on the faster side more
readily evident at a greater range of distances. This is seen by examining the results
obtained by 4 sets of inversions (Figure 2-11a to 2-11d), and is also reflected by plotting
77
(Figure 2-11e) the obtained depth and velocity of the shallowest layer on each side of
the fault. As shown in Figure 2-11e, these two defining parameters of the surface layers
remain largely unchanged on the slower side of the fault (triangles) until only stations
within 2 km are used, whereas on the faster side (circles) the low velocity surface layer
grows systematically. There is no apparent consistent variation in the properties of any
of the deeper layers with changes in proximity to the fault. The larger the fraction of
near fault stations that are used in the inversions, the more the images of the shallow
layers reflect the velocity structure of the near fault region, particularly as there is an
increase in the relative number of head wave arrivals. Since the results characterizing
the shallowest model layers are correlated clearly with the proximity of the employed
receivers from the fault, it is reasonable to assume that the imaged low velocity material
represents a damage fault zone layer like that seen in trapped wave studies.
The average velocity profile for the upper 3 km should be the most well
constrained as all the rays from all the events pass through it. However since this zone
resides above the shallowest seismicity (2-3 km), the internal variations within the
shallow low velocity layers can not be resolved, except to say that they do not extend
below 3 km depth. We note that we obtain good fits for the arrival time, polarities and
amplitude of the head and direct waves (Figure 2-14) using quarter-space velocities of
~3 km/s on the slow side and ~4.5 km/s on the fast side. These values are broadly
consistent with the near fault velocities at shallow depth obtained in the arrival time
inversions and also those from Thurber et al. (1997). In the tomographic images of
Thurber et al. (1997), the fault zone may be associated with a vertical structure in the
P/S wave velocity ratio. It is interesting to note that this vertical structure in the results
78
of Thurber et al. (1997) is also offset to the SW from the surface trace of the fault, in
agreement with the asymmetry of the shallow structure in the results of Figure 2-11.
The same asymmetry is also evident in the parameters used to obtain the
synthetic waveform fits of Figure 2-14 for the two shallow events with depths of ~3 km
at a number of near fault stations. The shallow events are in the depth range where the
structure may be well represented by two quarter spaces separated by a low velocity
layer. As indicated in Figure 2-14, the velocities required to fit the direct P arrival times
generally reduce with proximity to the fault. The gradient of the velocity reduction with
increasing proximity to the fault is particularly strong on the faster side. Lower
velocities on the fast side of the fault are also required to obtain the correct head wave
arrival times, implying that the P wave velocities at the material interface are lower
again than at stations within a few hundreds of meters SW of the fault. These
observations imply, like the changes in the results of Figures 2-11, that the shallow low
velocity layer is not centered on the SAF, but is rather asymmetrically located mostly
on the faster side of the fault.
A shallow asymmetric damage zone, with more damage on the faster side of the
fault, is an expected outcome of ruptures on an interface that separates different elastic
media (Ben-Zion and Shi, 2005). Such ruptures tend to propagate preferentially and/or
more vigorously in the direction of slip on the slower side of the fault (e.g., Weertman,
1980; Ben-Zion, 2001; Shi and Ben-Zion, 2006; Rubin and Ampuero, 2007), and they
produce damage-generating tensile radiation consistently on the faster side. Significant
damage generation is limited to the top few km of the crust (Ben-Zion and Shi, 2005),
since increasing normal stress suppresses the damage generation. Our inference that the
79
observed asymmetric shallow structure is associated with a preferred propagation
direction of earthquakes on the SAF south of Hollister is compatible with observations
of Rubin and Gillard (2000) of asymmetric along-strike distribution of aftershocks in
this section of the fault. We also note that similar asymmetric damage zones, which are
offset to the faster side of the fault, were recently observed in inversions of trapped
waves on the San Jacinto fault (Lewis et al., 2005) and detailed geological mapping in
the structure of several faults of the southern San Andreas system (Dor et al., 2006a,b).
The asymmetry of the shallow fault zone structures may be associated with other
features and mechanisms, such as local topography and various site effects. However, at
present the most consistent explanation for the available diverse observations appears to
be a preferred propagation direction of earthquake ruptures.
As noted by Ben-Zion et al. (1992) and in section 3, analyses of fault zone head
and trapped waves complement each other in high resolution imaging of fault zone
structures. The head waves can provide depth profiles of the overall velocity contrast
across the fault, while the trapped waves can be used to obtain the depth-averaged
lateral variations of the structure. The sensitivity of near-fault head waves to the internal
FZ structure can help to reduce some of the non-uniqueness in trapped-wave inversions.
Additional important information on the material and geometrical properties of fault
zones can be obtained using other potential signatures of damaged rocks such as
anisotropy and scattering (e.g., Boness & Zoback, 2004; Liu et al., 2004, 2005; Cochran
et al., 2006; Peng & Ben-Zion, 2005, 2006). A combined use of phases generated by
material interfaces and damaged FZ rocks, direct body waves, and accurate absolute
event locations, has the potential to provide unparalleled detailed images of fault zone
80
structures and their surrounding environment. Our results indicate that the inclusion of
fault zone head waves in inversions of seismic data for the structure of large fault zones
can increase considerably the resolution of the velocity contrast across the fault as a
function of depth. We also note that proper accounting for fault zone head waves can
improve the accuracy of focal mechanisms determined from data collected by near-fault
arrays (e.g., Ben-Zion, 1990; Kilb et al., 2006).
81
Chapter 2 References
Aki, K., and W.H.K. Lee, (1976). Determination of three-dimensional velocity
anomalies under a seismic array using first P arrival times from local earthquakes,
Part 1, A homogeneous initial model, J. Geophys. Res., 81, 4381-4399.
Andrews, D. and Y. Ben-Zion, (1997). Wrinkle-like slip pulse on a fault between
different materials, J. Geophys. Res., 102, 55-571.
Bakun, W.H., R.M. Stewart, C.G. Bufe and S.M. Marks, (1980). Implications of
seismicty for failure of a section of the San Andreas fault, Bull. Seim, Soc. Am.,
70, 185-201.
Bedrosain, P.A., M.J. Unsworth, G.D. Egbert and C.H. Thurber, (2004). Geophysical
images of the creeping segment of the San Andreas Fault: implications for the role
of crustal fluids in the earthquake process, Tectonophysics, 385, 137-158.
Ben-Zion, Y., (1989). The response of two joined quarter spaces to SH line sources
located at the material discontinuity interface, Geophys. J. Int., 98, 213-222.
Ben-Zion, Y., (1990). The response of two half spaces to point dislocations at the
material interface, Geophys. J. Int., 101, 507-528.
Ben-Zion, Y., (1998). Properties of seismic fault zone waves and their utility for
imaging low velocity structure, J. Geophys. Res., 103, 12,567-12,585.
Ben-Zion, Y., (2001). Dynamic ruptures in recent models of earthquake faults, J. Mech.
Phys. Solids, 49, 2209-2244.
Ben-Zion, Y. and K. Aki, (1990). Seismic radiation from an SH line source in a laterally
heterogeneous planer fault zone, Bull. Seism. Soc. Am., 80, 971-994.
Ben-Zion, Y. and P. Malin, (1991). San Andreas Fault zone head waves near Parkfield
California, Science, 251, 1592-1594.
Ben-Zion, Y. and D.J. Andrews, (1998). Properties and implications of dynamic rupture
along a material interface, Bull. Seism. Soc. Am., 88, 1085-1094.
Ben-Zion, Y. and Z. Shi, (2005). Dynamic rupture on a material interface with
spontaneous generation of plastic strain in the bulk, Earth Planet. Sci. Lett., 236,
486-496, DOI: 10.1016/j.epsl.2005.03.025.
Ben-Zion, Y., S. Katz and P. Leary, (1992). Joint inversion of fault zone head and direct
P arrivals for crustal structure near major faults, J. Geophys. Res., 97, 1943-1951.

82
Ben-Zion, Y., Z. Peng, D. Okaya, L. Seeber, J.G. Armbruster, N. Ozer, A.J. Michael, S.
Baris and M. Aktar, (2003). A shallow fault zone structure illuminated by trapped
waves in the Karadere-Duzce branch of the north Anatolian Fault, Western
Turkey, Geophys. J. Int., 152, 699-717.
Boness, N.L. and M.D. Zoback, (2004). Stress-induced seismic velocity anisotropy and
physical properties in the SAFOD Pilot Hole in Parkfield, CA, Geophys. Res.
Lett., 31, L15S17, doi:10.1029/2004GL019020.
Cochran, E.S., Y-G. Li and J.E. Vidale, (2006). Anisotropy in the Shallow Crust
Observed around the San Andreas Fault Before and After the 2004 M 6.0
Parkfield Earthquake, Bull. Seism. Soc. Am., 96, S364–S375, doi:
10.1785/0120050804.
Dibblee, T.W., (1974). Geologic map of the Salinas quadrangle, California, in Open-file
report 74-1021, U.S. Geological Survey.
Dor O., T.K. Rockwell and Y. Ben-Zion, (2006a). Geologic observations of damage
asymmetry in the structure of the San Jacinto, San Andreas and Punchbowl faults
in southern California: A possible indicator for preferred rupture propagation
direction, Pure Appl. Geophys., 163, 301-349, DOI 10.1007/s00024-005-0023-9.
Dor O., Y. Ben-Zion, T. K. Rockwell and J. Brune, (2006b). Pulverized Rocks in the
Mojave section of the San Andreas Fault Zone, Earth Planet. Sci. Lett., 245, 642-
654.
Eberhart-Phillips, D. and A.J. Michael, (1993). Three-Dimensional velocity structure,
seismicity and fault structure in the Parkfield region, central California, J.
Geophys. Res., 98, 15737-15758.
Ellsworth, W., (1975). Bear valley, California, earthquake sequence of February-March
1972, Bull. Seism, Soc. Am. 65, 483-506.
Igel, H., G. Jahnke and Y. Ben-Zion, (1997). Simulation of SH and P-SV wave
propagation in fault zones, Geophys. J. int., 128, 533-546.
Jahnke, G., H. Igel and Y. Ben-Zion, (2002). Three-dimensional calculations of fault
zone guided waves in various irregular structures, Geophys. J. Int., 151, 416-426.
Kilb, D., and J.L. Hardebeck, (2006). Fault Parameter Constraints Using Relocated
Earthquakes: A Validation of First Motion Focal Mechanism Data, Bull. Seism.
Soc. Am. 96, 1140-1158.
Korneev, V.A., R.M. Nadeau and T.V. McEvilly, (2003). Seismological studies at
Parkfield IX: Fault zone imaging using guided wave Attenuation, Bull. Seism.
Soc. Am., 93, 1415-1426.
83
Lewis, M.A, Z. Peng, Y. Ben-Zion and F. Veron, (2005). Shallow seismic trapping
structure in the San Jacinto Fault Zone, California, Geophys. J. Int., 162, 867–881,
doi:10.1111/j.1365-246X.2005.02684.x.
Li, Y.G., J.E. Vidale and K. Aki, (1997b). San Jacinto fault-zone guided waves: A
discrimination for recently active fault strands near Anza, California, J. Geophys.
Res., 102, 11,689-11,701.
Li, Y.G., P. Leary, K. Aki, and P.E. Malin, (1990). Seismic trapped modes in the
Oroville and San Andreas fault zones, Science, 249, 763-766, 1990.
Li, Y.G., K. Aki, D. Adams, A. Hasemi and W.H.K Lee, (1994). Seismic guided waves
trapped in the fault zone of the Landers, California, earthquake of 1992,
J.Geophys. Res., 99, 11,705-11,722.
Li, Y.G., W.L. Ellsworth, C.H. Thurber, P. Malin and K. Aki, (1997a). Fault-zone
guided waves from explosions in the San Andreas fault at Parkfield and Cienega
Valley, California, Bull. Seismol. Soc. Am., 87, 210-221.
Liu, Y., T.L. Teng and Y. Ben-Zion, (2004). Systematic analysis of shear wave splitting
in the aftershock region of the 1999 Chi-Chi earthquake: Evidence for shallow
anisotropic structure and lack of systematic temporal variations, Bull. Seism. Soc.
Am., 94, 2330-2347.
Liu, Y., T.L. Teng and Y. Ben-Zion, (2005). Near-surface seismic anisotropy,
attenuation and dispersion in the aftershock region of the 1999 Chi-Chi
earthquake, Geophys. J. Int., 160(2), 695-706.
McGuire, J., and Y. Ben-Zion, (2005). High-resolution imaging of the Bear Valley
section of the San Andreas Fault at seismogenic depths with fault-zone head
waves and relocated seismicity, Geophys. J. Int., 163, 152-164, doi:
10.1111/j.1365-246X.2005.02703.x.
Marra, F., R. Azzara, F. Bellucci, A. Caserta, G. Cultrera, G. Mele, B. Palombo, A.
Rovelli and E. Boschi, (2000). Large amplification of ground motion at rock sites
within a fault zone in Nocera Umbra (Central Italy), J. of Seis., 4 534-554.
Michael, A.J. and Y. Ben-Zion, (1998). Inverting fault zone head trapped waves with
genetic algorithm, EOS, Trans. Am. Geophys. Un., 81, F1145.
Peng, Z., Y. Ben-Zion, L. Zhu and A.J. Michael, (2003). Inference of a shallow fault
zone layer in the rupture zone of the 1992 Landers, California earthquake fro
locations of events generating trapped waves and travel time analysis, Geophys. J.
Int., 155, 1021-1041.

84
Peng, Z. and Y. Ben-Zion, (2005). Spatio-temporal variations of crustal anisotropy from
similar events in aftershocks of the 1999 M7.4 İzmit and M7.1 Düzce, Turkey,
earthquake sequences, Geophys. J. Int., 160(3), 1027-1043, doi: 10.1111/j.1365-
246X.2005.02569.x.
Peng, Z. and Y. Ben-Zion, (2006). Temporal changes of shallow seismic velocity
around the Karadere-Duzce branch of the north Anatolian fault and strong ground
motion, Pure Appl. Geophys., 163, 567-599, DOI 10.1007/s00024-005-0034-6.
Rovelli, A., A. Caserta, F. Marra and V. Ruggiero, (2002). Can seismic waves be
trapped inside an inactive fault zone? The case study of Nocera Umbra, central
Italy, Bull. Seis. Soc. Am., 92, 2217-2232.
Rubin, A.M., and D. Gillard, (2000). Aftershock asymmetry/rupture directivity among
central San Andreas fault micro earthquakes, J. Geophys. Res., 105, 19095-19110.
Rubin, A.M. and J.P. Ampuero, (2007). Aftershock asymmetry on a bimaterial
interface, J. Geophys. Res. 112, B05307, doi:10.1029/2006JB004337.
Scott, J. S., T.G. Masters and F.L. Vernon, (1994). 3-D velocity structure of the San
Jacinto fault zone near Anza, California -- I. P waves, Geophys. J. Int., 119, 611-
626.
Shi, Z. and Y. Ben-Zion, (2006). Dynamic rupture on a bimaterial interface governed by
slip-weakening friction, Geophys. J. Int., 165, doi: 10.1111/j.1365-
246X.2006.02853.x.
Thurber, C., S. Roecker, W. Ellsworth, Y. Chen, W. Lutter and R. Sessions, (1997).
Two-dimensional seismic image of the San Andreas fault in the northern Gabilan
range, central California: Evidence for fluids in the fault zone, Geophys. Res.
Lett., 24, 1591-1594
Waldhauser, F. and W. Ellsworth, (2000). A double-difference earthquake location
algorithm: method and application to the northern Hayward Fault, Bull. Seism.
Soc Am., 90, 1,330-1368.
Weertman, J., (1980). Unstable slippage across a fault that separates elastic media of
different elastic constants, J. Geophys. Res., 85, 1455-1461.
85
Chapter 3: Examination of Scaling Between Proposed Early Signals in
P Waveforms and Earthquake Magnitudes
Summary
We examine seismic waveforms generated by clusters of similar earthquakes
and additional events located on the Karadere-Duzce branch of the North Anatolian
fault to investigate possible scaling between proposed early signals in the P waveforms
and magnitudes. The employed signals can be divided into two classes. The first is
associated with Ellsworth-Beroza and Iio type signals that are assumed to reflect
signatures of seismic nucleation phases. The second are associated with Nakamura type
period and Wu-Zhao type maximum displacement in the early waveforms that are
related to measures of the corner frequency and local magnitude, respectively. The
impact on each signal of artifacts generated by the Finite Impulse Response (FIR) filter
is also examined. The use of repeating event clusters reduces considerably the
variability of path and site effects compared to previous regional and global studies.
This allows us to focus more accurately on variations in the potential scaling parameters
with the source magnitude. However, the magnitude range of the employed events is
limited to about 4 units. The FIR filter artifacts can generate spurious signals prior to
the P wave arrival, potentially affecting the Ellsworth-Beroza type signals, but are not a
significant influence on any of the other measurements. The results indicate that
candidates for the Ellsworth-Beroza signals, associated with identification of a weak
arrival before the main P phase, exist in less than 20% of the cases. The Iio type signals
can be determined by definition on all waveforms. Both of these potential signatures of
nucleation phases show little or no scaling with the final event size in the examined
86
data. In contrast, the Nakamura and Wu-Zhao type signals in the first few seconds,
which can also be determined by definition on all waveforms, appear to scale with the
final event magnitude, albeit with large scatter. The relatively large time windows used
in the analysis, compared to the rupture time of the employed events, implies that the
results provide little if any information on the physical process of the rupture initiation.
The scaling of the Nakamura and Wu-Zhao type early signals with the final event size
may be explained by statistical tendency of stronger initial rupture phases to propagate
larger distances.

3-1. Introduction
The possibility that the early portions of seismograms can be used to determine
the final event size is alluring because it implies some regular scaling in the physics of
earthquake rupture and the prospect of developing early-warning systems. Theoretical
considerations suggest the following three general types of cases. If ruptures occur on a
smooth frictional surface with homogeneous properties, the nucleation size (e.g.,
Dieterich, 1992; Ben-Zion and Rice, 1997; Lapusta and Rice, 2003) of all brittle
instabilities will be the same, and hence scaling between the nucleation process and
final event size is not expected. When ruptures occur on a smooth fault with
heterogeneous frictional properties, there will be a range of nucleation sizes. In such
cases, larger nucleation phases are likely to produce larger events, as long as the
growing ruptures do not encounter strong correlated heterogeneities like the free surface
or large segment boundaries (e.g., Hillers et al., 2006, 2007). Finally, if ruptures occur
on highly disordered fault zones with strong geometrical heterogeneities, the failure
87
process consists of cascading sub-events and the event size is influenced strongly by the
highly variable stress-strength conditions along the fault (e.g. Ben-Zion et al., 2003). In
the latter case, the final rupture size is not expected to be affected by the nucleation
process, but may still scale with the strength of the early sub-event at the hypocentral
region.
Iio (1995), Ellsworth and Beroza (1995) and Beroza and Ellsworth (1996)
suggested that properties of the initial pulse in seismograms that might reflect signatures
of nucleation phases scale with the final event size. These results have not been
supported by other studies (e.g., Anderson and Chen, 1995; Mori and Kanamori, 1996;
Abercrombie and Mori, 1996; Ishihara et al. 1992; Kilb and Gomberg 1999; Sato and
Mori 2006). However, recent observations by Olson and Allen (2005), Wu and Zhao
(2006) and Zollo et al. (2006), indicating that the final event magnitude may be at least
partially determined prior to the cessation of rupture, have re-awakened the search for
informative signals in the early portions of seismograms. Most of the previous studies
that search for scaling with M of features from the early P waveforms have used a
relatively small number of seismograms (e.g., N < 50) generated by events covering a
wide range of magnitudes (2 < M < 8). The employed earthquakes occurred in a large
region and/or many different tectonic domains. The inconsistent conclusions obtained in
the different studies may in part result from the small number of the employed
waveforms and the amalgamation of data from different domains.
In the present paper we analyze systematically a waveform dataset (N=190) for
several proposed signals in the early waveforms. The data were generated by clusters of
88
repeating earthquakes on the Karadere-Duzce branch of the North Anatolian Fault
(Figure 3-1) in the magnitude range 0 < M < 2.75.
0 5 10
km
1
2
3
5
7
9
10
17
19
22
30
31
FP
GE
LS
FI
MO VO
WF
BU
BV
CH
Izmit M 7.4
Dzce M 7.1
Almacik Block
Akyazi
Hendek

Adapazari
1943 M6.4
1957 M7.0
Depth range
0 ~ 5 km
5 ~ 10 km
10 ~ 15 km
> 15 km
0 50 100
km
Marmara
Sea
Istanbul
Black Sea
NAF
M 7.4
M 7.1
i
Eften Lake
00'
40° 40'
40° 50'
41°
30°30' 30° 40' 30° 50' 31° 00' 31° 10'
Dûzce

Figure 3-1: A location map of the study area (square box in the inset) around the Karadere-Dûzce branch
of the North Anatolian fault. The red and gold lines are the mapped ruptures of the 1999 Izmit and Dûzce
earthquakes, while the thin black lines show other mapped faults in the area. The circles give the average
location of 13 employed repeating earthquake clusters with depth shown by color and labeled by numbers
from Peng and Ben-Zion (2006). Additional used events are marked by diamonds. The employed 10
seismic stations labeled with 2 letter codes are marked by grey triangles and near by towns are denoted
with black squares.
For parts of the analysis this magnitude range is extended by using non-repeating events
from the same region. Choosing to primarily use clusters of similar events allows us to
minimize the effects of variations in path and site conditions, and thus focus on source
effects. The analysis examines several types of proposed signals in the early portion of
the P waveforms. The first type, related to the corner frequency, involves estimating the
predominant period (or frequency) in the early part of the seismogram (Nakamura,
1988; Nakamura and Tucker, 1988; Kanamori 2003, 2005; Allen and Kanamori 2003;
89
Olson and Allen 2005). The second type attempts to target signals that might reflect the
nucleation phases of earthquakes. These consist of the Iio method that measures the
slow initial rise of the P wave (Iio 1995), and the Ellsworth-Beroza method that
measures the duration between a low amplitude initial phase arrival P
1
and the more
impulsive main arrival P
2
(Ellsworth and Beroza 1995). The third type, related to
measures of the magnitude, follows the work of Wu and Zhao (2006) and Zollo et al.
(2006) in measuring the peak amplitude of displacement in the early part of the
seismogram.
In the next section we describe the general properties of the earthquakes,
stations and waveform data used in this study. In section 3 we describe the employed
analysis methods. The results are presented in section 4 and discussed in section 5. The
analysis indicates that the seismic nucleation signature proposed by Ellsworth and
Beroza (1995) is not a common observable feature of the analyzed events at the
employed sampling rate and hypocenter distances > 4 km. The early phases of the type
used by Iio (1995) and Ellsworth and Beroza (1995) show very limited scaling with the
final even size. The predominant period in the early waveforms determined by the
method of Olson and Allen (2005), as well as the peak amplitude of displacement (Wu
and Zhao, 2006), appear to scale with M within large fluctuations. Since the
measurements are done with larger time windows than the rupture duration of the used
small events, the results do not shed, on their own, light on the initiation process of the
ruptures.

90
3-2. Data
3-2.1 The seismic experiment
A week after the August 17, 1999, Mw7.4 Izmit earthquake, a 10-station
PASSCAL seismic network (Figure 3-1) was deployed along and around the Karadere-
Duzce branch of the North Anatolian fault in Turkey (Seeber et al., 2000; Ben-Zion et
al., 2003). All 10 stations had REFTEK recorders and three-component L22 short-
period sensors with a sampling frequency of 100 Hz. The Karadere-Duzce branch was
chosen for the deployment as uniquely for this region the rupture broke the surface in
bedrock, avoiding sedimentary layers that diffuse, attenuate and mask wave propagation
signals. Three months after the deployment, the November 12, 1999, Mw7.2 Duzce
earthquake started and propagated eastward from the Karadere-Duzce segment. As a
result of the aftershocks of the large Izmit and Duzce earthquakes, about 26,000 events
were recorded in the ~6 month duration of the seismic experiment.
Peng and Ben Zion (2005, 2006) performed cross correlations on waveforms of
~18,000 of the events that are within 20 km of the Izmit and Duzce rupture zones and
recorded by at least 3 stations. The analysis resulted in organization of 292 of the events
into 36 repeating earthquake clusters, each containing at least 5 events having median
waveform correlation coefficient of ≥ 0.95. We choose 13 of these clusters (circles in
Fig. 1) that contain the largest and smallest event sizes, the largest magnitude ranges
and most events. This leaves 1141 recorded waveforms generated by 142 events for the
analysis of variations of features in the early portions of the P waveforms with
magnitude. Figure 3-2 shows an example of the waveforms recorded at station MO for
cluster 5.
91
0 1 2 3 4 5

Time (sec)

Figure 3-2: Similar velocity waveforms generated by the repeating events in cluster 5 and recoded at
station MO. The seismograms have normalized amplitudes and are aligned on the P arrivals (black
dashed line) at 0 seconds.
The waveform shapes are highly similar in the P, S, and coda waves, implying that the
events have highly similar locations and focal mechanisms. However, the peak
amplitudes of the waveforms are different, indicating that the magnitude and possibly
other source properties (e.g., stress drop) of the different earthquakes are not the same.
The similarity of the waveforms generated by the repeating earthquake clusters,
and the related ability to accurately identify the same features in the waveforms of
events with different magnitudes, should enable a more accurate study of features that
might scale with the event size. Furthermore, any differences between features in the
early waveforms of events in a cluster, recorded at the same station, should be due to
their different size (or other source effects), as the propagation and site effects for these
92
events are nearly identical. One shortcoming of the analysis of waveforms associated
with repeating earthquake clusters, compared to previous works on the topic, is that the
maximum magnitude range in any single cluster is ~2. If the scaling with M of the
signals in the early P waveform is weak this may make its observation more difficult or
even impossible. To expand on the available magnitude range for some facets of the
analysis, larger non-cluster events from the same region are used. This adds additional
48 events (diamonds in Fig. 1) and 243 waveforms with M >2.8. Combining the results
of the different clusters and data from the larger events increases the magnitude range
but leads to some of the more usual mixing of source, path and site effects.
3-2.2 Site and Instrument effects
Complexity in the P waveforms can be generated by source, propagation, site
and recording instrument effects. Four of the stations (VO, LS, MO and BV) are located
tens to hundreds of meters from the rupture zones of two major earthquakes on a plate
boundary fault, and as such are expected to record seismic phases associated with fault
zones. Such seismic phases include fault zone trapped waves (Ben-Zion et al., 2003),
fault related anisotropy (Peng and Ben-Zion, 2004, 2005) and fault zone head waves
(e.g., Ben-Zion and Malin 1991; Lewis et al., 2007). Fault zone head waves are
refracted phases that propagate along sharp contrast of rock types between two sides of
a fault, similar to head wave refractions along the moho discontinuity and other material
interfaces. Fault zone head waves can be the first arriving phase for stations near the
fault on the side with the lower seismic velocity, and are characterized by emergent
low-amplitude phase followed by larger-amplitude sharper P body waves (Ben-Zion
93
1989). The head waves could impinge on this study because their properties are similar
to those expected from an accelerating preparatory or nucleation phase, and their
amplitude (and hence size of region where they exist) can scale with the earthquake
magnitude. Previous studies using our dataset did not attempt to identify head waves,
and it is not clear whether the Karadere-Duzce faults have significant material interfaces
that can generate such phases, but they could be present in the data
Another important class of possible interference with genuine source effects is
the type of filter used in the seismic recording systems. The commonly-used FIR filter
can generate precursory signals to impulsive arrivals that scale with M and could be
misidentified as nucleation phases (Scherbaum and Bouin, 1997). This effect is
illustrated in Figure 3-3, where we show several views of waveforms generated by all
the events in cluster 3 at station WF. Panels (a) and (b) show the original seismograms
with fairly obvious FIR artifacts between the two vertical lines. Panel (c) and (d) show
waveforms in which the FIR effects are removed (Scherbaum and Bouin, 1997;
Scherbaum, 1996), producing simple and very impulsive onsets of the P waves.
94
-1 -0.5 0 0.5 1 1.5
1
1.5
2
2.5
3
3.5
-0.2 -0.1 0 0.1 0.2 0.3 0.4
1.6
1.8
2
2.2
2.4
2.6
Time (sec)
Time (sec)
-1 -0.5 0 0.5 1 1.5
1
1.5
2
2.5
3
3.5

-0.2 -0.1 0 0.1 0.2 0.3 0.4
1.6
1.8
2
2.2
2.4
2.6
Time (sec)
Time (sec)
a) c)
b)
Event magnitude Event magnitude
d)
Figure 3-3: FIR artifacts. Normalized velocity seismograms generated by the events in cluster 3 and
recorded at station WF. The waveforms are aligned at the P wave arrivals (dashed vertical lines) and
plotted at their corresponding magnitude. The solid vertical lines mark times of earlier weak phases of the
type used by Ellsworth and Beroza (1995), when such can be identified on the seismograms. (a) Original
waveforms. (b) Enlarged waveforms around the P arrivals (c) Same as (a) with FIR filter artifacts
removed. (d) Same as (b) with FIR filter artifacts removed.
Scherbaum and Bouin (1997) showed that FIR artifacts are not always easily-
recognizable monochromatic oscillatory signals (as in Figure 3-3), but can be complex
and may be misidentified as the Ellsworth-Beroza type signals. Figures 3-4 and 3-5
show waveforms generated by an event recorded at a number of stations with and
without the FIR filter removed. For these records, the P
1
arrivals before the main P
2

arrivals are only slightly modified by the removal of the FIR filter (Figure 3-5).
Therefore filter artifacts alone cannot explain all the complexity in the early P
waveforms and the P
1
arrivals remain candidate signatures for source, path or site
effects. While Scherbaum and Bouin (1997) discussed the influence of the FIR filters on
the Ellsworth-Beroza type nucleation phase, spurious signals associated with the filter
could also affect quantities associated with the Iio-type and Nakamura-type signals.
Since the FIR filter should not affect the amplitude, however, there is no reason to
believe there would be any effect on the Wu-Zhao type signals. To carry out a thorough
analysis of the P waves onset, and to determine the FIR effects on various early signals,
95
we perform the majority of our analysis both with and without the acausal response of a
zero-phase FIR filter removed.
1
-1
0
1
-1
0
1
-1
0
1
-1
0
1
-1
0
1
-1
0
-1 -0.5 0 0.5 1 1.5 2 2.5 3
Time (sec)
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.02
-0.02
0
-0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3
Time (sec)
MO epdist:16.4 km baz:48.8
FP epdist:12.0 km baz:42.8
FI epdist:16.8 km baz:48.1
LS epdist:18.8 km baz:52.1
VO epdist:11.9 km baz:40.9
WF epdist:16.5 km baz:35.3
Figure 3-4: Velocity seismograms before the removal of the FIR filters, generated by an event in cluster
30 at 6 stations. The identification code, epicentral distance and back azimuth of the stations are written
in the rectangles on top. The vertical solid and dashed lines mark the picked P
1
and P
2
arrivals,
respectively.
96
-0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3
Time (sec)
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.2
-0.2
0
0.02
-0.02
0
-1 -0.5 0 0.5 1 1.5 2 2.5 3
Time (sec)
1
-1
0
1
-1
0
1
-1
0
1
-1
0
1
-1
0
1
-1
0
MO epdist:16.4 km baz:48.8
FP epdist:12.0 km baz:42.8
FI epdist:16.8 km baz:48.1
LS epdist:18.8 km baz:52.1
VO epdist:11.9 km baz:40.9
WF epdist:16.5 km baz:35.3
Figure 3-5: The same as Figure 3-4 but with the removal of the acausal response of the FIR filter. The
two arrivals P
1
and P
2
are still clearly apparent.
3-3. Methods
The Nakamura-type methods calculate the predominant period in the early
portion of the P wave, their original purpose being to provide fast estimates of
magnitude for earthquake early-warning systems. In one version, the predominant
period is defined as:

97
∫
∫
=
0
0
0
2
2
0
) (
) (
τ
τ
dt t u
dt t u
r
&
, (1a)

and

r
c
π
τ
2
= , (1b)
where τ
0
is set at three seconds and u is the ground displacement (Kanamori,
2004). A related procedure by Olson and Allen (2005) employs similar iterative
calculations for each point in the seismogram, and the maximum value after the P
arrival within four seconds is used as the final result. Here the predominant period at
time step i is:

i
i
D
X
p
i
π τ 2 = (2a)

where

2
1 i i i
x X X + =
−
α (2b)

and

98
2
1
i
i i
dt
dx
D D






+ =
−
α . (2c)

In (2b) and (2c) x
i
is the recorded ground velocity at time step i and α is a
smoothing constant, dependent on the sampling rate, set here at 0.99 as in Olson and
Allen (2005). The maximum value of τ
p
within 4 seconds of the P arrival, τ
pmax
, is the
factor (Figure 3-6) that is proposed to scale with M. The time after the P wave where
τ
pmax
occurs varies and is generally greater for larger events. For earthquakes within the
magnitude range of our study this is always much less than 4 seconds.
MO epdist: 2.8 km baz: 325 mag: 1.34
20
10
15
5
-0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4
Velocity
τ
p
(s)
Time (sec)
τ
pmax

Figure 3-6: A seismogram and the calculated values of τ
p
A vertical velocity seismogram (top) aligned
with the P wave at time 0, and the iteratively calculated value τ
p
using the method of Allen and Olson
(2005) (bottom). The maximum calculated value between 0.05 seconds and 4 seconds after the P wave
arrival is denoted τ
pmax
and used in the analysis. The station identification, the epicentral distance, back
azimuth, and magnitude of the event are shown on top.
99
The flexibility in the time window is presumably what allows this method to more
effectively scale with M compared with the Kanamori (2004) method that uses a fixed 3
second time window. For the examined small events, the time window over which τ
pmax

is calculated is far greater than the rupture duration. However, Allen and Olson (2005)
reported a continuous scaling with magnitude of their predominant period for events up
to M ~ 8. This means that for earthquakes larger than M ~ 5, the final event size is
determined from the predominant period before the rupture was arrested. It should be
noted, however, that the quality of the fit between the trend and the data in the study of
Allen and Olson (2005) is reduced for the higher magnitude events. Kanamori (2005)
has similar results with higher values of τ
c
correlated with increasing magnitudes for
events up to M ~ 8, but this is again based on a small number of large magnitude events.
If these results stand up as more data are analyzed, it would imply that M can be
established before the rupture has terminated and thus, at least for large earthquakes, the
initial part of the rupture has some control on the final event size.
Wu and Zhao (2006) and Zollo et al. (2006) observed relations between peak
amplitude in seismograms and M which, like the Nakamura-type methods, are proposed
for use in earthquake early-warning systems. In the procedure that we call the Wu-Zhao
method, the original velocity seismograms are integrated to displacement and then a
high pass filter at 0.075Hz is used to remove any low frequency drift that resulted from
the numerical integration. A peak displacement parameter Pd is chosen as the highest
absolute amplitude of the displacement within three seconds of the P wave arrival.
Having a dataset with known values of Pd, hypocentral distance r and magnitudes, a
simple regression is used to find the coefficients of the linear relation
100

log
10
(Pd) = A + BM + Clog
10
(r). (3)

Wu and Zhao (2006) and Zollo et al. (2006) suggested that the obtained
coefficients A, B and C can then be used with measured values of Pd and estimated r of
new events to estimate the final event magnitude M.
The initial onset of the P wave in a seismogram can be complex with more than
one distinct arrival. Umeda (1990) identified in a group of ten large earthquakes an
initial weak arrival P
1
, followed by a larger and more impulsive arrival P
2
, and
suggested that the time difference between these two scales with M. Ellsworth and
Beroza (1995, 1998) and Beroza and Ellsworth (1996) observed similar results using
more events covering a greater range of magnitudes. An event with such initial phases
can be seen in our data (Figure 3-7). As shown below, however, such early phases can
be identified only in a small subset of the data
101
-2 -1 0 1 2 3 4
-1
-0.5
0
0.5
1
x 10
4
BV epdist: 4.2 km baz: 35.2 mag: 2.31
-0.1 -0.05 0 0.05 0.1 0.15
-600
-400
-200
0
200
P
1
P
2
t
sip
Time (sec)
Time (sec)

Figure 3-7: An example seismogram with an initial low amplitude arrival P
1
(solid black line) followed
by a larger more impulsive arrival P
2
(dashed black line). The top panel shows 3 second seismogram,
while the lower panel shows expanded waveform around the picked phases. The red line is an extension
of the maximum gradient of the first motion used in the method of Iio (1995).
Iio (1995) proposed another way of distinguishing between different early
phases in the context of rupture models associated with slip-weakening friction and
expanding circular crack. The latter is expected to produce a ramp function at the P
wave onset, whereas the former is assumed by Iio (1995) to have a rise of the P velocity
pulse that becomes increasingly more gradual for larger events. Projecting the
maximum slope of the initial motion of the P wave back to its intercept with the base
line level (Figures 3-7 and3- 8) gives the arrival time as it would be if the P wave was a
ramp function. The time difference t
sip
between this intercept and the true first motion is
what we call the Iio signal. It should be noted that Ishihara et al. (1992) found in a study
of large earthquakes that the initial moment rates increase more rapidly for larger
102
earthquakes, which is essentially the opposite of the micro-earthquake observations of
Iio (1995).
t
sip
Figure 3-8: Calculation of t
sip
. An example seismogram with the maximum gradient of the first P wave
pulse projected to the horizontal axis as a ramp function. The resulting time difference t
sip
is used in the
analysis.
3-4. Results
The analysis is performed on the dataset of waveforms generated by 13 clusters
of repeating earthquakes and additional 48 events, both with and without the removal of
the FIR filter artifacts. The instrument response is not removed, limiting our ability to
resolve reliably source effects of events with M < 1 having corner frequencies higher
than ~50 Hz. The following seven measurements are made on each seismogram for
comparison with the final event magnitude. 1) The Ellsworth and Beroza (1995) signal
associated with the time difference between the arrival of P
1
and P
2
phases. 2) The Iio
(1995) signal based on the time between the P arrival and the time it would arrive if it
were a ramp function. 3) A modified Ellsworth-Beroza signal associated with the
absolute area under the seismogram between the P
1
and P
2
arrivals, motivated by the
possibility that the size rather than the duration of the process is important. 4) A
modified Iio signal associated with the absolute area under the seismogram between the
Iio time and the P
2
arrival. 5) The Kanamori (2004) signal involving the predominant
103
period τ
c
from equations (1a) and (1b). (6) The Allen and Olson (2005) signal based on
the predominant period τ
pmax
from equations (2a)-(2d). (7) The Wu and Zhao (2006)
signal using the peak ground displacement parameter Pd. For each of the events at each
of the stations we picked manually at least one (and sometime two) P wave arrival
times and the associated values of τ
pmax
, Pd, τ
c
and t
sip
were calculated. The Ellsworth-
Beroza P
1
arrival can only be identified on 328 seismograms out of the total of 1384. In
the other seismograms there is either no apparent complexity in the P wave onset, or the
pre-signal noise obscures small signals (if such exist) for a given event-station pair.
Therefore less data are contained in the plots that were associated with the Ellsworth-
Beroza type measurements.
Other than the Ellsworth-Beroza signal, the removal of the FIR filter affects
does not change significantly any of the other measures. For the ~330 seismograms for
which the Ellsworth-Beroza measurement was possible, the time between P
1
and P
2
in
analysis without removal of the FIR artifacts is highly variable, fluctuating between a
few hundredths and a few tenths of a second (Figure 3-9a). Within these variations there
is no clear scaling of the measurements with M over the magnitude range of the events.
The sum of the absolute amplitude of the seismogram between the P
1
and P
2
phases
(Figure 3-9c) also seems to have no clear scaling with the final event size. Removal of
the FIR artifacts (Figure 3-10a and 3-10c) reduces the number of points, but has no
other significant effect on the results. On the other hand, the Iio method does show a
tendency for the average value of t
sip
to increase with M (Figure 3-9b), and this feature
is unaffected by the removal of the FIR filter affects (Figure 3-10b). The Wu-Zhao
104
method is only performed on the uncorrected waveforms, as the amplitude measurement
should be unaffected by the filter artifact.
0
0.2
0.4
0.6
-0.02
0
0.02
0.04
0.06
0
5
10
0
100
200
0 1 2 3 4
0
0.2
0.4
0 1 2 3 4
0
0.2
0.4
Modified
Ellsworth-Berzoa
Modified
Iio
log(Area Under P
2
-P
1
) τ
c

Time Difference (sec)
P
1
-P
2
Area Under t
sip
Time Difference (sec)
t
sip
Magnitude
Ellsworth-Berzoa
Kanamori Olsen-Allen
Magnitude
τ
pmax

Iio
398 301 485
0.62
0.94 1.37
1.36 0.60
a) b)
c) d)
e) f)

Figure 3-9: Calculated signals from early in the P waveforms versus magnitudes. (a) The time difference
P
1
- P
2
. (b) The time difference t
sip
. (c) The area under the time between P
1
and P
2
. (d) The area under the
time difference t
sip
. (e) The calculated dominant period τ
c
. (f) The calculated dominant period τ
pmax
. The
small black symbols are measurements at individual stations while the larger black symbols are average
values at all stations for a given event. The circles and triangles correspond to events from the clusters,
while the squares and diamonds correspond to the larger non-cluster events. Some points that are
considerably outside the common range are indicated by arrows on top, at the corresponding magnitude
positions, with the measured values given by the numbers. The solid red lines give linear least-squares
fits to the mean values, while the dashed red lines connect running averages of 4 points above and below
the line fits.
105
0
0.2
0.4
0.6
-0.02
0
0.02
0.04
0.06
0
5
10
0
100
200
0 1 2 3 4
0
0.2
0.4
0 1 2 3 4
0
0.2
0.4
Modified
Ellsworth-Berzoa
Modified
Iio
log(Area Under P
2
-P
1
) τ
c

Time Difference (sec)
P
1
-P
2
Area Under t
sip
Time Difference (sec)
t
sip
Magnitude
Ellsworth-Berzoa Iio
Kanamori Olsen-Allen
Magnitude
τ
pmax

0.74
0.71 0.74 2.33 0.62 1.35
715
443
1036 652
251
a) b)
c) d)
e) f)

Figure 3-10: Similar to Figure 3-9 for data that has been processed to remove the FIR filter artifacts. The
solid red lines give linear least-squares fits to the mean values, while the dashed red lines connect running
averages of 4 points above and below the line fits. The dashed grey lines are the linear least-squares fits
of Figure 3-9.
Figure 3-11a shows the relationship between Pd, hypocentral distance, and M. To
extract the magnitude dependence of the measured values of Pd, equation (3) with the
constants determined through a best fit regression analysis, is used to correct the
measured values of Pd for distance effects by normalizing them to a reference distance.
Zollo et al. (2006) performed a similar procedure correcting the measured Pd values to
a reference distance of 10km, whereas in this study we correct the measurements to the
shortest hypocentral distance in the data (~6km). Figure 3-11b gives the corrected Pd,
after the linear correction for the distance effects is performed. The correction does not
seem to fully remove the affect of distance, as the points with the smallest hypocentral
distances are consistently low and outside the trend of the data. Some of the recorded
106
waveforms are known to be clipped at stations near the event hypocenter. While care
was taken not to use any obviously clipped waveforms in this analysis, the near-source
results with relatively low amplitude could perhaps be affected by clipped waveforms.
Magnitude
Magnitude
Corrected Pd
Log of hypocentral
distance (km)
a)
b)
10
1
10
1
10
0
10
1
10
2
10
3
10
4
10
5

0.5
1
1.5
2
2.5
3
3.5
0.5
1
1.5
2
2.5
3
3.5
Hypocentral distance (km)
Pd
0 0.5 1 1.5 2 2.5 3 3.5 4
10
1
10
2
10
3
10
4
10
5
10
6
10
7

0.9
1
1.1
1.2
1.3
1.4
1.5
1.6
a)
b)

Figure 3-11: Calculated corrected Pd. (a) The distribution of measured Pd against hypocentral
distance with colors representing the magnitudes. The lines give linear regression relations between Pd
and hypocentral distance using the magnitudes values written on the right. (b) Corrected values of Pd
using the regression the coefficients. See text for more explanation. The line gives the relation between
the corrected Pd and magnitude from the regression and the colors represent the log
10
of the hypocentral
distances.
107
Of the two techniques of measuring the predominant period, the Kanamori
(2004) method does not show clear scaling with M and has greatly varying values for
data with similar magnitudes. This is likely because the technique is ill-equipped to deal
with such small events because of the fixed 3 second time window, and is sensitive to
low frequencies which are not recorded reliably in our data. In contrast, the Allen and
Olson (2005) method with the recursive equations (2a)-(2d) appears to show significant
scaling with the final event size. The Wu-Zhao method (Figure 3-11b) shows the
clearest scaling with magnitude of the tested techniques. This may be related to the fact
that the Wu-Zhao method is similar to a measurement of the local magnitude. We note
that with the employed small events, both the Allen-Olson and Wu-Zhao measurements
are made in portions of the waveforms after the end of the rupture. Consequently, the
scaling of the corrected Pd and τ
pmax
with M do not indicate, on their own, any
dependence of the final event size on the initial part of the rupture.

3-5. Discussion
From earthquakes occurring in the 13 clusters we have used 142 events that
resulted in 1141 seismograms at the 10 stations of the temporary network. These were
augmented by 48 non-cluster events that added additional 243 seismograms. Prior to the
removal of the acausal response of the FIR filter, 142 of these events (75%) had one or
more stations for which a P
1
arrival could be picked. After the removal of the FIR filter
artifacts the number of events dropped to 128 (67%). Thus approximately 8% of the
candidates for the Ellsworth-Beroza type nucleation phases can be discounted as an
artifact of the processing and recording system. When the condition is imposed that for
108
a given event two or more stations record a P
1
arrival, the number of qualifying events
is 57 (30%) after the removal of the filter artifacts, compared to 74 (38%) before. When
at least three stations are required to record the P
1
arrival for each event, the number of
events that satisfy this condition is 33 (17%) after the removal of the FIR artifacts. In
summary, for a given event within our dataset there is a good chance that one station
will have a weak initial arrival prior to the main P wave. It is much less common to see
a P
1
phase across multiple stations in the network, something that occurs for less than
20% of the events. If the P
1
arrival reflects an important part of the earthquake source it
should be visible at most stations, barring problems of radiation patterns. It thus appears
that apart from 17% of the events, the identified P
1
arrivals can likely be attributed to
un-correlated noise, weak head wave phases or other near-fault waveform complexities.
Alternatively, if the P
1
phase is commonly generated by the source, the signal to noise
ratios may be too low and propagation distances too great for observing consistently
this phase.
Figure 3-12 shows the seismograms at station FP for all the events in cluster 30,
including the event producing the seismograms of Figures 3-4 and 3-5 that give
recordings at multiple stations. In Figure 3-4, a similar initial low amplitude P
1
arrival
occurs at most of the recording stations, so this event falls into the ~20% which may be
candidates for having signatures of nucleation phases. Given the approximately constant
time difference between P
1
and P
2
at the various stations, this feature could be the result
of source, path or site effects. Since the event used in Figure 3-4 belongs to a cluster of
repeating earthquakes, if the early arrival reflects either path or site effects one would
expect to see the same phase for all the events in the cluster. As this is not observed in
109
Figure 3-12, it seems reasonable to attribute the small initial arrival in Figure 3-4 to a
source effect of only that one event in the cluster. A signature of a nucleation phase
would be expected for all events, while Figure 3-12 clearly shows that for the other
events in the cluster, including one that is larger in magnitude, there is no apparent
complexity or initial phase in the early P waveform. What has been characterized here
as a P
1
arrival has the characteristics of a small foreshock in a similar location to the
following mainshock. Such small foreshocks may occur for 10-20% of the events in the
examined dataset. The amplitude and duration of these assumed foreshocks (Figure 3-
9a) have no apparent bearing on the magnitude of the event that they precede. For
events in the magnitude range of this study, with stations at hypocentral distances of ~4-
20 km, the candidates for nucleation signatures can thus be explained as either
processing artifacts, noise or foreshocks.
-1 -0.5 0 0.5 1 1.5
0
1
2
3
-0.2 -0.1 0 0.1 0.2 0.3 0.4
0.5
1
1.5
2
Time (sec)
Time (sec)
Event magnitude Event magnitude

Figure 3-12: Normalized velocity seismograms of the 5 similar earthquakes of cluster 30 recorded at
station FP. The waveforms are aligned at the P wave arrivals (dashed vertical lines) and plotted at their
corresponding magnitude. The top panel is the waveform from 1 second before to 1.5 seconds after the P
arrival, while the bottom panel is a closer view around the P wave arrival (at time 0). The solid vertical
line marks the single event in the cluster with a P
1
arrival (see Figures 4 and 5 for seismograms from that
event at other stations). None of the other events have clear P
1
arrival.
110
The Iio signals show some positive trend with magnitude in the least-squares fit
lines (Figure 3-9b and 3-10b). This, however, can be explained by a combination of the
increase in scatter with magnitude, and the fact that the t
sip
values have to be positive
such that the increased scatter is positively biased. This may lead to an increase in the
mean of the measured t
sip
values that is not necessarily a genuine feature. Furthermore
many of the measured times are on the order of or smaller than a single time interval
associated with the sampling rate (100 samples per second), and thus may not be
reliable. Instrumentation with higher sampling rates may be able to determine whether
or not a slow initial rise exists in the early portions of the P waveforms.
In our analysis, the Olson-Allen and Wu–Zhao methods of determining the
predominant period for the early portion of the waveform (τ
pmax
) and the peak
displacement amplitude (Pd), respectively, produce the only signals that show clear
scaling with M. These methods, however, provide no information about the rupture
process in our analysis with the small employed event sizes. With events below
magnitude 4-5, the time in the waveform when the value of τ
pmax
and Pd are calculated
(4 and 3 seconds, respectively) is after the end of the rupture. However, Olson and
Allen (2005) found some suggestion that the scaling with magnitude continues to events
larger than 5, and in that case the measure made in the seismogram is determined before
the rupture terminates. At the higher magnitude range of their study there are fewer data
points and the quality of the fit to the trend worsens (Rydelek and Horiuchi 2006).
However, if their trend is robust up to these larger events, it would imply that the final
event size is somewhat determined before the rupture is over. Lewis and Ben-Zion
(2008) found similar supporting evidence using a waveform dataset associated with
111
high sampling rate (up to 1.6 KHz), generated by very-near small events from deep
South Africa mines. In that study, consistent scaling of the final event size remains
apparent even after reducing the time windows used in the Wu-Zhao and Olson-Allen
methods to below the estimated rupture duration of the events. We note that the Olson-
Allen and Wu–Zhao methods are related to measures of the corner frequency and local
magnitude, respectively, and are thus associated with estimates of standard source
properties. The scaling of these signals with M in the studies of Olson and Allen (2005),
Wu and Zhao (2006), Zollo et al. (2006) and Lewis and Ben-Zion (2008) probably
reflect the fact that stronger initial failures are likely to propagate further and produce
larger final event sizes.

112
Chapter 3 References
Abercrombie, R. E. and J. Mori, (1996). Occurrence patterns of foreshocks to large
earthquakes in the western United States. Nature, 381, 303-307.
Allen, R.M. and H. Kanamori, (2003). The potential for earthquake early warning in
southern California. Science, 300, 786-789.
Anderson, J. and Q. Chen, (1995). Beginnings of earthquakes in the Mexican
subduction zone on strong-motion accelerograms. Bull. Seismol. Soc. Am., 85,
1107-1116.
Ben-Zion, Y., (1989). The response of two joined quarter spaces to SH line sources
located at the material discontinuity interface. Geophys. J. Int., 98, 213-222.
Ben-Zion, Y., M. Eneva, and Y. Liu, (2003). Large Earthquake Cycles and Intermittent
Criticality On Heterogeneous Faults Due To Evolving Stress and Seismicity, J.
Geophys. Res., 108, NO. B6, 2307, doi:10.1029/2002JB002121.
Ben-Zion, Y. and P. Malin, (1991). San Andreas Fault zone head waves near Parkfield
California. Science, 251, 1592-1594.
Ben-Zion, Y. and J. R. Rice, (1997). Dynamic simulations of slip on a smooth fault in
an elastic solid. J. Goephy. Res., 102 17,771-17,784.
Ben-Zion, Y., Z. Peng, D. Okaya, L. Seeber, J. G. Armbruster, N. Ozer, A. J. Michael,
S. Baris, and M. Aktar, (2003).A shallow fault zone structure illuminated by
trapped waves in the Karadere-Duzce branch of the North Anatolian Fault,
western Turkey. Geophys. J. Int., 152, 699-717.
Berzoa, G.C. and W.L. Ellsworth, (1996). Properties of the seismic nucleation phase.
Tectonophysics, 261, 209-227.
Dieterich, J. H., (1992). Earthquake nucleation on faults with rate and state dependent
strength, Tectonophysics, 211, 115-134.
Ellsworth, W.L. and G.C. Beroza, (1995). Seismic evidence for an earthquake
nucleation phase. Science, 268, 851-855.
Ellsworth, W.L. and G.C. Beroza, (1998). Observation of the seismic nucleation phase
in the Ridgecrest, California, earthquake sequence. Geophys. Res. Lett., 25, 401-
404
Hillers, G., Y. Ben-Zion and P.M. Mai, (2006). Seismicity on a fault controlled by rate-
and state-dependent friction with spatial variations of the critical slip distance. J.
Geophys. Res., 111, B01403, doi:10.1029/2005JB003859.
113
Hillers, G., P.M. Mai, Y. Ben-Zion and J-P. Ampuero, (2007). Statistical Properties of
Seismicity Along Fault Zones at Different Evolutionary Stages, Geophys. J. Int.,
169, 515–533, doi: 10.1111/j.1365-246X.2006.03275.x.
Iio, Y., (1995). Observations of the slow initial phase generated by micro earthquakes:
Implications for earthquake nucleation and propagation. J. Geophys. Res., 100
15,333-15,349.
Ishihara, Y., Y. Fukao, I. Yamada and H. Aoki, (1992).The rising slope of moment rate
functions: the 1989 earthquakes off east coast of Honshu. Geophys. Res. Lett., 19,
873-876.
Kanamori, H., (2004). The diversity of the physics of earthquakes. Proc. Jpn. Acad.
Ser., B 80 297-315.
Kanamori, H., (2005). Real-Time seismology and earthquake damage mitigation. Annu.
Rev. Earth Planet. Sci., 33 195-214.
Kilb, D. and J. Gomberg, (1999). The initial subevent of the 1994 Northridge,
California, earthquake: Is earthquake size predictable? Journal of Seismology, 3
409-420.
Lapusta, N. and J. R. Rice, (2003). Nucleation and early seismic propagation of small
and large events in a crustal earthquake model. J. Geoph. Res. 108, B4, 2205,
doi:10.1029/2001JB000793.
Lewis, M.A., Y. Ben-Zion, and J. McGuire, (2007). Imaging the deep structure of the
San Andreas Fault south of Hollister with joint analysis of fault-zone head and
direct P arrivals, Geophys. J. Int., 169, 1028-1042.
Lewis, M.A. and Y. Ben-Zion, (2008). Searching for scaling with magnitude of signals
in the early portion of the P waveform recorded in two deep South African mines,
ms. submitted.
Mori, J. and H. Kanamori, (1996). Initial rupture of earthquakes in the 1995 Ridgecrest,
California sequence. Geophys. Res. Lett., 23 2437-2440.
Nakamura, Y., (1998). On the urgent earthquake detection and alarm system
(UrEDAS). Presented at the Ninth world Conf. Earthquake. Eng., Tokyo.
Nakamura Y, and B.E. Tucker, (1988). Japan’s earthquake warning system: should it be
imported to California? Calif. Geol. 41, 33-40.
Olson, E.L. and R.M. Allen, (2005). The deterministic nature of earthquake rupture.
Nature, 438, 212-215.
114
Peng, Z. and Y. Ben-Zion, (2005). Spatio-temporal variations of crustal anisotropy from
similar events in aftershocks of the 1999 M7.4 Izmit and M7.1 Duzce, Turkey,
earthquake sequences. Geophys. J. Int., 160, 1027-1043.
Peng, Z. and Y. Ben-Zion, (2006). Temporal changes of shallow seismic velocity
around the Karadere-Duzce branch of the north Anatolian fault and strong ground
motion. Pure Appl. Geophys., 163, 567-600, DOI 10.1007/s00024-005-0034-6.
Rydelek, P. and S. Horiuchi, (2006). Is earthquake rupture deterministic? (reply),
Nature, 442, doi:10.1038/nayure04964.
Sato, K. and J. Mori, (2006).Scaling relationships of initiations for moderate to large
earthquakes. J. Geophys. Res., 111, B05306, doi:10.1029/2005JB003613.
Scherbaum, F. and M-P. Bouin, (1997). FIR filter effects and nucleation phases.
Geophys J. Int., 130 661-668.
Scherbaum, F., (1996). Of Poles and Zeros: fundamentals of Digital Seismology.
Kluwer, Dordrecht.
Umeda, Y., (1990). High amplitude seismic waves radiated from the bright spot of an
earthquake. Tectonophysics, 175, 81-92.
Wu, Y-M., L. Zhao, (2006). Magnitude estimation using the first three second P-wave
amplitude in earthquake early warning. Geophy. Res. Lett., 33,
doi:10.1029/2006GL026871.
Zollo, A., M. Lancieri and S. Nielsen, (2006). Earthquake magnitude estimation from
peak amplitudes of very early seismic signals on strong motion records. Geophy.
Res. Lett., 33, doi:10.1029/2006GL027795.
115
Chapter 4: Examination of Scaling Between Earthquake Magnitude
and Proposed Early Signals in P Waveforms From Very Near-source
Stations in a South African Gold Mine
Abstract
We analyze near-source seismograms, recorded in a deep South African mine
for scaling between signals early in P waveforms and the magnitude M of the events.
The data consist of recordings at 26 stations of 122 events with -1.5 < M < 2.5 and
hypocentral distances as low as ~100 m. We examine four signals, belonging to two
classes. Signals in the first class include Ellsworth-Beroza and Iio measurements
assumed to reflect signatures of seismic nucleation phases. Signals in the second class
consist of Nakamura type predominant period and Wu-Zhao corrected peak
displacement in the early waveforms, related to measures of the stress drop and local
magnitude. Candidates for the Ellsworth-Beroza signals, involving a weak arrival
before the main P phase, exist for only 20-30% of the waveforms. These phases do not
appear to be associated with a general component of the rupture initiation process. The
Iio signal increases with increasing M but may be attributed, at least partially, to
statistical effects. The predominant period and peak displacement have significant
correlations with M, which remain robust when the analysis time window is reduced
from 3 sec or before the S wave arrival to 0.01 sec. The latter is shorter than the
estimated rupture duration for events with M > ~0, implying that the earthquake size is
affected statistically by the initial rupture strength. If these correlations hold for larger
earthquakes, as suggested by other studies, the early predominant period and peak
displacement are useful for early warning systems.

116
4-1. Introduction
Observational studies on scaling relations between features in the early portions
of seismograms and the ultimate size of the generating events are motivated by two
general goals. The first involves efforts to constrain with seismic data the physics
associated with transitions (e.g., Rice, 1980; Dieterich, 1992; Ben-Zion, 2003) from
aseismic nucleation phases to dynamic rupture. Iio (1995) and Ellsworth and Beroza
(1995) suggested, using methods discussed in the next section, that properties of the
initial pulse in seismograms might reflect signatures of nucleation phases and scale with
the final event size. Some such scaling was observed in numerical simulations of
earthquakes on a fault governed by rate-state friction with heterogeneous critical slip
distance (Hillers et al., 2006, 2007). However, other observational studies did not
support the suggested scaling relations (e.g., Anderson and Chen, 1995; Mori and
Kanamori, 1996; Abercrombie and Mori, 1996; Ishihara et al. 1992; Kilb and Gomberg
1999; Sato and Mori 2006). The second general goal of studies on scaling between early
signals of seismograms and event sizes is associated with efforts to obtain rapid
estimates of magnitudes (e.g., Nakamura and Tucker, 1988; Kanamori, 2005) that can
be used as early-warning signals. Olson and Allen (2005), Wu and Zhao (2006) and
Zollo et al. (2006) presented evidence, discussed in more detail in the next section, that
the predominant period and maximum displacement in the early part of seismograms
scale with the magnitude M of the generating events.
The above observational studies have generally analyzed data generated by
earthquakes within large regions or sets of regions (sometimes using only a few tens of
events or less), leading to possible mixing of source, path, site and instrument effects
117
(compounded in some cases by small statistical samples). To focus on statistically
significant variations associated with source effects, Lewis and Ben-Zion (2007) used
seismograms generated by clusters of 190 repeating and near by earthquakes on the
Karadere-Duzce branch of the North Anatolian fault in Turkey (Peng and Ben-Zion,
2006). The results indicated little or no scaling for the Ellsworth and Beroza (1995) and
Iio (1995) nucleation signals, while moderate to clear scaling was found for the Olson-
Allen and Wu-Zhao type signals albeit within large scatter. This may be understood
intuitively since the Olson and Allen (2005) and Wu and Zhao (2006) signals are
related, respectively, to the stress drop and magnitude of the earthquakes.
Owing to data (sampling-rate) limitations, the analysis of Lewis and Ben-Zion
(2007) used relatively large time windows compared to the rupture times of the
employed events. Therefore, the results from that study can not be used to infer on
properties of the early or other sub-parts of the rupture process. In this work we perform
similar analysis using high sampling rate seismograms recorded within a deep South
African mine. The data were generated by events with hypocentral distances from ~0.1-
10 km, yet without the heterogeneity and damage that exist at these distances in a
regular fault zone setting. The rock volume between the events and recording stations
consists primarily of competent quartzite having little attenuation (Yamada et al., 2007).
Hence our study area is both smaller and more uniform compared to previous studies.
The high sampling rates and other properties of the observational data in this work
allow us to use time windows that are smaller than the estimated rupture durations of
the events.
118
1 1.2
8
-0.01 -0.005 0 0.005 0.01 0.015
0.2
0.25
0.3
0.35
0.4
Hypocentral distance (km)
Time (sec)
-1 -0.5 0 0.5 1
-2
-1.5
-1
-0.5
0
0.5
1
1.5
2
2.5
3
Normalized Time
Magnitude
a)
b)

Figure 4-1: Radial-component normalized velocity seismograms at a station K1 with 0 seconds
corresponding to the P wave arrival for each event. a) The seismograms are plotted at their station-
hypocentral distances from the event and the. seismograms from events within 450 m of the station are
shown at amplified scales. b) Example seismograms for events in the magnitude range between -1.6 and
2.4 at station K1.
119
The analysis employs seismograms generated by 122 events in the magnitude
range ~2.5 > M > ~-1.5 (Figure 4-1), recorded at up to 26 very near-source stations.
Only events that occur during the night are used in an attempt to avoid station triggers
directly related to the mining activity. The dataset is recorded at 26 stations, 2 of which
have fixed sampling rates of 1000 and 2000 samples per second. The sampling rate at
the other stations varies depending on the recorded amplitude in the range 400-10,000
samples per second. The total number of seismograms used in the study is 773, which is
sufficiently large to draw meaningful statistical conclusions.
In the next section we describe briefly the data and the employed analysis
methods. The results are presented in section 3 and discussed in section 4. Candidates
for the Ellsworth-Beroza signals, associated with identification of a weak arrival before
the main P phase, exist only in about a third of the data. The fraction is expected to be
reduced to about 20% in analysis that includes the removal of artifacts associated with
FIR filters (Scherbaum and Bouin, 1997; Lewis and Ben-Zion, 2007). The other signals
can be determined by definition on all waveforms. The results show some positive
scaling with M within considerable scatter of the identified early signals from all
methods. The predominant period in the early waveforms (Olson and Allen, 2005) and
the peak amplitude of displacement (Wu and Zhao, 2006; Zollo et al., 2006) have the
most significant scaling with M, albeit still within large fluctuations. Modifying the
techniques to use a smaller time window after the P wave arrival has only a minimal
impact on the observed scaling, even when the time window is smaller than the
estimated rupture duration for events with M > ~0.

120
4-2. Data and Methods
We use seismograms observed in the TauTona deep South African gold mine
during the quiet hours of no ore production on January 10
th
2006 when a relatively-large
M 2.4 event occurred. The dataset was recorded (Mendecki, 1997) by Tri-axial 4.5Hz
geophones with sampling rates that vary in the range 400-10000 Hz, and provided by
the ISS International Limited and Boettcher et al. (2007). The observed seismograms
are decimated to reduce the sampling rate, during which anti-alias filters may be applied
depending on the degree of decimation and the specific station set up. Figure 4-1 shows
the range of magnitudes and hypocentral distances of events recorded by station K1.
The employed events have hypocentral distances between ~60 m and ~10 km from the
seismometers, with many occurring within a few hundred of meters (Figure 4-1a). The
majority of the earthquakes have magnitudes smaller than 0.5 while the largest event
has a magnitude of 2.4 (Figure 4-1b). The seismograms of Figure 4-1 illustrate the high
quality and short hypocentral distance of the recorded data. Richardson and Jordan
(2002) obtained reliable corner frequencies of up to ~400 Hz from similar recordings,
and the ISS International Limited states that usable frequencies are between 3 and 2000
Hz. In the analysis discussed below we use the radial components of seismograms.
The initial onset of the P wave in a seismogram can be complex with more than
one distinct arrival. Figure 4-2a shows examples of seismograms generated by a M 1.5
event at 5 nearby stations, with the panels in Figure 4-2b showing enlarged views of the
P wave onset. The observations illustrate the existence of a variety of small amplitude
features with different character prior to the main P arrival at 0 seconds.
121
Hypocentral distance (km)
Time since P arrival (sec)
-0.04 -0.02 0 0.02 0.04 0.06 0.08 0.1 0.12
0.3
0.31
0.32
0.33
0.34
0.35
0.36
0.37
L1
J1
K1
M1
N1
-1
0
1
L1
-1
0
1
J1
-1
0
1
K1
-1
0
1
M1
0
1
N1
Time since P arrival (sec)
-1
-0.02 -0.01 0 0.01 0.02 0.03
a) b)

Figure 4-2: Radial-component normalized velocity seismograms at 5 stations for a magnitude 1.5 event
with the P wave arrivals at 0 seconds. a) The seismograms are plotted at their station-hypocentral distance
from the event. The codes of the different stations are marked above each record. b) Enlarged views of
the early portions of the P waves in (a).
Umeda (1990) identified in a group of ten large earthquakes an initial weak arrival P
1
,
followed by a larger and more impulsive arrival P
2
(Figure 4-3a), and argued that the
time difference between these two scales with the magnitude M. Ellsworth and Beroza
(1995, 1998) and Beroza and Ellsworth (1996) presented similar scaling results using
more events covering a greater range of magnitudes and proposed that the early weak
arrivals reflect signatures of different-sized nucleation phases. However, other
122
observational studies (e.g., Anderson and Chen, 1995; Mori and Kanamori, 1996;
Abercrombie and Mori, 1996; Kilb and Gomberg 1999; Sato and Mori 2006) found no
scaling and sometimes opposite scaling between the time intervals separating early P
arrivals and M.
-5
0
5
x 10
5

0.15 0.1

0.05 0 0.05 0.1 0.15 0.2
0
0.02
0.04
0.06
Time (sec)
Velocity
τ
p
(s)
P
1
P
2
τ
d
τ
p
max
t
sip
X1 Dist:7387m Mag:2.2
a)
b)
0
0.2
0.4
0.6
0.8
Absolute
Displacement
Pd
Pdt
c)

Figure 4-3: Obtaining signals (a) A radial velocity seismogram with picked P
1
(solid vertical line) and P
2

(dashed vertical line) arrivals marked. The gray lines denote the projection of the maximum gradient of
the main P wave arrival (P
2
) up to the base line level of the seismogram. The intercept of the projected
line and the horizontal occurs at some time t
sip
after the P arrival. (b) The parameter τ
p
calculated from the
seismogram in (a) using the method of Olson and Allen (2005). The circle denotes the maximum point
within a 4 second time window τ
pmax
, which occurs at the time τ
d
after the P arrival. (c) The seismogram
in (a) converted to its absolute displacement. The circle denotes the maximum displacement Pd within 3
seconds of the P arrival, which occurs at the time Pdt after the P arrival.
123
Scherbaum and Bouin (1997) and Lewis and Ben-Zion (2007) discuss artifacts
from the recording and digitization of seismograms that could be mistaken for
complexity in the P wave onset. The commonly-used FIR filter can generate
(Scherbaum and Bouin, 1997) precursory signals to impulsive arrivals that scale with M
and could be misidentified as nucleation phases. Lewis and Ben-Zion (2007) examined
waveforms with and without the FIR filter removed and found that filter artifacts alone
cannot explain all the complexity in the early P wave onsets. However, removing the
filter artifact resulted in a reduction of the weak P
1
arrivals from ~30% to ~20% of the
recorded waveforms. The filter effects on the other methods discussed below were
found to be minimal. Some of the ISS data acquisition units in the Tautona mine have
FIR filters. Nevertheless, in this study the FIR filters are not removed since the effects
of the filters on most of the employed signals are very small and the precise details of
each recording system are unknown to us.
Iio (1995) suggested another way of distinguishing between different P wave
onsets in the context of rupture models associated with slip-weakening friction and
expanding circular cracks. The latter is expected to produce a ramp function at the P
wave onset, whereas the former is assumed by Iio (1995) to produce a rise of the P
velocity pulse that becomes increasingly gradual for larger events. Projecting the
maximum slope of the initial motion of the P wave back to its intercept with the
baseline level (Figure 4-3a) gives the arrival time as if the P wave were a ramp
function. The time difference t
sip
between this intercept and the true first motion reflects
the gradualness of the onset. However, Ishihara et al. (1992) found in a study of large
124
earthquakes that the initial moment rate increases more rapidly for larger earthquakes,
which is essentially the opposite of the micro-earthquake observations of Iio (1995).
Olson and Allen (2005) calculated the predominant period in the early portion of
the P wave, building on previous early warning studies (e.g., Nakamura, 1988; Allen
and Kanamori, 2003). A value, τ
p
, is iteratively calculated (Figure 4-3b) for each point
in the seismogram as follows:

, 2
i
i
D
X
p
i
π τ = (1a)

where

2
1 i i i
x X X + =
−
α (1b)

and

2
1
i
i i
dt
dx
D D






+ =
−
α . (1c)

In (1b) and (1c), x
i
is the recorded ground velocity at time step i and α is a
smoothing constant, dependent on the sampling rate. The maximum value of τ
p
in a
given time interval (τ
pmax
) was proposed to scale with M. Here we calculate τ
pmax
within
a time window that varies from 6 samples after the main P arrival to 3 seconds after the
125
P arrival or before the S arrival. The time after the P wave arrival where τ
pmax
occurs
(Figure 4-3b) is labeled τ
d
. The value of τ
d
varies but is generally greater for larger
events.
Wu and Zhao (2006) and Zollo et al. (2006) observed relations between peak
amplitude in seismograms and M which, like the Olson and Allen (2005) method, are
proposed for use in earthquake early-warning systems. In this procedure the original
velocity seismograms are integrated to displacement and then high pass filtered at 0.075
Hz to remove any low frequency drift that resulted from the numerical integration. A
peak displacement Pd, is chosen as the highest absolute amplitude of the displacement
within 3 seconds of the P wave arrival or before the S wave arrival. This peak in
displacement amplitude occurs at some time Pdt after the P arrival (Figure 4-3c). Using
a data-set with known values of Pd, hypocentral distance r and M, a simple linear
regression is used to find the coefficients of the relation:

log
10
(Pd) = A + BM + Clog
10
(r). (2)

The obtained constants are used to correct the measured values of Pd at each
station for distance by normalizing them to a reference distance (Figure 4-4). We
generate the corrected peak ground displacement by normalizing Pd to the shortest
hypocentral distance in the dataset. Figure 4-4a shows how the peak displacement
decreases with hypocentral distance but increases with magnitude. Figure 4-4b indicates
that the corrected Pd values plot on a linear trend with M, supporting overall the
proposed scaling relation.
126
10 10
-2 -2
10 10
-1 -1
10 10
0
10 10
1
10 10
-4 -4
10 10
-2 -2
10 10
0
10 10
2
10 10
4
Hypocentral Distance (km) Hypocentral Distance (km)
Pd Pd

1

-0.5 -0.5
0
0.5 0.5
1
1.5 1.5
2
-1.5 -1.5
-1 -1
-0.5 -0.5
0
0.5 0.5
1
1.5 1.5
2
Magnitude
-2 -2 -1 -1 0 1 2 3
10 10
-2 -2
10 10
0
10 10
2
10 10
4
10 10
6
10 10
8
10 10
10 10
10 10
12 12
10 10
14 14
Magnitude Magnitude
Corrected Corrected Pd Pd

-2 -2
-1.5 -1.5
-1 -1
-0.5 -0.5
0
0.5 0.5
log
10
of
hypocentral
distance (km)
a)
b)

Figure 4-4: Calculating corrected P
d
. (a) The distribution of measured Pd against hypocentral distance
with shading representing the magnitudes. The lines give the linear regression relation between Pd,
hypocentral and magnitude obtained from fitting to this data, with the magnitude value for each line
written on the right. (b) Corrected values of Pd after the removal of the linear trend with distance shown
as lines in (a). See text for more explanation. The solid line gives the relation between the corrected Pd
and magnitude from the regression. The shading represents the log
10
of the hypocentral distances.
127
The methods of Iio (1995) and Ellsworth and Beroza (1995), seeking to find
evidence that reflect properties of nucleation phases, involve analysis of adjustable time
intervals that become naturally smaller for smaller event sizes. In contrast, the methods
of Wu and Zhao (2005) and Olson and Allen (2005) use a fixed time window within
which their measurements are made. Since the events in the mine are small (M < 2.5),
the typically used (3-4 sec) time windows are larger than the rupture durations. Olson
and Allen (2005) reported a continuous scaling with M of their predominant period for
events up to M ~ 8. With a time window of 4 seconds τ
pmax
is determined before the
rupture had arrested for earthquakes larger than M ~ 6.5. Based on these results, Olson
and Allen (2005) argued that the final rupture size is at least partially determined by the
rupture initiation process. This, however, has been the subject of an ongoing debate
(Rydelek and Horiuchi, 2006; Olson and Allen 2006; Rydelek et al., 2007; Zollo et al.,
2007), which is difficult to resolve because of large scatter in the data and insufficient
numbers of large events. In the following analyses we attempt to overcome this problem
by using the high sampling rate of the recorded data to reduce the employed time
windows below the rupture duration of the events with M > ~0. The 3 second time
window is reduced in steps to a minimum length of 0.01 seconds after the P arrival.
Figures 4-5 and 4-6 show how the measured peak values at station X1 change as the
time window is reduced for two events of magnitudes 2.4 and 0.9. The values of τ
pmax

and the maximum value of Pd change slightly in these cases. We find, however, that the
changes are overall small compared to the differences between values associated with
events of different magnitude, and that the changes for both small and large events
preserves the scaling.
128
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
0
0.2
0.4
0.6
0.8
1
0.02
0.04
0.06
Time (sec)
Velocity
τ
p
(s)
Absolute
Displacement
a)
b)
c)

Figure 4-5: Dependence of τ
pmax
with time window length. Changes in the measured τ
pmax
and corrected
peak ground displacement as the length of the employed time window is reduced. (a) A velocity
seismogram at station X1 from a magnitude 2.4 event progressively cut to the duration of six employed
time windows after the P arrival. The vertical dashed line represents the estimated rupture duration of this
event.(b) Calculations of τ
p
for the seismogram and estimates of τ
pmax
for the six different time windows
shown in (a). Circles and diamonds mark results associated, respectively, with time windows larger and
smaller than the estimated rupture duration. (c) Corresponding results with progressively smaller time
windows for the displacement and the value of Pd.
129
Time (sec)
Velocity
-0.05 0 0.05 0.1 0.15 0.2 0.25 0.3
0
0.2
0.4
0.6
0
0.005
0.01
0.015
0.02
τ
p
(s)
Absolute
Displacement
a)
b)
c)

Figure 4-6: Similar to Figure 4-5 for a smaller event with a magnitude of 0.9.
The values of τ
d
and Pdt at each time window length are compared to an
estimate of the rupture duration Tr:

,
7 . 0
0
s
R
v
r
T = (3a)
with
0
,
, 0
2 f
v
C r
s p
s p
π
= , (3b)

130
where
0
r is the rupture radius, f
0
the corner frequency, ν
p,s
denotes the wave velocity of
P or S waves and C
p,s
is a corresponding constant. The employed values for the
parameters, taken from Richardson and Jordan (2002), are C
p
= 2.01, C
s
= 1.32, ν
p
=
6100 m/s and ν
s
= 3650 m/s.
Alternatively, the rupture duration can be estimated using the scaling relation
between earthquake magnitudes and rupture length Lr:

, 510 . 1
9 . 1 33 . 0 −
=
M
Lr

(4)

with Lr in km. Dividing half the obtained Lr values by an assumed rupture velocity Vr =
0.9ν
s
gives a second independent estimate of the rupture duration
.
Equation (4) stems
from combing (Ben-Zion, 2008) the scaling relation P
0
= (16/7)∆εR
3
between the scalar
seismic potency P
0
of a circular crack with a radius R = L/2 sustaining a uniform strain
drop ∆ε assumed to be 10
−4
, and the empirical potency-magnitude scaling relation of
Ben-Zion and Zhu (2002) for M < 3.5 events log
10
P
0
= 1.0M – 4.7 with P
0
in km
2
-cm.
In the analysis of the next section, the rupture duration is estimated as the
average value derived from equations (3a) and (3b) for P waves and the value derived
from the scaling relation (4), assuming symmetric bilateral ruptures. If the ruptures are
predominately unilateral, as expected for major earthquakes on large faults that separate
different rock bodies (e.g., Weertman, 1980; Ben-Zion and Andrews, 1998; Ampuero
and Ben-Zion, 2008) and supported by observational analysis of large earthquakes
(Henry and Das 2001; McGuire et al., 2002), the rupture durations would be longer by
about a factor of 2 than the used estimates. The values of τ
pmax
and corrected Pd are
estimated from a range of different time windows, which are both larger and smaller
131
than the estimated rupture durations. This allows us to gauge the effects of the
employed time windows on the scaling of the measured quantities with M. If the scaling
with M persists after the reduction of the time windows below the rupture durations, it
implies that some aspect of the early part of the rupture influences the final event size.
4-3. Results
The techniques outlined in the last section are used to extract six signals from
the waveforms. 1) The Ellsworth and Beroza (1995) signal associated with the time
difference between the arrival of P
1
and P
2
phases (Figure 4-3a). 2) The Iio (1995)
signal t
sip
, which is the time between the observed P arrival and the arrival time if it
were a ramp function (Figure 4-3a). 3) A modified Ellsworth and Beroza (1995) signal
associated with the absolute area under the seismogram between the P
1
and P
2
arrivals,
motivated by the possibility (Lewis and Ben-Zion, 2007) that the size rather than the
duration of the process is important. 4) A corresponding modified Iio (1995) signal
associated with the absolute area under the seismogram within the time t
sip
. 5) The
Olson and Allen (2005) signal (Figure 4-3b) based on the predominant period τ
pmax

from equations (1a)-(1d). 6) The Wu and Zhao (2006) signal using the peak ground
displacement (Pd) corrected for distance (Figures 4-3c and 4-4).
For every event at each of the stations we manually picked at least one (and
sometimes two) P wave arrival times, and calculated the associated values of τ
pmax
, Pd,
and t
sip
. The P
1
arrival, as described in Ellsworth and Beroza (1995), can only be
identified at three or more stations for 38 out of 122 events (31%). As a result, less data
are contained in the scaling plots associated with this arrival (Figures 4-7a and 4-7c). In
the remaining seismograms there is either no apparent complexity in the P wave onset,
132
or the pre-signal noise obscures the small signals (if any exist) for a given event-station
pair.
-
0.5
0.0
0.5
1.0
log
10
hypo
distance (km)
10
-4
10
-3
10
-2
10
-1
10
-6
-5
10
10
-4
10
-3
10
-2
10
-1
10
0
10
1
10
2
10
3
10
4
10
0
10
1
10
2
10
3
10
4
10
5
10
6
-3
10
10
-2
10
-1
10
2
10
4
10
6
10
8
10
10
10
12
Time Difference (sec)
P
1
-P
2
Time Difference (sec)
Iio-P
2
Area Under Iio-P
2
Area Under P
2
-P
1
Predominant Period
τ
p
max
(sec)
Ellsworth-Beroza
Modified
Ellsworth-Beroza
Olson-Allen
Magnitude Magnitude
Iio
-2 -1 0 1 2 3
-2 -1
0 1 2 3
Modified
Iio
Wu/Zollo
Corrected Peak
Ground Displacement
a)
c)
b)
d)
e) f)
b = 0.392
r = 0.949
b = 0.686
r = 0.793
b = 0.488
r = 0.840
b = 1.430
r = 0.954
b = 0.610
r = 0.686
b = 0.255
r = 0.402

Figure 4-7: Plots of the 6 investigated parameters as a function of magnitude. The gray dots represent the
measurements at a given station-event pair and the triangles give the average of these values for each
event. For each method the log
10
of the data is fit with a least-squares line (solid black). The gradient (b)
and correlation coefficient (r) of the fit are written under the title. (a) The time difference between an
initial weak and secondary stronger arrival. (b) The time (t
sip
) between the actual and projected P arrival.
(c) The area under the seismogram between the times in (a). (d) The area under the seismogram between
the times in (b). (e) The parameter τ
pmax
against magnitude. (f) The distance-corrected maximum
displacement amplitude Pd. The shading of dots in this panel indicates the hypocentral distance between
each event and the various stations, with light denoting near and dark denoting distant.
The percentage of events with a P
1
arrival is higher than the ~22% before the removal
of the FIR filter in the study of Lewis and Ben-Zion (2007) where the hypocentral
133
distances were ~4-20 km. Based on the results of Lewis and Ben-Zion (2007), the
removal of the FIR filter artifacts from the data is expected to reduce the number of
identified early P
1
arrival by about 10%. We thus estimate that ~20% of the examined
South Africa mine events are associated with an initial sub-event. The higher value
compared with the study of Lewis and Ben-Zion (2007) may be related to the shorter
propagation distances and smaller seismic attenuation associated with the present study.
a)
-0.05 0 0.05 0.1 0.15 0.2
-0.5
0
0.5
1
1.5
2
2.5
3
3.5
G1
E1
F1
I1
A1
B1
S1
V1
D1
X1
Hypocentral distance (km)
Time since P arrival (sec)
-0.04 -0.02 0 0.02
0.2
0.4
0.6
G1
0.5
0.55
F1
1.2
1.4
1.6
I1
1.8
2
B1
2
2.5
S1
2.8
3
3.2
D1
Time since P arrival (sec)
b)

Figure 4-8: Variety of P wave onset style. Radial-component normalized velocity seismograms at 10
stations for a magnitude 2.4 event with the P wave arrivals at 0 seconds. a) The seismograms are plotted
at their station-hypocentral distance from the event. The codes of the different stations are marked above
each record. b) Enlarged views of the early portions of the P waves in from some stations in (a) with a
variety of features before the P wave arrival at 0 seconds.
134
Figure 4-8 illustrates the different characters of the P wave onsets that are
typically seen for a single event at different stations. Some stations (like F1 in this case)
show a sharp onset, while others (like S1 in this case) have a clear weak first arrival.
We find that any station can have initial P
1
weak arrivals for different events, without
clear relations to the event locations, implying that the early weak arrivals are probably
not generated by site effects. Typically, when P
1
and P
2
arrivals are generated by a
given event at different stations, the time differences between the arrivals do not
increase statistically with the source-station distances. The lack of statistical moveout
between the arrivals suggests that the phases are produced at the source rather than
being propagation effects. However, we note again that initial P
1
arrivals are not a clear
feature of all seismograms nor are they associated with all the events (Figure 4-7a).
All six methods show some tendency for the measured parameter to increase
with M over the range of event sizes used in this study (Figure 4-7). However, the
positive trends in the mean values are associated with large scatter around the mean for
each method. The corrected Pd and τ
pmax
have the best fit, strongest trend and highest
correlation coefficients. We note that the correction to Pd seems to have sometimes
over-corrected for the relation between distance and amplitude, as the points with the
smallest hypocentral distances tend to give lower values for a given magnitude (shading
in Figure 4-7f). This might imply that the relationship between the amplitude, distance
and M is more complex than the assumed linear relation of equation (2). The over-
correction could increase the trend with M, as only larger M events are recorded at
instruments at greater distances. Table 1 summarizes the standard deviation of the
differences between the predicted M values from the best fit lines and the actual
135
magnitudes for each tested method. This should provide an estimate of the size of errors
if the employed signals were used to estimate magnitudes.
To test if M can be determined prior to the cessation of rupture, we reduce the
time window in which the Olson and Allen (2005) and Wu and Zhao (2006)
measurements are obtained. If the event sizes are not controlled statistically by aspects
of the early rupture process, the trend between M and the parameters would disappear
(starting with the largest M events) as the time window is reduced below the rupture
durations. On the other hand, if some property of the early rupture process affects the
final event size, the positive trend would remain apparent when the time window is
reduced to below the event duration. The number of samples within the smallest
examined time window varies between 6 and 55 depending on the sampling rate.
However, the quantities Pd and τ
p
remain exactly the same regardless of the duration of
the time window. Only the time in which a maximum value is selected changes, and as
such the calculations are stable even though the number of data points in the time
window may be small. As illustrated in Figures 4-5 and 4-6, the values of τ
pmax
and the
maximum value of Pd change slightly when the time windows are reduced, but these
appear to be secondary effects.
Figure 4-9 shows the values of τ
pmax
using time windows that end 3, 0.05, 0.02
and 0.01 sec after the P arrival (left), and the differences between the time when the
measurements were made after the P arrival and the rupture durations estimated as
discussed in section 2 (right). As the time window gets smaller we observe only a slight
reduction in the gradient b of the least-squares fit to τ
pmax
, and the quality of the fit to
the data, represented by the Pearson correlation coefficient r, also becomes only
136
marginally smaller. The same holds for the corrected peak ground displacement (Figure
4-10).
Time Difference
τ
d
-Tr (sec)
Predominant Period
τ
p
max
(sec)
Time window before S wave
Time window before S wave
Magnitude Magnitude
Time window 0.05 sec
Time window 0.01 sec
Time window 0.02 sec
Time window 0.01 sec
Time window 0.02 sec
Time window 0.05 sec
Time Difference
τ
d
-Tr (sec)
Time Difference
τ
d
-Tr (sec)
Time Difference
τ
d
-Tr (sec)
Predominant Period
τ
p
max
(sec)
Predominant Period
τ
p
max
(sec)
Predominant Period
τ
p
max
(sec)
-2 -1 0 1 2 3 -2 -1 0123
-3
10
10
-2
10
-1
10
-3
10
-2
10
-1
10
-3
10
-2
10
-1
10
-3
10
-2
10
-1
-0.5
0.0
0.5
-0.05
0.00
0.05
-0.05
0.00
0.05
-0.05
0.00
0.05
b = 0.390
r = 0.939
b = 0.385
r = 0.937
b = 0.373
r = 0.934
b = 0.363
r = 0.931

Figure 4-9: The effects of time window change on τ
pmax
. (left) Scaling plots between the parameter τ
pmax

and M based on analysis with time windows that decrease from top to bottom. The other symbols and
notations are the same as in Figure 4-7e. (right) The time τ
d
when the τ
pmax
value was taken minus the
estimated rupture duration Tr for the time windows used in the left panels. When τ
d
- Tr becomes
negative the measure is made on the seismogram before the end of the rupture.
137
-0.5
0.0
0.5
1.0
10
2
10
4
10
6
10
8
10
10
10
12
-
0.5
0.0
0.5
1.0
10
2
10
4
10
6
10
8
10
10
10
12
-0.05
0.00
0.05
-0.05
0.00
0.05
10
2
10
4
10
6
10
8
10
10
10
12
10
2
10
4
10
6
10
8
10
10
10
12
-2 -1 0 1 2 3
-0.05
0.00
0.05
-2 -1 0 1 2 3
Time Difference
Pdt-Tr (sec)
Time window before S wave
Time window before S wave
Magnitude Magnitude
Time window 0.05 sec
Time window 0.01 sec
Time window 0.02 sec
Time window 0.01 sec
Time window 0.02 sec
Time window 0.05 sec
Time Difference
Pdt-Tr (sec)
Time Difference
Pdt-Tr (sec)
Time Difference
Pdt-Tr (sec)
log
10
hypo
distance (km)
Corrected Peak
Ground Displacement
Corrected Peak
Ground Displacement
Corrected Peak
Ground Displacement
Corrected Peak
Ground Displacement
b = 1.430
r = 0.954
b = 1.461
r = 0.952
b = 1.537
r = 0.943
b = 1.636
r = 0.930

Figure 4-10: The effects of time window change on corrected Pd. Similar to Figure 4-9 for analysis
associated with the corrected peak ground displacement. The symbols and notations are the same as in
Figure 4-7f.
At a time window of 0.01 seconds we begin to lose the ability to use the seismograms
with the lowest sampling rates because ~6 samples or more after the P arrival are
required for stable estimates. As such we cannot reduce the employed time interval
138
further or say when the scaling of τ
pmax
and the corrected peak ground displacement
with M break down, other than it is below 0.01 seconds. Figure 4-11 shows with a 0.01
second time window that many of the measurements of τ
pmax
are made using time
intervals lower than the estimated rupture duration. On average τ
d
is 6 times less than
the estimated duration of the event (Figure 4-11b). A corresponding plot for the
corrected peak ground displacement would give similar results. The observations shown
in Figures 4-9, 4-10 and 4-11 indicate that the reduction of the time window has
minimal effects on the scaling and correlation of the τ
pmax
and corrected Pd parameters
with M, even when the measurements are made with time windows shorter than the
estimated rupture duration of the events with M > ~0.
139
τ
d
-Tr (sec)
-0.05
0.00
0.05
10
-2
10
-1
10
0
10
1
10
2
10
3
Tr/τ
d
Difference: measurement
time and rupture duration
Ratio: rupture duration
measurement time
-2 -1 0 1 2 3
Magnitude
a)
b)

Figure 4-11: Comparisons between estimated rupture durations and time intervals used to calculate τ
pmax

when the time window is 0.01 seconds. The symbols are the same as in Figure 4-7. a) Similar to the
bottom right plot in Figure 4-9, except that the grey triangles represent positive values associated with
measurements made in a time window longer than the estimated rupture duration. b) The estimated
rupture duration (Tr) divided by the time after the P wave when the measurement was made (τ
d
). Black
triangles represent ratios larger than 1 associated with measurements that were made in a time shorter
than the estimated rupture duration.
4-4. Discussion
The transition from aseismic slip to dynamic rupture at a finite nucleation size
provided hopes that the precursory aseismic deformation associated with the nucleation
process may be observed geodetically (e.g., Scholz, 1988; Ohnaka, 2003). For such
efforts to be useful, the size of the nucleation zone should scale with the size of the
impending earthquake (to allow differentiation between the numerous small events and
the target large ones). Some such scaling is found in numerical simulations of
140
earthquakes on a fault governed by heterogeneous rate-state friction (Hillers et al.,
2006, 2007). However, careful analysis of geodetic measurements a few tens of km
from the hypocenters of several moderate and large earthquakes show no precursory
aseismic deformation at the available data resolution (e.g., Johnston and Linde, 2002;
Johnston et al., 2006). There were several efforts to observe signatures of the nucleation
phases on seismograms with various proposed signals (e.g., Umeda, 1990; Iio, 1995;
Ellsworth and Beroza, 1995). However, the results based on such signals have been
inconclusive (e.g., Anderson and Chen, 1995; Kilb and Gomberg 1999; Sato and Mori,
2006).
Nakamura and Tucker (1988), Olson and Allen (2005), Wu and Zhao (2006) and
others suggested scaling relations between early signals in the P waveforms and
magnitudes that can be useful for the development of early warning systems. The
proposed scaling relations are the subject of ongoing debate (Rydelek and Horiuchi,
2006; Olson and Allen 2006; Rydelek et al., 2007; Zollo et al., 2007). To clarify the
existence of various signals in the early P waveforms and their possible scaling with M,
we analyze systematically high-resolution seismic waveform data recorded in a deep
South African mine very close to the generating earthquakes.
Complexity in the P wave onset, like that described in Ellsworth and Beroza
(1995), could only be identified in approximately 31% of the waveforms. The lack of a
P
1
arrival could only be attributed to high pre-signal noise obscuring possible early
weak arrivals for some of the remaining 69% of seismograms. Based on the previous
work of Lewis and Ben-Zion (2007) we expect the percentage of P
1
arrivals to be
reduced by the removal of FIR filter artifacts to ~20%. The cause of the occasional
141
initial complexity in the seismograms could be a path effect relating to some local
structure, uncorrelated noise or a feature of the source. Lewis and Ben-Zion (2007)
suggested that such initial waveform complexity may be associated with foreshocks,
which occur near some events but have little or no bearing on their final size. Precisely
what the P
1
arrival represents remains unresolved; however, it is clear that it is not a
ubiquitous feature of the high-resolution seismic waveforms generated by the examined
mine events. We therefore conclude that the P
1
phase is not an important general early
part of the earthquake rupture process.
The t
sip
time proposed by Iio (1995) appears to show a tendency to increase with
increasing M (Figures 4-7b and 4-7d) for the employed data. However, the many of the
measured values of t
sip
are smaller than the sampling rate of even our high quality
seismograms. This raises doubts on whether the measured times have any significance.
Furthermore, this type of measurement is expected to have increasing scatter with
increasing magnitude because of attenuation, since only larger magnitude events are
recorded at more distant stations. Such purely statistical effect may explain much of the
increase in the mean value with M, as the scatter increases only in the positive direction.
The reduced hypocentral distance and attenuation, compared to the previous study of
Lewis and Ben-Zion (2007) on the North Anatolian fault, does not seem to have
allowed for significantly better detection and scaling of either the Ellsworth and Beroza
(1995) or the Iio (1995) proposed nucleation signatures.
The predominant period τ
pmax
(Olson and Allen 2005; Wurman et al., 2007) and
corrected peak ground displacement (Wu and Zhao 2006 and Zollo et al., 2006) both
show scaling with M, which is consistent with many previous studies. These results may
142
be expected since the predominant period and peak displacements in the early
waveforms are related to measures of the stress drop and magnitude of the events. Other
techniques of measuring a predominant period may be better and faster than the one
used here (e.g. Simons et al., 2006) or more suited to only larger event sizes (e.g.
Kanamori 2005). The combination of multiple techniques could also help to improve
earthquake early warning (Wu et al., 2007). The corrected peak ground displacement in
Figure 4-7f has the least scatter about the mean distribution. However we found that our
linear correction for amplitude with distance leaves a slight remaining relationship
between the corrected peak displacement and distance (shading in Figure 4-7f). If this is
generally the case, it would be problematic for a real time earthquake early warning
system when stations near the event report first.
Figure 4-12 shows the ratio between the estimated rupture duration and the time
windows used to determine the parameters in this work and 3 other studies that use
larger magnitude events. The rupture lengths of events with M < 6 in the other studies
are estimated using a magnitude-length scaling relation similar to equation (4) based on
a circular crack (Ben-Zion, 2008), while for events with M > 6 the empirical magnitude-
length scaling relation of Wells and Coppersmith (1994) is used. These estimates of
rupture length are divided by a representative rupture velocity of 3.0 km/s to obtain
rupture durations. As seen in Figure 4-12, the analysis of the largest events used in this
study is comparable to looking at M ~7-8 events of the other studies when the time
windows are scaled up.

143
10
-2
10
-1
10
0
10
1
Tr/Twin
-2 -1 0123 456 78
Magnitude
This study
Olson-Allen
Wu/Zhao
Zollo et al

Figure 4-12: A complication of estimated rupture durations divided by the time window in which the
measurements were made for the events (triangles) of this study, the 71 events (stars) from Olson and
Allen (2005), the 25 events (circles) from Wu and Zhao (2006), and the 12 events (squares) greater than
M6.5 from Zollo et al. (2006). The time windows used in the observations of the different studies are
0.01, 4, 3 and 2 seconds, respectively. Grey and black symbols denote rupture durations that are smaller
and larger than the employed measurement time windows, respectively.
The use of increasingly small time windows, so that the measurements of τ
pmax

and corrected peak displacement are made before the termination of rupture for events
above M ~ 0-0.5 in our study (Figures 4-9 to 4-11) and events above M ~ 6 in the other
studies summarized in Figure 4-12 indicates that there is essentially no effect on the
scaling of these signals with M. In other words, the initial parts of the P waveforms of
larger events have higher amplitude and lower frequency content than small ones. We
therefore conclude that the scaling of τ
pmax
and corrected peak displacement with M is
already established over a time in the seismogram after the initial P arrival that is less
than the duration of rupture. As the techniques associated with τ
pmax
and corrected peak
displacement measure different attributes of seismograms, and both show clear scaling
with M in a very short time after the P wave arrival, it appears that some information on
144
the final event size is contained in the initial part of the rupture, and carried to the
seismic stations while the events are still rupturing. This may be explained by a
statistical tendency of stronger initial rupture phases, with more slip at the onset, to
propagate farther and produce larger magnitude events (Heaton, 1990).
145
Chapter 4 References
Abercrombie, R. E. and J. Mori, (1996). Occurrence patterns of foreshocks to large
earthquakes in the western United States, Nature, 381, 303-307.
Ampuero, J.-P. and Y. Ben-Zion, (2008). Cracks, pulses and macroscopic asymmetry of
dynamic rupture on a bimaterial interface with velocity-weakening
friction, Geophys. J. Int., 173, 674-692, doi: 10.1111/j.1365-246X.2008.03736.x.
Anderson, J. and Q. Chen, (1995). Beginnings of earthquakes in the Mexican
subduction zone on strong-motion accelerograms, Bull. Seismol. Soc. Am., 85,
1107-1116.
Ben-Zion, Y., (2003). Appendix 2, Key Formulas in Earthquake Seismology, in
International Handbook of Earthquake and Engineering Seismology, eds. W. HK
Lee, H. Kanamori, P. C. Jennings, and C. Kisslinger, Part B, 1857-1875,
Academic Press.
Ben-Zion, Y., (2008). Collective Behavior of Earthquakes and Faults: Continuum-
Discrete Transitions, Progressive Evolutionary Changes and Different Dynamic
Regimes, Rev. Geophysics, in review.
Ben-Zion, Y. and D.J. Andrews, (1998). Properties and implications of dynamic rupture
along a material interface, Bull. seism. Soc. Am., 88, 1085-1094.
Ben-Zion Y. and L. Zhu, (2002). Potency-magnitude Scaling Relations for Southern
California Earthquakes with 1.0 < M
L
< 7.0, Geophys. J. Int., 148, F1-F5.
Beroza, G.C., and W.L. Ellsworth, (1996). Properties of the seismic nucleation phase,
Tectonophysics, 261, 209-227.
Boettcher, M. S., A. McGarr, et al., (2007). A Search For Minimum-Magnitude
Earthquakes In South African Gold Mines, Seism, Res, Lett., 78(2): 293.
Dieterich, J. H., (1992). Earthquake nucleation on faults with rate and state dependent
strength, Tectonophysics, 211, 115-134.
Ellsworth, W.L., and G.C. Beroza, (1995). Seismic evidence for an earthquake
nucleation phase, Science, 268, 851-855.
Ellsworth, W.L. and G.C. Beroza, (1998). Observation of the seismic nucleation phase
in the Ridgecrest, California, earthquake sequence, Geophys. Res. Lett., 25, 401-
404.
Heaton, T. H., (1990). Evidence for and implications of self-healing pulses of slip in
earthquake rupture, Physics of the Earth and Planetary Interiors, 64, 1 - 20.
146
Henry, C. and S. Das, (2001). Aftershock zones of large shallow earthquakes: Fault
dimensions, aftershock area expansion, and scaling relations, Geophys. J. Int.,147,
272-293.
Hillers, G., Y. Ben-Zion and P.M. Mai, (2006). Seismicity on a fault controlled by rate-
and state-dependent friction with spatial variations of the critical slip distance, J.
Geophys. Res., 111, B01403, doi:10.1029/2005JB003859.
Hillers, G., P.M. Mai, Y. Ben-Zion and J-P Ampuero, (2007). Statistical Properties of
Seismicity Along Fault Zones at Different Evolutionary Stages, Geophys. J. Int.,
169, 515–533, doi: 10.1111/j.1365-246X.2006.03275.x.
Iio, Y., (1995). Observations of the slow initial phase generated by micro earthquakes:
Implications for earthquake nucleation and propagation, J. Geophys. Res., 100,
15,333-15,349.
Ishihara, Y., Y. Fukao, I. Yamada and H. Aoki, (1992).The rising slope of moment rate
functions: the 1989 earthquakes off east coast of Honshu, Geophys. Res. Lett., 19,
873-876.
Johnston, M. J. S., and A. T. Linde, (2002). Implications of crustal strain during
conventional, slow and silent earthquakes, in International Handbook of
Earthquake and Engineering Seismology, eds. W. HK Lee, H. Kanamori, P. C.
Jennings, and C. Kisslinger, Part A, 589–605. Academic Press.
Johnston, M. J. S., R. D. Borcherdt, A. T. Linde, and M. T. Gladwin, (2006).
Continuous Borehole Strain and Pore Pressure in the Near Field of the 28
September 2004 M 6.0 Parkfield, California, Earthquake: Implications for
Nucleation, Fault Response, Earthquake Prediction, and Tremor, Bulletin of the
Seismological Society of America, 96, S56–S72, doi: 10.1785/0120050822.
Kanamori, H., (2005). Real-Time seismology and earthquake damage mitigation, Annu.
Rev. Earth Planet. Sci., 33, 195-214.
Kilb, D., and J. Gomberg, (1999). The initial subevent of the 1994 Northridge,
California, earthquake: Is earthquake size predictable?, Journal of Seismology, 3,
409-420.
Lewis, M.A., and Y. Ben-Zion, (2007). Examination of scaling between proposed early
signals in P waveforms and earthquake magnitudes, Geophys. J. Int., 171, 1258–
1268, doi: 10.1111/j.1365-246X.2007.03591.x.
McGuire, J.J., L. Zhao and T.H. Jordan, (2002). Predominance of unilateral rupture for
a global catalog of large earthquakes, Bull. seism. Soc. Am., 92, 3309-3317.
Mendecki, A., (1997). Seismic Monitoring in Mines. London, UK: Chapman and Hall.
147
Mori, J., and H. Kanamori, (1996). Initial rupture of earthquakes in the 1995 Ridgecrest,
California sequence, Geophys. Res. Lett., 23, 2437-2440.
Nakamura Y., and B.E. Tucker, (1988). Japan’s earthquake warning system: should it
be imported to California?, Calif. Geol. 41, 33-40
Ohnaka, M., (2003). A constitutive scaling law and a unified comprehension for
frictional slip failure, shear fracture of intact rock, and earthquake rupture, J.
Geophys. Res., 108(B2), 2080, doi:10.1029/2000JB000123.
Olson, E.L., and R.M. Allen, (2005). The deterministic nature of earthquake rupture.
Nature, 438, 212-215.
Olson, E.L., and R.M. Allen, (2006). Replying to P. Rydelek & S. Horiuchi Nature 442,
doi:10.1038/nayure04963. Nature, 442, doi:10.1038/nayure04964.
Peng, Z., and Y. Ben-Zion, (2006). Temporal changes of shallow seismic velocity
around the Karadere-Duzce branch of the north Anatolian fault and strong ground
motion, Pure Appl. Geophys., 163, 567-600.
Rice, J. R., (1980). The mechanics of Earthquake Rupture, in Physics of the Earth’s
Interior, eds. A. M. Dziewonski and E. Boschi, pp. 555-649, Italian Physical
Society / North Holland, Amsterdam.
Richardson, E., and T.H. Jordan, (2002). Seismicity in deep gold mines of South Africa:
Implications for tectonic earthquakes, Bull. Seis. Soc. Am., 92, 1766-1782.
Rydelek, P., and S. Horiuchi, (2006). Is earthquake rupture deterministic? (reply),
Nature, 442, doi:10.1038/nayure04963.
Rydelek, P., C. Wu and S. Horiuchi, (2007). Comments on “Earthquake magnitude
estimation from peak amplitude of very early seismic signals on strong motion
records”, Geophys. Res. Lett., 34, doi: 10.1029/2007GL029387.
Sato, K., and J. Mori, (2006). Scaling relationships of initiations for moderate to large
earthquakes, J. Geophys. Res., 111, B05306, doi:10.1029/2005JB003613.
Scherbaum, F., and M-P. Bouin, (1997). FIR filter effects and nucleation phases,
Geophys J. Int., 130, 661-668.
Scholz, C. H., (1988). The critical slip distance for seismic faulting, Nature, 336, 761-
63.
Simons, F. J., B.D.E. Dando and R. M. Allen, (2006). Automatic detection and rapid
determination of earthquake magnitude by wavelet multiscale analysis of the
primary arrival, Earth and Planetary Science Lett., 250, 214-223.
148
Umeda, Y., (1990). High amplitude seismic waves radiated from the bright spot of an
earthquake, Tectonophysics, 175, 81-92.
Wells, D.L. and K. J. Coppersmith, (1994). New Empirical Relationships among
Magnitude, Rupture Length, Rupture Width, Rupture Area, and Surface
Displacement, Bull. Seis. Soc. Am., 84, No. 4, 974-1002.
Weertman, J., (1980). Unstable slippage across a fault that separates elastic media of
different elastic constants, J. Geophys. Res., 85, 1455-1461.
Wu, Y-M., and L. Zhao, (2006). Magnitude estimation using the first three second P-
wave amplitude in earthquake early warning, Geophy. Res. Lett., 33,
doi:10.1029/2006GL026871.
Wu, Y-M., H. Kanamori, R.M. Allen and E. Hauksson, (2007). Determination of
earthquake early warning parameters, τ
c
and P
d
, for Southern California, Geophys.
J. Int., doi: 10.1111/j.1365-246X.2007.03430.x.
Wurman, G., R.M. Allen and P. Lombard, (2007). Toward earthquake early warning in
Northern California, J. Geoph. Res., 112, doi:10.1029/2006JB004830.
Yamada, T., J. J. Mori, S. Ide, R. E. Abercrombie, H. Kawakata, M. Nakatani, Y. Iio,
and H. Ogasawara, (2007). Stress drops and radiated seismic energies of
microearthquakes in a South African gold mine, J. Geophys. Res., 112, B03305,
doi:10.1029/2006JB004553.
Zollo, A., M. Lancieri, and S. Nielsen, (2006). Earthquake magnitude estimation from
peak amplitudes of very early seismic signals on strong motion records, Geophy.
Res. Lett., 33, doi:10.1029/2006GL027795.
Zollo, A., Lancieri, and S. Nielsen, (2007). Reply to comments by P.Rydelek et al. on
“Earthquake magnitude estimation from peak amplitude of very early seismic
signals on strong motion records”, Geophys. Res. Lett., 34, doi:
10.1029/2007GL030560.
149
Recapitulation
This thesis investigates the structure of fault zones and properties of the early
part of ruptures in large near-source seismic datasets. The velocity structure of major
faults is investigated using trapped waves (Chapter 1) and head waves (Chapter 2). How
early in a seismogram the final event size can be determined is explored using a number
of different previously proposed signals (Chapter 3 and 4).
The results demonstrate the ability of specific fault zone phases to provide
details on small scale features and variations in fault zone structure. Low velocity
damage zones were imaged around the San Jacinto and San Andreas faults using
trapped and head waves respectively. Both imaging techniques were able to further
resolve an apparent asymmetry in the location of the low velocity zone with respect to
the surface trace. This sense of asymmetry is consistent with that predicted by a
preferred rupture propagation direction based on the local velocity contrasts (Ben-Zion
and Shi, 2005).
This thesis also investigates which signals are robust for making rapid estimates
of earthquake magnitude and the effects on the different signals of FIR filter artifacts
introduced by the recording system. Detailed analysis making use of increasingly small
time windows reveals that two of the signals are capable of determining the final event
magnitude in a time smaller than the estimated rupture duration. This suggests that
certain aspects of the early rupture behavior may have some control on the final event
size.
150
References
Abercrombie, R. E., and J. Mori, (1996). Occurrence patterns of foreshocks to large
earthquakes in the western United States, Nature, 381, 303-307.
Aki, K., and W.H.K. Lee, (1976). Determination of three-dimensional velocity
anomalies under a seismic array using first P arrival times from local earthquakes,
Part 1, A homogeneous initial model, J. Geophys. Res., 81, 4381-4399.
Allen, R.M. and H. Kanamori, (2003). The potential for earthquake early warning in
southern California. Science, 300, 786-789.
Ampuero, J-P. and Y. Ben-Zion, (2008). Cracks, pulses and macroscopic asymmetry of
dynamic rupture on a bimaterial interface with velocity-weakening
friction, Geophys. J. Int., 173, 674-692, doi: 10.1111/j.1365-246X.2008.03736.x.
Anderson, J., and Q. Chen, (1995). Beginnings of earthquakes in the Mexican
subduction zone on strong-motion accelerograms, Bull. Seismol. Soc. Am., 85,
1107-1116.
Andrews, D.J. and Y. Ben-Zion, (1997). Wrinkle-like slip pulse on a fault between
different materials, J. Geophys. Res., 102, 55-571.
Bakun, W.H., R.M. Stewart, C.G. Bufe, and S.M. Marks, (1980). Implications of
seismicity for failure of a section of the San Andreas fault, Bull. Seim, Soc. Am.,
70, 185-201.
Bedrosain, P.A., M.J. Unsworth, G.D. Egbert and C.H. Thurber, (2004). Geophysical
images of the creeping segment of the San Andreas Fault: implications for the role
of crustal fluids in the earthquake process, Tectonophysics, 385, 137-158.
Ben-Zion, Y., (1989). The response of two joined quarter spaces to SH line sources
located at the material discontinuity interface. Geophys. J. Int., 98, 213-222.
Ben-Zion, Y., (1990). The response of two half spaces to point dislocations at the
material interface, Geophys. J. Int., 101, 507-528.
Ben-Zion, Y., (1998). Properties of seismic fault zone waves and their utility for
imaging low velocity structure, J. Geophys. Res., 103, 12,567-12,585.
Ben-Zion, Y., (2001). Dynamic ruptures in recent models of earthquake faults, J. Mech.
Phys. Solids, 49, 2209-2244.
Ben-Zion, Y., (2003). Appendix 2, Key Formulas in Earthquake Seismology, in
International Handbook of Earthquake and Engineering Seismology, eds. W. HK
Lee, H. Kanamori, P. C. Jennings, and C. Kisslinger, Part B, 1857-1875,
Academic Press.
151
Ben-Zion, Y., (2008). Collective Behavior of Earthquakes and Faults: Continuum-
Discrete Transitions, Progressive Evolutionary Changes and Different Dynamic
Regimes, Rev. Geophysics, in review.
Ben-Zion, Y. and K. Aki, (1990). Seismic radiation from an SH line source in a laterally
heterogeneous planer fault zone, Bull. Seism. Soc. Am., 80, 971-994.
Ben-Zion, Y. and P. Malin, (1991). San Andreas Fault zone head waves near Parkfield
California. Science, 251, 1592-1594.
Ben-Zion, Y. and J. R. Rice, (1997). Dynamic simulations of slip on a smooth fault in
an elastic solid. J. Goephys. Res., 102 17,771-17,784.
Ben-Zion, Y. and D.J. Andrews, (1998). Properties and implications of dynamic rupture
along a material interface, Bull. seism. Soc. Am., 88, 1085-1094.
Ben-Zion Y. and L. Zhu, (2002). Potency-magnitude Scaling Relations for Southern
California Earthquakes with 1.0 < M
L
< 7.0, Geophys. J. Int., 148, F1-F5.
Ben-Zion, Y. And C.G. Sammis, (2003). Characterization of fault zones, Pure Appl.
Geophys., 160, 677-715.
Ben-Zion, Y. and Z. Shi, (2005). Dynamic rupture on a material interface with
spontaneous generation of plastic strain in the bulk, Earth Planet. Sci. Lett., 236,
486-496, DOI: 10.1016/j.epsl.2005.03.025.
Ben-Zion, Y., S. Katz and P. Leary, (1992). Joint inversion of fault zone head and direct
P arrivals for crustal structure near major faults, J. Geophys. Res., 97, 1943-1951.
Ben-Zion, Y., M. Eneva, and Y. Liu, (2003). Large Earthquake Cycles and Intermittent
Criticality On Heterogeneous Faults Due To Evolving Stress and Seismicity, J.
Geophys. Res., 108, NO. B6, 2307, doi:10.1029/2002JB002121.
Ben-Zion, Y., Z. Peng, D. Okaya, L. Seeber, J. G. Armbruster, N. Ozer, A. J. Michael,
S. Baris, and M. Aktar, (2003). A shallow fault zone structure illuminated by
trapped waves in the Karadere-Duzce branch of the North Anatolian Fault,
western Turkey. Geophys. J. Int., 152, 699-717.
Berger, J., L.N. Baker, J.N. Brune, J.B. Fletcher, T.M. Hanks and F.L. Vernon, (1984).
The Anza array: a high dynamic-range, broad-band, digitally radio-telemetered
seismic array, Bull. Seism. Soc. Am., 74, 1469-1682.
Beroza, G.C., and W.L. Ellsworth, (1996). Properties of the seismic nucleation phase,
Tectonophysics, 261, 209-227.
Boettcher, M. S., A. McGarr, et al., (2007). A Search For Minimum-Magnitude
Earthquakes In South African Gold Mines, Seism, Res, Lett., 78(2): 293.
152
Boness, N.L., and M.D. Zoback, (2004). Stress-induced seismic velocity anisotropy and
physical properties in the SAFOD Pilot Hole in Parkfield, CA, Geophys. Res.
Lett., 31, L15S17, doi:10.1029/2004GL019020.
Crase, E., A. Pica, M. Noble, J. McDonald and A. Tarantola, (1990). Robust elastic
non-linear inversions: applications to real data, Geophysics, 55, 527-538.
Cochran, E.S., J.E. Vidale and Y-G. Li, (2003). Near-fault anisotropy following the
Hector Mine earthquake, J. Geophys. Res., 108, no. B9, 2436,
doi:10.1029/2002JB002352.
Cochran, E.S., Y-G. Li, and J.E. Vidale, (2006). Anisotropy in the Shallow Crust
Observed around the San Andreas Fault Before and After the 2004 M 6.0
Parkfield Earthquake, Bull. Seism. Soc. Am., 96, S364–S375, doi:
10.1785/0120050804.
Cormier, V. F. and P. Spudich, (1984). Amplification of ground motion and waveform
complexities in fault zones: examples from the San Andreas and the Calaveras
faults, Geophys. J. R. Astr. Soc., 79, 135-152.
Davis, P.M., J.M. Rubinstein, K.H. Liu, S.S. Gao, and L. Knopoff, (2000). Northridge
earthquake damage caused by geologic focusing of seismic waves, Science, 289,
1746-1750.
Dibblee, T.W., (1974). Geologic map of the Salinas quadrangle, California, in Open-file
report 74-1021, U.S. Geological Survey.
Dieterich, J. H. (1992). Earthquake nucleation on faults with rate and state dependent
strength, Tectonophysics, 211, 115-134.
Dor O., T.K. Rockwell, and Y. Ben-Zion, (2006a). Geologic observations of damage
asymmetry in the structure of the San Jacinto, San Andreas and Punchbowl faults
in southern California: A possible indicator for preferred rupture propagation
direction, Pure Appl. Geophys., 163, 301-349, DOI 10.1007/s00024-005-0023-9.
Dor O., Y. Ben-Zion, T. K. Rockwell and J. Brune, (2006b). Pulverized Rocks in the
Mojave section of the San Andreas Fault Zone, Earth Planet. Sci. Lett., 245, 642-
654.
Eberhart-Phillips, D. and A.J. Michael, (1993). Three-Dimensional velocity structure,
seismicity and fault structure in the Parkfield region, central California, J.
Geophys. Res., 98, 15737-15758.
Ellsworth, W., (1975). Bear valley, California, earthquake sequence of February-March
1972, Bull. Seism, Soc. Am. 65, 483-506.
153
Ellsworth, W.L., and G.C. Beroza, (1995). Seismic evidence for an earthquake
nucleation phase, Science, 268, 851-855.
Ellsworth, W.L. and G.C. Beroza, (1998). Observation of the seismic nucleation phase
in the Ridgecrest, California, earthquake sequence, Geophys. Res. Lett., 25, 401-
404.
Fletcher, J., L. Haar, T. Hanks, L. Baker, F. Vernon, J. Berger and J. Brune, (1987). The
digital array at Anza, California: Processing and initial interpretation of source
parameters, J. Geophys. Res., 92, 369-382.
Field, E.H., P.A. Johnson, I.A. Beresnev and Y. Zeng, (1997). Nonlinear ground-motion
amplification by sediments during the 1994 Northridge earthquake, Nature, 390,
599–602.
Fohrmann, M., G. Jahnke, H. Igel and Y. Ben-Zion, (2004). Guided waves from
sources out side faults: an indication of shallow trapping structure?, Pure appl,
Geophys., 161, 2125-2137.
Haberland, C. et al., (2003). Modeling of seismic guided waves at the Dead Sea
Transform, J. Geophys. Res., 108, no. B7, 2342, doi:10.1029/2002JB002309.
Hauksson, E., (2000). Crustal structure and seismicity distribution adjacent to the
Pacific and North America plate boundary in southern California, J. Geophys.
Res., 105, 13,875-13,903.
Heaton, T. H., (1990). Evidence for and implications of self-healing pulses of slip in
earthquake rupture, Physics of the Earth and Planetary Interiors, 64, 1 - 20.
Henry, C. and S. Das, (2001). Aftershock zones of large shallow earthquakes: Fault
dimensions, aftershock area expansion, and scaling relations, Geophys. J. Int.,
147, 272-293.
Hillers, G., Y. Ben-Zion and P.M. Mai, (2006). Seismicity on a fault controlled by rate-
and state-dependent friction with spatial variations of the critical slip distance, J.
Geophys. Res., 111, B01403, doi:10.1029/2005JB003859.
Hillers, G., P.M. Mai, Y. Ben-Zion and J-P Ampuero, (2007). Statistical Properties of
Seismicity Along Fault Zones at Different Evolutionary Stages, Geophys. J. Int.,
169, 515–533, doi: 10.1111/j.1365-246X.2006.03275.x.
Igel, H., G. Jahnke and Y. Ben-Zion, (1997). Simulation of SH and P-SV wave
propagation in fault zones, Geophys. J. int., 128, 533-546.
Igel, H., G. Jahnke and Y. Ben-Zion, (2002). Numerical simulation of fault zone guided
waves: accuracy and 3-D effects, Pure appl. Geophys., 159, 2067-2083.
154
Iio, Y., (1995). Observations of the slow initial phase generated by micro earthquakes:
Implications for earthquake nucleation and propagation, J. Geophys. Res., 100,
15,333-15,349.
Ishihara, Y., Y. Fukao, I. Yamada and H. Aoki, (1992). The rising slope of moment rate
functions: the 1989 earthquakes off east coast of Honshu, Geophys. Res. Lett., 19,
873-876.
Jahnke, G., H. Igel and Y. Ben-Zion, (2002). Three-dimensional calculations of fault
zone guided waves in various irregular structures, Geophys. J. Int., 151, 416-426.
Johnston, M.J.S. and A.T. Linde, (2002). Implications of crustal strain during
conventional, slow and silent earthquakes, in International Handbook of
Earthquake and Engineering Seismology, Vol. 81A, Academic Press, New York,
589–605.
Johnston, M. J. S., R. D. Borcherdt, A. T. Linde, and M. T. Gladwin (2006), Continuous
Borehole Strain and Pore Pressure in the Near Field of the 28 September 2004 M
6.0 Parkfield, California, Earthquake: Implications for Nucleation, Fault
Response, Earthquake Prediction, and Tremor, Bulletin of the Seismological
Society of America, 96, No. 4B, pp. S56–S72, doi: 10.1785/0120050822.
Kanamori, H. (2005), Real-Time seismology and earthquake damage mitigation, Annu.
Rev. Earth Planet. Sci., 33, 195-214.
Kanamori, H., (2004). The diversity of the physics of earthquakes. Proc. Jpn. Acad.
Ser., B 80 297-315.
Kilb, D., and Hardebeck, J.L., (2006). Fault Parameter Constraints Using Relocated
Earthquakes: A Validation of First Motion Focal Mechanism Data, Bull. Seism.
Soc. Am. 96, 1140-1158.
Kilb, D., and J. Gomberg, (1999). The initial subevent of the 1994 Northridge,
California, earthquake: Is earthquake size predictable?, Journal of Seismology, 3,
409-420.
Korneev, V.A., R.M. Nadeau and T.V. McEvilly, (2003). Seismological studies at
Parkfield IX: Fault zone imaging using guided wave Attenuation, Bull. Seism.
Soc. Am., 93, 1415-1426.
Lapusta, N. and J. R. Rice, (2003). Nucleation and early seismic propagation of small
and large events in a crustal earthquake model. J. Geoph. Res. 108, B4, 2205,
doi:10.1029/2001JB000793.
Lewis, M.A., Z. Peng, Y. Ben-Zion and F. Veron, (2005). Shallow seismic trapping
structure in the San Jacinto Fault Zone, California, Geophys. J. Int., 162, 867–881,
doi:10.1111/j.1365-246X.2005.02684.x.
155
Lewis, M.A., Y. Ben-Zion and J. McGuire, (2007). Imaging the deep structure of the
San Andreas Fault south of Hollister with joint analysis of fault-zone head and
direct P arrivals, Geophys. J. Int., 169, 1028-1042.
Lewis, M.A., and Y. Ben-Zion, (2007). Examination of scaling between proposed early
signals in P waveforms and earthquake magnitudes, Geophys. J. Int., 171, 1258–
1268, doi: 10.1111/j.1365-246X.2007.03591.x.
Lewis, M.A. and Y. Ben-Zion, (2008). Searching for scaling with magnitude of signals
in the early portion of the P waveform recorded in two deep South African mines,
ms. submitted.
Li, Y-G. and P.C. Leary, (1990). Fault zone trapped waves, Bull. Seis. Soc. Am., 80, 5
1245-1271.
Li, Y-G. and F.L. Vernon, (2001). Characterization of the San Jacinto fault zone near
Anza, California, by fault zone trapped waves, J. Geophys. Res., 106, 30,671-
30,688.
Li, Y-G., J.E. Vidale and K. Aki, (1997b). San Jacinto fault-zone guided waves: A
discrimination for recently active fault strands near Anza, California, J. Geophys.
Res., 102, 11,689-11,701.
Li, Y.G., P. Leary, K. Aki and P.E. Malin, (1990). Seismic trapped modes in the
Oroville and San Andreas fault zones, Science, 249, 763-766.
Li, Y.G., K. Aki, D. Adams, A. Hasemi and. W.H.K. Lee, (1994). Seismic guided
waves trapped in the fault zone of the Landers, California, earthquake of 1992,
J.Geophys. Res., 99, 11,705-11,722.
Li, Y.G., W.L. Ellsworth, C.H. Thurber, P. Malin and K. Aki, (1997a). Fault-zone
guided waves from explosions in the San Andreas fault at Parkfield and Cienega
Valley, California, Bull. Seismol. Soc. Am., 87, 210-221.
Liu, Y., T.L. Teng and Y. Ben-Zion, (2004). Systematic analysis of shear wave splitting
in the aftershock region of the 1999 Chi-Chi earthquake: Evidence for shallow
anisotropic structure and lack of systematic temporal variations, Bull. Seism. Soc.
Am., 94, 2330-2347.
Liu, Y., T.L. Teng and Y. Ben-Zion, (2005). Near-surface seismic anisotropy,
attenuation and dispersion in the aftershock region of the 1999 Chi-Chi
earthquake, Geophys. J. Int., 160(2), 695-706.
McGuire, J., and Y. Ben-Zion, (2005). High-resolution imaging of the Bear Valley
section of the San Andreas Fault at seismogenic depths with fault-zone head
waves and relocated seismicity, Geophys. J. Int., 163, 152-164, doi:
10.1111/j.1365-246X.2005.02703.x.
156
McGuire, J.J., L. Zhao and T.H. Jordan, (2002), Predominance of unilateral rupture for
a global catalog of large earthquakes, Bull. seism. Soc. Am., 92, 3309-3317.
Magistrale, H. and C. Sanders, (1995). P wave image of the Peninsular Ranges
batholith, southern California, Geophys. Res. Lett., 22, 2549-2552.
Marra, F., R. Azzara, F. Bellucci, A. Caserta, G. Cultrera, G. Mele, G., B. Palombo, A.
Rovelli and E. Boschi, (2000). Large amplification of ground motion at rock sites
within a fault zone in Nocera Umbra (Central Italy), J. of Seis., 4 534-554.
Mendecki, A. (1997), Seismic Monitoring in Mines. London, UK: Chapman and Hall.
Michael, A.J. and Y. Ben-Zion, (1998). Inverting fault zone head trapped waves with
genetic algorithm, EOS, Trans. Am. Geophys. Un., 81, F1145.
Michael A.J., S.L. Ross and H.D. Stenner, (2002). Displaced rocks, strong motion, and
the mechanics of shallow faulting associated with the 1999 Hector Mine,
California, earthquake, Bull. Seis. Soc. Am., 92 (4): 1561-1569.
Mizuno, T., K. Nishigami, H. Ito and Y. Kuwahara, (2004). Deep structure of the
Mozumi-Sukenbu fault, central Japan, estimated from the subsurface array
observation of the fault zone trapped waves, Geophys. J. Int., 159, 622-642.
Mori, J., (1993). Fault plane determinations for three small earthquakes along the San
Jacinto fault, California: search for cross faults, J. Geophys. Res., 98, 17, 711-17,
722.
Mori, J., and H. Kanamori, (1996). Initial rupture of earthquakes in the 1995 Ridgecrest,
California sequence, Geophys. Res. Lett., 23, 2437-2440.
Nakamura, Y., (1998). On the urgent earthquake detection and alarm system
(UrEDAS). Presented at the Ninth world Conf. Earthquake. Eng., Tokyo.
Nakamura Y. and B.E. Tucker, (1988). Japan’s earthquake warning system: should it be
imported to California?, Calif. Geol. 41, 33-40
Ohnaka, M. (2003). A constitutive scaling law and a unified comprehension for
frictional slip failure, shear fracture of intact rock, and earthquake rupture, J.
Geophys. Res., 108(B2), 2080, doi:10.1029/2000JB000123.
Olson, E.L. and R.M. Allen, (2005). The deterministic nature of earthquake rupture.
Nature, 438, 212-215.
Olson, E.L. and R.M. Allen, (2006). Replying to P. Rydelek & S. Horiuchi Nature 442,
doi:10.1038/nayure04963. Nature, 442, doi:10.1038/nayure04964.
157
Peng, Z. And Y. Ben-Zion, (2004). Systematic analysis of crustal anisotropy along the
Karadere-Düzce branch of the north Anatolian fault, Geophys. J. Int., 159, 253–
274.
Peng, Z. and Y. Ben-Zion, (2005). Spatio-temporal variations of crustal anisotropy from
similar events in aftershocks of the 1999 M7.4 Izmit and M7.1 Duzce, Turkey,
earthquake sequences. Geophys. J. Int., 160, 1027-1043.
Peng, Z. and Y. Ben-Zion, (2006). Temporal changes of shallow seismic velocity
around the Karadere-Duzce branch of the north Anatolian fault and strong ground
motion. Pure Appl. Geophys., 163, 567-600, DOI 10.1007/s00024-005-0034-6.
Peng, Z., Y. Ben-Zion, L. Zhu and A.J. Michael, (2003). Inference of a shallow fault
zone layer in the rupture zone of the 1992 Landers, California earthquake fro
locations of events generating trapped waves and travel time analysis, Geophys. J.
Int., 155, 1021-1041.
Rice, J. R., (1980). The mechanics of Earthquake Rupture, in Physics of the Earth’s
Interior, eds. A. M. Dziewonski and E. Boschi, pp. 555-649, Italian Physical
Society / North Holland, Amsterdam.
Richardson, E., and T.H. Jordan, (2002). Seismicity in deep gold mines of South Africa:
Implications for tectonic earthquakes, Bull. Seis. Soc. Am., 92, 1766-1782.
Rovelli, A., A. Caserta, F. Marra and V. Ruggiero, (2002). Can seismic waves be
trapped inside an inactive fault zone? The case study of Nocera Umbra, central
Italy, Bull. Seis. Soc. Am., 92, 2217-2232.
Rubin, A.M., and D. Gillard, (2000). Aftershock asymmetry/rupture directivity among
central San Andreas fault micro earthquakes, J. Geophys. Res., 105, 19095-19110.
Rubin, A. M. and J-P. Ampuero, (2007). Aftershock asymmetry on a bimaterial
interface, J. Geophys. Res., 112, B05307, doi:10.1029/2006JB004337.
Rydelek, P. and S. Horiuchi, (2006). Is earthquake rupture deterministic? (reply),
Nature, 442, doi:10.1038/nayure04963.
Rydelek, P., C. Wu and S. Horiuchi, (2007). Comments on “Earthquake magnitude
estimation from peak amplitude of very early seismic signals on strong motion
records”, Geophys. Res. Lett., 34, doi: 10.1029/2007GL029387.
Sanders, C.O. and H. Kanamori, (1984). A seismo-tectonic analysis of the Anza seismic
gap, San Jacinto fault zone, southern California, J. Geophys. Res., 89, 5,873–
5,890.
Sanders C.O. and H. Magistrale, (1997). Segmentation of the northern San Jacinto fault
zone, southern California, J. Geophys. Res., 102(B12), 27453-27467.
158
Sato, K., and J. Mori, (2006). Scaling relationships of initiations for moderate to large
earthquakes, J. Geophys. Res., 111, B05306, doi:10.1029/2005JB003613.
Scherbaum, F., and M-P. Bouin, (1997). FIR filter effects and nucleation phases,
Geophys J. Int., 130, 661-668.
Scherbaum, F., (1996). Of Poles and Zeros: fundamentals of Digital Seismology.
Kluwer, Dordrecht.
Scholz, C. H., (1988). The critical slip distance for seismic faulting, Nature, 336, 761-
63.
Scott, J. S., T.G. Masters and F.L. Vernon, (1994). 3-D velocity structure of the San
Jacinto fault zone near Anza, California -- I. P waves, Geophys. J. Int., 119, 611-
626.
Seeber, L., J.G. Armbruster, N. Ozer, M. Aktar, S. Baris, D. Okaya, Y. Ben-Zionand E.
Field, (2000). The 1999 Earthquake Sequence along the North Anatolia Transform
at the Juncture between the Two Main Ruptures, Barka et al. edit. in The 1999
Izmit and Duzce Earthquakes: preliminary results, Istanbul technical universit,
209-223.
Sharp, R.J., (1967). San Jacinto fault zone in the peninsular ranges of southern
California, Geol. Soc. Am. Bull., 78, 705-730.
Shapiro, N.M., M. Campillo, L. Stehly and M.H. Ritzwoller, (2005). High resolution
surface wave tomography from ambient seismic noise, Science, 307, 1615-1618.
Shi, Z. and Y. Ben-Zion, (2006). Dynamic rupture on a bimaterial interface governed by
slip-weakening friction, Geophys. J. Int., 165, doi: 10.1111/j.1365-
246X.2006.02853.x.
Simons, F. J., B.D.E. Dando and R. M. Allen, (2006). Automatic detection and rapid
determination of earthquake magnitude by wavelet multiscale analysis of the
primary arrival, Earth and Planetary Science Lett., 250, 214-223.
Spudich, P. and K.B. Olson, (2001). Fault zone amplified waves as a possible seismic
hazard along the Calaveras fault in central California, Geophys. Res. Lett., 28,
2533-2536.
Thurber, C., S. Roecker, W. Ellsworth, Y. Chen, W. Lutter and R. Sessions, (1997).
Two-dimensional seismic image of the San Andreas fault in the northern Gabilan
range, central California: Evidence for fluids in the fault zone, Geophys. Res.
Lett., 24, 1591-1594.
159
Thurber, C., S. Roecker, H. Zhang, S. Baher and W. Ellsworth, (2004). Fine-scale
structure of the San Andreas fault and location of the SAFOD target earthquakes,
Geophys. Res. Lett., 31, L12S02, doi 10.1029/2003GL019398.
Umeda, Y., (1990). High amplitude seismic waves radiated from the bright spot of an
earthquake, Tectonophysics, 175, 81-92.
Vidale, J.E., D.V. Helmberger and R.W. Clayton, (1985). Finite-difference
seismograms for SH waves, Bull. Seis. Soc. Am., 75, 1765-1782.
Wagner, G. S., (1998). Local wave propagation near the San Jacinto fault zone,
Southern California: observations from a three component seismic array, J.
Geophys. Res., 103, 7231-7246.
Waldhauser, F. and W. Ellsworth, (2000). A double-difference earthquake location
algorithm: method and application to the northern Hayward Fault, Bull. Seism.
Soc Am., 90, 1,330-1368.
Wells, D.L. and K.J. Coppersmith, (1994). New Empirical Relationships among
Magnitude, Rupture Length, Rupture Width, Rupture Area, and Surface
Displacement, Bull. Seis. Soc. Am., 84, No. 4, 974-1002.
Weertman, J., (1980). Unstable slippage across a fault that separates elastic media of
different elastic constants, J. Geophys. Res., 85, 1455-1461.
Wu, Y-M., and L. Zhao, (2006). Magnitude estimation using the first three second P-
wave amplitude in earthquake early warning, Geophy. Res. Lett., 33,
doi:10.1029/2006GL026871.
Wu, Y-M., H. Kanamori, R.M. Allen and E. Hauksson, (2007). Determination of
earthquake early warning parameters, τ
c
and P
d
, for Southern California, Geophys.
J. Int., doi: 10.1111/j.1365-246X.2007.03430.x.
Wurman, G., R.M. Allen and P. Lombard, (2007). Toward earthquake early warning in
Northern California, J. Geoph. Res., 112, doi:10.1029/2006JB004830.
Yamada, T., J. J. Mori, S. Ide, R. E. Abercrombie, H. Kawakata, M. Nakatani, Y. Iio,
and H. Ogasawara. (2007). Stress drops and radiated seismic energies of
microearthquakes in a South African gold mine, J. Geophys. Res., 112, B03305,
doi:10.1029/2006JB004553.
Zollo, A., M. Lancieri, and S. Nielsen, (2006). Earthquake magnitude estimation from
peak amplitudes of very early seismic signals on strong motion records, Geophy.
Res. Lett., 33, doi:10.1029/2006GL027795.

160
Zollo, A., Lancieri, and S. Nielsen, (2007). Reply to comments by P.Rydelek et al. on
“Earthquake magnitude estimation from peak amplitude of very early seismic
signals on strong motion records”, Geophys. Res. Lett., 34, doi:
10.1029/2007GL030560

Abstract (if available)

Abstract Seismic signals associated with near-fault waveforms are examined to determine fault zone structure and scaling of earthquake properties with event magnitude. The subsurface structure of faults is explored using fault zone head and/or trapped waves, while various signals from the early parts of seismograms are investigated to find out the extent to which they scale with magnitude.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Multi-scale imaging and monitoring of crustal and fault zone structures in southern California

PDF

Multi-scale imaging of major fault zones in Southern California

PDF

Heterogeneity of earthquake stress drops, focal mechanisms and active fault zones

PDF

High-resolution imaging and monitoring of fault zones and shallow structures: case studies in southern California and on Mars

PDF

Applying automated techniques to large seismic datasets for systematic analyses of phases, source, and structure

PDF

Multi-scale damage signatures across major strike-slip faults

PDF

Multi-scale imaging of the fault zone velocity structure: double-difference tomography, inversion of fault zone headwaves, and fault zone sensitivity kernels

PDF

Symmetry properties, pulverized rocks and damage architecture as signatures of earthquake ruptures

PDF

Detailed properties of seismic waveforms and earthquake sources in Southern California

PDF

Volumetric interactions between major ruptures and fault zones illuminated by small earthquake properties

PDF

Seismic velocity structure and spatiotemporal patterns of microseismicity in active fault zones of Southern California

PDF

Analysis of waveform and catalog data of aftershocks for properties of earthquakes and faults

PDF

Spatiotemporal variations of stress field in the San Jacinto Fault Zone and South Central Transverse Ranges

PDF

Dynamic rupture processes and seismic radiation in models of earthquake faults separating similar and dissimilar solids

PDF

Direct observation of fault zone structure and mechanics in three dimensions: a study of the SEMP fault system, Austria

PDF

Elements of seismic structures near major faults from the surface to the Moho

PDF

Stress-strain characterization of complex seismic sources

PDF

Reactivation of multiple faults in oilfields with injection and production

PDF

Detection and modeling of slow slip events as creep instabilities beneath major fault zones

PDF

Statistical analyses of ambient seismic noise spectra with applications to detecting imperfectly diffuse wave field and deriving attenuation and phase velocity information

Asset Metadata

Creator Lewis, Michael Antony (author)

Core Title Determining fault zone structure and examining earthquake early warning signals using large datasets of seismograms

School College of Letters, Arts and Sciences

Degree Doctor of Philosophy

Degree Program Geological Sciences

Publication Date 09/10/2008

Defense Date 05/05/2008

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag earthquake early warning,fault zones,head waves,OAI-PMH Harvest,trapped waves

Place Name Anza (city or populated place), California (states), faults: North Anatolian (geographic subject), faults: San Andreas (geographic subject), faults: San Jacinto (geographic subject), Hollister (city or populated place), mine sites: TauTona (geographic subject), South Africa (countries)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Ben-Zion, Yehuda (committee chair), Lee, Vincent W. (committee member), Sammis, Charles G. (committee member)

Creator Email malewis@usc.edu,omnifluff@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-m1595

Unique identifier UC1418098

Identifier etd-Lewis-2345 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-111252 (legacy record id),usctheses-m1595 (legacy record id)

Legacy Identifier etd-Lewis-2345.pdf

Dmrecord 111252

Document Type Dissertation

Rights Lewis, Michael Antony

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Repository Name Libraries, University of Southern California

Repository Location Los Angeles, California

Repository Email cisadmin@lib.usc.edu