Synthetic sensitivity analysis of high frequency radiation of 2011 Tohoku-Oki (M_W 9.0) earthquake

1.   Introduction

Frequency-dependent rupture process of the 2011 M_W 9.0 Tohoku-Oki earthquake (e.g., Lay et al., 2012) has been examined using multi-frequency (0.05-5 Hz) back-projection imaging (e.g., Koper et al. 2011a, b; Kiser and Ishii 2011; Wang and Mori 2011; Ishii 2011; Huang et al., 2012; Meng et al, 2012) and compressive sensing studies (e.g., Yao et al., 2011), inferred to correspond to depth-varying rupture properties. This interpretation was corroborated by analyzing seismic source spectra of interplate megathrust events in frequency bands up to 10 Hz. (Ye et al., 2013). Analyses of strong motion data also indicate that higher frequency (> 10 Hz) energy was radiated strongly at the down-dip region rather than at the up-dip region of the Tohoku-Oki megathrust fault (e.g., Ide et al. 2011; Iwakiri and Hoshiba 2012). Investigations of this phenomenon are important both for understanding of variations in basic source physics across the megathrust and the distance decay of strong shaking.

Even though the robust results of teleseismic back-projection imaging and source spectra analyses strongly indicate frequency-dependent rupture process using data with relative low frequency bands (< 10 Hz), the results of analyzing near-field data with higher frequency band could provide supplement about the intrinsic differences in rupture behavior across the megathrust, so as to help us improving dynamic source models. Therefore, it is significant to carefully examine whether using near-field data could solidly support the conclusion of back-projection. However, due to poor station coverage of near field observation (Fig. 1a) and geometry of the low angle interplate megathrust fault (Fig. 1c), other important factors, such as geometrical and anelastic attenuation as well as diffraction of seismic waves should also be taken into account in these analyses.

Figure  1.  a (Left) Strong motion records high-pass filtered above 10 Hz at nine stations distributed along the rupture. The maximum amplitude in centimeters per second per second is shown at the top-right of each waveform. (Right) Slip distribution. The S-wave arrival times from several source locations, denoted with colored stars and with time after the rupture initiation also shown, are plotted on the waveforms with length inversely proportional to the source-station distance, modified from (Ide et al. 2011). b Mesh of the rupture region, 16 × 8 model blocks of the central positions of 30 km × 30 km patches. i and j increase along strike and along dip separately. c Schematic cross-section of interplate megathrust fault

DownLoad: Full-Size Img PowerPoint

Because of the geometrical settings (Fig. 1c; Table 1) of the Japan Trench megathrust and the stations on shore, high frequency radiations from the shallower part of the rupture attenuate much more along the path compared with that from deeper parts (Fig. 1c; Table 1). Rupture inversion results, however, showed that slip of the Tohoku earthquake (Fig. 2), increases to maximum from the down-dip region to trench (e.g., Ide et al. 2011; e.g., Yue and Lay 2013). The tradeoff between attenuation and slip distribution makes it less certain about distribution pattern of high frequency radiation based only on simple frequency analysis of the data observed. It requires carefully quantitative analysis of the data taking into account also the wave propagation geometry and attenuation effects.

Table  1.  Geometry of the interplate megathrust fault

 | Show Table

DownLoad: CSV

Figure  2.  a Inverted rupture pattern on the fault plane for the Tohoku-Oki earthquake, modified from (Yue and Lay 2013). b Moment density distribution calculated from (a) and PREM model. Value of distance along dip increases from the trench to down-dip

DownLoad: Full-Size Img PowerPoint

Using a finite fault model (e.g., Yue and Lay 2013) and synthetic seismograms, we have performed a synthetic sensitivity analysis to investigate whether data recorded from stations along the east coast of Honshu island could be analyzed to resolve the possible high frequency (> 10 Hz) radiation from the up-dip region of the megathrust. By placing testing receivers in a 2D region from the island to the trench, we finally present the result of detectability of high frequency radiations (> 10 Hz) in the shallower parts of the rupture around the Tohoku region.

2.   Methods

