ALBERT

Description: High-melt areas of glaciers and ice sheets foster a rich spectrum of ambient seismicity. These signals not only shed light on source mechanisms (e.g. englacial fracturing, water flow, iceberg detachment, basal motion) but also carry information about seismic wave propagation within glacier ice. Here, we present two approaches to measure and potentially monitor phase velocities of high-frequency seismic waves (≥1 Hz) using naturally occurring glacier seismicity. These two approaches were developed for data recorded by on-ice seasonal seismic networks on the Greenland Ice Sheet and a Swiss Alpine glacier. The Greenland data set consists of continuous seismograms, dominated by long-term tremor-like signals of englacial water flow, whereas the Alpine data were collected in triggered mode producing 1–2 s long records that include fracture events within the ice (‘icequakes’). We use a matched-field processing technique to retrieve frequency-dependent phase velocity measurements for the Greenland data. In principle, this phase dispersion relationship can be inverted for ice sheet thickness and bed properties. For these Greenland data, inversion of the dispersion curve yields a bedrock depth of 541 m, which may be too small by as much as 35 per cent. We suggest that the discrepancy is due to lateral changes in ice sheet depth and bed properties beneath the network, which may cause unaccounted mixing of surface wave modes in the dispersion curve. The Swiss Alpine icequake records, on the other hand, allow for reconstruction of the impulse response between two seismometers. The direct and scattered wave fields from the vast numbers of icequake records (tens of thousands per month) can be used to measure small changes in englacial velocities and thus monitor structural changes within the ice.

Description: Automatic classification of local seismic events which are only recorded at single stations poses great challenges because of weak hypocentre constraints. This study investigates how single-station event clusters relate to geographic hypocentre regions and common source processes. Typical applications arise in local seismic networks where reliable ground truth by a dense temporal network precedes or follows a sparse (permanent) installation. The seismic signals for this study comprise a 3-month subset from a field campaign to map subduction below northern Chile (PISCO ’94). Due to favourable ground noise conditions in the Atacama desert, the data set contains an abundance of shallow and deeper earthquakes, and many quarry explosions. Often event signatures overlap, posing a challenge to any signal processing scheme. Pattern recognition must work on reduced seismograms to restrict parameter dimensionality. Continuous parameter extraction based on noise-adapted spectrograms was chosen instead of discrete representation by, for example, amplitudes, onset times or spectral ratios to ensure consideration of potentially hidden features. Visualization of the derived feature vectors for human inspection and template matching algorithms was hereby possible. Because event classes shall comprise earthquake regions regardless of magnitude, clustering based on amplitudes is prevented by proper normalization of feature vectors. Principal component analysis is applied to further reduce the number of features used to train a self-organizing map (SOM). The SOM will topologically arrange prototypes of each event class in a 2-D map. Overcoming the restrictions of this black-box approach, the arranged prototypes could be transformed back to spectrograms to allow for visualization and interpretation of event classes. The final step relates prototypes to ground-truth information, confirming the potential of automated, coarse-grain hypocentre clustering based on single-station seismograms. The approach was tested by a twofold cross-validation whereby multiple sets of feature vectors from half the events are compared by a one-nearest neighbour classifier in combination with an Euclidean distance measure resulting in an overall correct geographic separation rate of 80.5 per cent.

Description: Finite-fault earthquake source inversion is an ill-posed inverse problem leading to non-unique solutions. In addition, various fault parametrizations and input data may have been used by different researchers for the same earthquake. Such variability leads to large intra-event variability in the inferred rupture models. One way to understand this problem is to develop robust metrics to quantify model variability. We propose a Multi Dimensional Scaling (MDS) approach to compare rupture models quantitatively. We consider normalized squared and grey-scale metrics that reflect the variability in the location, intensity and geometry of the source parameters. We test the approach on two-dimensional random fields generated using a von Kármán autocorrelation function and varying its spectral parameters. The spread of points in the MDS solution indicates different levels of model variability. We observe that the normalized squared metric is insensitive to variability of spectral parameters, whereas the grey-scale metric is sensitive to small-scale changes in geometry. From this benchmark, we formulate a similarity scale to rank the rupture models. As case studies, we examine inverted models from the Source Inversion Validation (SIV) exercise and published models of the 2011 Mw 9.0 Tohoku earthquake, allowing us to test our approach for a case with a known reference model and one with an unknown true solution. The normalized squared and grey-scale metrics are respectively sensitive to the overall intensity and the extension of the three classes of slip (very large, large, and low). Additionally, we observe that a three-dimensional MDS configuration is preferable for models with large variability. We also find that the models for the Tohoku earthquake derived from tsunami data and their corresponding predictions cluster with a systematic deviation from other models. We demonstrate the stability of the MDS point-cloud using a number of realizations and jackknife tests, for both the random field and the case studies.

Description: Recently, various methods have been proposed and applied for earthquake source imaging, and theoretical relationships among the methods have been studied. In this study, we make a follow-up theoretical study to better understand the meanings of earthquake source imaging. For imaging problems, the point spread function (PSF) is used to describe the degree of blurring and degradation in an obtained image of a target object as a response of an imaging system. In this study, we formulate PSFs for earthquake source imaging. By calculating the PSFs, we find that waveform source inversion methods remove the effect of the PSF and are free from artefacts. However, the other source imaging methods are affected by the PSF and suffer from the effect of blurring and degradation due to the restricted distribution of receivers. Consequently, careful treatment of the effect is necessary when using the source imaging methods other than waveform inversions. Moreover, the PSF for source imaging is found to have a link with seismic interferometry with the help of the source–receiver reciprocity of Green's functions. In particular, the PSF can be related to Green's function for cases in which receivers are distributed so as to completely surround the sources. Furthermore, the PSF acts as a low-pass filter. Given these considerations, the PSF is quite useful for understanding the physical meaning of earthquake source imaging.

Description: We introduce the Complex-Step-Finite-Difference method (CSFDM) as a generalization of the well-known Finite-Difference method (FDM) for solving the acoustic and elastic wave equations. We have found a direct relationship between modelling the second-order wave equation by the FDM and the first-order wave equation by the CSFDM in 1-D, 2-D and 3-D acoustic media. We present the numerical methodology in order to apply the introduced CSFDM and show an example for wave propagation in simple homogeneous and heterogeneous models. The CSFDM may be implemented as an extension into pre-existing numerical techniques in order to obtain fourth- or sixth-order accurate results with compact three time-level stencils. We compare advantages of imposing various types of initial motion conditions of the CSFDM and demonstrate its higher-order accuracy under the same computational cost and dispersion–dissipation properties. The introduced method can be naturally extended to solve different partial differential equations arising in other fields of science and engineering.

Description: A discontinuous grid finite-difference (FD) method with non-uniform time step Runge–Kutta scheme on curvilinear collocated-grid is developed for seismic wave simulation. We introduce two transition zones: a spatial transition zone and a temporal transition zone, to exchange wavefield across the spatial and temporal discontinuous interfaces. A Gaussian filter is applied to suppress artificial numerical noise caused by down-sampling the wavefield from the finer grid to the coarser grid. We adapt the non-uniform time step Runge–Kutta scheme to a discontinuous grid FD method for further increasing the computational efficiency without losing the accuracy of time marching through the whole simulation region. When the topography is included in the modelling, we carry out the discontinuous grid method on a curvilinear collocated-grid to obtain a sufficiently accurate free-surface boundary condition implementation. Numerical tests show that the proposed method can sufficiently accurately simulate the seismic wave propagation on such grids and significantly reduce the computational resources consumption with respect to regular grids.

Description: SS precursor observations are a powerful tool to study the topography and character of transition zone discontinuities, especially in regions such as ocean basins where few seismic stations exist, precluding other high resolution approaches. Still, the available coverage is limited by the distribution of sources and stations, but also by the level of noise and by the fact that, in some distance ranges, interfering seismic phases mask the weak signal from the SS precursors. We introduce an array data processing tool, the local slant-stack filter, to address these challenges and clean up the otherwise noisy SS precursor record sections. We show that these filters are a powerful tool for extracting the weak yet coherent SS precursor signals while removing interfering seismic phases as well as random noise, yielding robust precursor traveltime measurements with spatial resolution higher than what can be achieved by the conventional common midpoint stacking method. The effectiveness of the filters are demonstrated by application to synthetic and real data. We systematically apply this filtering method to an SS precursor data set recorded by the U.S. Transportable Array that samples a vast region of the Pacific Ocean and its northwest margin, and present maps of 410 and 660 discontinuity topography. We discuss correlations observed between our discontinuity images and several fine-scale heterogeneities revealed by mantle shear wave tomography in the vicinity of Hawaii and the Pacific Superswell.

Description: The properties of the overburden transmission response in seismic depth imaging are of interest for the analysis of reflectivity illumination or blurring. An elastic transmission-operator reciprocity is derived here by a direct method exploiting a symmetry of the elastic impedance operator at the input and output levels. This impedance-operator symmetry is obtained in a companion paper by assuming the existence and completeness of the lateral modes of the elastic wave equation. The transmission operator, or wave propagator, is a surface-to-surface displacement-to-displacement operator and an appendix here explains heuristically why it does not immediately display its reciprocity property in the way that force-to-displacement Green functions do. Instead, the transmission operator must be augmented by some version of the impedance operator at the input and output levels. We may view this augmentation as re-establishing the connection to force-to-displacement reciprocity or as a separate reciprocity relation for displacements which have been normalized by input and output impedances. Similar statements apply to reflection operators for the region between the two levels. The main transmission-operator reciprocity is an exact property, but a second and approximate symmetry exists when waves that are evanescent in depth are neglected.

Description: This work is aimed at the automatic and fast characterization of the extended earthquake source, through the progressive measurement of the P -wave displacement amplitude along the recorded seismograms. We propose a straightforward methodology to quickly characterize the earthquake magnitude and the expected length of the rupture, and to provide an approximate estimate of the average stress drop to be used for Earthquake Early Warning and rapid response purposes. We test the methodology over a wide distance and magnitude range using a massive Japan earthquake, accelerogram data set. Our estimates of moment magnitude, source duration/length and stress drop are consistent with the ones obtained by using other techniques and analysing the whole seismic waveform. In particular, the retrieved source parameters follow a self-similar, constant stress-drop scaling (median value of stress drop = 0.71 MPa). For the M 9.0, 2011 Tohoku-Oki event, both magnitude and length are underestimated, due to limited, available P -wave time window (PTWs) and to the low-frequency cut-off of analysed data. We show that, in a simulated real-time mode, about 1–2 seconds would be required for the source parameter determination of M 4–5 events, 3–10 seconds for M 6–7 and 30–40 s for M 8–8.5. The proposed method can also provide a rapid evaluation of the average slip on the fault plane, which can be used as an additional discriminant for tsunami potential, associated to large magnitude earthquakes occurring offshore.

Description: We present a simple, fast, and robust method for automatic detection of P - and S -wave arrivals using a nearest neighbours-based approach. The nearest neighbour algorithm is one of the most popular time-series classification methods in the data mining community and has been applied to time-series problems in many different domains. Specifically, our method is based on the non-parametric time-series classification method developed by Nikolov. Instead of building a model by estimating parameters from the data, the method uses the data itself to define the model. Potential phase arrivals are identified based on their similarity to a set of reference data consisting of positive and negative sets, where the positive set contains examples of analyst identified P - or S -wave onsets and the negative set contains examples that do not contain P waves or S waves. Similarity is defined as the square of the Euclidean distance between vectors representing the scaled absolute values of the amplitudes of the observed signal and a given reference example in time windows of the same length. For both P waves and S waves, a single pass is done through the bandpassed data, producing a score function defined as the ratio of the sum of similarity to positive examples over the sum of similarity to negative examples for each window. A phase arrival is chosen as the centre position of the window that maximizes the score function. The method is tested on two local earthquake data sets, consisting of 98 known events from the Parkfield region in central California and 32 known events from the Alpine Fault region on the South Island of New Zealand. For P -wave picks, using a reference set containing two picks from the Parkfield data set, 98 per cent of Parkfield and 94 per cent of Alpine Fault picks are determined within 0.1 s of the analyst pick. For S -wave picks, 94 per cent and 91 per cent of picks are determined within 0.2 s of the analyst picks for the Parkfield and Alpine Fault data set, respectively. For the Parkfield data set, our method picks 3520 P -wave picks and 3577 S -wave picks out of 4232 station–event pairs. For the Alpine Fault data set, the method picks 282 P -wave picks and 311 S -wave picks out of a total of 344 station–event pairs. For our testing, we note that the vast majority of station–event pairs have analyst picks, although some analyst picks are excluded based on an accuracy assessment. Finally, our tests suggest that the method is portable, allowing the use of a reference set from one region on data from a different region using relatively few reference picks.

Description: We examine shear velocity anisotropy in the Yuha Desert, California using aftershocks of the 2010 M7.2 El Mayor-Cucapah earthquake. The Yuha Desert is underlain by a complex network of right- and left-lateral conjugate faults, some of which experienced triggered slip during the El Mayor-Cucapah earthquake. An automated method that implements multiple measurement windows and a range of bandpass filters is used to estimate the fast direction ( ) and delay time ( t ) of the split shear waves. We find an average oriented approximately north–south suggesting it is primarily controlled by the regional maximum compressive stress direction. However, the spatial variability in reveals that the fault structures that underlie the Yuha Desert also influence the measured splitting parameters. We infer that the northeast- and northwest-oriented reflect shear fabric subparallel to the conjugate fault structures. We do not observe a simple correlation between t and hypocentral distance. Instead, the observed spatial variation in t suggests that near-source variation in anisotropic strength may be equal to or more important than effects local to the station. No temporal variation in splitting parameters is observed during the 70-day period following the main shock. In this region of complex faulting, we observe a spatially variable pattern of anisotropy that is both stress- and structure-controlled. This study suggests that shear fabric can form even along short, discontinuous fault strands with minimal offset.

Description: The properties of the overburden transmission response are of particular interest for the analysis of reflectivity illumination or blurring in seismic depth imaging. The first step to showing a transmission-operator reciprocity property is to identify the symmetry of the so-called displacement-to-traction operators. The latter are analogous to Dirichlet-to-Neumann operators and they may also be called impedance operators. Their symmetry is deduced here after development of a formal spectral or modal theory of lateral wavefunctions in a laterally heterogeneous generally anisotropic elastic medium. The elastic lateral modes are displacement-traction 6-vectors and they are built from two auxiliary 3-vector lateral-mode bases. These auxiliary modes arise from Hermitian and anti-Hermitian operators, so they have familiar properties such as orthogonality. There is no assumption of down/up symmetry of the elasticity tensor, but basic assumptions are made about the existence and completeness of the elastic modes. A point-symmetry property appears and plays a central role. The 6-vector elastic modes have a symplectic orthogonality property, which facilitates the development of modal expansions for 6-vector functions of the lateral coordinates when completeness is assumed. While the elastic modal theory is consistent with the laterally homogeneous case, numerical work would provide confidence that it is correct in general. An appendix contains an introductory overview of acoustic lateral modes that were studied by other authors, given from the perspective of this new work. A distinction is drawn between unit normalization of scalar auxiliary modes and a separate energy-flux normalization of 2-vector acoustic modes. Neither is crucial to the form of acoustic pressure-to-velocity or impedance operators. This statement carries over to the elastic case for the 3-vector auxiliary- and 6-vector elastic-mode normalizations. The modal theory is used to construct the kernel of the elastic displacement-to-traction or impedance operator. Symmetry properties of this operator are then deduced, which is the main goal of this paper. The implications of elastic impedance-operator symmetry and the symplectic property for the transmission and reflection responses of finite regions are described in a companion paper.

Description: On 2014 April 1, a magnitude M w  8.1 interplate thrust earthquake ruptured a densely instrumented region of Iquique seismic gap in northern Chile. The abundant data sets near and around the rupture zone provide a unique opportunity to study the detailed source process of this megathrust earthquake. We retrieved the spatial and temporal distributions of slip during the main shock and one strong aftershock through a joint inversion of teleseismic records, GPS offsets and strong motion data. The main shock rupture initiated at a focal depth of about 25 km and propagated around the hypocentre. The peak slip amplitude in the model is ~6.5 m, located in the southeast of the hypocentre. The major slip patch is located around the hypocentre, spanning ~150 km along dip and ~160 km along strike. The associated static stress drop is ~3 MPa. Most of the seismic moment was released within 150 s. The total seismic moment of our preferred model is 1.72 x 10 21 N m, equivalent to M w  8.1. For the strong aftershock on 2014 April 3, the slip mainly occurred in a relatively compact area, and the major slip area surrounded the hypocentre with the peak amplitude of ~2.5 m. There is a secondary slip patch located downdip from the hypocentre with the peak slip of ~2.1 m. The total seismic moment is about 3.9 x 10 20 N m, equivalent to M w  7.7. Between the rupture areas of the main shock and the 2007 November 14 M w  7.7 Antofagasta, Chile earthquake, there is an earthquake vacant zone with a total length of about 150 km. Historically, if there is no big earthquake or obvious aseismic creep occurring in this area, it has a great potential of generating strong earthquakes with magnitude larger than M w 7.0 in the future.

Description: Tectonic tremor (TT) and low-frequency earthquakes (LFEs) have been found in the deeper crust of various tectonic environments globally in the last decade. The spatial-temporal behaviour of LFEs provides insight into deep fault zone processes. In this study, we examine recurrence times from a 12-yr catalogue of 88 LFE families with ~730 000 LFEs in the vicinity of the Parkfield section of the San Andreas Fault (SAF) in central California. We apply an automatic burst detection algorithm to the LFE recurrence times to identify the clustering behaviour of LFEs (LFE bursts) in each family. We find that the burst behaviours in the northern and southern LFE groups differ. Generally, the northern group has longer burst duration but fewer LFEs per burst, while the southern group has shorter burst duration but more LFEs per burst. The southern group LFE bursts are generally more correlated than the northern group, suggesting more coherent deep fault slip and relatively simpler deep fault structure beneath the locked section of SAF. We also found that the 2004 Parkfield earthquake clearly increased the number of LFEs per burst and average burst duration for both the northern and the southern groups, with a relatively larger effect on the northern group. This could be due to the weakness of northern part of the fault, or the northwesterly rupture direction of the Parkfield earthquake.

Description: We observe seasonal seismic wave speed changes (d v / v ) in the San Jacinto fault area and investigate several likely source mechanisms. Velocity variations are obtained from analysis of 6 yr data of vertical component seismic noise recorded by 10 surface and six borehole stations. We study the interrelation between d v / v records, frequency-dependent seismic noise properties, and nearby environmental data of wind speed, rain, ground water level, barometric pressure and atmospheric temperature. The results indicate peak-to-peak seasonal velocity variations of ~0.2 per cent in the 0.5–2 Hz frequency range, likely associated with genuine changes of rock properties rather than changes in the noise field. Phase measurements between d v / v and the various environmental data imply that the dominant source mechanism in the arid study area is thermoelastic strain induced by atmospheric temperature variations. The other considered environmental effects produce secondary variations that are superimposed on the thermal-based changes. More detailed work with longer data on the response of rocks to various known external loadings can help tracking the evolving stress and effective rheology at depth.

Description: We develop an algorithm for automatic identification of fault zone trapped waves in data recorded by seismic fault zone arrays. Automatic S picks are used to identify time windows in the seismograms for subsequent search for trapped waves. The algorithm calculates five features in each seismogram recorded by each station: predominant period, 1 s duration energy (representative of trapped waves), relative peak strength, arrival delay and 6 s duration energy (representative of the entire seismogram). These features are used collectively to identify stations in the array with seismograms that are statistical outliers. Applying the algorithm to large data sets allows for distinguishing genuine trapped waves from occasional localized site amplification in seismograms of other stations. The method is verified on a test data set recorded across the rupture zone of the 1992 Landers earthquake, for which trapped waves were previously identified manually, and is then applied to a larger data set with several thousand events recorded across the San Jacinto fault zone. The developed technique provides an important tool for systematic objective processing of large seismic waveform data sets recorded near fault zones.

Description: We studied the seismicity rate increase in the Kanto region around Tokyo following the 2011 Tohoku-Oki earthquake ( M w 9.0) to examine whether this increase was correlated with the static increases in the Coulomb failure function (CFF) of the Tohoku-Oki earthquake sequence. Because earthquakes in the Kanto region exhibit various focal mechanisms, the receiver faults for the CFF were assumed to be the focal mechanism solutions for nearly 19 000 earthquakes that previously occurred. Our results showed that the number of earthquakes for which the mechanism solutions had a positive CFF (~12 000) is much larger than those that had a negative CFF (~2000). Comparison of the CFF values for earthquakes before and after the Tohoku-Oki earthquake showed that the latter had more positive values; this supports the hypothesis that the coseismic stress change transferred from the Tohoku-Oki earthquake sequence is the major contributing factor to the increased seismicity rate in the Kanto region.

Description: We determined the first high-resolution P - and S -wave attenuation ( Qp and Qs ) tomography of the crust and upper mantle under the entire Nankai subduction zone from the Nankai Trough to the Japan Sea using a large number of high-quality t * data measured from P - and S -wave spectra of local earthquakes. The suboceanic earthquakes used in this study were relocated precisely using sP depth phases and ocean-bottom-seismometer data. The overall pattern of the obtained Q models is similar to that of velocity models of the study region. Our present results show that high- Q (i.e. weak attenuation) anomalies in the upper crust generally correspond to plutonic rocks widely exposed in the Nankai arc. Some of the low- Q (i.e. strong attenuation) anomalies in the upper crust along the Pacific coast are associated with the Cretaceous–Cenozoic accretionary wedge. Obvious low- Q anomalies exist in the crust under the active arc volcanoes. Most of the large inland crustal earthquakes are located in or around the low- Q zones in the crust. The subducting Philippine Sea slab is imaged clearly as a landward dipping high- Q zone. Prominent low- Q anomalies are revealed in the mantle wedge under the volcanic front and backarc area, which reflect the source zone of arc magmatism caused by slab dehydration and corner flow in the mantle wedge. Significant low- Q anomalies exist in the forearc mantle wedge, which reflects a highly hydrated and serpentinized forearc mantle wedge due to abundant fluids released from dehydration of the young and warm Philippine Sea slab.

Description: Differences between 3-D numerical predictions of earthquake ground motion in the Mygdonian basin near Thessaloniki, Greece, led us to define four canonical stringent models derived from the complex realistic 3-D model of the Mygdonian basin. Sediments atop an elastic bedrock are modelled in the 1D-sharp and 1D-smooth models using three homogeneous layers and smooth velocity distribution, respectively. The 2D-sharp and 2D-smooth models are extensions of the 1-D models to an asymmetric sedimentary valley. In all cases, 3-D wavefields include strongly dispersive surface waves in the sediments. We compared simulations by the Fourier pseudo-spectral method (FPSM), the Legendre spectral-element method (SEM) and two formulations of the finite-difference method (FDM-S and FDM-C) up to 4 Hz. The accuracy of individual solutions and level of agreement between solutions vary with type of seismic waves and depend on the smoothness of the velocity model. The level of accuracy is high for the body waves in all solutions. However, it strongly depends on the discrete representation of the material interfaces (at which material parameters change discontinuously) for the surface waves in the sharp models. An improper discrete representation of the interfaces can cause inaccurate numerical modelling of surface waves. For all the numerical methods considered, except SEM with mesh of elements following the interfaces, a proper implementation of interfaces requires definition of an effective medium consistent with the interface boundary conditions. An orthorhombic effective medium is shown to significantly improve accuracy and preserve the computational efficiency of modelling. The conclusions drawn from the analysis of the results of the canonical cases greatly help to explain differences between numerical predictions of ground motion in realistic models of the Mygdonian basin. We recommend that any numerical method and code that is intended for numerical prediction of earthquake ground motion should be verified through stringent models that would make it possible to test the most important aspects of accuracy.

Description: The M  ~ 7 1915 Fucino (Central Italy) earthquake represents one of the most destructive seismic events ever occurred in the Italian Peninsula. Several seismogenic faults have been proposed in the past decades as the source of the earthquake by means of different approaches and techniques that lead to a variety of speculations about the source mechanism and the fault location, often contrasting with one another. The 1915 earthquake produced a remarkable data set of 73 coseismic hydrological changes in the near and intermediate field that consist in variation of the flow of streams and springs, liquefaction, rise of water temperature and turbidity. In this paper, we study the coseismic water level changes induced by the 1915 earthquake in the near field to provide convincing clues on the geometry of the earthquake causative fault. We model the coseismic strain field induced by seventeen individual faults proposed through different approaches, and compare its pattern with the distribution of streamflow changes. We find: (i) clues on the most probable geometry of the earthquake causative fault. Best fits between modelled deformation and observed data are displayed by sources (derived by geological or seismological data) that share several distinctive features, as they are ~135°-striking, SW-dipping, 25–30-km-long normal faults located along the eastern side of the Fucino basin. These data point to the Serrone Fault and the Parasano Fault as the most likely causative structures and support the hypothesis that the coseismic ruptures observed in the field represented primary surface faulting. On the contrary, our calculations show that the Pescina Fault and the Ventrino Fault are secondary faults from the perspective of the hydrological response. Finally, one of the best scoring potential sources (from geological data) is a multifaulting system that considers the presence, in the central-western part of the basin, of fault splays synthetic and antithetic to the main seismogenic structures; therefore, we infer for these splays a possible active involvement in a 1915-like seismogenic process; (ii) evidence against a number of seismogenic structures that were previously associated with the earthquake. In particular, the plots of coseismic strain induced by sources uniquely derived by macroseismic or geodetic data prove to be inconsistent with the polarities of the hydrological signatures. Also, sources mainly characterized by reverse faulting and/or by right-lateral strike-slip component are discarded and (iii) as a final remark, we maintain that the study of the hydrological signatures of earthquake strain can offer an alternative tool in the investigation of the historical seismicity, to estimate the focal mechanism of major earthquakes capable of giving rise to a consistent data set of hydrological data.

Description: In this paper, we propose a weighted Runge–Kutta (RK) discontinuous Galerkin (WRKDG) method for wavefield modelling. For this method, we first transform the seismic wave equations in 2-D heterogeneous anisotropic media into a first-order hyperbolic system, and then combine the discontinuous Galerkin method (DGM) with a weighted RK time discretization. The time discretization is based on an implicit diagonal RK method and an explicit technique, which changes the implicit RK method into an explicit one. In addition, we introduce a weighting factor in the process. Linear and quadratic polynomials for spatial basis functions are typically employed. We investigate the properties of the method in great detail, including the stability criteria and numerical dispersion relations for solving the 2-D acoustic equations. Our analysis indicates that the stability condition for the WRKDG method is more relaxed compared with the classic total variation diminishing (TVD) RK discontinuous Galerkin (RKDG) method, resulting in a 1.7 times superiority for P 1 element and is about as efficient as TVD RKDG method for P 2 element in computational efficiency. We also demonstrate that the WRKDG method can suppress numerical dispersion more efficiently than the staggered-grid (SG) method on the same grid. The WRKDG method is applied to simulate the wavefields in a large velocity contrast model, a 2-D homogeneous transversely isotropic (TI) model, a fluid-filled fracture model, and a 2-D SEG/EAGE salt dome model. Regular rectangular and irregular triangular elements are used. The numerical results show that the WRKDG method can effectively suppress numerical dispersion and provide accurate information on the wavefield on a coarse mesh. Therefore, the method evidently reduces the scale of the problem and increases computational efficiency. In addition, promising numerical tests show that the WRKDG method combines well with split perfectly matched layer boundary conditions.

Description: This study considers the propagation of elastic waves in a fluid-saturated double-porosity solid. Constitutive relations are derived for the double porosity medium having the transverse isotropy with vertical axis of symmetry. The equations of motion are solved for the propagation of harmonic waves. Then, the amount of wave-induced fluid flow at any point is expressed in terms of normal strain components in the composite medium. For propagation in a plane, the equations of motion decouple into two independent systems. One of these represents the particle motion normal to the plane and identifies the SH-type motion for solid particles. The other system represents the in-plane motion of particles and governs the propagation of four coupled waves. Among these four waves, three are quasi-longitudinal ( qP 1 , qP 2 , qP 3 ) waves and one is a quasi-transverse ( qSV ) wave. A numerical example is solved to calculate the velocity and attenuation anisotropy of the five waves. Effect of local fluid flow is observed on the wave velocities as well as the polarization of the constituent particles in double porosity medium. Effects of wave-frequency, propagation direction, pore-fluid viscosity, permeability and radius of spherical inclusions are analysed on the velocities and attenuation.

Description: On 2010 January 18 and 22, two earthquakes of M W 5.3 and 5.2, respectively, occurred near the town of Efpalio on the western Gulf of Corinth. We performed a shear wave splitting analysis using the cross-correlation method and calculated V P /V S ratios for events that occurred in the epicentral area of the Efpalio earthquakes, between 2009 January and 2010 December. The data analysis revealed the presence of shear wave splitting in the study area, as well as variations of the splitting parameters and V P /V S ratios. The average values of time-delay, fast polarization direction and V P /V S ratio for the time period before the Efpalio earthquakes, were calculated at 2.9 ± 0.4 ms km –1 , 92° ± 10° and 1.76 ± 0.04, respectively, while after the occurrence of the earthquakes, including the aftershock sequence, they were calculated at 5.5 ± 0.5 ms km –1 , 82° ± 9° and 1.88 ± 0.04. A few months after the occurrence of the Efpalio earthquakes, the mentioned splitting parameters were calculated at 3.6 ± 0.4 ms km –1 and 83° ± 9°. V P /V S ratio exhibited a mean value of 1.87 ± 0.04. The mean fast polarization directions were in general consistent with the regional stress field, almost perpendicular to the direction of the extension of the Gulf of Corinth. The observed increase in the time-delays and V P /V S ratios after the Efpalio earthquakes indicates changes in the crustal properties, which possibly resulted from variations in the pre-existing microcrack system characteristics. We suggest that a migration of fluids in the form of overpressured liquids, which are likely originated from dehydration reactions within the crust, was triggered by the Efpalio earthquakes and caused the observed variations. The findings of this work are consistent with those of previous studies that have indicated the presence of fluids of crustal origin in the study area.

Description: Global parametric sensitivity analysis of a complete seismic risk scenario is performed by means of the response surface method based on sparse Chebyshev polynomial expansion. Hazard model, described via a mixed technique combining Green's function and stochastic methods, is used to simulate free-field seismic records. Vulnerability model, based on the period-dependent damage spectra, is constructed using principles of normalized hysteretic energy and displacement ductility. Approximation of the risk model by means of orthogonal polynomials allows to compute sensitivity indices for model parameters directly from the interpolation matrix. The Monte Carlo method runs on top of the meta-model and allows to build probability density functions of the damage index at target period. The proposed method for sensitivity analysis is useful for understanding the impact of aleatory uncertainties on the risk model outputs.

Description: In order to study fine scale structure of the Earth's deep interior, it is necessary to extract generally weak body wave phases from seismograms that interact with various discontinuities and heterogeneities. The recent deployment of large-scale dense arrays providing high-quality data, in combination with efficient seismic data processing techniques, may provide important and accurate observations over large portions of the globe poorly sampled until now. Major challenges are low signal-to-noise ratios (SNR) and interference with unwanted neighbouring phases. We address these problems by introducing scale-dependent slowness filters that preserve time-space resolution. We combine complex wavelet and slant-stack transforms to obtain the slant-stacklet transform. This is a redundant high-resolution directional wavelet transform with a direction (here slowness) resolution that can be adapted to the signal requirements. To illustrate this approach, we use this expansion to design coherence-driven filters that allow us to obtain clean PcP observations (a weak phase often hidden in the coda of the P wave), for events with magnitude M w  〉 5.4 and distances up to 80°. In this context, we then minimize a linear misfit between P and PcP waveforms to improve the quality of PcP–P traveltime measurements as compared to a standard cross-correlation method. This significantly increases both the quantity and the quality of PcP–P differential traveltime measurements available for the modelling of structure near the core–mantle boundary. The accuracy of our measurements is limited mainly by the highest frequencies of the signals used and the level of noise. We apply this methodology to two examples of high-quality data from dense arrays located in north America. While focusing here on body-wave separation, the tools we propose are more general and may contribute to enhancing seismic signal observations in global seismology in situations of low SNR and high signal interference.

Description: Azimuthal anisotropy of attenuation is a physical phenomenon related to the directional change of attenuation. This study examines the frequency properties and directional attenuation of SKS waves. The directional frequency-dependent characteristics of SKS waves are investigated in the frequency band of 0.02–0.5 Hz using data from 53 permanent seismic stations located throughout the northern Yangtze Craton, the southern North China Craton and adjacent areas. In addition to normal splitting behavior, the analysis reveals that many SKS splitting measurements exhibit a lemniscate shape, reflecting frequency differences along fast and slow polarization directions. Frequency analysis shows that spectral ratios between fast/slow components of the lemniscate-type splitting results fluctuate strongly in a higher frequency band of 0.2–0.5 Hz, and fluctuate less within the main frequency band of 0.02–0.2 Hz. For each station, the ratio of the peak amplitude of the fast/slow components can be represented as a cotangential function of event backazimuth multiplying with a constant = 0.42 ± 0.10. This transformation shows that the regional average angles consistently fall within the relatively narrow range of –46.5 ± 3° with respect to the north, suggesting that a regional tectonic controlling factor dictates the relatively uniform directional attenuation of SKS waves within the frequency band of 0.02–0.2 Hz. Further analysis is performed by projecting the SKS waves onto the components along and perpendicular to the regional average angles. The calculation also shows that, in the 0.02–0.2 Hz band, the relationship between amplitude ratio and event backazimuth matches a cotangential functions with the same best matching angles and constant a 〈 1. Synthetic calculations demonstrate that although different filters influence the splitting parameters, attenuation anisotropy cannot be explained by elastic anisotropic media, including multilayer anisotropy and anisotropy with a tilting symmetrical axis. This observed behavior of the SKS wave may arise from the combined effects of frequency-dependent attenuation anisotropy and small-scale heterogeneities in the crust and the upper mantle.

Description: Detection of low magnitude event is critical and challenging in seismology. We develop a new method, named the match and locate (M&L) method, for small event detection. The M&L method employs some template events and detects small events through stacking cross-correlograms between waveforms of the template events and potential small event signals in the continuous waveforms over multiple stations and components, but the stacking is performed after making relative traveltime corrections based on the relative locations of the template event and potential small event scanning through a 3-D region around the template. Compared to the current methods of small event detection, the M&L method places event detection to a lower magnitude level and extends the capability of detecting small events that have large distance separations from the template. The method has little dependence on the accuracy of the velocity models used, and, at the same time, provides high-precision location information of the detected small-magnitude events. We demonstrate the effectiveness of the M&L method and its advantage over the matched filter method using examples of scaled-down earthquakes occurring in the Japan Island and foreshock detection before the 2011 M w 9.0 Tohoku earthquake. In the foreshock detection, the M&L method detects four times more events (1427) than the templates and 9 per cent (134) more than the matched filter under the same detection threshold. Up to 41 per cent (580) of the detected events are not located at the template locations with the largest separation of 9.4 km. Based on the identified foreshocks, we observe five sequences of foreshock migration along the trench-parallel direction toward the epicentre of the M w 9.0 main shock.

