Real-Time Characterization of Finite Rupture and Its Implication for Earthquake Early Warning: Application of FinDer to Existing and Planned Stations in Southwest China

Introduction

Earthquake early warning (EEW) systems quickly detect earthquakes (possibly while still evolving) and alert users (e.g. the public, automated response applications, situational awareness users) prior to the onset of strong ground shaking at a given location (Allen et al., 2009; Satriano et al., 2011; Allen and Melgar, 2019). In the last few decades, EEW systems have been used successfully and are generally regarded as a promising tool for earthquake risk reduction (Strauss and Allen, 2016).

Many countries around the world have operational EEW systems; among them the systems in Japan, Mexico and Taiwan are probably best known and advanced (Espinosa-Aranda et al., 1995; Nakamura and Saita, 2007; Wu et al., 2007). Other countries and regions, such as Romania, Turkey, South Korea, Israel and the West Coast of the United States, have installed EEW systems, and these systems have provided valuable warnings to the public and automated emergency systems of key infrastructures (Clinton et al., 2016; Sheen et al., 2017; Suárez et al., 2018; Kohler et al., 2020; Kurzon et al., 2020). Some regions (e.g. Italy, Chinese mainland, Switzerland, Chile, Nicaragua, Costa Rica and El Salvador) are currently testing EEW systems within their seismic monitoring networks (Clinton et al., 2016; Leyton et al., 2018; Massin et al., 2019; and see next section for details of Chinese mainland), while others are only at the stage of discussing the need and feasibility of developing EEW systems (Bird et al., 2008; Stankiewicz et al., 2015; Romeu Petit et al., 2016; Sokos et al., 2016; Zuccolo et al., 2016; Ogweno et al., 2019). Recently, crowdsourcing (Hammon and Hippner, 2012) has been applied as a promising new technique to EEW, allowing the general public to record seismic data/records with micro-electro-mechanical systems (MEMS) in smart devices and to receive earthquake alerts in real-time on a global or regional scale (Minson et al., 2015). Some well-known examples are the Earthquake Network (Finazzi, 2016), MyShake (Kong et al., 2016), Quake Catcher Network (Cochran et al., 2009), Home Seismometer (Horiuchi et al., 2009), and Google’s EEW effort based on Android phones (Stogaitis et al., 2020).

EEW algorithms, such as the Virtual Seismologist (VS; Cua and Heaton, 2007), ElarmS (Chung et al., 2019) or OnSite (Wu and Kanamori, 2005), have demonstrated their usefulness for providing earthquake source parameters in near real-time. Although these algorithms are fast, they assume that seismic sources are a point in time and space, which is unsuitable for large magnitude earthquakes. Since ignoring the finiteness of earthquake ruptures, the point-source model fails to capture the true temporal-spatial variance of the radiation of high-frequency seismic energy in earthquakes with M > 6, and tends to saturate in earthquake magnitude estimation (Allen and Melgar, 2019). Furthermore, point-source model-based algorithms typically rely on early information from P-wave and/or rupture nucleation to predict the final size of the (evolving) fault rupture, even though earthquake rupture predictability remains controversial (Ide, 2019; Hutchison et al., 2020; Meier et al., 2020). Hence, a multi-dimensional characterization of the source model to represent the complexities of a large earthquake is necessary to improve EEW.

Several geodetic and seismic algorithms have been developed to estimate the extended size of fault ruptures in (near) real-time. For example, the G-larmS (Grapenthin et al., 2014a; 2014b), G-FAST (Crowell et al., 2016, 2018), BEFORES (Minson et al., 2014) and REGARD (Kawamoto et al., 2016, 2017) algorithms use high-rate GPS/GNSS (Global Positioning System/Global Navigation Satellite System) observations for real-time source inversion. Although usually providing an accurate characterization of seismic sources, these algorithms require dense network observations and tend to be computationally time-consuming. Moreover, GPS/GNSS-based algorithms are usually unable to provide earthquake alerts to sites that are close to the epicenter due to the sparsity of the existing monitoring networks (McGuire et al., 2021). A promising alternative is to adopt the seismic Finite-Fault Rupture Detector (FinDer) algorithm (Böse et al., 2012, 2015, 2018), which uses the spatial distribution of seismic high-frequency ground-motions to estimate the centroid, length (L) and strike (Θ) of an assumed line-source.

In this study, we explore the feasibility of operating FinDer in Chinese mainland, where a nationwide EEW projects are currently underway, including the installation of several thousands of low-cost MEMS sensors. We will begin this article with a review of Chinese EEW efforts. Then, using waveform playbacks of both real and simulated waveform data of recent seismicity, we will test FinDer for existing and planned stations, and discuss its performance with a view to possible future operation of FinDer in Chinese mainland.

Development of EEW in Chinese Mainland

In 1994, the first EEW system in Chinese mainland was deployed at the Daya Bay nuclear power plant; since then, Chinese mainland has been continuously developing and evolving EEW systems in order to improve their performance for earthquake risk mitigation (Li et al., 2004).

Top-Level Design

Around the turn of the 21st century, a series of scientific projects related to EEW and fast earthquake intensity reports, organized by the China Earthquake Administration (CEA), were launched (Li, 2014). Mainly due to constructions between 2001 and 2008, the density of broadband seismic and strong-motion networks in Chinese mainland has been greatly increased (Figure 1A), and all stations have been digitized (Zhou, 2006; Liu et al., 2008).

FIGURE 1. (A) Existing (as of 2018) and (B) planned (as of 2022) stations, deployed or upgraded as part of the National System for Fast Report of Intensities and Earthquake Early Warning of Chinese mainland. Red dots show the epicenters (from China Earthquake Networks Center, CENC) of the 2008 M_S 8.0 Wenchuan and 57 5.0 ≤ M_S < 7.0 earthquakes studied in this paper. White lines show the fault locations. (C) Magnitude distribution of earthquakes shown in (A) and (B). (D) The larger tectonic setting with the five key EEW regions: 1) north China, 2) central China north-south seismic belt, 3) southeast coast, 4) middle section of Tianshan Mountains of Xinjiang, and 5) Lhasa of Tibet.

City/Infrastructure Scale

In 1994, a small earthquake alert system with 12 strong-motion stations was deployed on-site around Daya Bay nuclear power plant. This system issued a warning when a trigger threshold of 0.01 g was exceeded. Following this, similar alert systems were installed at more than 10 additional nuclear power plants. In 2001, PetroChina Dalian Branch developed an alert system for petrochemical enterprises. In 2007, an earthquake monitoring and alarm system was developed for the Hebei-Nanjing gas pipeline. In the same year, as part of China's digital strong-motion network, a fast earthquake intensity report system was deployed in some Chinese cities, including Beijing, Tianjin, Lanzhou, Urumqi and Kunming (Zhou, 2006). For the Three Gorges water conservancy project, 15 seismic and 16 mobile stations were deployed to provide EEW and earthquake information service. In addition, the EEW system is functioning for more than 20 lines of high-speed railways such as the Beijing-Tianjin, the Beijing-Shanghai, the Chengdu-Chongqing, the Beijing-Shijiazhuang-Wuhan lines (Li, 2014).

Provincial Scale

To demonstrate the feasibility of EEW on provincial scale, China has implemented prototype test systems in three areas (Sichuan Earthquake Administration, 2015). In 2007, the Institute of Geophysics at CEA and the Department of Geosciences at National Taiwan University built an EEW prototype system for the Beijing Capital Region (Peng et al., 2011). This system was based on the Capital Circle Seismograph Network of China and included 94 broadband and 68 short-period stations, with an average interstation spacing of roughly 50 km. In 2014, six years after the Wenchuan M_S 8.0 earthquake, new system in the region was expanded to 130 broadband seismic and 80 strong-motion stations with an average interstation spacing of around 40 km for whole covered land, and around 15 km in urban area. The seismic network allows data to be transmitted to the data center in real-time. In March 2015, an earthquake fast report system was deployed in this region, which integrated additional data from 100 intensity sensors in the Tangshan region. The second provincial EEW demonstration system deployed in Fujian province was completed in September 2013 and consists of 128 stations with an average interstation spacing of around 30 km. The system accesses data from 16 stations of Taiwan in real-time and provides warning for the earthquakes in Taiwan (Zhang et al., 2016). The Lanzhou EEW demonstration system was completed in October 2014, and connected to 83 broadband seismic and 80 strong-motion stations from the provincial networks of Gansu, Qinghai and Ningxia. Additional EEW demonstration systems were built in the Sichuan-Yunnan region and the east part of Guangdong successively over the last years.

National Scale

In order to improve the capability of EEW, earthquake parameter and fast intensity reports, and earthquake seismological research on a nationwide scale, the National System for Fast Report of Intensities and Earthquake Early Warning project of Chinese mainland, led by the CEA, was launched and implemented in 2015 and 2018, respectively. A nationwide EEW system covering five key EEW zones (north China, central China north-south seismic belt, southeast coast, middle section of Tianshan Mountains of Xinjiang, and Lhasa of Tibet; Figure 1D) is currently being built, with an average interstation spacing of around 10–15 km. It is expected that by 2022, close to around 2,000 broadband stations (equipped with three-component broadband or very broadband seismometers and all with accelerometers), around 3,200 strong-motion stations (equipped with three-component accelerometers), and around 10,200 low-cost intensity sensors (equipped with MEMS; Peng et al., 2017) will be deployed or upgraded throughout Chinese mainland. The EEW system will build on the Jopens seismic monitoring platform developed by CEA. Some areas (e.g. Sichuan-Yunnan and Beijing-Tianjin-Hebei regions) that have completed most of the station deployment have taken the lead in testing EEW algorithms (China Earthquake Administration, 2020; Peng et al., 2020). The system successfully provided early warning services for the 2018 Yongqing M_S 4.3 earthquake of Hebei, the 2019 Hualian M_S 6.7 earthquake of Taiwan, and 2019 M_S 6.0 earthquake of Sichuan. At present, the project is actively supporting and guiding social entities to participate (Song et al., 2021), e.g. the Institute of Care-life at Chengdu.

The Finite-Fault Rupture Detector

The Finite-Fault Rupture Detector (FinDer) is an efficient algorithm to compute a line-source model of an ongoing earthquake fault rupture from real-time high-frequency seismic data using template matching (Böse et al., 2012, 2015, 2018). The resulting model has four independent parameters that characterize the seismic source: length, strike, centroid, and rupture directivity. In recent years, FinDer has been adopted in various EEW systems, including the U.S. West Coast ShakeAlert warning system (Given et al., 2018; Chung et al., 2019; Kohler et al., 2020). FinDer is also undergoing real-time testing in Central America (Nicaragua, Costa Rica, and El Salvador), Chile, and Switzerland (Böse et al., 2018). Integration of FinDer into the nationwide EEW system of Chinese mainland is underway. Furthermore, Böse et al. (2021) recently applied FinDer to calculate line-source models for large (M > 6) global earthquakes using felt reports.

FinDer Algorithm

FinDer determines earthquake line-source models by matching the spatial distribution of the recorded high-frequency strong ground-motion (usually acceleration) amplitudes with theoretical template maps, which are computed for different line-source lengths and magnitudes. The strike of the fault rupture, Θ, is determined by rotating the templates around various trial angles and calculating the respective misfit with the spatial distribution of recorded PGA. Combined with the optimum spatial position and orientation, the template with the smallest misfit and highest correlation with the recorded amplitudes is found from a combined grid-search and divide-and-conquer approach (Böse et al., 2018). FinDer is computationally highly efficient and allows updating outputs every second until peak shaking across the seismic network is reached.

Ground-Motion Prediction Equations

