# A wavelet-based detection and characterization of damped transient waves occurring in geophysical time-series: theory and application to the search for the translational oscillations of the inner core

Severine Rosat, Pascal Sailhac, Pascal Gegout

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

## ► To cite this version:

Severine Rosat, Pascal Sailhac, Pascal Gegout. A wavelet-based detection and characterization of damped transient waves occurring in geophysical time-series: theory and application to the search for the translational oscillations of the inner core. Geophysical Journal International, 2007, 171, pp.55-70.
10.1111/j.1365-246X.2007.03533.x. hal-00643195

HAL Id: hal-00643195

https://hal.science/hal-00643195v1

Submitted on 8 Nov 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_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# A wavelet-based detection and characterization of damped transient waves occurring in geophysical time-series: theory and application to the search for the translational oscillations of the inner core

Severine Rosat,<sup>1,2</sup> Pascal Sailhac<sup>2</sup> and Pascal Gegout<sup>2</sup>

Accepted 2007 June 21. Received 2007 March 28; in original form 2007 March 28

## SUMMARY

Geophysical time-series with non-stationary behaviour which contain transient excitations need specific tools to be analysed. Especially, it is difficult to detect and characterize transient excitations with attenuation by diffusion when their intensity is low compared to noise. This is a typical question in the search for the elusive translational motion of the inner core, the so-called Slichter triplet. This detection remains a non-trivial problem since the amplitude of the Slichter triplet is expected to be weak. Moreover, the characteristics (period, time of excitation and quality factor) of these low-intensity damped transient waves are poorly constrained. To improve the detection and characterization of such causal damped transient waves, we introduce a new continuous wavelet-based method using correlations with causal damped sinusoids that serve as the wavelet basis. This causal damped wavelet (CDW) method presents the advantage to give a direct estimate of the quality factor of the transient wave. After some synthetic tests, we present an application of the method to the time-varying gravity data recorded by superconducting gravimeters in the search for the Slichter triplet. No evidence of the Slichter triplet has been found in the 1-yr gravity data considered. However, we show that, thanks to the 3-D representation of our CDW method, it is possible to have a complete overview of the data set content and to constrain the detected events in terms of amplitude, quality factor, frequency and time.

**Key words:** causal damped transient waves, continuous wavelet transform, normal modes, Slichter triplet.

## 1 INTRODUCTION

Many time-series in geophysics exhibit non-stationary behaviours. While the series may contain dominant periodic signals, these signals can vary both in amplitude and in frequency over long periods of time. We consider the signal decomposition that consists of a sum of intermittent excitations with attenuation of the oscillations by diffusion. This kind of behaviour is present in signals found in seismology and electromagnetism for instance.

The example would be the gravity variations associated with the translational motion of the solid inner core, the so-called Slichter triplet (Slichter 1961). As the main restoring force is the Archimedean force, the period of the Slichter mode is directly linked to the value of the density jump at the inner core boundary (ICB). This density jump is $600 \text{ kg m}^{-3}$ for the PREM model (Dziewonski & Anderson 1981); while some recent estimations using normal-mode data (Masters & Gubbins 2003) or PKiKP/PcP amplitude ratios (Koper & Pyle 2004) have reached diverging conclusions for the density jump being greater or smaller, respectively, than that of PREM model. The observation of the Slichter modes would, there-

fore, constrain the density jump—as well as the viscosity—at the ICB. However, its detection has been a challenge since the first claim by Smylie (1992) of a probable, but widely disputed, observation based on time-varying gravity signals recorded by cryogenic gravimeters. Smylie *et al.* (1993) confirmed the previous detection using a product spectrum of four gravity records. The quest for the Slichter triplet began in fact with the venue of high resolution and low subseismic noise relative gravimeters based on the superconductivity. The records of these superconducting gravimeters (SGs) have been used by various authors to try to confirm Smylie’s claim, for instance, Jensen *et al.* (1995), Hinderer *et al.* (1995) using a product spectrum and cross-spectrum of SG data. Later, Courtier *et al.* (2000) proposed a stacking method, which enhances the harmonic degree-one modes and suggested a candidate for the Slichter modes very similar to the one detected by Smylie (1992). However, some more recent attempts by Rosat *et al.* (2003, 2004) and Guo *et al.* (2006, 2007) on less noisy SG data have not led to the confirmation of the previous candidate and recent theoretical developments (e.g. Rieutord 2002; Rogister 2003) concluded that Smylie’s observation was inconsistent with the predictions.

<footnote n="1"><sup>1</sup>Royal Observatory of Belgium, Brussels. E-mail: Severine.Rosat@oma.be</footnote>

<footnote n="2"><sup>2</sup>EOST (UMR CNRS-ULP 7516), Strasbourg, France</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Quality factors *Q*<sub>s</sub> of the Slichter modes as a function of its periods *T*<sub>s</sub> for three values of dynamic viscosity at the ICB. The upper curve is for a viscosity estimated in laboratory experiments and the lower curve corresponds to a viscosity estimated by Smylie & McMillan (2000) based on the claimed Slichter modes by Smylie (1992). The relation that links *Q*<sub>s</sub> to the period *T*<sub>s</sub> has been taken from the equations in Smylie & McMillan (1998). The limits for the Slichter periods as deduced from the estimated density jumps at the ICB are plotted in vertical dashed lines. The observed period by Smylie (1992) close to 3 hr is also represented. The synthetic cases we have considered are: (S1) *T*<sub>s</sub> = 4 hr and *Q*<sub>s</sub> = 10; (S2) *T*<sub>s</sub> = 9 hr and *Q*<sub>s</sub> = 10; (S3) *T*<sub>s</sub> = 4 hr and *Q*<sub>s</sub> = 1000; (S4) *T*<sub>s</sub> = 9 hr and *Q*<sub>s</sub> = 1000; (S5) *T*<sub>s</sub> = 6 hr and *Q*<sub>s</sub> = 10<sup>4</sup>.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Temporal and spectral representations of the Morlet wavelet with a frequency *f*<sub>0</sub>, respectively, from top to bottom of 0.1 and 1. The solid line corresponds to the real part and the imaginary part is plotted in dotted line.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Temporal and spectral representations of the mother wavelet CDW with a parameter $\alpha_0 = 0.1$ and a quality factor $Q_w$, respectively, from top to bottom of 10 and 60. The solid line corresponds to the real part and the imaginary part is plotted in dotted line.</figcaption></figure>

The main difficulty in the detection of these normal modes is that their mechanism of excitation is unknown. They could be randomly excited by some turbulent flow in the core or by a strong earthquake or by the atmosphere. Moreover, the damping rate and the frequency of the Slichter modes are subjected to debate. Crossley (1987), Crossley *et al.* (1991) and recently Rosat (2007) have studied the excitation of the Slichter modes by a strong earthquake. Their conclusion is that the surface amplitude of the Slichter modes after the past major seismic events is too small to be detected in time-varying gravity data recorded by SGs. Indeed, Rosat (2007) has demonstrated that a magnitude 9.7 vertical deep-slip seismic mechanism is required to excite the Slichter modes to a level detectable by SGs, that is, to the nanogal level (1 nGal = $10^{-2}$ nm s$^{-2}$). Other mechanisms of excitation have not yet been studied.