Description: Typical infrasound propagation paths extend into the middle and upper atmosphere before turning back to Earth near the stratopause or lower thermosphere. The modelling of infrasound propagation is complicated by several factors. Propagation of sound in the atmosphere becomes increasingly non-linear as the mean density of the atmosphere decreases. As a consequence, infrasound propagation, which can follow paths high into the atmosphere before turning back to earth, is intrinsically non-linear. Further, attenuation of sound also increases as the mean density of the atmosphere decreases. Linear propagation modelling predicts severe attenuation along thermospheric paths, in contradiction with observation. Finally, the currently available atmospheric specifications necessarily involve temporal and spatial interpolations which can lead to the omission of atmospheric fluctuations to which infrasound propagation is sensitive. In this paper, existing methods for modelling weakly non-linear propagation are extended to account for the moving inhomogeneous medium and are then used to study the modelling of waveforms produced by impulsive events. It is found that there is a substantial interplay between signal attenuation and non-linear distortion. Waveform steepening and shocking in the middle and upper atmosphere associated with higher harmonic generation is moderated by attenuation while attenuation in the upper atmosphere is moderated by period lengthening associated with low frequency generation. The modelling of observed data is then considered. Using traveltime data, corrections to the atmospheric specifications are obtained. The predicted waveform evolution, modelled using the developed non-linear propagation model and the corrected atmospheric specifications, matches the observations well.

Description: Scotland is a relatively aseismic region for the use of local earthquake tomography, but 40 yr of earthquakes recorded by a good and growing network make it possible. A careful selection is made from the earthquakes located by the British Geological Survey (BGS) over the last four decades to provide a data set maximising arrival time accuracy and ray path coverage of Scotland. A large number of 1-D velocity models with different layer geometries are considered and differentiated by employing quarry blasts as ground-truth events. Then, SIMULPS14 is used to produce a robust 3-D tomographic P -wave velocity model for Scotland. In areas of high resolution the model shows good agreement with previously published interpretations of seismic refraction and reflection experiments. However, the model shows relatively little lateral variation in seismic velocity except at shallow depths, where sedimentary basins such as the Midland Valley are apparent. At greater depths, higher velocities in the northwest parts of the model suggest that the thickness of crust increases towards the south and east. This observation is also in agreement with previous studies. Quarry blasts used as ground truth events and relocated with the preferred 3-D model are shown to be markedly more accurate than when located with the existing BGS 1-D velocity model.

Description: We present a new approach for understanding the origin and nature of seismic anomalies in the continental crust of the Northern Middle East. We have created detailed attenuation ( Q Lg ) and velocity ( V Lg ) models for the Northern Middle East based on the analysis of waveforms of regional seismic phase Lg from 3171 regional earthquakes recorded at 578 stations in Turkish and Iranian Plateaus and surrounding regions. The attenuation and velocity models are assumed to serve as proxies for the bulk average crustal shear wave attenuation ( Q β ) and velocities ( V s ). 31 232 reliable Lg spectra were collected and used to measure the two-station method (TSM) and reverse two-station/event method (RTM) Lg Q at 1 Hz ( Q 0 ) and its frequency-dependence factor ( ). The Lg Q 0 and values are measured over the individual TSM and RTM paths and are then used to perform an LSQR tomographic inversion for lateral variations in Q 0 and . The Lg Q 0 and models both correlate well with the major tectonic boundaries in the region. The tomographic models as well as the individual TSM and RTM measurements show lower values of Lg Q 0 over the Turkish-Anatolian Plateau (〈150) than those observed over the Iranian Plateau (150–400). Furthermore, we obtained the Lg group velocity model by inverting the time of the first arrival of the Lg waveform on each seismogram. Our Q measurements are strongly correlated with the measurements of Lg group velocity ( V Lg ) suggesting that the source of many of the low Q and velocity anomalies is likely the same. The regions where we see negative correlations are likely a result of Sn to Lg converted energy. Our results also have implications for the far field ground motions, suggesting that large earthquakes in eastern Iran could have a significant far field ground motions due to relatively low crustal attenuation within the Iranian plateau.

