# Impact of 3-D Earth structure on W-phase CMT parameters

Catalina Morales-Yáñez, Zacharie Duputel, Luis Rivera

► To cite this version:

Catalina Morales-Yáñez, Zacharie Duputel, Luis Rivera. Impact of 3-D Earth structure on W-phase CMT parameters. Geophysical Journal International, 2020, 223 (2), pp.1432-1445. 10.1093/gji/ggaa377. hal-03032581

HAL Id: hal-03032581

https://hal.science/hal-03032581

Submitted on 28 Aug 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Impact of 3-D Earth structure on W-phase CMT parameters

Catalina Morales-Yáñez, Zacharie Duputel and Luis Rivera

Institut de Physique du Globe de Strasbourg, UMR7516, Université de Strasbourg/EOST, CNRS, 67081 Strasbourg, France. E-mail: cmorales@unistra.fr

Accepted 2020 August 6. Received 2020 August 5; in original form 2020 January 31

## SUMMARY

We investigate the impact of unmodelled 3-D structural heterogeneity on inverted W-phase source parameters. We generate a large data set of synthetic seismograms accounting for the Earths 3-D structure for 250 earthquakes globally distributed. The W-phase algorithm is then used to invert for earthquake CMT parameters, assuming a spherical Earth model. The impact of lateral heterogeneity is assessed by comparing inverted source parameters with those used to compute the 3-D synthetics. Results show that the 3-D structure mainly affects centroid location while the effect on the other source parameters remains small. Centroid mislocations present clear geographical patterns. In particular, W-phase solutions for earthquakes in South America are on average biased 17 km to the east of the actual centroid locations. This effect is significantly reduced using an azimuthally well balanced distribution of seismological stations. Source parameters are generally more impacted by mantle heterogeneity while the scalar moment of shallow earthquakes seems to be mainly impacted by the crustal structure. Shallow earthquakes present a variability of $M_{r\theta}$ and $M_{r\phi}$ moment tensor elements, resulting both from the small amplitude and a larger uncertainty of the associated Green's functions.

**Key words:** Structure of the Earth; Inverse theory; Earthquake source observations; Surface waves and free oscillations; Wave propagation.

## 1 INTRODUCTION

The characterization of earthquake sources is based on different data types such as teleseismic body waves, surface waves, strong-motion waveforms, GNSS and InSAR data (e.g. Delouis et al. 2010; Polet & Thio 2011; Crowell et al. 2012; Ekström et al. 2012). Source inversions also rely on various representations such as centroid–moment–tensors (CMT; e.g. Dziewonski et al. 1981), multiple point sources (e.g. Tsai et al. 2005), linear and non-linear finite-fault parametrizations (e.g. Ide 2007). The resulting source models are affected by different types of uncertainties (e.g. Duputel et al. 2012b). A first source of uncertainty is the error induced by more or less imperfect measurements (e.g. ambient noise, incorrect instrument calibration, etc.). Another source of uncertainty, often overlooked, is associated with modeling errors (i.e. errors in model predictions) and in particular, prediction inaccuracies due to imperfections of the Earth model used for the inversion. While different strategies have been developed to quantify such uncertainties (Yagi & Fukahata 2011; Duputel et al. 2014; Hallo & Gallovič 2016), most applications rely on a simple 1-D (or spherical) Earth model (Hallo et al. 2017; Gombert et al. 2018). It is the case for W-phase CMT inversions for which 3-D effects are supposedly small but have not yet been fully quantified (Kanamori & Rivera 2008; Duputel et al. 2016). The W-phase corresponds to a long period signal (100–1000 s) that is conspicuous for large earthquakes between the P wave and the surface waves. This phase can be described as the superposition

of normal mode overtones, that have limited sensitivity to shallow lateral heterogeneity compared to fundamental mode surface waves (Kanamori 1993).

Lateral structural heterogeneity can significantly affect the propagation of seismic waves, which in return can impact source estimates if they are not properly accounted for. Different corrections can be applied to mitigate the effect of 3-D Earth structures. For example, Nakanishi & Kanamori (1982) proposed to use path-dependent Rayleigh wave phase velocities for moment tensor inversions. Another example is the Global CMT algorithm in which S362ANI is used as the 3-D model. The seismograms are calculated using 3-D structure and the path-average approximation (Dziewonski et al. 1984; Woodhouse & Dziewonski 1984). Despite these corrections, CMT solutions can still be affected by errors due to lateral heterogeneity (Ferreira & Woodhouse 2006). In particular, the comparison between interferometric synthetic aperture radar (InSAR) data and CMT solutions suggests large uncertainties in centroid locations (Ferreira et al. 2011). In the same line, Hjörleifsdóttir & Ekström (2010) found that even with phase corrections, Global CMT solutions are still contaminated by the Earth’s 3-D structure, causing an error close to 11 per cent in moment along with non-negligible regional biases in centroid locations (e.g. larger than 10–20 km in South America).

With the recent improvement of computing capabilities, several studies are now using more expensive numerical methods to better incorporate the effect of 3-D heterogeneity (Liu et al. 2004;
Duputel et al. 2016). This is particularly relevant to account for 3-D effects when using fundamental mode surface waves that are strongly affected by lateral heterogeneity near the surface (e.g. caused by oceans and continents; Dahlen & Tromp 1998). Some other phases, like the W-phase, are less impacted by shallow 3-D structures because they propagate deeper in the mantle where lateral heterogeneity is supposedly smaller. The inversion of W-phase proved to be very useful for rapid source characterization of large earthquakes (Duputel et al. 2012b) and has been implemented in various warning centres to quickly provide CMT solutions in near real time conditions (Hayes et al. 2009; Wang et al. 2017; Zhao et al. 2017; Riquelme et al. 2018). The robustness of this approach relies partially on the fact that the W-phase corresponds mainly to the superposition of normal-mode overtones at long period, which are not strongly influenced by the 3-D structure and can be efficiently synthesized in a spherical Earth model (Kanamori & Rivera 2008). Even if W-phase should be less affected by the 3-D structure, recent studies suggest that there might still be some impact (Duputel et al. 2016). However, such 3-D effects have not yet been quantified systematically. In this study, we assess the impact of 3-D Earth structures on W-phase CMT solutions. For this purpose, we evaluate the performance of the W-phase algorithm using a large data set of synthetic seismograms (~64 000 waveforms) computed for 3-D Earth models.

## 2 METHODOLOGY

To assess the effects of lateral structural heterogeneity on W-phase solutions, we compute a large database of 3-D synthetics for earthquakes at various locations and with different source parameters. From those synthetic seismograms, we then invert for CMT parameters using the (1-D) W-phase approach. The performance of the algorithm is then assessed by comparing inverted source parameters with those used to compute input 3-D synthetics.

### 2.1 Earthquake catalogue and 3-D synthetic database

We first define a set of earthquakes to be used in this study. To have a reasonably representative catalogue, we use a 5° × 5° grid at the surface of the Earth and extract CMT parameters of the largest events ($M_W \ge 6.5$) in each cell from the Global CMT (GCMT) catalogue between 1995 and 2017 (Dziewonski et al. 1981; Ekström et al. 2012). This results into a catalogue of 252 earthquakes whose distribution and focal mechanisms are shown in Fig. 1. We also design a realistic network composed of 254 seismological stations. These stations, presented in Fig. 2, mainly belong to networks IU, II, GT, IC, CU, G, MN, CN, GE, CI and BK (for more information see caption Fig. 2).

Using the earthquake catalogue and seismic network described above, we compute synthetic seismograms assuming two different 3-D global mantle models: S362ANI (Kustowski et al. 2008) and S40RTS (Ritsema et al. 2011), which are both completed with the CRUST2.0 crust structure (Bassin et al. 2000). We use the spectral element method code SPECFEM3D_GLOBE (Komatitsch et al. 2015) with a mesh-size ensuring accurate simulations for periods down to ~27 s. This resolution is compatible with the filter used for smaller magnitude earthquakes ($M_W = 6.5$) and allows to incorporate realistic features of the Earth. Simulations are conducted for each earthquake in our catalogue. Computing the complete catalogue takes ~40 hr using 48 Tesla K80 Nvidia GPUs on the University of Strasbourg HPC cluster. Using 24 GPUs, each earthquake

