Frequency Formula for a Pure Complex Tone in a DTFT

• Academia / Research
• DFT_FFT
• Tutorials

Introduction

This is an article to hopefully give a better understanding of the Discrete Fourier Transform (DFT) by deriving frequency formula for the closely related Discrete Time Fourier Transform (DTFT). The distinction between the two is the latter has a domain on the integers from negative infinity to positive infinity and can be evaluated for any frequency. The DFT has a finite domain of $N$ samples and is only evaluated with frequencies that are a whole number of cycles across the frame. This article considers a rectangular windowed version of the DTFT for a set of $N$ sample points. By evaluating just three bins of a pure complex tone signal, the frequency of the tone can be determined. Like the DFT formulas [1], these calculations will be most robust in the presence of noise when the evaluation frequency is close to the target tone's frequency.

Angles from Bin Numbers

When using a DFT, the evaluated frequencies correspond to whole number of cycles, and that cycle count defines the bin number, or bin index, within the produced spectrum. For the DTFT, these values can be fractional, but it does not change the correspondence to the angle step used in the calculations. It simply means $k$ can take non-integer values. Therefore, in this article my usual convention of subscripts for indices is replaced with function notation. $$ \beta(k) = k \cdot \frac{ 2\pi }{ N } \tag {1} $$ The frequency of the tone can also be put on the same scale with a similar conversion factor. $$ \alpha = f \cdot \frac{ 2\pi }{ N } \tag {2} $$ The value of $\alpha$ is the angle step size in radians between samples in the signal. This is also the same as the arc length along the unit circle.

Signal Values

For a pure complex tone, the signal value at each sample can be calculated using a standard formula. $$ S_n = M e^{ i ( \alpha n + \phi ) } = M e^{ i \phi } \left[e^{ i \alpha} \right]^n = M e^{ i \phi } \left[e^{ i \frac{ 2\pi }{ N } f} \right]^n = M e^{ i \phi } e^{ i 2\pi \frac{fn}{N} } \tag {3} $$ The signal values are simply points around the unit circle multiplied by the amplitude $(M)$. They start at $\phi$ distance around the perimeter, and then proceed step wise from there. There is much more on the math of this in my earlier articles [2],[3].

Normalized DFT/DTFT

The finite sampled version of the DTFT formula is exactly the same as the DFT. $$ Z(k) = \frac{ 1 }{ N } \sum_{n=0}^{N-1} S_n \cdot e^{ -i \beta(k) n } \tag {4} $$ Again, function notation is used since $k$ can now take on fractional values. The normalization is simply preferred to be consistent with my other articles. The normalization factor gets factored out as part of the calculation so there is no need to perform a normalization when working with libraries that return unnormalized spectrums.

Normalized DFT/DTFT of a Pure Tone

The definition of the signal can now be plugged into the definition of the DTFT. The amplitude and phase factor can be factored out of the summation and the result is a geometric series. $$ \begin{aligned} Z(k) &= \frac{ 1 }{ N } \sum_{n=0}^{N-1} M e^{ i ( \alpha n + \phi ) } \cdot e^{ -i \beta(k) n } \\ &= \frac{ M }{ N } e^{ i \phi } \sum_{n=0}^{N-1} \left( e^{ i ( \alpha - \beta(k) ) } \right)^n \\ \end{aligned} \tag {5} $$ Applying the geometric summation formula yields a closed form equation for the fractional bin index. $$ Z(k) = \frac{ M }{ N } e^{ i \phi } \cdot \frac{ 1 - e^{ i (\alpha - \beta(k) ) N } }{ 1 - e^{ i(\alpha - \beta(k) ) } } \tag {6} $$ Unlike the DFT case [1], since $k$ can't be assumed to be an integer, $\beta(k)N$ can't be assumed to be a whole multiple of $2\pi$ so it needs to be retained.

Value and Gap

The best results will occur when the bins nearest the frequency of the evaluated tones are used. Three bins are required, so a center oriented location with a plus or minus gap will be used. This makes for three equally spaced bins in the frequency domain. $$ k = v + g \tag {7} $$ Of course, these are readily tranlated onto the angle scale. $$ \beta(k) = \beta(v+g) = \beta(v) + \beta(g) \tag {8} $$ Plugging these into the bin value formula (6), then factoring on the value and gap values. $$ \begin{aligned} Z(v+g) &= \frac{ M }{ N } e^{ i \phi } \cdot \frac{ 1 - e^{ i (\alpha - \beta(v) - \beta(g) ) N } } { 1 - e^{ i(\alpha - \beta(v) - \beta(g) ) } } \\ &= \frac{ M }{ N } e^{ i \phi } \cdot \frac{ 1 - e^{ i (\alpha - \beta(v) ) N } e^{ -i \beta(g) N } } { 1 - e^{ i(\alpha - \beta(v) ) } e^{ -i \beta(g) } } \\ \end{aligned} \tag {9} $$ Substituting in zero or $-g$ for $g$ will give the other two values.

Twist Factors

The factored values, being on the unit circle, can be expressed as complex twist factor values. That is, values on the unit circle which will rotate any complex number you multiply it with by its angle. Their magnitude is one. $$ \begin{aligned} a &= e^{ i(\alpha - \beta(v) ) } = e^{ i \frac{ 2\pi }{ N } ( f - v ) } \\ b &= e^{ -i \beta(g) } = e^{ -i \frac{ 2\pi }{ N } g } \\ \end{aligned} \tag {10} $$ Plugging these in simplifies the bin value equation. $$ Z(v+g) = \frac{ M }{ N } e^{ i \phi } \cdot \frac{ 1 - a^N e^{ -i \beta(g) N } }{ 1 - a b } \tag {11} $$ Cross multiplying by the RHS denominator removes the quotient for simplicity. $$ \left( 1 - a b \right) Z(v+g) = \frac{ M }{ N } e^{ i \phi } \cdot \left( 1 - a^N e^{ -i \beta(g) N } \right) \tag {12} $$

System of Equations

The three values can now be used to build a system of equations. For convenience, three evenly space frequencies are considered. A center one, and plus or minus a gap. $$ \begin{aligned} \left( 1 - a / b \right) Z(v-g) &= \frac{ M }{ N } e^{ i \phi } \cdot \left( 1 - a^N e^{ i \beta(g) N } \right) \\ \left( 1 - a \right) Z(v) &= \frac{ M }{ N } e^{ i \phi } \cdot \left( 1 - a^N \right) \\ \left( 1 - a b \right) Z(v+g) &= \frac{ M }{ N } e^{ i \phi } \cdot \left( 1 - a^N e^{ -i \beta(g) N } \right) \\ \end{aligned} \tag {13} $$ These equations can be combined into a single vector equation. $$ \begin{bmatrix} Z(v-g) \\ Z(v) \\ Z(v+g) \end{bmatrix} - a \begin{bmatrix} 1/b & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & b \\ \end{bmatrix} \begin{bmatrix} Z(v-g) \\ Z(v) \\ Z(v+g) \end{bmatrix} = \frac{M}{N} e^{i\phi} \left( \begin{bmatrix} 1 \\ 1 \\ 1 \end{bmatrix} - a^N \begin{bmatrix} e^{ i \beta(g) N } \\ 1 \\ e^{ -i \beta(g) N } \end{bmatrix} \right) \tag {14} $$ The vectors and matrix can now be given names and the equation becomes more manageable. $$ \vec Z - a D \vec Z = \frac{M}{N} e^{i\phi} \left( \vec U - a^N \vec V \right) \tag {15} $$ At this point, the frequency is represented on the left hand side by $a$, and the amplitude and phase are on the right hand side.

Orthogonal Vector to U and V

Dotting the RHS with a vector that is orthogonal to $\vec U$ and $\vec V$ would zero it out. Since these are dimension three vectors, the standard cross product can be used to find such an orthogonal vector. $$ \begin{aligned} \vec K &= \vec U \times \vec V \\ &= \begin{bmatrix} \vec i & \vec j & \vec k \\ 1 & 1 & 1 \\ e^{ i \beta(g) N } & 1 & e^{ -i \beta(g) N } \\ \end{bmatrix} \\ &= \begin{bmatrix} e^{ -i \beta(g) N } - 1 & -e^{ -i \beta(g) N } + e^{ i \beta(g) N } & 1 - e^{ i \beta(g) N } \\ \end{bmatrix} \\ &= \left( e^{ i \beta(g) N/2 } - e^{ -i \beta(g) N/2 } \right) \begin{bmatrix} -e^{ -i \beta(g) N/2 } & e^{ i \beta(g) N/2 } + e^{ -i \beta(g) N/2 } & -e^{ i \beta(g) N/2 } \\ \end{bmatrix} \\ &= 2i \sin\left( \beta(g) N/2 \right) \vec W \end{aligned} \tag {16} $$ The scalar factor, which goes to zero with $g$, can be disregarded as the requirement is simply for an orthogonal vector. The vector portion can then be decomposed. $$ \begin{aligned} \vec W &= \begin{bmatrix} -e^{ -i \beta(g) N/2 } & e^{ i \beta(g) N/2 } + e^{ -i \beta(g) N/2 } & -e^{ i \beta(g) N/2 } \\ \end{bmatrix} \\ &= \cos( \beta(g) N/2 ) \begin{bmatrix} -1 & 2 & -1 \\ \end{bmatrix} + i \, \sin( \beta(g) N/2 ) \begin{bmatrix} 1 & 0 & -1 \\ \end{bmatrix} \\ \end{aligned} \tag {17} $$

Solve for the Frequency Twist Factor

At this point (15) can be dotted by $\vec W$ to zero out the RHS. $$ \vec W \vec Z - a \vec W D \vec Z = 0 \tag {18} $$ This leaves the twist factor representing the frequency on the LHS which is readily solved for. $$ a = \frac{\vec W \vec Z }{ \vec W D \vec Z } \tag {19} $$ All the terms on the RHS are known. The $\vec W$ and $D$ values comes from the bin selection, and the $\vec Z$ values come from the DTFT taken on the signal in question. This allows the value of $a$ to be calculated.

Solve for the Frequency

Recovering the frequency angle from the frequency twist factor is merely inverting the defintition. $$ \begin{aligned} \alpha - \beta(v) &= \arg( a ) \\ \alpha &= \beta(v) + \arg \left( \frac{\vec W \vec Z }{ [\vec W D] \vec Z } \right) \end{aligned} \tag {20} $$ Recovering the frequency from the frequency angle is just a scale change. $$ \begin{aligned} f &= v + \arg \left( a \right) \cdot \frac{N}{2\pi} \\ &= v + \operatorname{atan2} \left( Im[a], Re[a] \right) \cdot \frac{ N }{ 2\pi } \\ &= v + \arctan \left( \frac{Im[a]}{Re[a]} \right) \cdot \frac{ N }{ 2\pi } \\ \end{aligned} \tag {21} $$ Since $a$ is expected to be very close to one, the angle is going to be very small and there isn't any worry about using the regular inverse Tangent function, rather than atan2.

For a pure tone, in the absence of noise, this answer will be exact. In the case of the presence of white noise, there will be a slight bias and the answer will be very close to the Maximum Likelyhood Estimate (MLE). [Thanks, Royi Avatal]

Conclusion

Another exact analytic formula for tone parameter calculation has been presented. It is very similar to the prior formulas based on the pure DFT, and doesn't really offer any advantages in the limited testing I have done.

(Other articles are linked in my overview article. This one can be found in the "Complex Tones" section.)

References

[1] Dawg, C., "Three Bin Exact Frequency Formulas for a Pure Complex Tone in a DFT", dsprelated.com/showarticle/1043.php
[2] Dawg, C., "The Exponential Nature of the Complex Unit Circle", dsprelated.com/showarticle/754.php
[3] Dawg, C., "A Recipe for a Basic Trigonometry Table", dsprelated.com/showarticle/1471.php

Sign in