1 Introduction

The Japan Trench is one of the great earthquake generating regions in the world. Due to the fast convergence between the Pacific and Eurasian plates (DeMets et al. 1994) and strong interplate coupling (e.g., Suwa et al. 2006; Loveless and Meade 2011; Kanda et al. 2013), large (Mw7.0–8.2) earthquakes have occurred regularly in the past 150 years (e.g., Kanamori et al. 2006). On March 11, 2011, the Mw 9.0 Tohoku-Oki earthquake ruptured the central part of the Japan Trench (Fig. 1).

Fig. 1
figure 1

Focal mechanisms of 69 aftershocks plotted at the epicenters and superposed on a slip model (Wei et al. 2012). The aftershocks are separated into two groups: (blue) aftershocks likely related to thrust faulting on the plate interface and (green) aftershocks not consistent with thrust faulting on the plate interface. The red line is the Japan Trench and black lines are contours of the plate interface (Simons et al. 2011). The star indicates the epicenter of the mainshock

In addition to its high magnitude and devastating impact on communities in northern Honshu, the Tohoku-Oki earthquake had surprising seismo-tectonic characteristics (see Lay and Kanamori 2011 for a review). Geodetic (e.g., Ozawa et al. 2011; Simons et al. 2011) and seismic inversions (e.g., Yoshida et al. 2011) indicate that the co-seismic slip was concentrated up-dip of the hypocenter and had a peak value of about 50 m. The aftershocks are distributed over a 500 × 100 km2 area (e.g., Hayes 2011) and include a variety of thrust-, normal- and strike-slip-faulting events. The diversity of focal mechanisms was not observed before the mainshock (e.g., Nettles et al. 2011) and suggests a complex redistribution of stress during the mainshock slip. To link the aftershock sequence to coseismic stress changes during the mainshock, it is critical to quantify the location (both depth and epicenter) of the aftershocks since the mainshock stress changes are highly heterogeneous. In this study, we determine focal depths of 69 large aftershocks (Table 1) using teleseismic P waveform inversion, and we illustrate how variations in inferred location of an aftershock drastically change the manner in which the mainshock influenced the aftershocks.

Table 1 Focal depths and fault plane parameters of 69 earthquakes obtained from teleseismic waveform inversion

2 Hypocenters of Tohoku-Oki aftershocks

According to the Japan Meteorological Agency (JMA) catalog, about 740 aftershocks with magnitudes greater than 5.0 have occurred between March 11, 2011 and June 30, 2013. Nearly one third of them occurred within one day after the Tohoku-Oki mainshock. 7 aftershocks had magnitudes between 7.0 and 8.0, and 102 aftershocks had magnitudes between 6.0 and 7.0.

Although the JMA catalog is based on local and regional seismic data, it is likely that its estimates of Tohoku-Oki aftershock hypocenters are imprecise. All seismic stations are to the west and at distances of several hundred kilometers for most events. Incomplete azimuthal station coverage and relatively large epicentral distances may lead to errors in the focal depth determination of a few tens of kilometers due to trade-offs between location and origin time and the approximated wave speed structure. Global earthquake catalogs based on arrival times such as the preliminary determination of epicenters (PDE) catalog will have similar uncertainties. The global CMT (gCMT) locations are updated PDE locations using long-period waveforms. While epicentral relocations may suffer from the trade-offs with velocity structure, gCMT depths may be better constrained than PDE depths (Ekström et al. 2012).

The uncertainties are apparent when we compare hypocenters from the various catalogs. Figure 2 compares the aftershock locations of 69 Tohoku-Oki aftershocks as reported in the JMA, PDE, and gCMT catalogs and estimated in this study by P waveform analysis (see Sect. 2.1). The method we use does not constrain epicentral location. On average, the focal depth estimates differ by about 10 km but discrepancies larger than 20 km between our estimates and those in the JMA and PDE catalogs are common. Our estimates of focal depth agree best with the gCMT estimates because of the common exploitation of the secondary phases. The Tohoku-Oki aftershock epicenters in the PDE catalog are systematically to the west (by 24 km on average) with respect to the JMA and gCMT epicenters. Aftershock locations, determined using a long-term OBS network to the south of the source region (Shinohara et al. 2011), are consistent with those listed in the JMA catalog. Moreover, the JMA hypocenters of thrust events are closer to the plate interface than the PDE hypocenters. Therefore, we rely on the JMA catalog for the epicenters.