simulation is conducted in 20 min to obtain seismograms with a duration of 40 min. Beyond 3-D lateral heterogeneity, these computations include the effect of ellipticity, topography, rotation, attenuation and gravity (implemented using the Cowling approximation; Komatitsch & Tromp 2002a, b).

The synthetics are computed such that they correspond to a unit scalar moment and a step source time history (i.e. a Dirac delta moment rate function). To assess the effect of ambient noise for different magnitudes, we add actual seismic noise after scaling the 3-D synthetics to different earthquake sizes. We use noise records for each used channel on time-periods with no significant earthquake activity (on 1995/06/01, 2005/11/3 and 2015/06/09, depending on the seismic station). The final noisy synthetics $s(t, M_0)$ at a given station can then be written as:

$$s(t, M_0) = M_0 \times \hat{s}(t) \otimes f(t, M_0) \otimes I(t) + n(t), \quad (1)$$

where the symbol $\otimes$ is used to denote convolution. In this equation, $M_0$ is the scalar seismic moment, $\hat{s}(t)$ is the unit synthetic seismogram (i.e. computed for unit scalar moment and step time-history), $f(t, M_0)$ is the normalized moment rate function, $I(t)$ is the instrument response and $n(t)$ the raw noise record at the corresponding station. The source time function $f(t, M_0)$ is assumed to be an isosceles triangle function of unit area with a time-shift ($t_c$) and a half duration ($h_c$) defined as (Duputel et al. 2012b):

$$t_c(M_0) = h_c(M_0) = 2.6 \times 10^{-6} M_0^{1/3}. \quad (2)$$

In the above equations, $M_0$ is in N.m (i.e. $10^7$ dyne.cm) and $t_c$ is in seconds.

### 2.2 W-phase inversion

The W-phase inversion algorithm is based on a Green's functions database that is pre-computed for a 1-D Earth model. For consistency, we use the same code to compute 3-D synthetics and to create the 1-D Green's functions. We compute SPECFEM3D_GLOBE Green's functions using the 1-D Earth model STW105 (Kustowski et al. 2008), which is the reference model of S362ANI.

In general, for every source-station pair we have to compute 18 Green's functions (6 moment tensor components, 3 receiver orientations). However, for a spherical Earth model, only 10 Green's functions are required for every depth-distance pair (see Kanamori & Rivera 2008). We thus compute a database including those 10 elementary Green's functions for a range of epicentral distances and focal depths considered in the problem. Epicentral distances are discretized every 0.1° from 0° to 90°. The depth discretization is variable with 2, 5 and 10 km depth intervals, respectively, for depth ranges of 3.5–25.5 km, 25.5–50.5 km and 50–760 km.

The W-phase inversion algorithm consists of estimating the centroid moment tensor parameters (i.e. the moment tensor elements along with the centroid location in time and space). This relies on a grid-search approach to find the point source time and location that minimizes the RMS waveform misfit. For each explored point in time and space the inverse problem is linear for the moment tensor elements (Kanamori & Rivera 2008). As for standard W-phase implementations, we use stations within 5°–85° of epicentral distance and a time-window starting at the P-wave arrival, with a duration ($\delta t$) proportional to the epicentral distance ($\Delta$): $\delta t = 15[s/°]\Delta$. We also apply a data screening by iteratively removing channels having a relative RMS misfit larger than 500, 300 and 90 per cent of the data L2 norm (for more details, see Duputel et al. 2012b).
<table><caption>**Table 1.** Comparison of the impact of crustal and mantle lateral heterogeneity. We compare solutions obtained using a full 3-D structure with solutions retrieved from a 3-D crust with a 1-D mantle (3Dc–1Dm) or a 1-D crust with a 3-D mantle (1Dc–3Dm). $\langle\Delta M_{W3D} - \Delta M_{Wmodel}\rangle$ is the average difference of the magnitude discrepancy $\Delta M_{W3D}$ (measured for a 3-D model) and $\Delta M_{Wmodel}$ (where model is 3Dc–1Dm or 1Dc–3Dm). Similarly, $\langle\Delta_{3D} - \Delta_{model}\rangle$, $\langle\Delta\tau_{3D} - \Delta\tau_{model}\rangle$, $\langle\Delta h_{3D} - \Delta h_{model}\rangle$ and $\langle|\Delta\mathbf{x}_{3D}| - |\Delta\mathbf{x}_{model}|\rangle$ denote, respectively, the average beachball RMS difference, centroid time difference, depth difference and horizontal mislocations.</caption><thead><tr><th></th><th>3-D crust 1-D mantle</th><th>1-D crust 3-D mantle</th></tr></thead><tbody><tr><td>⟨ΔM<sub>W3D</sub> - ΔM<sub>Wmodel</sub>⟩</td><td>-0.00014</td><td>0.0041</td></tr><tr><td>⟨Δ<sub>3D</sub> - Δ<sub>model</sub>⟩</td><td>0.0068</td><td>0.0052</td></tr><tr><td>⟨Δτ<sub>3D</sub> - Δτ<sub>model</sub>⟩</td><td>0.012 s</td><td>-0.004 s</td></tr><tr><td>⟨Δh<sub>3D</sub> - Δh<sub>model</sub>⟩</td><td>2.69 km</td><td>0.44 km</td></tr><tr><td>⟨|Δx<sub>3D</sub>| - |Δx<sub>model</sub>|⟩</td><td>6.50 km</td><td>-3.34 km</td></tr></tbody></table>

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Earthquake catalogue. We use centroid moment tensor solutions from the GCMT database. Blue focal mechanisms represent earthquakes shallower than 50 km and red represents earthquakes deeper than 50 km.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Seismological network. We use a combination of 254 stations from the IRIS/IDA (II: Scripps Institution of Oceanography 1986), ASL/USGS (IU; CU; IC; GT; Albuquerque Seismological Laboratory (ASL)/USGS 1988, 2006, 1992, 1993), BDSN (BK; Northern California Earthquake Data Center 2014), CNSN (CN; Geological Survey of Canada 1980), GEOFON (GE; GEOFON Data Centre 1993), GEOSCOPE (G; Institut de Physique du Globe de Paris IPGP 1982), MedNet (MN; MedNet Project Partner Institutions 1990) and Southern California Seismic Network (CI; California Institute of Technology and United States Geological Survey Pasadena (SCSN) 1926).</figcaption></figure>
# 3 IMPACT OF 3-D STRUCTURE ON W-PHASE CMT PARAMETERS

To evaluate the impact of lateral heterogeneity, we compare the source parameters used to compute 3-D synthetic seismograms (input source parameters) with the retrieved W-phase CMT solutions (output source parameters).

## 3.1 Comparison between input and retrieved source parameters

We first focus on results obtained for $M_W = 7.5$ earthquakes. Before inversion, traces are filtered in the 150–500 s passband as proposed by Duputel et al. (2012b) for this magnitude. The comparison between input and retrieved source parameters is summarized in Fig. 3 for Earth models S362ANI and S40RTS. In order to have a complete view of the influence of 3-D structures on W-phase CMT parameters, we evaluate the performance of the algorithm using five different metrics.

First, we evaluate the magnitude difference, defined as:

$$ \Delta M_W = M_W^{\text{output}} - M_W^{\text{input}} = \frac{2}{3} \log_{10} \frac{M_0^{\text{output}}}{M_0^{\text{input}}}, \quad (3) $$

where $M_W^{\text{input}}$ and $M_W^{\text{output}}$ are, respectively, the input and retrieved moment magnitude. $M_0^{\text{input}}$ and $M_0^{\text{output}}$ are the corresponding scalar moments, calculated using $M_0 = \sqrt{\frac{M \cdot M}{2}}$ with M the moment tensor (i.e. following the definition of the total scalar moment by Silver & Jordan 1982). Figs 3(a) and (f) show that magnitude difference is tightly centred around 0, with 83 per cent of earthquakes for which $|\Delta M_W| < 0.03$. Largest positive magnitude differences are observed in Greece, the Caspian Sea and Nepal, while the largest negative values are obtained for earthquakes in Indonesia and Northern Japan.