To test the sensitivity of high frequency radiation, we parameterize the source region into 16 by 8 rectangle fault patches (Fig. 1b). Accelerograms of k-th station A_{(i, j)}^k(t) from source patch (i, j) can be expressed as the convolution of Green's function G_{(i, j)}^k(t), source time function M_{(i, j)}S(t) and instrumental response I(t):

where M_{(i, j)} denotes the moment of each fault patch, S(t) is the is the Heaviside function to ensure production of abundant high frequency energy. We simply assume the same S(t) of each fault patch in order to investigate only the propagation effects.

In this study we use the finite-fault slip model of (Yue and Lay 2013) derived from a joint inversion of high-rate geodesy and seismic data. This is based on the consideration that (1) geodetic and seismic data are jointly inverted in solution, yielding better slip estimation and resolution than other inversions using one type of data only. (2) This estimated slip distribution pattern is similar to that of (Ide et al. 2011), what we will have our results to compare with. To generate synthetic seismograms with sampling rate of 40 Hz and duration of about 200 s, we adopt the wavenumber integration method because of its higher efficiency and broader frequency bands than spectral element or finite difference method.

For simplicity and feasibility, we only consider layered homogeneous media following the PREM model between the source region and receivers. Due to the lack of data constraints, we do not consider lateral variations of Q structure. Although Q_S/Q_P generally ranges from 1 to 3 for frequencies higher than 1 Hz (e.g., Sato et al., 2013), after propagating from source to receiver, the dominant energy of seismic wave mostly depends on shear wave or surface wave, while the effect of Q_P is negligible. As a result, we simply assume a homogeneous Q structure where Q = Q_S = Q_P, independent of frequency. By varying Q, we obtain different Green's function G_{(i, j)}^k(Q) of media.

We then high pass filter the accelerograms above 10 Hz and integrate the k-th station record from the fault patch (i, j) into its L2 norm E_{(i, j)}^k, i.e.,

where T_b and T_e denotes the begin and end time of the record of an earthquake. According to our synthetic test, the E_{(i, j)}^k of vertical component is more sensitive to high frequency energy radiated from the shallower parts of the megathrust fault (Fig. 3) than the horizontal components. Plotting E_{(i, j)}^k onto the fault plane, we can examine at each station the observed high frequency energy distribution (ORHFED) pattern on the fault plane (Fig. 4).

Figure  3.  Synthetic test of received high frequency energy distribution with uniform slip of the megathrust fault when Q is infinite, station is set just above the southwest corner of the megathrust fault. a Vertical component, b Horizontal component, c Total. Apparently, the r_sum^k of vertical component is the highest so as to be the most sensitive component. Value of distance along dip increases from the trench to down-dip

DownLoad: Full-Size Img PowerPoint

Figure  4.  ORHFED patterns with different Q structures of station IWTH14 (a), IWTH27 (b), MYGH08 (c), IBRH14 (d), CHBH14 (e), north to south, respectively. E_{(i, j)}^k of a certain fault patch, which has been normalized from 0 to 1, represents the high frequency energy received by station k rather than the energy radiated from that fault patch. Value of distance along dip increases from the trench to down-dip

DownLoad: Full-Size Img PowerPoint

3.   Discussion

We calculate the ORHFED patterns for stations along the east coast of Honshu island shown in Fig. 1a, with the Q value ranges from 200 to 1, 000 (Fig. 4). Generally speaking, when Q < 600, received energy from the down-dip region of the megathrust dominates, while the received energy form the up-dip region area begins to be comparable with the energy from the down-dip region area as Q exceeding 600 for northern stations along the coast. The general trend indicates that, for the tradeoff between attenuation and slip distribution, the latter begins to dominate as Q increases. To better evaluate the tradeoff, we introduce a term r_sum^k, i.e.,