Description: Several studies have focused on the role of damage zone (DZ) on the hydromechanical behaviour of faults by assuming a fractured DZ (i.e. low stiffness/high permeability). Yet, this vision may not be valid in all geological settings, in particular, in high-porosity reservoirs as targeted by several underground exploitations. We investigate the impact of a high-stiff/low-permeable DZ on the shear reactivation of a blind, undetectable normal fault (1 km long, ≤10 m offset), with a 0.5 m thick low-porosity/permeability fault core during fluid injection into a high-porosity reservoir. The spatial distribution of effective properties (elastic moduli, Biot's coefficients and permeability) of DZ including deformation bands (DB; elliptic inclusions) and intact rock were derived using upscaling analytical expressions. The influence of DZ on the hydromechanical behaviour of the fault zone was numerically explored using 2-D plane-strain finite-element simulations within the framework of fully saturated isothermal porous media by accounting for an orthotropic elastic rheology. The numerical results showed that the presence of DB plays a protective role by reducing the potential for shear reactivation inside the fault core. On the other hand, they favour shear failure in the vicinity of the fault core (off-fault damage) by accelerating the decrease of the minimum principal effective stress while limiting the decrease of the maximum one. This behaviour is strongly enhanced by the fault-parallel DZ effective stiffness, but limited by the combined effect of fault-normal DZ effective permeability and of the Biot's coefficients. This can have implications for the location and size of aftershocks during fault reactivation.

Description: Full waveform inversion (FWI) aims to reconstruct high-resolution subsurface models from the full wavefield, which includes diving waves, post-critical reflections and short-spread reflections. Most successful applications of FWI are driven by the information carried by diving waves and post-critical reflections to build the long-to-intermediate wavelengths of the velocity structure. Alternative approaches, referred to as reflection waveform inversion (RWI), have been recently revisited to retrieve these long-to-intermediate wavelengths from short-spread reflections by using some prior knowledge of the reflectivity and a scale separation between the velocity macromodel and the reflectivity. This study presents a unified formalism of FWI, named as Joint FWI, whose aim is to efficiently combine the diving and reflected waves for velocity model building. The two key ingredients of Joint FWI are, on the data side, the explicit separation between the short-spread reflections and the wide-angle arrivals and, on the model side, the scale separation between the velocity macromodel and the short-scale impedance model. The velocity model and the impedance model are updated in an alternate way by Joint FWI and waveform inversion of the reflection data (least-squares migration), respectively. Starting from a crude velocity model, Joint FWI is applied to the streamer seismic data computed in the synthetic Valhall model. While the conventional FWI is stuck into a local minimum due to cycle skipping, Joint FWI succeeds in building a reliable velocity macromodel. Compared with RWI, the use of diving waves in Joint FWI improves the reconstruction of shallow velocities, which translates into an improved imaging at deeper depths. The smooth velocity model built by Joint FWI can be subsequently used as a reliable initial model for conventional FWI to increase the high-wavenumber content of the velocity model.

Description: The determination of near-surface attenuation for hard rock sites is an important issue in a wide range of seismological applications, particularly seismic hazard analysis. In this article we choose six hard to very-hard rock sites ( Vs 30 1030–3000 m s –1 ) and apply a range of analysis methods to measure the observed attenuation at distance based on a simple exponential decay model with whole-path attenuation operator r . The r values are subsequently decoupled from path attenuation ( Q ) so as to obtain estimates of near-surface attenuation ( 0 ). Five methods are employed to measure r which can be split into two groups: broad-band methods and high-frequency methods. Each of the applied methods has advantages and disadvantages, which are explored and discussed through the comparison of results from common data sets. In our first step we examine the variability of the individual measured r values. Some variation between methods is expected due to simplifications of source, path, and site effects. However, we find that significant differences arise between attenuation measured on individual recordings, depending on the method employed or the modelling decisions made during a particular approach. Some of the differences can be explained through site amplification effects: although usually weak at rock sites, amplification may still lead to bias of the measured r due to the chosen fitting frequency bandwidth, which often varies between methods. At some sites the observed high-frequency spectral shape was clearly different to the typical r attenuation model, with curved or bi-linear rather than linear decay at high frequencies. In addition to amplification effects this could be related to frequency-dependent attenuation effects [e.g. Q ( f )]: since the r model is implicitly frequency independent, r will in this case be dependent on the selected analysis bandwidth. In our second step, using the whole-path r data sets from the five approaches, we investigate the robustness of the near-surface attenuation parameter 0 and the influence of constraints, such as assuming a value for the regional crustal attenuation ( Q ). We do this by using a variety of fitting methods: least squares, absolute amplitude and regressions with and without fixing Q to an a priori value. We find that the value to which we fix Q strongly influences the near-surface attenuation term 0 . Differences in Q derived from the data at the six sites under investigation could not be reconciled with the average values found previously over the wider Swiss region. This led to starkly different 0 values, depending on whether we allowed for a data-driven Q , or whether we forced Q to be consistent with existing simulation models or ground motion prediction equations valid for the wider region. Considering all the possible approaches we found that the contribution to epistemic uncertainty for 0 determination at the six hard-rock sites in Switzerland could be represented by a normal distribution with standard deviation 0  = 0.0083 ± 0.0014 s.

Description: It was established over a decade ago that the normal modes of the Earth are continuously excited at times without large earthquakes, but the sources of the ‘seismic hum’ have remained unresolved. In addition to the normal modes of the Earth, we show spectral lines in seismic data with frequencies which correspond closely to normal modes of the Sun. Moreover, the widths of the low-frequency lines in the seismic spectra are similar to those of solar modes and much narrower than those of the Earth's normal mode peaks. These seismic lines are highly coherent with magnetic fields measured on both the Geostationary Operations Environmental Satellite (GOES)–10 satellite and the Advanced Composition Explorer (ACE) spacecraft located at L1, 1.5 million km sunward of Earth suggesting that the solar modes are transmitted to the Earth by the interplanetary magnetic field and solar wind. The solar modes are split by multiples of a cycle/day and, surprisingly, by the ‘quasi two-day’ mode and other frequencies. Both the phase of the coherences and slight frequency offsets between seismic and geomagnetic data at observatories exclude the possibility that these effects are simply spurious responses of the seismometers to the geomagnetic field. We emphasize data from low-noise seismic observatories: Black Forest (BFO), Piñon Flat (PFO), Eskdalemuir (ESK) and Obninsk (OBN). Horizontal components of seismic velocity show higher coherences with the external (ACE) magnetic field than do the vertical components. This effect appears to be larger near the seismic torsional, or T -mode, frequencies.

Description: Dynamic tilts (rotational motion around horizontal axes) change the projection of local gravity onto the horizontal components of seismometers. This causes sensitivity of these components to tilt, especially at low frequencies. We analyse the consequences of this effect onto moment tensor inversion for very long period (vlp) events in the near field of active volcanoes on the basis of synthetic examples using the station distribution of a real deployed seismic network and the topography of Mt. Merapi volcano (Java, Indonesia). The examples show that for periods in the vlp range of 10–30 s tilt can have a strong effect on the moment tensor inversion, although its effect on the horizontal seismograms is significant only for few stations. We show that tilts can be accurately computed using the spectral element method and include them in the Green's functions. The (simulated) tilts might be largely influenced by strain–tilt coupling (stc). However, due to the frequency dependence of the tilt contribution to the horizontal seismograms, only the largest tilt signals affect the source inversion in the vlp frequency range. As these are less sensitive to stc than the weaker signals, the effect of stc can likely be neglected in this application. In the converse argument, this is not necessarily true for longer periods, where the horizontal seismograms are dominated by the tilt signal and rotational sensors would be necessary to account for it. As these are not yet commercially available, this study underlines the necessity for the development of such instruments.

Description: Computationally efficient 3-D frequency-domain full waveform inversion (FWI) is applied to ocean-bottom cable data from the Valhall oil field in the visco-acoustic vertical transverse isotropic (VTI) approximation. Frequency-domain seismic modelling is performed with a parallel sparse direct solver on a limited number of computer nodes. A multiscale imaging is performed by successive inversions of single frequencies in the 3.5–10 Hz frequency band. The vertical wave speed is updated during FWI while density, quality factor Q P and anisotropic Thomsen's parameters and are kept fixed to their initial values. The final FWI model shows the resolution improvement that was achieved compared to the initial model that was built by reflection traveltime tomography. This FWI model shows a glacial channel system at 175 m depth, the footprint of drifting icebergs on the palaeo-seafloor at 500 m depth, a detailed view of a gas cloud at 1 km depth and the base cretaceous reflector at 3.5 km depth. The relevance of the FWI model is assessed by frequency-domain and time-domain seismic modelling and source wavelet estimation. The agreement between the modelled and recorded data in the frequency domain is excellent up to 10 Hz although amplitudes of modelled wavefields propagating across the gas cloud are overestimated. This might highlight the footprint of attenuation, whose absorption effects are underestimated by the homogeneous background Q P model ( Q P = 200). The match between recorded and modelled time-domain seismograms suggests that the inversion was not significantly hampered by cycle skipping. However, late arrivals in the synthetic seismograms, computed without attenuation and with a source wavelet estimated from short-offset early arrivals, arrive 40 ms earlier than the recorded seismograms. This might result from dispersion effects related to attenuation. The repeatability of the source wavelets inferred from data that are weighted by a linear gain with offset is dramatically improved when they are estimated in the FWI model rather than in the smooth initial model. The two source wavelets, estimated in the FWI model from data with and without offset gain, show a 40 ms time-shift, which is consistent with the previous analysis of the time-domain seismograms. The computational efficiency of our frequency-domain approach is assessed against a recent time-domain FWI case study performed in a similar geological environment. This analysis highlights the efficiency of the frequency-domain approach to process a large number of sources and receivers with limited computational resources, thanks to the efficiency of the substitution step performed by the direct solver. This efficiency can be further improved by using a block-low rank version of the multifrontal solver and by exploiting the sparsity of the source vectors during the substitution step. Future work will aim to update attenuation and density at the same time of the vertical wave speed.

Description: Inversions for the full slip distribution of earthquakes provide detailed models of earthquake sources, but stability and non-uniqueness of the inversions is a major concern. The problem is underdetermined in any realistic setting, and significantly different slip distributions may translate to fairly similar seismograms. In such circumstances, inverting for a single best model may become overly dependent on the details of the procedure. Instead, we propose to perform extended fault inversion trough falsification. We generate a representative set of heterogeneous slipmaps, compute their forward predictions, and falsify inappropriate trial models that do not reproduce the data within a reasonable level of mismodelling. The remainder of surviving trial models forms our set of coequal solutions. The solution set may contain only members with similar slip distributions, or else uncover some fundamental ambiguity such as, for example, different patterns of main slip patches. For a feasibility study, we use teleseismic body wave recordings from the 2012 September 5 Nicoya, Costa Rica earthquake, although the inversion strategy can be applied to any type of seismic, geodetic or tsunami data for which we can handle the forward problem. We generate 10 000 pseudo-random, heterogeneous slip distributions assuming a von Karman autocorrelation function, keeping the rake angle, rupture velocity and slip velocity function fixed. The slip distribution of the 2012 Nicoya earthquake turns out to be relatively well constrained from 50 teleseismic waveforms. Two hundred fifty-two slip models with normalized L1-fit within 5 per cent from the global minimum from our solution set. They consistently show a single dominant slip patch around the hypocentre. Uncertainties are related to the details of the slip maximum, including the amount of peak slip (2–3.5 m), as well as the characteristics of peripheral slip below 1 m. Synthetic tests suggest that slip patterns such as Nicoya may be a fortunate case, while it may be more difficult to unambiguously reconstruct more distributed slip from teleseismic data.

Description: Regional body-wave tomography, also called ACH tomography, is the inversion of relative traveltime residuals of teleseismic body waves measured at regional networks. We analyse the characteristics of the finite-frequency Fréchet kernels for P and S waves for this kind of tomography. Using a simplified geometry enables us to use the complete Green's function in the expression of the Fréchet kernels and analyse elements, which are usually neglected, like the importance of the near-field terms and the P -wave traveltime sensitivity to shear wave velocity variations. By comparing the kernels of the relative residuals and absolute ones, we show that relative residuals have a reduced sensitivity to heterogeneities of large dimensions, and that this reduction is a generalization of the fact that the average model is not recovered in ACH tomography. This sensitivity reduction affects equally short- and long-period residuals. We show in addition the presence of a sensitivity reduction at large depth for the long-period waves. Kernels and reflectivity impulse responses of the crust are used to analyse if crustal corrections should be made frequency-dependent in finite-frequency regional tomography. We find that in most cases the frequency dependence due to reverberations is substantial, and that in many realistic network configurations ray theory is unlikely to be well appropriate to compute crustal corrections for the long-period waves. We also find that the lateral dimensions of the crust affecting the traveltimes is frequency dependent and reaches, at long periods, 50 km for sedimentary basins and 100 km for Moho depth.

Description: A new method of fully nonlinear waveform fitting to measure interstation phase speeds and amplitude ratios is developed and applied to USArray. The Neighbourhood Algorithm is used as a global optimizer, which efficiently searches for model parameters that fit two observed waveforms on a common great-circle path by modulating the phase and amplitude terms of the fundamental-mode surface waves. We introduce the reliability parameter that represents how well the waveforms at two stations can be fitted in a time–frequency domain, which is used as a data selection criterion. The method is applied to observed waveforms of USArray for seismic events in the period from 2007 to 2010 with moment magnitude greater than 6.0. We collect a large number of phase speed data (about 75 000 for Rayleigh and 20 000 for Love) and amplitude ratio data (about 15 000 for Rayleigh waves) in a period range from 30 to 130 s. The majority of the interstation distances of measured dispersion data is less than 1000 km, which is much shorter than the typical average path-length of the conventional single-station measurements for source-receiver pairs. The phase speed models for Rayleigh and Love waves show good correlations on large scales with the recent tomographic maps derived from different approaches for phase speed mapping; for example, significant slow anomalies in volcanic regions in the western Unites States and fast anomalies in the cratonic region. Local-scale phase speed anomalies corresponding to the major tectonic features in the western United States, such as Snake River Plains, Basin and Range, Colorado Plateau and Rio Grande Rift have also been identified clearly in the phase speed models. The short-path information derived from our interstation measurements helps to increase the achievable horizontal resolution. We have also performed joint inversions for phase speed maps using the measured phase and amplitude ratio data of vertical component Rayleigh waves. These maps exhibit better recovery of phase speed perturbations, particularly where the strong lateral velocity gradient exists in which the effects of elastic focussing can be significant; that is, the Yellowstone hotspot, Snake River Plains, and Rio Grande Rift. The enhanced resolution of the phase speed models derived from the interstation phase and amplitude measurements will be of use for the better seismological constraint on the lithospheric structure, in combination with dense broad-band seismic arrays.

Description: Spatially heterogeneous inversions for the slip in major earthquakes are typically only available for modern, instrumentally recorded events. Stress reconstructions on active faults, which are potentially vital elements of earthquake hazard assessment, are usually based on one or two instrumental inversions and some semi-quantitative information concerning historical seismic activity. Here we develop a new Bayesian Monte Carlo inversion method for sparse, large uncertainty data, such as is available for many historical earthquakes. We use it to reconstruct the slip in two great historical earthquakes on the Sunda megathrust, from 200 yr of paleogeodetic data recorded in the stratigraphy of coral microatolls on the forearc islands. The technique is based on the stochastic forward modelling of many slip distributions, which are constrained by the observed fractal scaling in the slip field. A set of static elastic Green's functions is used to estimate the vertical displacement at the coral locations resulting from each trial slip distribution. The posterior expected value of the slip and its variance are obtained from the average of the trial slips, weighted by their ability to reproduce the coral displacements. We successfully test the method on a synthetic slip distribution, before applying it to reconstruct the slip in a recent, instrumentally recorded event. We then invert for the great 1797 and 1833 megathrust events under the Mentawai Islands west of Sumatra. When compared with the slip distributions in recent earthquakes, our results unambiguously indicate that the sequence is not consistent with the classic characteristic earthquake model; earthquakes tend to incrementally tessellate the active fault plane rather than repeatedly breaking a segment of it. The results indicate that homogeneous loading of a fault with heterogeneous initial stress is enough to explain the observations on the Mentawai segment of the Sunda megathrust; no time dependence of material properties, for example, is necessary. Furthermore, we speculate that even small amounts of nonlinearity in the rupture process would ensure that any tessellating sequence will not be repeated, even over very long times. This method could be adapted to incorporate other historical information, such as records of tsunami run up and Mercalli Intensity estimates. Finally, we suggest that reconstruction of stress based on long-term slip histories may prove useful in identifying the possible locations of future earthquakes.

Description: Mechanisms leading to the initiation and early-stage development of continental rifts remain enigmatic, in spite of numerous studies. Among the various rifting models, which were developed mostly based on studies of mature rifts, far-field stresses originating from plate interactions (passive rifting) and nearby active mantle upwelling (active rifting) are commonly used to explain rift dynamics. Situated atop of the hypothesized African Superplume, the incipient Okavango Rift Zone (ORZ) of northern Botswana is ideal to investigate the role of mantle plumes in rift initiation and development, as well as the interaction between the upper and lower mantle. The ORZ developed within the Neoproterozoic Damara belt between the Congo Craton to the northwest and the Kalahari Craton to the southeast. Mantle structure and thermal status beneath the ORZ are poorly known, mostly due to a complete paucity of broad-band seismic stations in the area. As a component of an interdisciplinary project funded by the United States National Science Foundation, a broad-band seismic array was deployed over a 2-yr period between mid-2012 and mid-2014 along a profile 756 km in length. Using P -to- S receiver functions (RFs) recorded by the stations, the 410 and 660 km discontinuities bordering the mantle transition zone (MTZ) are imaged for the first time. When a standard Earth model is used for the stacking of RFs, the apparent depths of both discontinuities beneath the Kalahari Craton are about 15 km shallower than those beneath the Congo Craton. Using teleseismic P - and S -wave traveltime residuals obtained by this study and lithospheric thickness estimated by previous studies, we conclude that the apparent shallowing is the result of a 100–150 km difference in the thickness of the lithosphere between the two cratons. Relative to the adjacent tectonically stable areas, no significant anomalies in the depth of the MTZ discontinuities or in teleseismic P - and S -wave traveltime residuals are found beneath the ORZ. These observations imply an absence of significant thermal anomalies in the MTZ and in the upper mantle beneath the incipient rift, ruling out the role of mantle plumes in the initiation of the ORZ. We propose that the initiation and development of the ORZ were the consequences of relative movements between the South African block and the rest of the African plate along a zone of lithospheric weakness between the Congo and Kalahari cratons. An area of thinner-than-normal MTZ is found at the SW corner of the study area. This anomaly, if confirmed by future studies, could suggest significant transferring of heat from the lower to the upper mantle.

Description: The Borborema province in NE Brazil is characterized by seismic sequences with small earthquakes that can last 10 yr or more. The seismicity in this region is concentrated in three main seismic zones. In this work, we investigate the stress field in one of these zones, the Acaraú Seismic Zone, which is located in the NW part of the Borborema province. This seismic zone exhibits earthquake sequences that contain repeated earthquakes with similar waveforms and a shallow depth. Using a local network, we investigated a seismic sequence close to the town of Santana do Acaraú from December 2009 to December 2010, and we present detailed results (velocity model, hypocentres and focal mechanism) from this network. In addition, we inverted seven focal mechanisms, including six that were used in previous studies, and determined the directions of the three main axes of the regional stress field. Selecting a very precise set of 12 earthquakes, we found an active seismic zone with a depth between 3.5 and 4.8 km and with a horizontal dimension of approximately 2.5 km in the NW–SE direction (azimuth of 118°) and a strike-slip focal mechanism. The new seismic fault and some of the previous seismic faults determined in previous studies occur near the continental-scale Transbrasiliano lineament, but they exhibit no direct relationship with that ancient structure. The stress field is characterized by NW–SE trending compression and NE–SW trending extension. This result suggests that the rheological contrast between the continental–oceanic crusts created flexural stresses with maximum horizontal compression parallel to the continental margin. This stress pattern occurs along the Potiguar basin and continues west as far as the Amazon fan along the Equatorial margin of Brazil. This stress field and related seismicity may be a characteristic of this type of passive margin that is generated during the transform shearing between the South America and Africa plates and that exhibits an abrupt oceanic–continent transition, steep continental slopes and high bathymetric gradients.

Description: Inversions for the full slip distribution of earthquakes provide detailed models of earthquake sources, but stability and non-uniqueness of the inversions is a major concern. The problem is underdetermined in any realistic setting, and significantly different slip distributions may translate to fairly similar seismograms. In such circumstances, inverting for a single best model may become overly dependent on the details of the procedure. Instead, we propose to perform extended fault inversion trough falsification. We generate a representative set of heterogeneous slipmaps, compute their forward predictions, and falsify inappropriate trial models that do not reproduce the data within a reasonable level of mismodelling. The remainder of surviving trial models forms our set of coequal solutions. The solution set may contain only members with similar slip distributions, or else uncover some fundamental ambiguity such as, for example, different patterns of main slip patches. For a feasibility study, we use teleseismic body wave recordings from the 2012 September 5 Nicoya, Costa Rica earthquake, although the inversion strategy can be applied to any type of seismic, geodetic or tsunami data for which we can handle the forward problem. We generate 10 000 pseudo-random, heterogeneous slip distributions assuming a von Karman autocorrelation function, keeping the rake angle, rupture velocity and slip velocity function fixed. The slip distribution of the 2012 Nicoya earthquake turns out to be relatively well constrained from 50 teleseismic waveforms. Two hundred fifty-two slip models with normalized L1-fit within 5 per cent from the global minimum from our solution set. They consistently show a single dominant slip patch around the hypocentre. Uncertainties are related to the details of the slip maximum, including the amount of peak slip (2–3.5 m), as well as the characteristics of peripheral slip below 1 m. Synthetic tests suggest that slip patterns such as Nicoya may be a fortunate case, while it may be more difficult to unambiguously reconstruct more distributed slip from teleseismic data.

Description: We image the P - and S -wave structure of the upper mantle in southwestern Scandinavia using a wavelet-based, multiscale parametrization and finite-frequency theory to model wave propagation. Relative traveltime residuals of direct P and S waves are measured in a high- and low-frequency band and are corrected for crustal structure using a detailed model for the study area. A range of resolution tests are used to find optimal damping values not only for variations in V P and V S separately, but also for perturbations in their ratio V P / V S . The tests show that features down to a size of 100 (150) km can be well resolved in the P (S) tomography. To ease comparison with previous studies we also perform ray-theoretical multiscale tomographies, and to test the degree of vertical smearing we evaluate different parametrizations in the vertical direction (wavelet-based multiscale and convolutional quelling). Our finite-frequency, multiscale images of variations in V P and V S confirm the existence of low velocities below southern Norway and Denmark and high velocities beneath the shield proper in Sweden, as seen in previous studies, but add more details to this simplified picture. The low velocities below southern Norway and Denmark are confined to a channel-like structure at about 100–200 km depth, and the lateral transition from low to high velocities follows zones of Carboniferous-Permian extension and magmatism very closely. A deeper low-velocity anomaly below central southern Norway emerges from the channel at 150 km depth and extends to a depth of 350 km. In the Swedish area we infer high-velocity anomalies in V P and V S , and negative anomalies in V P / V S that indicate a strongly depleted mantle. We propose that the episodic erosion and convective removal of an originally thick mantle lithosphere below southern Norway to its current thickness of about 100 km could have been a trigger for episodic uplift in the Mesozoic and Cenozoic.

Description: We analyse continuous recordings from 23 broadband seismic stations near Alberta, the southwestern sector of the Western Canada Sedimentary Basin. Noise correlation tomo-graphy based on vertical-component seismograms reveals below-average shear velocities at shallow and middle crustal depths in central Alberta, spanning across Proterozoic accreted terranes and Archean microcontinents. This observation likely results from extensive plate convergence and crustal melting during the Proterozoic eon. The overall correlation between the crustal velocities and presumed basement domains is lower than expected, however. In the lower crust, the main pattern of shear velocities is relatively concordant with the reported domain boundaries and key Precambrian structures appear to be intact. The shear velocities beneath the Loverna Block, the largest constituent of the Hearne craton, are 10 per cent higher than the regional average. This prominent northeast striking seismic anomaly is moderately correlated with the regional heat flow and potentially represents the remnant core of the Archean Hearne province. The associated high velocities extend into the western part of the Medicine Hat Block, a possible Archean microcontinent with a debatable origin, and contribute to a strong east–west structural gradient in the lower crust. The presence and the continuity of this anomalous structure imply extensive communications among the various basement domains in southern Alberta during the assembly of the North American continent.

Description: The distribution of interevent times between two successive earthquakes is generally modelled as Gamma distribution. However, significant deviations from Gamma distribution have been observed at small interevent times. These deviations can be easily modelled through the introduction of a new parameter in the interevent time distribution. The log-likelihood estimations of the model parameters performed for the data from four earthquake catalogues reveal that the proposed model fits the data well. Finally, the new parameter is interpreted in terms of the c parameter of the Omori law.

Description: Seismic tomography has revealed very high P -wave velocities, over 8.5 km s –1 , at shallow depths, 30–100 km, beneath New Zealand. Here we study fast, high-frequency arrivals at North and South Island stations that contain additional information about the crust and mantle structure. These arrivals, which are from earthquakes within or close to the land mass, have a characteristic high-frequency precursor followed by a lower frequency, larger amplitude, main phase. Precursors were seen on at least one station from 262 of 306 candidate events; the best-recorded 76 events were analysed for wave speed, frequency content and polarization. Time–distance plots are consistent with two phases travelling at 8.38 ± 0.03 and 6.93 ± 0.05 km s –1 . The precursor has typical frequencies 4–9 Hz, the second arrival 2–4 Hz. Polarizations are off-azimuth by 30° and steeper than predicted by ray tracing through a smooth 3-D tomographic model. These results are explained by propagation through a dipping layer of order 10 km thick with seismic velocity around 8.5 km s –1 ; it is too thin to propagate frequencies below 4 Hz and waves refract from it at a steep, out-of-plane angle, explaining the anomalous polarization. Ray paths cover a region coinciding with the subducted Hikurangi Plateau; the fast layer is interpreted as the lowest section of the plateau that has transformed to eclogite, which has the same fast seismic velocity that we observe. Unlike the fast, eclogitic layers identified in subduction zones such as the Kermadecs, this layer is shallower, at 30 km, than the eclogite transformation; we therefore propose that it formed at the base of the thick plateau prior to subduction.

Description: Dispersion analysis of Rayleigh waves is performed to assess the velocity of complex structures such as sedimentary basins. At short periods several modes of the Rayleigh waves are often exited. To perform a reliable inversion of the velocity structure an identification of these modes is thus required. We propose a novel method to identify the modes of surface waves. We use the spectral ratio of the ground velocity for the horizontal components over the vertical component (H/V) measured on seismic coda. We then compare the observed values with the theoretical H/V ratio for velocity models deduced from surface wave dispersion when assuming a particular mode. We first invert the Rayleigh wave measurements retrieved from ambient noise cross-correlation with the assumptions that (1) the fundamental mode and (2) the first overtone are excited. Then we use these different velocity models to predict theoretical spectral ratios of the ground velocity for the horizontal components over the vertical component (H/V). These H/V ratios are computed under the hypothesis of equipartition of a diffuse field in a layered medium. Finally we discriminate between fundamental and higher modes by comparing the theoretical H/V ratio with the H/V ratio measured on seismic coda. In an application, we reconstruct Rayleigh waves from cross-correlations of ambient seismic noise recorded at seven broad-band stations in the Valley of Mexico. For paths within the soft quaternary sediments basin, the maximum energy is observed at velocities higher than expected for the fundamental mode. We identify that the dominant mode is the first higher mode, which suggests the importance of higher modes as the main vectors of energy in such complex structures.

Description: Maps of Rayleigh wave group velocity are presented across Brazil at periods between 6 and 23 s from seismic noise tomography. This range is complementary to previous continent-scale tomographic studies, which generally focus on periods longer than 10 s. By cross-correlating ambient noise recorded at 53 stations from 1996 to 2012, we produce empirical Green's functions between about 550 pairs of stations. From these, group velocity dispersion curves are measured by applying a frequency–time analysis with phase-matched filtering. At periods shorter than 23 s, 200–300 measurements pass the selection criteria, which are based on the signal-to-noise ratio and seasonal variability. Further filtering is performed by rejecting measured velocities too incompatible with maps produced from overdamped tomographic inversions. The final group velocity maps have a spatial resolution typically better than 400–500 km within most of the eastern half of Brazil, and slightly degrading toward shorter and longer periods of the study. The dispersion maps depict a strong correlation with sediment thickness and large-scale geological features, such as the São Francisco craton, surrounding mobile belts and intracratonic basins. The observed trends in group velocity and preliminary 1-D shear velocity models inferred from a Markov chain Monte Carlo exploration at three selected locations, are consistent with previous local seismic studies of the crustal structure. However, comparison with velocities from continent-scale earthquake tomography reveals diverging trends in some regions that deserve further investigation. As broad-band station coverage is rapidly expanding in Brazil, a full inversion to retrieve the 3-D crustal and uppermost mantle structure from ambient noise tomography should become possible in the near future.

Description: We present a comprehensive study of thickness and composition of the crust; and the nature of crust–mantle boundary beneath Southern India using P -wave receiver function from 119 seismic stations. Data from distributed network of seismograph location encompass geological domains like mid to late Archean Dharwar craton, Archean and Proterozoic metamorphic terrains, Proterozoic basin, rifted margins and escarpments, and Deccan volcanics. Except for the mid to lower crust exhumed Archean terrains (of West Dharwar and Southern Granulite) all other geological domains have crustal thickness in the range 33–40 km. In the western Dharwar, crustal thickness increases from ~40 km in the north to over 50 km in the south. The Archean domain of granulite terrain is thicker (40–45 km) and more mafic compared to its counterpart in south deformed at 550 Ma. Most of the crustal blocks have low to moderate Vp / Vs (1.72–1.76) representing a felsic to intermediate composition. Exception to the above include Archean granulite terrain with high Vp / Vs (1.76–1.81) suggestive of more mafic crust beneath them. When accounted for the paleo burial depth of 15–25 km, the study suggests a possible Himalaya–Tibet like scenario beneath the mid-late Archean in southwestern Dharwar and north granulite terrain whose deeper crust has progressively densified. This led to a gradational crust–mantle transition that is otherwise sharp elsewhere. The study suggests a more homogenized and felsic nature of the Precambrian crust beneath the terrains formed after 2.6 Ga, possibly due to delamination of the mafic lower crust. Our study does not suggest any distinction between late Archean and Proterozoic crust. The Deccan volcanism at 65 Ma does not appear to have altered the crustal character beneath it and is similar to the adjoining late Archean east Dharwar craton. The western Ghat escarpment and the coastal plain formed due to separation of India from Madagascar are underlain by mafic lower crust.

Description: We have developed a novel method to detect and locate geophysical events that makes use of any sufficiently dense sensor network. This method is demonstrated using acoustic sensor data collected in 2013 at the USArray Transportable Array (TA). The algorithm applies Delaunay triangulation to divide the sensor network into a mesh of three-element arrays, called triads. Because infrasound waveforms are incoherent between the sensors within each triad, the data are transformed into envelopes, which are cross-correlated to find signals that satisfy a consistency criterion. The propagation azimuth, phase velocity and signal arrival time are computed for each signal. Triads with signals that are consistent with a single source are bundled as an event group. The ensemble of arrival times and azimuths of detected signals within each group are used to locate a common source in space and time. A total of 513 infrasonic stations that were active for part or all of 2013 were divided into over 2000 triads. Low (0.5–2 Hz) and high (2–8 Hz) catalogues of infrasonic events were created for the eastern USA. The low-frequency catalogue includes over 900 events and reveals several highly active source areas on land that correspond with coal mining regions. The high-frequency catalogue includes over 2000 events, with most occurring offshore. Although their cause is not certain, most events are clearly anthropogenic as almost all occur during regular working hours each week. The regions to which the TA is most sensitive vary seasonally, with the direction of reception dependent on the direction of zonal winds. The catalogue has also revealed large acoustic events that may provide useful insight into the nature of long-range infrasound propagation in the atmosphere.

Description: For dynamic rupture problems, numerical simulation methods, such as the finite-difference method, the finite-element method and the boundary integral element method, usually produce spurious high-frequency oscillations that are mainly generated by discontinuities in the friction law and poor resolution of the breakdown zone. Techniques have been developed to reduce the oscillations; for example, the application of a damping coefficient, the introduction of a Green's function with higher accuracy and the use of a high-frequency filter. Presently, the spectral element method (SEM) is an important method used to simulate strong ground motion because of its high precision in calculations and flexibility in gridding media. Its greatest advantage is that it applies the orthogonal property of Gauss–Lobatto–Legendre points to form a diagonal mass matrix and is thus suitable for parallel computation that greatly reduces the computational time. However, comparisons made in the SCEC/USGS Spontaneous Rupture Code Verification Project show that the SEM has larger high-frequency oscillations than some other numerical methods for dynamic rupture problems. In this paper, we propose a new time-marching scheme of the SEM that has the frequency response of suppressing high-frequency oscillations for the slip-weakening friction law. Computation in rupture problem illustrates that the scheme greatly reduces spurious high-frequency oscillations. Furthermore, in the Appendix of the paper we provide some formula derivation to distinguish our scheme from generalized velocity schemes.

Description: The problem of decomposing irregular data on the sphere into a set of spherical harmonics is common in many fields of geosciences where it is necessary to build a quantitative understanding of a globally varying field. For example, in global seismology, a compressional or shear wave speed that emerges from tomographic images is used to interpret current state and composition of the mantle, and in geomagnetism, secular variation of magnetic field intensity measured at the surface is studied to better understand the changes in the Earth's core. Optimization methods are widely used for spherical harmonic analysis of irregular data, but they typically do not treat the dependence of the uncertainty estimates on the imposed regularization. This can cause significant difficulties in interpretation, especially when the best-fit model requires more variables as a result of underestimating data noise. Here, with the above limitations in mind, the problem of spherical harmonic expansion of irregular data is treated within the hierarchical Bayesian framework. The hierarchical approach significantly simplifies the problem by removing the need for regularization terms and user-supplied noise estimates. The use of the corrected Akaike Information Criterion for picking the optimal maximum degree of spherical harmonic expansion and the resulting spherical harmonic analyses are first illustrated on a noisy synthetic data set. Subsequently, the method is applied to two global data sets sensitive to the Earth's inner core and lowermost mantle, consisting of PKPab-df and PcP-P differential traveltime residuals relative to a spherically symmetric Earth model. The posterior probability distributions for each spherical harmonic coefficient are calculated via Markov Chain Monte Carlo sampling; the uncertainty obtained for the coefficients thus reflects the noise present in the real data and the imperfections in the spherical harmonic expansion.

Description: The Ricker wavelet, which is often employed in seismic analysis, has a symmetrical form. Seismic wavelets observed from field data, however, are commonly asymmetric with respect to the time variation. In order to better represent seismic signals, asymmetrical wavelets are defined systematically as fractional derivatives of a Gaussian function in which the Ricker wavelet becomes just a special case with the integer derivative of order 2. The fractional value and a reference frequency are two key parameters in the generalization. Frequency characteristics, such as the central frequency, the bandwidth, the mean frequency and the deviation, may be expressed analytically in closed forms. In practice, once the statistical properties (the mean frequency and deviation) are numerically evaluated from the discrete Fourier spectra of seismic data, these analytical expressions can be used to uniquely determine the fractional value and the reference frequency, and subsequently to derive various frequency quantities needed for the wavelet analysis. It is demonstrated that field seismic signals, recorded at various depths in a vertical borehole, can be closely approximated by generalized wavelets, defined in terms of fractional values and reference frequencies.

Description: Normal-mode spectra may be used to investigate large-scale elastic and anelastic heterogeneity throughout the entire Earth. The relevant theory was developed a few decades ago, however—mainly due to computational limitations—several approximations are commonly employed, and thus far the full merits of the complete theory have not been taken advantage of. In this study, we present an exact algebraic form of the theory for an aspherical, anelastic and rotating Earth model in which either complex or real spherical harmonic bases are used. Physical dispersion is incorporated into the quadratic eigenvalue problem by expanding the logarithmic frequency term to second-order. Proper (re)normalization of modes in a 3-D Earth model is fully considered. Using a database of 41 earthquakes and more than 10 000 spectra containing 116 modes with frequencies less than 3 mHz, we carry out numerical experiments to quantitatively evaluate the accuracy of commonly used approximate mode synthetics. We confirm the importance of wideband coupling, that is, fully coupling all modes below a certain frequency. Neither narrowband coupling, in which nearby modes are grouped into isolated clusters, nor self-coupling, that is, incorporating coupling between singlets within the same multiplet, are sufficiently accurate approximations. Furthermore, we find that (1) effects of physical dispersion can be safely approximated based on either a fiducial frequency approximation or a quadratic approximation of the logarithmic dispersion associated with the absorption-band model; (2) neglecting the proper renormalization of the modes of a rotating, anelastic Earth model introduces only minor errors; (3) ignoring the frequency dependence of the Coriolis and kinematic matrices in a wideband coupling scheme can lead to ~6 per cent errors in mode spectra at the lowest frequencies; notable differences also occur between narrowband coupling and quasi-degenerate perturbation theory, which linearizes the eigenvalue problem as well.

Description: The region of West Bohemia/Vogtland in the Czech–German border area is well known for the repeated occurrence of earthquake swarms, CO 2 emanations and mofette fields. We present a local earthquake tomography study undertaken to image the Vp and Vp/Vs structure in the broader area of earthquake swarm activity. In comparison with previous investigations, more details of the near-surface geology, potential fluid pathways and features around and below the swarm focal zone could be revealed. In the uppermost crust, for the first time the Cheb basin and the Bublák/Hartoušov mofette fields were imaged as distinct anomalies of Vp and Vp/Vs. The well-pronounced low-Vp anomaly of the Cheb basin is not continuing into the Eger rift indicating a particular role of the basin within the rift system. A steep channel of increased Vp/Vs is interpreted as the pathway for fluids ascending from the earthquake swarm focal zone up to the Bublák/Hartoušov mofette fields. As a new feature, a mid-crustal body of high Vp and increased Vp/Vs is revealed just below and north of the earthquake swarm focal zone. It may represent a solidified intrusive body which emplaced prior or during the formation of the rift system. We speculate that enhanced fluid flow into the focal zone and triggering of earthquakes could be driven by the presence of the intrusive body if cooling is not fully completed. We consider the assumed intrusive structure as a heterogeneity leading to higher stress particularly at the junction of the rift system with the basin and prominent fault structures. This may additionally contribute to the triggering of earthquakes.

Description: We develop and apply an algorithm for deriving interstation seismic attenuation from cross-correlations of ambient noise recorded by linear arrays. Theoretical results on amplitude decay due to attenuation are used to form a linear least-square inversion for interstation Q R values of Rayleigh surface waves propagating along linear arrays having three or more stations. The noise wave field is assumed stationary within each day and the interstation distances should be greater than the employed wavelength. The inversion uses differences of logarithmic amplitude decay curves measured at different stations from cross-correlation functions within a given frequency band. The background attenuation between noise sources and receivers is effectively cancelled with this method. The site amplification factors are assumed constant (or following similar patterns) in the frequency band of interest. The inversion scheme is validated with synthetic tests using ambient noise generated by ray-theory-based calculations with heterogeneous attenuation and homogenous velocity structure. The interstation attenuation and phase velocity dispersion curves are inverted from cross-correlations of the synthetic data. The method is then applied to triplets of stations from the regional southern California seismic network crossing the Mojave section of the San Andreas fault, and a dense linear array crossing the southern San Jacinto Fault zone. Bootstrap technique is used to derive empirical mean and confidence interval for the obtained inverse Q values. The results for the regional stations yield Q R values around 25 for a frequency band 0.2–0.36 Hz. The results for the San Jacinto fault zone array give Q R values of about 6–30 for frequencies in the range 15–25 Hz.

Description: We present a new method to locate low-frequency earthquakes (LFEs) within tectonic tremor episodes based on time-reverse imaging techniques. The modified time-reverse imaging technique presented here is the first method that locates individual LFEs within tremor episodes within 5 km uncertainty without relying on high-amplitude P -wave arrivals and that produces similar hypocentral locations to methods that locate events by stacking hundreds of LFEs without having to assume event co-location. In contrast to classic time-reverse imaging algorithms, we implement a modification to the method that searches for phase coherence over a short time period rather than identifying the maximum amplitude of a superpositioned wavefield. The method is independent of amplitude and can help constrain event origin time. The method uses individual LFE origin times, but does not rely on a priori information on LFE templates and families. We apply the method to locate 34 individual LFEs within tremor episodes that occur between 2010 and 2011 on the San Andreas Fault, near Cholame, California. Individual LFE location accuracies range from 2.6 to 5 km horizontally and 4.8 km vertically. Other methods that have been able to locate individual LFEs with accuracy of less than 5 km have mainly used large-amplitude events where a P -phase arrival can be identified. The method described here has the potential to locate a larger number of individual low-amplitude events with only the S -phase arrival. Location accuracy is controlled by the velocity model resolution and the wavelength of the dominant energy of the signal. Location results are also dependent on the number of stations used and are negligibly correlated with other factors such as the maximum gap in azimuthal coverage, source–station distance and signal-to-noise ratio.

Description: Resolving the topography of the core–mantle boundary (CMB) and the structure and composition of the D '' region is key to improving our understanding of the interaction between the Earth's mantle and core. Observations of traveltimes and amplitudes of short-period teleseismic body waves sensitive to lowermost mantle provide essential constraints on the properties of this region. Major challenges are low signal-to-noise ratio of the target phases and interference with other mantle phases. In a previous paper (Part I), we introduced the slant-stacklet transform to enhance the signal of the core-reflected ( PcP ) phase and to isolate it from stronger signals in the coda of the P wave. Then we minimized a linear misfit between P and PcP waveforms to improve the quality of PcP–P traveltime difference measurements as compared to standard cross-correlation methods. This method significantly increases the quantity and the quality of PcP–P traveltime observations available for the modelling of structure near the CMB. Here we illustrate our approach in a series of regional studies of the CMB and D '' using PcP–P observations with unprecedented resolution from high-quality dense arrays located in North America and Japan for events with magnitude M w 〉5.4 and distances up to 80 $\deg$ . In this process, we carefully analyse various sources of errors and show that mantle heterogeneity is the most significant. We find and correct bias due to mantle heterogeneities that is as large as 1 s in traveltime, comparable to the largest lateral PcP–P traveltime variations observed. We illustrate the importance of accurate mantle corrections and the need for higher resolution mantle models for future studies. After optimal mantle corrections, the main signal left is relatively long wavelength in the regions sampled, except at the border of the Pacific large-low shear velocity province (LLSVP). We detect the northwest border of the Pacific LLSVP in the western Pacific from array observations in Japan, and observe higher than average P velocities, or depressed CMB, in Central America, and slightly lower than average P velocities under Alaska/western Canada.

Description: In order to characterize the subsurface structure of the Jakarta Basin, Indonesia, a dense portable seismic broad-band network was operated by The Australian National University (ANU) and the Indonesian Agency for Meteorology, Climatology and Geophysics (BMKG) between October 2013 and February 2014. Overall 96 locations were sampled through successive deployments of 52 seismic broad-band sensors at different parts of the city. Oceanic and anthropogenic noises were recorded as well as regional and teleseismic earthquakes. We apply regularized deconvolution to the recorded ambient noise of the vertical components of available station pairs, and over 3000 Green's functions were retrieved in total. Waveforms from interstation deconvolutions show clear arrivals of Rayleigh fundamental and higher order modes. The traveltimes that were extracted from group velocity filtering of fundamental mode Rayleigh wave arrivals, are used in a 2-stage Transdimensional Bayesian method to map shear wave structure of subsurface. The images of S wave speed show very low velocities and a thick basin covering most of the city with depths up to 1.5 km. These low seismic velocities and the thick basin beneath the city potentially cause seismic amplification during a subduction megathrust or other large earthquake close to the city of Jakarta.

Description: We develop an indirect boundary integral equation method (IBIEM) to solve the scattering of seismic waves by a 3-D layered alluvial basin. We adopt the dynamic Green's functions for concentrated loads for a layered half-space derived from the modified stiffness method. This new algorithm of Green's function can solve the near-source response efficiently and accurately, and also facilitates the meshless implementation of the IBIEM. The numerical accuracy and stability of the IBIEM are tested for a homogeneous, hemispherical alluvial basin, and a two-layered model. Based on the IBIEM, the effects of several important parameters, such as the incident frequency, the angle of incidence and the properties of the alluvial layers are investigated for incident plane P and SV waves, respectively. The results show that the local amplification effects of a 3-D layered alluvial basin on the ground motion are strikingly significant, and that the spatial variation of the displacement response is drastic. We also find that the thickness of the near-surface low-velocity alluvial layer has a pronounced influence on the frequency spectrum of ground motion within the basin. As for the thick low-velocity layer, the amplification effect on the displacement amplitude spectrum appears in a wide range of frequencies, with more resonant models in the same frequency range. As for the thin low-velocity layer, in contrast, the amplification effect is close to the homogeneous case but becomes more significant for high-frequency waves. The displacement amplification for a basin with a soft intermediate layer is larger than that of the homogeneous basin for the lower frequencies, but seems to be weakened for high-frequency waves. Additionally, the damping ratio of the alluvial layer can substantially reduce the displacement amplitude in the basin, especially in the range of resonant frequencies. Our results provide a better understanding of the 3-D wave focusing and basin-edge effect within 3-D layered alluvial basins.

Description: Robust event detection and picking is a prerequisite for reliable (micro-) seismic interpretations. Detection of weak events is a common challenge among various available event detection algorithms. In this paper we compare the performance of two event detection methods, the short-term average/long-term average (STA/LTA) method, which is the most commonly used technique in industry, and a newly introduced method that is based on the power spectral density (PSD) measurements. We have applied both techniques to a 1-hr long segment of the vertical component of some raw continuous data recorded at a borehole geophone in a hydraulic fracturing experiment. The PSD technique outperforms the STA/LTA technique by detecting a higher number of weak events while keeping the number of false alarms at a reasonable level. The time–frequency representations obtained through the PSD method can also help define a more suitable bandpass filter which is usually required for the STA/LTA method. The method offers thus much promise for automated event detection in industrial, local, regional and global seismological data sets.

Description: Several recent studies have demonstrated the importance of crustal corrections when inverting surface wave data to model lateral variations in mantle radial anisotropy. It has also been shown that the choice of the prior crustal model to correct the data can strongly influence the anisotropy model and potentially lead to different geodynamic interpretations. In comparing tomographic models of radial anisotropy obtained from different crustal corrections, these studies did not, however, determine quantitative model uncertainties. Nevertheless, mantle models resulting from different prior crustal corrections are statistically different only if the posterior model errors stemming from the non-uniqueness of the inverse problem are smaller than the effect of the crustal correction itself. Here, we applied a model space search approach to global fundamental and higher mode Rayleigh and Love wave phase velocity maps to determine reliable, quantitative model uncertainties on seismic velocities and radial anisotropy. The technique employed enabled us to describe the model space with a posterior probability density function, and therefore to test whether models obtained from different crustal corrections are statistically different. We thus assessed the significance of the choice of the crustal model by comparing the posterior model errors to the differences in mantle structure resulting from different crustal corrections. We tested prior crustal models CRUST2.0, CRUST1.0 and 3SMAC. Our study shows that the use of prior crustal corrections from different crustal models yields significant discrepancies in mantle velocities around 50 km depth and in radial anisotropy down to 100 km. The impact of the crustal correction on radial anisotropy can extend down to 250 km in some locations. We found that choosing 3SMAC instead of the other crustal models has a stronger influence on the mantle model, but that CRUST1.0 and CRUST2.0 yield statistically identical anisotropy models at all depths, except at a few grid cells. Importantly, the effect of the crustal model is most significant in continental regions and not so much beneath oceans, which has important consequences for determining the depth of continental roots. Our results therefore suggest that improving constraints on crustal structure in continents is essential for our understanding of continent formation. Our work also demonstrates that the prior crustal model does not significantly affect radial anisotropy and velocities at depths greater than 100 km. This implies that if geodynamic interpretations of radial anisotropy below 100 km depth were to account for tomographic model uncertainties, they would not depend on the choice of the prior crustal model. It is therefore important for geodynamicists and seismologists to work in concert and to put effort into determining quantitative tomographic model uncertainties before interpreting the results. Our results also caution against the use of 3SMAC to correct surface wave data for studies of the continental lithosphere and suggest that the solid Earth community would benefit from putting some efforts towards building a revised 3SMAC. The discrepancies between mantle models built based on 3SMAC crustal corrections and those based on CRUST1.0 or CRUST2.0 should also help shed light on the validity of the geodynamical assumptions made in the construction of models like 3SMAC.

Description: Stimulated by the recent advances in computational tools for the simulation of seismic wave propagation problems in realistic geological configurations, this paper presents a 3D physics-based numerical analysis of near-source ground motion during the M W 6.0 2012 May 29 earthquake in the Po Plain, Northern Italy. To reproduce with sufficient accuracy some of the most peculiar features of the near-source strong-motion records and of the spatial variability of damage distribution, this study required a sequence of investigations, starting from the analysis of a wide set of near-source records, to the calibration of an improved kinematic seismic source model, up to the development of a 3D numerical model of the portion of the Po Plain interested by the earthquake. The latter includes the basin geometry, characterized by sediment thickness sharply varying from few tens of metres to some kilometres. The spatial resolution of the numerical model is suitable to propagate frequencies up to about 1.5 Hz. Numerical simulations were performed using the open-source high-performance code SPEED, based on the Discontinuous Galerkin Spectral Elements method. The 3D numerical model, coupled with the updated slip distribution along the rupturing fault, proved successful to reproduce with good agreement, measured through quantitative goodness-of-fit criteria, the most relevant features of the observed ground motion. These include: (i) the large fault normal velocity peaks at the near-source stations driven by updip directivity effects; (ii) the small-scale variability at short distance from the source, resulting in the out-of-phase motion at stations separated by only 3 km distance; (iii) the propagation of prominent trains of surface waves, especially in the Northern direction; (iv) the map of earthquake-induced ground uplift with maximum values of about 10 cm, in substantial agreement with geodetic measurements and (v) the two-lobed pattern of the peak ground velocity map, well correlated with the distribution of macroseismic intensity.

Description: This paper introduces a calculation method for the effective elastic stiffness tensor matrix of the viscous-elastic TTI medium based on the Chapman theory. We then obtain the phase velocity formula and seismic wave polarization formula of the viscous-elastic TTI medium, by solving the Christoffel equation; solve the phase angle of reflection and transmission wave through the numerical method in accordance with the wave slowness ellipsoid; on the basis of this assumption, and assuming that qP , qS and SH waves occurred simultaneously at the viscous-elastic anisotropic interface, establish the sixth-order Zoeppritz equation in accordance with the boundary conditions; establish the models for the upper and lower media which are viscous-elastic HTI, TTI, etc., on the basis of the sixth-order Zoeppritz equation; and study the impact of fracture dip angle, azimuth angle and frequency on the reflection coefficient. From this we obtain the following conclusions: the reflection coefficient can identify the fracture strike and dip when any information pertaining to the media is unknown; dispersion phenomenon is obvious on the axial plane of symmetry and weakened in the plane vertical to the axial plane of symmetry; the vertical-incidence longitudinal wave can stimulate the qS wave when the dip angle is not 0° or 90° under the condition of coincidence between the symmetry planes of the upper and lower media; when the symmetry planes of the upper and lower media do not coincide and the dip angle is not 0° or 90°, then the vertical-incidence qP will stimulate the qS and SH waves at the same time; the dip angle can cause the reflection coefficient curve to have a more obvious dispersion phenomenon, while the included angle between the symmetry planes of the upper and lower media will weaken the dispersion except SH ; and the intercept of reflection coefficient is affected by the fracture dip and included angle between the symmetry planes of the upper and lower media.

Description: I perform a retrospective forecast experiment in the most rapid extensive continental rift worldwide, the western Corinth Gulf (wCG, Greece), aiming to predict shallow seismicity (depth 〈15 km) with magnitude M ≥ 3.0 for the time period between 1995 and 2013. I compare two short-term earthquake clustering models, based on epidemic-type aftershock sequence (ETAS) statistics, four physics-based (CRS) models, combining static stress change estimations and the rate-and-state laboratory law and one hybrid model. For the latter models, I incorporate the stress changes imparted from 31 earthquakes with magnitude M ≥ 4.5 at the extended area of wCG. Special attention is given on the 3-D representation of active faults, acting as potential receiver planes for the estimation of static stress changes. I use reference seismicity between 1990 and 1995, corresponding to the learning phase of physics-based models, and I evaluate the forecasts for six months following the 1995 M = 6.4 Aigio earthquake using log-likelihood performance metrics. For the ETAS realizations, I use seismic events with magnitude M ≥ 2.5 within daily update intervals to enhance their predictive power. For assessing the role of background seismicity, I implement a stochastic reconstruction (aka declustering) aiming to answer whether M 〉 4.5 earthquakes correspond to spontaneous events and identify, if possible, different triggering characteristics between aftershock sequences and swarm-type seismicity periods. I find that: (1) ETAS models outperform CRS models in most time intervals achieving very low rejection ratio R N = 6 per cent, when I test their efficiency to forecast the total number of events inside the study area, (2) the best rejection ratio for CRS models reaches R N = 17 per cent, when I use varying target depths and receiver plane geometry, (3) 75 per cent of the 1995 Aigio aftershocks that occurred within the first month can be explained by static stress changes, (4) highly variable performance on behalf of both statistical and physical models is suggested by large confidence intervals of information gain per earthquake and (5) generic ETAS models can adequately predict the temporal evolution of seismicity during swarms. Furthermore, stochastic reconstruction of seismicity makes possible the identification of different triggering processes between specific seismic crises (2001, 2003–04, 2006–07) and the 1995 aftershock sequence. I find that: (1) seismic events with M ≥ 5.0 are not a part of a preceding earthquake cascade, since they are characterized by high probability being a background event (average P back 〉 0.8) and (2) triggered seismicity within swarms is characterized by lower event productivity when compared with the corresponding value during aftershock sequences. I conclude that physics-based models contribute on the determination of the ‘new-normal’ seismicity rate at longer time intervals and that their joint implementation with statistical models is beneficial for future operational forecast systems.

Description: In this study, we obtain a set of 2-D global phase velocity and attenuation maps for Rayleigh waves between 5 and 25 mHz. Correcting the effect of focusing–defocusing is crucial in order to obtain reliable attenuation structure. Great circle linearized ray theory, which has been used to date, can give useful predictions of this effect if careful attention is paid to how the phase velocity model is smoothed. In contrast, predictions based on the 2-D finite-frequency kernels are quite robust in this frequency range and suggest that they are better suited as a basis for inversion. We use a large data set of Rayleigh wave phase and amplitude measurements to invert for the phase velocity, attenuation, source and receiver terms simultaneously. Our models provide 60–70 per cent variance reduction to the raw data though the source terms are the biggest contribution to the fit of the data. The attenuation maps show structures that correlate well with surface tectonics and the age progression trend of the attenuation is clearly seen in the ocean basins. We have also identified problematic stations and earthquake sources as a by-product of our data selection process.

Description: The crustal and upper mantle velocity structure in the northeastern Tibetan Plateau is obtained from joint analysis of receiver functions and Rayleigh wave dispersions. The resulting velocity model reveals a close correlation between the thick (〉60 km) crust and the presence of an intracrustal low-velocity zone beneath the Qiangtang and Songpan-Ganzi terranes as well as the northwestern Qilian orogen. However, the high V p / V s ratio of the crust is found only beneath the Qiangtang and Songpan-Ganzi terranes. The crustal low velocity zone does not appear in the west Qinling and southeastern Qilian orogens, which have a relatively thin (~50 km) crust, indicating that crustal channel flow is not the primary mechanism by which the northeastern Tibetan Plateau grows. A continuous low velocity zone from the mid-to-lower crust down to 160 km beneath the eastern Kunlun fault suggests an induced local mantle upwelling after partial detachment of the lithosphere.

Description: Using the spectral-element method, we explored the effect of topography of upper-mantle discontinuities on the traveltimes of SS precursors recorded on transverse component seismograms. The latter are routinely used to infer the topography of mantle transition zone discontinuities. The step from precursory traveltimes to topographic changes is mainly done using linearised ray theory, or sometimes using finite-frequency kernels. We simulated exact seismograms in 1-D and 3-D elastic models of the mantle. In a second simulation, we added topography to the discontinuities. We compared the waveforms obtained with and without topography by cross correlation of the SS precursors. Since we did not add noise, the precursors are visible in individual seismograms without the need of stacking. The resulting time anomalies were then converted into topographic variations and compared to the original topographic models. Based on the correlation between initial and inferred models, and provided that ray coverage is good, we found that linearised ray theory gives a relatively good idea on the location of the uplifts and depressions of the discontinuities. It seriously underestimates the amplitude of the topographic variations by a factor ranging between 2 and 7. Real data depend on the 3-D elastic structure and the topography. All studies to date correct for the 3-D elastic effects assuming that the traveltimes can be linearly decomposed into a structure and a discontinuity part. We found a strong non-linearity in this decomposition which cannot be modelled without a fully non-linear inversion for elastic structure and discontinuities simultaneously.

Description: We present a novel optimization approach to improve the convergence of interstation coda correlation functions towards the medium's empirical Green's function. For two stations recording a series of impulsive events in a multiply scattering medium, we explore the impact of coda window selection through a Markov Chain Monte Carlo scheme, with the aim of generating a gather of correlation functions that is the most coherent and symmetric over events, thus recovering intuitive elements of the interstation Green's function without any nonlinear post-processing techniques. This approach is tested here for a 2-D acoustic finite difference model, where a much improved correlation function is obtained, as well as for a database of small impulsive icequakes recorded on Erebus Volcano, Antarctica, where similar robust results are shown. The average coda solutions, as deduced from the posterior probability distributions of the optimization, are further representative of the scattering strength of the medium, with stronger scattering resulting in a slightly delayed overall coda sampling. The recovery of singly scattered arrivals in the coda of correlation functions are also shown to be possible through this approach, and surface wave reflections from outer craters on Erebus volcano were mapped in this fashion. We also note that, due to the improvement of correlation functions over subsequent events, this approach can further be used to improve the resolution of passive temporal monitoring.

Description: In this study we derive a spectral model describing the source, propagation and site characteristics of S waves recorded in central Italy. To this end, we compile and analyse a high-quality data set composed of more than 9000 acceleration and velocity waveforms in the local magnitude ( M l ) range 3.0–5.8 recorded at epicentral distances smaller than 120 km. The data set spans the time period from 2008 January 1 to 2013 May 31, and includes also the 2009 L'Aquila (moment magnitude M w 6.1, M l = 5.8) sequence. This data set is suitable for the application of data-driven approaches to derive the empirical functions for source, attenuation and site terms. Therefore, we apply a non-parametric inversion scheme to the acceleration Fourier spectra of the S waves of 261 earthquakes recorded at 129 stations. In a second step, with the aim of defining spectral models suitable for the implementation in numerical simulation codes, we represent the obtained non-parametric source and propagation terms by fitting standard parametric models. The frequency-dependent attenuation with distance r shows a complex trend that we parametrize in terms of geometrical spreading, anelastic attenuation and high-frequency decay parameter k. The geometrical spreading term is described by a piecewise linear model with crossover distances at 10 and 70 km: in the first segment, the spectral ordinates decay as 〈 tex – mathid = " IM 0001" 〉 r – 1.01 while in the second as 〈 tex – mathid = " IM 0002" 〉 r – 1.68 . Beyond 70 km, the attenuation decreases and the spectral amplitude attenuate as 〈 tex – mathid = " IM 0003" 〉 r – 0.64 . The quality factor Q ( f ) and the high-frequency attenuation parameter k , are 〈 tex – mathid = " IM 0004" 〉 Q ( f ) = 290 f 0.16 and k = 0.012 s, respectively, the latter being applied only for frequencies higher than 10 Hz. The source spectra are well described by 2 models, from which seismic moment and stress drops of 231 earthquakes are estimated. We calibrate a new regional relationship between seismic moment and local magnitude that improves the existing ones and extends the validity range to 3.0–5.8. We find a significant stress drop increase with seismic moment for events with M w larger than 3.75, with so-called scaling parameter  close to 1.5. We also observe that the overall offset of the stress-drop scaling is controlled by earthquake depth. We evaluate the performance of the proposed parametric models through the residual analysis of the Fourier spectra in the frequency range 0.5–25 Hz. The results show that the considered stress-drop scaling with magnitude and depth reduces, on average, the standard deviation by 18 per cent with respect to a constant stress-drop model. The overall quality of fit (standard deviation between 0.20 and 0.27, in the frequency range 1–20 Hz) indicates that the spectral model calibrated in this study can be used to predict ground motion in the L'Aquila region.

Description: This paper examines the factors contributing to the ‘sequencing’ of tsunami waves in the far field, that is, to the distribution of the maximum sea surface amplitude inside the dominant wave packet constituting the primary arrival at a distant harbour. Based on simple models of sources for which analytical solutions are available, we show that, as range is increased, the wave pattern evolves from a regime of maximum amplitude in the first oscillation to one of delayed maximum, where the largest amplitude takes place during a subsequent oscillation. In the case of the simple, instantaneous uplift of a circular disk at the surface of an ocean of constant depth, the critical distance for transition between those patterns scales as $r_0^3 / h^2$ where r 0 is the radius of the disk and h the depth of the ocean. This behaviour is explained from simple arguments based on a model where sequencing results from frequency dispersion in the primary wave packet, as the width of its spectrum around its dominant period T 0 becomes dispersed in time in an amount comparable to T 0 , the latter being controlled by a combination of source size and ocean depth. The general concepts in this model are confirmed in the case of more realistic sources for tsunami excitation by a finite-time deformation of the ocean floor, as well as in real-life simulations of tsunamis excited by large subduction events, for which we find that the influence of fault width on the distribution of sequencing is more important than that of fault length. Finally, simulation of the major events of Chile (2010) and Japan (2011) at large arrays of virtual gauges in the Pacific Basin correctly predicts the majority of the sequencing patterns observed on DART buoys during these events. By providing insight into the evolution with time of wave amplitudes inside primary wave packets for far field tsunamis generated by large earthquakes, our results stress the importance, for civil defense authorities, of issuing warning and evacuation orders of sufficient duration to avoid the hazard inherent in premature calls for all-clear.

Description: We have successfully estimated the full moment tensors of 22 local earthquakes with local magnitude ranging from 1.2 to 4.8 that occurred in the Hungarian part of the Pannonian basin between 1995 and 2014. We used a probabilistic waveform inversion procedure that takes into account the effects of the random noise contained in the seismograms, the uncertainty of the hypocentre determined from arrival times and the inaccurate knowledge of the velocity structure, while estimating the error affecting the derived focal parameters. The applied probabilistic approach maps the posterior probability density functions (PPDFs) for both the hypocentral coordinates and the moment tensor components. The final estimates are given by the maximum likelihood points of the PPDFs, while solution uncertainties are presented by histogram plots. The estimated uncertainties in the moment tensor components are plotted on the focal sphere in such a way, that the significance of the double couple (DC), the compensated linear vector dipole (CLVD) and the isotropic (ISO) parts of the source can be assessed. We have shown that the applied waveform inversion method is equally suitable to recover the source mechanism for low-magnitude events using short-period local waveforms as well as for moderate-size earthquakes using long-period seismograms. The non-DC components of the retrieved focal mechanisms are statistically insignificant for all the analysed earthquakes. The negligible amount of the ISO component implies the tectonic nature of the investigated events. The moment tensor solutions reported by other agencies for five of the M L 〉 4 earthquakes studied in this paper are very similar to those calculated by the applied waveform inversion algorithm. We have found only strike-slip and thrust faulting events, giving further support to the hypothesis that the Pannonian basin is currently experiencing a compressional regime of deformation. The orientations of the obtained focal mechanisms are in good agreement with the main stress pattern published for the Pannonian region. The azimuth of the subhorizontal P principal axis varies from about NNE-SSW in SW Hungary through NE-SW well inside the basin to around E-W in the NE part of the country. Most of the analysed earthquakes occurred on faults or subfaults differently oriented than the main fault system.

Description: Using data predominantly from the NECESS Array, but also incorporating surface wave data from surrounding networks, we present the results of a Bayesian Monte Carlo inversion of receiver functions, Rayleigh wave ellipticity (H/V ratio) and Rayleigh wave group and phase speeds from 8 to 80 s period for the 3-D shear velocity structure of the crust and uppermost mantle beneath Northeast China. We define the final model as the mean and standard deviation of the posterior distribution at each location on a 0.5° x 0.5° grid from the surface to 150 km depth. The primary scientific motivation is to investigate the expression of intracontinental volcanism across the region. The model lithosphere displays prominent features (middle and lower crustal velocity, Moho depth, lithospheric thickness) across the study area that coincide with the location of volcanoes, which are predominantly situated in two distinct volcanic regions, which we call the ‘Northeast China Lineated Quaternary Volcanic Zone’, found near the eastern margin of the Songliao Basin and extending to Changbaishan Volcano, and the ‘Northern and Southern Greater Xing'an Range Pleistocene Volcanic Zones’. There is a strong similarity between the lateral distribution of depth-integrated mantle velocity anomalies in our model with the teleseismic body wave model of Tang et al ., although the vertical distribution of anomalies differ.

Description: Recent geodetic studies along the San Jacinto Fault (SJF) in southern California revealed a shallower locking depth than the seismogenic depth outlined by microseismicity. This disagreement leads to speculations that creeping episodes drive seismicity in the lower part of the seismogenic zone. Whether deep creep occurs along the SJF holds key information on how fault slips during earthquake cycle and potential seismic hazard imposed to southern California. Here we apply a matched filter technique to 10  M  〉 4 earthquake sequences along the SJF since 2000 and obtain more complete earthquake catalogues. We then systematic investigate spatio-temporal evolutions of these aftershock sequences. We find anomalously large aftershock zones for earthquakes occurred below the geodetically inferred locking depth (i.e. 11–12 km), while aftershock zones of shallower main shocks are close to expectations from standard scaling relationships. Although we do not observe clear migration of aftershocks, most aftershock zones do expand systematically with logarithmic time since the main shock. All the evidences suggest that aftershocks near or below the locking depth are likely driven by deep creep following the main shock. The presence of a creeping zone below 11–12 km may have significant implications on the maximum sizes of events in this region.

Description: We develop and apply a novel technique to image ambient seismic noise sources. It is based on measurements of cross-correlation asymmetry defined as the logarithmic energy ratio of the causal and anticausal branches of the cross-correlation function. A possible application of this technique is to account for the distribution of noise sources, a problem which currently poses obstacles to noise-based surface wave dispersion analysis and waveform inversion. The particular asymmetry measurement used is independent of absolute noise correlation amplitudes. It is shown how it can be forward-modelled and related to the noise source power-spectral density using adjoint methods. Simplified sensitivity kernels allow us to rapidly image variations in the power-spectral density of noise sources. This imaging method correctly accounts for viscoelastic attenuation and is to first order insensitive to unmodelled Earth structure. Furthermore, it operates directly on noise correlation data sets. No additional processing is required, which makes the method fast and computationally inexpensive. We apply the method to three vertical-component cross-correlation data sets of different spatial and temporal scales. Processing is deliberately minimal so as to keep observations consistent with the imaging concept. In accord with previous studies, we image seasonally changing sources of the Earth's hum in the Atlantic, Pacific and the Southern Ocean. The sources of noise in the microseismic band recorded at stations in Switzerland are predominantly located in the Atlantic and show a clear dependence on both season and frequency. Our developments are intended as a step towards full 3-D inversions for the sources of ambient noise in various frequency bands, which may ultimately lead to improvements of noise-based structural imaging.

Description: Peak delay and envelope broadening of an S -wavelet with travel distance increasing are seen in short-period seismograms of small earthquakes. Those phenomena are results of scattering by random velocity inhomogeneities in the earth medium. As shown in sonic well-log data we may suppose that random velocity fluctuation has power-law spectra even in the seismic spectral range. As a simple mathematical model, we study how the envelope of a scalar wavelet varies in von Kármán-type random media, which have power-law spectra at large wavenumbers. Since the centre wavenumber of a wavelet is a unique scale in the power-law spectral range, using it as a reference, we divide the random media into the low-wavenumber spectral (long-scale) component and the high-wavenumber spectral (short-scale) component. For the wave propagation through the long-scale component of random media, we may apply the parabolic approximation to the wave equation. Using the Markov approximation, which is a stochastic extension of the phase screen method, we directly synthesize the energy density, which is the mean-square (MS) envelope of a wavelet in a given frequency band. The envelope duration increases according to the second power of travel distance. There is an additional factor, the wandering effect which increases the envelope duration according to the traveltime fluctuation. Wide angle scattering caused by the short-scale component of random media attenuates wave amplitude with travel distance increasing. We use the total scattering coefficient of the short-scale component as a measure of scattering attenuation per distance, which is well described by the Born approximation. Multiplying the exponential scattering attenuation factor by the MS envelope derived by the Markov approximation, we can synthesize the MS envelope reflecting all the spectral components of random media. When the random medium power spectra have a steep role-off at large wavenumbers, the envelope broadening is small and frequency independent, and scattering attenuation is weak. When the random medium power spectra have a small role-off, however, the envelope broadening is large and increases with frequency, and the scattering attenuation is strong and increases with frequency. The proposed synthesis of MS envelopes is fully analytic; therefore, it can be a theoretical basis for the evaluation of random heterogeneity of the earth medium from the analysis of seismogram envelopes.

Description: Characterization of seismic sources is an important aspect of seismology. Parameter uncertainties in such inversions are essential for estimating solution robustness, but are rarely available. We have developed a non-linear moment tensor inversion method in a probabilistic Bayesian framework that also accounts for noise in the data. The method is designed for point source inversion using waveform data of moderate-size earthquakes and explosions at regional distances. This probabilistic approach results in an ensemble of models, whose density is proportional to parameter probability distribution and quantifies parameter uncertainties. Furthermore, we invert for noise in the data, allowing it to determine the model complexity. We implement an empirical noise covariance matrix that accounts for interdependence of observational errors present in waveform data. After we demonstrate the feasibility of the approach on synthetic data, we apply it to a Long Valley Caldera, CA, earthquake with a well-documented anomalous (non-double-couple) radiation from previous studies. We confirm a statistically significant isotropic component in the source without a trade-off with the compensated linear vector dipoles component.

Description: During arc-continent collision, buoyant sections of sediments and rifted continental crust from a subducting plate will accrete to the forearc of the upper plate as long as this backstop remains intact. Deformation of the oceanic arc and forearc block may ultimately lead to accretion of these mafic rock units to the new orogen. The Taiwan mountain belt, which formed at ~6.5 Ma by oblique convergence between the Eurasian passive margin and the overriding Luzon arc in northern Taiwan, offers important insight in this process, since the collision is more advanced in the north than in the south. The incipient stage of arc-collision can be studied in southern Taiwan, while the northern portion of the orogen is presently undergoing collapse due to a flip in the subduction polarity between the Eurasian Plate and the Philippine Sea Plate. In this study, we seismically image the structure of the northern section of the mountain belt with a tomographic inversion. We present marine and land-based seismic refraction data, as well as local earthquake data, from transect T6 of the Taiwan Integrated Geodynamic Research (TAIGER) program across the Taiwan mountain belt and the adjacent Ryukyu arc. Our 2-D compressional seismic velocity model for this transect, which is based on a tomographic inversion of 10 213 P -wave arrival times, shows that the Eurasian crystalline continental crust thickens from ~24 km in the Taiwan Strait to ~40 km beneath the eastern Central Range of Taiwan. The detailed seismic velocity structure of the Taiwan mountain belt shows vertical continuity in the upper 15 km, which suggests that rocks are exhumed to the surface here from the middle crust in a near-vertical path. The continental crust of the westernmost Ryukyu arc is almost as thick (~40 km) as in the adjacent northern Central Range of Taiwan, and it appears to override the leading edge of the Philippine Sea Plate offshore northeastern Taiwan. If we assume that the western Ryukyu arc crust also thickened in the collision, then the mountain belt is wider and less thick in northern Taiwan than in central Taiwan (~50 km), which may be the result of post-collisional extension in the north.

Description: The Canterbury, New Zealand, earthquake sequence, which began in September 2010, occurred in a region of low crustal deformation and previously low seismicity. Because, the ensuing seismicity in the region is likely to remain above previous levels for many years, a hybrid operational earthquake forecasting model for Canterbury was developed to inform decisions on building standards and urban planning for the rebuilding of Christchurch. The model estimates occurrence probabilities for magnitudes M ≥ 5.0 in the Canterbury region for each of the next 50 yr. It combines two short-term, two medium-term and four long-term forecasting models. The weight accorded to each individual model in the operational hybrid was determined by an expert elicitation process. A retrospective test of the operational hybrid model and of an earlier informally developed hybrid model in the whole New Zealand region has been carried out. The individual and hybrid models were installed in the New Zealand Earthquake Forecast Testing Centre and used to make retrospective annual forecasts of earthquakes with magnitude M 〉 4.95 from 1986 on, for time-lags up to 25 yr. All models underpredict the number of earthquakes due to an abnormally large number of earthquakes in the testing period since 2008 compared to those in the learning period. However, the operational hybrid model is more informative than any of the individual time-varying models for nearly all time-lags. Its information gain relative to a reference model of least information decreases as the time-lag increases to become zero at a time-lag of about 20 yr. An optimal hybrid model with the same mathematical form as the operational hybrid model was computed for each time-lag from the 26-yr test period. The time-varying component of the optimal hybrid is dominated by the medium-term models for time-lags up to 12 yr and has hardly any impact on the optimal hybrid model for greater time-lags. The optimal hybrid model is considerably more informative than the operational hybrid model at long time-lags, but less so when the period of the Canterbury earthquakes is excluded from the tests. The results highlight the value of including medium-term models and a range of long-term models in operational forecasting. Based on the tests carried out here, the operational hybrid model is expected to outperform most of the individual models in the next 25 yr.

Description: On 2008 May 29, two magnitude M W ~ 6 earthquakes occurred on two adjacent N-S faults in the western South Iceland Seismic Zone. The first main shock was followed approximately 3 s later by the rupture on a parallel fault, about 5 km to the west. An intense aftershock sequence was mostly confined to the western fault and an E-W aligned zone, extending west of the main shock region into the Reykjanes oblique rift. In this study, a total of 325 well-constrained focal mechanisms were obtained using data from the permanent Icelandic SIL seismic network and a temporary network promptly installed in the source region following the main shocks, which allowed a high-resolution stress inversion in short time intervals during the aftershock period. More than 800 additional focal mechanisms for the time period 2001–2009, obtained from the permanent SIL network, were analysed to study stress changes associated with the main shocks. Results reveal a coseismic counter-clockwise rotation of the maximum horizontal stress of 11 ± 10° (95 per cent confidence level) in the main rupture region. From previous fault models obtained by inversion of geodetic data, we estimate a stress drop of about half of the background shear stress on the western fault. With a stress drop of 8–10 MPa, the pre-event shear stress is estimated to 16–20 MPa. The apparent weakness of the western fault may be caused by fault properties, pore fluid pressure and the vicinity of the fault to the western rift zone, but may also be due to the dynamic stress increase on the western fault by the rupture on the eastern fault. Further, a coseismic change of the stress regime—from normal faulting to strike-slip faulting—was observed at the northern end of the western fault. This change could be caused by stress heterogeneities, but may also be due to a southward shift in the location of the aftershocks as compared to prior events.

Description: Coda-wave interferometry is a technique which exploits tiny waveform changes in the coda to detect temporal variations of seismic properties in evolving media. Observed waveform changes are of two kinds: traveltime perturbations and distortion of seismograms. In the last 10 yr, various theories have been published to relate either background velocity changes to traveltime perturbations, or changes in the scattering properties of the medium to waveform decorrelation. These theories have been limited by assumptions pertaining to the scattering process itself—in particular isotropic scattering, or to the propagation regime—single-scattering and/or diffusion. In this manuscript, we unify and extend previous results from the literature using a radiative transfer approach. This theory allows us to incorporate the effect of anisotropic scattering and to cover a broad range of propagation regimes, including the contribution of coherent, singly scattered and multiply scattered waves. Using basic physical reasoning, we show that two different sensitivity kernels are required to describe traveltime perturbations and waveform decorrelation, respectively, a distinction which has not been well appreciated so far. Previous results from the literature are recovered as limiting cases of our general approach. To evaluate numerically the sensitivity functions, we introduce an improved version of a spectral technique known as the method of ‘rotated coordinate frames’, which allows global evaluation of the Green's function of the radiative transfer equation in a finite domain. The method is validated through direct pointwise comparison with Green's functions obtained by the Monte Carlo method. To illustrate the theory, we consider a series of scattering media displaying increasing levels of scattering anisotropy and discuss the impact on the traveltime and decorrelation kernels. We also consider the related problem of imaging variations of scattering properties based on intensity perturbations observed in the coda. The impact of anisotropy is particularly pronounced for the scattering and decorrelation sensitivity kernels, which probe spatial/temporal changes in the scattering properties of the medium. Compared to the isotropic case, scattering anisotropy strongly increases the sensitivity of coda waves in the vicinity of the single-scattering ellipse, which may have important implications for imaging applications. In addition to demonstrating the impact of non-isotropic scattering on the sensitivity kernels of coda waves, our work offers a practical solution to model this process accurately.

Description: This paper examines the factors contributing to the ‘sequencing’ of tsunami waves in the far field, that is, to the distribution of the maximum sea surface amplitude inside the dominant wave packet constituting the primary arrival at a distant harbour. Based on simple models of sources for which analytical solutions are available, we show that, as range is increased, the wave pattern evolves from a regime of maximum amplitude in the first oscillation to one of delayed maximum, where the largest amplitude takes place during a subsequent oscillation. In the case of the simple, instantaneous uplift of a circular disk at the surface of an ocean of constant depth, the critical distance for transition between those patterns scales as $r_0^3 / h^2$ where r 0 is the radius of the disk and h the depth of the ocean. This behaviour is explained from simple arguments based on a model where sequencing results from frequency dispersion in the primary wave packet, as the width of its spectrum around its dominant period T 0 becomes dispersed in time in an amount comparable to T 0 , the latter being controlled by a combination of source size and ocean depth. The general concepts in this model are confirmed in the case of more realistic sources for tsunami excitation by a finite-time deformation of the ocean floor, as well as in real-life simulations of tsunamis excited by large subduction events, for which we find that the influence of fault width on the distribution of sequencing is more important than that of fault length. Finally, simulation of the major events of Chile (2010) and Japan (2011) at large arrays of virtual gauges in the Pacific Basin correctly predicts the majority of the sequencing patterns observed on DART buoys during these events. By providing insight into the evolution with time of wave amplitudes inside primary wave packets for far field tsunamis generated by large earthquakes, our results stress the importance, for civil defense authorities, of issuing warning and evacuation orders of sufficient duration to avoid the hazard inherent in premature calls for all-clear.

Description: In this study we derive a spectral model describing the source, propagation and site characteristics of S waves recorded in central Italy. To this end, we compile and analyse a high-quality data set composed of more than 9000 acceleration and velocity waveforms in the local magnitude ( M l ) range 3.0–5.8 recorded at epicentral distances smaller than 120 km. The data set spans the time period from 2008 January 1 to 2013 May 31, and includes also the 2009 L'Aquila (moment magnitude M w 6.1, M l = 5.8) sequence. This data set is suitable for the application of data-driven approaches to derive the empirical functions for source, attenuation and site terms. Therefore, we apply a non-parametric inversion scheme to the acceleration Fourier spectra of the S waves of 261 earthquakes recorded at 129 stations. In a second step, with the aim of defining spectral models suitable for the implementation in numerical simulation codes, we represent the obtained non-parametric source and propagation terms by fitting standard parametric models. The frequency-dependent attenuation with distance r shows a complex trend that we parametrize in terms of geometrical spreading, anelastic attenuation and high-frequency decay parameter k. The geometrical spreading term is described by a piecewise linear model with crossover distances at 10 and 70 km: in the first segment, the spectral ordinates decay as 〈 tex – mathid = " IM 0001" 〉 r – 1.01 while in the second as 〈 tex – mathid = " IM 0002" 〉 r – 1.68 . Beyond 70 km, the attenuation decreases and the spectral amplitude attenuate as 〈 tex – mathid = " IM 0003" 〉 r – 0.64 . The quality factor Q ( f ) and the high-frequency attenuation parameter k , are 〈 tex – mathid = " IM 0004" 〉 Q ( f ) = 290 f 0.16 and k = 0.012 s, respectively, the latter being applied only for frequencies higher than 10 Hz. The source spectra are well described by 2 models, from which seismic moment and stress drops of 231 earthquakes are estimated. We calibrate a new regional relationship between seismic moment and local magnitude that improves the existing ones and extends the validity range to 3.0–5.8. We find a significant stress drop increase with seismic moment for events with M w larger than 3.75, with so-called scaling parameter  close to 1.5. We also observe that the overall offset of the stress-drop scaling is controlled by earthquake depth. We evaluate the performance of the proposed parametric models through the residual analysis of the Fourier spectra in the frequency range 0.5–25 Hz. The results show that the considered stress-drop scaling with magnitude and depth reduces, on average, the standard deviation by 18 per cent with respect to a constant stress-drop model. The overall quality of fit (standard deviation between 0.20 and 0.27, in the frequency range 1–20 Hz) indicates that the spectral model calibrated in this study can be used to predict ground motion in the L'Aquila region.

Description: We present a novel optimization approach to improve the convergence of interstation coda correlation functions towards the medium's empirical Green's function. For two stations recording a series of impulsive events in a multiply scattering medium, we explore the impact of coda window selection through a Markov Chain Monte Carlo scheme, with the aim of generating a gather of correlation functions that is the most coherent and symmetric over events, thus recovering intuitive elements of the interstation Green's function without any nonlinear post-processing techniques. This approach is tested here for a 2-D acoustic finite difference model, where a much improved correlation function is obtained, as well as for a database of small impulsive icequakes recorded on Erebus Volcano, Antarctica, where similar robust results are shown. The average coda solutions, as deduced from the posterior probability distributions of the optimization, are further representative of the scattering strength of the medium, with stronger scattering resulting in a slightly delayed overall coda sampling. The recovery of singly scattered arrivals in the coda of correlation functions are also shown to be possible through this approach, and surface wave reflections from outer craters on Erebus volcano were mapped in this fashion. We also note that, due to the improvement of correlation functions over subsequent events, this approach can further be used to improve the resolution of passive temporal monitoring.

Description: Direct imaging of simultaneous-source (or blended) data, without the need of deblending, requires a precise subsurface velocity model. In this paper, we focus on the velocity analysis of simultaneous-source data using the normal moveout-based velocity picking approach.We demonstrate that it is possible to obtain a precise velocity model directly from the blended data in the common-midpoint domain. The similarity-weighted semblance can help us obtain much better velocity spectrum with higher resolution and higher reliability compared with the traditional semblance. The similarity-weighted semblance enforces an inherent noise attenuation solely in the semblance calculation stage, thus it is not sensitive to the intense interference. We use both simulated synthetic and field data examples to demonstrate the performance of the similarity-weighted semblance in obtaining reliable subsurface velocity model for direct migration of simultaneous-source data. The migrated image of blended field data using prestack Kirchhoff time migration approach based on the picked velocity from the similarity-weighted semblance is very close to the migrated image of unblended data.

Description: There has been debate in the Collaboratory for the Study of Earthquake Predictability project about the most appropriate form of the likelihood function to use to evaluate earthquake forecasts in specified discrete space–time intervals, and also to evaluate the validity of the model itself. The debate includes whether the likelihood function should be discrete in nature, given that the forecasts are in discrete space–time bins, or continuous. If discrete, can different bins be assumed to be statistically independent, and is it satisfactory to assume that the forecasted count in each bin will have a Poisson distribution? In order to discuss these questions, we start with the most simple models (homogeneous Poisson), and progressively develop the model complexity to include self exciting point process models. For each, we compare the discrete and continuous time likelihoods. Examples are given where it is proven that the counts in discrete space–time bins are not Poisson. We argue that the form of the likelihood function is intrinsic to the given model, and the required forecast for some specified space–time region simply determines where the likelihood function should be evaluated. We show that continuous time point process models where the likelihood function is also defined in continuous space and time can easily produce forecasts over discrete space–time intervals.

Description: In order to know whether principal stress orientations in the source area rotated after the 2011 April 11 M w 6.6 Fukushima-Hamadori earthquake in NE Japan, we investigated detailed spatial distributions of stress orientations for both the pre- and post-main shock periods using a large amount of focal mechanism data. We applied stress tensor inversions to focal mechanism data from Japan's National Research Institute for Earth Science and Disaster Prevention's F-net broadband seismic network and the Japan Meteorological Agency (JMA). The 3 -axes estimated for the pre-main shock period are predominantly oriented WSW–ENE, and are relatively homogeneously in space. In contrast, the orientations of the 3 -axes show a significantly heterogeneous distribution in space for the post-main shock period. In the northern subarea of the focal region, the 3 -axes are oriented NW–SE. In the east and west portions of the central subarea, they are oriented NNW–SSE and WNW–ESE, respectively, almost perpendicular to each other. In the southern subarea, the 3 -axes are oriented WSW–ENE. On the whole, the 3 -axis orientations show concentric circle-like distribution surrounding the large slip area of the M w M w 6.6 main shock rupture. The change of principal stress axis orientations after the earthquake is not significant because of the sparse data set for the pre-main shock period. We calculated static stress changes from the M w 6.6 main shock and three M w 〉 5.5 earthquakes to compare with the observed stress axis orientations in the post-main shock period. The 3 -axis orientations of the calculated total static stress change show a concentric circle-like distribution surrounding the large slip area of the main shock, similar to that noted above. This observation strongly suggests that the spatially heterogeneous stress orientations in the post-main shock period were caused by the static stress change from the M w 6.6 main shock and other large earthquakes. In order to estimate the differential stress magnitude in the focal area, we calculated deviatoric stress tensors in the post-main shock period by assuming that they are the sum of the deviatoric stress tensors in the pre-main shock period and the static stress changes. Comparison of the calculated and observed stress tensors revealed differential stress magnitudes of 2–30 MPa that explain the observed stress orientations, considering the probable range of estimated stress ratios in the pre-main shock period.

Description: Deep focus earthquakes within the underthrust Indian lower crust beneath the Himalaya occur in very specific regions and have distinct source characteristics. The study of the source mechanisms of these earthquakes provides valuable constraints on the kinematics of deformation of the underthrust Indian Plate, and its influence on the active deformation of the overlying Himalayan wedge. One of the most significant regions of these deep focus earthquakes is beneath the Sikkim and Bhutan Himalaya. We study the source characteristics of the 2011 September 18 ( M w 6.9) deep focus Sikkim main shock and its major aftershocks using global, regional and local waveform data. We determined the focal mechanism of the main shock using moment tensor inversion of global P and SH waveforms, and ascertained the earthquake fault plane using rupture directivity from regional P -wave spectra. The main shock originated at 53 ± 4 km depth and ruptured at least 20 km thickness of the underthrust Indian lower crust. Faulting occurred on a near vertical dextral strike-slip fault oriented NW-SE (strike 127°, dip 81° and rake 167°), oblique to the local strike of the Himalayan arc. The rupture initiated from the SE end of the fault and propagated to the northwest. The main shock was followed by 20 small-to-moderate aftershocks ( m b  〉 3.0), which we relocated using phase arrival times. We computed the focal mechanisms of the larger ones ( m b  ≥ 3.5) using local waveform inversion. We find that all aftershocks originated SE of the main shock, between depths of 12 and 50 km, and have dominantly strike-slip mechanisms. Our results, combined with the source mechanisms of earthquakes from previous studies, reveals that the entire underthrust Indian crust is seismogenic and deforms by dextral strike-slip motion on oblique structures beneath the Sikkim and Bhutan Himalaya. These active oblique structures with transverse motion possibly mark the western boundary of the clock-wise rotating ‘microplates’ in northeast India observed from GPS geodesy.

Description: We investigate the validity of a proposed generalized Omori–Utsu law for the aftershock sequences for the Landers, Hector Mine, Northridge and Superstition Hills earthquakes, the four largest events in the southern California catalogue we analyse. This law unifies three of the most prominent empirical laws of statistical seismology—the Gutenberg–Richter law, the Omori–Utsu law, and a generalized version of Båth's law—in a formula casting the parameters in the Omori–Utsu law as a function of the lower magnitude cutoff m c for the aftershocks considered. By applying a recently established general procedure for identifying aftershocks, we confirm that the generalized Omori–Utsu law provides a good approximation for the observed rates overall. In particular, we provide convincing evidence that the characteristic time c is not constant but a genuine function of m c , which cannot be attributed to short-term aftershock incompleteness. However, the estimation of the specific parameters is somewhat sensitive to the aftershock selection method used. This includes c ( m c ), which has important implications for inferring the underlying stress field.

Description: Mantle plumes as well as ‘superplumes’ have been imaged in the lowermost mantle in tomographic studies. To investigate seismic resolution of deep mantle plume anomalies, we use a spectral element method (SEM) to simulate global seismic wave propagation in 3-D wave speed models and measure frequency-dependent P -, S -, Pdiff - and Sdiff -wave traveltime anomalies caused by plume structures in the lowermost mantle. We compare SEM time delay measurements with calculations based on ray theory and show that an anticorrelation between bulk sound wave speed and S -wave speed could be produced as an artifact. This is caused by different wavefront healing effects between P and S waves in thermal plume models. The differences in wave diffraction between the two types of waves depend on epicentral distance and wave frequency. We show that bulk-sound speed structure can not be recovered in ray-theoretical tomographic inversions when the lateral extent of the anomaly is smaller than the size of the Fresnel zone in the lowermost mantle. In addition, an anticorrelation between bulk sound speed and S -wave speed can be produced in ray-theoretical tomography when the size of the anomaly is less than ~2000 km; and, the artifacts become more pronounced as the lateral extent of the plume decreases. This indicates a chemical origin of ‘superplumes’ in the lowermost mantle may not be necessary to explain observed seismic traveltimes of core–mantle diffracted waves. The same set of Pdiff and Sdiff measurements are inverted using finite-frequency tomography based on Born sensitivity kernels. We show that wavefront healing effects can be accounted for in finite-frequency tomography to recover the true velocity model.

Description: Monitoring and understanding induced seismicity is critical in order to estimate and mitigate seismic risk related to numerous existing and emerging techniques for natural resource exploitation in the shallow-crust. State of the art approaches for guiding decision making, such as traffic light systems, rely heavily on data such as earthquake location and magnitude that are provided to them. In this context we document the monitoring of a deep geothermal energy project in St Gallen, Switzerland. We focus on the issues of earthquake magnitude, ground motion and macroseismic intensity which are important components of the seismic hazard associated to the project. We highlight the problems with attenuation corrections for magnitude estimation and site amplification that were observed when trying to apply practices used for monitoring regional seismicity to a small-scale monitoring network. Relying on the almost constant source-station distance for events in the geothermal ‘seismic cloud’ we developed a simple procedure, calibrated using several M L 〉 1.3 events, which allowed the unbiased calculation of M L using only stations of the local monitoring network. The approach determines station specific M L correction terms that account for both the bias of the attenuation correction in the near field and amplification at the site. Since the smallest events ( M L 〈 –1) were only observed on a single borehole instrument, a simple relation between the amplitude at the central borehole station of the monitoring network and M L was found. When compared against magnitudes computed over the whole network this single station approach was shown to provide robust estimates (±0.17 units) for the events down to M L = –1. The relation could then be used to estimate the magnitude of even smaller events ( M L 〈 –1) only recorded on the central borehole station. Using data from almost 2700 events in Switzerland, we then recalibrated the attenuation correction, extending its range of validity from a minimum source–station distance of 20 km down to 1 km. Based on this we could determine the component of the previously derived station specific M L corrections due to local amplification. We analysed ground-motion and detailed macroseismic reports resulting from the 2013 July 20 St Gallen M L = 3.5 ± 0.1 ( M w = 3.3–3.5 ± 0.1) ‘main shock’ and compared it to a similar M L = 3.4 ± 0.1 event ( M w = 3.2 ± 0.1) that occurred in 2006 at another deep geothermal project in Basel, Switzerland. Differences in ground motion amplitudes between the Basel and St Gallen events and to an extent, the associated macroseismic observations, were investigated in terms of the different source terms: M w for long-period motions and the source-corner frequency (related to the source rupture velocity and stress-drop) for short periods.

Description: Radial and azimuthal anisotropy in seismic wave speeds have long been observed using surface waves and are believed to be controlled by deformation within the Earth's crust and uppermost mantle. Although radial and azimuthal anisotropy reflect important aspects of anisotropic media, few studies have tried to interpret them jointly. We describe a method of inversion that interprets simultaneous observations of radial and azimuthal anisotropy under the assumption of a hexagonally symmetric elastic tensor with a tilted symmetry axis defined by dip and strike angles. We show that observations of radial anisotropy and the 2 component of azimuthal anisotropy for Rayleigh waves obtained using USArray data in the western United States can be fit well under this assumption. Our inferences occur within the framework of a Bayesian Monte Carlo inversion, which yields a posterior distribution that reflects both variances of and covariances between all model variables, and divide into theoretical and observational results. Principal theoretical results include the following: (1) There are two distinct groups of models (Group 1, Group 2) in the posterior distribution in which the strike angle of anisotropy in the crust (defined by the intersection of the foliation plane with Earth's surface) is approximately orthogonal between the two sets. (2) The Rayleigh wave fast axis directions are orthogonal to the strike angle in the geologically preferred group of models in which anisotropy is strongly non-elliptical. (3) The estimated dip angle may be interpreted in two ways: as a measure of the actual dip of the foliation of anisotropic material within the crust, or as a proxy for another non-geometric variable, most likely a measure of the deviation from hexagonal symmetry of the medium. The principal observational results include the following: (1) Inherent S -wave anisotropy ( ) is fairly homogeneous vertically across the crust, on average, and spatially across the western United States. (2) Averaging over the region of study and in depth, in the crust is approximately 4.1 ± 2 per cent. in the crust is approximately the same in the two groups of models. (3) Dip angles in the two groups of models show similar spatial variability and display geological coherence. (4) Tilting the symmetry axis of an anisotropic medium produces apparent radial and apparent azimuthal anisotropies that are both smaller in amplitude than the inherent anisotropy of the medium, which means that most previous studies have probably underestimated the strength of anisotropy.