To measure the impact of 3-D heterogeneity on focal mechanisms, we estimate the beachball RMS difference ($\Delta$), following Rivera & Kanamori (2014):

$$ \Delta = \frac{1}{2\sqrt{2}} (\mathbf{D} : \mathbf{D})^{1/2} \quad (4) $$

$\mathbf{D}$ is the difference between normalized moments tensor $\mathbf{D} = \hat{\mathbf{M}}^{\text{output}} - \hat{\mathbf{M}}^{\text{input}}$, with $\hat{\mathbf{M}} = \frac{\mathbf{M}}{M_0}$. With this definition, we have $\Delta = 0$ when the two moment tensors are equal and $\Delta = 1$ when they are opposite. Overall, beachball RMS difference (Figs 3b and g) are small with ~90 per cent of the events showing $\Delta < 0.1$.

We also evaluate the difference between centroid times,

$$ \Delta \tau = t_s^{\text{output}} - t_s^{\text{input}} \quad (5) $$

In Figs 3(c) and (h), we note a clear geographical consistency in $\Delta \tau$. In the Northeast Pacific Ocean, the output centroid time is generally later than the input, while in the Northwest Pacific and South Atlantic we see the opposite behavior. The standard deviation in $\Delta \tau$ is ~2 s for the entire data set, which correspond to ~6 per cent of the input rupture duration for a $M_W = 7.5$ event (cf. eq. 2). This level of uncertainty is acceptable given that the used data set is sampled at 1 sample per second.

Figs 3(d) and (i) show the centroid depth difference defined as

$$ \Delta h = h^{\text{output}} - h^{\text{input}} \quad (6) $$

The retrieved centroid depths are overall deeper than the input ones. For S362ANI the average depth bias is +3 km while it is 1 km for S40RTS. This bias remains relatively small given that the depth discretization of the W-phase Green’s function ranges from 2 to

10 km. The absolute depth difference is smaller than 10 km for 85 and 90 per cent of the events calculated with the S362ANI and S40RTS models, respectively.

Finally, we evaluate the horizontal centroid mislocations as

$$ \Delta \mathbf{x} = \mathbf{x}^{\text{output}} - \mathbf{x}^{\text{input}}, \quad (7) $$

where $\mathbf{x}^{\text{input}}$ and $\mathbf{x}^{\text{output}}$ are input and retrieved horizontal centroid locations. We estimate an average centroid mislocation $|\Delta \mathbf{x}| \sim 17$ km. Although such uncertainty is not negligible compared to the typical rupture length of $M_W = 7.5$ earthquakes, it is of the same order as the 0.1° source-receiver distance interval used in our Green’s function database. Results in Figs 3(e) and (j) also indicate that there are significant geographical biases on retrieved centroid locations. In particular, earthquakes in central and south America seem to be biased towards the East while they seem to be biased towards the North in Alaska, the Aleutians and Japan.

If we compare Figs 3(a)–(e) and (f)–(j), we note that the impact of lateral heterogeneity on W-phase solutions is globally consistent between both models S362ANI and S40RTS. The S40RTS model has a higher resolution (up to 40 harmonic degree) than the S362ANI model (18 harmonic degrees). The consistency between both models confirms that the W-phase is not sensitive to small details of the structure, yet it is sensitive to large-scale heterogeneity that is present in both models.

## 3.2 Long-period ambient seismic noise

To assess the influence of ambient seismic noise on the inversion results, we compare in Fig. 4 solutions obtained for $M_W = 7.5$ earthquakes with solutions obtained for $M_W = 6.5$ events for the Earth model S362ANI with Crust2.0. We recall that actual seismic noise has been added to synthetic seismograms (see eq. 1) such that the long-period signal-to-noise ratio is naturally larger for a $M_W = 7.5$ earthquake than for a $M_W = 6.5$ event. To avoid any bias due to the use of different stations, the stations set is fixed prior to inversion.

Figs 4(a) and (b) show that noise-free moment tensor solutions obtained for $M_W = 6.5$ are very similar to those obtained for $M_W = 7.5$. On the other hand, the solutions are significantly affected when ambient noise is added to the synthetics (cf. Figs 4d and e). More specifically, we observe that the dispersion of beachball RMS difference ($\Delta$) for $M_W = 6.5$ is more than two times larger than for $M_W = 7.5$ earthquakes. This difference results from long-period noise that is particularly prominent when analysing smaller magnitude earthquakes.

In order to limit the impact of ambient noise for moderate sized earthquakes, it has been previously proposed to use a magnitude dependent bandpass filtering (Hayes et al. 2009; Duputel et al. 2012b). Using the 100–250 s passband suggested by Duputel et al. (2012b) for $M_W = 6.5$ earthquakes (instead of 150–500 s for $M_W = 7.5$ events), we see in Fig. 4 f that the resulting solutions are significantly less affected by ambient seismic noise. A detailed comparison of input and inverted solutions for $M_W = 6.5$ is presented in Fig. S1 for Earth models S362ANI and S40RTS using the 100–250 s passband. Overall, we note a slightly larger dispersion for $M_W = 6.5$ than for $M_W = 7.5$ earthquakes (cf. Fig. 4). This likely results from the smaller signal-to-noise ratio of W-phase waveforms for moderate sized events and larger 3-D effects when using shorter period waveforms.
Earth model S362ANI

<figure><img src="image_5.png" /><figcaption>(b) Beachball RMS difference</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(c) Difference in centroid time</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(d) Depth difference</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(e) Horizontal mislocation</figcaption></figure>

<figure><img src="image_9.png" /></figure>

Earth model S40RTS

<figure><img src="image_10.png" /><figcaption>(g) Beachball RMS difference</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>(h) Difference in centroid time</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>(i) Depth difference</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>(j) Horizontal mislocation</figcaption></figure>

<figure><img src="image_14.png" /></figure>

**Figure 3.** Effects of 3-D heterogeneity for *M*<sub>W</sub> = 7.5 earthquakes. Maps show the following quantities from top to bottom: moment magnitude difference, beachball RMS difference, difference in time-shift, depth difference and horizontal mislocation. Panels from (a) to (e) correspond to the inversion of synthetic seismograms created with the Earth model S362ANI. Panels (f) to (j) correspond to the inversion of synthetic seismograms created with the Earth model S40RTS. In both cases, seismic noise has been added to synthetics prior to source inversion (see Section 2.1).

## 3.3 Relative influence of crustal and mantle heterogeneity

To further investigate the influence of the 3-D structure on W-phase solutions, we analyse separately the influence of the crust and the mantle. Using the procedure described in section 2, we compute synthetics seismograms with noise, scaled to *M*<sub>W</sub> = 7.5, with two different Earth models:

(1) A model with a 3-D crust and a 1-D mantle (hereafter noted 3Dc–1Dm), using CRUST2.0 and STW105.

(2) A model with a 1-D crust and a 3-D mantle (hereafter noted 1Dc–3Dm), using STW105 and S362ANI.

As before we compare input and estimated parameters using the five metrics defined in Section 3.1.

