We first define a complete set of orthogonal basis functions by δ ( n − m ) δ ( n − m ) for m = 0 , 1 , 2 , ⋯ , ∞ m = 0 , 1 , 2 , ⋯ , ∞ . The input x ( n ) x ( n ) is broken down into a set of inputs by taking an inner product of the input with each of the basis functions. This produces a set of input components, each of which is a single impulse weighted by a single value of the input sequence ( x ( n ) , δ ( n − m ) ) = x ( m ) δ ( n − m ) ( x ( n ) , δ ( n − m ) ) = x ( m ) δ ( n − m ) . Using the time invariant property of the system, δ ( n − m ) → h ( n − m ) δ ( n − m ) → h ( n − m ) and using the scaling property of a linear system, this gives an output of x ( m ) δ ( n − M ) → x ( m ) h ( n − m ) x ( m ) δ ( n − M ) → x ( m ) h ( n − m ) . We now calculate the output due to x ( n ) x ( n ) by adding outputs due to each of the resolved inputs using the superposition property of linear systems. This is illustrated by the following diagram:

If the system is linear but time varying, we denote the response to an impulse at n = m n = m by δ ( n − m ) → h ( n , m ) δ ( n − m ) → h ( n , m ) . In other words, each impulse response may be different depending on when the impulse is applied. From the development above, it is easy to see where the time-invariant property was used and to derive a convolution equation for a time-varying system as y ( n ) = h ( n , m ) * x ( n ) = ∑ m = − ∞ ∞ h ( n , m ) x ( m ) . y ( n ) = h ( n , m ) * x ( n ) = ∑ m = − ∞ ∞ h ( n , m ) x ( m ) . Unfortunately, relaxing the linear constraint destroys the basic structure of the convolution sum and does not result in anything of this form that is useful.

By a change of variables, one can easily show that the convolution sum can also be written y ( n ) = h ( n ) * x ( n ) = ∑ m = − ∞ ∞ h ( m ) x ( n − m ) . y ( n ) = h ( n ) * x ( n ) = ∑ m = − ∞ ∞ h ( m ) x ( n − m ) . If the system is causal, h ( n ) = 0 h ( n ) = 0 for n < 0 n < 0 and the upper limit on the summation in (2.2) becomes m = n m = n . If the input signal is causal, the lower limit on the summation becomes zero. The form of the convolution sum for a linear, time-invariant, causal discrete-time system with a causal input is y ( n ) = h ( n ) * x ( n ) = ∑ m = 0 n h ( n − m ) x ( m ) y ( n ) = h ( n ) * x ( n ) = ∑ m = 0 n h ( n − m ) x ( m ) or, showing the operations commute y ( n ) = h ( n ) * x ( n ) = ∑ m = 0 n h ( m ) x ( n − m ) . y ( n ) = h ( n ) * x ( n ) = ∑ m = 0 n h ( m ) x ( n − m ) .

Convolution is used analytically analyze linear systems and it can also be used to calculate the output of a system by only knowing its impulse response. This is a very powerful tool because it does not require any detailed knowledge of the system itself. It only uses one experimentally obtainable response. However, this summation cannot only be used to analyze or calculate the response of a given system, it can be an implementation of the system. This summation can be implemented in hardware or programmed on a computer and become the signal processor.

Some of the properties and characteristics of convolution and of the systems it represents can be better described by a matrix formulation than by the summation notation. The first L L values of the discrete-time convolution defined above can be written as a matrix operator on a vector of inputs to give a vector of the output values.

If the input sequence x x is of length N N and the operator signal h h is of length M M , the output is of length L = N + M − 1 L = N + M − 1 . This is shown for N = 4 N = 4 and M = 3 M = 3 by the rectangular matrix operation [ y 0 y 1 y 2 y 3 y 4 y 5 ] = [ h 0 0 0 0 h 1 h 0 0 0 h 2 h 1 h 0 0 0 h 2 h 1 h 0 0 0 h 2 h 1 0 0 0 h 2 ] [ x 0 x 1 x 2 x 3 ] [ y 0 y 1 y 2 y 3 y 4 y 5 ] = [ h 0 0 0 0 h 1 h 0 0 0 h 2 h 1 h 0 0 0 h 2 h 1 h 0 0 0 h 2 h 1 0 0 0 h 2 ] [ x 0 x 1 x 2 x 3 ] It is clear that if the system is causal ( h ( n ) = 0 h ( n ) = 0 for n < 0 n < 0 ), the H H matrix is lower triangular. It is also easy to see that the system being time-invariant is equivalent to the matrix being Toeplitz [1]. This formulation makes it obvious that if a certain output were desired from a length 4 input, only 4 of the 6 values could be specified and the other 2 would be controlled by them.

Although the formulation of constructing the matrix from the impulse response of the system and having it operate on the input vector seems most natural, the matrix could have been formulated from the input and the vector would have been the impulse response. Indeed, this might the appropriate formulation if one were specifying the input and output and designing the system.

