The DFT Magnitude of a Real-valued Cosine Sequence

This blog may seem a bit trivial to some readers here but, then again, it might be of some value to DSP beginners. It presents a mathematical proof of what is the magnitude of an N-point discrete Fourier transform (DFT) when the DFT's input is a real-valued sinusoidal sequence.

To be specific, if we perform an N-point DFT on N real-valued time-domain samples of a discrete cosine wave, having exactly integer k cycles over N time samples, the peak magnitude of the cosine wave's positive-frequency spectral component will be

where A is the peak amplitude of the discrete cosine sequence. Graphically, Eq. (1) is shown in Figure 1.

I presented the equivalent of Eq. (1) in the 3rd edition of my DSP textbook. However, I recently received an e-mail from a reader who asked me:

"(I can see that) |X(k)|=AN/2 indeed correct, but
I cannot see mathematically where it comes from!
Such a simple, intuitive thing, but believe it or
not, even after having gone through your book,
I just can't figure out where this is derived from!
Sometimes the obvious things are made difficult by
the reader himself!

After reading that e-mail I thought, "Didn't I give the derivation of Eq. (1) in my book?" The answer is, shame on me, "No." Out of curiosity I checked four other DSP textbooks from my bookshelf and they chose not to provide Eq. (1) at all. How strange! I believe Eq. (1) is important and later in this blog I explain why. In any case, below is the mathematical derivation of Eq. (1).

Deriving Equation (1)

We start with a discrete time-domain cosine sequence x(n) as:

where A is the peak value of the cosine wave, k is the integer number of complete sinusoidal cycles occurring in the N samples, and variable n is the time index. We can show the desired N-point DFT of x(n), X(m), as:

where m is the frequency-domain index.

OK. Remembering Euler's relationship of cos(α) = (e^jα + e^-jα)/2, we can write X(m) as:

The left summation in Eq. (4) represents a positive-frequency spectral component, and the right summation represents a negative-frequency spectral component. In deriving Eq. (1) we need concentrate only on the positivefrequency summation which is:

Equation (5) is a geometric series of the form:

where variable q = -j2π(m-k)/N is a simple dummy variable we’ll use to keep the next equation simple. For m ≠ k Eq. (6) can be written in closed-form as:

(Equation (7) is a very important equation in the field of DSP, however it's typically presented as a simple unjustified fact in DSP books. If you want to see the full derivation of Eq. (7), see reference [1] or reference [2].)

Substituting q = -j2π(m-k)/N into the right side of Eq. (7), we rewrite Eq. (5) as:

OK, we have what we want. Equation (8) is a closed-form expression for the positive-frequency DFT of a real-valued input cosine sequence. (We could perform the algebraic acrobatics to convert Eq. (8) into a familiar sin(x)/x form, but we need not do that here.)

With the original DFT input being exactly integer k cycles of a cosine sequence, to verify Eq. (1) we evaluate Eq. (8) when the frequency index is m = k. Doing that we run into trouble because spectral sample X(m=k) is:

My original solution to the Eq. (9) problem was to use l'Hôpital's Rule to evaluate Eq. (8) when m = k.

That solution, while valid, was overly complicated and unnecessary. As DSP guru Prof. Dilip Sarwate [3] showed me, a simpler solution is to merely evaluate Eq. (6) for q = 0 (m = k in Eq. (5)) as

Why Equation (1) is Important

Beyond the facts that knowing Eq. (1) will make you live longer, improve your sex life, and reduce global warming, Eq. (1) is important when computing fast Fourier transforms (FFTs) on fixed-point binary hardware. When you're computing FFT magnitude samples using binary two's-complement number representation, you want to avoid numerical overflow errors. So this means, for example, if you're using 16-bit hardware you want no computed FFT magnitude sample value to be greater than 2¹⁵.

Exploring this idea, let's assume that the input data to an N-point FFT, computed using B-bit hardware (B = 16 in the above paragraph), is b-bits in word width. That means the maximum positive peak amplitude value of an input cosine wave is A_max = 2^b-1. From our Eq. (1) we can write:

However, to avoid overflow using two's-complement arithmetic:

Rearranging Eq. (14) we have:

Sorry for all the algebra rigmarole. What Eq. (14) is telling is, based on Eq. (1), there is a maximum-sized FFT you can perform on fixed-point hardware while avoiding numerical overflow. For example, on B = 16-bit hardware and b = 8-bit input samples, the largest sized FFT you can compute to determine spectral magnitude values while avoid numerical overflow is

That Nmax = 512 value is shown in Figure 2 above the b = 8 x-axis value.

So this means if we perform 1024- or 2014-point FFTs, with a 16-bit machine accepting 8-bit input samples, it's possible that numerical overflow will occur.

Now, during a System Design Review meeting, if someone suggests performing 2048-point FFTs using 16-bit hardware on 8-bit data, you can pipe up with Eq. (16), derived from Eq. (1), to warn your colleagues of potential numerical overflow. And then tell them that some sort of FFT-internal data-checking and conditional truncation, or multi-word processing, will be necessary to perform those 2048-point FFTs. (Think of the increased software development cost. Yikes!) Or else you could just recommend using floating-point hardware.

A Few Final Thoughts

To determine the negative-frequency X(N-k) DFT sample's magnitude, we keep in mind that due to the circular nature of the DFT, X(N-k) = X*(k) = AN/2, where the '*' symbol means complex conjugate. (See reference [4].)

If our N-point DFT's input, in Eq. (3), had been a sine wave sequence,

the above derivation method, using Euler's relationship of sin(α) = (e^jα - e^-jα)/j2, would produce the same positive-frequency result of X(k) = AN/2.

One last thought from me, and it's a criticism. I stated that I couldn't find a derivation of Eq. (1) in the DSP books on my bookshelf. Well, the closest match I found was in a popular college DSP textbook. In that book the authors have the homework problem:

"Compute the DFT of x(n) = cos(2πnk_o/N), where 0 ≤ n ≤ N-1"

In their book's Solution Manual (I own a legal copy), the authors' solution consisted of nothing more than the following two lines:

Now you tell me, ...how can anyone expect a first-semester DSP student to go from Eq. (17), by inspection in one step, to Eq. (18)? I found no discussion in that book's DFT chapter that would show the student how to go from Eq. (17) to the delta-function expression in Eq. (18). It's no wonder why DSP can be so painful to learn in college classrooms!

References
[1] R. Lyons, “Understanding Digital Signal Processing”, 3rd Edition, Prentice Hall Publishing, 2011, Appendix B.
[2] http://en.wikipedia.org/wiki/Geometric_series
[3] http://en.wikipedia.org/wiki/L%27H%C3%B4pital%27s_rule
[4] J. Proakis and D. Manolakis, Digital Signal Processing- Principles, Algorithms, and Applications, 3rd Ed., Prentice Hall, 1996, pp. 413.

Sign in