Bottom     Previous     Contents

Annex E - Interpolation and Decimation.

Often, time domain processing involves changing the sampling frequency of signals during various processing stages. This is sometimes called 'multi-rate processing', and exploits the fact that any signal only needs a sampling frequency equal to its band width (for real signals we include the negative frequencies in this band width). In particular, if the band width is reduced we can reduce the number crunching burden on subsequent processing stages by a corresponding reduction the sampling frequency. It would be useful to have a way to do this in the frequency domain.

Interpolation.

Conceptually, in the time domain interpolation by a factor M (M is an integer) involves 'stuffing' (M-1) zeros between the original signal samples. In the frequency domain this introduces (M-1) aliased images of the original signal. A filter is then applied to remove the unwanted images, leaving us with an interpolated version of the original signal. (Of course this isn't the way things are actually done in practice, its more usual to use separate filters to generate each interpolated point and then interleave the results. But the maths is the same, either way.)

Suppose we wanted to perform the equivalent operation in the frequency domain. We have the DFT (FFT) of an N point sequence f(n), and we want to generate the DFT (FFT) of an MN point sequence g(n) which is obtained by stuffing (M-1) zeros between the samples of f(n).

So..

equation

equation

The DFT of g(n) is (ignoring the zeros):

equation

F(k) is periodic (with period N), so the DFT of g is simply obtained by catenating M identical (scaled) copies of DFT of f.

We could use this as a practical means of interpolating (by multiplying G(k) by the DFT of the interpolation filter implulse response and taking the inverse transform). In practice this probably isn't the fastest way because it requires an MN point inverse FFT. Implementing M separate filters will only require M N point inverse FFTs, which is slightly faster.

Interpolating Periodic Sequences.

The above result can be used to show how to interpolate if f(n) is known to be periodic (with period N). In this case we can (conceptually) implement an ideal 'brick wall' interpolating filter, simply by padding F(k) with (M-1)N trailing zeros and taking the MN point inverse FFT of the result. Again, in practice it may be more efficient to break this up into M N point inverse FFTs.

The interpolated result is (n= 0 .. N-1, m= 0 .. M-1):

equation

This inverse DFT can be evaluated for each m separately.

Decimation.

Suppose we have the DFT F(k) of an MN point sequence, but we only want N samples from the inverse DFT which are separated by M samples, the remaining (M-1)N samples being discarded. One way to achieve this would be to evaluate the full MN point DFT and discard the samples we don't want, but this wouldn't be very efficient. The method used for mixed radix algorithms provides us a clue how to do better.

The full time domain sequence is given by the inverse DFT (n= 0 .. N-1, m= 0 .. M-1):

equation

Assume (for the sake of simplicity) that we're only interested in values for m=0:

equation

equation

equation

So we can see that the DFT of the decimated sequence is given by:

equation

Fractional Re-Sampling.

Now that we know how to interpolate and decimate by integral factors, we can easily generalise this to change the sample rate of an N point sequence by any factor which is a rational number (I/D), where I and D are mutually prime integers. Basically we interpolate by factor I to generate the DFT of an NI point sequence and then decimate by D, to get the DFT of an (NI/D) point sequence. Obviously this requires that N is a multiple of D (N = MD).

Let the original sequence be f(n) (n= 0 .. MD-1), its DFT is:

equation

This is interpolated by factor I by (conceptually) catenating I copies of F and multiplying by the IMD point DFT of an appropriate interpolating filter H(k):

equation

The result is then decimated by D, to generate the IM (=NI/D) point DFT of the result:

equation

This is the expression that a routine performing this operation will evaluate (for k= 0 .. IM-1):

equation

The 2D equivalent of this operation is one way of achieving arbitrary scaling of images.

©1999 - Engineering Productivity Tools Ltd.


    Top