Finite-Fault Rupture Detector (
FinDer
):
Going Real-Time in Californian
ShakeAlert
Warning System
by M. Böse, C. Felizardo, and T. H. Heaton
INTRODUCTION
Rapid detection of local and regional earthquakes and issu-
ance of fast alerts for impending shaking is considered ben-
eficial to save lives, reduce losses, and shorten recovery times
after destructive events (
Allen
et al.
, 2009
). Over the last two
decades, several countries have built operational earthquake
early warning (
EEW
) systems, including Japan (
Hoshiba
et al.
,
2008
), Mexico (
Espinosa-Aranda
et al.
,1995
), Romania
(
M
ă
rmureanu
et al.
, 2011
), Turkey (
Erdik
et al.
, 2003
), Tai-
wan (
Hsiao
et al.
,2011
), and China (
Peng
et al.
, 2011
). Other
countries, such as the United States (
Böse, Allen,
et al.
,2013
),
Italy (
Satriano
et al.
,2011
), and Switzerland (
Behr
et al.
,
2015
), are currently developing systems or evaluating algo-
rithms in their seismic real-time networks.
Over the past eight years, scientists at the California Institute
of Technology (Caltech), the University of California
–
Berkeley,
the University of Southern California, the University of Washing-
ton, the U.S. Geological Survey (USGS), and the Swiss Federal
Institute of Technology (ETH Zurich, Switzerland) developed the
U.S. west-coast-wide
ShakeAlert
EEW
demonstration system for
California, Oregon, and Washington (
Böse, Allen,
et al.
,2013
).
Real-time alerts are currently being shared with around 300 indi-
viduals, companies, and emergency response organizations to
gather feedback about system performance, to educate potential
end users about
EEW
, and to identify needs and applications of a
future operational warning system.
To quickly determine earthquake magnitudes and (point-
source) hypocenter locations, the Californian
ShakeAlert
system
processes and interprets real-time waveform data streams from sev-
eral hundred seismic broadband velocity and strong-motion accel-
eration sensors as part of the California Integrated Seismic
Network (CISN). Designed as a hybrid system,
ShakeAlert
com-
bines the estimates from one single-sensor and two network-based
EEW
algorithms that run in parallel, including
τ
c
P
d
Onsite
(
Wu
et al.
,2007
;
Böse, Hauksson, Solanki, Kanamori, Wu,
et al.
,2009
;
Böse, Hauksson, Solanki, Kanamori, and Heaton, 2009
),
ElarmS
(
Allen, 2007
;
Kuyuk
et al.
, 2014
), and the
Virtual Seismologist
(
Cua
et al.
, 2009
;
Behr
et al.
, 2015
). Once an earthquake has been de-
tected, each algorithm starts indep
endently sending reports to the
ShakeAlert
Decision Module, which
—
based on these reports
—
de-
termines the most probable event parameters (location, magnitude,
and origin time) and uncertainties, which are then continuously
updated with progressing time as new data arrive. Using these esti-
mated source parameters, the
ShakeAlert
UserDisplay software,
which is installed at the end user
’
s site and subscribes to the
Shake-
Alert
system, predicts and displays the arrival and intensity of ex-
pected peak shaking (
Böse, Allen,
et al.
, 2013
).
Encouraged by its promising real-time performance during
the recent moderate-size 2014
M
5.1 La Habra and
M
6.0 South
Napa earthquakes in southern and northern California, the
Shake-
Alert
Science team decided in 2014 to include the Finite-fault
rupture Detector algorithm,
FinDer
(
Böse
et al.
,2012
), as a fourth
seismic algorithm in the
ShakeAlert
systemtoenabletheuseof
rupture-to-site distances in grou
nd-motion predictions and hence
improve
EEW
performance during large-magnitude (
M
≥
6
:
0
)
earthquakes. Based on a rapid near- or far-source classification
and comparison with precalculated templates,
FinDer
provides
rapid estimates of the centroid position (
lat
centroid
,
lon
centroid
),
length
L
, and strike
Θ
of an ongoing fault rupture, assuming
alinesource(Fig.
1
). The integration of
FinDer
into
ShakeAlert
was finalized in April 2015. Currently,
FinDer
is operated as
a MathWorks MATLAB (
www.mathworks.com/products/
matlab
, last accessed September 2015) stand-alone code with
a C++ waveform-processing module installed at three CISN da-
tacenters at Caltech/USGS-Pasadena, UC Berkeley, and USGS-
Menlo Park. Translation of the
FinDer
algorithm to C++ by
USGS-Menlo Park, ETH Zurich, and Caltech is currently
under way.
We will start this article with a brief review of the
FinDer
algorithm developed by
Böse
et al.
(2012)
, followed by a sum-
mary and demonstration of recent improvements to the algo-
rithm, including error estimates, usage of generic and fault-
specific templates, and extension
to subduction-zone earthquakes.
We will demonstrate and evaluate the real-time and off-line per-
formance of
FinDer
during the recent 2014
M
5.1 La Habra and
M
6.0 South Napa (California) earthquakes, as well as the 2010
M
7.2 El Mayor
–
Cucapah (Mexico) and 2011
M
9.0 Tohoku-
Oki (Japan) earthquakes.
REVIEW AND EXTENSION OF THE
FinDer
ALGORITHM
The
FinDer
algorithm uses image-recognition techniques to
detect and model finite-fault ruptures using spatial images
of observed ground motions and theoretical templates modeled
from empirical ground-motion prediction equations (
GMPE
s).
1692 Seismological Research Letters Volume 86, Number 6 November/December 2015
doi: 10.1785/0220150154
As described by
Böse
et al.
(2012)
and illustrated in Figure
1
,
FinDer
processing involves three main steps: (1) near- or far-
source classification of spatiall
y interpolated high-frequency ob-
servations, which gives a binary image
I
x; y
;(2)applicationof
matching by correlation (e.g.,
Gonzales
et al.
, 2004
)todetermine
the correlation
R
x; y
j
L;
Θ
of template
T
x; y
j
L;
Θ
and
I
x; y
at any spatial location
x; y
; the geographic coordinates
lat
centroid
and
lon
centroid
of the rupture centroid come from the template
spatial location that corresponds to the maximum spatial integral
of the given template cross correlation with the data; (3) minimi-
zation of the misfit between image
I
x; y
and templates
T
x; y
j
L;
Θ
gives
L
,
Θ
,and
lat
centroid
and
lon
centroid
of the as-
sumed line source. The output parameters are updated whenever
a new set of ground-motion amplitudes becomes available (for
instance, every second). This way,
FinDer
keeps track of a rupture
evolving from a small point source to a large finite-rupture event.
Images
I
x; y
(Fig.
1
, left) are obtained from map-
projected and phase-independent (
Böse
et al.
, 2012
) interpo-
lated logarithmic peak ground acceleration (
PGA
) amplitudes.
We define
PGA
as the largest absolute value recorded on any of
the three ground-motion components over a constant time win-
dow (e.g., 120 s), corresponding to the rupture duration of the
largest expected earthquake in a given region. Because the large
high-frequency motions are typically observed at only small rup-
ture-to-site distances (e.g.,
Yamada
et al.
, 2007
), a ground-motion
threshold,
gm_threshold
, is applied for near- or far-source clas-
sification, where image pixels with
PGA
obs
≥
gm
threshold
are
set to 1 and otherwise are set to 0. For the examples shown in
this article, we set
gm
threshold
70 cm
=
s
2
. We are using the
term
“
near source
”
not in a strictly physical sense, but rather to
characterize the spatial zone around the fault rupture with
PGA
values exceeding some arbitrary threshold.
observed
(high-frequency)
ground motion
amplitudes
interpolated
ground motion
map
binary map
“image”
centroid, L
,
θ
, likelihood
2D finite-fault
rupture
near/far source classification
generic
symmetric
templates
generic
asymmetric
templates
(for subduction
zones)
“Matching by
correlation” ,
minimization
of misfit
Observations
Models (“
templates”
)
fault-specific
templates
length
L
strike
Θ
Θ
▴
Figure 1.
Scheme of
FinDer
algorithm. Processing involves three steps: (1) near- or far-source classification of spatially interpolated
high-frequency observations, which gives a binary image (left); (2) application of matching by correlation to determine the spatial cor-
relation of the image and a template (right) for a given rupture length
L
and strike
Θ
; the location of maximum correlation gives the rupture
centroid. Templates are generated from ground-motion prediction equations (
GMPE
s). Currently, we use three sets of templates (right):
generic symmetric, generic asymmetric (for subduction zones), and fault specific; (3) minimization of the misfit between image and
templates gives
L
,
Θ
, centroid, and likelihood
L
of the assumed line source. Estimates are updated regularly as new data become
available (e.g., every second).
Seismological Research Letters Volume 86, Number 6 November/December 2015 1693
Templates
T
x; y
j
L;
Θ
(Fig.
1
, right) are precalculated
binary images generated for various rupture lengths
L
and
strikes
Θ
from empirical
GMPE
s, in which pixels close to
the rupture (
PGA
GMPE
≥
gm
threshold
) are set to 1 and all
other pixels (
PGA
GMPE
<gm
threshold
) to 0. The templates
thus consist of equidistant buffer zones around the line source,
bordered by two parallel lines and two semicircles, which result
in quasi-elliptical shapes (Fig.
1
, right). The relationship be-
tween the cutoff distance and the rupture length
L
(and thus
magnitude
M
) is nonlinear due to the nonlinear nature of the
underlying
GMPE
s for
PGA
, which saturate at large magnitudes
and close distances.
Because of the assumed symmetry of the generic templates,
the results of
Θ
and
Θ
180
° are equivalent. It is assumed that
near-source
PGA
s are not affected by rupture directivity, as dis-
cussed by
Spudich and Chiou (2008)
. For the application of
FinDer
to earthquakes in California and northern Mexico as de-
scribedinthisarticle,weuse
GMPE
sby
Cua and Heaton (2009)
and empirical rupture length
–
magnitude relations by
Wells
and Coppersmith (1994)
; for application to the Tohoku-Oki
subduction-zone earthquake (Japan), we apply relations after
Si
and Midorawaka (2000)
and
Blaser
et al.
(2010)
.
Modeled templates
T
x; y
j
L;
Θ
are used as spatial filters to
determine the correlation
R
x; y
j
L;
Θ
with an observed image
I
x; y
at each spatial location of
T
in
I
(the spatial extent of
I
is
much larger than of
T
). This method is known as matching by
correlation (e.g.,
Gonzales
et al.
,2004
;Fig.
1
,middle).Toincrease
computational speed, we implemented the algorithm in the Fou-
rier domain by taking advantage of the correlation theorem:
EQ-TARGET;temp:intralink-;df1;311;418
I
x; y
⋆
T
x; y
j
L;
Θ
⇔
~
I
k
x
;k
y
~
T
k
x
;k
y
j
L;
Θ
;
1
which relates the spatial correlation to the product of the Fourier
transforms
~
I
k
x
;k
y
and
~
T
k
x
;k
y
j
L;
Θ
in the wavenumber do-
main
k
x
;k
y
, in which the star denotes correlation and the asterisk
denotesthecomplexconjugate.Thisapproachiscomputationally
highly efficient and could be easily parallelized for additional speed,
though this has not been done yet in the current
FinDer
imple-
mentation (
Böse
et al.
, 2012
).
Although
FinDer
was originally developed for generic sym-
metric templates assuming simple line sources with symmetric
quasi-elliptical ground-shaking distributions around the finite-
fault rupture (
Böse
et al.
,2012
), we recently added a new set
of fault-specific templates to provide more accurate rupture es-
timates along curved fault segments, such as along the Big Bend
section of the San Andreas fault (
SAF
) in southern California. In
addition, we created a third set of generic asymmetric templates,
particularly designed to meet the requirements of subduction-
zone environments, where seismic observations are typically
limited to the hanging wall side of the fault rupture (Fig.
1
,
right; shaded area corresponds to downward-dipping fault
plane). Table
1
gives a summary of the three currently used
sets of templates.
The finite-fault parameters (
L
,
Θ
, and centroid position)
are determined in a two-step procedure. In the first step,
matching by correlation is used to determine the optimum spa-
tial position
x
′
;y
′
of a given template
T
in image
I
. The po-
sition
x
′
;y
′
refers to the lower left corner of the template,
assuming that the origin (
0
,
0
) of both
T
and
I
are in the lower
Table 1
Description of Template Sets Used in This Study
Generic Symmetric
Templates
Fault-Specific Templates*
Generic Asymmetric Templates
(Subduction Zones)
Range
L
:5
–
300 km,
Δ
L
5km
L
:15
–
350 km,
Δ
L
5km
L
:10
–
450 km,
Δ
L
10 km
Θ
:0°
–
179°,
ΔΘ
1
°
Θ
: follows fault-surface-trace of SAF
for various rupture start and end
positions, that is
Θ
is not fixed in a
given template; three overlapping
subsets such as southern SAF, middle
SAF, and northern SAF
Θ
:0
–
357°,
ΔΘ
3
°
Size
77
×
77
cells
(each cell:
5
×
5km
)
56
×
123
(southern SAF)
155
×
155
cells
(each cell:
10
×
10 km
)
78
×
78
(middle SAF)
67
×
78
(northern SAF) (each cell:
5
×
5km
)
Total number
10,800
1022
2700
456 (southern SAF)
378 (middle SAF)
188 (northern SAF)
GMPE
†
Cua and Heaton (2009)
Cua and Heaton (2009)
Si and Midorikawa (2000)
L
–
M
relation
Wells and Coppersmith (1994) Wells and Coppersmith (1994)
Blaser
et al.
(2010)
*SAF, San Andreas fault.
†
GMPE, ground-motion prediction equation.
1694 Seismological Research Letters Volume 86, Number 6 November/December 2015
left corner. The rupture centroid position (
lat
centroid
,
lon
centroid
)
is determined from the template midpoint at the position
x
′
;y
′
. In the second step, we minimize the normalized sum
of squared errors to find the optimum template and thus
L
and
Θ
. We modify the original error function defined in
Böse
et al.
(2012)
and define the new misfit
E
for a given image and a
generic or fault-specific template for length
L
and strike
Θ
as
EQ-TARGET;temp:intralink-;df2a;41;661
E
L;
Θ
P
x
′
;y
′
I
x
x
′
;y
y
′
−
T
x
′
;y
′
j
L;
Θ
2
P
x
′
;y
′
σ
2
m
x
′
;y
′
I
x
x
′
;y
y
′
T
x
′
;y
′
j
L;
Θ
;
2a
in which the summation is done over
x
′
0
...
w
−
1
,
y
′
0
...
h
−
1
for template dimensions
w
×
h
. The weighting
term
σ
2
m
x
′
;y
′
in equation
(2a)
accounts for the model uncer-
tainties of the
GMPE
s(
σ
GMPE
) used for template generation,
which we define as
EQ-TARGET;temp:intralink-;df2b;41;503
σ
2
m
x
′
;y
′
8
<
:
R
−
R
min
=
R
gm
threshold
−
R
min
1if
R
min
≤
R
≤
R
gm
threshold
2
−
R
−
R
gm
threshold
=
R
max
−
R
gm
threshold
if
R
gm
threshold
≤
R
≤
R
max
1otherwise
;
2b
in which
R
gm
threshold
is the rupture-to-cell distance at which
PGA
mean
gm
threshold
(e.g.,
gm
threshold
70 cm
=
s
2
),
R
min
is the rupture-to-cell distance at which
PGA
mean
−
σ
GMPE
gm
threshold
,and
R
max
is the rupture-to-cell distance at which
PGA
mean
σ
GMPE
gm
threshold
. Equation
(2b)
describes a
simple bilinear function that increases linearly from 1 to 2 for
distances between
R
min
and
R
mean
and then decreases linearly
back to 1 for distances between
R
mean
and
R
max
. For distances
smaller than
R
min
or larger than
R
max
σ
2
m
x
′
;y
′
is 1.
The template (and thus
L
and
Θ
) that shows the best
agreement with image
I
x; y
is found through minimization
of the misfit function in equation
(2a)
,thatis
E
L;
Θ
→
min
.
For the generic templates, we apply the simplex search algo-
rithm by
Lagarias
et al.
(1998)
with upper and lower bounds
for
L
and
Θ
, which is a direct search algorithm that does not
require numerical or analytic gradients. The convergence rate
toward the minimum is excellent and usually requires 30
–
50
iterations only. For fault-specific templates, we consider all pos-
sible templates (and thus all
L
) for a given region (northern,
middle, and southern
SAF
; Table
1
); and, for each of these tem-
plates, we determine
E
and the optimum solution. If its mini-
mum misfit is smaller than for the generic templates (for a
simple 1D line source), the fault-specific template is preferred,
and otherwise we use the generic one. At the moment, we test
both generic symmetric and fault-specific templates in California
and asymmetric templates only in subduction-zone environ-
ments (Japan). In our current MATLAB implementation,
the total optimization procedure takes
∼
1
–
2s
on an i686
@2826.177 MHz machine.
As a new feature of
FinDer
, we recently started to estimate
uncertainties for the best solutions of rupture length
~
L
and strike
~
Θ
by approximation of the multivariate likelihood function
L
.For
generic templates, we keep the centroid position and rupture
length of the best solution (
e
lat
centroid
;
f
lon
centroid
;
~
L
) fixed and vary
the strike
Θ
from 0° to 179° to determine the estimated univari-
ate likelihood function for strike as
EQ-TARGET;temp:intralink-;df3;311;655
L
Θ
j
~
L; I
x; y
1
σ
d
2
π
p
e
−
0
:
5
E
σ
2
d
:
3
Equation
(3)
is a probability density function that de-
scribes how well the predictions from a given set of model
parameters
Θ
;L
—
and thus from a given template
T
—
fit
the observed data image
I
. The standard deviation
σ
d
describes
the combined effects of assumed errors in the observations (in-
cluding effects of interpolation) and the assumed prediction
errors in the forward model, that is the
GMPE
s used for tem-
plate generation. We experimented with various
σ
d
and found
that
σ
d
0
:
1
gives reasonable results that are in good agree-
ment with geodetic inversions (
Minson, Murray,
et al.
, 2014
).
In the next step, we keep the centroid position and strike of
the best solution (
e
lat
centroid
;
f
lon
centroid
;
~
θ
)fixedandwalkthrough
all
L
(5
–
300 km) to determine
L
L
j
~
Θ
;I
x; y
, similar to equa-
tion
(3)
. For fault-specific templates, we determine the misfit for
each template and thus for each
L
, and can easily apply equa-
tion
(3)
to calculate the likelihood function.
Although a simultaneous inversion of the likelihood func-
tions for
L
,
Θ
, and the centroid position would be preferred,
such computations are unrealistic under the strict time con-
straints in
EEW
. Also, assuming that
L
and
Θ
are largely uncor-
related appears valid for at least moderate-to-large earthquake
ruptures.
FinDer
PERFORMANCE
We will now demonstrate the performance of
FinDer
with the
new features described above for four recent earthquakes.
Three of them are crustal strike-slip, normal, and oblique fault-
ing events that occurred recently in California (2014
M
5.1 La
Habra and
M
6.0 South Napa earthquakes) and in northern
Mexico (2010
M
7.2 El Mayor
–
Cucapah); the fourth event is
the 2011
M
9.0 Tohoku-Oki (Japan) subduction-zone earth-
quake. The latter is particularly challenging because station
coverage is along one side of the fault rupture only. USGS
ShakeMaps
(
Wald
et al.
, 1999
) and finite-source solutions
Seismological Research Letters Volume 86, Number 6 November/December 2015 1695
(
Mai and Thingbaijam, 2014
) for the four earthquakes are
shown in Figure
2
.
Real-Time Performance
The
M
5.1 La Habra (California) earthquake on 28 March 2014
was felt widely throughout Orange, Los Angeles (LA), Ventura,
Riverside, and San Bernardino Counties with a maximum ob-
served instrumental intensity of VII close to the epicenter, lo-
cated 1 km east of La Habra (Table
2
). The Global Centroid
Moment Tensor (Global CMT) moment tensor (see
Data
and Resources
) shows oblique faulting, with a northward-dipping
plane that approximately aligns with the Puente Hills thrust fault.
The strike of the rupture is determined as
Θ
obs
232
°, which,
due to symmetry in our 2D approach, is equivalent to
Θ
obs
52
°(Table
2
). The
ShakeAlert
EEW
demonstration system de-
tected the La Habra earthquake within 4.3 s of the event origin
and provided 4 s of warning to Caltech in Pasadena, around
30 km from the epicenter (see
Data and Resources
).
During the La Habra earthquake,
FinDer
was running in
real-time test mode at Caltech and getting waveform streams
from 420 CISN strong-motion stations. However,
FinDer
was
not connected to the
ShakeAlert
system. Although processing
speed was not yet optimized,
FinDer
detected the earthquake
within 11.5 s of the event origin and sent out a sequence of
three internal reports with updated event information. The fi-
nite-fault solutions in these reports show an excellent agreement
with the CMT solution (Table
2
), as well as with the aftershock
distribution of the La Habra earthquake (
L
obs
≈
10 km
;Fig.
3
).
▴
Figure 2.
U.S. Geological Survey (USGS)
ShakeMaps
showing the spatial distribution of observed instrumental intensities (modified
Mercalli intensity [
MMI
]) in the (a)
M
5.1 La Habra, (b)
M
6.0 South Napa, (c)
M
7.2 El Mayor
–
Cucapah, and (d)
M
9.0 Tohoku-Oki earth-
quakes. Stars show epicenters, gray lines show 1D line source (b and c) or 2D rupture plane (c and d), respectively.
1696 Seismological Research Letters Volume 86, Number 6 November/December 2015
The 2014
M
6.0 South Napa earthquake on 24 August at
10:20:44 UTC occurred close to the West Napa fault. The
Global CMT solution shows a strike-slip event with
Θ
obs
157
°; from the aftershock distribution and empirical
relations by
Wells and Coppersmith (1994)
, we estimate
L
obs
≈
15
–
20 km
.
ShakeAlert
detected this event within
5.1 s (
Grapenthin
et al.
, 2014
). The initial event location
and moment magnitude were off by 3 km and 0.9 magnitude
units (estimated magnitude was 5.1) compared with the Ad-
vanced National Seismic Network catalog.
FinDer
, again oper-
ating in real-time test mode but not connected to the
ShakeAlert
system, detected the event within 16 s. The fairly sparse distri-
bution of strong-motion stations around the Napa earthquake
led to larger uncertainties in the
FinDer
-estimated rupture strike,
compared, for example, with the La Habra earthquake. Although
the final strike estimate is off by around
∼
30
°
–
40
°, the CMT
solution of
Θ
obs
157
° is within the
FinDer
real-time-estimated
uncertainty range (Fig.
3
).
Off-Line Performance
Following the La Habra and South Napa earthquakes, we started
to improve the
FinDer
algorithm and code in terms of robust-
ness and processing speed. To test our new implementation, we
replayed archived waveforms from historic and simulated earth-
quakes using the
EarthWorm Tankplayer
software. In this
article, we present the off-line-results from the La Habra and
Napa earthquakes, as well as from the 2010
M
7.2 El Mayor
–
Cucapah and 2011
M
9.0 Tohoku earthquakes.
Figure
4
shows the
FinDer
-estimated rupture length
L
,
strike
Θ
, and magnitude
M
for the four events. Magnitudes
were estimated from empirical L
–
M relations by
Wells and
Coppersmith (1994)
and
Blaser
et al.
(2010)
. The shaded areas
show the estimated uncertainties of the rupture parameters;
they refer to the standard deviation estimated from the like-
lihood functions
L
, assuming Gaussian distributions; if
L
is
non-Gaussian (e.g., likelihood function for strike at an early
stage of the La Habra earthquake, Fig.
5
), only the value with
maximum likelihood is shown. Figure
5
displays three ran-
domly picked screenshots for each event to illustrate the tem-
poral rupture evolution.
The off-line-predicted rupture parameters for the La
Habra and South Napa earthquakes (Figs.
4
and
5
) are quite
similar to the real-time results (Fig.
2
). Again, because of the
sparseness of seismic instrumentation around the South Napa
earthquake, uncertainties in the
FinDer
output are quite large.
The
M
7.2 El Mayor
–
Cucapahearthquakeon4April2010
occurred in northern Baja California in an area with a
high level of historical seismicity, around 60 km south of the
Mexico
–
United States border (
Hauksson
et al.
, 2011
). The event
exhibited complex faulting, p
ossibly starting as a smaller
M
∼
6
normal-faulting earthquake, followed
∼
15 s
later by the normal
or strike-slip faulting mainshock with
Θ
obs
132
°(
Wei
et al.
,
2011
). The aftershock zone extends from the southern end of the
ElsinorefaultzoneinsouthernCaliforniaalmosttothenorthern
tip of the Gulf of California, over a length of
L
obs
120 km
(
Hauksson
et al.
,2011
).
AtthetimeoftheElMayor
–
Cucapah earthquake,
Shake-
Alert
did not yet exist. However, all three point-source algorithms
that are contributing to the current system were running in real-
time test mode:
τ
c
P
d
Onsite
detected the event within 16 s,
Virtual Seismologist
within 28 s of origin. Both algorithms esti-
mated the event as
M
∼
6
, that is, they were obviously confused
bythesmallerforeshockandunderestimatedthesizeofthemain-
shock.
ElarmS
did not process this event, because it was outside
of the CISN network. At the time of the El Mayor
–
Cucapah
earthquake, none of the Centro de Investigación Científica y de
Educació n Superior de Ensenada (CICESE) Baja California,
Mexico Seismic Network (Red Sismica del Noroeste de Mexico
[RESNOM]) stations deployed in northern Mexico was stream-
ing real-time waveform data to CISN, so only Californian sta-
Table 2
Earthquake Source Parameters
Earthquake
Origin Time (UTC) (yyyy/
mm/dd hh:mm:ss)
Hypocenter (latitude/
longitude [º], depth [km])
Moment
Magnitude (M)
Rupture Length
L
obs
(km)
Rupture
Strike
Θ
obs
(°)
La Habra
(oblique)
2014/03/29 04:09:42
33.932/
−
117.917, 4.8
5.1
10*
52
†
,
‖
South Napa
(strike slip)
2014/08/24 10:20:44
38.215/
−
122.312, 11.3
6.0
15
–
20*
157
†
,
‖
El Mayor
–
Cucapah
(complex)
2010/04/04 22:40:42
32.259/
−
115.287, 10
7.2 (starting as
M
∼
6
:
3
‡
)
120
‡
132
‡
,
‖
Tohoku-Oki
(megathrust)
2011/03/11 05:46:23
38.297/142.272, 30
9.0
380
§
202
§
*Aftershock distribution.
†
Global Centroid Moment Tensor (CMT) catalog (see
Data and Resources
).
‡
Wei
et al.
(2011)
.
§
Lay
et al.
(2011)
.
‖
Because of our 2D approximation of the rupture,
Θ
obs
and
Θ
obs
−
180
° are equivalent when using symmetric templates.
Seismological Research Letters Volume 86, Number 6 November/December 2015 1697
tions were available in real time. Only recently, eight CICESE
stations were added to the CISN live streams and are now being
processed routinely by
ShakeAlert
.
For our
FinDer
off-line-test, we complemented the CISN
waveform dataset with the data from eight CICESE stations
(Fig.
5
;see
Data and Resources
). Although the near-source sta-
tion coverage is still quite poor, both the length and strike of the
northwestern part of the bilateral rupture can be well recovered
(Figs.
2
,
4
,and
5
). The rupture centroid appears to be slightly
shifted toward the east but is within the uncertainty range of our
current map resolution (5 km). The final rupture length is esti-
mated as 70 km, which corresponds to
M
∼
7
:
1
(
Wells and
Coppersmith, 1994
) and is thus in good agreement with the
observed moment magnitude. The output of
FinDer
seems
to be largely unaffected by the smaller
M
∼
6
foreshock and
inherent complexities of this event.
The last earthquake that we will analyze here is the 2011
M
9.0 Tohoku-Oki (Japan) megathrust earthquake, which
caused tremendous numbers of casualties and damage in Iwate,
Miyagi, and Fukushima, particularly from a Pacific-wide tsu-
nami.
Simons
et al.
(2011)
estimate that peak displacements of
around 60 m occurred along the central section of the fault
plane (Fig.
2
). Kinematic inversions by
Minson, Simons,
et al.
(2014)
suggest that the rupture propagated with a low average
velocity of about
1
:
2km
=
s
and that most of the rupture oc-
curred within
∼
120
–
125 s
.
FinDer Real-Time Performance:
2014 La Habra:
(right)
L
obs
=~ 10 km
L
est
= 10
±
5 km
M
obs
= 5.1
M
est
= 5.8
±
0.5
Θ
obs
=~ 52
o
Θ
est
= 64
±
8
o
2014 South Napa:
(bottom)
L
obs
=~ 15 km
L
est
=~ 10
±
5 km
M
obs
= 6.0
M
est
= 5.8
±
0.5
Θ
obs
=~ 157
o
Θ
est
=~ 119
±
34
o
(a)
(b)
(c)
FinDer
+34
o
▴
Figure 3.
FinDer
real-time output for the 2014 (a)
M
5.1 La Habra and (b, c)
M
6.0 South Napa earthquakes. First reports were available
within 11.5 s (La Habra) and 16 s (South Napa) of event origin. Maps and data block show results from the final
FinDer
report. Triangles and
squares, respectively, mark the California Integrated Seismic Network (CISN) strong-motion stations at which peak ground accelerations
are
PGA
obs
<
70 cm
=
s
2
and
PGA
obs
≥
70 cm
=
s
2
. (b) Because of sparse instrumentation, the strike
Θ
of the Napa earthquake is not well
constrained. Though the final strike estimate is off by
∼
40
°, the Global Centroid Moment Tensor (CMT) solution (
Θ
obs
157
°) is within the
FinDer
real-time-estimated uncertainty range. Solutions in (b) and (c) are equivalent, because they are producing almost the same misfit
between observed and predicted high-frequency ground motions. For both earthquakes, there is an excellent agreement between the
FinDer
-predicted line source (red), CMT solutions (Table
2
), and distribution of aftershocks (yellow circles).
1698 Seismological Research Letters Volume 86, Number 6 November/December 2015
The EEW system operated by the Japanese Meteorological
Agency (JMA;
Hoshiba
et al.
,2008
) detected the event within
∼
25 s
from event origin and successfully released a warning to
the Japanese public. However, because the event magnitude was
underestimated and source finiteness not considered, ground
motions were largely underestimated, in particular in the Kanto
district (
Hoshiba
et al.
,2011
).
For our
FinDer
off-line test, we use strong-motion records
of the Tohoku-Oki earthquake from 273 K-NET stations op-
erated by the Japanese National Research Institute for Earth
Science and Disaster Prevention (NIED; see
Data and Resour-
ces
). Because the seismic observations in our test are limited to
onshore stations and are thus restricted to one side of the fault
rupture (Fig.
5
, bottom), we use our generic asymmetric tem-
plate set (Fig.
1
,right;Table
1
). Constrained by the geometry
of the Japan trench, we assume a lower and upper strike limit of
90° and 270°, respectively, that is a westward-orientated dip.
In general,
FinDer
results (Fig.
4
) become available a bit
later than in the JMA
EEW
system, mainly because we are con-
sidering K-NET stations only. Although the strike stays quite
stable after 40 s, we observe a steady increase in rupture length
(and thus magnitude) up to 270 km and
M
8.5 within 160 s
after rupture nucleation, which seems to be consistent with the
source time function of the Tohoku-Oki earthquake (
Minson,
Simons,
et al.
,2014
). The rupture centroid, in general, is not
well constrained for a low-dipping fault plane and is also affected
by boundary conditions of the ground-motion interpolation
caused by the one-sided station coverage. In the test shown
in this article (Figs
4
and
5
), we set image pixels at
≥
10 km
distance from the closest recording station to 0 to avoid
artifacts from interpolation.
Using KiK-net stations in addition to K-NET stations has
no significant impact on the
FinDer
performance. The apparent
underestimation of rupture length (
L
est
≈
270 km
)isthuslikely
not due to station density, but mainly a result of using the
GMPE
by
Si and Midorikawa (2000)
, which predicts quite large shaking
far from the fault rupture. The
FinDer
solution is consistent
with the observed ground-motion distribution and applied
GMPE
and, as such, is useful to predict future shaking at larger
distances (using the same
GMPE
) as needed for
EEW
. For illus-
tration, in Figure
6
we compare the temporal evolution of ob-
served and predicted JMA intensities in Tokyo. The predicted
values are computed from
Midorikawa
et al.
(1999)
and
Si and
Midorikawa (2000)
using the distance to the
FinDer
-estimated
finite-fault rupture. Lead times between predicted and observed
intensity levels range up to 50 s.
Minson, Simons,
et al.
(2014)
report that high-frequency
motions in the Tohoku-Oki earthquake were predominantly
radiated down-dip of the regions of largest fault slip, particu-
larly south of Sendai. Because
FinDer
uses high-frequency mo-
tions to determine rupture dimensions, the lack of correlation
between high- and low-frequency motions could potentially
lead to an underestimation of rupture dimensions, too.
DISCUSSION
Aside from the estimated finite-fault rupture parameters, there
are two further aspects that make
FinDer
distinct from other
current
EEW
algorithms: (1)
FinDer
searches for suspicious
ground-motion patterns by looking at observed ground-motion
images as a whole. This means that the
FinDer
output is a true
network solution that combines spatially distributed observa-
tions instead of averaging over them, and (2)
FinDer
canbeop-
erated in two modes: either as a complete stand-alone algorithm
or by being triggered by another algorithm or module. In the
first mode,
FinDer
scans real-time ground-motion amplitudes
continuously for suspicious patterns, and whenever the correla-
tion between the observed ground motions (in space and time)
and the theoretical templates exceeds some predefined threshold
(e.g.,
R>
0
:
9
), it starts sending out reports with the computed
source parameters, completely autonomous from any other
(point-source) algorithm. In the second mode,
FinDer
requires
a trigger from another algorithm to get started and can then be
used to verify a potential detection and compute source
parameters.
Traditional picker and (pick-to-event) associators are de-
signed to provide high-fidelity locations. If there is a complex
sequence, such as a seismic swarm or increased aftershock ac-
tivity, traditional associators can become confused, because it
becomes unclear which picks go with which event. This was,
for instance, the case during the intense aftershock activity fol-
lowing the
M
9.0 Tohoku-Oki (Japan) earthquake, causing nu-
merous false and missed alerts in the JMA
EEW
system
(
Hoshiba
et al.
, 2011
). In
EEW
, we do not really care about
▴
Figure 4.
Temporal evolution of
FinDer
-predicted (a) rupture
length
L
, (b) moment magnitude
M
, and (c) rupture strike
Θ
in
off-line test runs using archived waveforms of the
M
5.1 La Habra,
M
6.0 South Napa,
M
7.2 El Mayor
–
Cucapah, and
M
9.0 Tohoku-
Oki earthquakes (compare with Table
2
). The shaded areas show
68% confidence intervals estimated from the likelihood functions.
Snapshots are shown in Figure
5
.
Seismological Research Letters Volume 86, Number 6 November/December 2015 1699
▴
Figure 5.
Three randomly chosen times after nucleation taken during off-line test runs of
FinDer
for the (a)
M
5.1 La Habra, (b)
M
6.0
South Napa, (c)
M
7.2 El Mayor
–
Cucapah, and (d)
M
9.0 Tohoku-Oki earthquakes (compare with Fig.
2
). The red lines show
FinDer
-
estimated ruptures; near-source (
PGA
obs
≥
70 cm
=
s
2
) and far-source (
PGA
obs
<
70 cm
=
s
2
) classified stations are shown by squares
and triangles. Templates, which best fit near-source areas (cyan), are shown in gray. The likelihood functions quantifying reliability
of estimated rupture length
L
and strike
Θ
are shown below each map.
1700 Seismological Research Letters Volume 86, Number 6 November/December 2015
these details, we just want to know the general time and place
where the shaking originates.
FinDer
is therefore a simpler
and more robust approach to asso
ciation. These characteris-
tics make
FinDer
likely also well suited for applications in
highly noisy environments, as is typically observed in dense
low-cost seismic and geodetic sensor networks, such as Qua-
keCatcher (
Cochran
et al.
, 2009
), the Community Seismic
Network (
Clayton
et al.
,2011
), or future crowd-sourced
networks based on consumer smart devices (e.g.,
Minson
et al.
,
2015
).
It is important to point out that the
FinDer
-estimated line
sources depend on the
GMPE
used for template generation.
Usually,
FinDer
results agree well with the observed fault rup-
tures; however, there is a chance that source dimensions are
under- or overestimated (as could be seen in the case of the
Tohoku-Oki earthquake; Figs.
4
and
5
). Still, at any given time
after rupture nucleation, the solution is consistent with the cur-
rent spatial seismic ground-motion observations and the applied
GMPE
, and this is exactly what is needed for
EEW
(Fig.
6
). To
minimize the effects of shortcomings in the
GMPE
sonthe
FinDer
result, we build our templates from ground-motion
thresholds rather than complete
GMPE
s. The resulting binary
images and templates are also less affected by uncertainties
caused by the spatial interpolation of ground-motion ampli-
tudes, particularly in sparsely instrumented regions.
To account for increasing uncertainties in spatially inter-
polated ground-motion values with increasing distances from
the recording stations, we recently explored the usefulness of
down-weighting image pixels that are far from observational
sites, similar to what is done in the USGS
ShakeMaps
(
Worden
et al.
, 2010
). However, we found that for large
gm_thresholds
and application to dense station networks as in California, the
effect of down-weighting is quite negligible. However, for sub-
duction-zone events with one-sided station coverage, the ben-
efit of down-weighting becomes more important.
Böse
et al.
(2012)
found that consideration of site effects of seismic
ground motions has no significant impact on the
FinDer
solution.
It is also important to keep in mind that
FinDer
does not
predict the future evolution of the fault rupture. Still,
FinDer
-
estimated rupture dimensions allow prediction of future shaking
at larger distances from the fault (as needed for
EEW
), which are
more accurate for large earthquakes compared with what a
point-source algorithm can achieve, because rupture-to-site dis-
tances can be taken into account. An instructive example is given
by the
M
7.8 ShakeOut scenario earthquake (
Graves
et al.
,2008
)
with a 300-km-long fault rupture along the
SAF
, starting at
Bombay Beach and heading toward LA. Using a point-source
approximation of the earthquake, shaking intensities in LA
are predicted as light to moderate (modified Mercalli intensity
[
MMI
]IV
–
V), whereas shaking in LA predicted with the
FinDer
finite-fault solution is severe (
MMI
VIII;
Böse
et al.
,2014
). This
value is much more consistent with seismic ground-motion sim-
ulations by
Graves
et al.
(2008)
, which predict severe to violent
shaking (
MMI
VIII
–
IX). In a big (finite-source) earthquake,
much larger areas are affected by damaging ground-shaking com-
pared with a small-to-moderate event; thus, many more people
could benefit from
EEW
with warning times of several tens of
seconds (
Heaton, 1985
).
Although
FinDer
itself does not provide predictions of fu-
ture rupture, its output can be used for recognition of the fault
along which rupture is occurring. This information, along with
observed slip amplitudes and known fault characteristics, has po-
tential to provide estimates of future rupture evolution (
Böse
and Heaton, 2010
). Also, as demonstrated by
Böse
et al.
(2014)
using the Southern California Earthquake Center Cyber-
Shake dataset, the
FinDer
output can help to determine the di-
rection in which a fault rupture is propagating and thus can help
to enhance ground-motion predictions with consideration of di-
rectivity effects.
Last but not least,
FinDer
source dimensions and uncer-
tainty estimates can constrain Global Positioning System
(
GPS
)-based inversions of fault slip and magnitudes without
saturation in large earthquakes. With the development of the
FinDer-GPSlip
and
FinDer-BEFORES
algorithms,
Böse,
Heaton, and Hudnut (2013)
and
Minson, Böse,
et al.
(2014)
provide the first seismic-geodetic approaches to
EEW
that are
consistent with both seismic and geodetic observations at any
time after rupture nucleation. Traditionally, finite-fault-slip
models are inverted for a known fixed-fault geometry. In
an
EEW
setting, however, rupture geometry is unknown
apri-
ori
, and thus, it becomes necessary to simultaneously solve for
the fault geometry and slip, which is a nonlinear and compu-
tationally expensive inverse problem (
Minson, Murray,
et al.
,
2014
). Using both seismic and geodetic real-time observations
helps constrain both source geometry and slip and can replace
simplified assumptions such as of San Andreas parallel rup-
ture strike, which is typical in
GPS
inversions (e.g.,
Grape-
nthin
et al.
,2014
).
0
50
100
150
200
250
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Seconds after origin
JMA intensity in Tokyo
observed
predicted
▴
Figure 6.
Temporal evolution of observed and
FinDer
-pre-
dicted Japanese Meteorological Agency (JMA) intensity in Tokyo
during the 2011
M
9.0 Tohoku-Oki earthquake.
Seismological Research Letters Volume 86, Number 6 November/December 2015 1701
CONCLUSIONS AND OUTLOOK
Finite-fault rupture extent and azimuth are crucial for accu-
rately predicting seismic ground motions in large earthquakes
with moment magnitudes
M
≥
6
:
0
where source finiteness,
rupture-to-site distances, and
rupture directivity control
ground shaking. Thus, detecting and modeling finite-fault
ruptures in real time are essential to
EEW
and rapid emer-
gency response. Following a period of extensive real-time and
off-line testing, the finite-fault rupture detector algorithm,
FinDer
(
Böse
et al.
,2012
), was recently successfully integrated
in the California-wide
ShakeAlert
EEW
demonstration sys-
tem. Since April 2015,
FinDer
has been actively contributing
to real-time
ShakeAlert
reports automatically sent to several
hundreds of test users in California, thus complementing the
three previous point-source algorithms (
τ
c
P
d
Onsite
,
ElarmS
,
and
Virtual Seismologist
) and improving
ShakeAlert
perfor-
mance during large earthquakes.
Currently,
FinDer
is operated as a MATLAB stand-alone
code with a C++ waveform-processing module installed at
three CISN datacenters (Caltech/USGS-Pasadena, University
of California
–
Berkeley, and USGS-Menlo Park) and is contin-
uously scanning real-time waveform data streams from around
420 strong-motion stations in California for
PGA
patterns
indicative of earthquakes. In a joint effort of USGS-Menlo
Park, ETH Zurich, and Caltech, we have recently started to
translate the current
FinDer
code to C++. The goal is to ob-
tain a faster and more flexible implementation that allows eas-
ier maintenance and better integration into various seismic
processing systems, including
EarthWorm/AQMS
(as currently
used in
ShakeAlert
) and
SeisComp3
(as currently dominantly
used in Europe;
Olivieri and Clinton, 2012
). The new C++
code will benefit from existing and widely tested free software
libraries for computer vision (
Open Source Computer Vision
;
http://opencv.org
, last accessed September 2015) and geo-
graphic mapping (
Generic Mapping Tools
;
Wessel
et al.
,2013
).
During the next two years, we plan to increase the robust-
ness of the current code and develop the second-generation
FinDer 2
algorithm with important new features: the current
FinDer
code supports utilization of a single ground-motion
threshold only, which in the current
ShakeAlert
installation is
set to
70 cm
=
s
2
. This value is usually exceeded only in moder-
ate-to-large-magnitude earthquakes, for instance within 20 km
from the fault rupture for
M
∼
6
:
0
or within 40 km for
M
∼
7
:
0
(
Cua and Heaton, 2009
). A high threshold provides stable de-
tection of large events at the cost of reduced detection speed and
missed detections of smaller earthquakes.
FinDer 2
will combine multiple ground-motion thresh-
olds: starting with small values that are typically observed in
small earthquakes or shortly after rupture nucleation,
FinDer
2
will progressively increase these thresholds to allow detecting
earthquakes with
M
>
3
:
5
with a gradual refinement of finite-
fault parameter estimates in large events. This refinement will
also benefit from the new error estimates in
FinDer
as de-
scribed in this article.
DATA AND RESOURCES
Seismic waveforms used in this study were downloaded from
California Integrated Seismic Network (CISN) (
http://
www.cisn.org
, last accessed July 2015), NIED (
http://www.
kik.bosai.go.jp
, last accessed July 2015), and Red Sismica del
Noroeste de Mexico (RESNOM)/Centro de Investigación
Científica y de Educació n Superior de Ensenada (CICESE)
(
http://resnom.cicese.mx/sitio/
, last accessed July 2015). A
description of the
ShakeAlert
performance during the 2014
La Habra earthquake was obtained from
http://earthquake.
usgs.gov/earthquakes/eventpage/ci15481673#general_
summary
(last accessed August 2015). The Global Centroid
Moment Tensor (CMT) project database was searched using
www.globalcmt.org/CMTsearch.html
(last accessed August
2015). USGS
ShakeMaps
in Figure
2
were downloaded from
http://earthquake.usgs.gov/earthquakes/shakemap/
(last
accessed August 2015). Maps in Figures
3
and
5
were made using
Generic Mapping Tools
(GMT;
Wessel
et al.
,2013
).
ACKNOWLEDGMENTS
The authors would like to thank Sarah Minson, Deborah
Smith, Egill Hauksson, Men-Andrien Meier, Yannik Behr,
Carlo Cauzzi, John Clinton, and two anonymous reviewers
and editors for discussions on
FinDer
and/or proofreading
and reviewing of the manuscript. This work was funded by
the U.S. Geological Survey and the Gordon and Betty Moore
Foundation.
REFERENCES
Allen, R. M. (2007). The ElarmS earthquake early warning methodology
and its application across California, in
Earthquake Early Warning
Systems
, P. Gasparini, G. Manfredi, and J. Zschau (Editors),
Springer, Berlin, Germany, 21
–
44, ISBN-13 978-3-540-72240-3.
Allen, R. M., P. Gasparini, O. Kamigaichi, and M. Böse (2009). The sta-
tus of earthquake warning around the world: An introductory over-
view,
Seismol. Res. Lett.
80,
no. 5, 682
–
693.
Behr, Y., J. F. Clinton, P. Kästli, C. Cauzzi, R. Racine, and M.-A. Meier
(2015). Anatomy of an earthquake early warning (EEW) alert: Pre-
dicting time delays for an end-to-end EEWsystem,
Seismol. Res. Lett.
86,
no. 3, 830
–
840, doi:
10.1785/0220140179
.
Blaser, L., F. Krüger, M. Ohrnberger, and F. Scherbaum (2010). Scaling
relations of earthquake source parameter estimates with special fo-
cus on subduction environment,
Bull. Seismol. Soc. Am.
100,
2914
–
2926.
Böse, M., and T. H. Heaton (2010). Probabilistic prediction of
rupture length, slip and seismic ground motions for an ongoing
rupture: Implications for early warning for large earthquakes,
Geophys. J. Int.
183,
no. 2, 1014
–
1030, doi:
10.1111/j.1365-
246X.2010.04774.x
.
Böse, M., R. Allen, H. Brown, G. Cua, M. Fischer, E. Hauksson, T.
Heaton,M.Hellweg,M.Liukis,D.Neuhauser,
et al.
(2013). CISN
ShakeAlert
—
An earthquake early warning demonstration system for
California, in
Early Warning for Geological Disasters
—
Scientific Methods
and Current Practice
,F.WenzelandJ.Zschau(Editors),Springer,Berlin,
Germany, ISBN: 978-3-642-12232-3.
Böse, M., R. Graves, D. Gill, S. Callaghan, and P. Maechling (2014).
CyberShake-derived ground-motion prediction models for the Los
1702 Seismological Research Letters Volume 86, Number 6 November/December 2015