More studies have been done concerning the possible decay rate of the inner core free oscillation. Three sources of damping have been considered previously: the anelastic deformations of the inner core, the viscous dissipation and the magnetic dissipation. Crossley *et al.* (1991) have shown that the damping due to the seismic anelasticity of the inner core and mantle has a $Q$ value of the order of 5000 with a corresponding damping time of 400 d. The damping due to the outer core viscosity has been formulated by Smylie & McMillan (1998) and also by Rieutord (2002). The estimates of the dynamic viscosity of the outer core range from $1.6 \times 10^{-2}$ Pa s using laboratory experiments (Rutter *et al.* 2002), giving a $Q$ value of the order of $10^7$ to $1.2 \times 10^{11}$ Pa s (Smylie & McMillan 2000) using the claimed Slichter modes by Courtier *et al.* (2000), giving a $Q$ value less than 10 for a Slichter period of 5.5 hr. Smylie (1992) measured a quality factor of the resonance between 100 and 400, for a damping time of 6 d, and direct measurements on fig. 5 of Courtier *et al.* (2000) shows even higher values, while Eckman numbers found by Smylie & McMillan (2000) imply quality factors less

than 10 (Rieutord 2002). Mathews & Guo (2005) have proposed an upper limit of $1.7 \times 10^5$ Pa s using nutation data corresponding to a $Q$ value of 5000 (cf. also Guo *et al.* 2007). The magnetic damping of the inner core oscillation has been studied by Buffett & Goertz (1995). They have shown that the $Q$ value should be between $5.8 \times 10^5$ and 2200 for a magnetic field ranging from 0.0005 to 0.001 T corresponding to an e-folding time of the oscillation as long as 108 yr or as short as 150 d for a nominal period of 5.2 hr.

The damping rate of the Slichter modes is, therefore, poorly constrained. The period of the Slichter modes covers also a wide range of possibilities. The predicted PREM period is 5.42 hr for an ICB density jump of 600 kg m$^{-3}$, while observations of the density jump at the ICB give a value close to 4 hr for a density jump of 800 kg m$^{-3}$ (Masters & Gubbins 2003) or a value close to 9 hr for a density jump of 300 kg m$^{-3}$ (Koper & Pyle 2004). The claimed Slichter modes by Smylie (1992) gives a period closer to 3 hr while recent theoretical computations predict a period between 4 and 6 hr (Rieutord 2002, Rogister 2003). The predicted values for the quality factor and the period of the Slichter triplet are summarized in the graph Fig. 1.

The detectability of the Slichter modes depends on the magnitude of the response of the inner core to the geodynamic excitation process and on its decay rate. As the Slichter modes have not yet been undoubtedly detected, we can suppose that: (1) either they are not excited at all; (2) or they are continuously excited but their amplitudes are too small so that they are hidden by the noise (3) or they are quickly damped out after each excitation and in this case the usual spectral approaches are inefficient.

Lets follow the later hypothesis and assume that previous attempts to detect the Slichter modes based on spectral methods were unfruitful, because a possible strong attenuation by dissipation was neglected. Indeed, a frequency Fourier analysis reaches its limit in the case of a non-often repeated excitation and a strong
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Causal damped wavelet transform and damping identification method applied to the synthetic cases (S1) ($f_s = 7 \times 10^{-5}$ Hz, $Q_s = 10$), (S2) ($f_s = 3 \times 10^{-5}$ Hz, $Q_s = 10$), (S3) ($f_s = 7 \times 10^{-5}$ Hz, $Q_s = 1000$), (S4) ($f_s = 3 \times 10^{-5}$ Hz, $Q_s = 1000$) and (S5) ($f_s = 4.6 \times 10^{-5}$ Hz, $Q_s = 10^4$) excited at $t_s = 41.67$ d (synthetic cases in Fig. 3) and injected in a white noise (SNR = 10). From top to bottom are, respectively, plotted the time-series, the energy $|W_x|^2$ of the CDW transform, the phase, and $|W_x|^2$ versus ($t, Q_w$) at the estimated frequency $f_e$. For the synthetic case (S1), we have also plotted (b) $|W_x|^2$ versus (f, $Q_w$) at the estimated time $t_e = t_s = 41.7$ d and (c) $|W_x|^2$ and phase versus (t, f) with $Q_w = Q_s = 10$. The contour lines corresponding to the levels 99, 95 and 90 per cent of the maximum are plotted.</figcaption></figure>

attenuation of the oscillation. Time-frequency or wavelet representations (Morlet *et al.* 1982) are needed in that case for the detection and the characterization of the oscillation (time of excitation, frequency and decay rate). In acoustics, it is well known that damped free oscillations can be detected and characterized using the so-called asymptotic wavelet analysis for stationary phase approximation (see e.g. Saracco *et al.* 1990; Delprat *et al.* 1992; Carmona *et al.* 1997; Ruzzene *et al.* 1997; Staszewski 1997; Le & Argoul 2004). However, this technique is not effective for analysing very noisy data (Shibata

& Takagiwa 1996); furthermore the gravity variations concerned are not just damped free oscillations with time symmetries that resemble to Morlet wavelets but they can show some rapid jumps followed by damped sinusoids. Therefore, we propose to use a continuous wavelet method based upon a special mother wavelet built from a causal damped sinusoid with zero values at negative time. This new mother wavelet has a shape very similar to a causal exponentially damped oscillation, analogue to the so-called free induction decay (FID) signal in nuclear magnetic resonance (NMR) spectroscopy.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** (Continued.)</figcaption></figure>

In NMR, signals are consequences of strong artificial magnetic excitations and the stationary phase approximation on Morlet or Gabor representations works well (e.g. Antoine *et al.* 2000; Antoine & Coron 2001). To approach the analysis of very noisy signals, we rather consider families of causal damped wavelet (CDW). Thus the stationary phase approximation is not necessary. The signal is mapped into a 3-D representation of correlation coefficients versus time, frequency and quality factor, allowing to directly estimate the quality factor of the detected transient signals.

In the first part of this paper, the theory and the methodology to detect and characterize a damping oscillation are developed and illustrated using synthetic data. In a second part we show a first application of the CDW to actual SG data in the search for the Slichter triplet.

## 2 THEORY AND METHODOLOGY

Wavelet transforms enable to detect and characterize transient waves with damping that occur in a dispersive continuum. In this work, we consider that the transient events to be found are damped sinusoids; there are three target parameters to be estimated for each event: its starting time, its pulsation $\omega_0$ related to the frequency $f_0$, and its damping rate $\alpha_0$ related to the quality factor $Q = \omega_0/2\alpha_0$.

