Exact Near Instantaneous Frequency Formulas Best at Peaks (Part 2)

• Basics
• Tips and Tricks
• Tutorials

Introduction

This is an article that is a continuation of a digression from trying to give a better understanding of the Discrete Fourier Transform (DFT). It is recommended that my previous article "Exact Near Instantaneous Frequency Formulas Best at Peaks (Part 1)"[1] be read first as many sections of this article are directly dependent upon it.

A second family of formulas for calculating the frequency of a single pure tone in a short interval in the time domain is presented. It is very similar to the first family presented in the previous article. The formulas also work for complex tones. The formulas are still exact, meaning no approximations were taken in the derivation. The advantage of the second family over the first is that it offers a more robust solution over the same sampling interval. It also requires marginally more calculations to achieve this. Like the first family, this one should be evaluated at peaks of real tones and the formulas become indeterminate and noise sensitive at zero crossings. In the complex case, there is no such limitation.

Pure Tone Signal Definitions

This is the equation used in all my real tone blog articles to describe a discrete single pure tone: $$ S_n = M \cdot \cos( \alpha n + \phi ) \tag {1} $$ This is the corresponding equation for a pure complex tone. $$ S_n = M \cdot e^{i(\alpha n + \phi )} \tag {2} $$ In both types, the parameter being solved for with these families of formulas is $\alpha$, which is the one pertaining to the frequency.

One of the Family

Here is a member of the second family of formulas: $$ \begin{aligned} P_{n,1} &= S_{n+2} + S_{n-2} \\ P_{n,2} &= S_{n+4} + S_{n-4} \\ P_{n,3} &= S_{n+6} + S_{n-6} \\ P_{n,4} &= S_{n+8} + S_{n-8} \end{aligned} \tag {3} $$ $$ \alpha = \frac{1}{2} \cos^{-1}\left( \frac{ 30 S_n + 26 P_{n,1} + 16 P_{n,2} + 6 P_{n,3} + P_{n,4} } { 40 S_n + 30 P_{n,1} + 12 P_{n,2} + 2 P_{n,3} } \right) \tag {4} $$ Like in the previous article, this one is the $d=2,k=4$ case. The difference is that this family uses all the neighbor pairs in both the numerator and denominator. This makes the formulas more robust, that is noise resistant, though it does make for a few more calculations.

See Appendix A for a list of the formulas up to $k=9$.

Neighbor Pairs of Signal Values

The two neighbor sample values, on either side of the sample center by the same interval, are added together to form a neighbor pair sum. For both types of signals the formula is the same. In this article, the sample spacing will be incorporated from the start. $$ P_{n,m} = S_{n+md} + S_{n-md} = 2 S_n \cos[ (\alpha d) m ] \tag {5} $$ $$ S_n \cos[ (\alpha d) m ] = \frac{P_{n,m}}{2} \tag {6} $$

Powers of a Cosine Binomial

Instead of taking the powers of the cosine function, this family is based on taking the powers of a binomial containing the cosine. The calculations are a little more complicated and are shown in Appendix B. $$ \left[1 + \cos(\theta)\right]^0 = 1 \tag {7} $$ $$ \left[1 + \cos(\theta)\right]^1 = 1 + \cos(\theta) \tag {8} $$ $$ \left[1 + \cos(\theta)\right]^2 = \frac{3}{2} + 2 \cos(\theta) + \frac{1}{2} \cos(2\theta) \tag {9} $$ $$ \left[1 + \cos(\theta)\right]^3 = \frac{5}{2} + \frac{15}{4} \cos(\theta) + \frac{3}{2} \cos(2\theta) + \frac{1}{4} \cos(3\theta) \tag {10} $$ $$ \left[1 + \cos(\theta)\right]^4 = \frac{35}{8} + 7 \cos(\theta) + \frac{7}{2} \cos(2\theta) + \cos(3\theta) + \frac{1}{8} \cos(4\theta) \tag {11} $$ The main difference is that all the power terms appear instead of skipping every other one like the first family does.

Signal Calculations

The signal calculation takes the general form: $$ W_{n,k} = S_n \left[1 + \cos[ (\alpha d) m ]\right]^k \tag {12} $$ Where $n$ is the center sample index and $k$ is the degree of the exponentiation.

The first few examples can be found from the cosine power series in the previous section. $$ W_{n,0} = S_n \tag {13} $$ $$ W_{n,1} = S_n + S_n \cos[ (\alpha d) m ] = S_n + \frac{1}{2} P_{n,1} \tag {14} $$ $$ \begin{aligned} W_{n,2} &= \frac{3}{2} S_n + 2 S_n \cos[ (\alpha d) m ] + \frac{1}{2} S_n \cos(2\alpha) \\ &= \frac{3}{2} S_n + P_{n,1} + \frac{1}{4} P_{n,2} \end{aligned} \tag {15} $$ $$ W_{n,3} = \frac{5}{2} S_n + \frac{15}{8} P_{n,1} + \frac{3}{4} P_{n,2} + \frac{1}{8} P_{n,3} \tag {16} $$ $$ W_{n,4} = \frac{35}{8} S_n + \frac{7}{2} P_{n,1} + \frac{7}{4} P_{n,2} + \frac{1}{2} P_{n,3} + \frac{1}{16} P_{n,4} \tag {17} $$ The coefficients for the first nine formulas can be found in Appendix C.

Frequency Calculation

Finding the frequency value is similar to how it was done with the first family. Taking the quotient of two adjacent signal calculations yields an expression for the frequency term. $$ q = \frac{ W_{n,k} }{ W_{n,k-1} } = \frac{S_n \left[1 + \cos[ (\alpha d) m ]\right]^{k}}{S_n \left[1 + \cos[ (\alpha d) m ]\right]^{k-1}} =1 + \cos[ (\alpha d) m ] \tag {18} $$ Solving for the frequency is then very straighforward: $$ \alpha = \cos^{-1}(q-1) = \cos^{-1}\left(\frac{ W_{n,k} }{ W_{n,k-1} }-1\right) \tag {19} $$

Base Case

The base case occurs when $k$ is set to zero and $d$ is set to one. $$ q = \frac{ W_{n,1} }{ W_{n,0} } = \frac{S_n + \frac{1}{2} P_{n,1}}{S_n} = 1 + \frac{ P_{n,1} }{2 S_n} \tag {20} $$ Applying a little algebra shows that the base case of this family is the same as the base case of the first family. $$ \alpha = \cos^{-1}\left( \frac{ P_{n,1} }{2 S_n} \right) = \cos^{-1}\left( \frac{ S_{n+1}+S_{n-1} }{2 S_n} \right) \tag {21} $$

Better Signal Value

Similar as what was done in the previous article, (12) can be rearranged to give a formula for an estimated value of the signal at a given sample point. $$ G_n = \frac{ W_{n,k} }{ \left[1 + \cos(\alpha d)\right]^k } = \frac{ W_{n,k} }{ q^k } \tag {22} $$ In the noiseless case, this guess will match the true signal value. In a noisy case it is likely to be more accurate.

Simplifying the Formula

The fractions in the coefficients can be eliminated by multiply the numerator and denominator of the $q$ calculation by the appropriate power of two. $$ q = \frac{ W_{n,k} }{ W_{n,k-1} } = \frac{ 2^k W_{n,k} }{ 2^k W_{n,k-1} } \tag {23} $$ Subtracting the one from $q$ to get $r$ can also be precalculated, resulting is simpler formulas for the frequency calculation. $$ r = q - 1 = = \frac{ 2^k W_{n,k} }{ 2^k W_{n,k-1} } - 1 = \frac{ 2^k W_{n,k} - 2^k W_{n,k-1} }{ 2^k W_{n,k-1} } \tag {24} $$ The actual calculation for the frequency also is simplified. $$ \alpha = \frac{1}{d} \cos^{-1}(r) \tag {25} $$ The coefficients for the first nine formulas in the family can be found in Appendix A.

Numerical Example

Here is a sample calculation of the formula given in (3). Since it is a noiseless case, the accuracy is good to the precision of the numerical variables, which is higher than the values displayed.

    alpha =  0.0626894

   S[144] =  2.6701126
   S[145] =  2.7086362
   S[146] =  2.7365186
   S[147] =  2.7536500
   S[148] =  2.7599633
   S[149] =  2.7554336
   S[150] =  2.7400787
   S[151] =  2.7139589
   S[152] =  2.6771768

      Pn1 =  5.5090836
      Pn2 =  5.4765972
      Pn3 =  5.4225951
      Pn4 =  5.3472894

      Vn3 = 22.0147123
      Vn4 = 43.9861803

        q =  1.9980357

  alpha_n =  0.0626894
       Gn =  2.7599633

The Larger Family

Replacing the one in (12) with a variable $x$ makes the two families of equations covered in this and my previous article merely two points on a much larger set. $$ W_{n,k}(x) = S_n \left[x + \cos[ (\alpha d) m ]\right]^k \tag {26} $$ The first family from the previous article corresponds to $x=0$ and the second family, covered in this article, corresponds to $x=1$. The effect of varying x is changing the weighting of the sample points around the center.

People who are good with filters (and I'm not really one) will recognize that the numerator and denominator of the frequency calculations are equivalent to finite impulse response (FIR) filters. I have steered away from employing this terminology because the derivation of these formulas is inherently center based and FIRs are generally defined relative to the last sample. Therefore calling them FIRs may lead to more confusion than benefit.

Limited Applicability

The formulas in this article and the previous article are really only applicable to single tone low noise signals. This makes their usefulness rather limited except for those specific applications. The mathematics are pretty neat though and in those few cases where these formulas are useful, they are probably going to be the best ones. I have not done any extensive testing to see what settings are optimal for any given situation. One reason is the limited applicability. Another reason is because of the wide range of parameter settings would make the testing large and cumbersome, and the results would correspondingly be large and cumbersome also. Having said that, it is still on my todo list, but readers are welcome to beat me to it.

Conclusion

A second family of exact time domain frequency formulas for single tones has been derived. It turns out just to be one of infinitely many such families. A lot more testing is needed to determine what the optimal parameter values are for any given situation.

Reference

[1] Dawg, Cedron, Exact Near Instantaneous Frequency Formulas Best at Peaks (Part 1)

Appendix A. Simplified Formulas

These are the coefficients for the first nine simplified formulas given by (23). The numerator and denominator need to be dotted with the vector $ \begin{bmatrix} S_n & P_{n,1} & P_{n,2} & P_{n,3} & ... \end{bmatrix} $ in order to be evaluated. $$ 1:\frac{ \begin{bmatrix} 0 & 1 \end{bmatrix} }{ \begin{bmatrix} 2 \end{bmatrix} } \tag {27} $$ $$ 2:\frac{ \begin{bmatrix} 2 & 2 & 1 \end{bmatrix} }{ \begin{bmatrix} 4 & 2 \end{bmatrix} } \tag {28} $$ $$ 3:\frac{ \begin{bmatrix} 8 & 7 & 4 & 1 \end{bmatrix} }{ \begin{bmatrix} 12 & 8 & 2 \end{bmatrix} } \tag {29} $$ $$ 4:\frac{ \begin{bmatrix} 30 & 26 & 16 & 6 & 1 \end{bmatrix} }{ \begin{bmatrix} 40 & 30 & 12 & 2 \end{bmatrix} } \tag {30} $$ $$ 5:\frac{ \begin{bmatrix} 112 & 98 & 64 & 29 & 8 & 1 \end{bmatrix} }{ \begin{bmatrix} 140 & 112 & 56 & 16 & 2 \end{bmatrix} } \tag {31} $$ $$ 6:\frac{ \begin{bmatrix} 420 & 372 & 255 & 130 & 46 & 10 & 1 \end{bmatrix} }{ \begin{bmatrix} 504 & 420 & 240 & 90 & 20 & 2 \end{bmatrix} } \tag {32} $$ $$ 7:\frac{ \begin{bmatrix} 1584 & 1419 & 1012 & 561 & 232 & 67 & 12 & 1 \end{bmatrix} }{ \begin{bmatrix} 1848 & 1584 & 990 & 440 & 132 & 24 & 2 \end{bmatrix} } \tag {33} $$ $$ 8:\frac{ \begin{bmatrix} 6006 & 5434 & 4004 & 2366 & 1092 & 378 & 92 & 14 & 1 \end{bmatrix} }{ \begin{bmatrix} 6864 & 6006 & 4004 & 2002 & 728 & 182 & 28 & 2 \end{bmatrix} } \tag {34} $$ $$ 9:\frac{ \begin{bmatrix} 22880 & 20878 & 15808 & 9828 & 4928 & 1940 & 576 & 121 & 16 & 1 \end{bmatrix} }{ \begin{bmatrix} 25740 & 22880 & 16016 & 8736 & 3640 & 1120 & 240 & 32 & 2 \end{bmatrix} } \tag {35} $$

Appendix B. Cosine Binomial Calculations

These are the raw calculations for the first few cosine binomial expansions. Each equation is then tested for accuracy by plugging in $0$, $\pi$, and $\pi/2$ to ensure the results are $2^k$, $0$, and $1$ respectively. $$ \begin{aligned} \left[1 + \cos(\theta)\right]^2 &= 1 + 2 \cos(\theta) + \cos^2(\theta) \\ &= 1 + 2 \cos(\theta) + \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^2 \\ &= 1 + 2 \cos(\theta) + \frac{e^{i2\theta}+2+e^{-i2\theta}}{4} \\ &= \frac{3}{2} + 2 \cos(\theta) + \frac{1}{2} \cos(2\theta) \end{aligned} \tag {36} $$ $$ \begin{aligned} \theta = 0 \quad &\Rightarrow \quad \frac{3}{2} + 2 + \frac{1}{2} = 4 \\ \theta = \pi \quad &\Rightarrow \quad \frac{3}{2} - 2 + \frac{1}{2} = 0 \\ \theta = \frac{\pi}{2} \quad &\Rightarrow \quad \frac{3}{2} - \frac{1}{2} = 1 \end{aligned} \tag {37} $$ $$ \begin{aligned} \left[1 + \cos(\theta)\right]^3 &= 1 + 3 \cos(\theta) + 3 \cos^2(\theta) +\cos^3(\theta) \\ &= 1 + 3 \frac{e^{i\theta}+e^{-i\theta}}{2} + 3 \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^2 + \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^3 \\ &= 1 + \frac{3}{2} ( e^{i\theta}+e^{-i\theta} ) + \frac{3}{4} ( e^{i2\theta}+2+e^{-i2\theta} ) \\ &+ \frac{1}{8} ( e^{i3\theta}+3e^{i\theta}+3e^{-i\theta}+e^{-i3\theta} ) \\ &= \frac{5}{2} + \frac{15}{8} ( e^{i\theta}+e^{-i\theta} ) + \frac{3}{4} ( e^{i2\theta}+e^{-i2\theta} ) + \frac{1}{8} ( e^{i3\theta}+e^{-i3\theta} ) \\ &= \frac{5}{2} + \frac{15}{4} \cos(\theta) + \frac{3}{2} \cos(2\theta) + \frac{1}{4} \cos(3\theta) \end{aligned} \tag {38} $$ $$ \begin{aligned} \theta = 0 \quad &\Rightarrow \quad \frac{5}{2} + \frac{15}{4} + \frac{3}{2} + \frac{1}{4} = 8 \\ \theta = \pi \quad &\Rightarrow \quad \frac{5}{2} - \frac{15}{4} + \frac{3}{2} - \frac{1}{4} = 0 \\ \theta = \frac{\pi}{2} \quad &\Rightarrow \quad \frac{5}{2} -\frac{3}{2} = 1 \end{aligned} \tag {39} $$ $$ \begin{aligned} \left[1 + \cos(\theta)\right]^4 &= 1 + 4 \cos(\theta) + 6 \cos^2(\theta) + 4 \cos^3(\theta) + \cos^4(\theta) \\ &= 1 + 4 \frac{e^{i\theta}+e^{-i\theta}}{2} + 6 \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^2 \\ &+ 4 \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^3 + \left[ \frac{e^{i\theta}+e^{-i\theta}}{2}\right]^4 \\ &= 1 + 2 ( e^{i\theta}+e^{-i\theta} ) + \frac{3}{2} ( e^{i2\theta}+2+e^{-i2\theta} ) \\ &+ \frac{1}{2} ( e^{i3\theta}+3e^{i\theta}+3e^{-i\theta}+e^{-i3\theta} ) \\ &+ \frac{1}{16} ( e^{i4\theta}+4e^{i2\theta}+6+4e^{-i2\theta}+e^{-i4\theta} ) \\ &= \frac{35}{8} + \frac{7}{2} ( e^{i\theta}+e^{-i\theta} ) + \frac{7}{4} ( e^{i2\theta}+e^{-i2\theta} ) \\ &+ \frac{1}{2} ( e^{i3\theta}+e^{-i3\theta} ) + \frac{1}{16} ( e^{i4\theta}+e^{-i4\theta} ) \\ &= \frac{35}{8} + 7 \cos(\theta) + \frac{7}{2} \cos(2\theta) + \cos(3\theta) + \frac{1}{8} \cos(4\theta) \end{aligned} \tag {40} $$ $$ \begin{aligned} \theta = 0 \quad &\Rightarrow \quad \frac{35}{8} + 7 + \frac{7}{2} + 1 + \frac{1}{8} = 16 \\ \theta = \pi \quad &\Rightarrow \quad \frac{35}{8} - 7 + \frac{7}{2} - 1 + \frac{1}{8} = 0 \\ \theta = \frac{\pi}{2} \quad &\Rightarrow \quad \frac{35}{8} - \frac{7}{2} + \frac{1}{8} = 1 \end{aligned} \tag {41} $$ The general form looks like this: $$ \begin{aligned} \left[1 + \cos(\theta)\right]^n &= \sum_{k=0}^{n} { \binom{n}{k} \cos^k(\theta) } \\ &= \sum_{k=0}^{n} { \binom{n}{k} \left[ \frac{ e^{i\theta} + e^{-i\theta} }{ 2 } \right]^k } \\ &= \sum_{k=0}^{n} { \binom{n}{k} \frac{ 1 }{ 2^k } \sum_{j=0}^{k} { \binom{k}{j} e^{i\theta (k-j)} e^{-i\theta j } } } \\ &= \sum_{k=0}^{n} { \sum_{j=0}^{k} { \frac{ 1 }{ 2^k } \binom{n}{k} \binom{k}{j} e^{i\theta (k-2j)} } } \\ &= \sum_{k=0}^{n} { \sum_{j=0}^{k} { \frac{ 1 }{ 2^k } \frac{n!}{ (n-k)! \cdot j! \cdot (k-j)! } e^{i\theta (k-2j)} } } \end{aligned} \tag {42} $$

Appendix C. Calculation Coefficients

These are the coefficients for the first nine signal calculation. When used in the frequency calculation it is better to use the simplified coefficient values found in Appendix A.

1:	1	$\frac{1}{2}$
2:	$\frac{3}{2}$	1	$\frac{1}{4}$
3:	$\frac{5}{2}$	$\frac{15}{8}$	$\frac{3}{4}$	$\frac{1}{8}$
4:	$\frac{35}{8}$	$\frac{7}{2}$	$\frac{7}{4}$	$\frac{1}{2}$	$\frac{1}{16}$
5:	$\frac{63}{8}$	$\frac{105}{16}$	$\frac{15}{4}$	$\frac{45}{32}$	$\frac{5}{16}$	$\frac{1}{32}$
6:	$\frac{231}{16}$	$\frac{99}{8}$	$\frac{495}{64}$	$\frac{55}{16}$	$\frac{33}{32}$	$\frac{3}{16}$	$\frac{1}{64}$
7:	$\frac{429}{16}$	$\frac{3003}{128}$	$\frac{1001}{64}$	$\frac{1001}{128}$	$\frac{91}{32}$	$\frac{91}{128}$	$\frac{7}{64}$	$\frac{1}{128}$
8:	$\frac{6435}{128}$	$\frac{715}{16}$	$\frac{1001}{32}$	$\frac{273}{16}$	$\frac{455}{64}$	$\frac{35}{16}$	$\frac{15}{32}$	$\frac{1}{16}$	$\frac{1}{256}$
9:	$\frac{12155}{128}$	$\frac{21879}{256}$	$\frac{1989}{32}$	$\frac{4641}{128}$	$\frac{1071}{64}$	$\frac{765}{128}$	$\frac{51}{32}$	$\frac{153}{512}$	$\frac{9}{256}$	$\frac{1}{512}$

Sign in