Digital Signal Processing

By Steven W. Smith, Ph.D.

- 1: The Breadth and Depth of DSP
- 2: Statistics, Probability and Noise
- 3: ADC and DAC
- 4: DSP Software
- 5: Linear Systems
- 6: Convolution
- 7: Properties of Convolution
- 8: The Discrete Fourier Transform
- 9: Applications of the DFT
- 10: Fourier Transform Properties
- 11: Fourier Transform Pairs
- 12: The Fast Fourier Transform
- 13: Continuous Signal Processing
- 14: Introduction to Digital Filters
- 15: Moving Average Filters
- 16: Windowed-Sinc Filters
- 17: Custom Filters
- 18: FFT Convolution
- 19: Recursive Filters
- 20: Chebyshev Filters
- 21: Filter Comparison
- 22: Audio Processing
- 23: Image Formation & Display
- 24: Linear Image Processing
- 25: Special Imaging Techniques
- 26: Neural Networks (and more!)
- 27: Data Compression
- 28: Digital Signal Processors
- 29: Getting Started with DSPs
- 30: Complex Numbers
- 31: The Complex Fourier Transform
- 32: The Laplace Transform
- 33: The z-Transform
- 34: Explaining Benford's Law

Your laser printer will thank you!

Convolution via the Frequency Domain

Suppose that you despise convolution. What are you going to do if given an
input signal and impulse response, and need to find the resulting output signal?
Figure 9-8 provides an answer: transform the two signals into the frequency
domain, multiply them, and then transform the result back into the time domain.
This replaces one convolution with two DFTs, a multiplication, and an Inverse
DFT. Even though the intermediate steps are very different, the output is
*identical* to the standard convolution algorithm.

Does anyone hate convolution enough to go to this trouble? The answer is yes.
Convolution is avoided for two reasons. First, convolution is *mathematically*
difficult to deal with. For instance, suppose you are given a system's impulse
response, and its output signal. How do you calculate what the input signal is?
This is called deconvolution, and is virtually impossible to understand in the
time domain. However, deconvolution can be carried out in the frequency
domain as a simple *division*, the inverse operation of multiplication. The
frequency domain becomes attractive whenever the complexity of the Fourier
Transform is less than the complexity of the convolution. This isn't a matter of
which you like better; it is a matter of which you hate less.

The second reason for avoiding convolution is *computation speed*. For example,
suppose you design a digital filter with a kernel (impulse response) containing
512 samples. Using a 200 MHz personal computer with floating point numbers,
each sample in the output signal requires about one millisecond to calculate,
using the standard convolution algorithm. In other words, the throughput of the
system is only about 1,000 samples per second. This is 40 times too slow for
high-fidelity audio, and 10,000 times too slow for television quality video!

The standard convolution algorithm is slow because of the large number of
multiplications and additions that must be calculated. Unfortunately, simply
bringing the problem into the frequency domain via the DFT doesn't help at all.
Just as many calculations are required to calculate the DFTs, as are required to
directly calculate the convolution. A breakthrough was made in the problem in
the early 1960s when the *Fast Fourier Transform* (FFT) was developed. The
FFT is a clever algorithm for rapidly calculating the DFT. Using the FFT,
convolution by multiplication in the frequency domain can be hundreds of times
faster than conventional convolution. Problems that take hours of calculation
time are reduced to only minutes. This is why people get excited about the FFT,
and processing signals in the frequency domain. The FFT will be presented in
Chapter 12, and the method of FFT convolution in Chapter 18. For now, focus
on how signals are convolved by frequency domain multiplication.

To start, we need to define how to multiply one frequency domain signal by
another, i.e., what it means to write: *X*[*f*] × *H*[*f*] = *Y*[*f*]. In polar form, the magnitudes are multiplied: *MagY*[*f*] = *MagX*[*f*] ? *MagH*[*f*], and the phases are
added: *PhaseY*[*f*] = *PhaseX*[*f*] ? *PhaseH*[*f*]. To understand this, imagine a
cosine wave entering a system with some amplitude and phase. Likewise, the
output signal is also a cosine wave with some amplitude and phase. The polar
form of the frequency response directly describes how the two amplitudes are
related and how the two phases are related.

