Improved finite-source inversion through joint measurements of rotational and translational ground motions: A theoretical study

With the prospects of seismic equipment being able to measure rotational ground motions in a wide frequency and amplitude range in the near future we engage in the question how this type of ground motion observation can be used to solve the seismic inverse problem. In this paper, we focus on the question, whether finite source inversion can benefit from additional observations of rotational motion. Keeping the overall number of traces constant, we compare observations from a surface seismic network with 44 3-component translational sensors (classic seismometers) with those obtained with 22 6-component 5 sensors (with additional 3-component rotational motions). Synthetic seismograms are calculated for known finite-source properties. The corresponding inverse problem is posed in a probabilistic way using the Shannon information content as measure how the observations constrain the seismic source properties. We minimize the influence of the source receiver geometry around the fault by statistically analyzing six-component (three velocity and three rotation rate) inversions with a random distribution of receivers. The results show that with the 6-C subnetworks the source properties are not only equally well recovered 10 (even that would be benefitial because of the substantially reduced logistics installing half the sensors) but statistically some source properties are almost always better resolved. We assume that this can be attributed to the fact that the (in particular vertical) gradient information is contained in the additional motion components. We compare these effects for strike-slip and normal-faulting type sources and confirm that the increase in inversion quality for kinematic source parameters is even higher for the normal fault. This indicates that the inversion benefits from the additional information provided by the horizontal rota15 tion rates, i.e. information about the vertical displacement gradient.

to be infinitesimal, the displacement of a particle at position u(x + δx) can be described with the following equation (Aki and Richards, 2002): In Equation 1, is the symmetric strain tensor, therefore containing six instead of nine independent components, and ω = ( 1 2 ∇ × u) are the rotations around each of the three translational axes of u. Therefore the full description of a dis-5 placement around a point x needs three translational, three rotational and six strain components. The rotation components are mathematically related to the space derivatives of the translational motion components u x , u y , u z : Although rotational ground motions are caused by earthquake source processes, measurements thereof have not been carried out until the 1990s (Nigbor, 1994;Takeo, 1998). Nevertheless, the possible benefit has been considered decades before - 10 in the first edition of Aki andRichards, 2002, which was published in 1980. Even then, the methods that were used mostly consisted of several three-component seismometers and were not able to be used for seismology because of their low sensitivity (Suryanto, 2006). Since then, many scientists commended the benefits of studying rotational ground motions for earth sciences (Twiss et al., 1993;Spudich et al., 1995;Takeo and Ito, 1997;Igel et al., 2007;Fichtner and Igel, 2009), physics (DeSalvo, 2009;Lantz et al., 2009) and also for engineering applications (Trifunac, 2009). In the last years, new instruments, e.g. ring-15 laser technology (Schreiber et al., 2009;Velikoseltsev et al., 2012) or adapted gyroscopes (Bernauer et al., 2012), that are capable of measuring rotational motions with a much higher precision, have been developed and possibly bring new insight into seismological applications.
With the potential availability of seismometers measuring rotational ground motions, the goal of this study is to test the effect of incorporating rotational ground motion data into finite source inversions for different fault mechanisms. Following the 20 approach by (Bernauer et al., 2014), we employ a Bayesian, i.e. probabilistic, inversion algorithm and compare the inversions with different data types and different earthquake source mechanisms. In the last decade there have been several studies showing that probabilistic methods are well suited for ill-posed inversion problems of finite source earthquakes, because they overcome the drawbacks of regularization techniques such as local minima ((Monelli et al., 2009;Fichtner and Tkalčić, 2010)).
We directly compare the results with and without incorporation of rotational ground motion data for a pure strike-slip 25 earthquake and for a dip-slip earthquake, i.e. the strike also propagates in vertical direction. Because rotations around the horizontal axes (i.e. ω x and ω y in Equation (2)) carry information about vertical gradients, which is not possible to resolve with common arrays on the Earth's surface (Bernauer et al., 2014), we expect to see an even greater benefit for inversions for the dip-slip scenario due to the higher energy in the horizontal rotation components. To better assess the quality differences we keep the number of seismograms used in the inversion the same for all scenarios.
Because using three translational and three rotational ground motion recordings doubles the amount of seismograms, we half the number of seismic stations when using the so-called 6C (six-component) data, as compared to only three translational seismograms (3C/three-component data). In the following chapter, we present the Bayesian inversion scheme. Next, we explain the parametrization of the synthetic earthquakes, the spatial distribution of the seismic stations and define the different inversion 5 scenarios. We will then present the inversion results with and without incorporation of the rotational ground motion data for the two different earthquake mechanisms.