The results are presented in Figs 5 and S3. To evaluate the relative impact of crust and mantle heterogeneity, we also show in Table 1 the average discrepancy in magnitude difference (Δ*M*<sub>W</sub>) beachball RMS difference (Δτ), centroid time difference (Δτ), depth difference (Δh) and horizontal mislocations (Δx). If the average difference for the 3Dc–1Dm model is closer to zero than for the
<figure><img src="image_15.png" /><figcaption>**Figure 4.** Beachball RMS difference (Δ) after W-phase inversion in different frequency passband with and without ambient seismic noise. Panels (a) and (d) correspond to $M_W = 7.5$ earthquakes inverted in the 150–500 s passband. Panels (b) and (e) are for the same passband but for $M_W = 6.5$. (c) and (f) are for $M_W = 6.5$ with the 100–250 passband.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 5.** Effects of crustal and mantle heterogeneity on W-phase solutions obtained for $M_W = 7.5$ earthquakes. Maps show from top to bottom: magnitude difference, difference in time-shift and horizontal mislocation. Figures from (a) to (c) correspond to the inversion of synthetic seismograms created with a 3-D crust (CRUST2.0) and a 1-D mantle (STW105). Figures (d) to (f) correspond to the inversion of synthetic seismograms created with a 1-D crust (STW105) and a 3-D mantle (S362ANI). As described in the main text, ambient noise has been added to synthetics prior to source inversion.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 6.** Impact of station distribution on horizontal mislocation. Maps show the horizontal mislocation, after inversion of noisy synthetic seismograms created with the Earth model S362ANI, computed for the real global seismological network (a) and for an azimuthally well balanced array (b).</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>**Figure 7.** Effect of earthquake depth on retrieved source parameters. (a) Magnitude difference, (b) beachball RMS difference and (c) centroid time difference is shown for as a function of centroid depth.</figcaption></figure>

1Dc–3Dm model, we can postulate that the evaluated parameter is more sensitive to crustal heterogeneity than to the mantle, and vise versa. Table 1 suggests that the crust has a larger impact on the magnitude difference than mantle. This is consistent with Figs 5(a) and (d) that shows a larger dispersion of the magnitude difference with the Earth model 3Dc–1Dm than 1Dc–3Dm. Table 1 indicates that the other parameters seem to be mainly impacted by mantle heterogeneity (consistently with Figs 5 and S3). This is particularly true for the centroid location, for which mantle heterogeneity induces mislocations that are much larger than those induced by the crust (see Figs 5c and f). Interestingly, we see in Figs 6(b)–(c) and (e)–(f) that 3-D heterogeneity in the crust and in the mantle results into anti-correlated effects in terms of centroid time and horizontal mislocation (this is further discussed in Section 4.2).

# 4 DISCUSSION

## 4.1 Centroid depth and location uncertainties

As shown in Fig. 3, lateral heterogeneity induces non-negligible uncertainties both in terms of centroid depth and centroid location. We note a centroid depth bias of about +3 km to larger depths (Figs 3d, i and Figs S1d, i). A similar effect was noted by Hjörleifsdóttir & Ekström (2010), who reported a bias of +6 km in GCMT solutions due to 3-D heterogeneity. Similar to GCMT solutions, our results indicate negligible trade-offs between centroid times and centroid depths. While GCMT depths seem to be mainly impacted by local velocity structures at the source and receivers, Figs S3(b) and 3(d)

suggest that W-phase depths are more affected by upper mantle heterogeneity, where the W-phase travel. The measured W-phase depths bias (ranging from 1.2 to 4.5 km in Figs 3d, i and Fig. S1) remains however small relative to the value reported by Hjörleifsdóttir & Ekström (2010) and is of the same order as the depth discretization of the W-phase Green’s function database (ranging from 2 to 10 km, see Section 2.2). The comparison between this work and Hjörleifsdóttir & Ekström (2010), can be problematic because the GCMT algorithm uses different phases (body, mantle and surface waves) and filter passbands (between 40 and 350 s for $M_W = 6.0$) for the inversion. In addition, Hjörleifsdóttir & Ekström (2010) focuses mainly on moderate sized earthquakes ($M_W = 5.5$ and $M_W = 6.0$) while we focus on larger earthquakes ($M_W = 6.5$ and $M_W = 7.5$). Nevertheless, the procedure to assess the effect of 3-D heterogeneity on CMT solutions is similar to our study.

As pointed out in Section 3.1, there are clear geographical patterns in centroid mislocations. Centroid locations in the Circum-Pacific belt seem to be biased outward of the Pacific Ocean while events in the middle east are typically shifted 10–20 km to the north. Such geographical patterns are also reported for GCMT locations that are biased towards West in South America (Hjörleifsdóttir & Ekström 2010; Weston et al. 2011, 2012). Interestingly, GCMT solutions in the western part of the Americas are biased in the direction opposite to W-phase solutions. This is consistent with actual observations showing that GCMT locations are systematically located at the west of W-phase centroids in this region (Duputel et al. 2012a). Such difference is likely related to the fact that GCMT incorporates 3-D corrections (SH8/U4L8 Earth model before 2011 and S362ANI
<figure><img src="image_19.png" /><figcaption>**Figure 8.** Moment tensor resolution matrix as a function of depth for an earthquake in Nepal (lat = 27.91°N, lon = 85.33°E). Panels (a)–(e) correspond to the resolution matrix at different depths. Panel (f) shows the variation in depth of the diagonal elements of the resolution matrix, horizontal lines show the Moho and middle crust depth for STW105 (1-D model) and CRUST2.0 (3-D model).</figcaption></figure>

since 2011, Dziewonski *et al.* 1992) while W-phase solutions are computed at much longer period assuming a spherical model.

Unmodelled 3-D heterogeneity can cause errors in earthquake locations that can be amplified by an uneven distribution of seismological stations with limited path coverage (Bondár *et al.* 2004; Bai *et al.* 2006; Bondár & Storchak 2011). To evaluate the effect of station distribution on W-phase estimates, we compute another set of 3-D synthetic seismograms assuming an azimuthally well balanced array around each earthquake. For each event in our database, we assume a circular array made of stations distributed at 15°, 30°, 45°, 60° and 75° of epicentral distances (see Fig. S2). Fig. 6 compares centroid mislocations obtained using the actual global seismological network and an azimuthally well balanced array of stations. It can

readily be observed that the use of a well balanced array yields to significantly smaller mislocations. Geographical patterns are in fact similar to that observed for the real global network but the average centroid mislocation is reduced from 17 to 11 km. This suggests that largest centroid mislocations are actually due to a combined effect of lateral velocity variations and limited path coverage (e.g. with paths only crossing the Pacific Ocean).

## 4.2 Variation of uncertainty as a function of depth

The impact of 3-D heterogeneity on W-phase solutions changes dramatically with the source depth. Fig. 7 shows that uncertainties
<figure><img src="image_20.png" /><figcaption>**Figure 9.** Comparison of Green's function waveforms for 1-D and 3-D Earth models. These seismograms are computed for an earthquake in Nepal at a depth of 35.5 km for the station BFO (II, 60° from the epicentre). Green's functions are for each moment tensor element computed assuming $M_{ij} = 10^{28}$ dyne.cm (where i, j are moment tensor indexes). Waveforms are filtered in the 150–500 s bandpass using a causal 4th order butterworth filter. The two red dots indicate the time-window of the W-phase. The blue waveform corresponds to the 1-D Earth model STW105 and the red waveform corresponds to the 3-D Earth model S362ANI with Crust2.0.</figcaption></figure>

on the magnitude, focal mechanism and centroid time are significantly larger at shallow depth. Magnitude differences larger than 0.1, beachball RMS differences larger than 0.06 and centroid time differences larger than 3 s are only found for earthquakes shallower than 50 km. To analyse these results, we explore how the moment tensor resolution evolves as a function of depth. Let us first write the W-phase moment tensor solution obtained assuming a 1-D (spherical) Earth model as:

$$ \hat{\mathbf{m}} = \mathbf{G}_{1D}^{-g} \mathbf{d}. \qquad (8) $$

