Digital Signal Processing Lab Manual: Subject Code: ECE 3161

Download as pdf or txt
Download as pdf or txt
You are on page 1of 45

DIGITAL SIGNAL PROCESSING

LAB MANUAL
Subject Code: ECE 3161
For

5th Semester B. Tech. ECE

DEPARTMENT OF
ELECTRONICS & COMMUNICATION ENGINNEERING

Prepared by Approved by

Staff, ECE Department, M.I.T., Manipal.

Date: 01-12-2020

1
Instructions to the Students:

1. All the students are required to come prepared for the experiments to be done in the
lab.
2. Students should try to analyze and understand the solved problems and then try to
solve the unsolved problems of the experiment in the lab.
3. Maintaining an observation copy is compulsory for all, where the programs (codes)
and results of the unsolved-problems (exercise) are to be noted down.
4. Students have to get their results verified and observation copies checked by the
instructor before leaving the lab.
5. Maintain a separate folder in the computer you use, where you save all the programs
you do in the lab.
6. Use of external storage media during lab is not allowed.
7. Maintain the timings and the discipline of the lab.
8. Students will be evaluated in every lab based on individual performance, observation
copy and involvement in the lab.

Evaluation plan

 In-semester Assessment Marks : 60% (60 Marks)


 Continuous evaluation component (lab sessions 1 to 6): 10 marks each.
 Assessment is based on conduction of each experiment, exercise problems,
answering the questions related to the experiment.
 End semester assessment: 40 % (40 marks)
 write up (coding, analysis): 20 marks
 Concept/Theory: 10 marks
 Viva-Voce: 10 marks

2
Lab Session 1: Working with signals in time domain

Objective:

To study signals and systems in time-domain and understand the


concepts of sampling, convolution, correlation and impulse response.

Signal generation and sampling

The signals we use in the real world (such as voice signal) are "analog"
signals. To process these signals in computing devices, we need to
convert them in to "digital" form. While an analog signal is continuous
in both time and amplitude, a digital signal is discrete in both time and
amplitude. To convert a signal from continuous time to discrete time,
sampling is used. The value of the signal is taken at certain intervals
in time to get discrete samples of the signal.

The Sampling Theorem states that a signal can be exactly reproduced


if it is sampled at a frequency Fs, where Fs is equal or greater than twice
the maximum frequency in the signal (Nyquist rate). The sampling
frequency or sampling rate, Fs, is the average number of samples
obtained in one second (samples per second), thus Fs = 1/T. As the
sampling frequency decreases, the signal separation also decreases.
When the sampling frequency drops below the Nyquist rate, the
frequencies will crossover and cause aliasing. In MATLAB signals can
be simulated using functions where the duration and the sampling time
can be specified.
Lab Session 2: Analysis in Frequency Domain
Lab Session 3: Analysis in Z-Domain
Lab Session 4 : IIR Filter Design - 1

Objective: To design and study analog infinite-duration impulse response (IIR) filters.

Theory: Digital IIR filters are mainly designed from analog filters (classical techniques).
The widely used IIR filters are Butterworth, Chebyshev l (Type I and II), and Elliptic.

Butterworth filters:

The magnitude response low pass Butterworth filter is

| ( Ω)| =

Here N is the order of the filter and Ωc = 3dB cut-off frequency (rad/sec).

A typical magnitude response plot is shown below

Figure 4.1 Magnitude response of Butterworth LPF

In the above figure, Rp and Rs represent the pass band and stop band attenuations,
respectively, and Ω and Ω are the corresponding edge frequencies in rad/sec.

Chebyshev: A Chebyshev filter with ripples in pass-band is categorized as type-1,


whereas, the one with ripples in the stop-band comes under type-2. The figures shown
below are examples of type-I filters.
Figure 4.2 Magnitude response of type-1 Chebyshev filter with N odd

Figure 4.3 Magnitude response of type-1 Chebyshev filter with N even Elliptic:
They have ripples in both pass band and stop-band.
Figure 4.4 Magnitude response of elliptic filter

Implementation in MATLAB:

The following functions are available in MATLAB to design an analog low pass
filters.

1. [N,Wc] = buttord(wp,ws,Rp,Rs,'s') finds the minimum order N and cutoff


frequency Wc for an analog low pass Butterworth filter where wp and ws are the
edge frequencies in radians per second.
2. [b,a] = butter(N,Wc,'s') returns the transfer function coefficients of an Nth order
analog lowpass Butterworth filter with cutoff angular frequency Wc.
3. buttap(N) designs a unity cutoff LPF (these are referred to as prototype filters)
4. Similarly,

cheb1ord/cheby1,

cheb2ord/cheby2,