Description: Scattering and refractions that occur in the heterogenous near-surface beneath seismic stations can strongly affect the relative amplitudes recorded by three-component seismometers. Using data from Tungurahua volcano we have developed a procedure to correct these ‘ site effects’. We show that seismic noise signals store site information, and then use their normalized spectral amplitudes as site frequency response functions. The process does not require a reference station (as per the S -wave and coda methods) or assume that the vertical amplitude is constant (the H/V component ratio method). Correcting the site effects has three consequences on data analysis: (1) improvement of the seismic source location and its energy estimation; (2) identification of a strong influence on the volcanic acoustic seismic ratio (VASR) and (3) decoupling the air wave impact on the ground caused by explosions or eruption jets. We show how site effect corrections improve the analysis of an eruption jet on 2006 July 14–15, appearing two periods of strong acoustic energy release and a progressive increase of the seismic energy, reaching the maximum before finishing the eruption.

Description: We propose a method to invert surface wave dispersion data directly for 3-D variations of shear wave speed, that is, without the intermediate step of phase or group velocity maps, using frequency-dependent ray tracing and a wavelet-based sparsity-constrained tomographic inversion. A fast marching method is used to compute, at each period, surface wave traveltimes and ray paths between sources and receivers. This avoids the assumption of great-circle propagation that is used in most surface wave tomographic studies, but which is not appropriate in complex media. To simplify the problem we consider quasi-stratified media with smoothly varying seismic properties. We represent the 3-D shear wave speed model by means of 1-D profiles beneath grid points, which are determined from all dispersion data simultaneously using a wavelet-based sparsity-constrained tomographic method. The wavelet coefficients of the wave speed model are estimated with an iteratively reweighted least squares algorithm, and upon iteration the surface wave ray paths and the data sensitivity matrix are updated using the newly obtained wave speed model. To demonstrate its feasibility, we apply the method to determine the 3-D shallow crustal shear wave speed variations in the Taipei basin of Taiwan using short period interstation Rayleigh wave phase velocity dispersion measurements extracted from the ambient noise cross-correlation method. The results are consistent with previous studies and reveal strong shallow crustal heterogeneity that correlates with surface geology.