Ground-motion prediction equations (GMPEs) are an essential element of FinDer. They are used to generate templates that FinDer matches with spatial distributions of observed peak ground acceleration (PGA) amplitudes. In most installations, FinDer adopts the GMPEs developed by Cua and Heaton (2009) that were derived from strong-motion records in California (M 2–M 7.3) and the Next Generation Attenuation (NGA) strong-motion dataset (M 2–M 8). Real-time and offline tests of FinDer in the U.S. West Coast, Switzerland, Italy, Japan, and China (Böse et al., 2012, 2015, 2018, 2021; Li et al., 2020a) have confirmed that the templates are globally applicable to crustal earthquakes. In previous work (Li et al., 2020a), we compared the GMPE of Cua and Heaton (2009) with the regional relationship used in the fifth-generation hazard zoning map in southwestern China, as well as with PGA values, that were observed during the 2008 M_S 8.0 Wenchuan, 2013 M_S 7.0 Lushan and 2017 M_S 7.0 Jiuzhaigou earthquakes. We found that the GMPEs and PGA values are in good agreement, so we use the same templates based on Cua and Heaton (2009) here in this study.

FinDer Magnitude

For small earthquakes (M < 5.5) and at the early stage of rupture in larger events, the FinDer magnitude, M_FD, is determined from the regression of P- and S-wave acceleration amplitudes using relations of Cua and Heaton (2009). For larger events (M ≥ 5.5) M_FD is computed from empirical rupture length-magnitude relationships (typically Wells and Coppersmith, 1994) using the length of the FinDer determined line-source. FinDer line-source models were introduced to characterize the temporal-spatial distribution of high-frequency ground-motions (PGA) rather than to provide an accurate source characterization (Li et al., 2020a). These high-frequency ground-motions are of primary concern in practice and therefore are critically important to real-time seismology. Because of its high-frequency physical background, M_FD is more consistent with the energy magnitude (Picozzi et al., 2017) rather than the long-period M_W. In general, M_FD can be thought of as a scaling factor that quantifies the level and temporal-spatial distribution of high-frequency ground-motions (Böse et al., 2018; Li et al., 2020b). M_FD does not necessarily show great agreement with other magnitude scales, especially with those that are sensitive to longer periods such as M_W. Compared to the magnitude scales that are based on seismic displacement amplitudes, M_FD has the strong advantage that it does not saturate in large earthquakes. This is because FinDer does not use absolute amplitudes, but spatial patterns of ground-motions in order to derive source dimensions which scale with magnitude.

Study Region and Available Datasets

The Sichuan-Yunnan Region

We select the Chinese Sichuan and Yunnan provinces (98.5^oE–106.0^oE, 22.0^oN–34.0^oN) as our study area (Figure 1). This region is exposed to a very high earthquake risk and has been hit by a number of large and damaging earthquakes in the past, including the 1833 Songming M 8 earthquake with about 6,700 fatalities, the 1879 Wudu M 8 earthquake with about 22,000 fatalities, and the 2008 Wenchuan M_S 8.0 earthquake with about 90,000 fatalities (Department of Earthquake Disaster Prevention of the State Seismological Bureau, 1995; Department of Earthquake Disaster Prevention of the China Earthquake administration, 1999). Since January 1, 2009, roughly 530,000 earthquakes have been recorded in this region and their waveforms have been archived by the China Strong Motion Networks Center (CSMNC) and the China Earthquake Networks Center (CENC; Zheng et al., 2010) at CEA.

Before 2008, 105 broadband seismic and 383 strong-motion stations were deployed in the study region (Figure 1A). These stations had interstation distances of around 30 km along the line from Chengdu to Kunming and of 15–20 km in some key areas (Figure 2A). The median interstation distance in the current networks is around 60 km; in 30% of the study area the interstation distance is less than 35 km (Figure 2C). We computed these values from the average distance of a given site to its closest three broadband and/or strong-motion stations (Kuyuk and Allen, 2013) using a mesh of 0.1 × 0.1 resolution.

FIGURE 2. Maps showing the interstation distances for the (A) existing, and (B) planned network shown in Figure 1. Areas colored in dark green have the densest station coverage. For any given site (computed on a grid of 0.1 × 0.1 resolution), the interstation distance is calculated as the average distance to the three closest seismic stations. In order to compare with networks in Japan, California and Beijing Capital Region, we apply the same color scale range as in Kuyuk and Allen (2013) and Li et al. (2016). The histogram in (C) shows the variability in the interstation distances across the Sichuan-Yunnan region. The interstation distance is currently around 60 km for half of the study region and will be reduced to around 30 km for three quarters of the whole region by 2022.

By 2022, 434 broadband seismic stations (equipped with three-component seismometers and accelerometers), 488 strong-motion stations (equipped with three-component accelerometers) and 1,609 low-cost intensity sensors (equipped with MEMS) will be added or upgraded throughout the study region through the construction of the National System for Fast Report of Intensities and Earthquake Early Warning (Figure 1B). With these installations, three quarters of Sichuan and Yunnan region is expected to reach a mean interstation distance of about 30 km; 70 and 40% of the area will have interstation distances of less than 30 and 15 km, respectively. In large parts of Chengdu and Kunming urban areas, the interstation distance will be even 10–15 km only (Figure 2C).

Waveform Data

In this study, we use seismic waveform records of the 2008 Wenchuan M_S 8.0 earthquake and of 57 earthquakes with 5.0 ≤ M_S ≤ 7.0 and depth of H ≤ 35 km that occurred in the study region between January 1, 2009 and August 4, 2019 (Figure 1 and Table 1). Based on their moment magnitude, M_W, we divide these 58 earthquakes into two groups: for earthquakes in the first group (with M_W < 6) we use both broadband and strong-motion records; for the second group, consisting of six earthquakes with M_W ≥ 6, we use strong-motion records only. Throughout this paper, we will treat the earthquakes in the first group as point-source events, and those in the second group as finite-source events (with rupture lengths of 10–300 km).

TABLE 1. Source parameters of the 2008 Wenchuan M_S 8.0 and 57M_S ≥ 5.0 earthquakes (2009-01-01–2019-08-04) in the study region (98.5^oE–106.0^oE, 22.0^oN–34.0^oN) from China Earthquake Networks Center (CENC) and the Global Centroid-Moment-Tensor (GCMT).

