Abstract

A new methodology to automatically determine the source process and near-source strong ground motions using regional and local distance data was used to study the October 16, 1999 Hector Mine earthquake (M_W7.1). Broadband displacement data were used to determine the seismic moment tensor and resolve the fault plane ambiguity independent of aftershock location and surface faulting information. The identified NNW striking fault plane was used to invert for the distribution of fault slip over a planar surface. The slip model was used to predict near-fault peak ground velocity, which was found to agree well with the TriNet ShakeMap, and also with near-fault observations that were not available at the time of our initial analysis. Our method if automated is capable of providing near-source strong ground motion information within 30 minutes of the earthquake origin time even in cases where there are few near-fault recording instruments.

Introduction

It is possible to rapidly determine the time, location and magnitude of an earthquake, and while this information is of value in characterizing the extent of an earthquake emergency a regional assessment of the extent of strong shaking is more relevant to emergency response. In southern California the TriNet ShakeMap (Wald et al., 1999) provides near-realtime characterization of strong ground motions. Where there is very dense coverage, this approach provides a ‘ground truth’ assessment of shaking level. We have developed a regional distance finite fault methodology that makes use of near-realtime broadband data to determine the gross properties of the earthquake source process. The derived source parameters are then used to simulate near-fault strong ground motions. In areas where there is very sparse local coverage (if any at all) the model provides strong ground motion estimates that include the influences of earthquake-specific source finiteness and rupture directivity. In areas where there are adequate numbers of local instruments the method may be integrated with direct observations of strong ground motions using the approach developed by Wald et al. (1999), and can provide needed redundancy in the event of seismic telemetry or other network failure that interrupts the retrieval of near-source strong-motion information. A seamless, regional scale, shakemap method is needed by the emergency response community since much of the supporting transportation, water and power infrastructure, as well as at-risk rural communities, extend beyond areas where dense instrumentation is likely to be deployed. Our approach is being developed for the monitoring situation that exists in sparsely instrumented areas of northern California and for the monitoring situation anticipated for the Advanced National Seismic System (ANSS, e.g. USGS, 1999). We demonstrate the technique using data recorded for the recent October 16, 1999 M_W 7.1 Hector Mine, California, earthquake (02:46:44UTC; 34.594 N 116.271 W; depth=6 km) using broadband data from TriNet (Figure 1).

Figure 1. Map showing regional faults, the surface trace of the fault (bold line), and TriNet stations (triangles). Stations used in the quasi-automatic slip inversion are denoted by open solid triangles. The HEC and BC3 were used to determine a revised slip model. Fault plane solution obtained by inversion of low-frequency waves at the BDSN station KCC (A) is compared to that for an inversion, which included TriNet stations (filled triangles).

Finite Source Methodology

We use three component broadband data that may be obtained from either on-scale velocity sensors or co-located strong-motion instruments, such as is available from the Berkeley Digital Seismic Network (BDSN), TERRAscope, and TriNet, to invert the representation theorem (equation 1) for finite source parameters (e.g. Hartzell and Heaton, 1983). The observed seismogram is related to the spatio-temporal integration of slip distributed over the discretized plane where,

, (1)

and U_n is the n^th component of observed displacement, and m, , and u are the rigidity, a fault orientation unit vector and the fault slip, respectively. G_ni,j is the Green’s function that describes the wave propagation from each point on the fault to the recording station. x is a vector describing the relative location of the source and receiver and x and t are the spatial and temporal variables of integration. The indices n, i, and j refer to the component of displacement and the geographical directions north, east and up, respectively. The quantity is equivalent to , which describes the spatial and temporal seismic moment tensor.

The finite fault inverse method is initiated from seismic moment tensor input, where the moment tensor is used to define the orientation of the two possible double-couple nodal planes and the moment magnitude is used to scale the source model dimension using the relations of Wells and Coppersmith (1994),

A=10^{-3.49+0.91*Mw}; L=10^{-2.44+0.59*Mw}; W=A/L; (2)