Description: In order to know whether principal stress orientations in the source area rotated after the 2011 April 11 M w 6.6 Fukushima-Hamadori earthquake in NE Japan, we investigated detailed spatial distributions of stress orientations for both the pre- and post-main shock periods using a large amount of focal mechanism data. We applied stress tensor inversions to focal mechanism data from Japan's National Research Institute for Earth Science and Disaster Prevention's F-net broadband seismic network and the Japan Meteorological Agency (JMA). The 3 -axes estimated for the pre-main shock period are predominantly oriented WSW–ENE, and are relatively homogeneously in space. In contrast, the orientations of the 3 -axes show a significantly heterogeneous distribution in space for the post-main shock period. In the northern subarea of the focal region, the 3 -axes are oriented NW–SE. In the east and west portions of the central subarea, they are oriented NNW–SSE and WNW–ESE, respectively, almost perpendicular to each other. In the southern subarea, the 3 -axes are oriented WSW–ENE. On the whole, the 3 -axis orientations show concentric circle-like distribution surrounding the large slip area of the M w M w 6.6 main shock rupture. The change of principal stress axis orientations after the earthquake is not significant because of the sparse data set for the pre-main shock period. We calculated static stress changes from the M w 6.6 main shock and three M w 〉 5.5 earthquakes to compare with the observed stress axis orientations in the post-main shock period. The 3 -axis orientations of the calculated total static stress change show a concentric circle-like distribution surrounding the large slip area of the main shock, similar to that noted above. This observation strongly suggests that the spatially heterogeneous stress orientations in the post-main shock period were caused by the static stress change from the M w 6.6 main shock and other large earthquakes. In order to estimate the differential stress magnitude in the focal area, we calculated deviatoric stress tensors in the post-main shock period by assuming that they are the sum of the deviatoric stress tensors in the pre-main shock period and the static stress changes. Comparison of the calculated and observed stress tensors revealed differential stress magnitudes of 2–30 MPa that explain the observed stress orientations, considering the probable range of estimated stress ratios in the pre-main shock period.

Description: A new effective model is presented to compute horizontal-to-vertical spectral ratios (HVSR) relative to ambient vibrations, under the assumption that these are originated by a distribution of spatially correlated random surface sources. The major novelty of this model lies in the description of both ground displacement and sources as stochastic fields defined on the Earth's surface, stationary in time and homogeneous in space. In this frame, the power spectral density of the displacement stochastic field can be written as a function of the power spectral density of the force stochastic field and of the subsoil properties, through the relevant Green's function. Spatial correlation between ambient vibration sources is shown to be a necessary condition to warrant convergence of the integrals defining the frequency power spectra of the displacement field that make up the HVSR curve. Furthermore, it is shown that this HVSR curve may be significantly affected by the effective range of the force-field correlation on the Earth's surface. This formalization reduces computational efforts with respect to the previous version of the model based on distributed surface sources and may provide synthetic HVSR-curve patterns that are in line with those given by that computationally more troublesome version, as well as with those deduced under the assumption that the ambient vibrations constitute a diffuse wavefield.

Description: We present an extension of the nodal discontinuous Galerkin method for elastic wave propagation to high interpolation orders and arbitrary heterogeneous media. The high-order lagrangian interpolation is based on a set of nodes with excellent interpolation properties in the standard triangular element. In order to take into account highly variable geological media, another set of suitable quadrature points is used where the physical and mechanical properties of the medium are defined. We implement the methodology in a 2-D discontinuous Galerkin solver. First, a convergence study confirms the hp -convergence of the method in a smoothly varying elastic medium. Then, we show the advantages of the present methodology, compared to the classical one with constant properties within the elements, in terms of the complexity of the mesh generation process by analysing the seismic amplification of a soft layer over an elastic half-space. Finally, to verify the proposed methodology in a more complex and realistic configuration, we compare the simulation results with the ones obtained by the spectral element method for a sedimentary basin with a realistic gradient velocity profile. Satisfactory results are obtained even for the case where the computational mesh does not honour the strong impedance contrast between the basin bottom and the bedrock.