When frequency domain multiplication is carried out in *rectangular form* there
are cross terms between the real and imaginary parts. For example, a sine wave
entering the system can produce both cosine and sine waves in the output. To
multiply frequency domain signals in rectangular notation:

Focus on understanding multiplication using *polar notation*, and the idea of
cosine waves passing through the system. Then simply accept that these more
elaborate equations result when the same operations are carried out in
rectangular form. For instance, let's look at the *division* of one frequency
domain signal by another. In polar form, the division of frequency domain
signals is achieved by the inverse operations we used for multiplication. To
calculate: *H*[*f*] = *Y*[*f*]/*X*[*f*], divide the magnitudes and subtract the phases, i.e., *MagH*[*f*] = *MagY*[*f*]/*MagX*[*f*], *PhaseH*[*f*] = *PhaseY*[*f*] - *PhaseX*[*f*]. In rectangular form this becomes:

Now back to frequency domain convolution. You may have noticed that we
cheated slightly in Fig. 9-8. Remember, the convolution of an *N* point signal
with an *M* point impulse response results in an *N*+*M*-1 point output signal. We
cheated by making the last part of the input signal all *zeros* to allow this
expansion to occur. Specifically, (a) contains 453 nonzero samples, and (b)
contains 60 nonzero samples. This means the convolution of the two, shown in
(c), can fit comfortably in the 512 points provided.

Now consider the more general case in Fig. 9-9. The input signal, (a), is 256
points long, while the impulse response, (b), contains 51 nonzero points. This
makes the convolution of the two signals 306 samples long, as shown in (c).
The problem is, if we use frequency domain multiplication to perform the
convolution, there are only 256 samples *allowed* in the output signal. In other
words, 256 point DFTs are used to move (a) and (b) into the frequency domain.
After the multiplication, a 256 point Inverse DFT is used to find the output
signal. How do you squeeze 306 values of the correct signal into the 256 points
provided by the frequency domain algorithm? The answer is, you can't! The
256 points end up being a distorted version of the correct signal. This process
is called circular convolution. It is important because you want to *avoid* it.

To understand circular convolution, remember that an *N* point DFT views the
time domain as being an infinitely long periodic signal, with *N* samples per
period. Figure (d) shows three periods of how the DFT views the output signal
in this example. Since *N* = 256, each period consists of 256 points: 0-255, 256-511, and 512-767. Frequency domain convolution tries to place the 306 point
*correct output signal*, shown in (c), into each of these 256 point periods. This
results in 49 of the samples being pushed into the neighboring period to the
right, where they overlap with the samples that are legitimately there. These
overlapping sections add, resulting in each of the periods appearing as shown
in (e), the *circular convolution*.

Once the nature of circular convolution is understood, it is quite easy to avoid.
Simply pad each of the signals being convolved with enough zeros to allow the output signal room to handle the *N*+*M*-1 points in the correct
convolution. For example, the signals in (a) and (b) could be padded with zeros
to make them 512 points long, allowing the use of 512 point DFTs. After the
frequency domain convolution, the output signal would consist of 306 nonzero
samples, plus 206 samples with a value of zero. Chapter 18 explains this
procedure in detail.

Why is it called *circular* convolution? Look back at Fig. 9-9d and examine the
center period, samples 256 to 511. Since all of the periods are the same, the
portion of the signal that flows out of this period to the *right*, is the same that
flows into this period from the *left*. If you only consider a single period, such
as in (e), it *appears* that the right side of the signal is somehow *connected* to the
left side. Imagine a snake biting its own tail; sample 255 is located next to
sample 0, just as sample 100 is located next to sample 101. When a portion of
the signal exits to the right, it magically reappears on the left. In other words,
the *N* point time domain behaves as if it were *circular*.

In the last chapter we posed the question: does it really matter if the DFT's time
domain is viewed as being *N* points, rather than an infinitely long periodic
signal of period *N*? Circular convolution is an example where it *does* matter.
If the time domain signal is understood to be *periodic*, the distortion
encountered in circular convolution can be simply explained as the signal
expanding from one period to the next. In comparison, a rather bizarre
conclusion is reached if only *N* points of the time domain are considered. That
is, frequency domain convolution acts as if the time domain is somehow
wrapping into a circular ring with sample 0 being positioned next to sample *N*-1.