ellipord/ellip functions can be used to design Chebyshev lowpass (Type 1 and


2), and Elliptic lowpass respectively.

Note: Refer to MATLAB documentation for more information.


Example-1: Design an analog Butterworth low-pass prototype filter of order 5. Use it
to get an unnormalized Butterworth filter with cutoff frequency of 0.2π rad/sec.
Solution:
clc; clear all; close all;
N = 5;

% given order of the filter

omegac = 0.2*pi;

% Given cutoff frequency

[z, p, k ] = buttap(N);

% Getting the prototype


p1 = p* omegac;

k1 = k* omegac ^N;

B = real(poly(z));

b0 = k1;

b = k1*B;
a = real(poly(p1));
Solution:
b=

0.0979

a= 1.0000 2.033 2.0671 1.2988 0.5044 0.0979


Example-2: Design an analog Butterworth low-pass filter, given the following
specifications and also plot the frequency response. Rp = 1dB, Rs = 40dB, fp = 100Hz
and fs = 500Hz

Solution:

clc;
clear all;
close all;
rp = 1;

%Pass-band attenuation
rs = 40;

%Stop-band attenuation
wp = 2*pi*100;

%Frequency should be in rad/sec


ws = 2*pi*500;
[N, wc] = buttord(wp,ws,rp,rs,'s');
%Getting the order N and the 3dB cut-off freq
[b, a] = butter(N,wc,'s');

%Getting the numerator & denominator polynomial coeffs


[H, W] = freqs(b,a); %Complex freq.response of analog filter
plot(W/(2*pi),20*log10(abs(H))), grid on;
title('Frequency Response');
xlabel(' Frequency in Hz '),ylabel(' Magnitude in dB ');
disp(' the order and cut off frequency in Hz are '),disp( N),disp(
wc/(2*pi));
Solution:
Order: 4
Cut-off frequency(Hz): 158.119
Frequency in Hz

Figure 4.5 Frequency response

Exercise

1. Given the following filter specifications: Passband edge: 200Hz; Rp=


1dB; Stopband edge: 700Hz; Rs= 50 dB.
a. Design an analog lowpass Butterworth filter
b. Design an analog lowpass Chebyshev Type 1 filter
c. Design an analog lowpass Chebyshev Type 2 filter
d. Design an analog lowpass Elliptic filter

Plot the frequency response for each case and compare them.

2. Given the following filter specifications: Passband edge: 1kHz; Rp= 1dB;
Stopband edge: 200Hz; Rs= 50 dB.
a. Design an analog high pass Butterworth filter
b. Design an analog high pass Chebyshev Type 1 filter

Note: Refer “butter” function to design high pass and other filters.
Lab session 5: IIR Filter Design – 2

Objective: To study analog-to-digital filter transformation methods for IIR Filter


design.

Theory: The analog filters designed in the previous experiment can be transformed into
digital filter using impulse invariance or bilinear transformation (BLT) methods.

Implementation details: The following functions are available in MATLAB:

1. [bz, az] = impinvar(b, a, Fs) creates a digital filter with numerator and
denominator coefficients bz and az, respectively, whose impulse response is
equal to the impulse response of the analog filter with coefficients b and a, scaled
by 1/Fs (Fs = sampling frequency).
2. [numd, dend] = bilinear(num, den, Fs) convert an s-domain transfer function
given by num and den to a discrete equivalent using BLT method.

Refer to MATLAB documentation for more information.

Note: The functions buttord and butter can also be used to design a digital filter directly.
For example, the following function returns the order N and the cutoff frequency of
digital filter with specifications, wp, ws, rp, rs.

[N, wc] = buttord(wp, ws, rp, rs);

Here rp and rs must be in dB & wp and ws must be normalized to π or half the sampling
frequency. (e.g. If fp is 100Hz, then wp = 2*fp/Fs). Therefore, the range for wp and ws
is always 0 to 1.

3. [b, a]=butter(N, wc) - By default the function butter(N, wc) uses


bilinear transformation with frequency prewarping.

4. Similarly,

cheb1ord/cheby1,

cheb2ord/cheby2,
ellipord/ellip functions can be used to design Chebyshev lowpass (Type 1 and
2), and Elliptic digital IIR filters respectively.

Refer to MATLAB documentation for more information.

Example-1: Design a low-pass digital filter using an analog Butterworth filter to satisfy
= 0.2 , = 0.3 , = 1 , = 15 .
Use the impulse invariance technique to digitize (use impinvar function).

Solution:

clc, clear all , close all;


% Digital Filter specification
wp = 0.2*pi; ws = 0.3*pi; rp = 1; rs = 15;

% Analog design
T = 1;

% assumed
wa_p = wp/T;
wa_s = ws/T;

% Butterworth filter design