In this equation, $\hat{\mathbf{m}}$ is the inverted moment tensor, $\mathbf{d}$ is the data vector (i.e. concatenated W-phase waveforms) and $\mathbf{G}_{1D}^{-g}$ is the least-squares generalized inverse operator computed assuming a 1-D Earth model (i.e. $\mathbf{G}_{1D}^{-g} = (\mathbf{G}_{1D}^T \mathbf{G}_{1D})^{-1} \mathbf{G}_{1D}^T$, where $\mathbf{G}_{1D}$ is the 1-D Green's function matrix; Kanamori & Rivera 2008). The data vector can be represented as

$$ \mathbf{d} = \mathbf{G}_{3D}\mathbf{m} + \mathbf{\epsilon}, \qquad (9) $$

where $\mathbf{G}_{3D}$ are Green's functions accounting for 3-D heterogeneity, $\mathbf{m}$ is the 'true' moment tensor (used to compute SPECFEM3D_GLOBE synthetics) and $\mathbf{\epsilon}$ is the ambient seismic noise in the data. We can then rewrite eq. (8):

$$ \tilde{\mathbf{m}} = \mathbf{R}_{3D}\mathbf{m} + \mathbf{G}_{1D}^{-g}\mathbf{\epsilon}, \qquad (10) $$

<figure><img src="image_21.png" /><figcaption>**Figure 10.** Comparison between input and output moment tensor solutions (for $M_W = 7.5$ earthquakes and the S362ANI Earth model). Blue dots correspond to shallow earthquakes (depth < 50 km) and red dots correspond to event depths larger than 50 km. We see a larger dispersion for $M_{r\theta}$ and $M_{r\phi}$ components for shallow earthquakes.</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>**Figure 11.** Comparison between the scalar moment ratio and the ratio of $\sin(2\delta)$ ($\delta$ is the dip angle) for dip-slip earthquakes in our catalogue. These results correspond to solutions obtained for $M_W = 7.5$ and the S362ANI model with a rake angle $80° < |\lambda| < 100°$. Earthquakes at depth shallower than 50 km are colored in blue while deeper events are colored in red.</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>**Figure 12.** Relative difference between Green's functions computed for 1-D and 3-D Earth models for an earthquake in Nepal (same event as in Figs 8 and 9). The relative difference is defined as $D = \sum_{i=1}^{N} \int |s_i^{3D}(t) - s_i^{1D}(t)|dt / \sum_{i=1}^{N} \int |s_i^{3D}(t)|dt$, where $s_i^{1D}(t)$ and $s_i^{3D}(t)$ are, respectively, 1-D and 3-D seismograms, $i$ is the channel index and $N$ is the total number of channels in our data set ($N \sim 360$). We note larger difference for Green's functions calculated for $M_{r\theta}$ and $M_{r\phi}$. To a lesser extent, we also note that differences is generally larger for $M_{rr}$ than for the other moment tensor elements.</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>**Figure 13.** Centroid time differences in Eurasia, for different Earth models and $M_W = 7.5$ earthquakes. (a) histograms of centroid time difference for the 1-D crust and 1-D mantle (red), model S362ANI (black) and the 1-D crust and 1-D mantle (blue) and (b), (c) and (d) show the corresponding distribution of centroid time differences.</figcaption></figure>

where the resolution matrix $\mathbf{R}_{3D}$ is defined as:

$$ \mathbf{R}_{3D} = \mathbf{G}_{1D}^{-g} \mathbf{G}_{3D}. \qquad (11) $$

The matrix $\mathbf{R}_{3D}$ provides a direct measure of the impact of 3-D structure on moment tensor inversion results. From eq. (10) the difference between inverted and input moment tensor can then be

written as:

$$ \tilde{\mathbf{m}} - \mathbf{m} = (\mathbf{R}_{3D} - \mathbf{I})\mathbf{m} + \mathbf{G}_{1D}^{-g} \epsilon \qquad (12) $$

where $\mathbf{I}$ is the identity matrix. If the resolution matrix is not the identity $\mathbf{I}$, then the error in the estimated moment tensor will result both from the impact of 3-D heterogeneity (i.e. the term $(\mathbf{R}_{3D} - \mathbf{I})\mathbf{m}$) and from the propagation of noise on the estimated solution (i.e. the term $\mathbf{G}_{1D}^{-g} \epsilon$). We can see that when the values of $\mathbf{R}_{3D}$ are close to $\mathbf{I}$, the moment tensor components are not significantly affected by 3-D heterogeneity (error in model estimates will then only be caused by ambient noise $\epsilon$).

Fig. 8 illustrates the variation of $\mathbf{R}_{3D}$ with depth for an earthquake in Nepal (lat = 27.9°N, lon = 85.3°E), where the crust is thicker than the 1-D model used to compute the Green's functions. We note that off-diagonal elements of $\mathbf{R}_{3D}$ are close to zero, indicating that elements on the diagonal can roughly be interpreted as scaling factors between true and inverted moment tensor components. Results show that most diagonal elements of $\mathbf{R}_{3D}$ remain close to one, with the exception of those associated with $M_{r\theta}$ and $M_{r\phi}$ moment tensor components. The elements of $\mathbf{R}_{3D}$ corresponding to $M_{r\theta}$ and $M_{r\phi}$ show large variations as function of depth. Other examples provided in the online supplementary for earthquakes in Chile and Turkey also show larger variations for $M_{r\theta}$ and $M_{r\phi}$ components (see Figs S5 and S6). In Nepal, from the surface to the middle crust of the 1-D model they are smaller than one (Fig. 8a). These elements then increase to reach ~2 in the vicinity of the 1-D Moho and the 3-D middle crust (Fig. 8c). Finally, they decrease and stay close to one at depths larger than the 3-D Moho (Fig 8e). This depth variation of $\mathbf{R}_{3D}$ indicates that $|M_{r\theta}|$ and $|M_{r\phi}|$ will systematically be overestimated below the 1-D middle crust in Nepal. This is consistent with inversion results in the region presented in Fig. 3, showing that the scalar moment is overestimated at depth larger than 15 km (i.e. below the 1-D middle crust). Fig. 9 compares 1-D and 3-D Green's functions computed at a centroid depth of 35.5 km. It clearly appears that 1-D Green's functions for $M_{r\theta}$ and $M_{r\phi}$ are smaller that Green's functions computed in a 3-D Earth, which explains the overestimation of $M_{r\theta}$ and $M_{r\phi}$ and the corresponding large values in $\mathbf{R}_{3D}$.

More generally, if we compare inverted moment tensor solutions with the input values (as shown in Fig. 10), we observe larger dispersion on $M_{r\theta}$ and $M_{r\phi}$ elements, for shallower events. Such, larger dispersion on $M_{r\theta}$ and $M_{r\phi}$ can also be observed when comparing W-phase and GCMT solutions (Duputel et al. 2012b). For shallow dip-slip earthquakes, the larger uncertainty on $M_{r\theta}$ and $M_{r\phi}$ results in a moment-dip trade-off, where different solutions corresponding to a constant value of $M_0\sin 2\delta$ (where $\delta$ is the dip angle) will be associated with the same data misfit (e.g. Bukchin 2006). We should therefore have the relation

$$ \frac{M_0^{\text{output}}}{M_0^{\text{input}}} = \frac{\sin(2\delta_{\text{input}})}{\sin(2\delta_{\text{output}})}. \qquad (13) $$

Fig. 11 compares the moment ratio $M_0^{\text{output}}/M_0^{\text{input}}$ with the right side of eq. (13), for earthquakes with dip-slip mechanisms in our catalogue (i.e. with a rake angle 80° < $|\lambda|$ < 100°). Although there is some variability induced by ambient noise and 3-D heterogeneity, the estimated scalar moments and dip angles for earthquakes deeper than 50 km are relatively close to the input values. On the other hand, shallow earthquakes clearly depict an additional variability that can largely be attributed to the $M_0\sin 2\delta$ trade-off. There is some scatter around the trade-off curve in eq. (13). This can be induced by observational and modeling errors and can also be related to deviations from a pure dip-slip mechanism.
<table><caption>**Table 2.** Overall impact of 1-D heterogeneity on source parameters. For the two studied 1-D Earth models and for both $M_W = 7.5$ and $M_W = 6.5$ earthquakes, we show the mean and standard deviation of the magnitude difference ($\Delta M_W$), the beachball RMS difference ($\Delta$), the centroid time difference ($\Delta \tau$), the centroid depth difference ($\Delta h$) and the centroid mislocation ($|\Delta x|$).</caption><thead><tr><th>Model</th><th></th><th>ΔM<sub>W</sub></th><th>Δ</th><th>Δτ [s]</th><th>Δh [km]</th><th>|Δx| [km]</th></tr></thead><tbody><tr><td>S362ANI</td><td>(M<sub>W</sub> = 7.5)</td><td>0 ± 0.03</td><td>0.05 ± 0.04</td><td>-0.4 ± 1.9</td><td>3.1 ± 6.5</td><td>17.0</td></tr><tr><td>S40RTS</td><td>(M<sub>W</sub> = 7.5)</td><td>0 ± 0.03</td><td>0.05 ± 0.04</td><td>0 ± 2</td><td>1.2 ± 5.5</td><td>17.3</td></tr><tr><td>S362ANI</td><td>(M<sub>W</sub> = 6.5)</td><td>0 ± 0.03</td><td>0.06 ± 0.05</td><td>0.1 ± 2</td><td>4.5 ± 5.7</td><td>16.5</td></tr><tr><td>S40RTS</td><td>(M<sub>W</sub> = 6.5)</td><td>0 ± 0.03</td><td>0.06 ± 0.05</td><td>0.6 ± 2.1</td><td>2 ± 5.4</td><td>17.3</td></tr></tbody></table>

The larger uncertainty on $M_{r\theta}$ and $M_{r\phi}$ estimates at long period is usually interpreted as a direct consequence of the low amplitude of the associated Green's functions close to the free surface (e.g. Dziewonski et al. 1981; Kanamori & Given 1982). Moreover, we note in Fig. 9 that Green's functions associated with $M_{r\theta}$ and $M_{r\phi}$ are also more affected by the 3-D structure than the other components. This is confirmed in Figs 12 and S7 showing that, for earthquakes at different locations, the relative impact of lateral heterogeneity is generally larger for $M_{r\theta}$ and $M_{r\phi}$. It shows that 3-D heterogeneity in the source region is a significant contributor to the observed differences for $M_{r\theta}$ and $M_{r\phi}$ components of shallow earthquakes.

## 4.3 Anticorrelation between crustal and mantle effects

In Section 3.4, we found that 3-D heterogeneity in the crust and in the mantle seems to have an anticorrelated impact on spatio-temporal centroid coordinates.

For example, we see that 3-D crustal heterogeneity induces positive centroid time delays in Eurasia and South America, while mantle heterogeneity induces negative delays in the same regions (see Figs 5b and e). In the same way, the crust induces negative delays at the interface between Australian, Pacific and Antarctic plates, while mantle heterogeneity induces positive delays. Although such anticorrelation is less clear for centroid locations, we see that the crust induces a westward bias in South America, while the mantle biases locations towards the east in the same area (see Figs 5c and f).

Fig. 13 presents a detailed analysis of the relative effect of crust and mantle heterogeneity on centroid times estimated in Eurasia, where the differences of 3Dc–1Dm (3-D crust and 1-D mantle) and 1Dc–3Dm (1-D crust and 1-D mantle) are quite clear. Fig. 13(a) shows the distribution of centroid time differences ($\Delta \tau$) for 3Dc–1Dm, 1Dc–3Dm and the full 1-D model assuming $M_W = 7.5$ earthquakes (W-phase solutions obtained using a 150–500 s bandpass filter). It clearly illustrates the anticorrelation between crustal and mantle effects on centroid times that is also visible in Figs 13(b) and (d). On the other hand, we note that values of $\Delta \tau$ for the full 1-D model are smaller, suggesting that mantle and crustal effects compensates. Such anticorrelation is also found for $M_W = 6.5$ earthquakes (using a 100–250 s bandpass filter) in Fig. S8 that also show a larger impact of crustal heterogeneity on centroid times (resulting into larger $\Delta \tau$ values). The larger crustal effects for $M_W = 6.5$ earthquakes is most likely due to the use of shorter period waveforms that are more sensitive to shallow heterogeneity.

The 1-D Global mantle models used in this study are obtained by assuming a specific crustal model (the CRUST2.0 model). The crust structure is then linked to the mantle as these models rely on a fixed crustal model to fit the observations. We think that this could explain partly the observed anticorrelation of centroid times estimates in Fig. 13. Such mutual cancellation of crustal and mantle

effects has also been reported in another context by Koelemeijer (2014).

## 5 CONCLUSIONS

We quantify the impact of lateral structural heterogeneity on W-phase CMT solutions by simulating synthetic earthquakes in 1-D Earth models and comparing the retrieved source parameters with the input values. To assess the influence of long-period noise, we add actual ambient noise to synthetics assuming two earthquake magnitudes ($M_W = 6.5$ and $M_W = 7.5$). The results obtained for 1-D Earth models S362ANI and S40RTS and moment magnitudes $M_W = 6.5$ and 7.5 are summarized in Table 2.

On average, the impact of 1-D heterogeneity is relatively small on retrieved W-phase solutions. The resulting discrepancy on moment magnitude estimates is smaller than 0.03 for 80 per cent of earthquakes. Beachball RMS differences are smaller than 0.1 for ~90 per cent of earthquakes. The impact on centroid time is of the same order as the sampling period. Overall, lateral heterogeneity mainly impacts centroid locations, especially if the station coverage is limited. This is particularly visible in the Circum-Pacific belt for which many earthquakes have paths that only cross the Pacific Ocean for large ranges of source to station azimuths.

Results show that the impact of 1-D heterogeneity is consistent between Earth models S362ANI and S40RTS, suggesting that the W-phase is mainly sensitive to large-scale heterogeneity that exists in both models. Although uncertainties estimated for $M_W = 6.5$ earthquakes are often larger than for $M_W = 7.5$ events, shifting the passband towards higher frequencies clearly mitigates the impact of long period ambient noise. The remaining larger dispersion observed for $M_W = 6.5$ events likely results from larger 1-D effects when using shorter period waveforms.

Results also indicate that most source parameters seems to be mostly affected by mantle heterogeneity, where most of the W-phase propagating energy is confined. However, crustal heterogeneity seems to have a larger impact on scalar moment estimates for shallow earthquakes. In general, the uncertainty on scalar seismic moment, focal mechanism and centroid time increases significantly at shallow depth. This can be attributed to the difficulty to constrain $M_{r\theta}$ and $M_{r\phi}$ moment tensor elements at long period due to the small amplitude of the associated Green's functions for shallow earthquakes. In addition, our tests indicate that Green's functions for $M_{r\theta}$ and $M_{r\phi}$ are more impacted by lateral heterogeneity than the other components.

Although the impact of lateral heterogeneity on W-phase solutions is relatively moderate on average, source models can clearly be improved by incorporating 1-D Green's functions (e.g. Ferreira et al. 2011; Duputel et al. 2016; Hejrani et al. 2017). While it is still challenging for global real-time applications, the improvement of
computational capabilities makes rapid 1-D CMT inversions already possible at regional scale (Wang & Zhan 2019). As lateral heterogeneity can have non-negligible impact on earthquake locations, such improvements could be instrumental for rapid earthquake response and tsunami warning purposes.

## ACKNOWLEDGEMENTS

The authors would like to acknowledge the Editor, Ana Ferreira, an anonymous reviewer and Göran Ekström for their comments. We also would like to thank Stephen Beller, Hiroo Kanamori, Javier Ruiz and Sophie Lambotte for helpful discussion at various stages of this work. Some of the computations were performed on the facilities of the High Performance Computing Center of the University of Strasbourg. This work was supported by the Initiative d'Excellence (IDEX) funding framework (Université de Strasbourg). This project has also received funding from the European Research Council (ERC, under the European Union's Horizon 2020 research and innovation programme under grant agreement No. 805256) and from Agence Nationale de la Recherche (project ANR-17-ERC3-0010). We wish to thank the network and station operators for their commitment to collect high-quality seismic data. Some figures were made using the Generic Mapping Tools (GMT) software version 5.3.1 (Wessel et al. 2013).