Description: The 2003 Bam, Iran, earthquake ( M w  = 6.6) was recorded by the BAM accelerometer station. Since the causative fault was located just below the city, the accelerometer recorded the main shock, a foreshock and several local aftershocks. To study the scenario of rupturing, we simulated all components of the observed main shock waveform via the empirical Green's function method. 28 selected aftershocks and the single foreshock are used to simulate the main shock in the frequency range of 0.5–5 Hz. Since the events were very close to the station, some small events may not have similar path effects to the main shock. Therefore, it is essential to employ some appropriate changes to the waveforms to alleviate path difference effects. The starting point of the rupture is identified in the centre of the strong motion generation area and is located approximately 5 km south of the BAM station and in depth of about 7 km. The horizontal simulated components imply that the main shock was located west of the BAM station. In contrast, significant variation in the ratio of amplitudes in EW and NS components may be used to discuss the possibility of dissimilarity in the focal mechanism of the small events. Most aftershocks with similar mechanisms to the main shock, that is similar EW/NS maximum amplitude ratio, have capacity to simulate certain peaks of their horizontal components. However, some small events with different mechanisms are only able to simulate the peaks of up to one horizontal component. Some changes were applied to the empirical Green's function method to incorporate two small events by using a combined fault model. While the two aftershocks have different mechanisms, some combinations may improve simulations. The rupture initiating point at the middle of the fault plane and improved simulations by combination of two fault surfaces with different focal mechanisms may suggest a bilateral rupture and combination of two focal mechanisms for the main shock of the Bam earthquake.