Probabilistic Inversion Scheme
Although computation time is drastically increased when applying probabilistic inversion schemes, it avoids drawbacks such as 10 regularization or falling into local minima in the iterations during the process of minimizing the misfit between model and data.
We use the basic ideas of bayesian inversion for a finite fault developed by (Tarantola, 2005) and (Mosegaard and Tarantola, 1995): For each inversion, a certain number of parameters is stored in the vector m, given by In Equation (4), κ is again a normalization constant, χ l (m) denots the measure of misfit between the predicted model data and the observed data and s l represents the noise level. 20 The prior pdf ρ(m) is set to be constant in all model parameters within a certain interval of interest. Since all model parameters are independent, the combined prior pdf can be obtained by simply multiplying the probability densities for all parameters.
In order to calculate the posterior pdf, we let the probabilistic inversion do random walks through the entire model space. We approximate the posterior pdf for all 26 parameters by applying the Metropolis algorithm to 1, 000, 000 test models. Once the entire model space is sampled, one can illustrate it by marginal density distributions of each of the parameters. The shape of 25 the these distributions gives a discrete approximation of the posterior pdf.
To quantify the posterior pdf compared to the prior pdf, i.e. the relative information gain, we use the so-called Shannon's measure of information gain: The relative information gain I(ρ; σ) with respect to the posterior pdf ρ and the prior pdf σ is defined in (Tarantola, 2005) as Since the base of the logarithm in Equation 5 is 2, the unit of I is termed a bit. One major advantage of this method is the possible comparison between results of very different experiments by quantification of the inversion quality. For further details 5 on the method we refer to the publications of (Bernauer et al., 2014) and (Tarantola, 2005).

Forward Model
The We illustrate receiver and fault setup in Figure 1. 44 Receivers, marked as grey, and red triangles, are setup around the source with a spacing of 0.4 • latitude and 0.25 • longitude to guarantee a good coverage in all directions from the fault. Since we want to decrease the number of receivers by a factor of 2 when incorporating rotational ground motion data to keep the amount of data used in the inversion the same, we will limit our recording to the 22 stations marked with red triangles, when using both 15 data types. The fault, illustrated as a black line, has a constant latitude of 0 • and the epicenter of the earthquake is marked with a yellow star. At this point we want to emphasize that due to the theoretical nature of this study the geographical location of fault and receivers are in no way related to real data of past earthquakes or geographical areas.
We brake up the fault plane into 8 · 3 = 24 quadratic subfaults. Each subfault is filled by 8 · 8 = 64 moment tensor point sources with a regular spacing of 500 meters between the point sources. The side-length of each subfault is therefore about 4km. We . Note that there are two parameters (i.e. rupture velocity and rise time) that are not linearly related to displacement or velocity, contrary to the 24 slip amplitudes. We therefore invert for 26 kinematic source parameters.
The parametrization of the slip on the fault is illustrated in Figure 3. We show heterogeneous slip amplitude as varying colors for each subfault, ranging from 0.6 meters in all four corners up to 3.2 meters in the eastern center part of the fault, 5 i.e. one of the 4 subfaults closest to the hypocenter (yellow star). Because rupture velocity (2700m/s) and rise time (0.8s) are homogeneous across the fault, it is sufficient to state the single values, which are valid for the entire fault. We add Gaussian noise to all synthetic seismograms to render our synthetic study more realistic. Because we calculate velocity and rotation rate seismograms individually, the noise on both data types is not correlated. We set the noise level to 10% of the maximum amplitude observed in all recordings for each respective data type for all scenarios. This ensures that a potential increase or 10 decrease in inversion quality does not result from different noise levels.
We compare the inversions for two different fault mechanisms. For Scenario I we choose a pure left-lateral strike-slip event (rake parallel in strike direction). For Scenario II, its dip is changed from 90 • to 45 • and the sliding now occurs perpendicular to strike direction (rake λ = −90 • ). The fault shows pure dip-slip behavior and fractures also in vertical direction ( Figure 1).

15
If we compute the total energy in all three rotational components measured at all stations, the rotation rates in the horizontal components make up a very large part for the dip-slip fault ( 57%), whereas the strike-slip fault features the highest rotation rate ratio in the vertical components( 40% for rotations around both horizontal axes), but horizontal rotation only contributes 20% of the total energy. This makes sense, since we created our strike-slip earthquake to solely rupture in horizontal direction.
To better assess the relevance of this study we check if the synthetic absolute rotation rates are in a range where they can 20 be measured by state-of-the-art measurement devices. We therefore calculate the maximum amplitude for all 44 receivers for each of the three rotation rate components. They are shown in Figure 4.
The logarithmic scale for the rotation rates on the x-axes helps to get a sense for the amplitudes' order of magnitude. The  (Table 1). This indicates that shallower subfaults are resolved better in the inversion and rotation rates are likely to improve these subfaults more than deeper subfaults. It seems like the bottom row of the subfaults does not show any improvement by the incorporation of rotational ground motion data. This agrees with several studies (e.g. (Bernauer et al., 2014;Cotton and Campillo, 1995;Semmane et al., 2005;Monelli and Mai, 2008;Mendoza and Hartzell, 1989)). The total cumulative information gain over all 24 subfaults is with 14.98bit about 9.86% larger 20 for the 6C inversion than for the 3C inversion with 13.63bit.
We illustrate the posterior pdf's for the two parameters that are not linearly related to displacement, i.e. rupture velocity and rise time, in Figure 6. We use the same labelling as in Figure 5. We cannot resolve the exact true values very well, but if one compares the information gain between the 3C inversion and the 6C inversion, there is a clear increase when incorporating rotation rate data into the inversion scheme. There is a 45.02% increase in information gain for inverting for rupture velocity 25 and a 17.77% increase for inverting for rise time. This is again in agreement with the results from (Bernauer et al., 2014) where the true values are also not captured and the increase in information gain is around the same order for both values.
In Figure 7 we illustrate the results for a dip-slip fault in the same manner as Figure 5. The decrease of inversion quality with depth is clearly visible again. A comparison between the cumulative information gain per row leads us to the assumption that slip amplitudes in the shallowest subfault layer seem well resolved whereas the deeper layers tend to be poorly resolved in the bit(6C)). We use these numbers to calculate the increase in information gain for the 6C inversion for each row (Table 1). The highest increase can be seen in the top row. The bottom row features only a little less increase than the row above. There is a much larger difference in the middle and the bottom row for a strike-slip fault. This information tells us that for a normal fault, we are able to improve the inversion quality for the parameters of the deeper parts of the fault much better than in Scenario I by incorporating rotation ground motion data. The 6C inversion posterior pdf also seems to capture the true slip values quite well, especially on the sides and corners of the fault, compared to the results for a strike-slip fault. The total cumulative increase in information gain for all 24 subfaults is with a value of 32.13% more than three times higher than in Scenario I. The results for rupture velocity and rise time for a dip-slip fault are presented in Figure 8. Although we see a clear improvement between 3C and 6C posterior pdf's, both parameters are not captured correctly by the inversions. There is, nevertheless, an improvement 5 in resolving the true value with the 6C inversions. Especially for the rise time, the maximum of the posterior pdf is closer to the true value than for the 3C inversion.
The increase in information gain for rupture velocity (56.02%) and for rise time (28.47%) is also higher for Scenario II than for Scenario I.
The 6C inversions in the previous chapters are carried out with a manually predefined station selection.

10
The unequal distribution of energy ratios of rotation rate components might lead to the assumption that there is a higher information gain for the 6C inversion in Scenario II, because more receivers are selected north and south of the fault than in east and west direction. These receivers feature a higher energy ratio of the horizontal rotation components. To minimize this factor we run 200 6C inversions with 22 random stations each inversion, for both faulting types, and we show the results in Figure 9.
The results for a strike-slip fault are shown in the left column and for a dip-slip fault in the right column. We look on the slip 15 for the top, middle and bottom row individually to better assess the inversion quality dependence with depth for the different scenarios. Therefore there are 5 histograms for the information gain for the 200 6C inversions. We mark the information gain value for the 3C inversion for each scenario with a red vertical line in the histogram plots. Slip (in gray), rupture velocity (in blue) and rise time (in green) are all colored differently to distinguish between the three types of parameters.
(1) For the top layer of the subfaults, the information gain for a dip-slip fault (right) reaches the highest values, which is also 20 the result from the previous inversion. There are some receiver selections that even increase the information gain compared to the results from the inversion with regular receiver spacing.
(2) The middle subfault layer already features a smaller variance in the information gain distribution for both faulting types.
We also see a decrease in information gain increase, compared to the top layer, which is to be expected from the results before.
Additionally, there is a higher information gain increase for a dip-slip fault (right).