First, we briefly recall the basic theory (e.g. Morlet *et al.* 1982; Holschneider 1995; Mallat 1998) of the continuous wavelet transform and give the expression for the CDW; then we present an application of the CDW method to synthetic transient oscillations that mimic the Slichter mode. It is compared with the stationary phase approximation approach using Morlet wavelets. Finally, we give the
<figure><img src="image_8.png" /><figcaption>Figure 4. (Continued.)</figcaption></figure>

step by step procedure to use our CDW method in the detection and characterization of any damping oscillation.

## 2.1 Continuous wavelet transform

We define the continuous wavelet transform (CWT) of the temporal signal $x$ as the inner-product:

$$W_x(b, a) = \int_{-\infty}^{+\infty} x(t)h_{b,a}^*(t) dt, \qquad (1)$$

where $b$ and $a$ are the translation (time) and dilatation (scale) parameters, respectively. The function $h^*$ is the complex conjugate of $h$, and $h_{b,a}$ are the ‘daughters’ obtained by translation and dilatation of the mother wavelet $h_0$:

$$h_{b,a}(t) = \frac{1}{a} h_0 \left( \frac{t-b}{a} \right). \qquad (2)$$

To define a ‘mother’ wavelet $h_0$ which is translated and dilated to generate the wavelet basis, one has to check the admissibility condition which specifies that the normalization factor $C_g$, must be finite:

$$C_g = \int_0^{+\infty} \frac{|\hat{h}_0(\omega)|^2}{\omega} d\omega < \infty. \qquad (3)$$

In other words, to be admissible as a wavelet, a function must have zero mean and be localized both in time and in frequency space (Farge 1992). Different mother wavelets obey the admissibility condition (e.g. Ricker or Mexican hat, Cauchy wavelet). For historical reasons and comparisons with the windowed Fourier transform, one often uses the Morlet wavelet. It is symmetric in time and its damping is controlled by the width of a Gaussian envelope (Morlet et al. 1982; Grossmann et al. 1989). The expression for the Morlet wavelet
<figure><img src="image_9.png" /><figcaption>Figure 4. (Continued.)</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 5. Wavelet Transform of a damped sinusoid of frequency $f_s = 3 \times 10^{-5}$ Hz and with a quality factor $Q_s = 10$ (synthetic case S2), corresponding to a characteristic time of 1.23 d. The signal has been injected in a white noise so that the signal-to-noise ratio is 10. The wavelet analysis has been performed using (a) a CDW family with a quality factor $Q_w$ equal to 100, (b) a Morlet wavelet family with a frequency $10^{-4}$ Hz. The temporal signal, the phase and the energy $|W_x|^2$ of the wavelet transform, and the $Q_w$ determination graph, respectively, $\ln\{\frac{|W_x|\{b,a(b)\}|^2}{|h_0(\omega_e)|^2}\}$ (eq. B1) for Morlet wavelet and $|W_x|^2$ versus $(t, Q_w)$ at the frequency $f_e$ for CDW.</figcaption></figure>

in Fourier domain reads:

$$ \hat{h}_0(\omega) = e^{-(\omega-\omega_0)^2/2} + \eta(\omega), \quad (4) $$

where $\eta$ is a small corrective term that guarantees the integral in eq. (3) converges and that the expression (4) vanishes for $\omega \le 0$. Thus the Morlet wavelet obeys the admissibility condition (Holschneider 1995). Time and frequency representations of a Morlet wavelet of pulsation $\omega_0 = 2\pi$ are shown in Fig. 2.

## 2.2 The causal damped wavelet

Since the signal to be detected, the Slichter modes, may have a shape more similar to a causal damped oscillation than to the function used in Gabor or Morlet analysis, let us consider correlations with a new wavelet family based on the causal damped oscillation representation. Causal damped oscillations can be evaluated
<figure><img src="image_11.png" /><figcaption>Figure 6. Error of the estimated quality factor $Q_e$ with respect to $Q_s$ as a function of the temporal ridge limits $t_i$ and $t_f$ used to fit a slope in the Morlet damping identification method (eq. B1).</figcaption></figure>

using the following Lorentz function in the Fourier domain:

$$ \hat{g}_0(\omega) = \frac{1}{2[\alpha_0 + i(\omega - \omega_0)]}, \quad (5) $$

where $\omega_0 = 2\pi f_0$ is the pulsation related to the frequency $f_0$ of the damped oscillation, and $\alpha_0 = \omega_0/2Q$ is the damping rate with the corresponding quality factor $Q$.

We define a continuous wavelet representation in a way similar to that classically used for the Morlet wavelet (Grossmann et al. 1989). Instead of the sinusoid with a Gaussian window, we consider a causal sinusoid with an exponential decay in order to obtain a basis dedicated to the representation and detection of causal damped sinusoids.

For the previous function (in eq. 5), the normalizing factor $C_g$ is defined by

$$ C_g = \int_0^{+\infty} \frac{1}{\alpha_0^2 \omega + \omega(\omega - \omega_0)^2} d\omega. \quad (6) $$

Since the integral in the above expression is diverging for $\omega \to 0$ (and converging to zero for $\omega \to +\infty$), $C_g$ is infinite; so the function $g_0$ is not admissible. This lack of admissibility of $g_0$ is due to the presence of long-period terms. In order to use the causal damping as a mother wavelet, a function with fast decay in period is added. That supplementary term is similar to the small corrective term $\eta$ of the Morlet wavelet (eq. 4); it removes the low-frequency components causing the divergence of $C_g$.

Thus, we consider another mother wavelet $h_0$ expressed in the Fourier domain by:

$$ \hat{h}_0(\omega) = \frac{1}{\alpha_0 + i(\omega - \omega_0)} - \frac{1}{\alpha_0 - i\omega_0(1 + \frac{\omega}{\omega_1})}. \quad (7) $$

where $\alpha_0 = \omega_0/2Q$ with $\omega_0$ the mother pulsation and $Q$ the quality factor of the corresponding wave. The pulsation $\omega_1$ is a normalizing factor which controls the amplitude and quality factor of the corrective term that makes the wavelet to be admissible. Taking the inverse

Fourier transform in eq. (7) leads to the following expression for the CDW:

$$ h_0(t) = H(t) \exp(-\alpha_0 + i\omega_0)t + \frac{\omega_1}{\omega_0} H(-t) \exp\left(\frac{\alpha_0\omega_1}{\omega_0} - i\omega_0\right) t, \quad (8) $$

where $H$ refers to the Heaviside step function. The temporal and spectral representations of the CDW $h_0$ are plotted in Fig. 3 for different values of $Q$ and $\omega_0$. In these plots and the applications hereafter, the pulsation $\omega_1$ has been fixed to 1 Hz and the pulsation $\omega_0$ to $10^{-6}$ Hz. This choice of $\omega_0$ much smaller than $\omega_1$ leads to a correcting term that damps much faster as the quality factor $Q_1$ associated with the correcting term is related to the wavelet quality factor $Q$ via the relation $Q_1 = \frac{\omega_0}{\omega_1}Q$ and its maximum amplitude is equal to the ratio $\omega_1/\omega_0$.

The corresponding normalization factor $C_h$ is finite, and is equal to (see Appendix A):

$$ C_h = \left( \frac{2Q}{\omega_0} \right)^2 \frac{\left( \frac{\omega_0^2}{\omega_1^2} - 1 \right) \ln\left( \frac{\omega_0}{\omega_1} \right) + 4\left( 1 + \frac{\omega_0}{\omega_1} \right)^2 Q \tan^{-1}(2Q)}{\left( 1 - \frac{\omega_0}{\omega_1} \right)^2 + 4Q^2 \left( 1 + \frac{\omega_0}{\omega_1} \right)^2}. \quad (9) $$

A family of CDW is defined for a given $Q = Q_w$, so we can consider many families for various $Q_w$. Indeed, dilation and translation does not change the quality factor of the damped sinusoid wavelets. In the following we will see, thanks to numerical simulations, that the choice of $Q_w$ is primordial for the analysis. In particular, we will see that the correlation between the CDW and a damped oscillating signal is the largest when $Q_w$ is equal to the quality factor $Q_e$ of the signal. Hence we can estimate directly the quality factor of the signal from maxima of the CDW transform.
<figure><img src="image_12.png" /><figcaption>**Figure 7.** CDW analysis with $Q_w = 100$ of the stacked 1-min time-varying gravity residuals at the five SG sites Cantley, Canberra, Matsushiro, Strasbourg and Sutherland for the year 2001. The temporal signal and the FFT amplitude spectra in the subseismic frequency range (0.03–0.08) mHz and in the seismic frequency band (0.4–0.7) mHz, are, respectively, plotted for (a) the equatorial prograde, (b) the axial and (c) the equatorial retrograde sequences resulting from the stacking process.</figcaption></figure>

## 2.3 Causal damped wavelet analysis of synthetic Slichter signals

The analysis of synthetic signals is necessary to understand the behaviour of the wavelet decomposition on the basis of damped oscillations: particularly this helps to understand the influence of the quality factor $Q_w$ of the wavelet. We consider synthetic data sets that mimic the translational motion of the inner core. These events are considered to obey the following shape: $x(t) = H(t - t_0)A(t - t_0)\cos\omega_s(t - t_0)$, where $H$ is the Heaviside step function and $A$ is the amplitude with exponential decay related to the pulsation $\omega_s$ via the quality factor $Q_s$, and $t_s$ is the starting time of the transient signal.

In the following, we analyse a damped sinusoid excited once at the time $t_s = 41.67$ d (1000 hr) at the frequency $f_s$ with a quality factor $Q_s$, a damping time $\tau_s$ corresponding to the five synthetic cases indicated in Fig. 1. These synthetic cases are: (S1) $f_s = 7 \times 10^{-5}$ Hz ($T_s = 4$ hr), $Q_s = 10$ ($\tau_s = 0.53$ d); (S2) $f_s = 3 \times 10^{-5}$ Hz ($T_s = 9$ hr), $Q_s = 10$ ($\tau_s = 1.23$ d); (S3) $f_s = 7 \times 10^{-5}$ Hz ($T_s = 4$ hr), $Q_s = 1000$ ($\tau_s = 53$ d); (S4) $f_s = 3 \times 10^{-5}$ Hz ($T_s = 9$ hr), $Q_s = 1000$ ($\tau_s = 123$ d) and (S5) $f_s = 4.6 \times 10^{-5}$ Hz ($T_s = 6$ hr), $Q_s = 10^4$ ($\tau_s = 801$ d). These values of frequency and quality factor sample the possible parameters predicted for the Slichter modes based on theoretical computations and on indirect observations of the density jump at the ICB (cf. Section 1). We do not consider the

cases of higher $Q$ values, as it would correspond to non-damping oscillations in the timescales of actual gravity records. The transient signal has been injected in a white noise so that the signal-to-noise ratio is 10. This large amplitude is used to better illustrate our CDW method, while, in actual records, the Slichter triplet's amplitude is expected to be hardly above the noise level.

The CDW analysis results of the five synthetic time-series (S1)–(S5) are plotted in Fig. 4. For every case, from top to bottom we have represented the temporal signal, the energy $|W_x|^2$ and the phase of the CDW transform, and $|W_x|^2$ versus $(t, Q_w)$ at the estimated frequency $f_e$ of the local maximum of $|W_x|^2$.

Note that we have considered the case of a single transient wave, while the Slichter mode in a realistic rotating elliptical Earth model should be modelled by three transient oscillations (hence so-called Slichter triplet), which are largely split in frequency.

## 2.4 Damping identification with the causal damped wavelet

The complete characterization of a damped transient wave requires the estimation of the quality factor of the geophysical process associated with the excitation and then damping of the wave. Intuitively, we can expect that, when the parameter $Q_w$ of the CDW is equal to the quality factor $Q_s$ of the signal analysed, the correla-
<figure><img src="image_13.png" /><figcaption>**Figure 8.** CDW transform [||*W*<sub>*x*</sub>||<sup>2</sup> versus (*t*,*Q*<sub>w</sub>)] of the axial series plotted in Fig. 7 (b), (a) around the seismic mode 0S<sub>3</sub> and (b) around the seismic mode 0S<sub>4</sub> with *Q*<sub>w</sub> = 500.</figcaption></figure>

tion between the CDW and the signal is stronger. So the idea is to compute the maximum of the energy |*W*<sub>*x*</sub>|^2 of the CDW transform at the local maximum for different *Q*<sub>w</sub>, that is, for different CDW families, and the *Q*<sub>w</sub> for the maximum of these maxima give the estimated value *Q*<sub>e</sub> of *Q*<sub>s</sub>. The energy |*W*<sub>*x*</sub>|^2 for the previous synthetic cases (S1)–(S5) are plotted versus the wavelet quality factor *Q*<sub>w</sub> and the time *t* at *f*<sub>e</sub> in Fig. 4 (lowest graphs). Similar maps of |*W*<sub>*x*</sub>|^2 versus (*f*,*Q*<sub>w</sub>) can be plotted at the estimated time *t*<sub>e</sub> (example for the synthetic case S1 in Fig. 4b). The energy |*W*<sub>*x*</sub>|^2 and the phase of the CDW transform of the synthetic case (S1) using *Q*<sub>w</sub> = *Q*<sub>s</sub> = 10 are represented in Fig. 4(c). Note the much higher values of |*W*<sub>*x*</sub>|^2 in this case. From the *Q*<sub>w</sub>−*t* maps of |*W*<sub>*x*</sub>|^2, we obtain a good approximation of the quality factor of the transient wave for the five cases (S1)–(S5) considered.

This damping determination method is only applicable to the CDW. Some studies have already treated the problem of damping measurement in timescale using Gabor or Morlet decompositions (e.g. Delprat et al. 1992; Staszewski 1997; Slavic et al. 2003). The method is described in Appendix B. It works also with the CDW but the introduction of the CDW family aims in fact to estimate the damping rate of the signal directly. We will see in the next section that the CDW damping estimation method is more stable in the

presence of noise than the damping identification method used with the Morlet wavelet, for instance.

## 2.5 Comparison with the Morlet wavelet

For a comparison between the Morlet wavelet and the CDW, we use the synthetic signal (S2), that is to say a damping signal excited at the time $t_s = 1000$ hr = 41.67 d oscillating at the frequency $f_s = 3 \times 10^{-5}$ Hz with a quality factor $Q_s = 10$. The signal-to-noise ratio is 10. The Morlet wavelet and the CDW transforms are plotted in Fig. 5 with, from top to bottom, the temporal signal, the phase, the energy $|W_x|^2$, the logarithm of $(|W_x|/|\hat{n}_0|)$ (eq. B1) at the frequency $f_e = 3 \times 10^{-5}$ Hz for the Morlet wavelet (based on the damping identification method described in Appendix B) and $|W_x|^2$ versus $(t, Q_w)$ for CDW at the frequency $f_e = 3 \times 10^{-5}$ Hz. The damping identification method with Morlet gives a good estimate of $Q_s$, however the choice of the ridge (local maximum) limits $t_i$ and $t_f$ in time can perturb considerably the estimation of $Q_s$ as the slope does not clearly appear. To illustrate this, we vary the ridge limits $t_i$ and $t_f$ and compute, for different pairs $(t_i, t_f)$, the error
<figure><img src="image_14.png" /><figcaption>**Figure 9.** $|W_x|^2$ versus $(t, Q_w)$ at the estimated frequency $f_e$ (a) for the three detected singlets $\alpha$, $\beta$, and $\gamma$ of $0S_3$, (b) for the detected singlet $\delta$ of $0S_4$, corresponding to the maxima seen in Fig. 8(a) and (b).</figcaption></figure>

$|Q_s - Q_e|/Q_s$. The resulting error map is plotted in Fig. 6. Note that the error can be larger than 100 per cent.

Finally, the damping identification method described in Appendix B can give very large errors, because small errors on the ridge limits result in changes in the slope, which lead to values of $Q_e$ largely different. So the damping identification method using CDW enables a more direct and more stable estimation of $Q_s$ in the noisy case.

## 2.6 Methodology for the causal damped wavelet application

As we have introduced a parameter $Q$ in the wavelet function, we can define various wavelet families for different values of $Q = Q_w$. We propose a methodology for the characterization of any damping signal, with frequency $f_s$ and quality factor $Q_s$ using the CDW. However, before analysing the wavelet transform to extract the parameters of the damped oscillation, the signal must first be detected. The detection step can be performed by any time-frequency or wavelet representation: the CDW family is not superior to the Morlet wavelet for instance for the detection. In the case of using the CDW, the choice of $Q_w$ for the detection is arbitrary; however a better frequency resolution is obtained for higher $Q_w$. Once the presence of a harmonic signal has been identified, the CDW can be applied to the data set limited in the time and frequency ranges that surround the detected signal. In the first stage, $Q_w$ is arbitrary—except if we already have some *a priori* knowledge of the decay rate of the signal detected. From the energy $|W_x|^2$ (squared module) and the phase, the time of excitation and the frequency can be estimated. The max-

imum of the squared modulus of the CDW transform $|W_x|^2$ is found at a time $t$ equal to $t_e$ and a frequency $f$ equal to $f_e$. A map of $|W_x|^2$ can be plotted as a function of $Q_w$ and $t$ at the frequency $f_e$. From the maximum of this map, the quality factor $Q_e$ can be evaluated. By iteration on $Q_w$ and based on the maximum of correlation of the CDW and the signal at the estimated frequency $f_e$, a more precise value of $Q_e$ can be obtained.

We have developed a new CDW family that works well for the detection and characterization of damped transient waves that mimic the Slichter mode. In the next section, we will apply our CDW analysis method to actual time-varying SG gravity data in the search for the surface gravity effect of the Slichter triplet.

## 3 APPLICATION TO ACTUAL GRAVITY RECORDS

The previous searches for the Slichter modes were based on different stacking and spectral methods applied to long SG records (usually longer than 1 yr), which implicitly supposed that the Slichter signal was continuously present or statistically well represented in the data series, equivalently meaning that its quality factor $Q$ was very large or if small, its excitation was often repeated. As these precedent studies have not led to the undisputable detection of the Slichter triplet, we must consider the case that the Slichter modes have a small $Q$ (of the order of hundred or less), as proposed by Smylie (1992).

In view of all the possible periods for the Slichter modes, the CDW method must be applied on a wide frequency range. In the following,
<figure><img src="image_15.png" /><figcaption>**Figure 10.** $|W_x|^2$ versus $(t,f)$ with $Q_w = 100$ in the subseismic frequency band (0.05–0.08 mHz) of the three sequences (a) equatorial prograde, (b) axial and (c) equatorial retrograde, resulting from the stacked 1-min time-varying gravity residuals at the five SG sites Cantley, Canberra, Matsushiro, Strasbourg and Sutherland for the year 2001.</figcaption></figure>

we consider that the decay rate of the Slichter modes is fast (small $Q$). We suppose also that the Slichter modes have been excited to a level detectable with SGs, meaning with amplitude larger than $10^{-2}$ nm s$^{-2}$ (1 nGal).

We apply the CDW analysis method to the same data set used in Rosat et al. (2004, 2006) in their searches for the Slichter modes using stacking and scanning methods. The data set are the stacked 1-min time-varying gravity residuals obtained from the records at the five SG sites Cantley (Canada), Canberra (Australia), Matsushiro (Japan), Strasbourg (France) and Sutherland (South Africa) during the year 2001. The choice of these stations was based on their low-noise level in the frequency range considered (periods between 1h and 8h) and their geographical distribution. The gravity records were corrected for the polar motion effect, for the local tides and for the local atmospheric pressure effect via a local barometric admittance, in order to reduce the noise level. The residuals obtained were band-passed filtered with cut-off periods of 15 min and 12 hr to remove any remaining long-period signal as well as high-frequency disturbances. As the atmospheric pressure effect was corrected using a fitted admittance, the diurnal thermal waves were still present in the data. To remove these diurnal harmonics, a Lorentz function was fitted to the spectral peaks and these solar tides resonances were then put to the background level using the same method as Smylie

et al. (1993) (Rosat et al. 2006). The same process was applied to the tidal waves remaining in the residuals, as well as to the waves identified as being non-linear oceanic tides. These gravity residuals were then stacked. The stacking method used was the one proposed by Courtier et al. (2000) that takes advantage of the harmonic degree-one pattern of the Slichter modes and leads to three sequences called equatorial prograde, axial and equatorial retrograde corresponding to the three oscillations representing the Slichter modes in a realistic rotating non-spherical Earth. In a worry to save some computation time, we decimate the stacked data to 10 min. That reduces considerably the size of our data sets without disturbing the signals in the frequency band considered.

We analyse with the CDW the three resulting sequences: equatorial prograde, axial and equatorial retrograde (cf. Rosat et al. 2004, 2006) which should, respectively, enhance the $m = -1, 0$ and 1 singlet of the Slichter modes. The temporal signals and the corresponding FFT amplitude spectra, in the subseismic and in the seismic frequency bands, are plotted in Fig. 7. We perform the CDW analysis of these three series with $Q_w = 100$ in two frequency ranges: (1) the seismic band (0.3–0.8) mHz, where we detect two events in the band (0.45–0.66) mHz (Fig. 7b); (2) the subseismic band (0.03–0.08) mHz, where the Slichter modes frequencies are predicted, but we focus in the frequency band (0.05–0.08) mHz to avoid the
<figure><img src="image_16.png" /><figcaption>**Figure 11.** 3-D-scatter representation and projection in two planes of the local maxima of $|W_x|^2$ as a function of the frequency, $Q_w$ and time for the three stacked sequences. The colour of the scatters is proportional to the amplitude $|W_x|^2$. These local maxima have been selected according to a threshold of 60 per cent of the global maximum of $|W_x|^2$. Note that Fig. 10 is a continuous representation and a projection of this 3-D representation in the plane $Q_w = 100$ for a threshold of 0 per cent.</figcaption></figure>

quart-diurnal tidal band and because this band contains both the previous claim by Courtier *et al.* (2000) and the possible candidates detected by Rosat *et al.* (2006).

## 3.1 CDW analysis in the seismic band (0.45–0.66) mHz

Two main maxima of correlation around the day 175 and close to the frequencies 0.5 and 0.65 mHz are detected. They correspond to the well-known seismic normal modes $_0S_3$ and $_0S_4$ that have been excited by the $M_w = 8.4$ Peruvian earthquake on 2001 June 23, respectively, at the frequencies 0.47 and 0.65 mHz with PREM quality factors of 417.5 and 373.2 (Masters & Widmer 1995).

To better distinguish these detected signals, we perform two CDW analyses in smaller frequency ranges, respectively, between 0.45 and 0.49 mHz (Fig. 8a) for $_0S_3$ and between the frequencies 0.62 and 0.66 mHz for $_0S_4$ (Fig. 8b). Because of the Earth's rotation and ellipticity, the seismic normal modes are split into $2l + 1$ singlets, where $l$ is the spherical harmonic degree and $m$ the harmonic order taking the values $-l$ to $+l$. In order to distinguish the singlets of $_0S_3$ and the singlets of $_0S_4$, we have used a CDW family of $Q_w = 500$ in the Fig. 8. Taking smaller $Q_w$-values does not resolve the singlets in frequency. For the three singlets of $_0S_3$ (noted $\alpha, \beta$ and $\gamma$ in the figures) and one singlet of $_0S_4$ (noted $\delta$ in the figures) that appear in our stacked SG data, we have plotted $|W_x|^2$ versus $(t, Q_w)$ in Fig. 9. The quality factors for $_0S_3$ are underestimated (Fig. 9a). The reason may be the small amplitude of $_0S_3$ compared to $_0S_4$, the

poorly resolved singlets of $_0S_3$ and the stacking process that aims at enhancing the degree one above other harmonic degree modes. For $_0S_4$, we obtain $f_e = 0.06432$ mHz and $Q_e = 495 \pm 90$ for the singlet $m = 0$. The uncertainties on $Q_e$ are large but they could be considerably reduced by concentrating the CDW analysis in the seismic frequency band and by using 1-min sampled data instead of 10-min sampled data. As the purpose of this work is not to analyse the seismic normal modes, we will not develop further their CDW analysis here. Besides, other methods are more appropriate.

## 3.2 CDW analysis in the subseismic band (0.05–0.08) mHz

We have seen that the CDW performs well for the detection of the seismic modes present in the SG data. We perform again the CDW analysis on the three sequences of Fig. 7, but we focus this time the analysis in the subseismic frequency band, between 0.05 and 0.08 mHz. This frequency band excludes the quart-diurnal tidal group where many non-linear tides are imperfectly modelled. Note that $|W_x|^2$ is smaller than in the seismic band (Fig. 10). Some events appear in the CDW transform of the axial sequence. They are due to the presence of an unexplained spike in the stacked time-series.

If we compute the CDW transform of the three sequences for different $Q_w$-CDW families, we can make a 3-D scatter representations of the energy $|W_x|^2$ versus $(tf, Q_w)$, where for clarity we select only the scatters corresponding to an amplitude larger than 60 per cent of the global maximum of $|W_x|^2$. These 3-D scatter plots are
<figure><img src="image_17.png" /><figcaption>**Figure 11.** (Continued.)</figcaption></figure>

represented in Fig. 11. Note that no event with $Q_e > 40$ has been detected in the frequency range 0.05–0.08 mHz. The most significant events (with the largest correlation with the CDW) correspond to the dark-red scatters where $Q_e \approx 5$, at the day 250. They correspond to the spike present in the data at that time and which put some energy at almost all the frequencies. No other significant event has been detected.

The peaks attributed to the Slichter triplet by Courtier *et al.* (2000) do not appear in our data set, as well as the predicted frequencies of the Slichter triplet (e.g. Rieutord 2002; Rogister 2003). The candidates detected by Guo *et al.* (2007) that have passed their test for at least two data sets have not been detected by our CDW analysis. The possible candidates detected by Rosat *et al.* (2006) based on the same stacked data sets have not been detected by our CDW method. So we may conclude that they were much probably a noise effect, or due to the presence of the spike, the origin of which is unknown.

Finally, we have not found any clear evidence of the Slichter triplet in the CDW transform of the stacked 1-yr SG records, which would undoubtedly confirm their detection and the previous claim. However, our CDW method enables to put some constraints on the possible events by associating them with a physical parameter (the quality factor) associated to an excitation process. Hence it is possible to distinguish a realistic damped harmonic oscillation from a much probable noise artefact.

# 4 CONCLUSION AND PERSPECTIVES

We have developed and tested a wavelet-based method to detect and characterize damped transient waves occurring in geophysical

time-series. The wavelet analysis is constrained by a parameter $Q_w$ defining a CDW family. The quality factor of the geophysical wave is directly determined by the maximum of correlation with various CDW families.

This new wavelet-based technique has been designed for the detection and characterization of damped transient waves; it can be an efficient way to detect small amplitude signals like those caused by the Slichter translational motion of the inner core, which is one of the largest challenges remaining in fundamental geophysics. The search for the free oscillation of the inner core using this wavelet-based method has not led yet to their detection and has not confirmed the previous claimed observations. However, the CDW method constrains the detection in amplitude, time, frequency and quality factor, associated to a physical excitation process.

Nevertheless, further work is needed to detect the Slichter modes, particularly the application of the CDW method to more data sets. Similarly, we suggest that this wavelet method would be useful in the detection of other transient small amplitude electromagnetic and seismic signals associated to tectonic or volcanic activity.

# ACKNOWLEDGMENTS

This contribution was partly supported by CNRS (Intérieur de la Terre), the Japan Society for the Promotion of Science (JSPS) and PRODEX (Contract of the Belgian Federal Science Policy). The authors would like to thank the GGP station managers for their data sets.
<figure><img src="image_18.png" /><figcaption>**Figure 11.** (Continued.)</figcaption></figure>

## REFERENCES

Antoine, J.P. & Coron, A., 2001. Time-frequency and time-scale approach to magnetic resonance spectroscopy, *J. Comput. Methods Sci. Eng. (JCMSE)*, **1**, 327–352.

Antoine, J.P., Coron, A. & Dereppe, J.M., 2000. Water peak suppression: time-frequency vs. time-scale approach, *J. Magnetic Resonance*, **144**, 189–194.

Buffett, B.A. & Goertz, D.E., 1995. Magnetic damping of the translational oscillations of the inner core, *Geophys. J. Int.*, **120**, 103–110.

Carmona, R.A., Hwang, W.L. & Torresani, B., 1997. Characterization of signals by the ridges of their wavelet transform, *IEEE Trans. Signal Proc.*, **45**, 2586–2590.

Courtier, N. et al., 2000. Global superconducting gravimeter observations and the search for the translational modes of the inner core, *Phys. Earth Planet. Int.*, **117**, 3–20.

Crossley, D., 1987. The excitation of core modes by Earthquakes, in *Structure and Dynamics of Earth's Deep Interior, Geophys. Monogr., 46/IUGG Ser.*, Vol. 1, pp. 41–50, eds D.E. Smylie & R. Hide, AGU, Washington, DC.

Crossley, D.J., 1992. Eigensolutions and seismic excitation of the Slichter mode triplet for a fully rotating Earth model, *EOS* **73**(43), 60.

Crossley, D., Hinderer, J. & Legros, H., 1991. On the excitation, detection and damping of core modes, *Phys. Earth Planet. Int.*, **68**, 97–116.

Dahlen, F.A. & Sailor, R.V., 1979. Rotational and elliptical splitting of the free oscillations of the Earth, *Geophys. J. R. Astron. Soc.*, **58**, 609–623.

Delprat, N., Escudié, B., Guillemain, P., Kronland-Martinet, R., Tchamitchian, Ph. & Torrésani, B., 1992. Asymptotic wavelet and Gabor analysis: extraction of instantaneous frequencies, *IEEE Trans. Inform. Theory*, **38**, 644–664.

Dziewonski, A.M. & Anderson D.L., 1981. Preliminary Reference Earth Model (PREM), *Phys. Earth Planet. Int.*, **25**, 297–356.

Farge, M., 1992. Wavelet transforms and their applications to turbulence, *Annu. Rev. Fluid. Mech.*, **24**, 395–457.

Gilbert, F. & Dziewonski, A., 1975. An application of normal mode theory to the retrieval of structural parameters and source mechanisms from seismic spectra. *Phil. Trans. R. Soc. Lond.*, **A278**, 187–269.

Grossmann, A., Kronland-Martinet, R. & Morlet, J., 1989. Reading and understanding continuous wavelet transforms, in *Wavelets: Time-Frequency Methods and Phase Space*, pp. 2–20, eds J. Combes, A. Grossmann & Ph. Tchamitchian, Springer-Verlag, Marseille, France.

Guo, J.Y., Dierks, O., Neumeyer, J. & Shum, C.K., 2006. Weighting algorithms to stack superconducting gravimeter data for the potential detection of the Slichter modes, *J. Geodyn.*, **41**, 326–333.

Guo, J.Y., Dierks, O., Neumeyer, J. & Shum, C.K., 2007. A search for the Slichter modes in superconducting gravimeter records using a new method, *Geophys. J. Int.*, **168**, 507–517.

Hinderer, J., Crossley, D. & Jensen, O., 1995. A search for the Slichter triplet in superconducting gravimeter data, *Phys. Earth Planet. Int.*, **90**, 183–195.

Holschneider, M., 1995. *Wavelets, an Analysis Tool*, Clarendon Press, Oxford, 423p.

Jensen, O.G., Hinderer, J. & Crossley, D.J., 1995. Noise limitations in the core-mode band of superconducting gravimeter data, *Phys. Earth Planet. Int.*, **90**, 169–181.

Koper, K.D. & Pyle, M.L., 2004. Observations of PKiKP/PcP amplitude ratios and implications for Earth structure at the boundaries of the liquid core, *J. Geophys. Res.*, **109**, B03301.

Le, T.-P. & Argoul, P., 2004. Continuous wavelet transform for modal identification using free decay response, *J. Sound Vib.*, **277**, 73–100.

Mallat, S., 1998. *A Wavelet Tour of Signal Processing*, Academic Press, San Diego, 637p.

Masters, G. & Gubbins, D., 2003. On the resolution of the density within the Earth, *Phys. Earth Planet. Int.*, **140**, 159–167.

Masters, T. G., Widmer, R., 1995. Free oscillations: frequencies and attenuations, in *Global Earth Physics, A Handbook of Physical Constants*, pp. 104–125, Am. Geophys. Union, Washington DC.
Mathews, P.M. & Guo, J.Y., 2005. Visco-electromagnetic coupling in precession-nutation theory, *J. Geophys. Res.*, **110**(B2), B02402. doi: 10.1029/2003JB002915.

Morlet, J., Arens, G., Fourgeau, E. & Giard, D., 1982. Wave propagation and sampling theory part I and II, *Geophysics*, **47**, 203–236.

Rieutord, M., 2002. Slichter modes of the Earth revisited, *Phys. Earth Planet. Int.*, **131**, 269–278.

Rogister, Y., 2003. Splitting of seismic-free oscillations and of the Slichter triplet using the normal mode theory of a rotating, ellipsoidal earth, *Phys. Earth Planet. Int.*, **140**, 169–182.

Rosat, S., 2007. Optimal seismic source mechanisms to excite the Slichter mode, *J. Geodesy*, in press.

Rosat, S., Hinderer, J., Crossley, D. & Rivera, L., 2003. The search for the Slichter mode: comparison of noise levels of superconducting gravimeters and investigation of a stacking method, *Phys. Earth Planet. Int.*, **140**, 183–202.

Rosat, S., Hinderer, J., Crossley, D. & Boy, J.P., 2004. Performance of superconducting gravimeters from long-period seismology to tides, *J. Geodyn.*, **38**(3–5), 461–476.

Rosat, S., Rogister, Y., Crossley, D. & Hinderer, J., 2006. A search for the Slichter Triplet with Superconducting Gravimeters: impact of the density jump at the inner core boundary, *J. Geodyn.*, **41**, 296–306.

Rutter, M.D., Secco, R.A., Uchida, T., Liu, H., Wang, Y., Rivers, M.L. & Sutton, S.R., 2002. Towards evaluating the viscosity of the Earth's outer core: an experimental high pressure study of liquid Fe-S (8.5 wt.% S), *Geophys. Res. Lett.*, **29**(8), 1217. doi: 10.1029/2001GL014392.

Ruzzene, M., Fasana, A., Garibaldi, L. & Piombo, B., 1997. Natural frequencies and dampings identification using wavelet transform: application to real data, *Mech. Sys. Signal Proc.*, **11**, 207–218.

Saracco, G., Guillemain, P. & Kronland-Martinet, R., 1990. Characterization of elastic shells by the use of the wavelet transform, *IEEE-Ultrasonics*, **2**, 881–885.

Shibata, R. & Takagiwa, M., 1996. Consistency of frequency estimates based on the wavelet transform, *J. Time Series*, **18**, 641–662.

Slavic, J., Simonovski, I. & Boltezar, M., 2003. Damping identification using a continuous wavelet transform: application to real data, *J. Sound Vibrat.*, **262**(2), 291–307.

Slichter, L.B., 1961. The fundamental free mode of the Earth's inner core, *Proc. Natl. Acad. Sci., USA*, **47**, 186–190.

Smylie, D.E., 1992. The inner core translational triplet and the density near earth's center, *Science*, **255**, 1678–1682.

Smylie, D.E. & McMillan, D.G., 1998. Viscous and rotational splitting of the translational oscillations of Earth's solid inner core, *Phys. Earth Planet. Int.*, **106**, 1–18.

Smylie, D.E. & McMillan, D.G., 2000. The inner core as a dynamic viscometer, *Phys. Earth Planet. Int.*, **117**, 71–79.

Smylie, D.E., Hinderer, J., Richter, B. & Ducarme, B., 1993. The product spectra of gravity and barometric pressure in Europe, *Phys. Earth Planet. Int.*, **80**, 135–157.

Staszewski, W.J., 1997. Identification of damping in MDOF systems using time-scale decomposition, *J. Sound Vibrat.*, **203**, 283–305.

Tchamitchian, P. & Torresani, B., 1992. Ridge and skeleton extraction from the wavelet transform, in *Wavelets and Their Applications*, pp. 123–151, ed M.B. Ruskai, Jones and Bartlett, Boston.

## APPENDIX A: COMPUTATION OF THE NORMALIZATION FACTOR

The normalization factor is defined as:

$$ Ch = \int_{0}^{+\infty} \frac{|\hat{h}_{0}(\omega)|^2}{\omega} d\omega = \int_{0}^{+\infty} \frac{\left| \frac{1}{\omega_0 + i(\omega - \omega_0)} - \frac{1}{\omega_0 - i\omega_0(1 + \frac{\omega}{\omega_1})} \right|^2}{\omega} d\omega \\ = \int_{0}^{+\infty} \frac{\omega(1 + \frac{\omega_0}{\omega_1})^2}{[\alpha_0^2 + (\omega - \omega_0)^2][\alpha_0^2 + \omega_0^2(1 + \frac{\omega}{\omega_1})^2]} $$

For the calculation of this integral, we have used a decomposition of rational functions into a sum of partial fractions then we have integrated. However, since there is problems of convergence at both the upper and lower limits of the integral, we rather decided to perform the integration between $\epsilon$ and $1/\epsilon$, then simplifying for $\epsilon \to 0$ and eventually putting $\epsilon = 0$ to obtain the limit integral.

Thus, let us consider the following integral ($\epsilon > 0$):

$$ I_{\epsilon} = \int_{\epsilon}^{1/\epsilon} \frac{\omega(1 + \frac{\omega_0}{\omega_1})^2}{[\alpha_0^2 + (\omega - \omega_0)^2][\alpha_0^2 + \omega_0^2(1 + \frac{\omega}{\omega_1})^2]} d\omega. $$

This integral can be decomposed into the rational functions:

$$ I_{\epsilon} = \int_{\epsilon}^{1/\epsilon} \frac{A + B(\omega - \omega_0)}{\alpha_0^2 + (\omega - \omega_0)^2} d\omega + \int_{\epsilon}^{1/\epsilon} \frac{C + D(1 + \frac{\omega}{\omega_1})}{\alpha_0^2 + \omega_0^2(1 + \frac{\omega}{\omega_1})^2} d\omega, $$

where A, B, C and D are constant expressed as a function of $\omega_0$ and $\alpha_0$. Finally, the integration and putting $\epsilon = 0$ and $\alpha_0 = \omega_0/2Q$ yield to the expression for the normalization factor:

$$ C_h = \lim_{\epsilon \to 0} I_{\epsilon} = \left( \frac{2Q}{\omega_0} \right)^2 \frac{\left( \frac{\omega_0^2 - 1}{\omega_1} \right) \ln\left(\frac{\omega_0}{\omega_1}\right) + 4\left(1 + \frac{\omega_0}{\omega_1}\right)^2 Q \tan^{-1}(2Q)}{\left(1 - \frac{\omega_0}{\omega_1}\right)^2 + 4Q^2 \left(1 + \frac{\omega_0}{\omega_1}\right)^2}, $$

where $\omega_0$ is the mother pulsation and $Q$ the quality factor of the corresponding wave.

## APPENDIX B: DAMPING IDENTIFICATION METHOD

The method is based on a wavelet reconstruction formula for asymptotic signals where the stationary phase approximation is assumed: locally, the signal pulsation is constant while its amplitude is slowly modulating.

First step of the method is the selection of the ridges (e.g. Tchamitchian & Torresani 1992), that is, the lines of local maxima. We use the wavelet transform modulus maxima: this yields lines of coordinates $[b, a(b)]$ with wavelet modulus $|W_x[b, a(b)]|$.

Second step is to consider one line $[b, a(b)]$ to find the scale $a_e = a(b)$ for which the wavelet modulus is maximum; this scale $a_e$ then provides the corresponding pulsation $\omega_e = \omega_0/a_e$, where $\omega_0$ is the pulsation of the mother wavelet.

Third step is a linear regression in a semi-log plot of wavelet modulus maxima versus time $b$, used to determine $Q_e$ from the slope value. This is done by fitting the following expression:

$$ \ln \left( \frac{|W_x[b, a(b)]|^2}{|\hat{h}_0(\omega_e)|^2} \right) = 2 \ln A_0 - \frac{\omega_0 b}{Q_e}, \quad (B1) $$

where $|\hat{h}_0(\omega_e)|^2$ is the energy of the wavelet at the pulsation $\omega_e = 2\pi f_e$, and $A_0$ is a constant. For the CDW,

$$ |\hat{h}_0(\omega)|^2 = \frac{\omega^2 \left(1 + \frac{\omega_0}{\omega_1}\right)^2}{[\alpha_0^2 + (\omega - \omega_0)^2][\alpha_0^2 + \omega_0^2(1 + \frac{\omega}{\omega_1})^2]} $$

$Q_e$ is the resulting estimate for the quality factor $Q_s$ of the analysed signal.