DSPRelated.com
Blogs

Exact Near Instantaneous Frequency Formulas Best at Zero Crossings

Cedron DawgJuly 20, 2017

Introduction

This is an article that is the last of my digression from trying to give a better understanding of the Discrete Fourier Transform (DFT). It is along the lines of the last two.

In those articles, I presented exact formulas for calculating the frequency of a pure tone signal as instantaneously as possible in the time domain. Although the formulas work for both real and complex signals (something that does not happen with frequency domain formulas), for real signals they are only applicable near the peak values of the signal and perform poorly, or not at all, at zero crossings. This article fills that gap by providing formulas that work best at zero crossings and perform poorly, or not at all, at peaks.

You should first read my two previous on "Exact Near Instantaneous Frequency Formulas Best at Peaks"[1][2] to understand the context of this article. Just like the previous formulas, these also work with complex signals.

Pure Tone Signal Definition

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} $$

Differences of Real Neighbor Pairs

Instead of summing neighbor pairs, this approach takes the difference. This is what makes them work best at zero crossings and perform poorly at peaks.

These are the definitions of the neighbor pair signal values for a real valued signal. $$ \begin{aligned} S_{n+m} &= M \cdot \cos[ \alpha (n+m) + \phi ] \\ &= M \cdot \cos[ (\alpha n + \phi) + (\alpha m) ] \\ &= M \cdot \left[ \cos( \alpha n + \phi ) \cos( \alpha m ) - \sin( \alpha n + \phi ) \sin( \alpha m ) \right] \end{aligned} \tag {3} $$ $$ S_{n-m} = M \cdot \left[ \cos( \alpha n + \phi ) \cos( \alpha m ) + \sin( \alpha n + \phi ) \sin( \alpha m ) \right] \tag {4} $$ Taking their difference simplifies the equation as two of the terms cancel. Unlike the sum version in the previous articles, the expression for the center signal value does not reappear. $$ \begin{aligned} D_{n,m} &= S_{n+m} - S_{n-m} \\ &= -2M \sin( \alpha n + \phi ) \sin( \alpha m ) \end{aligned} \tag {5} $$

Ratios of Real Neighbor Differences

When cosines are raised to a power, as in my two previous articles, the result is a linear expression of cosines of multiples of the angle. This is not true for sines. When sines are raised to a power, the result is a mix of sines and cosines of multiples of the angle. Therefore the same approach that was used in the previous articles can't be used here in the same manner.

Instead, a different trick is employed by exploiting the properties of the double angle formula for sines. $$ \begin{aligned} D_{n,2m} &= -2M \sin( \alpha n + \phi ) \sin( \alpha 2m ) \\ &= -2M \sin( \alpha n + \phi ) 2 \sin( \alpha m ) \cos( \alpha m ) \end{aligned} \tag {6} $$ By applying the double angle formula on $ D_{n,2m} $, a copy of $ D_{n,m} $ emerges which can be factored out by taking the quotient of the two. $$ \begin{aligned} Q_{n,m} &= \frac{ D_{n,2m} }{ D_{n,m} } \\ &= \frac{ -2M \sin( \alpha n + \phi ) 2 \sin( \alpha m ) \cos( \alpha m ) }{ -2M \sin( \alpha n + \phi ) \sin( \alpha m ) } \\ &= 2 \cos( \alpha m ) \end{aligned} \tag {7} $$ This results in an expression which is a simple cosine.

Differences of Complex Neighbor Pairs

The same approach can be taken using the definition of a pure complex tone. $$ \begin{aligned} S_{n+m} &= M \cdot e^{ i ( \alpha (n+m) + \phi ) } \\ &= M \cdot e^{ i ( \alpha n + \phi ) } e^{ i \alpha m } \\ \end{aligned} \tag {8} $$ $$ S_{n-m} = M \cdot e^{ i ( \alpha n + \phi ) } e^{ -i \alpha m } \tag {9} $$ Taking the difference, and applying Euler's Equation[3], shows that the difference is an imaginary multiple of the signal value. $$ \begin{aligned} D_{n,m} &= S_{n+m} - S_{n-m} \\ &= M \cdot e^{ i ( \alpha n + \phi ) } \left( e^{ i \alpha m } - e^{ -i \alpha m } \right) \\ &= S_n \cdot 2i \cdot \sin( \alpha m ) \end{aligned} \tag {10} $$ The result is quite different than the real values case. Because the signal value reappears in the formula, the value of $\alpha$ can be directly calculated from this formula. $$ \alpha_{n,m} = \frac{1}{m} \sin^{-1} \left( \frac{ S_{n+m} - S_{n-m} }{ 2 i S_{n} } \right) \tag {11} $$ This formula only works with complex tones. The problem with this formula is that it only works up to half the Nyquist frequency. It is analogous to the base case formulas in my previous two blogs which works in either case: $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left( \frac{ S_{n+m} + S_{n-m} }{ 2 S_{n} } \right) \tag {12} $$ The problem with this formula in the complex case is that it does not distinguish between positive and negative values of $\alpha$. A tone with a negative frequency will get a positive frequency as the answer. This problem can be avoided by using (11) in conjunction or instead.

Ratios of Complex Neighbor Differences

The same ratio approach as used in the real case can be used to find the quotient value in the complex case. $$ \begin{aligned} D_{n,2m} &= S_n \cdot 2i \cdot \sin( \alpha 2m ) \\ &= S_n \cdot 2i \cdot 2 \sin( \alpha m ) \cos( \alpha m ) \end{aligned} \tag {13} $$ $$ \begin{aligned} Q_{n,m} &= \frac{ D_{n,2m} }{ D_{n,m} } \\ &= \frac{ S_n \cdot 2i \cdot 2 \sin( \alpha m ) \cos( \alpha m ) }{ S_n \cdot 2i \cdot \sin( \alpha m ) } \\ &= 2 \cos( \alpha m ) \end{aligned} \tag {14} $$ Surprisingly, or perhaps not, the result is identical to the real case.

Simple Frequency Calculation

The value for the frequency term $ \alpha $ can now be solved for from both (7) for the real case and (14) for the complex case. $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1}\left( \frac{ Q_{n,m} }{ 2 } \right) \tag {15} $$ Plugging in the definition of $ Q_{n,m} $ from (5) and (7), or (10) and (14), gives the result in terms of signal values: $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( \frac{ S_{n+2m} - S_{n-2m} }{ S_{n+m} - S_{n-m} } \right) \right] \tag {16} $$ Since the frequency term ($\alpha$) can be calculated at any point along the signal, and with any spacing. Any weighted average, with the weights adding up to one, can be used to calculate $\alpha$ over an interval in order to mitigate the effects of noise.

Comparison to Turner's Formula

In Rick Lyons' article[4], which prompted my digression from the frequency domain to the time domain, he cites a formula by Turner for a four point calculation of frequency for a pure real tone. This is equation (2) from Lyons' article translated into my nomenclature and generalized with a sample spacing parameter ($m$) built in: $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( \frac{ S_{n+2m} - S_{n-m} }{ S_{n+m} - S_{n} } - 1 \right) \right] \tag {17} $$ The similarity between (16) and (17) is striking. They are both four sample formulas. The difference is the width of the stance on the signal. The samples in (16) are centered around sample $n$ and are five samples wide. The samples in (17) are centered around sample $n+\frac{1}{2}$ and are four samples wide. Thus in terms of being more instantaneous, Turner's formula is superior.

Like the formulas derived in this article, Turner's formula will work best at zero crossings and perform poorly at peaks for real signals. Also, like the formulas in this article, Turner's 4 Sample formula will also work on complex tones, with the positive/negative caveat mentioned above. This is not mentioned in Turner's derivation or Lyons' article. The proof of this is shown in Appendix A.

Like (16), (17) can be calculate at various points along an interval with various spacing and used in a weighted average to get a noise mitigated result.

Turner's 3 Sample formula, equation (1) in Lyons' article, is the same as (12) when converted to arbitrary sample spacing. This is the base case for the formulas that are best at the peak.

Cosine Binomial Frequency Calculation

The quotient of the differences approach has the advantage that the amplitude term has already been cancelled out. Therefore the direct calculation approach as shown above is possible. However, since (7) and (14) can be used to find the cosine values of multiples of $\alpha$, the same approach as used in the previous two articles can be used as well.

The two families of formulas derived in the previous two articles can be represented by this formula with $x$ taking the value of zero or one. $$ q = \frac{\left[x + \cos(\alpha)\right]^{k}}{\left[x + \cos(\alpha)\right]^{k-1}} =x + \cos(\alpha) \tag {18} $$ When the binomial expressions are multiplied out the result is an expression of cosine values of multiples of $\alpha$. The resulting coefficients for $x=0$ and $x=1$ can be found in the previous articles. For other values of $x$, the coefficient values can be pre-computed for implementation. Once $q$ has been evaluated, the value of $\alpha$ can be found by this equation: $$ \alpha = \cos^{-1}(q-x) \tag {19} $$ Once again, since it is an inverse cosine being evaluated, negative $\alpha$ values in complex signal cases will seem to be positive.

Since the Turner formulas can also be used to find the cosine value of $\alpha$ multiples from signal values, they can also be used in the binomial expansion evaluations.

Conclusion

Although they work for complex signals, these formulas were derived to provide near instaneous frequency readings for real pure tone signals. Two set of formulas were derived in the previous articles for when the signal is near peaks and the formulas in this article are for when the signal is near zero crossings. The application of these formulas is limited to the use on relatively noiseless pure tones. They will be most useful in applications where low latency is critical. When latency is not critical, frequency domain solutions, such as evaluating two or three bins of a DFT and using corresponding formulas will yield superior results. If anybody actually ever uses these formulas, I would love to hear about it in the comments.

Reference

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

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

[3] Dawg, Cedron, The Exponential Nature of the Complex Unit Circle

[4] Lyons, Rick, Sinusoidal Frequency Estimation Based on Time-Domain Samples

Appendix A: Proof of Turner's Formula in Complex Pure Tones

Here is the proof that Turner's 4 Sample Formula works for a complex pure tone as well as a real one as long as the tone has a positive frequency.

Start by plugging in the signal definition (2) into the formula (17): $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( \frac{ M \cdot e^{i[\alpha (n+2m) + \phi ]} - M \cdot e^{i[\alpha (n-m) + \phi ]} }{ M \cdot e^{i[\alpha (n+m) + \phi ]} - M \cdot e^{i[\alpha n + \phi ]} } - 1 \right) \right] \tag {20} $$ Every term in the quotient can have the signal formula factored out. $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( \frac{ e^{i\alpha 2m } - e^{-i\alpha m } }{ e^{i\alpha m } - 1 } - 1 \right) \right] \tag {21} $$ The numerator and denominator can be divided by $ e^{i\alpha \frac{1}{2}m } $ to make the exponents symmetric about zero. $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( \frac{ e^{i\alpha \frac{3}{2}m } - e^{-i\alpha \frac{3}{2}m } }{ e^{i\alpha \frac{1}{2}m } - e^{-i\alpha \frac{1}{2}m } } - 1 \right) \right] \tag {22} $$ The quotient can now be recognized as being the difference of two cubes divided by a difference. $$ \frac{ a^3 - b^3 }{ a - b } = a^2 + ab + b^2 \tag {23} $$ This simplifies the quotient considerably. $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ \frac{ 1 }{ 2 } \left( e^{i\alpha m } + 1 + e^{-i\alpha m } - 1 \right) \right] \tag {24} $$ The ones cancel, and applying Euler's Equation[3] to the other terms leaves a cosine expression. $$ \alpha_{n,m} = \frac{1}{m} \cos^{-1} \left[ cos( \alpha m ) \right] \tag {25} $$ For positive values of $\alpha$, less than $\pi$, the inverse cosine cancels the cosine. Then the $m$'s cancel and this remains: $$ \alpha_{n,m} = \alpha \tag {26} $$ For any sample point ($n$), for any spacing ($m$), the answer the formula will give is the frequency term of the signal.

Quite Easily Done.



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.

Please login (on the right) if you already have an account on this platform.

Otherwise, please use this form to register (free) an join one of the largest online community for Electrical/Embedded/DSP/FPGA/ML engineers: