Article Volume 14, Number 7 8 July 2013 doi: 10.1002/ggge.20144 ISSN: 1525-2027
A new seismogeodetic approach applied to GPS and accelerometer observations of the 2012 Brawley seismic swarm: Implications for earthquake early warning Jianghui Geng, Yehuda Bock, Diego Melgar, Brendan W. Crowell, and Jennifer S. Haase Cecil H. and Ida M. Green Institute of Geophysics and Planetary Physics, Scripps Institution of Oceanography, University of California, San Diego, La Jolla, California, USA ([email protected]
 The 26 August 2012 Brawley seismic swarm of hundreds of events ranging from M1.4 to M5.5 in the Salton Trough, California provides a unique data set to investigate a new seismogeodetic approach that combines Global Positioning System (GPS) and accelerometer observations to estimate displacement and velocity waveforms. First in simulated real-time mode, we analyzed 1–5 Hz GPS data collected by 17 stations fully encircling the swarm zone at nearsource distances up to about 40 km using precise point positioning with ambiguity resolution (PPPAR). We used a reference network of North American GPS stations well outside the region of deformation to estimate fractional-cycle biases and satellite clock parameters, which were then combined with ultrarapid orbits from the International GNSS Service to estimate positions during the Brawley seismic swarm. Next, we estimated seismogeodetic displacements and velocities from GPS phase and pseudorange observations and 100–200 Hz accelerations collected at three pairs of GPS and seismic stations in close proximity using a new tightly coupled Kalman ﬁlter approach as an extension of the PPP-AR process. We can clearly discern body waves in the velocity waveforms, including P-wave arrivals not detectable with the GPS-only approach for earthquake magnitudes as low as Mw 4.6 and signiﬁcant static offsets for magnitudes as low as Mw 5.4. Our study shows that GPS networks upgraded with strong motion accelerometers can provide new information for improved understanding of the earthquake rupture process and be of critical value in creating a robust early warning system for any earthquake of societal signiﬁcance. Components: 10,952 words, 8 ﬁgures, 4 tables. Keywords: GPS geodesy; seismology; seismogeodesy; tightly coupled Kalman ﬁlter; earthquake early warning; precise point positioning. Index Term: 1200 Geodesy and Gravity: Seismic cycle related deformations, Instruments and techniques. Received 11 February 2013; Revised 8 April 2013; Accepted 10 April 2013; Published 8 July 2013. Geng, J., Y. Bock, D. Melgar, B. W. Crowell, and J. S. Haase (2013), A new seismogeodetic approach applied to GPS and accelerometer observations of the 2012 Brawley seismic swarm: Implications for earthquake early warning, Geochem. Geophys. Geosyst., 14, 2124–2142, doi:10.1002/ggge.20144.
© 2013. American Geophysical Union. All Rights Reserved.
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
1. Introduction  Global Positioning System (GPS) networks are able to observe crustal deformation throughout the entire earthquake cycle from slow interseismic slip to strong coseismic motions. For a large seismic event, high-rate GPS can provide rapid estimates of broadband displacements, including static offsets and dynamic motions of arbitrarily large magnitude [e.g., Nikolaidis et al., 2001; Larson et al., 2003; Bock et al., 2004; Larson, 2009]. High-rate GPS-derived displacements can quickly estimate earthquake magnitude for tsunami warnings [Blewitt et al., 2006], produce centroid moment tensor solutions [Melgar et al., 2012], model ﬁnite fault slip [Crowell et al., 2012; Ohta et al., 2012; Wright et al., 2012], and track seismic wave ﬁelds [Grapenthin and Freymueller, 2011].  It is recognized that high-rate GPS can also play an important role in earthquake early warning (EEW) by providing estimates of permanent displacement within minutes of initiation [e.g., Crowell et al., 2009; Allen et al., 2011]. This is especially valuable close to the source for large (>M7) events where broadband seismometers clip and accelerometer data cannot be objectively integrated to produce reliable displacements in real time [Boore and Bommer, 2005; Emore et al., 2007; Melgar et al., 2013]. Typical EEW systems [e.g., Gasparini et al., 2007; Allen et al., 2009b] depend on conventional seismic instruments, and employ P-wave detection to predict the arrival and intensity of destructive S and surface waves [Heaton, 1985; Nakamura, 1988; Allen and Kanamori, 2003]. However, algorithms only based on seismic data tend to saturate; it is difﬁcult to distinguish an event of magnitude 7 from a larger magnitude of 8 or 9 [Wu and Zhao, 2006; Brown et al., 2009, 2011]. Although GPS excels in providing critical estimates of static offsets, GPS-derived dynamic motions by themselves are not accurate enough to identify millimeter-level or even smaller amplitude P-waves. Furthermore, P-wave arrivals have most of their energy in the vertical direction, making it more difﬁcult for GPS because of the signiﬁcantly less precise vertical component [e.g., Bock et al., 2000]. To be able to detect P-wave arrivals and address the problem of magnitude saturation, Bock et al.  applied a multirate Kalman ﬁlter [Smyth and Wu, 2006] to combine high-rate GPS (1–50 Hz) displacements and accelerometer (100–200 Hz) data in near real time to estimate 3-D seismogeodetic waveforms with millimeterlevel or better precision in displacement and 1
mm/s or better in seismic velocity. They demonstrated the capability of measuring P-wave arrivals for near-source stations deployed in southern California during the 2010 Mw 7.2 El Mayor-Cucapah earthquake in northern Baja California. Thus, seismogeodesy improves on both seismic-only and GPS-only methods, by providing the full spectrum of seismic motions from the detection of P-wave arrivals to the estimation of static displacements.  Because of the relatively small magnitudes of the earthquakes and the excellent distribution of nearby GPS stations, the 26 August 2012 Brawley swarm provides a unique data set to examine the lower bound on the sensitivity of seismic waveforms estimated from high-rate GPS with and without strong-motion accelerometer data. Several hundred events were recorded by the California Integrated Seismic Network (http:// www.cisn.org/). The swarm started at about 15:30 UTC with six events of M < 2.0 and three M2.5 events occurring within a few minutes. The largest events and the focus of this study (Table 1) occurred at 19:20.04.5 UTC (Mw ¼ 4.6) (‘‘Event 1’’), 19:31:22.9 (Mw ¼ 5.4) followed immediately at 19:33:00.8 (Mw ¼ 4.9) (‘‘Event 2––doublet’’), 20:57:58.2 (Mw ¼ 5.5) (‘‘Event 3’’), and 23:33:25.1 (Mw ¼ 4.6) (‘‘Event 4’’) (USGS/NEIC PDE-Q database, htp://earthquake.usgs.gov/earthquakes/eqarchives). The earthquakes occurred on a northeast striking fault zone located about 6 km north of the northwest end of the Imperial fault [http://www.scsn.org/2012Brawley.html], an area that has a history of seismic swarms including one earthquake with a maximum magnitude of M5.1 in 2005 at Obsidian Buttes [Lohman and McGuire, 2007], and another in June 2008 in the same general location as the 2012 event. Chen and Shearer  summarized the history of swarm data in the Imperial Valley and characterized the migration patterns and earthquake mechanisms since 1981. Initial 3-D earthquake relocations and moment tensor inversions [Hauksson et al., 2012] using the double-difference method [Hauksson and Shearer, 2005] indicated a left-lateral strikeslip motion on near vertical fault planes for Events Table 1. Location and Magnitude of the Events Analyzeda Event
1 4.6 2—doublet 5.4 4.9 3 5.5 4 4.6
33.0380 33.0193 33.0210 33.0243 33.0390
Longitude 115.5553 115.5632 115.5540 115.5495 115.5260
Depth (km) Time (UTC) 5 5 14.5 9 11
19:20:04.5 19:31:22.9 19:33:00.8 20:57:58.2 23:33:25.1
Source is the SCEC database (http://www.data.scec.org).
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
Figure 1. Brawley seismic swarm (red dots) surrounded by real-time continuous GPS stations (blue diamonds) and available continuous strong motion stations (open yellow circles) for the period of the swarm. The focal mechanisms are those computed by the SCEC for the four events considered in this study (Table 1). Coseismic displacements and 95% conﬁdence ellipses are from the 24 h SOPAC/JPL combination.
1–3, with a normal faulting component for Event 2, and predominantly normal motion for Event 4 (Figure 1).  We present a novel seismogeodetic analysis method to analyze the earthquake swarm, based on precise point positioning with ambiguity resolution (PPP-AR) [Ge et al., 2008; Geng et al., 2012], supplemented by accelerometer data at locations where GPS receivers and strong motion accelerometers are in close proximity. Precise point positioning methods are attractive because current relative network positioning approaches become more cumbersome as the number of GPS stations to be processed increases from hundreds to thousands at active plate boundaries. The seismogeodetic approach applies a tightly coupled Kalman ﬁlter to GPS and accelerometer data at the observation level as an extension of the GPS-only
PPP-AR approach to estimate displacement and velocity waveforms. This one-step approach differs from the two-step approach presented by Bock et al.  to apply a multirate Kalman ﬁlter to previously estimated GPS displacements and raw accelerometer data.  We discuss the applicability of seismogeodetic PPP-AR to EEW systems, both in the determination of static offsets and the detection of P-wave arrivals, in light of our analysis of the 2012 Brawley seismic swarm data.
2. Theory  In this section we describe the technical details of the seismogeodetic PPP-AR approach and summarize its advantages compared to current GPSonly and loosely coupled seismogeodetic methods. 2126
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
2.1. Ambiguity Resolution for a Single GPS Station  The two main approaches to GPS analysis can be classiﬁed as relative network positioning [e.g., Dong and Bock, 1989; Blewitt, 1989] and precise point positioning [Zumberge et al., 1997]. In either case, observations to GPS satellites from a ground receiver consist of phase and pseudorange measurements at two radio frequencies (‘‘L1’’ at 1575.42 MHz and ‘‘L2’’ at 1227.60 MHz). The phase measurements have integer-cycle phase ambiguities, which are the total number of cycles from the satellite to the receiver. It is critical for precise real-time positioning applications to be able to resolve the phase ambiguities [e.g., Bock et al., 2000].  Baseline vectors between stations in a network are estimated as part of relative network positioning. Common-mode errors due to clock and hardware biases completely cancel in doubly differenced (between satellites and between stations) phase observations, thereby revealing the integer nature of the phase ambiguities. Likewise, common-mode atmospheric errors due to tropospheric and ionospheric refraction are reduced as the distances between stations shorten. Absolute position estimates are then derived by ﬁxing (or tightly constraining) the true-of-date coordinates of one or more reference stations within the network to precise a priori values with respect to a global terrestrial reference frame.  Ambiguity resolution for a single GPS station is difﬁcult because undifferenced phase ambiguity estimates contain noninteger biases, which originate in receiver and satellite hardware. To recover the integer properties of undifferenced ambiguities in PPP analysis, the fractional cycle parts of the noninteger biases can be estimated using a network of reference stations [e.g., Ge et al. 2008; Geng et al., 2012] outside the region of active deformation. Then the fractional-cycle bias (FCB) estimates allow clients to attempt ambiguity resolution for a single station using the PPP-AR method [e.g., Geng et al., 2011]. This approach also requires precise satellite clock and orbit information. In the following section, we describe a tightly coupled Kalman ﬁlter that extends PPP-AR analysis by adding very high-rate strong motion accelerometer data.
2.2. Tightly Coupled Seismogeodetic Filter  For seismogeodesy, we have available very high-rate collocated accelerometer data in addition
to high-rate GPS data. Bock et al.  presented a loosely coupled multirate Kalman ﬁlter that optimally combines GPS displacement and accelerometer data in a two-step process that can be implemented in real time. In the ﬁrst step, the GPS phase and pseudorange data are analyzed to estimate station displacements––this can be done using either relative network positioning or PPPAR. In the second step, the GPS displacements are combined with the accelerometer data. There is no feedback between the two steps so we call this a loosely coupled Kalman ﬁlter. Here we present a tightly coupled Kalman ﬁlter that operates on the raw GPS and accelerometer data in a single step. This formulation is applicable to the estimation of FCBs as part of the PPP analysis of the reference GPS network outside of the zone of deformation, as well as for individual PPP-AR clients within the seismically active region.  Without loss of generality, we assume that the integer-cycle ‘‘wide-lane’’ (the difference between L1 and L2) ambiguities have been resolved through a linear combination of phase and pseudorange measurements [e.g., Teunissen and Kleusberg, 1998]. Since the wide-lane wavelength is about 86.2 cm compared to the 19 and 24 cm wavelengths of the L1 and L2 phases, respectively, this is straightforward even for reference networks of global extent. The wide-lane ambiguities are then applied to the analysis of ionosphere-free carrier phase observations (which have noninteger ambiguities), leaving integer-cycle ‘‘narrow-lane’’ ambiguities with a wavelength of 10.7 cm [e.g., Dong and Bock, 1989]. Then, the narrowlane carrier-phase measurement for reference station i (i ¼ 1, . . . ,r) to satellite j (j ¼ 1, . . . ,s) is given by Lji ¼ ji þ cti þ Mij Ti þ Nij þ Bi Bj þ "ji ;
where ji denotes the geometric distance between station i and satellite j; c is the speed of light in vacuum, ti is the receiver clock error, Ti is the zenith tropospheric delay and Mij is the mapping function, ¼ c=ðf1 þ f2 Þ is the narrow-lane wavelength where f1 and f2 are L1 and L2 frequencies, respectively, Nij denotes the narrow-lane integer ambiguity, Bi and Bj denote the fractional part of noninteger receiver- and satellite-speciﬁc hardj ware biases, respectively, and2 "i denotes the ranj dom error with "i N 0; L . Multipath effects, higher-order ionospheric effects, etc. are ignored for brevity. As part of the GPS analysis we obtain a real-valued estimate for Nij þ Bi Bj ; Bi is 2127
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
eliminated by differencing between satellites leaving the integer part Nij and the fractional part Bj to be estimated. Once the Bj estimates are collected from the r reference stations, we can derive the FCB estimate for satellite j (j ¼ 1, . . . ,s) by
B ji =r;
have achieved the PPP-AR solution for an individual station (client).  For a tightly coupled GPS/accelerometer solution, after correcting the raw accelerometer data for gain, we model the accelerometer data at station l and epoch k as Akl ¼ akl þ bkl þ "kl ;
where Bj ji is the estimate of Bj at station i. Increasing the number of stations in the reference network will improve the reliability and accuracy of Bj [Geng et al., 2012]. Because the result is averaged over many stations, it is not necessary for the stations to be the same as the PPP client stations, only that the analysis in the PPP clients be consistent with the larger reference network. j
 The FCB products for all satellites B (j ¼ 1, . . . ,s) determined in the ﬁrst step are then distributed to the individual PPP clients, where single-station ambiguity resolution can be attempted. Similar to equation (1), at a particular client l (to distinguish from subscript i used for reference stations) at epoch k, we linearize the narrow-lane carrier-phase observation Ljkl ¼ Djkl pkl þ ctkl þ Mklj Tkl þ Nlj þ Bl Bj þ "jl ;
where Akl is the corrected accelerometer measurement, akl is the true acceleration, bkl is an acceleration bias which is estimated as a random walk parameter to accommodate slowly time-varying changes, especially during earthquakes when the bias can change signiﬁcantly, and "kl is the random error with "kl N 0; 2A . Other errors due to, for example, instrument tilts are reasonably presumed to be minimal during the 26 August 2012 Brawley swarm and ignored in this study. Equation (5) is combined with equation (3) for the measurement update. The transition equation for the state vector in equation (4) takes the form of 2
3 2 xkl 1 4 vkl 5 ¼ 4 0 1 akl 0 0
3 32 xk1;l 2 =2 54 vk1;l 5; ak1;l 1
where is the sampling interval of the accelerometer data. Equation (6) is used for the Kalman ﬁlter time update.
2.3. Advantages of the Tightly Coupled Kalman Filter
and 8 T > < pkl ¼"½ xkl vkl akl # j xkl xk ; j > 0 0 : Dkl ¼ jkl
where denotes the increment of a parameter estimate, pkl is the state vector comprising the position increment xkl , the velocity vkl and the acceleration akl , and Djkl is the design matrix in which xkl and xjk are the station and satellite posi tions, respectively, and "jl N 0; 2L . The satelj lite bias Bj is assigned the value B , and the receiver Bl is assimilated into the receiver clock estimate, or equivalently can be removed by between-satellite differencing. As a result, the estimate for Nlj þ Bl Bj is reduced to an estimate for Nlj where the integer property has been recovered. Ambiguity resolution for a single station can then be attempted. Once successfully ﬁxing Nlj to integers for all visible satellites, we
 The key difference between the tightly coupled Kalman ﬁlter presented here and the loosely coupled ﬁlter in Bock et al.  is that in the tightly coupled case the accelerations are used as additional data to resolve ambiguities. The accelerometer data are applied as tight constraints on the position variation between epochs. This single-step process improves cycle-slip repair for GPS carrier-phase data and rapid ambiguity resolution after GPS outages [Grejner-Brzezinska et al., 1998]. This is conﬁrmed with the Brawley swarm data as discussed in the supporting information (section 1 and Figure S1).1  As with the loosely coupled ﬁlter the tightly coupled ﬁlter also minimizes step functions possibly introduced by tilt in the accelerometer observations. Furthermore, biases, bkl, in the 1 Additional supporting information may be found in the online version of this article.
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
accelerometer data are estimated along with other parameters of interest. Hence, no pre-event mean needs to be eliminated from the acceleration data before starting the Kalman ﬁlter. In addition, bkl is estimated at each epoch as a random walk parameter to mitigate possible drift in the accelerometer data due to translations (indistinguishable from rotations) [Trifunac and Todorovska, 2001] and temperature changes.
2.4. Advantages of GPS PPP-AR Over Relative Network Positioning  There are several advantages to the PPP-AR approach. First, it is highly efﬁcient because the satellite clocks and FCBs are estimated only once for positioning any number of clients. The data processing at each client is independent and does not affect any other clients. Relative network positioning is complicated by the need to assign baselines, overlapping Delaunay triangles [Crowell et al., 2009], or overlapping subnetworks [Bock et al., 2011]. This is a critical difference as one is faced with the challenge of analyzing hundreds to thousands of stations in real time. Furthermore, intermittent station dropouts complicate relative network positioning. Therefore, PPP-AR can be efﬁciently applied to large GPS networks deployed over a wide area such as around the Circum-Paciﬁc Seismic Belt or at isolated stations in remote areas. Second, the PPP-AR approach compared to
relative network positioning does not require a local reference station, which might be displaced during a large event. Instead, PPP-AR requires a continental- or global-scale reference network well outside the zone of expected deformation, yet still has the same satellites visible as the client stations. Finally, if undifferenced ambiguities can be successfully ﬁxed the positioning accuracy is comparable to that of relative network positioning [e.g., Bertiger et al., 2010; Geng et al., 2010].
3. Data Analysis 3.1. An Operational Real-Time PPP-AR System at SOPAC  At the time of the 2012 Brawley swarm and as part of a prototype EEW system for the Western U.S., we had already implemented at the Scripps Orbit and Permanent Array Center (SOPAC) an operational PPP-AR service center for estimating satellite clocks and FCBs. These parameters are intended for distribution to PPP clients in real time whether at a centralized processing facility, at a remote computer with internet access to data from a speciﬁc station, at a local processor at the remote station or within the GPS receiver itself. With predicted ultrarapid satellite orbits from the international GNSS service (IGS), the SOPAC service center generates satellite clocks every second for
Figure 2. Distribution of high-rate GPS stations used by the SOPAC PPP-AR system. Solid green circles denote the 46 stations used for clock estimation whereas solid red triangles denote the 48 stations used for FCB determination. The solid black star shows the location of the 26 August 2012 Brawley seismic swarm. 2129
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
each visible GPS satellite using 46 reference stations, and FCBs every 5 s using 48 stations across North America (Figure 2—Some of the stations are overlapping). The 1 Hz reference station data are collected from IGS and UNAVCO’s Plate Boundary Observatory (UNAVCO/PBO) servers. The reference network for satellite clock estimation is chosen to be of continental scale to reduce errors; another option is to use satellite clocks from existing sources based on a global distribution of stations. On the other hand, the FCB network is chosen to be as close to the area of interest as possible, while staying sufﬁciently outside the region of expected deformation. The reference stations for the SOPAC service are chosen to be further than 200 km from the primary zones of tectonic deformation in California, Oregon and Washington to avoid contamination of the satellite clocks and FCBs during a large seismic event in that region. For the FCB determination, we choose a reference network that is close to the western U.S. coast (rather than a single reference station that is chosen for relative network positioning). Should one of the FCB network stations be subject to dynamic motions, the impact would be minimized by the averaging in equation (2). We note that for the western U.S. the FCB determination is challenged by the sparseness of real-time off-shore reference stations in the Paciﬁc. The reference station positions are ﬁxed to true-of-date estimates with respect to ITRF2008 [Altamimi et al., 2011] produced through an (1–2 week) extrapolation of combined modeled time series based on SOPAC’s and Jet Propulsion Laboratory’s (JPL)’s routine weekly analysis of 24 h, 30 s sampled data from a global and regional set of continuous GPS stations (http://sopac.ucsd.edu/processing/coordinates/sector.shtml) and made available through the GPS Explorer data portal (http://geoapp.ucsd.edu/). Real-time SOPAC PPP-AR results for the highrate station GLRS in southern California (Figure 1) from 11 August to 5 September 2012 show that ambiguity resolution in PPP improves the rootmean-square (RMS) difference between the position estimates and ground truth from 20, 30, and 61 mm without ambiguity resolution to 15, 12, and 40 mm for the North, East, and Up components, respectively.
3.2. GPS PPP-AR Analysis  For the 26 August 2012 Brawley earthquake swarm, we used SOPAC’s PPP-AR system to process high-rate GPS data from 16 PBO stations and 1 Southern California Integrated GPS Network
(SCIGN) station in this region (Figure 1). The stations stream 1 Hz GPS data, the normal operational setting for most real-time GPS networks. For this study, we requested after the fact that UNAVCO/PBO download 5 Hz data from the receiver buffers at the 16 PBO stations. We then processed the 5 Hz data in a simulated real-time mode using SOPAC’s 1 Hz satellite clocks, 5 s FCBs, and the predicted IGS ultrarapid orbits published for that period. The 1 Hz data at USGS station BOMG were also processed in this way. Further details on PPP-AR estimation are provided in section 2 of the supporting information.  In practice, the transmission of high-rate GPS data from the reference network (or a client station) to server has a typical latency of 0.4–1.0 s. The satellite clocks and FCB parameters are continuously estimated and made immediately available to clients. The PPP-AR processing at individual stations (clients) can then be performed at each epoch with a delay about 1 s after data arrival at the data analysis center.
3.3. Tightly Coupled Kalman Filter Analysis  There are only a few collocated GPS and seismic stations in southern California (Figure 1). For this study we identiﬁed three suitable collocations. The ﬁrst pair consists of GPS station P506 sampling at 5 Hz and accelerometer site WLA sampling at 200 Hz, separated by 2.6 km and approximately 8 km from the source. This site pair was the only collocation available within the Brawley Seismic Zone. Emore et al.  demonstrated good agreement of 1 Hz GPS and strong-motion data with instrument separations of up to 4 km. The second pair is GPS station P494 sampling at 5 Hz and strong motion site WES sampling at 100 Hz, separated by about 80 m and 34 km from the largest event. The third and most closely spaced pair (within 10 m) is BOMG (1 Hz) and strong motion site BOM (200 Hz), about 42 km from the swarm events. The strong motion sensors are observatory grade EpiSensor accelerometers and the broadband sensors are STS2 instruments, both on 24 bit Quanterra data loggers. For each pair, we combined the 1–5 Hz PPP-AR derived displacements with 100–200 Hz acceleration measurements. We also considered pair P493/ NP.286 (Figure 1) but the strong-motion site NP.286 is located inside the second story of a building and its record is complicated by the building response. 2130
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
Table 2. Standard Deviation (1, mm) of the Differences Between the 1–5 Hz GPS-Derived Ambiguity Fixed PPP and Instantaneous Position Estimates, and the True-of-Date SOPAC/JPL Positions Over 5 h Preceding the Brawley Swarm Instantaneous Relative Network Positioning
P506 P499 P495 P502 P498 P503 P507 P497 P501 P744 P510 P493 P508 P509 P496 P494 BOMG Mean
8 8 8 9 9 9 9 7 10 6 11 7 9 9 9 8 11 9
7 6 8 6 5 6 7 5 6 4 8 4 7 7 5 5 10 6
40 35 34 25 22 21 37 24 29 30 33 19 25 30 24 25 41 29
13 13 13 13 12 13 15 12 13 12 14 12 14 13 13 13 N/A 13
16 15 15 15 15 14 16 15 14 14 17 14 16 17 16 15
54 53 53 54 58 53 51 56 56 58 53 48 56 57 65 65
 As a ﬁrst step the raw accelerometer data were corrected for gain. In the tightly coupled Kalman ﬁltering, the precision of the raw accelerometer data from site WLA was taken to be 10 mm/s2 for all three channels to account for its relatively long distance from station P506 (alternatively, we could have applied a time shift to the accelerometer data). In contrast, we applied 0.1 mm/s2 for sites WES and BOM, which are closer to their GPS counterparts. The process noise of the accelerometer bias parameters is presumed to be 0.001 mm/s2.5. As a comparison, in the loosely coupled approach of Bock et al.  the pre-event means are removed from the accelerometer data. Then, the two ﬁlter parameters for the periods of shaking, the system variance q and the measurement variance r, are determined from 60 s of pre-event noise on each channel of the accelerometer and GPS, respectively. Also, a near-real-time smoother with 5 s lag is applied to the loosely coupled ﬁlter.
4. Results 4.1. PPP-AR Analysis Using Only GPS Data  We evaluated the PPP-AR approach, ﬁrst using only GPS data in the context of real-time analysis for EEW. The results are of interest in the case that only GPS data are available in real time
near the source of a large earthquake. Furthermore, the GPS-only approach can be used for rapid and reliable estimates of earthquake magnitude, fault geometry and rupture characteristics for large events often faster than traditional seismic methods on their own [e.g., Melgar et al., 2012; Crowell et al., 2012; Ohta et al., 2012; Wright et al., 2012].  The RMS of the differences between the 1 or 5 Hz GPS-derived ambiguity-ﬁxed position estimates and the true-of-date SOPAC/JPL positions are 31, 14, and 32 mm on average while the standard deviations (1) are 9, 6, and 29 mm in the North, East and Up directions, respectively, based on 5 h of data for all sites before the start of Brawley swarm activity. In contrast, the RMS for the ambiguity-ﬂoat positions are 29, 20, and 57 mm, and the standard deviations (1) are 9, 17, and 50 mm, respectively. This is consistent with earlier results that indicate that ambiguity resolution improves the positioning accuracy of real-time PPP, especially in East and Up components [e.g., Geng et al., 2011]. In addition, for the 5 h data span, we compared PPP-AR to instantaneous relative positioning [Bock et al., 2000, 2011] (Table 2). We were able to assign for the latter a reference station just outside the zone of deformation because the earthquake magnitudes were so small. We can see that PPP-AR outperforms instantaneous relative positioning by 31%, 60%, and 48% in the standard deviations for the North, East and Up components, respectively. This improvement is attributed to the forward (smoothing) Kalman ﬁlter used in PPP-AR, whereas instantaneous relative positioning processes each epoch of data independently [Bock et al., 2000], and to improved single-epoch ambiguity resolution in the PPP-AR process. Other relative positioning methods that use Kalman ﬁlter estimation should provide the same level of precision as PPP-AR [e.g., Bertiger et al., 2010; Geng et al., 2010], but are more complicated to apply in real-time scenarios. In the remainder of this paper, we take the individual standard deviation (1) from the real-time PPPAR processing as the precision of the GPS position time series (Table 2).  We show the displacement waveforms for the 17 GPS stations during all four events in the supporting information (Figures S2–S5). From the horizontal components of the three stations P506, P499, and P495 that are nearest to the epicenter ( 3 earthquake at near-ﬁeld distances (10 km) and M > 5 at regional distances (100 km) (Robert Clayton, personal communication). We can infer that these types of accelerometers directly collocated at the GPS stations would have been adequate to monitor the larger events of the 2012 Brawley swarms.
6. Conclusions  We have established that single-epoch GPSonly PPP-AR is a highly efﬁcient and viable
method to analyze high-rate GPS data, with the same precision as relative network positioning methods but without the requirement to choose a local reference station. Rather, positions are independently estimated with respect to an external reference frame realized through the true-of-date global coordinates of the reference network stations and estimated satellite orbits and clocks. FCBs estimated from the reference network allow for ambiguity resolution at each station within the zone of deformation. This methodology is particularly useful for continuous GPS networks of hundreds to thousands of stations deployed within zones of tectonic deformation and seismic risk because the computation time scales linearly as the number of stations increases. However, like relative network positioning, it is limited in precision to about 20 mm for detection of real-time coseismic displacements and is unable to detect Pwave arrivals.  We have shown the lower bound of sensitivity for GPS-only and combined GPS/accelerometer seismogeodetic waveforms using data collected from four near-source events of the 26 August 2012 Brawley seismic swarm, ranging in magnitude from Mw 4.6 to 5.5. During periods of seismic shaking the tightly coupled seismogeodetic PPPAR analysis provides continuous unclipped estimates of broadband displacements in all three components that clearly discern the arrival of Pwaves. Seismogeodetic velocities are estimated throughout the period of seismic shaking with a precision of about 0.1–1 mm/s. Together, seismogeodetic broadband displacements and velocities from near-source stations can form the basis for EEW systems for any magnitude earthquake of consequence, which is more robust against false alarms than seismic-only or GPS-only approaches. Furthermore, seismogeodesy can be used to better understand the physics of the earthquake rupture process.
Acknowledgments  We thank the two reviewers and the Associate Editor for providing useful comments. Melinda Squibb and Anne Sullivan handled the real-time GPS data. Peng Fang was responsible for the 24 h GPS analysis at SOPAC, and Angelyn Moore for the analysis at JPL. GPS 5 Hz data for the 2012 Brawley swarm were provided by the Plate Boundary Observatory (PBO) operated by UNAVCO for EarthScope (http:// www.earthscope.org) and supported by NSF grant EAR0323309. Real-time 1 Hz GPS data over North America were collected from the International GNSS Service and PBO. 2140
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
Accelerometer data were obtained from the California Integrated Seismic Network, and University of California Santa Barbara (courtesy of Jamie Steidl). Seismic event information is from the Southern California Earthquake Center (SCEC) data archive. This paper was funded by NASA grants AIST11 NNX09AI67G, ROSES NNX12AK24G, 06-MEaSUREs06-0085, and SCEC award 12083.
References Allen, R. M., and H. Kanamori (2003), The potential for earthquake early warning in southern California, Science, 300, 786–789. Allen, R. M., and A. Ziv (2011), Application of real-time GPS to earthquake early warning, Geophys. Res. Lett., 38, L16310, doi:10.1029/2011GL047947. Allen, R. M., P. Gasparini, O. Kamigaichi, and M. Bose (2009a), The Status of earthquake early Warning around the world: An introductory overview, Seismol. Res. Lett., 80, 682–693, doi:10.1785/gssrl.80.5.682. Allen, R. M., P. Gasparini, and O. Kamigaichi (Eds.) (2009b), Earthquake early warning, special issue, Seismol. Res. Lett., 80(5), 682–782. Altamimi, Z., X. Collilieux, and L. Metivier (2011), ITRF2008: An improved solution of the international terrestrial reference frame, J. Geod., 85(8), 457–473. Bertiger, W., S. D. Desai, B. Haines, N. Harvey, A. W. Moore, S. Owen, and J. P. Weiss (2010), Single receiver phase ambiguity resolution with GPS data, J. Geod., 84(5), 327–337. Blewitt, G. (1989). Carrier phase ambiguity resolution for the global positioning system applied to geodetic baselines up to 2000 km, J. Geophys. Res., 94, 10,187–10,203. Blewitt, G., C. Kreemer, W. C. Hammond, H.-P. Plag, S. Stein, and E. Okal (2006), Rapid determination of earthquake magnitude using GPS for tsunami warning systems, Geophys. Res. Lett., 33, L11309, doi:10.1029/2006G L026145. Bock, Y., R. Nikolaidis, P. J. de Jonge, and M. Bevis (2000), Instantaneous geodetic positioning at medium distances with the global positioning system, J. Geophys. Res., 105, 28,223–28,254. Bock, Y., L. Prawirodirdjo, and T. I. Melbourne (2004), Detection of arbitrarily large dynamic ground motions with a dense high-rate GPS network, Geophys. Res. Lett., 31, L06604, doi:10.1029/2003GL019150. Bock, Y., D. Melgar, and B. W. Crowell (2011), Real-time strong-motion broadband displacements from collocated GPS and accelerometers, Bull. Seismol. Soc. Am., 101(5), 2904–2925, doi:10.1785/0120110007. Boore, D. M., and J. J. Bommer (2005), Processing of strongmotion accelerograms : Needs, options, and consequences, Soil Dyn. Earthquake Eng., 25, 93–115. Böse, M., T. Heaton, and E. Haukkson (2012), Rapid estimation of earthquake source and ground-motion parameters for earthquake early warning using data from a single threecomponent broadband or accelerometer sensor, Bull. Seismol. Soc. Am., 102(2), 738–750, doi:10.1785/0120 110152. Brown, H., R. M. Allen, and V. Grasso (2009), Testing ElarmS in Japan, Seismol. Res. Lett., 80, 727–739. Brown, H. M., R. M. Allen, P. Hellweg, O. Khainovski, D. Neuhauser, and A. Souf (2011), Development of the ElarmS methodology for earthquake early warning: Realtime appli-
cation in California and ofﬂine testing in Japan, Soil Dyn. Earthquake Eng., 31, 188–200. Chen, X., and P. M. Shearer (2011), Comprehensive analysis of earthquake source spectra and swarms in the Salton Trough, California, J. Geophys. Res., 116, B09309, doi:10.1029/2011JB008263. Choi, K., A. Bilich, K. Larson, and P. Axelrad (2004), Modiﬁed sidereal ﬁltering: Implications for high-rate GPS positioning, Geophys. Res. Lett., 31, L22608, doi:10.1029/ 2004GL021621. Crowell B. W. (2013), Using GPS to rapidly detect and model earthquakes and transient deformation events, PhD thesis, Univ. of Calif., San Diego, Calif. Crowell, B. W., Y. Bock, and M. B. Squibb (2009), Demonstration of earthquake early warning using total displacement waveforms from real-time GPS networks, Seismol. Res. Lett., 80(5), 772–782, doi:10.1785/gssrl.80.5.772. Crowell, B. W., Y. Bock, and D. Melgar (2012), Real-time inversion of GPS data for ﬁnite fault modeling and rapid hazard assessment, Geophys. Res. Lett., 39, L09305, doi:10.1029/2012GL051318. Dong, D. N., and Y. Bock (1989), Global positioning system network analysis with phase ambiguity resolution applied to crustal deformation studies in California, J. Geophys. Res., 94, 3949–3966. Ellsworth, W. L., and G. C. Beroza (1995), Seismic evidence for an earthquake nucleation phase, Science, 268(5212), 851–855. Emore, G., J. Haase, K. Choi, K. M. Larson, and A. Yamagiwa (2007), Recovering absolute seismic displacements through combined use of 1-Hz GPS and strong motion accelerometers, Bull. Seismol. Soc. Am., 97(2), 357–378, doi:10.1785/ 0120060153. Gasparini, P., G. Manfredi, and J. Zschau (Eds.) (2007), Earthquake Early Warning Systems, Springer, New York. Ge, M., G. Gendt, M. Rothacher, C. Shi, and J. Liu (2008), Resolution of GPS carrier-phase ambiguities in precise point positioning (PPP) with daily observations, J. Geod., 82(7), 389–399. Geng, J., F. N. Teferle, X. Meng, and A. H. Dodson (2010), Kinematic precise point positioning at remote marine platforms. GPS Solutions, 14(4), 343–350. Geng, J., F. N. Teferle, X. Meng, and A. H. Dodson (2011), Towards PPP-RTK: Ambiguity resolution in real-time precise point positioning, Adv. Space Res., 47(10), 1664–1673. Geng, J., C. Shi, M. Ge, A. H. Dodson, Y. Lou, Q. Zhao, and J. Liu (2012), Improving the estimation of fractional-cycle biases for ambiguity resolution in precise point positioning, J. Geod., 86(8), 579–589. Genrich, J. F., and Y. Bock (1992), Rapid resolution of crustal motion at short ranges with the global positioning system, J. Geophys. Res., 97, 3261–3269. Genrich, J. F., and Y. Bock (2006), Instantaneous geodetic positioning with 10–50 Hz GPS measurements: Noise characteristics and implications for monitoring networks, J. Geophys. Res., 111, B03403, doi:10.1029/2005JB003617. Grapenthin, R., and J. Y. Freymueller (2011), The dynamics of a seismic wave ﬁeld: Animation and analysis of kinematic GPS data recorded during the 2011 Tohoku-oki earthquake, Japan, Geophys. Res. Lett., 38, L18308, doi:10.1029/ 2011GL048405. Grejner-Brzezinska, D. A., R. Da, and C. Toth (1998), GPS error modeling and OTF ambiguity resolution for highaccuracy GPS/INS integrated system, J. Geod., 72(11), 626–638. 2141
GENG ET AL : SEISMOGEODETIC EARTHQUAKE EARLY WARNING
Hauksson, E., and P. Shearer (2005), Southern California hypocenter relocation with waveform cross-correlation: Part 1. Results using the double-difference method, Bull. Seismol. Soc. Am., 95, 896–903. Hauksson, E., W. Yang, and P. M. Shearer (2012), Waveform relocated earthquake catalog for southern California (1981 to June 2012), Bull. Seismol. Soc. Am., 102(5), 2239–2244. Heaton, T. H. (1985), A model for a seismic computerized alert network, Science, 228, 987–990. Kanamori, H. (2005), Real-time seismology and earthquake damage mitigation, Annu. Rev. Earth Planet. Sci., 33, 195–214, doi:10.1146/annurev.earth.33.092203.122626. Langbein, J., and Y. Bock (2004), High-rate real-time GPS network at Parkﬁeld: Utility for detecting fault slip and seismic displacements, Geophys. Res. Lett., 31, L15S20, doi:10.1029/2003GL019408. Larson, K. (2009), GPS seismology, J. Geod., 83(3):227–233. Larson, K., P. Bodin, and J. Gomberg (2003), Using 1-Hz GPS data to measure deformations caused by the Denali fault earthquake, Science, 300, 1421–1424. Lohman, R. B., and J. J. McGuire (2007), Earthquake swarms driven by aseismic creep in the Salton Trough, California, J. Geophys. Res., 112, B04405, doi:10.1029/2006JB004596. Melgar, D., Y. Bock, and B. W. Crowell (2012), Real-time centroid moment tensor determination for large earthquakes from local and regional displacement records, Geophys. J. Int., 188(2), 703–718, doi:10.1111/j.1365-246X.2011. 05297.x. Melgar, D., Y. Bock, D. Sanchez, and B. W. Crowell (2013), On robust and reliable automated baseline corrections for strong motion seismology, J. Geophys. Res., 118, doi: 10.1029/2012JB009937. Nakamura, Y. (1988), On the urgent earthquake detection and alarm system (UrEDAS), paper presented a Proceedings of the 9th World Conference on Earthquake Engineering, Int. Assoc.for Earthquake Eng., Tokyo, Kyoto, 2–9 August. Nikolaidis R., Y. Bock, P. J. de Jonge, P. Shearer, D. C. Agnew, and M. Van Domselaar (2001), Seismic wave observations with the global positioning system, J. Geophys. Res., 106, 21,897–21,916. Ohta, Y., et al. (2012), Quasi real-time fault model estimation for near-ﬁeld tsunami forecasting based on RTK-GPS analy-
sis: Application to the 2011 Tohoku-Oki earthquake (9.0 Mw), J. Geophys. Res., 117, B02311, doi:10.1029/2011 JB008750. Olsen, E. L., and R. M. Allen (2005), The deterministic nature of earthquake rupture, Nature, 438, 212–215, doi:10.1038/ nature04214. Olsen, E. L., and R. M. Allen (2006), Olsen and Allen Reply, Replying to: P. Rydelek and S. Horiuchi, Nature, 442, E6, doi:10.1038/nature04963. Rios, J. A., and E. White (2002), Fusion ﬁlter algorithm enhancements for a MEMS GPS/IMU, ION GPS 2002, paper presented at the 2002 National Technical Meeting of The Institute of Navigation, Institute of Navigation, San Diego, Calif. Rydelek, P., and S. Horiuchi (2006), Is earthquake rupture deterministic? Arising from: E. L. Olson and R. M. Allen, Nature, 438, 212–215 (2005), Nature, 442, E5. Smyth, A., and M. Wu (2006), Multirate Kalman ﬁltering for the data fusion of displacement and acceleration response measurements in dynamic system monitoring, Mech. Syst. Signal Process., 21, 706–723. Teunissen, P., and A. Kleusberg (Eds.) (1998), GPS for Geodesy, 2nd ed., 650 pp., Springer, Berlin. Trifunac, M. D., and M. I. Todorovska (2001), A note on the useable dynamic range of accelerographs recording translation, Soil Dyn. Earthquake Eng., 21, 275–286, doi:10.1016/ S0267-7261(01)00014-8. Wright T. J., N. Houlie, M. Hildyard, and T. Iwabuchi (2012), Realtime, reliable magnitudes for large earthquakes from 1 Hz GPS precise point positioning: The 2011 Tohoku-Oki (Japan) earthquake, Geophys. Res. Lett. 39, L12302, doi:10.1029/2012GL051894. Wu, Y.-M., and L. Zhao (2006), Magnitude estimation using the ﬁrst three seconds P-wave amplitude in earthquake early warning, Geophys. Res. Lett. 33, L16312, doi:10.1029/ 2006GL026871. Zumberge, J. F., M. B. Heﬂin, D. C. Jefferson, M. M. Watkins, and F. H. Webb (1997), Precise point positioning for the efﬁcient and robust analysis of GPS data from large networks, J. Geophys. Res., 102(B3), 5005–5017, doi:10.1029/ 96JB03860.