Observed Waveforms

For each of the earthquakes in our dataset, we exclude the recordings of obviously malfunctioning broadband/strong-motion sensors and those with PGA smaller than 3 cm/s². For baseline correction we subtract from each trace the mean noise amplitude taken over a time window of several seconds before the P-wave onset.

The broadband seismic records provided by the Data Management Centre of China National Seismic Network at Institute of Geophysics of CEA start from the origin time (OT) of the earthquakes, so there is no need for time corrections. We only correct for the sensor gain and differentiate the time series to obtain acceleration. A large number of high-gain broadband seismic waveforms recorded during the large earthquakes, however, are clipped. Therefore, we only use broadband data for earthquakes in the first group with M_W < 6. In real-time operation, FinDer uses both strong-motion and broadband station data, but discards sensors with clipped amplitudes. In smaller and distant earthquakes, (unclipped) broadband records can serve as an important supplement to strong-motion records if station density is sparse. The broadband seismic three-component traces for each event in the first group are provided in the Supplementary Material.

The strong-motion records are obtained from CSMNC. These records start at a certain time (usually 20 s) before the P-wave arrival. To recover the absolute time, we take the source parameters in Table 1 and align the theoretical P-wave arrival (assuming V_P = 6.0 km/s) with the P-arrival in each record. FinDer uses the temporal-spatial distribution of PGA, including the information that certain stations have not yet recorded strong shaking at a given time. Therefore, we fill pre-event data gaps with random noise amplitudes that we determine from the recorded noise before the P-wave. As this is only for testing FinDer, we are not concerned about the details of this noise, but are primarily interested in simulating realistic amplitudes and a representative temporal evolution of PGA. The simulated noise before P-wave arrivals does not impact the FinDer results. All strong-motion three-component traces for the studied earthquakes are shown in the Supplementary Material.

Seismic waveform data recorded by low-cost intensity sensors have not yet been systematically archived. As described in the next section, we simulate these records from an empirical waveform envelope approach. In order to make these simulations as realistic as possible, we add observed background noise. The typical noise level of a low-cost intensity sensor is shown in Figure 3 in the time domain and in Supplementary Figure 1 in the frequency domain. The waveforms were originally in velocity. To convert them to acceleration we applied the same differentiation method as used for the broadband records. The maximum pulse acceleration of the horizontal (north-south and east-west) and the vertical components in Figure 3 is ±1 cm/s² and ±8 cm/s², respectively. The background noise level for all three-component traces is ±0.5 cm/s². The background noise level is obviously higher than that recorded by broadband seismic (usually about ±0.01 cm/s²) and strong-motion (usually about ±0.05 cm/s²) stations.

FIGURE 3. Three-component noise time series recorded at a low-cost intensity sensor (23.67^oN, 116.64^oE). The records start on February 19, 2019, and are 50 h long. The noise is largest on the vertical component, likely due to human activities. See Supplementary Figure S1 for frequency domain plot.

Simulated Waveforms

To simulate the acceleration time series and the temporal evolution of PGA at sites where broadband, strong-motion and low-cost intensity sensor stations are expected to be deployed in the future (Figure 1B), we compute waveform envelopes developed by Cua and Heaton (2009) for the first group of earthquakes (M_W < 6), and combine stacked waveform envelopes (Yamada and Heaton, 2008) with stochastic simulations (Böse et al., 2012) for the second group (M_W ≥ 6). For the second group, the fault rupture is divided into a number of smaller sub-faults, each radiating P- and S-waves once the rupture front arrives. The radiated waves from each sub-fault are computed from a stochastic time series assuming a simple Brune source model (Brune, 1970), random phase, and a characteristic waveform envelope (Cua and Heaton, 2009). Based on the reference line-source model in Supplementary Table 1 (which is derived from the 1 day-long aftershock distribution), we model the Wenchuan earthquake with 25 sub-faults of 10 km length (around M 6); for the other five earthquakes we use sub-fault of 5 km length (around M 5.5). We assume a constant rupture speed of 2.8 km/s and use the crustal velocity model of Wang et al. (2003) (Supplementary Table 2) for the Sichuan-Yunnan region to model the onsets of the P- and S-wave envelopes.

In order to make our simulated waveforms as realistic as possible (in terms of both the envelope shape and peak amplitudes), we compare them to the observed broadband seismic and strong-motion data at the same station (see Supplementary Material) and conduct the following corrections: first, we use the observed broadband records (that have reliable time measurements) to calibrate the phase arrivals in the simulations. Second, we use both strong-motion and (unclipped) broadband records to calibrate peak amplitudes (PGA).

Figure 4 shows the peak values (PGA) and temporal evolution of amplitudes for both the observed and simulated series at each station before and after the calibration. This figure compares the times of when the series reach PGA (diamonds) and first reach 4.6 cm/s² (rectangles), respectively. We select 4.6 cm/s², because this is the FinDer trigger threshold that we will use later (see FinDer Thresholds). Results for other thresholds are shown in Supplementary Figure 2. We determine the following corrections: 1) using the P-wave arrivals of observed broadband records as benchmarks for our calibration, we move the phase arrivals in the simulations forward by 2 s (that is the onsets of the original simulations tend to be 2 s late); we apply the same time correction also to the simulated data at strong-motion stations, and 2) using the peak amplitudes of the observed strong-motion and broadband records as a reference, we amplify the amplitudes in the synthetics by a factor of 2 for the strong-motion stations; amplitude corrections for simulations of broadband records are not required.

FIGURE 4. Comparison of observed and simulated waveforms for the 58 studied earthquakes in terms of (A) peak ground acceleration (PGA) amplitudes and, (B) the time when they first reach PGA and 4.6 cm/s² (which we will later use as a trigger threshold for FinDer), respectively, before correction. Inf means that the three-component records of the corresponding station never reach this threshold. (C) and (D) show the same comparison after amplitude and time correction.