Fig. 2
figure 2

Comparison of hypocenters of 69 aftershocks from various catalogs. Histograms in (a), (b), and (c) show depth difference with respect to number of events. d, e, and f are plots of event ID versus difference in depth, longitude, and latitude, respectively

2.1 P waveform analysis

gCMT source mechanism solutions are based on relatively long-period waveforms including 40-s long-period body-waves. These data include a variety of seismic signals which is ideal for constraining the point-source moment tensor or focal mechanism. However, focal depth is best constrained by the modeling of broadband waveforms for relatively large earthquakes. Following well-established approaches (e.g., Langston and Helmberger 1975; Christensen and Ruff 1985; McCaffrey and Nábělek 1986), we estimate focal depths by modeling the broadband recording of P waves and trailing depth phases pP and sP recorded at teleseismic distances. For shallow (<30 km) earthquakes, P, pP, and sP signals may interfere due to the finite source duration. Therefore, we employ the waveform inversion method, developed by Kikuchi and Kanamori (1982), which involves the matching of complete P waveforms to a synthetic waveform. We use the t* operator with a value of 1 s to model teleseismic P wave attenuation.

We examine vertical-component recordings at teleseismic distances between 30° and 95°. We require that the aftershocks have relatively simple source time functions, and select waveform data in which the P, pP, and sP phases are recorded with high signal-to-noise ratios because the time delays between sP and P and between pP and P are the primary constraints on focal depth. In computing synthetic waveforms, we ignore spatial finiteness and parameterize the moment rate functions using four overlapping triangles, each with a width of 2–3 s wide. Therefore, source time durations can be up to 8 s long.

We calculate synthetic waveforms using the Jeffreys–Bullen model for the receiver regions and a two-layer velocity and density model for the crust and mantle structure of the Japan Trench region, based on wide-angle reflection and refraction studies in the region (Miura et al. 2005). In layer 1 (the crust), VP = 6.6 km/s, VS = 3.8 km/s, and ρ = 2.87 kg/m3. In layer 2 (the mantle) VP = 8.0 km/s, VS = 4.6 km/s, and ρ = 3.30 kg/m3. The interface between layers 1 and 2 is at 60 km depth. Obviously, a two-layered model does not represent the complex plate boundary structure of Japan Trench region. However, we have found by experimentation that models with multiple crustal layers (though still 1D) have negligible effects on the teleseismic P wave synthetics. The relative times between P, pP, and sP are primarily controlled by the average wave speeds. In addition, we have found that the presence of an ocean is not important since most of the pP and sP energy reflects off the sea floor. Therefore, the focal depths determined in this study are depths below the sea floor.

The accuracy of earthquake depth estimates is controlled by numerous factors, including the effects of the 3D crust and mantle structure, the point-source approximation, source parameters trade-offs, and data quality. It is beyond the scope of this paper to explore these effects in detail. Changes in the average wave speed of 5 %, representing the laterally heterogeneous seismic structure, translates into errors in focal depth of 2 km. We estimate that the focal depths are accurate to within ±4 km given by our experimentation with the 1D velocity model and visual inspection of waveform fits.

In the analysis, we estimate the focal depth by finding the optimal waveform match between band-pass filtered (0.03–1 Hz) recordings, and synthetics that includes P and 30 s of its coda (Fig. 3). We apply the analysis to every 0.5 km depths between 2 and 60 km. The strike and dip of the fault plane are obtained from the gCMT catalog and are held fixed because P waveform data alone poorly constrain the source mechanism. We invert the waveforms to obtain the optimal focal depth, fault rake, and the moment rate function. The focal depths listed in Table 1 are determined by overall waveform fit and by the match to the pP-P and sP-P travel time differences.