## REFERENCES

* Albuquerque Seismological Laboratory (ASL)/USGS, 1988. Global Seismograph Network (GSN - IRIS/USGS). International Federation of Digital Seismograph Networks.

* Albuquerque Seismological Laboratory (ASL)/USGS, 1992. New China Digital Seismograph Network. International Federation of Digital Seismograph Networks.

* Albuquerque Seismological Laboratory (ASL)/USGS, 1993. Global Telemetered Seismograph Network (USAF/USGS). International Federation of Digital Seismograph Networks.

* Albuquerque Seismological Laboratory (ASL)/USGS, 2006. Caribbean USGS Network. International Federation of Digital Seismograph Networks.

* Bai, L., Wu, Z., Zhang, T. & Kawasaki, I., 2006. The effect of distribution of stations upon location error: statistical tests based on the double-difference earthquake location algorithm and the bootstrap method, *Earth, Planets Space*. **58**(2), e9–e12, doi.org/10.1186/BF03353364.

* Bassin, C., Laske, G. & Masters, G., 2000. The current limits of resolution for surface wave tomography in North America, *EOS, Trans. Am. Geophys. Un.*, **81**(F897).

* Bondár, I. & Storchak, D., 2011. Improved location procedures at the international seismological centre, *J. geophys. Int.*, **186**(3), 1220–1244, doi.org/10.1111/j.1365-246X.2011.05107.x.