Adopting these corrections, we simulate the time series at sites where new stations are expected to be deployed in the future (Figure 1B). We add Gaussian white noise with power spectral density (PSD) of −120 cm²/(s⁴·Hz) and −100 cm²/(s⁴·Hz) to broadband seismic and strong-motion stations, respectively. Additionally, we use a series of random time windows to intercept the observed MEMS data in Figure 3 and add them as background noise to the simulated data in low-cost intensity sensors. In this process, we do not remove the large background ground-motion spikes, in order to make our simulations as realistic as possible. As an example, Figure 5 shows the observed and simulated vertical envelopes/waveforms for two events with M_W ≥ 6 and M_W < 6, respectively. The simulated envelopes/waveforms for the other events are shown in the Supplementary Material.

FIGURE 5. Observed and simulated waveforms (vertical component) for the (A) 2008 Wenchuan M_W 7.9 earthquake and (B) a smaller M_W 4.9 earthquake, corresponding to events no.1 and 48 in Table 1. The Wenchuan earthquake is modelled from multiple stacked envelopes and a stochastic time series; the smaller event from a single waveform envelope. See main text for further details.

Results

With the processing described above, we obtain two datasets: dataset 1 contains the observed data recorded at the existing stations (Figure 1A); dataset 2 contains both the data from dataset 1 plus the simulations at future stations (Figure 1B). In this section, we use these two datasets to simulate and analyze the off-line performance of FinDer in Chinese mainland. In order to mimic the temporal evolution of the FinDer estimated line-source models, we determine in this retrospective study the peak absolute ground-motion amplitudes (taken over all three waveform/envelope components) in time windows of 1 s and feed the logarithmic value of these amplitudes into FinDer (Böse et al., 2018). However, before running the waveform playbacks, we first need to configure the FinDer trigger thresholds, which depend on station density.

FinDer Trigger Thresholds

As described in Böse et al. (2018), FinDer uses a cascade of increasing PGA thresholds for its binary template matching. The smallest PGA threshold (and the minimum number of stations at which this threshold needs to be exceeded) controls the triggering of the algorithm. Typically, these parameters are set as 2 cm/s² and 3, respectively (Böse et al., 2018). However, it is necessary to tune these values in order to account for the respective station density and noise characteristics of the network in which FinDer is applied. In this study, we use the 2008 Wenchuan earthquake to define the trigger parameters.

We are testing in Figure 6 four PGA trigger thresholds for the Wenchuan earthquake (2.0, 4.6, 10.5, and 23.2 [cm/s²]), corresponding to the average peak acceleration caused by an earthquake of M = 2.5, 3.0, 3.5, and 4.0 at 5 km distance will produce (Cua and Heaton, 2009). In parallel, we vary the minimum number of triggered stations from 1 to 4. Details are given in Supplementary Table 3. For some configurations, e.g. for a trigger threshold of 2 cm/s², FinDer creates multiple solutions (split events), which is not desired.

FIGURE 6. Testing FinDer trigger parameters for different PGA thresholds (top to bottom: 2.0, 4.6, 10.5, and 23.2 cm/s²) and minimum numbers of stations (left to right: 1–4), at which this threshold needs to be exceeded, using the 2008 M_S 8.0 Wenchuan earthquake. The maps show the final FinDer line-source models (black lines) for the respective trigger configuration. In general, the trigger parameters depend on the density and noise conditions of the network in which FinDer is applied. For some configurations, e.g. for a trigger threshold of 2 cm/s², FinDer creates multiple solutions (split events). As a compromise between accuracy and timeliness of the line-source models, we chose in this study a trigger threshold of 4.6 cm/s² and a minimum number of three stations.

As a compromise between the accuracy and timeliness (the lower the threshold and the fewer stations are required, and the faster the first FinDer model is calculated) we chose for this study a trigger threshold of 4.6 cm/s² and a minimum number of three stations at which this threshold needs to be exceeded. With this setting, FinDer triggers 12 s after the OT of the Wenchuan earthquake (neglecting data latencies). The final FinDer line-source model is estimated as L = 249 km long (corresponding M_FD 7.9) and with the strike of 55^o. These values are in good agreement with the results of field surveys following the Wenchuan earthquake that identified two separate parallel surface ruptures of L = 240 km and L = 90 km at a strike of N45^oE (Xu et al., 2008). The FinDer line-source model characterizes the entire fault rupture including those segments with no surface exposure. Compared with the results of Li et al. (2020a), who tested FinDer with a larger trigger threshold (20 cm/s²), the length of the line-source model is underestimated by about 50 km (corresponding to 0.1 magnitude units). This difference is caused by sparse station density around the Wenchuan earthquake as discussed in Li et al. (2020b).

To determine the FinDer trigger parameters for this study (Figure 6), we used the Wenchuan records in Dataset 1, which contains only observed data, i.e. no simulations at future stations. In general, we expect that the trigger threshold can be chosen to be lower, the denser the network. The choice of the threshold, however, should be discussed in the context of the whole system: higher thresholds come at the cost of the system being triggered later, which means less available warning time. In practice, one can only trade-off the threshold (i.e. the system trigger time) and the accuracy. There is currently no strict empirical relationship between threshold and station density. For a better comparison of the performance of FinDer in a sparse and dense sensor network, we will adopt for this study the same trigger parameters for both datasets.

Relationship Between FinDer Trigger Thresholds and Station Density

In this section, we investigate the relationship between the FinDer trigger thresholds and station density. Using the trigger configuration optimized for the Wenchuan earthquake (4.6 cm/s² at three or more stations) we plot in Figure 7 all earthquakes that trigger FinDer and those that do not. A first-order function to separate triggered and non-triggered events can be determined from a representative GMPE. This is illustrated by the blue line in Figure 7A, which is determined by transforming the GMPE of Cua and Heaton (2009) into an explicit function of the epicentral distance, R_epi, and M_W, R_epi = log₁₀ (M^a + b) + c. For our trigger threshold of 4.6 cm/s² and 0 ≤ R_epi ≤ 120 km, we find a = 1.91, b = 5.93, and c = 1.64. As shown in Figure 7A (blue line), this relation, however, does not provide a satisfying boundary between triggered and non-triggered events. This is expected, because our data is plotted as of function of D_epi, that is the average distance between the epicenter to the closest three broadband seismic and/or strong-motion stations, and not of R_epi. Keeping the same functional form and coefficients a and b, we find that FinDer has a high probability to trigger, if D_epi ≤ log₁₀ (M^a + b)+c, where a = 1.91, b = 5.93, and c = 2.34 for M = M_W ≥ 4.8, and c = 2.49 for M = M_S ≥ 5.0, respectively (red lines in Figures 7A,B). Note, however, that these functions provide a first-order approximation only, because there is in general considerable variability in the seismic ground-motions and the two datasets (observed data v.s. observed + simulated data) cannot fully constrain these relations as we will discuss later.

FIGURE 7. Testing FinDer trigger settings (4.6 cm/s² at neighboring three stations) for all 58 earthquakes as a function of interstation distance, D_epi, and (A) moment magnitude, M_W (here of GCMT) and (B) surface-wave magnitude, M_S (here of CENC), respectively. D_epi is defined as the average distance between the epicenter to the three closest stations. The blue line shows the GMPE of Cua and Heaton (2009) for PGA = 4.6 cm/s². The two red lines illustrate the boundary condition: FinDer is expected to trigger, if D_epi ≤ log₁₀ (M^a + b)+c, where a = 1.91, b = 5.93, and c = 2.34 for M = M_W ≥ 4.8, and c = 2.49 for M = M_S ≥ 5.0, respectively. The datasets used here, however, cannot fully constrain these relationships, and our boundary conditions provide a rough estimate only.

When calculating D_epi in Figure 7, we also included stations that did not work properly. This means that our results consider the possible failure of stations. This likely explains the two data outliers in Figure 7: a M_W 4.8 for which FinDer triggers even though station density is sparse (Figure 7A), and a M_S 5.1 event that FinDer misses (Figure 7B). Of course, in general there is also a significant level of variability in the radiation of high-frequency motions.

Based on the results in Figure 7, we can roughly determine a magnitude completeness, M_FDc, which defines the lowest magnitude at which FinDer is likely to trigger, given a certain interstation distance, D_epi. For one-third of the Sichuan-Yunnan region, where the current interstation distance is around 30 km, the smallest detectable events are M_FDc ≈ M_S 5.3. With the planned network densification with interstation distances of around 20 km in half of Sichuan-Yunnan after 2020, we expect M_FDc ≈ M_S 5.0 and M_FDc ≈ M_S 4.4 in some urban areas.

Playback Results

Figures 8–11 and Supplementary Table 4 illustrate and summarize the FinDer playback results for all 58 earthquakes. Additional information is given in the Supplementary Material. Figure 8A shows the magnitude difference between the FinDer estimated final magnitude, M_FD, and the reported M_W from GCMT and M_S from CENC, respectively. The average difference is 0.73 for M_W and 0.53 for M_S in dataset 1, and 0.64 for M_W and 0.47 for M_S in dataset 2. Overall, M_FD tends to be more consistent with M_S, which is in line with the high-frequency physical background of M_FD as discussed in FinDer Magnitude. The comparison of the two datasets suggests that the future networks will be able to better constrain earthquake magnitudes once more stations are deployed.

FIGURE 8. (A) Difference between the final FinDer estimated magnitude, M_FD, and the reported M_W/M_S. (B) Distance between the estimated final epicenter (for M_W ≥ 6 only)/centroid (for all events) and the reported epicenter/centroid. The solid and dashed lines show the average values for M_W and M_S for dataset 1 (black) and dataset 2 (gray), respectively.

FIGURE 9. Temporal evolution of FinDer magnitude, M_FD, from event origin time (OT) for all 58 earthquakes (neglecting data latencies).

FIGURE 10. Temporal evolution of the distance between the FinDer centroid and the reported epicenter (for M_W < 6)/centroid (for M_W ≥ 6) from event origin time (OT) for all 58 earthquakes (neglecting data latencies).

FIGURE 11. Comparison of the strike of the final FinDer line-source and the GCMT focal mechanism for all 58 earthquakes. The absence of FinDer strike means that FinDer has not been triggered.

Figure 8B illustrates the distance between the estimated final centroid as the midpoint of the line-source and the reported epicenter (for M_W < 6)/centroid (for M_W ≥ 6). The FinDer centroid defines the mid-point of the FinDer line-source. Therefore, Figure 8B also shows the distance between the estimated and reported epicenters of the six larger earthquakes. It should be noted that 22 out of 25 earthquakes that do not generate FinDer triggers in dataset 1, produce triggers for the extended dataset 2, which includes additional simulations at future stations. Except for the significant difference between the two datasets for the two major earthquakes, the results for the other events are similar (Figure 8B). This can be explained by the fact that for the 33 earthquakes that triggered FinDer in dataset 1, the additional stations in dataset 2 have little impact on the epicenter constraints. This implies that the FinDer estimated epicenter/centroid is stable even when the network is quite sparse (provided that the azimuthal gap is not too large; see Böse et al. (2021)). In Figure 8B, the average value of the epicenter/centroid difference (18 km) obtained from dataset 2 is smaller compared with that from dataset 1 (24 km). This means that the 25 earthquakes that triggered FinDer only for the denser network in dataset 2 all have smaller epicenter errors.

Figures 9–10 show the temporal evolution of M_FD for all 58 earthquakes, as well as the distance between the FinDer epicenter/centroid and reported ones from OT, respectively. For the Wenchuan earthquake, FinDer yields its final line-source around 120 s after OT. For the other five large earthquakes (M 6.0–6.9), the final results are obtained within 30–40 s after OT. For the smaller earthquakes (M < 6), the final values are usually reached within 10–20 s. In general, FinDer provides 5–10 s faster and more accurate estimates of the final rupture for the earthquakes in dataset 2, which simulates a denser network. Figure 11 shows the final FinDer strike from dataset 1 with the background of the focal mechanism from GCMT. The line-source model calculated by FinDer immediately after the earthquake can help to quickly determine the underlying seismogenic fault if prior local tectonic/geological information is available, and that is of great significance in the likelihood of strong aftershocks, time-dependent seismic hazard assessment, and even operational earthquake forecasting (OEF). Note, however, that FinDer does not require any prior information on fault locations.

Implications for EEW

The goal of EEW is to estimate the earthquake’s damage potential and to issue an alarm of expected ground-motions to areas before the strong shaking starts. The primary concern of EEW is the time delay between triggering and OT as this time is directly related to the EEW blind-zone (or no-warning zone). Figure 12 illustrates the relationship between D_epi and the time when FinDer triggers relative to OT, T_tri. Based on a series of evenly distributed modelled point-sources and stations, Kuyuk and Allen (2013) determined a semi-quantitative empirical curve that describes the relationship between the minimum number of triggered stations, the average interstation distance and the radius of the blind-zone. In Figure 12 we regress a similar curve for T_tri as the function of D_epi and find that

T_tri can be converted into the radius of the no-warning zone by simply multiplying the S-wave velocity (e.g. V_S = 3.5 km/s).

FIGURE 12. Relationship between the interstation distance, D_epi, defined as the average distance between the epicenter to the three closest stations, and FinDer trigger time from OT, T_trig.

The regression curve in Figure 12 has two characteristics: 1) for D_epi > 15 km, T_tri (and thus the radius of the no-warning zone) grows exponentially with D_epi, which means that network densification has a very strong impact on reducing T_tri; 2) for D_epi < 15 km, T_tri saturates and further network densification has little impact on T_tri. However, it should be noted that this does not imply that a densification beyond 15 km is useless. Given the actual operation of a seismic network, it is common for a subset of stations to malfunction. In a threshold triggered system (e.g. based on FinDer), stations near the epicenter that are not functioning properly have the same effect as an increase in D_epi. Densifying networks thus effectively controls the growth of no-warning areas (Li et al., 2016).

The on-going network densification in Sichuan-Yunnan from the current mean value of D_epi = 60 km (Figure 2A) to D_epi = 30 km (Figure 2B) will significantly reduce T_tri (and thus no-warning areas) across the whole region after 2020: based on Eq. 1 we expect T_tri to decrease from currently around 8–9 s for one-third of the region to around 4–5 s for half, and even 3–4 s in the key urban areas. These estimates do not include data latencies. Preliminary tests for newly-built and upgraded stations suggest that data latencies of 1 s are a realistic expectation (Peng et al., 2020; Wenhui Huang, written communication, 2021).

Finally, we investigate how much warning time FinDer could provide for different levels of ground-shaking (again neglecting data latencies). We doing this by analyzing the potential performance on a per-station basis. We define the warning time (at each station) as the time interval between the first prediction (which is not necessarily in the first report) that ground-motions will exceed 10 cm/s² using the GMPEs developed by Cua and Heaton (2009) and the actual first exceedance of this level. If the observed and predicted PGA values both exceed that threshold, the station is considered to have hit the warning (i.e. true positive, TP). The occurrence of peak shaking usually occurs later, so the true warning times might actually be a bit longer.

Figure 13 shows the TP ratios of stations with successful alerting to the total number of stations for different intensity bins and for different ranges of warning time (≥1 s, ≥ 5 s, ≥ 10 s, ≥ 15 s, and ≥ 20 s). The results are calculated for the Chinese intensity scale, but the difference to the MMI scale is small (Hu, 2009). If the set of broadband seismic and strong-motion data used in this study had been available in real-time and without (or with short) delay, 40–70% of sites experiencing intensities of V–VIII and 20% experiencing IX-X could have been issued a warning with 5–10 s for the current station network (Figure 13A). Once the full network is deployed, these warning times are expected to increase: as shown in Figure 13B the TP ratios are expected to increase by 5–20% for 5–20 s of warning for V and VIII, and by almost 20% for IX for all warning time intervals. Adapting the FinDer trigger setting to account for the denser network, will probably lead to additional improvements in the new network.

FIGURE 13. Percentage of true positives (TP) for different levels of intensity (both Chinese and MMI scale) and different ranges of warning time using FinDer models for the (A) current network (dataset 1) and (B) future network (dataset 2). Numbers in brackets on x-axis give station counts for each intensity bin.

Discussion

The triggering of an EEW algorithm obviously depends on the earthquake magnitude and the network density near the epicenter. We found that FinDer has a high probability to trigger if the network density around the epicenter, D_epi, is less than log₁₀ (M^a + b) + c, where a = 1.91, b = 5.93, and c = 2.34 for M = M_W ≥ 4.8, and c = 2.49 for M = M_S ≥ 5.0, respectively. The relationship is unlikely linear. Otherwise, it would suggest that FinDer will not trigger for any earthquake below magnitude 4.5, regardless of station density. Real-time and offline tests of FinDer in other parts of the world, however, have shown that earthquakes as small as M 2.5 can be detected (Böse et al., 2018). Our dataset is probably insufficient to fully constrain the boundary condition for FinDer triggering.

Peng et al. (2020) recently analysed the alert performance of a hybrid demonstration EEW system installed in the same region as studied here. This system is equipped with MEMS-based sensors and broadband seismic stations with low-latency data transmission of 0.5–1.5 s, and determines the magnitude by the P_d scaling algorithm. For nearly half of the 130 M_L 3.0–5.1 earthquakes that occurred in this region between 2017 and 2018 the demonstration system triggered and sent first alerts within 15 s from OT. Using a different dataset and a larger magnitude range (M 5.0–8.0), our results in Figure 12 suggest that FinDer could detect nearly two-thirds of the events within 12 s from OT (neglecting latencies).

