Main Content

This example shows how to design filters for decimation and interpolation. Typically lowpass filters are used for decimation and for interpolation. When decimating, lowpass filters are used to reduce the bandwidth of a signal prior to reducing the sampling rate. This is done to minimize aliasing due to the reduction in the sampling rate. When interpolating, lowpass filters are used to remove spectral images from the low-rate signal. For general notes on lowpass filter design see the example on Designing Low Pass FIR Filters

To start consider changing the rate of a signal by a factor of 2. Halfband filters are an efficient way of doing this. Halfband FIR filters are implemented in `dsp.FIRHalfbandInterpolator`

and `dsp.FIRHalfbandDecimator`

. Their IIR counterparts, `dsp.IIRHalfbandInterpolator`

and `dsp.IIRHalfbandDecimator`

can be an even more efficient way of interpolating/decimating by 2.

Fs = 1e6; hbInterp = dsp.FIRHalfbandInterpolator('TransitionWidth',Fs/10,... 'SampleRate',Fs); fvtool(hbInterp) % Notice gain of 2 (6 dB) in the passband hbDecim = dsp.FIRHalfbandDecimator('TransitionWidth',Fs/10,... 'SampleRate',Fs); fvtool(hbDecim)

The sampling rate Fs refers to the input signal. In the case of interpolation, the filter retains most of the spectrum from 0 to Fs/2 while attenuating spectral images. For decimation, the filter passes about half of the band, that is 0 to Fs/4, and attenuates the other half in order to minimize aliasing. The amount of attenuation can be set to any desired value for both interpolation and decimation. If unspecified, it defaults to 80 dB.

Halfband filters can be cascaded for efficient multistage rate conversion. For example, interpolating/decimating by 8 can be done by cascading 3 halfband interpolators/decimators. Rather than designing the 3 filters by hand, `dsp.SampleRateConverter`

will design all 3 filters in a very efficient way. Note that `dsp.SampleRateConverter`

can be used for designs that go beyond using halfband filters. For more details on this, see Efficient Sample Rate Conversion Between Arbitrary Factors.

src = dsp.SampleRateConverter('InputSampleRate',10e3,... 'OutputSampleRate',8*10e3,'Bandwidth',9e3); info(src)

ans = 'Overall Interpolation Factor : 8 Overall Decimation Factor : 1 Number of Filters : 3 Multiplications per Input Sample: 96.000000 Number of Coefficients : 66 Filters: Filter 1: dsp.FIRInterpolator - Interpolation Factor: 2 Filter 2: dsp.FIRInterpolator - Interpolation Factor: 2 Filter 3: dsp.FIRInterpolator - Interpolation Factor: 2 '

Sometimes it is desirable to design a filter for changing the rate by a rational factor regardless of the actual sampling frequencies involved. `designMultirateFIR(L,M)`

designs an FIR filter for interpolation by an integer factor L and decimation by an integer factor M. `designMultirateFIR`

returns filter coefficients. These coefficients are to be used with `dsp.FIRDecimator`

(L=1), `dsp.FIRInterpolator`

(M=1), and `dsp.FIRRateConverter`

(general case).

```
L = 5;
M = 1;
b = designMultirateFIR(L,M);
firInterp = dsp.FIRInterpolator(L,b);
fvtool(firInterp) % Notice gain of L in the passband
```

Optional inputs to `designMultirateFIR`

allow for steeper transitions and better stopband attenuation. For a steeper transition, specify a larger half polyphase length (the default is 12).

hpl = 20; b2 = designMultirateFIR(L,M,hpl); firInterp2 = dsp.FIRInterpolator(L,b2); h = fvtool(firInterp,firInterp2); legend(h, 'Polyphase length = 24','Polyphase length = 40')

When decimating, the bandwidth of a signal is reduced to an appropriate value so that minimal aliasing occurs when reducing the sampling rate. Suppose a signal that occupies the full Nyquist interval (i.e. has been critically sampled) has a sampling rate of 800 Hz. The signal's energy extends up to 400 Hz. If we'd like to reduce the sampling rate by a factor of 4 to 200 Hz, significant aliasing will occur unless the bandwidth of the signal is also reduced by a factor of 4. Ideally, a perfect lowpass filter with a cutoff at 100 Hz would be used. In practice, several things will occur: The signal's components between 0 and 100 Hz will be slightly distorted by the passband ripple of a non-ideal lowpass filter; there will be some aliasing due to the finite stopband attenuation of the filter; the filter will have a transition band which will distort the signal in such band. The amount of distortion introduced by each of these effects can be controlled by designing an appropriate filter. In general, to obtain a better filter, a higher filter order will be required.

Let's start by designing a simple lowpass decimator with a decimation factor of 4.

M = 4; % Decimation factor Fp = 80; % Passband-edge frequency Fst = 100; % Stopband-edge frequency Ap = 0.1; % Passband peak-to-peak ripple Ast = 80; % Minimum stopband attenuation Fs = 800; % Sampling frequency fdDecim = fdesign.decimator(M,'lowpass',Fp,Fst,Ap,Ast,Fs) %#ok

fdDecim = decimator with properties: MultirateType: 'Decimator' Response: 'Lowpass' DecimationFactor: 4 Specification: 'Fp,Fst,Ap,Ast' Description: {4x1 cell} NormalizedFrequency: 0 Fs: 800 Fs_in: 800 Fs_out: 200 Fpass: 80 Fstop: 100 Apass: 0.1000 Astop: 80

The specifications for the filter determine that a transition band of 20 Hz is acceptable between 80 and 100 Hz and that the minimum attenuation for out of band components is 80 dB. Also that the maximum distortion for the components of interest is 0.05 dB (half the peak-to-peak passband ripple). An equiripple filter that meets these specs can be easily obtained as follows:

decimFilter = design(fdDecim,'equiripple', 'SystemObject', true); measure(decimFilter) specScope1 = dsp.SpectrumAnalyzer(... 'PlotAsTwoSidedSpectrum', false, ... 'SpectralAverages', 50, 'OverlapPercent', 50, ... 'Title', 'Decimator with equiripple lowpass filter',... 'YLimits', [-50, 0], 'SampleRate', Fs/M*2); for k = 1:1000 inputSig = randn(500,1); % Input decimatedSig = decimFilter(inputSig); % Decimator specScope1(upsample(decimatedSig,2)); % Spectrum % The upsampling is done to increase X-limits of SpectrumAnalyzer % beyond (1/M)*Fs/2 for better visualization end release(specScope1);

ans = Sample Rate : 800 Hz Passband Edge : 80 Hz 3-dB Point : 85.621 Hz 6-dB Point : 87.8492 Hz Stopband Edge : 100 Hz Passband Ripple : 0.092414 dB Stopband Atten. : 80.3135 dB Transition Width : 20 Hz

It is clear from the measurements that the design meets the specs.

Nyquist filters are attractive for decimation and interpolation due to the fact that a 1/M fraction of the number of coefficients is zero. The band of the Nyquist filter is typically set to be equal to the decimation factor, this centers the cutoff frequency at (1/M)*Fs/2. Note that halfband filters are Nyquist filters for the case M=2. Note also that `designMultirateFIR`

designs Nyquist filters as well.

In this example, the transition band is centered around 400/M = 100 Hz.

TW = 20; % Transition width of 20 Hz fdNyqDecim = fdesign.decimator(M,'nyquist',M,TW,Ast,Fs) %#ok

fdNyqDecim = decimator with properties: MultirateType: 'Decimator' Response: 'Nyquist' DecimationFactor: 4 Specification: 'TW,Ast' Description: {2x1 cell} Band: 4 NormalizedFrequency: 0 Fs: 800 Fs_in: 800 Fs_out: 200 TransitionWidth: 20 Astop: 80

A Kaiser window design can be obtained in a straightforward manner.

nyqDecim = design(fdNyqDecim,'kaiserwin','SystemObject', true); specScope2 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ... 'SpectralAverages', 50, 'OverlapPercent', 50, ... 'Title', 'Decimator with Nyquist filter',... 'YLimits', [-50, 0],... 'SampleRate', Fs/M*2); for k = 1:1000 inputSig = randn(500,1); % Input decimatedSig = nyqDecim(inputSig); % Decimator specScope2(upsample(decimatedSig,2)); % Spectrum % The upsampling is done to increase X-limits of SpectrumAnalyzer % beyond (1/M)*Fs/2 for better visualization end release(specScope2);

Suppose the signal to be filtered has a flat spectrum. Once filtered, it acquires the spectral shape of the filter. After reducing the sampling rate, this spectrum is repeated with replicas centered around multiples of the new lower sampling frequency. An illustration of the spectrum of the decimated signal can be found from:

NFFT = 4096; [H,f] = freqz(nyqDecim,NFFT,'whole',Fs); figure; plot(f-Fs/2,20*log10(abs(fftshift(H)))) grid on hold on plot(f-Fs/M,20*log10(abs(fftshift(H))),'r-') plot(f-Fs/2-Fs/M,20*log10(abs(fftshift(H))),'k-') legend('Baseband spectrum', ... 'First positive replica', 'First negative replica') title('Aliasing with Nyquist filter'); fig = gcf; fig.Color = 'White'; hold off

Note that the replicas overlap somewhat, so aliasing is introduced. However, the aliasing only occurs in the transition band. That is, significant energy (above the prescribed 80 dB) from the first replica only aliases into the baseband between 90 and 100 Hz. Since the filter was transitioning in this region anyway, the signal has been distorted in that band and aliasing there is not important.

On the other hand, notice that although we have used the same transition width as with the lowpass design from above, we actually retain a wider usable band (90 Hz rather than 80) when comparing this Nyquist design with the original lowpass design. To illustrate this, let's follow the same procedure to plot the spectrum of the decimated signal when the lowpass design from above is used