Fig. 3
figure 3

Record section of vertical-component displacement seismograms for event 34 showing P (at time = 0) and the depth phase pP and sP. Waveforms are band-pass-filtered from 0.03 to 1 Hz. The phases -pP and -sP are near-source surface reflections above the hypocenter and have similar slowness and propagation paths through the mantle as P (see Bai et al. 2006 for the identification of the depth phases). Red dotted lines show the window used for the waveform inversions

Figure 4 shows examples of waveform inversion results for aftershocks 1, 6, and 31. Event 1 is a normal-faulting earthquake located close to the plate interface. The dip angle of the west-dipping nodal plane is similar to the dip of the mainshock. Event 6 is a reverse-faulting earthquake with a mechanism similar to the mainshock. However, this event is much shallower than the plate interface and the discrepancy between the gCMT’s and this study’s focal depth estimate is the largest among all 69 aftershocks of Table 1. Event 31 is a thrust-faulting event on the plate interface. This event is deep with clearly separated pP and sP phases from P.

Fig. 4
figure 4

Modeling result for events 1 (upper panel), event 6 (middle panel), and event 31 (lower panel). The waveforms are (from top to bottom) the raw data, the synthetics at the preferred depth, and the synthetics at the depth listed in the gCMT catalog. The numbers in parentheses following the depth are least-squares misfit between observed and synthetic seismograms

3 Results

We analyze available waveform data for 69 aftershocks with moment magnitudes between 5.4 and 7.0 using the global seismic network waveform archived from the Incorporated Research Institutes for Seismology (IRIS) Data Management Center (DMC). Waveform data for many of the aftershocks on March 11 were either of poor quality or not available at the IRIS/DMC. We have therefore chosen aftershocks beginning on March 12. The focal depths (ranging from 2 to 57 km) (Fig. 5) and source parameters of the 69 earthquakes investigated here are given in Table 1.

Fig. 5
figure 5

Focal depths of the 69 aftershocks expressed as height above the slab interface. The aftershocks are separated into two groups: shallow-angle thrust events (blue) and non-thrust events (green). Focal mechanisms are plotted using a lower hemisphere projection. See the main text for the definition

Figure 1 shows the focal mechanisms on a contour map of the slab interface. We relate 35 “thrust events” (blue beach-balls) to thrust faulting on the plate interface between 10 and 25 km depth (18 events) and between 35 and 55 km depth (17 events). The dip angle for these two groups increases from about 15° to about 25° with depth, in agreement with the steepening of the slab as mapped by wide-angle seismic imaging (e.g., Ito et al. 2004; Miura et al. 2005). Almost all of these thrust events are outside the region of large mainshock slip (e.g., Ozawa et al. 2011; Simons et al. 2011; Yoshida et al. 2011) and none are up-dip from the mainshock. The thrust events are all within 10 km from the slab interface model, except for the event 61, which has a large difference in longitude between the JMA and PDE catalogs.

We define 34 “non-thrust” events as aftershocks that appear inconsistent with thrust faulting on the plate interface (green beach-balls). These events have normal or strike-slip faulting mechanisms, or they have reverse faulting mechanisms on or near the slab interface with anomalous slip directions. Four events (16, 37, 44, and 62) are normal-faulting outer-rise events that are located seaward of the Japan Trench. There is a remarkable subset of events located within 100 km southeast from the mainshock, which have been used to argue for a total stress release (Hasegawa et al. 2011) or a dynamic overshoot (Ide et al. 2011) in the mainshock. The non-thrust events vary from 60 km above to 40 km below the slab interface. Some (e.g., event 1 and 9) are on the slab interface.