FinDer determines seismic line-source models assuming that the dislocations are evenly distributed along this model. Such a model can be understood as the integral of the Haskell model (Haskell, 1964, 1966, 1969) in time and the width of the rupture surface. This assumption is theoretically reasonable. Using a list of candidate source models (homogeneous slip model, triangular slip model, k-square model, slip tip taper model, and restricted stochastic source model), Li et al. (2020b) applied the Akaike information criterion (AIC) to identify the model that is most consistent with inverted slip distributions for the 2011 Van (Turkey) M_W 7.1 earthquake. The authors found that the homogeneous slip model with the smallest degree-of-freedom is only inferior to the k-square model and the triangular slip model. This suggests that the homogeneous slip model performs reasonably well.

The recently proposed extension of FinDer, FinDerS (+), assumes an elliptical slip model described by two independent parameters: the maximum slip and a skewness factor that controls the asymmetry of the slip profile (Hutchison et al., 2020; Böse et al., 2021). This model has similar features and the same degree-of-freedom like the triangular slip model, so this could be an important improvement, even though a systematic evaluation of source models for a representative set of earthquakes is needed, such as provided through the Source Inversion Validation (SIV) project and US Geological Survey finite-source earthquake database.

Combing different (point-source/finite-fault) algorithms allows EEW systems to deliver fast and accurate alerts for a large range of magnitudes. However, this combination in an actual EEW system is challenging. In the US ShakeAlert system, for instance, the Solution Aggregator (SA) module associates the source parameters from the EPIC and the FinDer algorithms if locations are within 100 km and OTs within 30 s, respectively. The combined source parameters are then determined as the average of the estimated parameters weighted by the normalized uncertainties of each algorithm (Kohler et al., 2020). In the future, both the event association and solution combination will likely move to the ground-motion space (Minson et al., 2017).

Conclusion

Rapid finite-rupture information, such as provided by FinDer, is important for EEW (and rapid response), in particular in large earthquakes (Böse et al., 2012; Li et al., 2020b; Böse et al., 2021). In this study, we tested the performance of the FinDer algorithm (Böse et al., 2018) for 58 5.0 ≤ M_S ≤ 8.0 earthquakes in the Sichuan-Yunnan region of China. We used playbacks of waveforms recorded by the high-gain seismic network stations and simulated waveforms/envelopes for planned future station installations. Based on the 2008 Wenchuan earthquake, we determined the optimum FinDer trigger setting as 4.6 cm/s² at neighboring three stations as a compromise between speed and accuracy of alerts. Overall, the FinDer estimated line-source models agree well with the catalog source parameters and focal mechanisms, but they can be computed faster compared to traditional approaches. The planned densification of seismic networks in Chinese mainland in the near future is expected to further reduce magnitude uncertainties and allow to better constrain epicenter and centroid. In another study (Böse et al., 2021) we recently developed an empirical equation that relates FinDer performance (in terms of length, strike, and location of the line-source) to the azimuthal gap between observing stations. While calibrated for the moment magnitude M_W, we compared FinDer magnitudes in this study also to the surface wave magnitude M_S, which is more commonly used in China. We found that the overall difference between M_FD and M_S is smaller than between M_FD and M_W, which is consistent with the high-frequency physical background of M_FD.

We estimated the magnitude of completeness for FinDer, M_FDc, as 5.3 for one-third of the Sichuan-Yunnan region with interstation distance of around 30 km, and as 4.4 in some urban areas after 2022. The FinDer trigger time, T_tri, may be reduced from currently around 8–9 s for one-third of the Sichuan-Yunnan region with D_epi = 30 km to around 4–5 s for half, and to 3–4 s in the key urban areas (plus data latencies). If the waveform data used in this study had been available in real-time, 40–70% of sites experiencing intensities of V–VIII and 20% experiencing IX-X could have been issued a warning with 5–10 s, before the arrival of the S-wave.

This study provides a useful reference for the planned future installation of FinDer in the nationwide EEW system of Chinese mainland. With the increasing number of broadband seismic, strong-motion and low-cost intensity stations with low latency being deployed in the near-future, FinDer has the potential to quickly report finite-source parameters (e.g. length, strike, rupture direction, rupture velocity etc.), and is expected to significantly improve both EEW and rapid loss estimates in China (Li et al., 2020a).

Data Availability Statement

Strong-motion waveforms used in this study were provided by the China Strong Motion Network Center (CSMNC) at the Institute of Engineering Mechanics (China Earthquake Administration). Broadband seismic records were obtained from the Data Management Centre of China National Seismic Network at Institute of Geophysics (IGP, SEISDMC, doi: 10.11998/SeisDmc/SN) at China Earthquake Administration. The earthquake catalog was downloaded from China Earthquake Networks Center (CENC; http://10.5.160.18/console/exit.action, last accessed: December 2020). The moment magnitude and focal mechanism were extracted from Global Centroid-Moment-Tensor (GCMT; https://www.globalcmt.org/, last accessed: January 2021). The fault data is taken from GMT Chinese community (https://gmt-china.org/data/, last accessed: January 2021).

Author Contributions

JL: computed FinDer models, conducted analysis, wrote manuscript; MB: gave idea and guided research, wrote manuscript; YF: helped with discussions and writing manuscript; CY: provided the planned station information, helped with discussing the Chinese mainland’s EEW project.

Funding

This project is supported by the Guangdong Basic and Applied Basic Research Foundation (Grant No. 2020A1515110844). JL’s visit at ETH Zürich is funded by the China Scholarship Council (201804190032) and the Swiss Seismological Service (SED).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors would like to thank Liangqiong Lou (China Strong Motion Networks Center), Xiaoxiao Song (Hubei Earthquake Agency), Xiaojing Ma (Guangdong Earthquake Agency), Tianyang Gao and Yan Zhang (Institute of Geophysics, China Earthquake Administration) for helping us to apply for and download the earthquake data. The authors would also like to thank Wenhui Huang (Shenzhen Academy of Disaster Prevention and Reduction), Chaoyong Peng (Institute of Geophysics, China Earthquake Administration), John F. Clinton (Swiss Federal Institute of Technology Zürich), Frédérick Massin (Swiss Federal Institute of Technology Zürich) for related scientific discussions, as well as the topic editor Masumi Yamada and two reviewers for their feedback and suggestions.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/feart.2021.699560/full#supplementary-material

References