[H,f] = freqz(decimFilter,NFFT,'whole',Fs); figure; plot(f-Fs/2,20*log10(abs(fftshift(H)))) grid on hold on plot(f-Fs/M,20*log10(abs(fftshift(H))),'r-') plot(f-Fs/2-Fs/M,20*log10(abs(fftshift(H))),'k-') legend('Baseband spectrum', ... 'First positive replica', 'First negative replica') title('Aliasing with lowpass filter'); fig = gcf; fig.Color = 'White'; hold off

In this case, there is no significant overlap (above 80 dB) between replicas, however because the transition region started at 80 Hz, the resulting decimated signal has a smaller usable bandwidth.

When interpolating a signal, the baseband response of the signal should be left as unaltered as possible. Interpolation is obtained by removing spectral replicas when the sampling rate is increased.

Suppose we have a signal sampled at 48 Hz. If it is critically sampled, there is significant energy in the signal up to 24 Hz. If we want to interpolate by a factor of 4, we would ideally design a lowpass filter running at 192 Hz with a cutoff at 24 Hz. As with decimation, in practice an acceptable transition width needs to be incorporated into the design of the lowpass filter used for interpolation along with passband ripple and a finite stopband attenuation. For example, consider the following specs:

L = 4; % Interpolation factor Fp = 22; % Passband-edge frequency Fst = 24; % Stopband-edge frequency Ap = 0.1; % Passband peak-to-peak ripple Ast = 80; % Minimum stopband attenuation Fs = 48; % Sampling frequency fdInterp = fdesign.interpolator(L,'lowpass',Fp,Fst,Ap,Ast,Fs*L) %#ok

fdInterp = interpolator with properties: MultirateType: 'Interpolator' Response: 'Lowpass' InterpolationFactor: 4 Specification: 'Fp,Fst,Ap,Ast' Description: {4x1 cell} NormalizedFrequency: 0 Fs: 192 Fs_in: 48 Fs_out: 192 Fpass: 22 Fstop: 24 Apass: 0.1000 Astop: 80

An equiripple design that meets the specs can be found in the same manner as with decimators

EQRInterp = design(fdInterp,'equiripple','SystemObject', true); specScope4 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ... 'SpectralAverages', 50, 'OverlapPercent', 50, ... 'Title', 'Interpolator with equiripple lowpass filter', ... 'SampleRate', Fs*L); for k = 1:1000 inputSig = randn(500,1); % Input intrpSig = EQRInterp(inputSig); % Interpolator specScope4(intrpSig); % Spectrum end release(specScope4);

Notice that the filter has a gain of 6 dBm. In general interpolators will have a gain equal to the interpolation factor. This is needed for the signal being interpolated to maintain the same range after interpolation. For example,

release(EQRInterp); sine = dsp.SineWave('Frequency', 18, 'SampleRate', Fs, ... 'SamplesPerFrame', 100); intrpSig = EQRInterp(sine()); arrayPlot = dsp.ArrayPlot('YLimits', [-2, 2], ... 'Title', 'Sine wave interpolated'); arrayPlot(intrpSig(200:300)) % Plot the output

Note that although the filter has a gain of 4, the interpolated signal has the same amplitude as the original signal.

Similar to the decimation case, Nyquist filters are attractive for interpolation purposes. Moreover, given that there is a coefficient equal to zero every L samples, the use of Nyquist filters ensures that the samples from the input signal are retained unaltered at the output. This is not the case for other lowpass filters when used for interpolation (on the other hand, distortion may be minimal in other filters, so this is not necessarily a huge deal).

TW = 2; fdNyqIntrp = fdesign.interpolator(L,'nyquist',L,TW,Ast,Fs*L) %#ok nyqInterp = design(fdNyqIntrp,'kaiserwin', 'SystemObject', true); specScope5 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ... 'SpectralAverages', 30, 'OverlapPercent', 50, ... 'Title', 'Interpolator with Nyquist filter',... 'SampleRate', Fs*L); for k = 1:1000 inputSig = randn(500,1); % Input intrpSig = nyqInterp(inputSig); % Decimator specScope5(intrpSig); % Spectrum end release(specScope5);

fdNyqIntrp = interpolator with properties: MultirateType: 'Interpolator' Response: 'Nyquist' InterpolationFactor: 4 Specification: 'TW,Ast' Description: {2x1 cell} Band: 4 NormalizedFrequency: 0 Fs: 192 Fs_in: 48 Fs_out: 192 TransitionWidth: 2 Astop: 80

In an analogous manner to decimation, when used for interpolation, Nyquist filters allow some degree of imaging. That is, some frequencies above the cutoff frequency are not attenuated by the value of Ast. However, this occurs only in the transition band of the filter. On the other hand, once again a wider portion of the baseband of the original signal is maintained intact when compared to a lowpass filter with stopband-edge at the ideal cutoff frequency when both filters have the same transition width.