The locations of several of the 69 aftershocks have been reported previously. Huang and Zhao (2013) have relocated large earthquakes (M > 6) from March 9, 2011 to December 31, 2011 using P and S arrival times and three different velocity models. The focal depths of the 36 events in common with Table 1 differ ~7 km on average with a maximum value of 19 km for event 27. Nakajima et al. (2011) estimated event 27 to be ~18 km deeper than the plate interface (~10 km deeper than the plate interface in this study) and argues this to be seismic evidence for reactivation of a buried hydrated fault in the Pacific slab. Lay et al. (2013) estimated the event 65 consisted of three subevents, with the first one at 52.5 km depth (44.5 km in this study) and regarded it as indicating compressional stress accumulating at the uppermost mantle of the Pacific plate around the outer-rise regions.

4 Relationship between mainshock stress change and aftershocks

Several of the normal aftershocks were located west of the trench in an area that is expected to be in horizontal, trench perpendicular compression (e.g., Ide et al. 2011; Hasegawa et al. 2011; Kato et al. 2011; Imanishi et al. 2012). The state of stress that leads to an earthquake can be estimated from a focal mechanism (e.g., Gephart and Forsyth 1984; Arnold and Townend 2007), as well as from more detailed observations or models of coseismic slip (e.g., Angelier 1979; Medina Luna and Hetland 2013), assuming that the direction that the fault slipped during the aftershock (i.e., the slip rake on either of the nodal planes) is parallel with the direction of maximum shear stress on the fault plane prior to the earthquake rupture. In the case of an aftershock sequence, the stresses that lead to the aftershocks might be due to both the background stress and the change in stress during the mainshock. For brevity, we refer to the change in stress during the Tohoku-Oki mainshock as the mainshock stress change (MSC). Due to the large and concentrated coseismic fault slip in the Tohoku-Oki mainshock (e.g., Simons et al. 2011; Wei et al. 2012), the MSC could be on order of 10 MPa near the mainshock, although the magnitude and orientations of the mainshock stress changes are highly heterogeneous. Due to the heterogeneities in MSC, interpretations of aftershock focal mechanisms in context of the Tohoku-Oki mainshock rely on accurate estimates of aftershock hypocenters. Slight changes in the aftershock location relative to the coseismic slip in the mainshock, could result in significant changes in the MSC, and thus different interpretations of the aftershock focal mechanism in relation to the stress change in the mainshock. For example, to the west of the Tohoku-Oki earthquake, the MSC results in potentially significant unloading of the trench-perpendicular compressive stresses at shallower depths, while resulting in potentially significant loading of the trench-perpendicular compressive stresses at greater depths.

In the Tohoku-Oki aftershock sequence, there were several aftershocks with unexpected mechanisms, including several normal-faulting aftershocks located near the mainshock epicenter. Ide et al. (2011) recognized that the dip of one of the nodal planes of several of these aftershocks is similar to the dip of the megathrust. They postulated that these aftershocks are evidence of dynamic overshoot during the mainshock. This implies that the accumulated reverse-sense shear stress on the megathrust was not only completely relieved during the mainshock, but that the final shear stress on the megathrust was normal-sense immediately following the mainshock. Assuming that the stress drop in the aftershocks is about 1–10 MPa (e.g., Kanamori and Anderson 1975), the dynamic overshoot would need to be similarly large to cause normal coseismic slip on the megathrust. This amount of dynamic overshoot could be a substantial fraction of the total stress drop in the mainshock, even assuming that the stress drop in the Tohoku-Oki earthquake may have been as large as 50 MPa (Kanda et al. 2013). Such large normal-sense shear would also likely lead to normal-sense afterslip on the megathrust, which has not been observed to date (e.g., Johnson et al. 2012; Uchida and Matsuzawa 2013).

In contrast to the dynamic overshoot interpretation of Ide et al. (2011), Hasegawa et al. (2011) concluded that the normal-faulting aftershocks near the Tohoku-Oki mainshock merely indicate near total stress release in the mainshock. This conclusion was based on the two-dimensional theory developed by Hardebeck and Hauksson (2001), which assumes that the aftershock occurs in close proximity to the region of maximum stress release in the mainshock and that the slip occurred along a fault parallel to the fault that ruptured in the mainshock. Because the normal aftershocks are indicative of trench perpendicular extension after the mainshock, whereas the stress prior to the mainshock was most likely in trench perpendicular compression, this large rotation of the stress field indicated near complete stress drop in the mainshock (Hasegawa et al. 2011). We find that many of these normal mechanism aftershocks near the Tohoku-Oki epicenter are above or below the megathrust by more than our estimated depth uncertainty of ±4 km (Fig. 5). For instance, events 1 and 2 are both normal-mechanism aftershocks to the south of the mainshock, and we infer these two events to be up to 10 km off of the megathrust (Fig. 5). We also note that the shallower of the two nodal planes was not necessarily the fault plane, and it is equally likely that the steeper of the two nodal planes in these normal aftershocks was the fault plane.

For the aftershocks near the main coseismic slip region of the Tohoku-Oki mainshock, the direction in which the mainshock loads either of the nodal planes (the plausible fault surfaces) depends strongly on the precise location of that aftershock. We illustrate this point for event 1 by computing the shear stress change on each of the two nodal planes due to the MSC, assuming various centroid locations of this aftershock (Fig. 6). We vary the epicenter of the aftershock between the PDE and JMA determined epicenters, where for event 1, the JMA epicenter is further offshore compared to the PDE epicenter (Fig. 2). We depict the sense of shear stress change by plotting focal mechanisms for the expected aftershock mechanisms under the assumptions that (1) only the MSC leads to the fault failure and (2) the fault slips in the direction of the maximum shear stress. We assume the gCMT’s nodal plane strike and dip, and determine the rake by computing the direction of the maximum shear stress change on each of the nodal planes due to MSC. In other words, the calculated rake is the direction in which the Tohoku-Oki earthquake loaded each of the nodal planes. We then assume, a pure double-couple to compute the orientation and rake of the second nodal plane. To calculate the MSC, we use the median coseismic slip model of Simons et al. (2011) and use the dislocation algorithm of Meade (2007) in a Poisson half-space. More recent coseismic slip models of the Tohoku-Oki earthquake include larger coseismic slip closer to the trench than in the Simons et al. (2011) model (e.g., Wei et al. 2012), and different coseismic slip models would predict different MSC at the trial locations of the aftershock shown in Fig. 6. The aftershock, we use to illustrate the high spatial variability of the MSC is about 50 km to the south of the mainshock, and the degree of variability would be similar to what is presented here for other coseismic slip models. The sense that each of the nodal planes was loaded during the mainshock is highly dependent on the precise location of the aftershock (Fig. 6). In general, assuming the coseismic slip model of Simons et al. (2011), the Tohoku-Oki mainshock loaded the nodal planes in a reverse-sense at depths significantly shallower than the preferred depth we determine here. At the preferred depth we determine, as well as deeper depths, the sense that the mainshock loaded the nodal planes depends on epicentral location (Fig. 6). For locations closer to the PDE determined epicenter, the nodal planes are loaded in an oblique-normal to strike-slip sense. For locations closer to the JMA epicenter, the nodal planes are loaded in a more normal-sense. Particularly, the steeper dipping nodal plane is loaded in almost the same direction as the gCMT determined slip rake for a broad range of depths near our preferred depth and closer to the JMA epicenter (Fig. 6).

Fig. 6
figure 6

Depiction of the stress change during the Tohoku-Oki mainshock at various locations corresponding to probable centroid locations of aftershock 1. The left and right columns of both panels are assuming epicenters as determined by PDE and JMA, respectively, while the middle two columns are assuming epicentral locations evenly spaced between the PDE and JMA determined epicenters. The stress change is depicted using a focal mechanism, where the strike and dip of each of the indicated nodal plane is from the gCMT solution and the slip rake on that nodal plane is calculated as the direction of the maximum shear stress change on that fault during the mainshock (i.e., the focal mechanism if only the coseismic stress change caused fault slip). Centroid locations are shown relative to the preferred depth determined in this study, and with epicenters varying between the PDE and JMA determined epicenters. The calculated direction of the maximum shear stress is indicated below each beach-ball, and color of the beach-ball corresponds to how close that direction is to the slip rake on that nodal plane in the gCMT solution (blue corresponds to those cases where the mainshock loaded the nodal plane in the same direction that it slipped, and red corresponds to those cases where the mainshock loaded the nodal plane opposite to the direction that it slipped)

5 Discussion and conclusions

All of the 69 aftershocks are shallower than 60 km, which is the lower limit of the seismogenic zone along the subducting plate in this region (Igarashi et al. 2001). The 35 thrust events are within 10 km of the slab interface model as modeled by Kanda et al. (2013) and can be interpreted as thrust faulting on the megathrust. Although we estimate uncertainties in focal depth determination to be ±4 km, it is possible that the slight offset is due to mislocated epicenters or variations in the megathrust geometry that are not accounted for in the larger scale megathrust models (Zhan et al. 2012).

The non-thrust aftershocks are up to 60 km from the megathrust, which is significantly more than errors in the focal depth or epicenters. It is likely that events 27, 35, 39, 45, 60, 64, and 65 occurred within the subducting plate and that events 5, 6, 11, 28, 29, 30, 52, 54, and 55 occurred within the overriding plate.

Earthquakes near the trench are usually shallow and associated with compressional or extensional stress regimes in the upper part of the slab due to the slab bending (e.g., Christensen and Ruff 1983). Three aftershocks (event 39, 64, and 65) near the trench were relatively deep. Event 65 had a high-angle reverse fault and located ~38 km below the plate interface. Events like event 65 may play an important role in the hydration of the uppermost mantle of the Pacific plate prior to the subduction.

The change in stress during a large earthquake can result in changes to both the shear and normal stresses on nearby faults (e.g., King et al. 1994). Assuming that faults slip in the direction of the maximum shear stress resolved on the fault plane (e.g., Angelier 1979) and that the MSC is large enough to affect the mechanisms of the aftershocks, there are three interpretations of the stresses that lead to an aftershock in light of the MSC. (In these interpretations we do not focus on the change in fault normal stress due to the MSC, rather we assume that the changes in fault normal stress are conducive to fault failure in the aftershocks and we only focus on the affect of shear stress changes). The first is that the change in shear stress on the fault due to the MSC is roughly collinear with the shear stress on the fault prior to the mainshock (case 1 in Fig. 7). In this case, the mainshock effectively pushes the fault closer to failure, for which slip occurs in the same sense that it had been loaded in due to the background tectonic stresses. This is unlikely for the non-thrust aftershocks. Except for the outer-rise aftershocks, all non-thrust aftershocks are located in regions where the stresses prior to the mainshock were likely near-trench-perpendicular compression (e.g., Hasegawa et al. 2011; Kanda et al. 2013). The second case is when the change in shear stress on the aftershock fault plane is much larger than the shear stresses that were accumulated on that fault prior to the mainshock (case 2 in Fig. 7). In this case, the MSC largely determines the sense of fault failure during the aftershock. This case may be true for shallow aftershocks involving failure on steeply dipping faults, where if prior to the mainshock the region was in near horizontal compression for which there would be little shear stress resolved on that fault plane prior to the mainshock. The third case is when the magnitude of shear stress change during the mainshock is comparable to the magnitude of shear stress on the fault prior to the mainshock (case 3 in Fig. 7). In this case, the fault will fail in a direction that is neither the direction in which the fault was loaded prior to the mainshock nor the direction of the shear stress change during the mainshock. For instance, if a fault was loaded in a reverse-slip sense prior to the mainshock, but not at a level in which it would fail, and the shear stress change during the mainshock was in a strike-slip sense, assuming that the magnitude of the resulting shear stress was sufficient for failure, the fault would then slip in an oblique-sense. This last case is a likely scenario for the Tohoku-Oki aftershock sequence. Precise locations of the aftershocks are therefore crucial to interpret those aftershocks in light of the stress changes during the Tohoku-Oki mainshock.

Fig. 7
figure 7

Cartoons illustrating three cases of how shear stress changes on a fault due to the mainshock might affect coseismic slip direction in an aftershock