25
( 3) The deepest layer of subfaults shows the same information gain for 3C and 6C inversions for a strike-slip fault (left). The Figure shows that including rotational ground motion does not contribute to improving the results for the slip on the bottom row. This is different for a dip-slip fault (right). Every inversion increases the information gain for slip amplitudes in the bottom layer. We see some inversions in which the information gain doubles compared to the 3C inversion. This major difference in the two scenarios is one of the key aspects of our results. The ability to invert for the slip amplitudes at higher depth is something 30 that we only achieve with the incorporation of rotational ground motion rate data into the inversion. The statistical evaluation of 200 inversions with varying receiver consolidates this.
(4) The inversions for rupture velocity result in a higher information gain for a dip-slip fault, even with random station distribution. Although inversions for rupture velocity for a strike-slip event (Scenario I) already benefit from rotation rate data, the improvement is even higher for a dip-slip (Scenario II) event. We compare the results for both types of faulting in Figure 9 and see that for every row (all five parameters) the increase in 5 information gain between 3C and 6C inversions is higher for a pure dip-slip event. The random selection of the receivers used in the inversions does not undo this result.
For the dip-slip mechanism, the inversion quality depends much more on the spatial distribution of stations (larger spreading in Figure 9). This means, that much more care needs to be taken when choosing locations with respect to the fault geometry for rotational sensors compared to translational sensors.

Discussion
With the range of measured rotational ground motions between 10 −7 and 10 −4 rad/s we are able to draw practical conclusions from this synthetic study since effects in this range should be observable soon (Schreiber et al., 2009).
We show a summary of the improvement of the information gain between the different scenarios and parameters in percent in Table 1. Subfaults are analysed by the three horizontal layers (varying depth). The largest difference in increase in information 15 gain for the slip amplitudes is in the deepest subfault layers. Rupture velocity and rise time benefit from rotation rate data for both types of scenarios, but the improvement is larger for the dip-slip event. The results for slip, rupture velocity and rise time are also on the same order as the results in (Bernauer et al., 2014) for a similar strike-slip event.
Rotation data seem to be most beneficial for kinematic parameters of sources that fracture at least partly in vertical direction.
Additionally, the increase in information gain for rupture velocity and rise time suggests possible improvements in inversions 20 for other seismological research areas. These parameters are important for rheology and friction laws used in dynamic rupture simulations (Tinti et al., 2009). Directivity measurements could also benefit from high-quality rupture velocity inversions (Somerville et al., 1997). Rotation rate data improved the inversions for the rise time as a parameter of the source time function.
This helps when working with far-field ground displacements which are related to the derivative of the source time function.
A successful decrease in non-uniqueness in probabilistic finite source inversion suggests manufacturing portable rotation sen-25 sors with the mentioned sensitivities. The noise level, which was set to be equal to the noise level for translational ground motions sensors of 10%, should not exceed this to successfully improve inversion results (Bernauer et al., 2014).

Conclusion
Overall, non-uniqueness, which is an eminent problem in inversion theory, has been significantly reduced by the incorporation of synthetic rotational ground motion data.

30
The results from this study clearly suggest that seismological studies can benefit from the ability of seismometers to also measure rotational ground motions in addition to the conventional three translational components. This applies not only to finite source inversions, which are a crucial part of seismological research, but probably to all work that relies on data recorded by seismic receivers for treating/inverting real earthquakes and suffers problems such as non-uniqueness due to the non-linearity of the problem.
We successfully showed how rotation rate measurements can reduce the degree of non-uniqueness in probabilistic finite source 5 inversion. This is most prominent for an earthquake which also fractures in vertical direction, due to the higher energy in the horizontal rotation components. It is very likely that the vertical displacement gradient, included in these components, leads to the observed decrease in non-uniqueness since this information can not be obtained with translation recordings on the Earth's surface.
Acknowledgements. We would like to thank Lion Krischer for the support in professional code optimization. The authors acknowledge funding from the ERC advanced project ROMY (http://romy-erc.eu/). Further information necessary to replicate and build upon the reported research is available via correspondance by email.