* Bondár, I., Myers, S.C., Engdahl, E.R. & Bergman, E.A., 2004. Epicentre accuracy based on seismic network criteria, *J. geophys. Int.*, **156**(3), 483–496, doi.org/10.1111/j.1365-246X.2004.02070.x.

* Bukchin, B.G., 2006. Specific features of surface wave radiation by a shallow source, *Izvest., Phys. Solid Earth*, **42**(8), 712–717, doi.org/10.1134/S1069351306080088.

* California Institute of Technology and United States Geological Survey Pasadena (SCSN), 1926. Southern California Seismic Network. International Federation of Digital Seismograph Networks.

* Crowell, B.W., Bock, Y. & Melgar, D., 2012. Real-time inversion of GPS data for finite fault modeling and rapid hazard assessment, *Geophys. Res. Lett.*, **39**(9), doi.org/10.1029/2012GL051318.

* Dahlen, F. & Tromp, J., 1998. *Theoretical Global Seismology*, Princeton Univ. Press.

* Delouis, B., Nocquet, J.-M. & Vallée, M., 2010. Slip distribution of the February 27, 2010 Mw = 8.8 Maule earthquake, central chile, from static

and high-rate GPS, InSAR, and broadband teleseismic data, *Geophys. Res. Lett.*, **37**(17), doi.org/10.1029/2010GL043899.

* Duputel, Z., Rivera, L., Fukahata, Y. & Kanamori, H., 2012a. Uncertainty estimations for seismic source inversions, *J. geophys. Int.*, **190**(2), 1243–1256, doi.org/10.1111/j.1365-246X.2012.05554.x.

* Duputel, Z., Rivera, L., Kanamori, H. & Hayes, G., 2012b. W phase source inversion for moderate to large earthquakes (1990–2010), *J. geophys. Int.*, **189**(2), 1125–1147, doi.org/10.1111/j.1365-246X.2012.05419.x.

* Duputel, Z., Agram, P.S., Simons, M., Minson, S.E. & Beck, J.L., 2014. Accounting for prediction uncertainty when inferring subsurface fault slip, *J. geophys. Int.*, **197**(1), 464–482, doi.org/10.1093/gji/ggt517.

* Duputel, Z., Vergne, J., Rivera, L., Wittlinger, G., Farra, V. & Hetényi, G., 2016. The 2015 Gorkha earthquake: a large event illuminating the main Himalayan thrust fault, *Geophys. Res. Lett.*, **43**(6), 2517–2525, doi.org/10.1002/2016GL068083.

* Dziewonski, A., Franzen, J. & Woodhouse, J., 1984. Centroid-moment tensor solutions for January-March, 1984, *Phys. Earth planet. Inter.*, **34**(4), 209–219.

* Dziewonski, A.M., Chou, T.-A. & Woodhouse, J.H., 1981. Determination of earthquake source parameters from waveform data for studies of global and regional seismicity, *J. geophys. Res.*, **86**(B4), 2825–2852, doi.org/10.1029/JB086iB04p02825.

* Dziewonski, A.M., Ekström, G. & Salganik, M.P., 1992. Centroid-moment tensor solutions for July–September 1991, *Phys. Earth planet. Inter.*, **72**(1–2), 1–11.

* Ekström, G., Nettles, M. & Dziewoski, A., 2012. The global CMT project 2004–2010: centroid-moment tensors for 13,017 earthquakes, *Phys. Earth planet. Inter.*, **200–201**, 1–9.

* Ferreira, A.M.G. & Woodhouse, J.H., 2006. Long-period seismic source inversions using global tomographic models, *J. geophys. Int.*, **166**(3), 1178–1192, doi.org/10.1111/j.1365-246X.2006.03003.x.

* Ferreira, A. M.G., Weston, J. & Funning, G.J., 2011. Global compilation of interferometric synthetic aperture radar earthquake source models: 2. Effects of 3-D earth structure, *J. geophys. Res.*, **116**(B8), doi.org/10.1029/2010JB008132.

* GEOFON Data Centre, 1993. Geofon seismic network. deutsches geoforschungszentrum gfz.

* Geological Survey of Canada, 1980. Canadian National Seismograph Network (CNSN). International Federation of Digital Seismograph Networks.

* Gombert, B., Duputel, Z., Jolivet, R., Simons, M., Jiang, J., Liang, C., Fielding, E. & Rivera, L., 2018. Strain budget of the Ecuador–Colombia subduction zone: a stochastic view, *Earth planet. Sci. Lett.*, **498**, 288–299, doi.org/10.1016/j.epsl.2018.06.046.

* Hallo, M. & Gallovič, F., 2016. Fast and cheap approximation of green function uncertainty for waveform-based earthquake source inversions, *J. geophys. Int.*, **207**(2), 1012–1029, doi.org/10.1093/gji/ggw320.

* Hallo, M., Asano, K. & Gallovič, F., 2017. Bayesian inference and interpretation of centroid moment tensors of the 2016 Kumamoto earthquake sequence, Kyushu, Japan, *Earth, Planets Space*, **69**(1), 134.

* Hayes, G.P., Rivera, L. & Kanamori, H., 2009. Source inversion of the W-phase: real-time implementation and extension to low magnitudes, *Seismol. Res. Lett.*, **80**(5), 817.

* Hejrani, B., Tkalčić, H. & Fichtner, A., 2017. Centroid moment tensor catalogue using a 3-D continental scale earth model: application to earthquakes in Papua New Guinea and the Solomon Islands, *J. geophys. Res.*, **122**(7), 5517–5543.

* Hjörleifsdóttir, V. & Ekström, G., 2010. Effects of three-dimensional earth structure on cmt earthquake parameters, *Phys. Earth planet. Inter.*, **179**(3), 178–190.

* Ide, S., 2007. 4.07 - slip inversion, in *Treatise on Geophysics*, pp. 193–223, ed. Schubert, G., Elsevier.

* Institut de Physique du Globe de Paris (IPGP) and Ecole et Observatoire des Sciences de la Terre de Strasbourg (EOST), 1982. Geoscope, French Global Network of Broad Band Seismic Stations. Institut de Physique du Globe de Paris (IPGP).

* Kanamori, H., 1993. W phase, *Geophys. Res. Lett.*, **20**(16), 1691–1694.
Kanamori, H. & Given, J.W., 1982. Use of long-period surface waves for rapid determination of earthquake source parameters 2. Preliminary determination of source mechanisms of large earthquakes (ms $>= 6.5$) in 1980, *Phys. Earth planet. Inter.*, **30**(2), 260–268.

Kanamori, H. & Rivera, L., 2008. Source inversion of W-phase: speeding up seismic tsunami warning, *J. geophys. Int.*, **175**(1), 222–238, doi.org/10.1111/j.1365-246X.2008.03887.x.

Koelemeijer, P., 2014, Normal mode studies of long wavelength structures in Earth's lowermost mantle, *PhD thesis*, University of Cambridge.

Komatitsch, D. & Tromp, J., 2002a. Spectral-element simulations of global seismic wave propagation – I. Validation, *J. geophys. Int.*, **149**(2), 390–412.

Komatitsch, D. & Tromp, J., 2002b. Spectral-element simulations of global seismic wave propagation – II. Three-dimensional models, oceans, rotation and self-gravitation, *J. geophys. Int.*, **150**(1), 303–318.