The basic convolution defined in (Equation 1), derived in (Equation 2), and given in matrix form in (Equation 3) relates the input to the output for linear systems. This is the form of convolution that is related to multiplication of the DTFT and z-transform of signals. However, it is cyclic convolution that is fundamentally related to the DFT and that will be efficiently calculated by the fast Fourier transform (FFT) developed in Part III of these notes. Matrix formulation of length-L cyclic convolution is given by [ y 0 y 1 y 2 ⋮ y L − 1 ] = [ h 0 h L − 1 h L − 2 ⋯ h 1 h 1 h 0 h L − 1 h 2 h 2 h 1 h 0 h 3 ⋮ ⋮ h L − 1 ⋯ h 0 ] [ x 0 x 1 x 2 ⋮ x L − 1 ] [ y 0 y 1 y 2 ⋮ y L − 1 ] = [ h 0 h L − 1 h L − 2 ⋯ h 1 h 1 h 0 h L − 1 h 2 h 2 h 1 h 0 h 3 ⋮ ⋮ h L − 1 ⋯ h 0 ] [ x 0 x 1 x 2 ⋮ x L − 1 ]

This matrix description makes it clear that the matrix operator is always square and the three signals, x ( n ) x ( n ) , h ( n ) h ( n ) , and y ( n ) y ( n ) , are necessarily of the same length.

There are several useful conclusions that can be drawn from linear algebra [1]. The eigenvalues of the non-cyclic are all the same since the eigenvalues of a lower triangular matrix are simply the values on the diagonal.

Although it is less obvious, the eigenvalues of the cyclic convolution matrix are the N N values of the DFT of h ( n ) h ( n ) and the eigenvectors are the basis functions of the DFT which are the column vectors of the DFT matrix. The eigenvectors are completely controlled by the structure of H H being a cyclic convolution matrix and are not at all a function of the values of h ( n ) h ( n ) . The DFT matrix equation from (3.10) is given by X = F x and Y = F y X = F x and Y = F y where X X is the length-N vector of the DFT values, H H is the matrix operator for the DFT, and x x is the length-N vector of the signal x ( n ) x ( n ) values. The same is true for the comparable terms in y y .

The matrix form of the length-N cyclic convolution in (3.10) is written y = H x y = H x

Taking the DFT both sides and using the IDFT on x x gives F y = Y = F H x = F H F − 1 X F y = Y = F H x = F H F − 1 X If we define the diagonal matrix H d H d as an L L by L L matrix with the values of the DFT of h ( n ) h ( n ) on its diagonal, the convolution property of the DFT becomes Y = H d X Y = H d X This implies H d = F H F − 1 and H = F − 1 H d F H d = F H F − 1 and H = F − 1 H d F which is the basis of the earlier statement that the eigenvalues of the cyclic convolution matrix are the values of the DFT of h ( n ) h ( n ) and the eigenvectors are the orthogonal columns of F F . The DFT matrix diagonalizes the cyclic convolution matrix. This is probably the most concise statement of the relation of the DFT to convolution and to linear systems.

An important practical question is how one calculates the non-cyclic convolution needed by system analysis using the cyclic convolution of the DFT. The answer is easy to see using the matrix description of H H . The length of the output of non-cyclic convolution is N + M − 1 N + M − 1 . If N − 1 N − 1 zeros are appended to the end of h ( n ) h ( n ) and M − 1 M − 1 zeros are appended to the end of x ( n ) x ( n ) , the cyclic convolution of these two augmented signals will produce exactly the same N + M − 1 N + M − 1 values as non-cyclic convolution would. This is illustrated for the example considered before. [ y 0 y 1 y 2 y 3 y 4 y 5 ] = [ h 0 0 0 0 h 2 h 1 h 1 h 0 0 0 0 h 2 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 ] [ x 0 x 1 x 2 x 3 0 0 ] [ y 0 y 1 y 2 y 3 y 4 y 5 ] = [ h 0 0 0 0 h 2 h 1 h 1 h 0 0 0 0 h 2 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 0 0 0 0 h 2 h 1 h 0 ] [ x 0 x 1 x 2 x 3 0 0 ] Just enough zeros were appended so that the nonzero terms in the upper right-hand corner of H H are multiplied by the zeros in the lower part of x x and, therefore, do not contribute to y y . This does require convolving longer signals but the output is exactly what we want and we calculated it with the DFT-compatible cyclic convolution. Note that more zeros could have been appended to h h and x x and the first N + M − 1 N + M − 1 terms of the output would have been the same only more calculations would have been necessary. This is sometimes done in order to use forms of the FFT that require that the length be a power of two.

If fewer zeros or none had been appended to h h and x x , the nonzero terms in the upper right-hand corner of H H , which are the ``tail" of h ( n ) h ( n ) , would have added the values that would have been at the end of the non-cyclic output of y ( n ) y ( n ) to the values at the beginning. This is a natural part of cyclic convolution but is destructive if non-cyclic convolution is desired and is called aliasing or folding for obvious reasons. Aliasing is a phenomenon that occurs in several arenas of DSP and the matrix formulation makes it easy to understand.