where A is the fault area, L is the fault length, and W is the fault width. Because we do not know the rupture plane a priori nor whether the earthquake is unilateral or bilateral, we double the source dimension and add a 20% safety factor when constructing the fault model space. This is done to ensure that the fault model for each of the two possible nodal planes is large enough to accommodate unilateral rupture in either direction. The fault model is discretized and the slip for each discrete sub-fault is solved for.

The integration over dΣ may be transformed to an integration in time by considering a model of the rupture evolution. In our quasi-realtime application we simplify the problem in which a sub-fault is allowed to slip only once at the arrival time of the rupture front. The slip rise time is also assumed to be constant and is estimated from the total rupture duration T=L/V_r, where V_r is the rupture velocity initially assumed to be 80% of the shear wave velocity. The dislocation rise time is taken to be 0.10*T after the observation of Heaton (1990). Alternatively, the empirical relationship, T_R=4.37*10^-7*M₀^1/3 (rewritten for SI units), of Somerville et al. (1999) may be used. To stabilize the inversion slip positivity and smoothing constraints are applied (e.g. Hartzell and Heaton, 1983). Although the source parameterization is simpler than what is commonly applied to strong-motion data sets (e.g. Wald and Heaton, 1994), it is appropriate given the regional distances of the stations used and our objective to use the method to obtain fault dimension and the gross slip distribution in a timeframe useful for emergency response.

Any waveform inverse method for source parameters is only as good as the Green’s functions that are assumed to represent regional wave propagation characteristics. We use the SoCal model to calculate complete Green’s functions because it has been shown that the model is effective in describing regional wave propagation in southern California (Dreger and Helmberger, 1993). The application of this method requires that calibrated Green’s functions are available.

Finite Source Modeling

Our method begins with an estimate of the seismic moment tensor. Figure 1 shows the seismic moment tensor solution obtained by inverting three-component long-period data (0.02 to 0.05 Hz) from the Berkeley Digital Seismic Network station KCC, located 407 km NW of the event. We use this preliminary solution (strike=343^o, rake=175^o, dip=70^o, Mo=4.01*10¹⁹ Nm, M_W=7.0) because it was the one disseminated via email 33 minutes after the event and it is the moment tensor solution that we used in our quasi-realtime implementation of the finite source inversion the morning of the earthquake. This solution agrees well with first motions and other regional and teleseismic moment tensor solutions. It is very similar to the moment tensor solution (strike=330^o, rake=180^o, dip=82^o, Mo=2.07*10¹⁹ Nm) we obtained using the five stations in Figure 1 including KCC. Based on the initial moment magnitude of 7.0 the source dimension (L) was determined to be 48 km and the total dimension of our source model, 115 km. Given L and a rupture velocity of 2.91 km./s (80% of the shear wave velocity at hypocentral depth) the total duration was estimated to be T=16.8s, and a rise time of 1.7s was obtained. If Somerville et al. (1999) is used a rise time of 1.5s is obtained. The factor of two lower scalar seismic moment of the 6 station moment tensor inversion results in a 90 km fault model, which is large enough to accommodate the slip (below), and rise time of 1.2s. The down-stream processing depends upon a moment tensor solution for input and an incorrect moment tensor will produce unreliable results. If the event cannot be approximated by a point-source for f≤0.05 Hz (~M_W7.5) then alternative input such as aftershocks, and teleseismic moment tensor solutions may be used.

The first stage of the finite source processing is the identification of the causative fault plane. To accomplish this we perform a series of line source calculations using broadband (0.02 to 5.0 Hz) displacement data for the two double-couple nodal planes of the moment tensor inversion. In the line source calculations the spatial integration over dS in equation 1 is collapsed to an integration over fault length, dL. Figure 2 compares the variance reduction versus rupture velocity for both the NNW-striking and ENE-striking planes. The NNW-striking plane provides a superior fit to the data, and there is a peak in the variance reduction curve for a rupture velocity of 2.6 km/s. Additionally, the line source results indicate that most of the slip in the event was located SSE of the assumed hypocenter. The calculation of 20 line source models (10 rupture velocities for each plane) required 4.2 minutes of processing time on a 200MHz Sun Ultra 1 Workstation.

Figure 2. Percent variance reduction as a function of rupture velocity.

Figure 3. Fault slip. The entire inverse model space is shown and the hypocenter plotted at the center of the model (white circle).

The second stage inverts the broadband displacement data for planar slip assuming the fault orientation and rupture velocity determined in stage one. For the Hector Mine earthquake the finite fault inversion consisted of 440 sub-faults on the 110x25 km² fault model and required 5.0 minutes of processing time. Figure 3 shows the distribution of slip that was determined using the 5 three-component stations that were available in realtime, and by applying a 3s delay in rupture to account for the emergent onset of the waveforms. The fit to the displacement data is found to be quite good (Figure 4). Most of the slip in the Hector Mine earthquake was located 5-33km SSE of the hypocenter (Figure 3), although the asperity 11-18 km to the northwest indicates a bilateral rupture. There is a deeper patch of slip located approximately 3 km SSE of the hypocenter that fits large lapse time signal in the data, but whose depth is likely an artifact of the simplifying assumptions. The average and peak slip were found to be 138 and 826 cm respectively. The total seismic moment was determined to be 6.85*10¹⁹ Nm equaling a moment magnitude of 7.2. Aftershocks confirm the bilateral nature of the source, and the dimension of the slip model compares well with the observed 40 km of surface faulting (Scientists from the USGS, SCEC and CDMG, 2000).

Figure 4. Synthetic displacement seismograms (dashed) for the model shown in Figure 3 are compared to observations (solid). Both the data and synthetics have been bandpass filtered between 0.02 to 5.0 Hz.

Prediction of Near Source Ground Motions

The final stage of processing solves the forward representation problem (equation 1) using the derived fault slip distribution to determine near-fault strong ground motions. The SoCal model (Dreger and Helmberger, 1993) was used in the forward simulations. A map of simulated PGV (Figure 5A) is obtained from the geometric mean of bandpass filtered (0.02 to 3.0 Hz) synthetic seismograms. Recall that 0.02 to 5.0 Hz data was used to determine the slip distribution. The predicted PGV are found to agree well with observed values, even though the observations are not band limited and the predictions are. In particular, station HEC has a peak value of 37 cm/s and the station lies between the 20 and 40 cm/s contours of predicted PGV even though it was not used to derive the slip model. The area of greatest predicted shaking, has a peak value of 97.6 cm/s and is in an area SE of the epicenter.

Figure 5. (a) Simulated PGV obtained from the slip map in Figure 3 (contours begin at 20 cm/s and are plotted at intervals of 20 cm/s) is compared to observations (values shown in cm/s). The peak velocity is 97.6 cm/s. The boxes identify the four closest stations whose data was available in realtime, and the surface fault trace is shown in bold. (b) Simulated PGV (circles) is combined with direct observations following the interpolation and site correction method of Wald et al. (1999). Contour intervals are the same as in (a). Small scale variation along contours is an artifact of the site corrections. The peak velocity is 99.5 cm/s.

It is possible to combine the model predictions with direct observations (Figure 5B) to produce a more comprehensive shakemap using the procedure outlined in Wald et al. (1999). The circles (Figure 5B) show the locations where model predictions were used to interpolate the map between observed values. It is satisfying that the 20 and 40 cm/s contours only required slight modification from the model predicted shakemap (Figure 5A).

The allowed source process in the applied method has been simplified considerably. In the quasi-realtime results presented here we applied a 3s delay to account for the emergent onset of the waveforms, and we found that the PGV map is qualitatively similar to that in Figure 5A however the peak in PGV is further to the SE. It has since been found that the Hector Mine earthquake occurred on a curved fault surface (Scientists of the USGS, SCEC and CDMG, 2000). We have examined the spatio-temporal complexity of the event using a 3-segment model to approximate fault curvature, variable slip rake angle, and multiple time windows to allow for variable dislocation rise time. A detailed discussion of the results is beyond the scope of this letter and will be the focus future work (but can be viewed by following this link). We have found, however, that when these complexities are considered the slip tends to focus near the hypocenter and is bilateral in nature. The dislocation rise time varies and can be as large as 6s. Forward predictions of PGV based on the revised model are qualitatively very similar to that given in Figure 5A, except that the peak values are approximately 30% lower. Thus the estimated PGV shakemaps for three-levels of approximation describe an elongate, near-fault region of elevated ground motions and the patterns agree well with the observed near-fault observations.

Conclusions

We have used near-regional distance displacement seismograms to constrain the seismic moment tensor, the causative fault plane, and the slip distribution for the October 16, 1999 Hector Mine earthquake, and have shown that it is possible to predict near-fault strong shaking from the derived source parameters. Although the software was not operating in an automated fashion at the time of the earthquake, it was possible to simulate a realtime application. The determination of the best fault plane and slip distribution were obtained the morning of the earthquake, and were subsequently verified by the observed fault offsets and aftershock locations. This method may be applied in both densely and sparsely instrumented regions. In densely instrumented regions it may provide needed redundancy in the event of telemetry or other failure of the primary realtime strong-motion network. In sparsely instrumented areas it may be the only viable approach to estimate the level of near-fault strong-ground motions. Finally, the computer processing times of the individual components of the method as applied to the Hector Mine earthquake indicate that the method may be applied in near-realtime and provide strong shaking information within the 30-minute post-earthquake time frame.

Acknowledgements

We are grateful to the TriNet organizations for making the data used in this study available in a timely fashion. David Oppenheimer and an anonymous reviewer provided constructive comments that greatly improved the manuscript. The methodology used in this paper was conceived using funding from the USGS NEHRP program (USDI 1434-HQ-GR-02742) and was developed as part of a PGE/PEER Phase II contract (PGE-009566). This is contribution number 00-x of the UC Berkeley Seismological Laboratory.

References

Dreger, D. S., and D. V. Helmberger (1993), Determination of Source Parameters at Regional Distances with Single Station or Sparse Network Data, Journ. Geophys. Res., 98, 8107-8125.

Hartzell, S. H., and T. H. Heaton (1983). Inversion of strong ground motion and teleseismic waveform data for the fault rupture history of the 1979 Imperial Valley, California, earthquake, Bull. Seism. Soc. Am., 73, 1553-1583.

Heaton, T. H. (1990) Evidence for and implications of self-healing pulses of slip in earthquake rupture, Phys. Earth and Planet. Int., 64, 1-20.

Scientists from the USGS, SCEC and CDMG (2000). Preliminary report on the 16 October 1999 M7.1 Hector Mine, California, earthquake, Seism. Res. Lett., 71, 11-23.

Somerville, P. G., N. F. Smith, R. W. Graves, and N. A. Abrahamson (1997) Modification of empirical strong ground motion attenuation relations to include the amplitude and duration effects of rupture directivity, Seism. Res. Lett., 68, 199-222.

Somerville, P. G., K. Irikura, R. Graves, S. Sawada, D. Wald, N. Abrahamson, Y. Iwasaki, T. Kagawa, N. Smith, and A. Kowada (1999). Characterizing crustal earthquake slip models for the prediction of strong ground motion, Seism. Res. Lett., 70, 59-80.

U.S. Geologic Survey (1999). An assessment of seismic monitoring in the United States: Requirement for an Advanced National Seismic System, USGS Circular, 1188, 55p.

Wald, D. J., V. Quitoriano, T. H. Heaton, H. Kanamori, C. W. Scrivner, and C. B. Worden (1999) TriNet ``ShakeMaps'': Rapid Generation of Instrumental Ground Motion and Intensity Maps for Earthquakes in Southern California, Earthquake Spectra, 15, 537-555.

Wald, D. J., and T. H. Heaton (1994). Spatial and temporal distribution of slip for the 1992 Landers, California, earthquake, Bull. Seism. Soc. Am., 84, 668-691.

Wells, D. L., and K. J. Coppersmith (1994) New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement, Bull. Seism. Soc. Am., 84, 974-1002.