21EEL66 – DIGITAL SIGNAL PROCESSING LABORATORY

LABORATORY MANUAL
For the academic year 2024-25

Name:………………………

USN:………………………..

Department of Electrical & Electronics Engineering

AMC ENGINEERING COLLEGE, BANGALORE – 83
 COURSE OBJECTIVES

1. To help students in developing software skills

2. To explain the the use of MATLAB/Sci Lab/Python software in conducting the experiments
of signal processing laboratory ,evaluating the DFT and IDFT of given sequence
3. To explain generation of different types of signals both in continuous and discrete time
domains
4. To explain verification of linear and circular convolution of given sequences
5. To explain evaluating the DFT and IDFT of given sequence
6. To design and implementation of IIR and FIR Filters for given frequency specifications and
realize them.

COURSE OUTCOMES

1. Analyze the sampling theorem

2. Provide solution for a given difference equation
3. Evaluation impulse response of a system
4. Analyze the performance of convolution of a given sequences
5. Analyze the computation of DFT and IDFT of a given sequence
6. Design the IIR and FIR filters for the given specifications

(i)
 DIGITAL SIGNAL PROCESSING LABORATORY

Subject Code : 21EEL66 IA Marks : 50

No. of Practical Hrs. / Week : 03 Exam Hours : 03
Total No. of Practical Hrs : 42 Exam Marks : 100

SN Experiment Page
No
1 Generation of Different signals in both Continuous and discrete time domains. 7
2 Verification of Sampling Theorem both in time and Frequency Domains 14
3 To perform Basic operations on given sequences-Signal folding, evaluation of even 18
and odd
4 Evaluation of Impulse Response of a System 25
5 Evaluation of Linear and Circular Convolution of given sequences 27
6 Computation of N-Point DFT by using Definite equations. 33
7 Evaluation of Circular Convolution of two sequences using DFT and IDFT approach 35
8 Solution of a difference Equation 37
9 Computation of N-Point DFT and IDFT by FFT Method 41
10 Design and implementation of IIR filter to meet the given specification 44
11 Design and implementation of FIR filter to meet the given specification 59

(ii)
 Digital Signal Processing Laboratory -21EEL66

INTRODUCTION TO MATLAB

MATLAB, which stands for MATrix LABoratory, is a state-of-the-art mathematical software

package for high performance numerical computation and visualization provides an interactive
environment with hundreds of built in functions for technical computation, graphics and
animation and is used extensively in both academia and industry. It is an interactive program for
numerical computation and data visualization, which along with its programming capabilities
provides a very useful tool for almost all areas of science and engineering.

At its core ,MATLAB is essentially a set (a “toolbox”) of routines (called “mfiles” or “mex
files”) that sit on your computer and a window that allows you to create new variables with
names (e.g. voltage and time) and process those variables with any of those routines (e.g. plot
voltage against time, find the largest voltage, etc). It also allows you to put a list of your
processing requests together in a file and save that combined list with a name so that you can run
all of those commands in the same order at some later time. Furthermore, it allows you to run
such lists of commands such that you pass in data.

MATLAB Windows:
MATLAB works with through these basic windows

COMMAND WINDOW
This is the main window ,it is characterized by MATLAB command prompt >> when you launch
the application program MATLAB puts you in this window all commands including those for
user-written programs ,are typed in this window at the MATLAB prompt.

THE CURRENT DIRECTORY WINDOW

The Current Directory window displays a current directory with a listing of its contents. There is
navigation capability for resetting the current directory to any directory among those set in the
path. This window is useful for finding the location of particular files and scripts so that they can
be edited, moved, renamed, deleted, etc. The default current directory is the Work subdirectory
of the original MATLAB installation directory.

Dept of EEE, AMCEC Page 1

 Digital Signal Processing Laboratory -21EEL66

THE COMMAND HISTORY WINDOW

The Command History window, at the lower left in the default desktop, contains a log of
commands that have been executed within the Command window. This is a convenient feature
for tracking when developing or debugging programs or to confirm that commands were
executed in a particular sequence during a multistep calculation from the command line.

GRAPHICS WINDOW
The output of all graphics commands typed in the command window are flushed to the graphics
or figure window, a separate gray window with white background color the user can create as
many windows as the system memory will allow.

EDIT WINDOW
This is where you write edit, create and save your own programs in files called M files.

Input-output
MATLAB supports interactive computation taking the input from the screen and flushing, the
output to the screen. In addition it can read input files and write output files.

Data Type
The fundamental data distinct data objects- integers, real numbers, matrices, character strings,
structures and cells. There is no need to declare variables as real or complex, MATLAB
automatically sets the variable to be real.

Dimensioning
Dimensioning is automatic in MAT required for vectors or arrays .we can find the dimensions of
an existing matrix or a vector with the size and length commands.

Where to work in MATLAB?

All programs and commands can be entered either in the
a) Command window b) As an M file using MATLAB editor

Note: Save all M files in a folder in the current directory. Otherwise you have to locate the file
during compiling.

Dept of EEE, AMCEC Page 2

 Digital Signal Processing Laboratory -21EEL66

Typing quit in the command prompt>> quit, will close MATLAB Development Environment.
For any clarification regarding plot etc, which are built in functions type help topic i.e. help plot

Basic Instructions in MATLAB

1. T = 0: 1:10 This instruction indicates a vector T which as initial value 0 and final value
10 with an increment of 1 in MATLAB is the array.
T = [0 1 2 3 4 5 6 7 8 9 10]

2. F= 20: 1: 100
F = [20 21 22 23 24 ……… 100]

3. T= 0:1/ pi: 1
T= [0, 0.3183, 0.6366, 0.9549]

4. zeros (1, 3)

The above instruction creates a vector of one row and three columns whose values are
zero Output= [0 0 0]

5. Transpose a vector
Suppose T= [1 2 3],
Then transpose T‟= 1
2
3
6. Matrix Operation

a) If a = [1 2 3] b = [4 5 6]
a.*b = [4 10 18]
b) If v = [0:2:8]
v = [0 2 4 6 8]
v (1:3)
ans [0 2 4]
v (1:2:4)
ans [ 0 4]
c) A = [1 2 3; 3 4 5; 6 7 8]
A=
1 23
3 45
6 78

Dept of EEE, AMCEC Page 3

 Digital Signal Processing Laboratory -21EEL66

A(2,3)
ans 5
A(1:2,2:3)
ans =
23
45
A(:,2)
ans =
2
4
7
A(3,:)
ans =
678

Operations on vector and matrices in MATLAB

MATLAB utilizes the following arithmetic operators;

+ Addition
- Subtraction
* Multiplication
/ Division
^ Power Operator
„ transpose

Relational operators in MATLAB

1) Syntax of the for loop is shown below
for k = array
commands
end
The commands between for and end statements are executed for all values stored in the array.
2) Syntax for the if loop
if expression
commands

Dept of EEE, AMCEC Page 4

 Digital Signal Processing Laboratory -21EEL66

end
This construction is used if there is one alternative only.
Two alternatives requires the following construction
if expression
commands (evaluated if expression is true)
else
commands (evaluated if expression is false)
end

3) Syntax of the switch-case construction is

switch expression (scalar or string)
case value1 (executes if expression evaluates to value1)
commands
case value2 (executes if expression evaluates to value2)
commands
...
otherwise
statements
end
Switch compares the input expression to each case value. Once the match is found it executes the
associated commands.

Basic Functions in MATLAB

1) Plot Syntax: plot (x,y)

Plots vector y versus vector x. If x or y is a matrix, then the vector is plotted versus the rows or
columns of the matrix.
2) Stem Syntax: stem(Y)
Discrete sequence or "stem" plot. Stem (Y) plots the data sequence Y as stems from the x axis
terminated with circles for the data value. If Y is a matrix then each column is plotted as a
separate series.

Dept of EEE, AMCEC Page 5

 Digital Signal Processing Laboratory -21EEL66

3) Subplot Syntax: Subplot (2,2,1)

This function divides the figure window into rows and columns. Subplot (2 2 1) divides
the figure window into 2 rows and 2 columns 1 represent number of the figure.

4) Disp Syntax: disp(X)

Description: disp(X) displays an array, without printing the array name. If X contains a
text string, the string is displayed. Another way to display an array on the screen is to type its
name, but this prints a leading "X=," which is not always desirable. Note that disp does not
display empty arrays.

5) xlabel Syntax: xlabel('string')

Description: xlabel('string') labels the x-axis of the current axes.

6) ylabel Syntax : ylabel('string')

Description: ylabel('string') labels the y-axis of the current axes.
7) Title Syntax : title('string')
Description: title('string') outputs the string at the top and in the center of the current
axes.
8) grid on Syntax : grid on
Description: grid on adds major grid lines to the current axes.

Dept of EEE, AMCEC Page 6

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO1
Generation Of Different Signals In Both Continuous And Discrete Time Domains

Aim: To generate basic signals and display it

THEORY:-
One of the more useful functions in the study of linear systems is the "unit impulse function."
An ideal impulse function is a function that is zero everywhere but at the origin, where it is infinitely high.
However, the area of the impulse is finite. This is, at first hard to visualize but we can do so by using the
graphs shown below.

Key Concept: Sifting Property of the Impulse

If b>a, then

Example: Another integral problem

Assume a<b, and evaluate the integral

Dept of EEE, AMCEC Page 7

 Digital Signal Processing Laboratory -21EEL66

Solution:
We now that the impulse is zero except at t=0 so

Unit Step Function

The unit step function and the impulse function are considered to be fundamental functions in engineering,
and it is strongly recommended that the reader becomes very familiar with both of these functions.
The unit step function, also known as the Heaviside function, is defined as such:

Sometimes, u(0) is given other values, usually either 0 or 1. For many applications, it is irrelevant what the
value at zero is. u(0) is generally written as undefined.
Derivative
The unit step function is level in all places except for a discontinuity at t = 0. For this reason, the derivative
of the unit step function is 0 at all points t, except where t = 0. Where t = 0, the derivative of the unit step
function is infinite. The derivative of a unit step function is called an impulse function. The impulse function
will be described in more detail next.
Integral
The integral of a unit step function is computed as such:

Dept of EEE, AMCEC Page 8

 Digital Signal Processing Laboratory -21EEL66

Dept of EEE, AMCEC Page 9

 Digital Signal Processing Laboratory -21EEL66

MATLAB Program:

%Unit step function

subplot(3,4,1);
t=-10:10;
y=[zeros(1,10),ones(1,11)];
stem(t,y);
title('unit step function');
xlabel('t(sec)');
ylabel('u(t)');
%Impulse function
subplot(3,4,2);
t=-10:10;
y=[zeros(1,10),ones(1,1),zeros(1,10)];
stem(t,y);
title('unit impulse function');
xlabel('t(sec)');

Dept of EEE, AMCEC Page 10

 Digital Signal Processing Laboratory -21EEL66

ylabel('i(t)');
%Ramp function
subplot(3,4,3);
t=0:10;

y=[t];
stem(t,y)
title('ramp function');
xlabel('t(sec)');
ylabel('r(t)');
%Sine function
subplot(3,4,4);
t=linspace(-10,10);
y=sin(t);
plot(t,y);
title('sine function');
xlabel('t(sec)');
ylabel('s(t)');
%Cosine function
subplot(3,4,5);
t=linspace(-10,10);
y=1+cos(t);
plot(t,y);
title('cosine function');
xlabel('t(sec)');
ylabel('c(t)');
%Falling exponential
subplot(3,4,6);
t=linspace(0,10);
y=exp(-t);
plot(t,y);
title('falling exponential function');
xlabel('t(sec)');

ylabel('e(t)');

Dept of EEE, AMCEC Page 11

 Digital Signal Processing Laboratory -21EEL66

%Rising exponential
subplot(3,4,7);
t=linspace(-10,10);
y=exp(0.1*t);
plot(t,y);
title('rising exponential function');
xlabel('t(sec)');
ylabel('e(t)');
%Rectangular function
subplot(3,4,8);
t=-10:10;
y=[zeros(1,9),ones(1,3),zeros(1,9)];
stem(t,y);
title('rectagular function');
xlabel('t(sec)');
ylabel('r(t)');
%Sinc function
subplot(3,4,9);
t=linspace(-10,10);
y=sinc(2*3.14*0.1*t);
plot(t,y);
title('sinc function');
xlabel('t(sec)');
ylabel('sinc(t)');
%Triangular wave
subplot(3,4,10);
t=-10:10;
output=((1-2*abs(t))/20);

plot(t,output);
title('triangular wave');
xlabel('t(sec)');
ylabel('x(t)');
%Complex exponential
subplot(3,4,11);

Dept of EEE, AMCEC Page 12

 Digital Signal Processing Laboratory -21EEL66

t=-10:10;
c=complex(3,2);
d=exp(c*t);
plot(t,abs(d));
disp(angle(d));
title('complex exponential wave');
xlabel('t(sec)');

ylabel('x(t)');

Result

Dept of EEE, AMCEC Page 13

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO 2

Verification Of Sampling Theorem

Aim: To verify Sampling theorem for a signal of given frequency.

Theory:
 Sampling is a process of converting a continuous time signal (analog signal) x(t) into a
discrete time signal x[n], which is represented as a sequence of numbers. (A/D converter)
 Converting back x[n] into analog (resulting in x’(t)) is the process of reconstruction. (D/A
converter)
 For x’(t) to be exactly the same as x(t), sampling theorem in the generation of x(n) from
x(t) is used. The sampling frequency fs determine the spacing between samples.
 Aliasing-A high frequency signal is converted to a lower frequency, results due to under
sampling. Though it is undesirable in ADCs, it finds practical applications in stroboscope
and sampling oscilloscopes.
 It states that exact reconstruction of a continuous time base band signal from its samples
is possible, if the signal is band limited and the sampling frequency is greater than twice
the signal bandwidth
 Nyquist ratefs=2*fd
 Under sampling fs<2*fd
 Over sampling fs>2*fd

MATLAB Implementation:

Step 1: MATLAB can generate only discrete time signals. For an approximate analog signal xt,
choose the spacing between the samples to be very small (≈0), say 50µs = 0.00005. Next choose
the time duration, say xt exists for 0.05seconds.(tfinal in program) (for low frequency say <1000
Hz choose 0.05 secs, for higher choose 0.01 secs or lesser as appropriate).
Now begin with the vector that represents the time base-
t = 0:0.00005:0.05;

Dept of EEE, AMCEC Page 14

 Digital Signal Processing Laboratory -21EEL66

The colon (:) operator in MATLAB creates a vector, in the above case a time vector running
from 0 to 0.05 in steps of 0.00005. The semicolon (;) tells MATLAB not to display the result.
Given t, the analog signal xt of frequency fd is generated as (cos(ωt)=cos(2πft))-
Xt = cos(2*pi*fd*t); pi is recognized as 3.14 by MATLAB.
Step 2: To illustrate oversampling condition, choose sampling frequency fs0=2.2*fd. For this
sampling rate T0=1/fs0, generate the time vector as n1 = 0:T0:0.05; & over sampled discrete
time signal x1 = cos(2*pi*fd*n1);
[Alternately let n1 be a fixed number of samples, say n = 0:10; & x1 = cos(2*pi*n*fd/fs0);]
The discrete signal is converted back to analog signal (reconstructed) using the MATLAB plot
function (FOH). In one plot both the analog signal (approximate one in blue color) and the
reconstructed signal (in red) are plotted for comparison
Step 3: Repeat step 2 for different sampling frequencies, i.e., fs=1.3*fd & fs=2*fd for under
sampling and Nyquist sampling conditions.

MATLAB PROGRAM:

tfinal=0.05;
t=0:0.00005:tfinal;
fd=input('Enter analog freuency');
xt=cos(2*pi*fd*t); %define analog signal for comparison
%Condition for under sampling
fs1=1.3*fd; %fs1<2*fd
n1=0:1/fs1:tfinal; %define the time vector
xn1=cos(2*pi*n1*fd); %Generate the under sampled signal
subplot (3,1,1); %plot the analog & sampled signals
plot(t,xt,'b',n1,xn1,'r*-');
title('under sampling plot');
xlabel('time');
ylabel('amplitude')

Dept of EEE, AMCEC Page 15

 Digital Signal Processing Laboratory -21EEL66

%Condition for Nyquist plot

fs2=2*fd;
n2=0:1/fs2:tfinal;
xn2=cos(2*pi*fd*n2);
subplot(3,1,2);
plot(t,xt,'b',n2,xn2,'r*-');
title('Nyquist plot');
xlabel('time');
ylabel('amplitude')
%Condition for oversampling
fs3=5*fd;
n3=0:1/fs3:tfinal;
xn3=cos(2*pi*fd*n3);
subplot(3,1,3);
plot(t,xt,'b',n3,xn3,'r*-');
title('Oversampling plot');
xlabel('time');
ylabel('amplitude');
legend('analog’,’ discrete');

Result:
Enter analog frequency 200
The plots are shown in Fig.1.1

Output:

Dept of EEE, AMCEC Page 16

 Digital Signal Processing Laboratory -21EEL66

Fig 1.1 : Plots of a sampled cosine wave of 200Hz

Dept of EEE, AMCEC Page 17

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO 3
Basic Operations On Given Signals & Sequences

AIM: To perform the basic operations on given signals

(a)Signal Folding
(b)Evaluation of Even and Odd
THEORY: Basic Operations on Signals
Time shifting: y(t)=x(t-T)The effect that a time shift has on the appearance of a signal If T is
a positive number, the time shifted signal, x (t -T ) gets shifted to the right, otherwise it gets shifted
left.

Dept of EEE, AMCEC Page 18

 Digital Signal Processing Laboratory -21EEL66

PROCEDURE:
 Open MATLAB
 Open new M-File
 Type the Program
 Save in Current Directory
 Compile and Run the program
 For the Output see Command Window/Figure Window
MATLAB PROGRAM:
(a)Signal Folding
clc;
clear all;
close all;
t=0:0.001:2;
Dept of EEE, AMCEC Page 19
 Digital Signal Processing Laboratory -21EEL66

s=sin(2*pi*5*t);
m=length(s);
n=[-m:m];
y=[0,zeros(1,m),s];
subplot(2,1,1);
plot(n,y,'g');
xlabel('time');
ylabel('amplitude');
title('original signal');
y1=[fliplr(s),0,zeros(1,m)];
subplot(2,1,2);
plot(n,y1,'r');
xlabel('time');
ylabel('amplitude');
title('folded signal');

Dept of EEE, AMCEC Page 20

 Digital Signal Processing Laboratory -21EEL66

(b)Evaluation of Even &Odd

clc;
clear all;
close all;
h=input('enter no.of samples');
m=(h-1)/2;
n=-m:m;
x=input('enter sample values');
subplot(4,1,1);
stem(n,x,'g');
xlabel('time');
ylabel('amplitude');
title('original sequence');
xmir=fliplr(x);
subplot(4,1,2);
stem(n,xmir,'r');
xlabel('time');
ylabel('amplitude');
title('folded sequence');
xeven=(x+xmir)/2;
subplot(4,1,3);
stem(n,xeven,'r');
xlabel('time');
ylabel('amplitude');
title('even part of sequence');
xodd=(x-xmir)/2;
subplot(4,1,4);
stem(n,xodd,'g');
xlabel('time');
ylabel('amplitude');
title('odd part of sequence');

OUTPUT:
Enter the number of samples 5
Enter the sample Values[1,2,3,4,5]

Dept of EEE, AMCEC Page 21

 Digital Signal Processing Laboratory -21EEL66

RESULT:
In this Experiment various operations on signals and evaluation of Odd and Even sequences
performed and output obtained.

Dept of EEE, AMCEC Page 22

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO 4

Impulse Response Of A Given System

Aim: To find the impulse response h(n) of the given LTI system whose response y(n) to an
input x(n) is given.

Theory:

In signal processing, the impulse response or impulse response function (IRF), of a dynamic
system is its output when presented with a brief input signal, called an IMPULSE. More
generally, an impulse response is the reaction of any dynamic system in response to some
external change. In both cases, the impulse response describes the reaction of the system as a
function of time.

Fig.4.1 A LTI system

 A discrete time LTI system (also called digital filters) as shown in Fig.2.1 is represented
by
o A linear constant coefficient difference equation, for example,
y[n]  a1 y[n 1]  a2 y[n  2]  b0 x[n]  b1 x[n 1]  b2 x[n  2];
o A system function H(z) (obtained by applying Z transform to the difference
Y (z) b  b z 1  b z 2
equation). H (z)   0 1 2
X (z) 1  a1 z 1  2a z 2

o A frequency response function obtained by applying DTFT on the impulse

response h[n] (or by replacing z with ejΩ in H(z)) to get
Y (e j ) b  b e j  b e2 j
H (e )  X (e j )  0
j  1 2

1  a1 e j  a 2e2 j

Dept of EEE, AMCEC Page 23

 Digital Signal Processing Laboratory -21EEL66

 Given the difference equation or H(z), the impulse response of the LTI system is found
using filter or impz MATLAB functions.

 Given only the input sequence x[n] and the output sequence y[n], we can find the impulse
function h[n] by using the inverse operation deconv. (The conv operation convolves 2
sequences x[n] and h[n] to obtain the output y[n]. If both y[n] and x[n] are finite then the
impulse response is finite (FIR). The deconvolution operation is valid only if the LTI
system is „invertible‟

MATLAB Implementation:

h = deconv(y,x): y the output sequence should be of longer length than x. DECONV - Deconvolution
and polynomial division. [Q, R] = DECONV(B, A) deconvolves vector A out of vector B. The result
is returned in vector Q and the remainder in vector R such that B = conv(A,Q) + R. If A and B are
vectors of polynomial coefficients, deconvolution is equivalent to polynomial division. The result of
dividing B by A is quotient Q and remainder R. (Remember conv is polynomial multiplication).
Note: Valid y and x should be given, else the conv verification will result in a slightly different y
(because of the remainder generated by deconv).
h=impz(b,a,N), returns the impulse response h[n], where b and a are coefficient arrays obtained from
difference equation, N = length of the response required. (refer expt.7). For the given difference
equation b=[1]; a=[1 -1 0.9].
IMPZ Impulse response of digital filter. [H,T] = IMPZ(B, A) computes the impulse response of the
filter B/A choosing the number of samples for you, and returns the response in column vector H and a
vector of times (or sample intervals) in T (T = [0 1 2 ...]').
[H, T] = IMPZ(B, A, N) computes N samples of the impulse response.

Dept of EEE, AMCEC Page 24

 Digital Signal Processing Laboratory -21EEL66

MATLAB PROGRAM:

N=input('Length of impulse response required=')

b=input('coefficient of x[n]=') %Enter x(n) Coefficients
a=input('coefficient of y[n]=') %Enter y(n) Coefficients
h=impz(b,a,N) %Impulse Response
%Plot the Impulse Response Sequence
n=0:N-1 %Time Vector for Plotting
subplot(3,1,1)
stem(b)
xlabel('time')
ylabel('amplitude')
title('input1')
subplot(3,1,2)
stem(a)
xlabel('time')
ylabel('amplitude')
title('input2')
subplot(3,1,3)
stem(n,h)
xlabel('n')
ylabel('h(n)')
title('impulse response')

RESULT:
Length of impulse response required=5
Coefficient of x[n] = [2 -1]
Coefficient of y[n] = [1 3 2]

The impulse response is h = [2 -7 17 -37 77]

Dept of EEE, AMCEC Page 25

 Digital Signal Processing Laboratory -21EEL66

Fig 4.2 : Impulse Response

Dept of EEE, AMCEC Page 26

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO 5

Evaluation of Linear Convolution and Circular Convolution

Aim: To obtain convolution of two finite duration sequences.

Theory:
 The output y[n] of a LTI (linear time invariant) system can be obtained by convolving the
input x[n] with the system‟s impulse response h[n].
 
 The convolution sum is y[n]  x[n]  h[n]  x[k]h[n  k]  x[n  k]h[k]
k  k 

 x[n] and h[n] can be both finite or infinite duration sequences.

 h[n]  0.9  u[n] ), we can
n
Even if one (or both) of the sequences is infinite (say,
analytically evaluate the convolution formula to get a functional form (closed form
solution) of y[n].
 If both the sequences are of finite duration, then we can use the MATLAB function
„conv‟ to evaluate the convolution sum to obtain the output y[n]. Convolution is
implemented as polynomial multiplication (refer MATLAB help).
 The length of y[n] = xlength + hlength -1.
 The conv function assumes that the two sequences begin at n=0 and is invoked by
y=conv(x,h).
 Even if one of the sequences begin at other values of n, say n=-3,or n=2; then we need to
provide a beginning and end point to y[n] which are ybegin = xbegin + hbegin and
yend = xend + hend respectively.

Dept of EEE, AMCEC Page 27

 Digital Signal Processing Laboratory -21EEL66

5(a) LINEAR CONVOLUTION MATLAB PROGRAM:

clc;

x1=input('enter the first sequence');

subplot(3,1,1);

stem(x1);

xlabel('time');

ylabel('amplitude');

title('plot of first sequence');

x2=input('enter the second sequence');

subplot(3,1,2);

stem(x2);

xlabel('time')

ylabel('amplitude');

title('plot of second sequence');

f=conv(x1,x2);

disp('output of linear convolution

is'); disp(f);

xlabel('time index

n');

ylabel('amplitudef');

subplot(3,1,3);

stem(f);

title('linear conv of sequence');

Dept of EEE, AMCEC Page 28

 Digital Signal Processing Laboratory -21EEL66

RESULT:

Enter the first sequence [1 2 3 4]

Enter the second sequence [5 6 7 8]

The result of convolved sequence is [5 16 34 60 61 52 32]

The plots are shown in Fig.3.1

Fig 4.1A : Linear Convolution

Dept of EEE, AMCEC Page 29

 Digital Signal Processing Laboratory -21EEL66

5(b) CIRCULAR CONVOLUTION PROGRAM

a. Circular convolution using the convolution summation formula
x1(n)=[2,1,2,1] and x2(n)=[1, 2 ,3,4]
clc;
clear
all;
a=input('Enter the sequence
x1(n)='); b=input('Enter the
sequence x2(n)='); n1=length(a);
n2=length(b);
N=max(n1,n2)
;
x=[a zeros(1,(N-
n1))]; for i=1:N
k=i;
for j=1:n2
H(i,j)=x(k)*b(j
); k=k-1;
if(k==0)
k=N;
end
end
end
y=zeros(1,N);
M=H';
for j=1:N
for i=1:n2
y(j)=M(i,j)+y(j);

Dept of EEE, AMCEC Page 30

 Digital Signal Processing Laboratory -21EEL66

end
end
disp('The output sequence is y(n)=');

disp(y);
subplot(3,1,1)
;
stem(a,'filled');
title('x1(n)');
xlabel('n--->');
ylabel('Amplitude');
subplot(3,1,2);
stem(b,'filled');
title('x2(n)');
xlabel('n--->');
ylabel('Amplitude');
subplot(3,1,3);
stem(y,'filled');
title('y(n)');
xlabel('n--->');
ylabel('Amplitude');

Dept of EEE, AMCEC Page 31

 Digital Signal Processing Laboratory -21EEL66

Output:
Enter the sequence x1(n)=[2 1 2 1]
Enter the sequence x2(n)=[1 2 3 4]
The output sequence is y(n)=
14 16 14 16

Fig 5B) Circular convolutions using the convolution summation formula

Dept of EEE, AMCEC Page 32

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO: 6

N – POINT DFT

Aim: To compute N – point DFT of a given sequence and to plot the magnitude and phase
spectrum.

Algorithm:

1) Value of N in N point DFT analysis.

2) Enter the input sequence

3) Calculate the DFT using the defining equation of DFT.

4) Calculate and plot the magnitude and phase of DFT

MATLAB PROGRAM:

N=input('enter the value of N=');

x=input('enter the sequence for which DFT is to be
calculate='); n=[0:1:N-1];
k=[0:1:N-1];
WN=exp(-
1j*2*pi/N); nk=n'*k;
WNnk=WN.^nk;
Xk=x*WNnk;
magX=abs(Xk)
phaseX=angle(Xk)*180/
pi figure(1);
subplot(2,1,1);
plot(k,magX);
subplot(2,1,2);
plot(k,phaseX)
;

Dept of EEE, AMCEC Page 33

 Digital Signal Processing Laboratory -21EEL66

Result:

enter the value of N=4

enter the sequence for which DFT is to be calculate=[1 2 3 4]
magX =
10.0000 2.8284 2.0000 2.8284
phaseX =
0 135.0000 -180.0000 -135.0000

Fig 6.1: Magnitude and phase spectrum for N=4

Dept of EEE, AMCEC Page 34

 Digital Signal Processing Laboratory -21EEL66

Experiment No: 7

Evaluation of Circular Convolution by DFT and IDFT method

Aim: To compute circular convolution by DFT and IDFT method.

Algorithm:

1) Enter the sequence x(n) and h(n).

2) Compute the DFT of both the sequences using the command fft.

3) Multiply both the DFT‟s.

4) Compute the IDFT of the previous result using the command ifft

5) Plot x(n), h(n),y(n).

MATLAB PROGRAM

%Circular convolution using DFT and IDFT method

x1= input('First sequence = ');
x2= input('Second sequence =
');
N=max(length(x1),length(x2));
%obtain
DFTs X1=
fft(x1,N);
X2= fft(x2,N);
%Perform vector
multiplication y= X1.*X2;
%ifft to get
y[n] yn=
ifft(y,N);
disp('circular convolution of x1 &x2 is

Dept of EEE, AMCEC Page 35

 Digital Signal Processing Laboratory -21EEL66

yn='); disp(yn);
%plot
stem(yn);
title('Circular convolution output y(n)');

Result:

First sequence = [1 2 3 4]

Second sequence = [1 2 3 4]

circular convolution of x1 &x2 is yn=

26 28 26 20

Fig 7.2: Circular convolution using DFT and IDFT

Dept of EEE, AMCEC Page 36

 Digital Signal Processing Laboratory -21EEL66

Experiment No 8

Solution of a Given Difference Equation

Aim: To find the solution of a given difference equation.

Algorithm:

1) Enter the coefficients of x and y from the given difference equation.

2) Use command filter filters the input data x using a rational transfer function defined by the
numerator and denominator coefficients b and a.

3) The command filtic assumes that the input x is 0 in the past. (only for finding forced response)

4) Plot x(n) and y(n).

MATLAB PROGRAM:

a. Find the output of y(n) -0.5 y(n-1) = x(n) with x(n)=u(n)-u(n-4)

clc;

clear all;

b=input('enter the coefficients of x');

a=input('enter the coefficients of y');

n=[-5:50];

x=[(n>=0)]-[(n>4)];

figure(1);

subplot(2,1,1);

stem(n,x,'filled');

title('input sequence x(n)');

Dept of EEE, AMCEC Page 37

 Digital Signal Processing Laboratory -21EEL66

xlabel('n--->');

ylabel('x');

subplot(2,1,2);

y=filter(b,a,x);

stem(n,y,'filled');

title('output sequence y(n)');

xlabel('n--->');

ylabel('y');

OUTPUT:

enter the coefficients of x[1 -1]

enter the coefficients of y[1 -0.5]

Dept of EEE, AMCEC Page 38

 Digital Signal Processing Laboratory -21EEL66

b. Find the output of y(n) -0.25 y(n-1) -0.125 y(n-2)= x(n) with x(n)=u(n)-u(n10) and initial
condition y(-1)=1,y(-2)= -2

clc;

clear all;

b=input('enter the coefficients of x=');

a=input('enter the coefficients of y=');

c=input('enter the initial conditions=');

n=[-5:50];

ic=filtic(b,a,c);

x=[(n>=0)]-[(n>10)];

y=filter(b,a,x,ic);

subplot(2,1,1);

stem(n,x,'filled');

title('input sequence x(n)');

xlabel('n');

ylabel('x');

subplot(2,1,2);

stem(n,y,'filled');

title('output sequence y(n)');

xlabel('n');

ylabel('y')

Dept of EEE, AMCEC Page 39

 Digital Signal Processing Laboratory -21EEL66

OUTPUT:

enter the coefficients of x=1

enter the coefficients of y=[1 -0.25 -0.125]

enter the initial conditions=[1-2]

Dept of EEE, AMCEC Page 40

 Digital Signal Processing Laboratory -21EEL66

Experiment No: 9

Computation of N-point DFT and IDFT by FFT Method

Aim: To calculate DFT and IDFT by FFT using MATLAB Algorithm:

1) Get the input sequence
2) Find the FFT of the input sequence using MATLAB function.
3) Find the IFFT of the input sequence using MATLAB function.
4) Display the above outputs using stem function.

MATLAB PROGRAM:

clc;

clear all;

xn=input('Enter the input sequence');

N=length(xn);

n=0:1:N-1;

k=0:1:N-1;

Xk=fft(xn,N);

disp('DFT of the given sequence is');

disp(Xk)

disp('The magnitude sequence is');

MagX=abs(Xk)

disp('The phase sequence is');

PhaseX=angle(Xk)*180/pi

disp('IDFT of the sequence Xk is');

Dept of EEE, AMCEC Page 41

 Digital Signal Processing Laboratory -21EEL66

xn=ifft(Xk)

subplot(2,2,1);

stem(n,xn,'filled');

title('input sequence');

ylabel('xn');

xlabel('n');

subplot(2,2,2);

stem(n,xn,'filled');

title('IDFT Sequence');

ylabel('xn');

xlabel('n');

subplot(2,2,3);

stem(k,MagX,'filled');

title('Magnitude Spectrum');

ylabel('Mag Xk');

xlabel('K');

subplot(2,2,4);

stem(k,PhaseX,'filled');

title('Phase Spectrum');

ylabel('PhaseXk');

xlabel('K');

Dept of EEE, AMCEC Page 42

 Digital Signal Processing Laboratory -21EEL66

Output:

Enter the input sequence[4 3 2 1]

DFT of the given sequence is
10.0000 2.0000 - 2.0000i 2.0000 2.0000 + 2.0000i
The magnitude sequence is
MagX = 10.0000 2.8284 2.0000 2.8284
The phase sequence is
PhaseX = 0 -45 0 45

IDFT of the sequence Xk is

xn = 4 3 2 1

Dept of EEE, AMCEC Page 43

 Digital Signal Processing Laboratory -21EEL66

EXPERIMENT NO : 10

DESIGN OF IIR FILTERS

Aim: To design and implement IIR filters to meet the given specifications.

Algorithm:

1) Enter the pass band ripple (rp) and stop band ripple (rs).

2) Enter the pass band frequency (wp) and stop band frequency (ws).

3) Get the sampling frequency (fs).

4) Calculate normalized pass band frequency, and normalized stop band frequency w1 and w2
respectively.

w1 = 2 * wp /fs

w2 = 2 * ws /fs

5) Make use of the following function to calculate order of filter

Butterworth filter order [n,wn]=buttord(w1,w2,rp,rs)

Chebyshev filter order [n,wn]=cheb1ord(w1,w2,rp,rs)

6) Design an nth order digital lowpass Butterworth or Chebyshev filter using the

following statements.

Butterworth filter [b, a]=butter (n, wn)

Chebyshev filter [b,a]=cheby1(n,0.5,wn)

OR

Design an nth order digital high pass Butterworth or Chebyshev filter using the following

statement.

Dept of EEE, AMCEC Page 44

 Digital Signal Processing Laboratory -21EEL66

Butterworth filter [b,a]=butter (n, wn,‟high‟)

Chebyshev filter [b,a]=cheby1 (n, 0.5, wn,'high')

7) Find the digital frequency response of the filter by using „freqz( )‟ function

[H,w]=freqz(b,a,512,fs)

8) Calculate the magnitude of the frequency response in decibels (dB)

mag=20*log10 (abs (H))

9) Plot the magnitude response [magnitude in dB Vs normalized frequency (Hz]]

10) Calculate the phase response using an = angle (H)

11) Plot the phase response [phase in radians Vs normalized frequency (Hz)].

MATLAB PROGRAM

Design of Butterworth filters(Low Pass, High Pass, Bandpass and Bandstop)

Using Bilinear Transformation:

%using bilinear transformation LPF and HPF butterworth

clc; clear all; close all;
warning off;
disp('enter the IIR filter design specifications');
rp=input('enter the passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter the passband freq');
ws=input('enter the stopband freq');
fs=input('enter the sampling freq');
w1=2*wp/fs; %normalized pass band frequency
w2=2*ws/fs; %normalized stop band frequency
[n,wn]=buttord(w1,w2,rp,rs); % Find the order n and cutoff frequency

Dept of EEE, AMCEC Page 45

 Digital Signal Processing Laboratory -21EEL66

ch=input('give type of filter 1:LPF,2:HPF');

switch ch
case 1
disp('Frequency response of Butterworth IIR LPF is:');
[b,a]=butter(n,wn); % Find the filter coefficient of LPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs); % to get the transfer func of the filter
mag=20*log10(abs(H));
phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');
title('Magnitude Response of the desired Butterworh LPF');
subplot(212);
plot(w,phase);
grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Butterworh LPF');
case 2
disp('Frequency response of IIR Butterworth HPF is:');
[b,a]=butter(n,wn,'high'); % Find the filter coefficients of HPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs); % to get the transfer func of the filter
mag=20*log10(abs(H));
phase=angle(H);

Dept of EEE, AMCEC Page 46

 Digital Signal Processing Laboratory -21EEL66

subplot(211);
plot(w,mag);
grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');

Output:

enter the IIR filter design specifications

enter the passband ripple = 0.5

enter the stopband ripple = 50

enter the passband freq = 1200

enter the stopband freq = 2400

enter the sampling freq = 10000

give type of filter 1:LPF,2:HPF1

Frequency response of Butterworth IIR LPF is:

Fig 9.1: Magnitude and phase response of Butterworth IIR Low Pass Filter

Dept of EEE, AMCEC Page 47

 Digital Signal Processing Laboratory -21EEL66

Frequency response of Butterworth IIR HPF is:

Fig 9.2: Magnitude and phase response of Butterworth IIR High Pass Filter

% using bilnear tranformation for BPF and BSF butterworth

clc; clear all;
warning off;
disp('enter the IIR filter design specifications');
rp=input('enter the passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter the passband freq');
ws=input('enter the stopband freq');
fs=input('enter the sampling freq');
w1=2*wp/fs;%normalized pass band frequency
w2=2*ws/fs;%normalized stop band frequency
[n]=buttord(w1,w2,rp,rs);% Find the order n and cutoff frequency
wn=[w1 w2];

Dept of EEE, AMCEC Page 48

 Digital Signal Processing Laboratory -21EEL66

ch=input('give type of filter 1:BPF,2:BSF');

switch ch
case 1
disp('Frequency response of Butterworth IIR BPF is:');
[b,a]=butter(n,wn,'bandpass');%Find the filter coefficient ofBPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer func of the filter
mag=20*log10(abs(H));

phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');
title('Magnitude Response of the desired Butterworh BPF');
subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Butterworh BPF');
case 2
disp('Frequency response of IIR Butterworth BSF is:');
[b,a]=butter(n,wn,'stop'); % Find the filter coefficients of BSF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer func of the filter
mag=20*log10(abs(H));
phase=angle(H);
subplot(211);

Dept of EEE, AMCEC Page 49

 Digital Signal Processing Laboratory -21EEL66

plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');
title('Magnitude Response of the desired Butterworh BSF');
subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Butterworh BSF');
end
Output:
enter the IIR filter design specifications
enter the passband ripple 0.3
enter the stopband ripple 40
enter the passband freq 1500
enter the stopband freq 2000
enter the sampling freq 9000
give type of filter 1:BPF,2:BSF1
Frequency response of Butterworth IIR BPF is:

Fig 9.3: Magnitude and phase response of Buttrworth IIR Band Pass Filter

Dept of EEE, AMCEC Page 50

 Digital Signal Processing Laboratory -21EEL66

enter the IIR filter design specifications

enter the passband ripple 0.4
enter the stopband ripple 46
enter the passband freq 1100
enter the stopband freq 2200
enter the sampling freq 6000
give type of filter 1:BPF,2:BSF2
Frequency response of IIR Butterworth BSF is:

Fig 9.4: Magnitude and phase response of Butterworth IIR Band Stop Filter

Note: For design of filter by impulse invariance method use the command impinvar
instead of bilinear.

Dept of EEE, AMCEC Page 51

 Digital Signal Processing Laboratory -21EEL66

Design of Chebyshev filters(Low Pass, High Pass, Bandpass and Bandstop)

Using Bilinear Transformation:
%using bilinear transformation LPF and HPF chebyshev
clc; clear all; close all;
warning off;
disp('enter the IIR filter design specifications');
rp=input('enter the passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter the passband freq');
ws=input('enter the stopband freq');
fs=input('enter the sampling freq');
w1=2*wp/fs;%normalized pass band frequency
w2=2*ws/fs;%normalized stop band frequency
[n,wn]=cheb1ord(w1,w2,rp,rs);% Find the order n and cutoff
frequency
ch=input('give type of filter 1:LPF,2:HPF');

switch ch
case 1
disp('Frequency response of Chebyshev IIR LPF is:');
[b,a]=cheby1(n,rp,wn); % Find the filter coefficient of LPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer function of the
filter
mag=20*log10(abs(H));
phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');

Dept of EEE, AMCEC Page 52

 Digital Signal Processing Laboratory -21EEL66

title('Magnitude Response of the desired Chebyshev LPF');

subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Chebyshev LPF');
case 2
disp('Frequency response of IIR Chebyshev HPF is:');
[b,a]=cheby1(n,rp,wn,'high'); % Find the filter coefficie of HPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer func of the filter
mag=20*log10(abs(H));
phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');
title('Magnitude Response of the desired Chebyshev HPF');
subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Chebyshev HPF');
end
Output:
enter the IIR filter design specifications
enter the passband ripple 0.2
enter the stopband ripple 45

Dept of EEE, AMCEC Page 53

 Digital Signal Processing Laboratory -21EEL66

enter the passband freq 1300

enter the stopband freq 1500
enter the sampling freq 10000
give type of filter 1:LPF,2:HPF1
Frequency response of Chebyshev IIR LPF is:

Fig 9.5: Magnitude and phase response of Chebyshev IIR Low Pass Filter

enter the IIR filter design specifications

enter the passband ripple 0.3
enter the stopband ripple 60
enter the passband freq 1500
enter the stopband freq 2000
enter the sampling freq 9000
give type of filter 1:LPF,2:HPF2
Frequency response of IIR Chebyshev HPF is:

Dept of EEE, AMCEC Page 54

 Digital Signal Processing Laboratory -21EEL66

Fig 9.6: Magnitude and phase response of Chebyshev IIR High Pass Filter

%using bilnear tranformation for BPF and BSF chebyshev

clc; clear all; close all;
warning off;
disp('enter the IIR filter design specifications');
rp=input('enter the passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter the passband freq');
ws=input('enter the stopband freq');
fs=input('enter the sampling freq');
w1=2*wp/fs;%normalized pass band frequency
w2=2*ws/fs;%normalized stop band frequency
[n]=cheb1ord(w1,w2,rp,rs);% Find the order n & cutoff frequency
wn=[w1 w2];
ch=input('give type of filter 1:BPF,2:BSF');
switch ch

Dept of EEE, AMCEC Page 55

 Digital Signal Processing Laboratory -21EEL66

case 1
disp('Frequency response of Chebyshev IIR BPF is:');
[b,a]=cheby1(n,rp,wn,'bandpass'); % Find the filter coeff of BPF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer func of the filter
mag=20*log10(abs(H));

phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');
xlabel('--> Normalized frequency in Hz');
title('Magnitude Response of the desired Chebyshev BPF');
subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Chebyshev BPF');
case 2
disp('Frequency response of IIR Chebyshev BSF is:');
[b,a]=cheby1(n,rp,wn,'stop'); % Find the filter coeffici of BSF
z=tf(b,a)
[num,den]=bilinear(b,a,fs);
t=tf(num,den,1/fs)
[H,w]=freqz(b,a,512,fs);% to get the transfer func of the filter
mag=20*log10(abs(H));
phase=angle(H);
subplot(211);
plot(w,mag);grid on;
ylabel('--> Magnitude in dB');

Dept of EEE, AMCEC Page 56

 Digital Signal Processing Laboratory -21EEL66

xlabel('--> Normalized frequency in Hz');

title('Magnitude Response of the desired Chebyshev BSF');
subplot(212);
plot(w,phase);grid on;
ylabel('--> Phase in radians');
xlabel('--> Normalized frequency in Hz');
title('Phase Response of the desired Chebyshev BSF');
end
Output:
enter the IIR filter design specifications
enter the passband ripple 0.4
enter the stopband ripple 35
enter the passband freq 2000
enter the stopband freq 2500
enter the sampling freq 10000

give type of filter 1:BPF,2:BSF1

Frequency response of Chebyshev IIR BPF is:

Fig 9.7: Magnitude and phase response of Chebyshev IIR Band Pass Filter

Dept of EEE, AMCEC Page 57

 Digital Signal Processing Laboratory -21EEL66

enter the IIR filter design specifications

enter the passband ripple 0.25
enter the stopband ripple 40
enter the passband freq 2500
enter the stopband freq 2750
enter the sampling freq 7000
give type of filter 1:BPF,2:BSF2
Frequency response of IIR Chebyshev BSF is:

Fig 9.8: Magnitude and phase response of Chebyshev IIR Band Stop Filter

Dept of EEE, AMCEC Page 58

 Digital Signal Processing Laboratory -21EEL66

Experiment No 11

DESIGN OF FIR FILTERS

Aim: To design and implement FIR filters using different windows to meet the given
specifications.
MATLAB PROGRAM
%FIR Low Pass/High pass filter design using Rectangular/Bartlett Window

clc; clear all; close all;

rp=input('enter passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter passband freq');
ws=input('enter stopband freq');
fs=input('enter sampling freq ');
w1=2*wp/fs;
w2=2*ws/fs;
num=-20*log10(sqrt(rp*rs))-13;
dem=14.6*(ws-wp)/fs;
n=ceil(num/dem);
n1=n+1;
if(rem(n,2)~=0)
n1=n; n=n-1;
end
c=input('enter your choice of window function 1. rectangular 2.
bartlett');
if(c==1)
y=rectwin(n1);
disp('Rectangular window filter response');
end
if (c==2)
Dept of EEE, AMCEC Page 59
 Digital Signal Processing Laboratory -21EEL66

y=bartlett(n1);
disp('bartlett window filter response');
end
ch=input('give type of filter 1:LPF,2:HPF');
switch ch
case 1
b=fir1(n,w1,y);
[h,o]=freqz(b,1,256);
m=20*log10(abs(h));
plot(o/pi,m);
title('LPF');
xlabel('(a) Normalized frequency-->');
ylabel('Gain in dB-->');

case 2
b=fir1(n,w1,'high',y);
[h,o]=freqz(b,1,256);
m=20*log10(abs(h));
plot(o/pi,m);
title('HPF');
xlabel('(b) Normalized frequency-->');
ylabel('Gain in dB-->');
end
Output:
enter passband ripple 0.05
enter the stopband ripple 0.04
enter passband freq 1500
enter stopband freq 2000
enter sampling freq 9000
enter your choice of window function 1. rectangular 2. bartlett1
Rectangular window filter response

Dept of EEE, AMCEC Page 60

 Digital Signal Processing Laboratory -21EEL66

give type of filter 1:LPF,2:HPF1

Fig 10.1:LPF Using Rectangular Window

Fig 10.2:HPF Using Rectangular Window

Dept of EEE, AMCEC Page 61

 Digital Signal Processing Laboratory -21EEL66

Fig 10.3: LPF Using Bartlett Window

Fig 10.4: HPF Using Bartlett Window

Dept of EEE, AMCEC Page 62

 Digital Signal Processing Laboratory -21EEL66

%FIR Band Pass/BandStop filter design using Rectangular/Bartlett window

clc; clear all; close all;
rp=input('enter passband ripple');
rs=input('enter the stopband ripple');
wp=input('enter passband freq');
ws=input('enter stopband freq');
fs=input('enter sampling freq ');
w1=2*wp/fs;
w2=2*ws/fs;
num=-20*log10(sqrt(rp*rs))-13;
dem=14.6*(ws-wp)/fs;
n=ceil(num/dem);
n1=n+1;
if(rem(n,2)~=0)
n1=n; n=n-1;
end
c=input('enter your choice of window function 1. rectangular 2.
bartlett');
if(c==1)
y=rectwin(n1);
disp('Rectangular window filter response');
end
if (c==2)
y=bartlett(n1);
disp('bartlett window filter response');
end
ch=input('give type of filter 1:BPF,2:BSF');
wn=[w1 w2];
switch ch
case 1
b=fir1(n,wn,y);

Dept of EEE, AMCEC Page 63

 Digital Signal Processing Laboratory -21EEL66

[h,o]=freqz(b,1,256);
m=20*log10(abs(h));
plot(o/pi,m);
title('BPF');
xlabel('(a) Normalized frequency-->');
ylabel('Gain in dB-->');
case 2
b=fir1(n,wn,'stop',y);
[h,o]=freqz(b,1,256);
m=20*log10(abs(h));
plot(o/pi,m);
title('BSF');
xlabel('(b) Normalized frequency-->');
ylabel('Gain in dB-->');
end

Output:
enter passband ripple 0.05
enter the stopband ripple 0.04
enter passband freq 1500
enter stopband freq 2000
enter sampling freq 9000
enter your choice of window function 1. rectangular 2. bartlett1
Rectangular window filter response
give type of filter 1:BPF,2:BSF1

Dept of EEE, AMCEC Page 64

 Digital Signal Processing Laboratory -21EEL66

Fig 10.5: BPF Using Rectangular Window

Fig 10.6: BSF Using Rectangular Window

Dept of EEE, AMCEC Page 65

 Digital Signal Processing Laboratory -21EEL66

Fig 10.7: BPF Using Bartlett Window

Fig 10.8: BSF Using Bartlett Window

Note: Repeat the above program for Blackman, Hammimg and Hannning Window.

Dept of EEE, AMCEC Page 66