[Nb, wcb] = buttord(wa_p,wa_s, rp, rs, 's');
[bb_s,ab_s] = butter(Nb, wcb, 's');
disp('For Butterworth, the order and cut off frequency are '),
disp( Nb),disp(wcb);
[Hb,wb] = freqs(bb_s,ab_s);

% Analog freq. response

% Transformation from analog to digital


[bb_z, ab_z] = impinvar(bb_s,ab_s,1/T);
[Hbd, wbd] = freqz(bb_z,ab_z);

% Digital freq. response figure(1);


plot(wbd/pi,20*log10(abs(Hbd))),grid on;
ylabel('Magnitude in dB' ) ;
xlabel('Frequency in pi units’) ;

Figure 5.1 Frequency response

Example-2: Simulate a signal with 50Hz and 200 Hz frequency. Filter the signal through
a 3rd order digital Butterworth LPF with 3-dB cutoff frequency 100 Hz. Plot the spectrum
of the input and the filtered signal.
Solution:

clear all; clc ; close all

Fs = 1000; % Sampling rate


t = 0:1/(Fs):0.5;
wc = 2*100/Fs

%Normalized 3 dB cutoff freq. for the filter


x = sin(2*pi*t*50) + sin(2*pi*t*200);
[b,a] = butter(3,wc);

%Getting the numerator & denominator polynomial coeffs


[H,W] = freqz(b,a,Fs);

%Complex freq. response of analog filter


figure(1);
k=W*Fs/(2*pi);

% defining frequency vector for plotting


plot(k,20*log10(abs(H))), grid on; title('Frequency Response');
xlabel(' Frequency in Hz '),ylabel(' Magnitude in dB ');

% Filter testing
yb = filter(b,a,x) ;
% Getting spectrum and plot
N = 512;
w = [0:N/2 - 1]* (Fs/N) ;

%defining frequency vector for plotting


X = fft(x,N);
Yb = fft(yb,N);
figure(2);
subplot(211), plot(w,abs(X(1:N/2)));
title('Spect rum of input signal ' ) ;
subplot(212), plot(w,abs(Yb(1:N/2)));
title('Spectrum of filtered signal ' ) ;
Figure 5.2 Frequency response

Figure 5.3 Input and output for example-2


Exercise:

1. Transform

into a digital filter H(z) using the impulse invariance technique in


which T = 0.1s. First find the solution analytically and then verify
your result with the help of a MATLAB program. Also plot the
magnitude responses of the analog and the resulting digital filter and
comment on the result.

2. Design a low-pass digital filter using a


i)

Chebyshev-I analog filter ii)

Chebyshev-II analog filter iii)

Elliptic analog filter


to satisfy = 0.2 , = 0.3 , =1 , = 15 .
Use the impulse invariance technique.

3. Transform the system function given in Question 1 above using the bilinear
transformation assuming T = 1s.

4. Repeat Question 2 using bilinear transformation.

5. Certain digital low-pass filter has the following specifications: Pass-band


frequency = 0.2π rad, Pass-band ripple = 1dB

Stop-band frequency = 0.3π rad, Stop-band attenuation = 20db.

Write a MATLAB program to design this filter using Butterworth,


Chebyshev and elliptic methods. Plot the frequency response in each
case and compare the performance (Use functions like butter,
cheby1…etc).
6. Design using MATLAB, an analog Butterworth low-pass filter that has
maximum of 1dB attenuation at 30rad/sec and at least 30dB attenuation at
40rad/sec. Digitize this using bilinear and impulse invariant transformation.
Use suitable sampling rate. Plot the frequency response in each case.

Note: Refer MATLAB functions to design high pass and other filters.
Lab session 6: FIR Filter Design

Objective: To design and study finite-duration impulse response (FIR) filters.

Theory: A typical FIR (low-pass) filter response is as shown below

Pass band ripple in dB is = −20 log10[(1 − )/(1 + )]

Stop band attenuation in dB is = −20 log10

Pass band deviation in dB = −20 log10(1 + )

For an Nth order filter, there are N+1 coefficients or filter length is N+1.

Important design methods are window based design and frequency sampling design.

Window based Design: The commonly used windows are hamming window and Kaiser
Window.

Hamming window,

The normalized transition width is δw, (transition width divided by sampling frequency or
(ws-wp) divided by 2π) and filter length are related by δw=3.3/ , and the designed filter
will have a stop band attenuation of 53 dB.

Kaiser window,

Where I0[.] is modified zero order Bessel function and β is a parameter that depends on
N.

The filter length M is given by


where

The other windows are rectangular, hanning, Bartlett, Blackman and raised cosine.

Frequency Sampling Design: The desired frequency response is sampled at discrete


frequency points 2πk/N to get frequency domain points from which filter taps are
obtained.

For linear phase filter, [ ] = [ ] ∠ [ ] where,

[ ] = sampled magnitude value

2 −1
− ( ), = 0,1, … 2
∠ [ ]={ 2 ( − ) −1 and
{ }, = +1
2 ,…, +1

where, , and if N is even, lower integer of is taken for k.

To get filter coefficients, N-point inverse FFT is taken.

Figure 6.1 Frequency Sampling Design example


Note: The length M can be estimated as follows
For LPF/HPF:

, where D and F are functions of and given by

= [ 1(log10 )2 + 2 log10 + 3] log10 + [ 4(log10 )2 + 5 log10 + 6]

= 11.01217 + 0.51244 (log10 − log10 ) 1

= 5.309 × 10−3, 2= 7.114 × 10−2, 3= −4.761 × 10−1,

4= −2.66 × 10−3, 5= −0.5941, 6= −0.4278.

and is the normalized transition width.

For BPF:

, where C and G are functions of and given by

= [ 1(log10 )2 + 2 log10 + 3] log10 + [ 4(log10 )2

+ 5 log10 + 6]

= −14.61 log10( / ) − 16.9

1= 0.01201, 2= 0.09664, 3= −0.51325, 4= 0.00203, 5= −0.5705, 6=


−0.44314.

Example-1: Design, using window based approach, a low-pass digital FIR filter whose
desired specifications are: wp = 0.2π rad, ws = 0.3π rad, Rp = 0.2dB, As = 50dB.

a) Hamming window

Solution:
clc, clear all, close all;
wp = 0.2*pi; ws = 0.3*pi;
N = ceil(6.6*pi/(ws-wp));

% estimates order of the filter


b = fir1(N, 0.2); [H,W] = freqz(b,1,512);
plot(W/pi,20*log10(abs(H))), grid on, title('FIR filter using Hamming window');
xlabel('Normalized Frequency'), ylabel('Gain in dB');

% this filter provides about 53dB attenuation in stop band


% and pass band ripples are less than 0.25dB as expected.

Solution:

Normalized Frequency

Figure 6.2 FIR filter frequency response


Example-2: Design the filter specified in problem no.1 using frequency sampling design
technique assuming linear phase.
Solution:

clc, clear all, close all;


M = 20;
% assumed
alpha = (M-1)/2;
Hrs = [1,1,1,zeros(1,15),1,1];
k1 = 0: floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)*k1/M, alpha*(2*pi)*(M-k2)/M];
H = Hrs.*exp(j*angH); h = real(ifft(H,M));
[H,W]=freqz([h],1,1024);
plot(W/pi,20*log10(abs(H))), zoom on, grid on;
title('FIR LPF using frequency sampling design');
xlabel('Normalized Frequency'), ylabel('Gain in dB');
% Note that stop band attenuation is not very good, is about 15dB only
Solution:

Normalized Frequency
Figure 6.3 FIR filter frequency response
Exercise:

1. Design, using Kaiser window based approach, a low-pass digital FIR filter whose
desired specifications are: wp = 0.2π rad, ws = 0.3π rad, Rp = 0.2dB, As = 50dB.

2. Design and plot the frequency response of a linear phase band-pass filter whose
specifications are as follows: Lower pass band edge = 400Hz, upper pass-band
edge = 500Hz, lower stop- band edge = 300Hz and upper stop-band edge = 600Hz,
Rp = 0.5dB and As = 40dB. The sampling frequency is 2 kHz.
a. using hamming window
b. using Kaiser window

3. Design the same filter in sample problem 2 with M = 25 and compare them.

4. The attenuation in stop band was observed to be not adequate in the above
problem. To improve, a few non-zero points are considered in the transition band.
These values are normally estimated by optimization techniques. In the above
problem, use Hrs(k) = [ 1, 1, 1, 0.5, 0.1 , zeros(1,11), 0.1,0.5,1,1] Compare the
results and find out whether there is any improvement. If there is, then at what
cost?
5. Design the filter specified in problem 2 using frequency sampling designing
technique.
6. Certain signal has 10Hz, 30Hz and 50Hz components. Simulate this signal. Design
suitable FIR filter using the following methods to select only 30Hz component.
Filter the signal using the filters, and plot the spectrum of the input and the filtered
signal. Hence compare the performance of the filters. a) Window method
b) Frequency sampling designing
Reference
1. Sanjit K. Mitra ,“Digital Signal Processing Laboratory Using Matlab”,
Mcgraw-Hill College, 2005.
2. D.G. Manolakis and Vinay K. Ingle, “Applied Digital Signal Processing”,
Cambridge University Press, 2012.
3. MATLAB help

You might also like