Komatitsch, D., et al., 2015. SPECFEM3D globe v7.0.0 [software].

Kustowski, B., Ekström, G. & Dziewoski, A.M., 2008. Anisotropic shear-wave velocity structure of the earth's mantle: a global model, *J. geophys. Res.*, **113**(B6), doi.org/10.1029/2007JB005169.

Liu, Q., Polet, J., Komatitsch, D. & Tromp, J., 2004. Spectral-element moment tensor inversions for earthquakes in southern California, *Bull. seism. Soc. Am.*, **94**(5), 1748.

MedNet Project Partner Institutions, 1990. Mediterranean Very Broadband Seismographic Network (mednet).

Nakanishi, I. & Kanamori, H., 1982. Effects of lateral heterogeneity and source process time on the linear moment tensor inversion of long-period rayleigh waves, *Bull. seism. Soc. Am.*, **72**(6A), 2063.

Northern California Earthquake Data Center, 2014. Berkeley Digital Seismic Network (BDSN) [data set].

Polet, J. & Thio, H.K., 2011. Rapid calculation of a centroid moment tensor and waveheight predictions around the north Pacific for the 2011 off the Pacific coast of Tohoku earthquake, *Earth, Planet Space, 63*, doi.org/10.5047/eps.2011.05.005.

Riquelme, S., Medina, M., Bravo, F., Barrientos, S., Campos, J. & Cisternas, A., 2018. W-phase real-time implementation and network expansion from 2012 to 2017: the experience in Chile, *Seismol. Res. Lett.*, **89**(6), 2237.

Ritsema, J., Deuss, A., van Heijst, H.J. & Woodhouse, J.H., 2011. S40rts: a degree-40 shear-velocity model for the mantle from new Rayleigh wave dispersion, teleseismic traveltime and normal-mode splitting function measurements, *J. geophys. Int.*, **184**(3), 1223–1236.

Rivera, L. & Kanamori, H., 2014. Diagnosing source geometrical complexity of large earthquakes, *Pure appl. Geophys.*, **171**(10), 2819–2840.

Scripps Institution of Oceanography, 1986. Global Seismograph Network (GSN - IRIS/USGS). International Federation of Digital Seismograph Networks.

Silver, P.G. & Jordan, T.H., 1982. Optimal estimation of scalar seismic moment, *Geophys. J. Int.*, **70**(3), 755–787.

Tsai, V.C., Nettles, M., Ekström, G. & Dziewonski, A.M., 2005. Multiple CMT source analysis of the 2004 Sumatra earthquake, *Geophys. Res. Lett.*, **32**(17), L17304.

Wang, D., et al., 2017. Evaluation of W phase CMT based PTWC real-time tsunami forecast model using DART observations: events of the last decade, *AGU Fall Meeting Abstracts*, , abstract #NH21D-02.

Wang, X. & Zhan, Z., 2019. Moving from 1-D to 3-D velocity model: automated waveform-based earthquake moment tensor inversion in the Los Angeles region, *Geophys. J. Int.*, **220**(1), 218–234, doi.org/10.1093/gji/ggz435.

Wessel, P., Smith, W.H.F., Scharroo, R., Luis, J. & Wobbe, F., 2013. Generic mapping tools: improved version released, *EOS, Trans. Am. Geophys. Un.*, **94**(45), 409–410.

Weston, J., Ferreira, A. M.G. & Funning, G.J., 2011. Global compilation of interferometric synthetic aperture radar earthquake source models: 1. Comparisons with seismic catalogs, *J. geophys. Res.*, **116**(B8), doi.org/10.1029/2010JB008131.

Weston, J., Ferreira, A.M. & Funning, G.J., 2012. Systematic comparisons of earthquake source models determined using InSAR and seismic data, *Tectonophysics*, **532-535**, 61–81, doi.org/10.1016/j.tecto.2012.02.001.

Woodhouse, J.H. & Dziewonski, A.M., 1984. Mapping the upper mantle: three-dimensional modeling of earth structure by inversion of seismic waveforms, *J. geophys. Res.*, **89**(B7), 5953–5986, doi.org/10.1029/JB089iB07p05953.

Yagi, Y. & Fukahata, Y., 2011. Introduction of uncertainty of green's function into waveform inversion for seismic source processes, *J. geophys. Int.*, **186**(2), 711–720, doi.org/10.1111/j.1365-246X.2011.05043.x.

Zhao, X., Duputel, Z. & Yao, Z., 2017. Regional w-phase source inversion for moderate to large earthquakes in china and neighboring areas, *J. geophys. Res.*, **122**(12), 10 052–10 068, doi.org/10.1002/2017JB014950.

## SUPPORTING INFORMATION

Supplementary data are available at GJI online.

**Figure S1** Effects of 1-D heterogeneity for $M_W = 6.5$ earthquakes. Maps show the following quantities from top to bottom: moment magnitude difference, beachball RMS difference, difference in time-shift, depth difference and horizontal mislocation. Panels from (a) to (e) correspond to the inversion of synthetic seismograms created with the Earth model S362ANI. Panels (f) to (j) correspond to the inversion of synthetic seismograms created with the Earth model S40RTS. In both cases, seismic noise has been added to synthetics prior to source inversion (see Section 2.1 of the main text).

**Figure S2** Example of circular array for an earthquake in the middle of the Pacific. The inset on the lower left-hand corner shows the station distribution from the actual seismological network that would have been used for that earthquake.

**Figure S3** Effects of crustal and mantle heterogeneity on W-phase solutions obtained for $M_W = 7.5$ earthquakes. Top maps indicate beachball RMS difference ($\Delta$) and bottom map show the depth difference ($\Delta h$). Panels (a) and (b) correspond to the inversion of synthetic seismograms created with a 1-D crust (CRUST2.0) and a 1-D mantle (STW105). Panels (c) and (d) correspond to the results of the inversion of synthetic seismograms created with a 1-D crust (STW105) and a 1-D mantle (S362ANI). As described in the main text, ambient noise has been added to synthetics prior to source inversion.

**Figure S4** Effects of crustal and mantle heterogeneity on W-phase solutions obtained for $M_W = 6.5$ earthquakes. Maps show from top to bottom: moment magnitude difference, beachball RMS difference, difference in time-shift, depth difference and horizontal mislocation. Panels from (a) to (e) correspond to the inversion of synthetic seismograms created with a 1-D crust (CRUST2.0) and a 1-D mantle (STW105). Panels (f) to (j) correspond to the inversion of synthetic seismograms created with a 1-D crust (STW105) and a 1-D mantle (S362ANI). As described in the main text, ambient noise has been added to synthetics prior to source inversion.

**Figure S5** Moment tensor resolution matrix as a function of depth for an earthquake in Chile (lat = 40.64°S, lon = 75.37°W). Figures from (a) to (e) correspond to the resolution matrix at different depths. Panel (f) shows the variation in depth of the value of the diagonal components of the resolution matrix, horizontal lines show the Moho and middle crust depth for STW105 (1-D model) and crust 2.0 (1-D model).

**Figure S5** Moment tensor resolution matrix as a function of depth an earthquake in Turkey (lat = 41.01°N, lon = 29.97°E). Panels from (a) to (e) correspond to the resolution matrix at different depths. Panel (f) shows the variation in depth of the value of the diagonal components of the resolution matrix, horizontal lines show the Moho and middle crust depth for STW105 (1-D model) and crust 2.0 (1-D model).
**Figure S6** Relative difference between Green's functions computed for 1-D and 1-D Earth models for earthquakes in (a) Chile and (b) Turkey (same events as in Figs S5 and S6). The relative difference is defined as in Fig. 12 of the main text.

**Figure S7** Centroid time differences in Eurasia, for different Earth models and $M_W = 6.5$ earthquakes. (a) histograms of centroid time difference for the 1-D crust and 1-D mantle (red), model S362ANI

(black) and the 1-D crust and 1-D mantle (blue) and (b), (c) and (d) show the corresponding distribution of centroid time differences.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper.