# Fitting a Damped Sine Wave

A damped sine wave is described by

$$ x_{(k)} = A \cdot e^{\alpha \cdot k} \cdot cos(\omega \cdot k + p)

\tag{1}$$

with frequency $\omega$ , phase p , initial amplitude A and damping constant $\alpha$ . The $x_{(k)}$ are the samples of the function at equally spaced points in time.

With $x_{(k)}$ given, one often has to find the unknown parameters of the function. This can be achieved for instance with nonlinear approximation or with DFT – methods.

I present a method to find the unknown parameters which is based on linear algebra only.

The equation (1) is the solution to the difference equation

$$x_{(k+2)} = 2 \cdot \sqrt b \cdot cos(\omega) \cdot x_{(k+1)} - b

\cdot x_{(k)} \tag{2}$$

This difference equation is a linear equation system for the unknown constants b and $2 \cdot \sqrt b \cdot cos(\omega) $.

The overdetermined linear system

$$

\left(

\begin{matrix}

x_{(k+2)} \\

x_{(k+3)} \\

. \\

.

\end{matrix}

\right)

=

\left(

\begin{matrix}

x_{(k+1)} & x_{(k )} \\

x_{(k+2)} & x_{(k+1)} \\

. & .\\

. & .

\end{matrix}

\right)

\cdot

\left(

\begin{matrix}

2 \cdot \sqrt b \cdot cos(\omega) \\

-b

\end{matrix}

\right)

\tag{3}$$is solved using linear regression analysis.

There is a pitfall in this calculation: For high oversampling rates, i.e. small resulting $\omega$, the numerical results with additional noise are poor. Performance is best near half the Nyquist rate, which results in approximately 4 samples per sine wave. The $x_{(k)}$ can be rearranged accordingly without losing samples or accuracy, such as taking every 8th sample for the $x_{(k)},x_{(k+1)}, ... $ The resultant $\omega$ has to be corrected by this integer factor.

The $b$ is related to $\alpha$ by

$$e^\alpha\ = \sqrt b \tag{4}$$

For the derivation of this relation see the appendix.

So the $\omega$, $\alpha$ and b are known.

The amplitude A and the phase $p$ are still unknown and calculated by

'mixing with a quadrature carrier' :

The system of equations

$$

\left(

\begin{matrix}

x_{(0)} \\

x_{(1)} \\

. \\

.

\end{matrix}

\right)

=

\left(

\begin{matrix}

e^{\alpha \cdot 0} \cdot \sin{(\omega \cdot 0)} & e^{\alpha \cdot 0}

\cdot \cos{(\omega \cdot 0)} \\

e^{\alpha \cdot 1} \cdot \sin{(\omega \cdot 1)} & e^{\alpha \cdot 1}

\cdot \cos{(\omega \cdot 1)} \\

. & .\\

. & .

\end{matrix}

\right)

\cdot

\left(

\begin{matrix}

u \\

v

\end{matrix}

\right)

\tag{5}$$

is solved for u and v. With the u, v interpreted as the complex number $u+iv$ we obtain the amplitude and phase

$$A=|u+iv|$$

$$p=\sphericalangle{(u-iv)}$$

In real setups the samples $x_{(k)}$ suffer from additional noise. It turns out that the results for frequency $\omega$ and phase $p$ are very stable in noise. That is not the case for amplitude A nor for the damping factor $\alpha$. These parameters are adjusted with a gradient iteration. Amplitude and damping are slightly varied until the best fit to the sampled data is achieved.

I published a corresponding Matlab script in the 'Mathworks file exchange' system. You can find it by searching for 'matlab fit damped sine wave'. Download is free and you do not need a registration. I hope that this algorithm is useful and download numbers will rise.

Cheers

Detlef

Appendix

Replacing $x_{(k)}$ in (2) with its definition in (1) yields

$$ \begin{align}

& A e^{\alpha \cdot (k+2)} \cos {(\omega (k+2) + p)} = \\

& \\

& 2 \sqrt{b} \cos{(\omega) A e^{\alpha(k+1)}} \cdot \cos{(\omega (k+1)+p)} \\

& - b \cdot A \cdot e^{\alpha \cdot k} \cdot \cos{(\omega k + p)} \\

\end{align}

$$

Division by $A \cdot e^{\alpha \cdot k}$ and application of addition theorems for cos-functions gives

$$\begin{align}

& e^{2 \alpha} [\cos{(\omega k + p)} \cos{(2 \omega)}-\sin{(\omega k + p)} \sin{(2 \omega)} ] = \\ \\

& 2 e^{\alpha}\sqrt{b}\cos{\omega} [\cos{(\omega k + p)} \cos{(\omega)}-\sin{(\omega k + p)} \sin{(\omega)}] \\

& - b \cos{(\omega k+p)} \\

\end{align}

$$

Rearranging terms

$$\begin{align}

0=& & \cos{(\omega k + p)} [e^{2 \cdot \alpha} \cos{(2 \omega)}-2e^\alpha \sqrt b \cos ^2 {\omega} +b ]- \\

& & \sin{(\omega k + p)} [e^{2 \cdot \alpha} \sin{(2 \omega)}-2e^\alpha \sqrt b \sin {\omega}\sin {\omega} ]= \\ \\

0= & & \cos{(\omega k + p)} [e^{2 \cdot \alpha} \cos{(2 \omega)}-2e^\alpha \sqrt b \cos ^2 {\omega} +b ]- \\

& & \sin{(\omega k + p)}\sin{(2 \omega)} \cdot [e^{2 \alpha } -e^{ \alpha }\sqrt b ] \\

\end{align}

$$

The last equation can only hold for any k if the expressions in the squared brackets equal 0:

$$ 0= e^{2 \alpha } -e^{ \alpha }\sqrt b

$$

$$ \Rightarrow e^{\alpha } = \sqrt b

$$

- Comments
- Write a Comment Select to add a comment

Kaz

$\omega$ is a normalized frequency, it is related to the sample frequency. So for a given signal, if the sample rate gets higher, the $\omega$ is going down. From the linear system of equations the $\cos \omega$ term is calculated. This function is flat for small $\omega$, a small change in the $\cos \omega$ will result in a big change for $\omega$. The situation ist best near $\pi / 2$, because the $\cos \omega$ function ist steepest there, and this sample rate is corresponding to Nyquist / 2 . And yes, higher oversampling rate gives better results, you use all samples, but you have to rearrange the samples. Take a look at the Matlab-script.

Cheers

Detlef

Regards,

[-Rick-]

thanks for the hint. For Google 'matlab fit damped sine wave', the MathWork Exchange link is the first hit, at least for my Google results. I hope, that your efforts to understand my blog were crowned with success. I'd be more than glad to give some additional clarifications.

Cheers

Detlef

To post reply to a comment, click on the 'reply' button attached to each comment. To post a new comment (not a reply to a comment) check out the 'Write a Comment' tab at the top of the comments.

Registering will allow you to participate to the forums on ALL the related sites and give you access to all pdf downloads.