where r_sum^k denotes the ratio between the sum of energy E_{(i, j)}^k from up-dip region to that of down-dip region. As Q increases, r_sum^k increases steadily (Fig. 5). When r_sum^k is less than the cutoff ratio of 0.1, the total energy received by k-station from the up-dip region is at least an order of magnitude less than that from the down-dip region. In other words, we propose that when r_sum^k is less than the cutoff ratio of 0.1, the high frequency energy from the up-dip region is unresolvable in a complicate megathrust source rupture process because seismic phases overlap in time domain. After examining of the r_sum^k curves of stations shown in Fig. 1a, we conclude that whether these stations along the coastline can resolve the high frequency energy radiated from up-dip region depends on the Q value in our model. Generally speaking, when Q is lower than 500, the high frequency energy radiated from the shallower parts is unresolvable using data from stations located along the east coastline of the Honshu island due to the dominance of attenuation effect, irrespective of the amount of high frequency energy radiated from the maximum slip in the up-dip region area.

Figure  5.  The relation between r_sum^k and Q value. IWTH14 to CHBH14 are results of different stations from north to south, respectively. The black dashed line denotes the cut ratio 0.1

DownLoad: Full-Size Img PowerPoint

Technically, non-geometrical seismic wave attenuation is usually caused by two mechanisms, scattering due to lateral heterogeneity and intrinsic attenuation. Unlike intrinsic attenuation which is mainly associated with anelastic processes, scattering redistributes energy instead of causing energy loss per cycle (e.g., Sato et al., 2013). In this study, we use wavenumber integration in a layered homogeneous media to generate Green's function without taking scattering into consideration. Therefore, for a certain value of Q, we may overestimate E_{(i, j)}^k because the L2 norms of the seismic records decrease considerably after shifting energy from direct arrivals back to the coda, and so does the r_sum^k. Since only anelastic attenuation is considered in our model, it actually represents less attenuation effect in real media.

Q_S^-1 generally takes the value on the order of 0.01 around 1 Hz and on the order of 0.001 at 20 Hz following the frequency-dependent relationship of Q_S^-1 ∝ f^-n, where the exponent n ranges from 0.0 to 1.0 and varies regionally, with depth in the crust (e.g., Sato et al., 2013; Romanowicz and Mitchell 2007). Q_S dominates in our model because of the amplitude of S wave or surface wave is much larger than P wave. Especially in calculating the L2 norms of the records, the E_{(i, j)}^k affected by Q_P is negligible.

Hong (2010) investigated the Lg attenuation in east Asia with both continental and oceanic environments. Their results suggest significantly low quality factor at 1 Hz of 100 or less in the Sea of Japan including our study area. As a result, attenuation of shear wave is strong in crust. The distributed cracks caused by the high seismicity of the hanging wall result in highly intrinsic friction and scattering, which probably is the main reason for high attenuation. Together with the general attenuation behavior of the crust, we would expect significant low Q_S (probably less than 500) in 10 to 20 Hz in this region even through the average frequency dependency parameter n is still unclear. As a result, the stations along the east coast of Honshu Island probably could not resolve the high frequency (> 10 Hz) radiation from the up-dip region due to the dominance of attenuation.

For further discussion, we calculate as the receiver is placed in a 2D region from the island to the trench, where . In the plot regions with lower than -1, the high frequency energy radiated from the up-dip region of the megathrust fault cannot be resolved (Fig. 6). As Q increases, the area with lower than -1 shrinks. Received high frequency energy from the shallower parts is not resolvable by analyzing data recorded by stations on Honshu island when Q < 500. As Q increases, the high frequency energy radiated from the shallower parts becomes resolvable at greater regions for the receiver station to be located at, moving from the east coast toward the west coast on the island.

Figure  6.  Distribution of testing stations. b-k The maps for different Q values, where denotes the logarithm of r_sum^k. If the region has a lower than -1, station in this region could not resolve the high frequency energy radiated form the shallower parts of the megathrust fault. The write contours denote the slip distribution of Tohoku-Oki earthquake inverted by (Yue and Lay 2013)

DownLoad: Full-Size Img PowerPoint

Spectra of moderate earthquakes from various depths of the megathrust fault are usually analyzed to study the regional difference of the source radiation. (Iwakiri and Hoshiba 2012) analyze the frequency-dependent trend of spectral ratio (up to 20 Hz) has been analyzed to indicate the possible different source time functions of interplate earthquakes located from the up-dip region to down-dip region. However, the main factors in affecting the wave spectra include both regional difference in the source radiation and the path effects. (e.g., Iwakiri and Hoshiba 2012). (Ye et al. 2013) use empirical Green's function (EGF) method to suppress propagation effects. The fairly reliable source spectra can only be extracted up to 5 Hz by carefully deconvoluting the EGF events because of the clearly increasing noise levels. In addition, the energy of surface waves could be much stronger rather than body waves for shallower events (Fig. 7). While calculating the spectral ratio of an event from the down-dip region to one from the up-dip region, surface wave could decrease the value in relative low frequency band compared to high frequency band, which could result in the biased estimation that deeper events radiate more strongly at high frequency. Here, we set six hypothetical testing events into two groups from up-dip region to down-dip region of the megathrust fault (Table 2). With a radius of 200, 300 or 400 km, we let station evenly distribute along a circle with azimuths range from 202° to 22° at a 10° interval. We then sum up the spectrums of signals at the stations for a window of 100-s-long and calculate the ratio of the spectra of event pair E4 & E1, E5 & E2, and E6 & E3. Figure 8 indicates the spectral ratios with different direct distances and Q structures in frequency band from 0.2 to 20 Hz. Despite the same source time function for each event pair, the spectral ratios increase considerably with frequency, indicating that lower spectral ratios in low frequency mainly result from the dramatic amplitude difference of surface wave vibrated with sources in different focal depths. Overall, while using spectral ratios to analyze the source properties, frequency-dependent trend of spectral ratio could be explained by the combination of difference in source time functions, propagation effects and the strong bias introduced by surface waves. Extracting reliable broadband source spectra (up to 20 Hz) of the interplate events requires high quality near-fault records or analysis of the high-resolution attenuation structure for the high-frequency range (> 10 Hz). One of the most challenging missions of seismic observation is to overcome the difficulty of maintaining an ocean bottom seismic array from the island to the trench.

Figure  7.  Synthetic accelerograms from event pair Event 1 & Event 4 (Table 2). The ratio of maximum amplitude of Event 1 to that of Event 4 is 6.84. Both accelerograms are calculated with Q = 500, epicentral distance = 300 km and azimuth = 252°, respectively. P and S arrivals are marked by light gray arrows

DownLoad: Full-Size Img PowerPoint

Table  2.  Parameters of synthetic events and test stations

 | Show Table

DownLoad: CSV

Figure  8.  Spectral ratios of testing events from deeper parts and those from shallower parts of the megathrust fault. a Ratios of events pairs E4 & E1, E5 & E2, and E6 & E3 (Table 2.), direct distance from event to station is 200 km, spectra are calculated from 50 or 100-s-long signals. b and c are same as (a), but direct distance from event to station are 300 and 400 km. Black lines denote model with Q = 200, whereas the gray lines demote model with Q = 500

DownLoad: Full-Size Img PowerPoint

According to our result (Fig. 6), hypothetical stations close to the trench are extremely sensitive to the high frequency radiated from the up-dip region. OBS (ocean bottom seismometers) data near the trench have significant advantages in analyzing the depth-varying rupture process. It could be verified by broadband OBS data in the future that whether the up-dip region area could generate as much high frequency (> 10 Hz) radiation as down-dip region does.

4.   Conclusions

The results of our synthetic tests indicate that the possible high frequency (> 10 Hz) energy radiated from the up-dip region could hardly be resolved by stations because of the high attenuation of this region, whereas the OBS stations located near the megathrust fault can detect the up-dip region generating high frequency radiation. In addition, when analyzing spectral ratios of moderate earthquakes from various depths of the megathrust fault, one should be careful about the bias introduced by the strong low frequency surface wave of earthquakes from shallower parts of the fault, which could result in erroneous interpretation that deeper events radiate more strongly at high frequency.

Acknowledgments

This work made use of GMT, SAC. and COMPUTER PROGRAMS IN SEISMOLOGY software. We thank Professor Z. Shen, Y. Zhou, J. Ning, Z. Ge, Haiming Zhang and two anonymous reviewers for providing thoughtful discussions and comments. This work is supported by NSFC grants (Nos. 91128210 and 40821062).