DEPARTMENT OF ELECTRONICS AND

COMMUNICATION ENGINEERING

gec, BHARUCH

Digital Signal Processing

Subject Code: 3171003

B.E. IV (7th Sem.) (EC)
 CERTIFICATE
This is to certify that Mr./Ms………………….

……………………Of B.E …........... Semester

......…. Enrollment No …………….……Branch

……………….. has been found satisfactory in the

continuous internal evaluation of the laboratory,

practical and term work in the subject Digital Signal

Processing for the academic year 20… – 20…

Date: Sign. of teacher

2
 INDEX
Sr. No. Title Date Sign

1 Introduction to SCILAB.

Discrete Integrator and Differentiator

2

Impulse Response & Total Response

3

Application of Correlation
4

Z-transform of finite sequences

5

Pole-zero plot on z-plane

6

Frequency response of filters

7

Linear and Circular convolution using DFT-IDFT

8

9 Time-domain and frequency domain characteristics of

different windows

10 FIR filter design

11 IIR Filter design

3
EXPERIMENT 1 Date ___________

Introduction to SCILAB
Scilab is a free and open-source, cross-platform numerical computational package and
a high-level, numerically oriented programming language. It can be used for signal
processing, statistical analysis, image enhancement, fluid dynamics simulations, numerical
optimization, and modelling, simulation of explicit and implicit dynamical systems and (if the
corresponding toolbox is installed) symbolic manipulations.
Scilab is a high-level, numerically oriented programming language. The language
provides an interpreted programming environment, with matrices as the main data type. By
using matrix-based computation, dynamic typing, and automatic memory management, many
numerical problems may be expressed in a reduced number of code lines, as compared to
similar solutions using traditional languages, such as Python, C, or C++. This allows users to
rapidly construct models for a range of mathematical problems. While the language provides
simple matrix operations such as multiplication, the Scilab package also provides a library of
high-level operations such as correlation and complex multidimensional arithmetic.
Scilab also includes a free package called Xcos for modelling and simulation of
explicit and implicit dynamical systems, including both continuous and discrete sub-systems.
Xcos is the open source equivalent to Simulink from MathWorks.
As the syntax of Scilab is similar to MATLAB, Scilab includes a source code
translator for assisting the conversion of code from MATLAB to Scilab. Scilab is available
free of cost under an open source license. Due to the open-source nature of the software,
some user contributions have been integrated into the main program.

Install SCILAB
Scilab is numerical computation software that anybody can freely download. Available
under Windows, Linux and Mac OS X, Scilab can be downloaded at the following address:
http://www.scilab.org/

Become familiar with Scilab

The useful workspace in Scilab consists of several windows:
● The console for making calculations,
● The editor for writing programs,
● The graphics windows for displaying graphics,
● The embedded help.
After double-clicking the icon to launch Scilab, the Scilab environment by default consists
of the following docked windows – console, files and variables browsers, and command
history

4
In the console after the prompt “ --> “, just type a command and press the Enter (Windows
and Linux) or Return key (Mac OS X) on the keyboard to obtain the corresponding result.
--> 57/4
ans =
14.25
--> (2+9)^5
ans =
161051.
It is possible to come back at any moment with the keyboard's arrow keys ← ↑ → ↓ or with
the mouse. The left and right keys are used to change the instructions and the up and down
keys are used to come back on a previously executed command.

Particular numbers
%e and %pi represent respectively e and π:

--> %e

%e = 2.7182818

--> %pi

%pi = 3.1415927

%i represents the i of complexes in input and is displayed i in output:

--> 2+3*%i

ans = 2. + 3.i

For not displaying the results In adding a semicolon “ ; “ at the end of a command line, the
calculation is done but the result is not displayed.

5
-->(1+sqrt(5))/2;

--> (1+sqrt(5))/2

ans = 1.618034

The menu bar

The menus listed below are particularly useful. Applications

● The command history allows you to find all the commands from previous sessions to
the current session.
● The variables browser allows you to find all variables previously used during the
current session.

Edit

Preferences (in the Scilab menu under Mac OS X) allow you to set and customize colours,
fonts and font size in the console and in the editor, which is very useful for screen projection.
Clicking on Clear Console clears the entire content of the console. In this case, the command
history is still available and calculations made during the session remain in memory.
Commands that have been erased are still available through the keyboard’s arrow keys.

Control

To interrupt a running program, you can:

● Type pause in the program or click on Control> Interrupt in the menu bar (Ctrl X
under Windows and Linux or Command X under Mac OS X), if the program is
already running. In all cases, the prompt “ --> “ will turn into “ -1-> “, then into “
-2-> “…, if the operation is repeated.
● To return to the time before the program interruption, type resume in the console or
click on Control > Resume.
● To quit for good a calculation without any possibility of return, type abort in the
console or click on Control > Abort in the menu bar.

The editor
Typing directly into the console has two disadvantages: it is not possible to save the
commands and it is not easy to edit multiple lines of instruction. The editor is the appropriate
tool to run multiple instructions.

Opening the editor

To open the editor from the console, click on the first icon in the toolbar or on
Applications > SciNotes in the menu bar. The editor opens with a default file named
“Untitled 1“.

Writing in the editor

6
Typing in the editor is like in any word processor. In the text editor, opening and closing
parentheses, end loops, function and test commands are added automatically. However, these
features can be disabled in Options > Auto-completion on the menu, by clicking on the two
below en tries enabled by default:

● (,[,…
● if, function,…

While in principle each instruction should be entered on a separate line, it is possible to type
multiple statements on the same line separating each statement with a semicolon “ ; “. A
space at the start of the line called indentation is automatic when a loop or a test is started. In
the following example, we calculate 10 terms of the sequence (𝑢𝑛) defined by:

{𝑢1 = 1 𝑢𝑛+1 = 2𝑢𝑛 − 3

Scilab keywords
● backslash — (\) left matrix division: Exact or least square solution
● brackets [,;] — Concatenation. Recipients of an assignment. Results of a function
● colon (:) — Ranging operator. Addresses all elements along an array dimension or of
a list.
● comma — (,) comma; instruction, argument separator
● comments — (// or /*...*/) comments
● comparison — comparison, relational operators
● dollar — ($) last index
● dot — (.) symbol
● equal — (=) assignment, comparison, equal sign
● getscilabkeywords — returns a list with all Scilab keywords.
● hat — (^) exponentiation
● greater — (>) greater than comparison
● minus — (-) subtraction operator. Sign change

7
 ● not — (~) logical not
● parentheses — ( ) left and right parenthesis
● per cent — (%) special character
● plus (+) — Numerical addition. Text concatenation (glueing)
● quote — (') transpose operator, string delimiter
● semicolon — (;) ending expression and row separator
● slash — (/) right divisions. System's feedback. Comments
● star — (*) multiplication operator
● symbols — Scilab operator names
● tilde — (~) logical not

Control flow
● abort — interrupt evaluation.
● break — keyword to interrupt loops
● case — keyword used in statement "select"
● continue — keyword to pass control to the next iteration of a loop
● do — language keyword for loops
● else — keyword in if-then-else and select-case-then-else
● elseif — keyword in if-then-else
● end — end keyword
● for — keyword entering a non-conditional loop
● halt — stop execution
● if — keyword for conditional execution
● pause — temporarily pauses the running execution, and allows instructions in the
console.
● resume — return or resume execution and copy some local variables
● return — return or resume execution and copy some local variables
● select — select keyword
● then — a keyword in control flows 'if' and 'select'
● catch — beginning of catch block in try-catch control instruction
● while — Opens a block of instructions iterated on a heading condition

Signal Processing
● Convolution - Correlation
o conv — discrete 1-D convolution.
o conv2 — discrete 2-D convolution.
o convol — convolution
o convol2d — discrete 2-D convolution, using FFT.
o corr — correlation, covariance
o hank — covariance to hankel matrix
o xcorr — Computes discrete auto or cross correlation
● Filters
o analpf — create analog low-pass filter

8
o buttmag — Power transmission of a Butterworth filter
o casc — cascade realization of filter from coefficients (utility function)
o cheb1mag — response of Chebyshev type 1 filter
o cheb2mag — response of type 2 Chebyshev filter
o ell1mag — magnitude of elliptic filter
o eqfir — minimax approximation of FIR filter
o eqiir — Design of iir filters
o faurre — filter computation by simple Faurre algorithm
o ffilt — coefficients of FIR low-pass
o filt_sinc — samples of sinc function
o filter — filters a data sequence using a digital filter
o find_freq — parameter compatibility for elliptic filter design
o frmag — magnitude of FIR and IIR filters
o fsfirlin — design of FIR, linear phase filters, frequency sampling technique
o group — group delay for digital filter
o hilbert — Discrete-time analytic signal computation of a real signal using
Hilbert transform
o iir — iir digital filter
o iirgroup — group delay Lp IIR filter optimization
o iirlp — Lp IIR filter optimization
o kalm — Kalman update
o lev — Yule-Walker equations (Levinson's algorithm)
o levin — Toeplitz system solver by Levinson algorithm (multidimensional)
o lindquist — Lindquist's algorithm
o remez — Remez exchange algorithm for the weighted chebyshev
approximation of a continuous function with a sum of cosines.
o remezb — Minimax approximation of magnitude response
o sgolay — Savitzky-Golay Filter Design
o sgolaydiff — Numerical derivatives computation using Savitzky-Golay filter.
o sgolayfilt — Filter signal using Savitzky-Golay Filter.
o srfaur — square-root algorithm
o srkf — square root Kalman filter
o sskf — steady-state Kalman filter
o syredi — Design of iir filters, syredi code interface
o system — observation update
o trans — low-pass to other filter transform
o wfir — linear-phase FIR filters
o wfir_gui — Graphical user interface that can be used to interactively design
wfir filters
o wiener — Wiener estimate
o wigner — 'time-frequency' wigner spectrum
o window — compute symmetric window of various type
o yulewalk — least-square filter design
o zpbutt — Butterworth analog filter

9
 o zpch1 — Chebyshev analog filter
o zpch2 — Chebyshev analog filter
o zpell — lowpass elliptic filter
● Identification
o frfit — frequency response fit
o lattn — recursive solution of normal equations
o lattp — Identification of MA part of a vector ARMA process
o mrfit — frequency response fit
o phc — Markovian representation
o rpem — Recursive Prediction-Error Minimization estimation
● Spectral estimation
o cepstrum — cepstrum calculation
o cspect — two sided cross-spectral estimate between 2 discrete time signals
using the correlation method
o czt — chirp z-transform algorithm
o mese — maximum entropy spectral estimation
o pspect — two sided cross-spectral estimate between 2 discrete time signals
using the Welch's average periodogram method.
● Transforms
o idct — Inverse discrete cosine transform.
o idst — Inverse discrete sine transform.
o ifft — Inverse fast Fourier transform.
o fft2 — two-dimension fast Fourier transform
o fftshift — rearranges the FFT output, moving the zero frequency to the centre
of the spectrum
o hilb — FIR approximation to a Hilbert transform filter
o ifftshift — an inverse of fftshift
● bilt — bilinear or biquadratic transform SISO system given by a zero/poles
representation
● detrend — remove constant, linear or piecewise linear trend from a vector
● intdec — Changes sampling rate of a signal
● unwrap — unwrap a Y(x) profile or a Z(x,y) surface. Unfold a Y(x) profile
● xcov — Computes discrete auto or cross covariance

Marks/Grade Teacher’s Signature

10
EXPERIMENT 2 Date ___________

Discrete Integrator and Differentiator

AIM : (a) Find out the output of Discrete Integrator of input x(n) =
cos(2πn/20)u(n) using Accumulator.

(b) Find out the discrete differentiator of x(n) =

cos(2πn/20)u(n) using the Forward difference system and
backward difference system.

Part-(a) Accumulator

Program:

clc; // Clear command window

close; // Close all figures
clear; // Clear all variables
N=50; //Number of Samples
n=0:N-1;
x=cos(2*%pi*n/20); // Input signal
s=0;
for i=1:N
y(i)=s+x(i);
s=y(i);
end
subplot(211),plot2d3(n,x);xlabel("n"); ylabel("x(n)"); title("Input Signal");
subplot(212),plot2d3(n,y);xlabel("n"); ylabel("y(n)"); title("Output Signal");

11
Plots:

Part-(b) Backward Difference System

Program:

clc; // Clear command window

close; // Close all figures
clear; // Clear all variables
N=50;
n=0:N-1;
x=cos(2*%pi*n/20);
x1=[x 0];
x2=[0 x];
y=x1-x2;
m=0:N;

12
subplot(211),plot2d3(n,x);xlabel('n');ylabel('x(n)');title("Input Signal");
subplot(212),plot2d3(m,y);xlabel('n');ylabel('y(n)');title("Output Sinal");

Plots:

Part-(b) Forward Difference System

Program:

clc; // Clear command window

close; // Close all figures
clear; // Clear all variables
N=50;
n=0:N-1;
x=cos(2*%pi*n/20);
x1=[x 0];
x2=[0 x];

13
y=x1-x2;
m=-1:N-1;
subplot(211),plot2d3(n,x);xlabel('n');ylabel('x(n)');title("Input Signal");
subplot(212),plot2d3(m,y);xlabel('n');ylabel('y(n)');title("Output Sinal");
Plots:

CONCLUSION :

Marks/Grade Teacher’s Signature

14
EXPERIMENT 3 Date ___________

Impulse Response & Total Response

AIM : (a) Find the convolution of x1(n) = {1,2,2,1} and x2(n)={1,-1,2}.

↑ ↑

(b) Determine the impulse response and total response of the following
LTI systems described by their constant coefficient difference equation.
Assume y(-1)=1 and y(-2)=0 for all systems

(i) y(n) – ½ y(n-1) = x(n) with x(n) = 10cos(πn/4)u(n)

(ii) y(n) – 0.9 y(n-1) + 0.81 y(n-2) = x(n) with x(n) = u(n)

(iii) y(n)-4y(n-1)-4y(n-2) = x(n)+3x(n-1) with x(n) = n (1/3)n u(n)

Part (a) Convolution

Program:

clc;
clear;
close;
x1=[1 2 2 1]; //sequence x1 and its index n1
n1=-2:1;
x2=[1 -1 2]; //sequence x2 and its index n2
n2=-2:0;
y=conv(x1,x2); // Convolution
// lower index of y=lower index of x1 + lower index of x2
// higher index of y= higher index of x1 + higher index of x2
m=min(n1)+min(n2):max(n1)+max(n2); // index of y
plot2d3(m,y); xlabel("n");ylabel("y(n)");title("Output of Convolution");

15
Plots:

Part (b)(i) Impulse response of y(n) – ½ y(n-1) = x(n)

Program:

clc;
close;
clear;
N=20;
n=0:N-1;
x1=[1 zeros(1,N-1)]; //Impulse signal
yint1=1; // y(-1)=0
xint1=0;
y=[yint1 zeros(1,N)];
x=[xint1 x1];

16
for i=1:N
y(i+1)=(1/2)*y(i)+x(i+1); //y(n) – ½ y(n-1) = x(n)
end
m=-1:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
subplot(2,1,2),plot2d3(m,y);xlabel("n");ylabel("y(n)");title("Output Signal");
Plots:

Part (b)(i) Total response of y(n) – ½ y(n-1) = x(n) with x(n) =

10cos(πn/4)u(n)

Program:

clc;
close;
clear;
N=20;

17
n=0:N-1;
x1=10*cos(n*%pi/4);//input signal
yint1=1; // y(-1)=0
xint1=0;
y=[yint1 zeros(1,N)];
x=[xint1 x1];
for i=1:N
y(i+1)=(1/2)*y(i)+x(i+1); //y(n) – ½ y(n-1) = x(n)
end
m=-1:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
subplot(2,1,2),plot2d3(m,y);xlabel("n");ylabel("y(n)");title("Output Signal");
Plots:

18
 Part (b)(ii) Impulse response of y(n) – 0.9 y(n-1) + 0.81 y(n-2) = x(n)

Program:

clc;
close;
clear;
N=20;
n=0:N-1;
x1=[1 zeros(1,N-1)]; //Impulse signal
yint1=1;
yint2=0;
xint1=0;
xint2=0;
y=[yint2 yint1 zeros(1,N)];
x=[xint2 xint1 x1];
for i=1:N
y(i+2)=0.9*y(i+1)-0.81*y(i)+x(i+2);//y(n) – 0.9 y(n-1) + 0.81 y(n-2) = x(n)
end
m=-2:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
subplot(2,1,2),plot2d3(m,y);xlabel("n");ylabel("y(n)");title("Output Signal");

19
Plots:

Part (b)(ii) Total response of y(n) – 0.9 y(n-1) + 0.81 y(n-2) = x(n) with x(n)
= u(n)

Program:

clc;close;clear;
N=20;
n=0:N-1;
x1=ones(1,N); //Input signal x(n)=u(n)
yint1=1;
yint2=0;
xint1=0;
xint2=0;
y=[yint2 yint1 zeros(1,N)];
x=[xint2 xint1 x1];

20
for i=1:N
y(i+2)=0.9*y(i+1)-0.81*y(i)+x(i+2); //y(n) – 0.9 y(n-1) + 0.81 y(n-2) = x(n)
end
m=-2:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
Plot:

Part (b)(iii) Impulse response of y(n)-4y(n-1)-4y(n-2) = x(n)+3x(n-1)

Program:

clc;close;clear;
N=20;
n=0:N-1;
x1=[1 zeros(1,N-1)]; //Impulse Signal
yint1=1;

21
yint2=0;
xint2=0;
xint1=0;
y=[yint2 yint1 zeros(1,N)];
x=[xint2 xint1 x1];
for i=1:N
//Equation y(n)-4y(n-1)-4y(n-2)=x(n)+3x(n-1)
y(i+2)=4*y(i+1)+4*y(i)+x(i+2)+3*x(i+1);
end
m=-2:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
subplot(2,1,2),plot2d3(m,y);xlabel("n");ylabel("y(n)");title("Output Signal");
Plot:

22
Part (b)(iii) Total response of y(n)-4y(n-1)-4y(n-2) = x(n)+3x(n-1) with input
x(n) = n (1/3)n u(n)

Program:

clc;close;clear;
N=20;
n=0:N-1;
x1=n.*(1/3)^n.*ones(1,N);//[1 zeros(1,N-1)];//n.*(1/3)^n.*ones(1,N);
yint1=1;
yint2=0;
xint2=0;
xint1=0;
y=[yint2 yint1 zeros(1,N)];
x=[xint2 xint1 x1];
for i=1:N
y(i+2)=4*y(i+1)+4*y(i)+x(i+2)+3*x(i+1); //y(n)-4y(n-1)-4y(n-2) =
x(n)+3x(n-1)
end
m=-2:N-1;
subplot(2,1,1),plot2d3(m,x);xlabel("n");ylabel("x(n)");title("Input Signal");
subplot(2,1,2),plot2d3(m,y);xlabel("n");ylabel("y(n)");title("Output Signal");

23
Plot:

CONCLUSION :

Marks/Grade Teacher’s Signature

24
EXPERIMENT 4 Date ___________

Application of Correlation
AIM: (a) For time delay estimation in radar. Let.
x(n) = {0,1,2,3,2,1,0}
↑
be a transmitted signal and y(n) be the received y(n) = α x(n‐D) +
w(n) with α=0.5.
(i) Plot ryx(l) and rwx(l). How can we determine from the plot that
target is present?
(ii) Explain how we can measure the delay D by computing the
cross-correlation ryx(l).
(b) Let one periodic sequence x(n) = cos(nπ/4) is corrupted by noise as
in equation y(n) = x(n) + w(n) estimated the period of noisy signal
using auto-correlation.

Part (a) Time Delay Estimation in Radar

Program:

clc;clear;close;
x=[0 1 2 3 2 1 0];// original signal x
n=-3:3;
D=10; // Delay betweens Tx and Rx signal
alpha=0.5; //Attenuation of signal x in channel
y1=alpha*[zeros(1,D) x]; // y(n-D) time delay D
N=length(y1);
w=0.05*rand(1,N); //Noise signal (max value 0.05)
y=y1+w;
/*cross correlation of ouput signal y with original x when
object is present*/
ryx=xcorr(y,x);
/* cross correlation of output signal w (noise) with
original signal x when object is absent*/
rwx=xcorr(w,x);////cross correlation

25
l=-(N-1):(N-1);
subplot(211), plot2d3(l,ryx);
xlabel("l");ylabel("ryx(l)");title("When object is present");
subplot(212), plot2d3(l,rwx);
xlabel("l");ylabel("rwx(l)");title("When object is absent");
Plots:

Part (a) Period Estimation in Noisy periodic signal

Program:

clc;
clear;
close;
N=20;
n=0:N-1;

26
L=8; //period of signal
f=1/L;
x=cos(2*%pi*f*n); //periodic signal with period L
w=0.05*rand(1,length(x)); //Noise signal (max value 0.05)
y=x+w;
ryy=xcorr(y,y); //Auto-correlation of y
l=-(N-1):(N-1);
plot2d3(l,ryy);xlabel("l"); ylabel("ryy(l)"); title("Auto-correlation of y(n)");
Plots:

CONCLUSION :

Marks/Grade Teacher’s Signature

27
EXPERIMENT 5 Date ___________

Z-transform

AIM : Find out the z-transform of following finite sequence and write

the statement of ROC

(i) x1(n)={1, 0, 3, -1, -2}
↑
(ii) x2(n)={-3,-2,-1,0,1}
↑
(iii) x3(n)={2,-1,3,2,1,0,2,3}
↑
(iv) x4(n) = 𝛅(n)
(v) x5(n) = 𝛅(n-5)
(vi) x6(n) = 𝛅(n+5)

Program:

clc;
clear;
z=%z;
x1=[1 0 3 -1 -2];
n1=0:4;
x1z=sum(x1.*z^(-n1));
disp("ztransform of x1 is ",x1z);
disp("ROC is entire z-plane except z=0");
x2=[-3 -2 -1 0 1];
n2=-4:0;
x2z=sum(x2.*z^(-n2));
disp("ztransform of x2 is ",x2z);
disp("ROC is entire z-plane except z=Inf");
x3=[2 -1 3 2 1 0 2 3];
n3=-2:5;

28
x3z=sum(x3.*z^(-n3));
disp("ztransform of x3 is ",x3z);
disp("ROC is entire z-plane except z=0 and z=Inf");
x4=[0 1 0];
n4=-1:1;
x4z=sum(x4.*z^(-n4));
disp("ztransform of x3 is ",x4z);
disp("ROC is entire z-plane");
x5=[zeros(1,5) 1];
n5=0:5;
x5z=sum(x5.*z^(-n5));
disp("ztransform of x5 is ",x5z);
disp("ROC is entire z-plane except z=0");
x6=[1 zeros(1,5) ];
n6=-5:0;
x6z=sum(x6.*z^(-n6));
disp("ztransform of x5 is ",x6z);
disp("ROC is entire z-plane except z=Inf");
Results:

"ztransform of x1 is "
-2 -z +3z² +z⁴
--------------
z⁴
"ROC is entire z-plane except z=0"

"ztransform of x2 is "
1 -z² -2z³ -3z⁴
"ROC is entire z-plane except z=Inf"

29
"ztransform of x3 is "
3 +2z +z³ +2z⁴ +3z⁵ -z⁶ +2z⁷
----------------------------
z⁵
"ROC is entire z-plane except z=0 and z=Inf"

"ztransform of x3 is "
1
-
1
"ROC is entire z-plane"

"ztransform of x5 is "
1
--
z⁵
"ROC is entire z-plane except z=0"

"ztransform of x5 is "
z⁵
"ROC is entire z-plane except z=Inf"

CONCLUSION :

Marks/Grade Teacher’s Signature

30
EXPERIMENT 6 Date ___________

Pole-zero plot on z-plane

AIM :Find out poles and zeros of the following systems and display

them on the z-plane. If all the systems are causal, check their
stability.
−1
1+3𝑧
(i) H(z) = −1 −2
1+3𝑧 +2𝑧
−2
1+2𝑧
(ii) H(z) = −2
1+𝑧

1
(iii) H(z) = −1 −1 2
(1−2𝑧 ) (1−𝑧 )
−1 −2
1+0.5𝑧 +𝑧
(iv) H(z) = −1 −2
1+0.4𝑧 +0.4𝑧
−1
1+3𝑧
Part (i) H(z) = −1 −2
1+3𝑧 +2𝑧

Program:
clear;clc ;close;
z=%z;
HZ=(1+3*z^-1)/(1+3*z^-1+2*z^-2);
Zk=roots(HZ(2));
Pk=roots(HZ(3));
plzr(HZ);
disp(Zk,"Zeros of filter are");
disp(Pk,"Poles of filter are");
Results:
Zeros of filter are

- 3.

Poles of filter are

- 2.

- 1.

31
Plots:

−2
1+2𝑧
Part (ii) H(z) = −2
1+𝑧

Program:
clear;clc ;close;
z=%z;
HZ=(1+2*z^-2)/(1+z^-2);
Zk=roots(HZ(2));
Pk=roots(HZ(3));
plzr(HZ);
disp(Zk,"Zeros of filter are");
disp(Pk,"Poles of filter are");
Results:

Zeros of filter are

1.4142136i

- 1.4142136i

Poles of filter are

32
 -i

Plots:

1
Part (iii) H(z) = −1 −1 2
(1−2𝑧 ) (1−𝑧 )

Program:
clear;clc ;close;
z=%z;
HZ=1/((1-2*z^-1)*(1-z^-1)^2);
Zk=roots(HZ(2));
Pk=roots(HZ(3));
plzr(HZ);
disp(Zk,"Zeros of filter are");
disp(Pk,"Poles of filter are");

33
Results:

Zeros of filter are

Poles of filter are

2.

1.0000000

1.

Plots:

34
 −1 −2
1+0.5𝑧 +𝑧
Part (iv) X(z) = −1 −2
1+0.4𝑧 +0.4𝑧

Program:
clear;clc ;close;
z=%z
HZ=(1+0.5*z^-1+z^-2)/(1+0.4*z^-1+0.4*z^-2);
Zk=roots(HZ(2));
Pk=roots(HZ(3));
plzr(HZ);
disp(Zk,"Zeros of filter are");
disp(Pk,"Poles of filter are");
Results:
Zeros of filter are
- 0.25 + 0.9682458i
- 0.25 - 0.9682458i
Poles of filter are
- 0.2 + 0.6i
- 0.2 - 0.6i
Plots:

35
CONCLUSION :

Marks/Grade Teacher’s Signature

36
EXPERIMENT 7 Date ___________

Frequency response of filters

AIM : Draw magnitude response and phase response of the systems

having transfer function H(ω) given below:

(i) 3-tap Moving Average system

( )
−𝑗ω
1−𝑒
(ii) H(ω) = −𝑗ω
1+0.9𝑒

( )
−𝑗2ω
1−𝑒
(iii) H(ω) = −𝑗2ω
1+0.9𝑒
−1 −2
1+6𝑧 +𝑧
(iv) 𝐻(𝑧) = −1 −2 −1
(1−2𝑧 +2𝑧 )(1−0.5𝑧 )
−1 −2
0.9382(1−1.618𝑧 +𝑧 )
(v) 𝐻(𝑧) = −1 −2
(1−1.5162𝑧 +0.8782𝑧 )

Identify the type of each filter.

Part (i) 3-tap Moving Average system

Program:

clear;
clc ;
close ;
w=-%pi:0.01:%pi;
j=%i;
H=1/3*(1+exp(-j*w)+exp(-2*j*w));
//caluculation of Phase and Magnitude of H
[phase_H,m]=phasemag(H);
Hm=abs(H);
subplot(2,1,1);
plot(w,Hm);
xlabel('Frequency in Radians')

37
ylabel('abs(Hm)');
title('MAGNITUDE RESPONSE');
subplot(2,1,2);
plot(w,phase_H);
xlabel('Frequency in Radians');
ylabel('<(H)');
title('PHASE RESPONSE');
Plots:

( )
−𝑗ω
1−𝑒
Part (ii) H(ω) = −𝑗ω
1+0.9𝑒

Program:

clear;clc ;close ;
w=-%pi:0.01:%pi;
j=%i;

38
H=(1-exp(-j*w))./(1+0.9*exp(-j*w));
//caluculation of Phase and Magnitude of H
[phase_H,m]=phasemag(H);
Hm=abs(H);
subplot(2,1,1);
plot2d(w,Hm);
xlabel('Frequency in Radians')
ylabel('abs(Hm)');
title('MAGNITUDE RESPONSE');
subplot(2,1,2);
plot2d(w,phase_H);
xlabel('Frequency in Radians');
ylabel('<(H)');
title('PHASE RESPONSE');
Plots:

39
 ( )
−𝑗2ω
1−𝑒
Part (iii) H(ω) = −𝑗2ω
1+0.9𝑒

Program:
clear;
clc;
close ;
w=-%pi:0.01:%pi;
j=%i;
H=(1-exp(-2*j*w))./(1+0.9*exp(-2*j*w));
[phase_H,m]=phasemag(H);
Hm=abs(H);
subplot(2,1,1);
plot2d(w,Hm);
xlabel('Frequency in Radians')
ylabel('abs(Hm)');
title('MAGNITUDE RESPONSE');
subplot(2,1,2);
plot2d(w,phase_H);
xlabel('Frequency in Radians');
ylabel('<(H)');
title('PHASE RESPONSE');

40
Plot:

−1 −2
1+6𝑧 +𝑧
Part (iv) 𝐻(𝑧) = −1 −2 −1
(1−2𝑧 +2𝑧 )(1−0.5𝑧 )

Program:
clear;
clc ;
close ;
w=-%pi:0.01:%pi;
j=%i;
H=(1+6*exp(-j*w)+exp(-j*2*w))./((1-2*exp(-j*w)+2*exp(-2*j*w)).*(1-0.5*ex
p(-j*w)));
//caluculation of Phase and Magnitude of H
[phase_H,m]=phasemag(H);
Hm=abs(H);
subplot(2,1,1);

41
plot2d(w,Hm);
xlabel('Frequency in Radians')
ylabel('abs(Hm)');
title('MAGNITUDE RESPONSE');
subplot(2,1,2);
plot2d(w,phase_H);
xlabel('Frequency in Radians');
ylabel('<(H)');
title('PHASE RESPONSE');
Plot:

42
 −1 −2
0.9382(1−1.618𝑧 +𝑧 )
Part (iv) 𝐻(𝑧) = −1 −2
(1−1.5162𝑧 +0.8782𝑧 )

Program:
clear;
clc ;
close ;
w=-%pi:0.01:%pi;
j=%i;
H=0.932*(1-1.618*exp(-j*w)+exp(-j*2*w))./((1-1.5162*exp(-j*w)+0.8782*exp
(-2*j*w)));
//caluculation of Phase and Magnitude of H
[phase_H,m]=phasemag(H);
Hm=abs(H);
subplot(2,1,1);
plot2d(w,Hm);
xlabel('Frequency in Radians')
ylabel('abs(Hm)');
title('MAGNITUDE RESPONSE');
subplot(2,1,2);
plot2d(w,phase_H);
xlabel('Frequency in Radians');
ylabel('<(H)');
title('PHASE RESPONSE');

43
Plot:

CONCLUSION :

Marks/Grade Teacher’s Signature

44
EXPERIMENT 8 Date ___________

Discrete Fourier Transform

AIM : (a) Determine DFT of the sequence x(n)= {1,-1,-2,3,-1} and find

out IDFT of the resulting sequence.

(b)Find circular convolution of two sequences x1(n)={2,1,2,1} and
x2(n)={1,2,3,4} using DFT-IDFT method
(c) Find linear convolution of two sequences x1(n)= {2,1,2,1} and
x2(n)={1,2,3,4} using DFT-IDFT method

Part(a) DFT of the sequence x1(n)= {1,-1,-2,3,-1}

Program:

clear;
clc;
close;
x = [1 -1 -2 3 -1];
//Computing DFT
X = fft(x,-1);
disp("DFT of sequence x is",X);
//Computing IDFT
x_inv = fft(X,1);
disp("IDFT of above sequence is:",x_inv);

Results:

"DFT of sequence x is"

column 1 to 3
0. + 0.i -0.427051 + 2.9389263i 2.927051 - 4.7552826i
column 4 to 5

45
 2.927051 + 4.7552826i -0.427051 - 2.9389263i

"IDFT of above sequence is:"

1. -1. -2. 3. -1.

Part(b) circular convolution of two sequences x1(n)={2,1,2,1} and

x2(n)={1,2,3,4}
Program:

clear;
clc;
close;
x1 = [2 1 2 1];
x2 = [1 2 3 4];
//Computing DFT
X1 = fft(x1,-1);
disp("DFT of x1",X1);
X2 = fft(x2,-1);
disp("DFT of x2",X2);
//Multiplication of 2 DFTs
X3 = X1.*X2
disp("Multiplication of 2 DFTs",X3)
//inverse DFT of X3
x3 =fft(X3,1);
disp("Circular Convolution of x1 and x2",x3);
Results:

"DFT of x1"
6. 0. 2. 0.

46
 "DFT of x2"
10. + 0.i -2. + 2.i -2. + 0.i -2. - 2.i

"Multiplication of 2 DFTs"
60. + 0.i 0. + 0.i -4. + 0.i 0. + 0.i

"Circular Convolution of x1 and x2"

14. 16. 14. 16.
Part(b) linear convolution of two sequences x1(n)={2,1,2,1} and
x2(n)={1,2,3,4}
Program:

clear;
clc;
close;
x1 = [2 1 2 1];
x2 = [1 2 3 4];
N1 = length(x1)
N2 = length(x2)
//Length of linear convolution
N = N1+N2-1;
//Padding zeros to Make Length of x1 and 'x2
//Equal to length of linear convolution
x1p = [x1,zeros(1,N-N1)]
x2p = [x2,zeros(1,N-N2)]
//Computing DFT
X1 = fft(x1p,-1)
X2 = fft(x2p,-1)
//Multiplication of 2 DFTs
Y = X1.*X2
47
//Linear Convolution Result
y =fft(Y,1);
disp("Linear Convolution using DFT-IDFT is",y);
yc=conv(x1,x2);
disp("Linear Convolution using command conv",yc);

Results:

"Linear Convolution using DFT-IDFT is"

2. 5. 10. 16. 12. 11. 4.

"Linear Convolution using command conv"

2. 5. 10. 16. 12. 11. 4.

CONCLUSION :

Marks/Grade Teacher’s Signature

48
EXPERIMENT 9 Date ___________

Characteristic of Windows

AIM : Draw time domain and frequency domain characteristics

of different windows.

Program:

clc;
clear;
close;
N = 31;
n=0:N-1;
wre = window('re',N); // Rectangular window
wbr = window('tr',N); // Bartlett (triangle) window
whm = window('hm',N); // Hamming window
whn = window('hn',N); // Hanning window
wkr = window('kr', N, 6); // Kaiser window
subplot(121);
plot(n, [wre; wbr; whm; whn; wkr]);
xlabel("n");
ylabel("w(n)");
title("Time domain Characteristic");
L = 256; // L-point frequency response
[Wre,fr] = frmag(wre, L);
[Wbr,fr] = frmag(wbr, L);
[Whm,fr] = frmag(whm, L);
[Whn,fr] = frmag(whn, L);
[Wkr,fr] = frmag(wkr, L);
subplot(122);

49
plot(fr, 20*log10([Wre; Wbr; Whm; Whn; Wkr]))
xlabel("Normalized frequency");
ylabel("Magnitude (dB)");
legend(["Rectangular"; "Bartlett"; "Hamming"; "Hanning"; "Kaiser"]);
title("Frequency domain Characteristic");

Plots:

CONCLUSION :

Marks/Grade Teacher’s Signature

50
EXPERIMENT 10 Date ___________

FIR filter design

AIM : (a) Design 50th order digital FIR linear phase lowpass filter of
π
ωc= 4
using Rectangular, Bartlett, Hamming and Hanning
window

(b) Design the 50th order and 100th order digital FIR linear
phase bandpass filter using rectangular window, having lower
π 3π
cutoff frequency ωcl = 4 and higher cutoff frequency ωch = 8

Part(a)
π
FIR linear phase lowpass filter of ωc= 4
using Rectangular

Program:

clc;clear;close;
wc=%pi/4;//cutoff frequency rad/sec
fc=wc/(2*%pi);//normalized cutoff frequency
ftype='lp'//string : 'lp','hp','bp','sb' (filter type)
N=50;// order of filter
n=0:N-1;
wtype="re";
wp=[0 0];//2-vector of window parameters. Kaiser window fpar(1)>0 fpar(2)=0.
//for other types of window this parameters are [0 0]
[h,hm,fr]=wfir(ftype,N,[fc 0],wtype,wp);
subplot(1,2,1);
//plot frequency response of filter
plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);

51
//plot impulse response of filter
plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");
Plots:

π
FIR linear phase lowpass filter of ωc= 4
using Bartlett

Program:

clc;clear;close;
wc=%pi/4;//cutoff frequency rad/sec
fc=wc/(2*%pi);//normalized cutoff frequency
ftype='lp'//string : 'lp','hp','bp','sb' (filter type)
N=50;// order of filter
n=0:N-1;
wtype="tr";

52
wp=[0 0];//2-vector of window parameters. Kaiser window fpar(1)>0 fpar(2)=0.
//for other types of window this parameters are [0 0]
[h,hm,fr]=wfir(ftype,N,[fc 0],wtype,wp);
subplot(1,2,1);
plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);
plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");
Plots:

π
FIR linear phase lowpass filter of ωc= 4
using Hamming

Program:

clc;clear;close;

53
wc=%pi/4;//cutoff frequency rad/sec
fc=wc/(2*%pi);//normalized cutoff frequency
ftype='lp'//string : 'lp','hp','bp','sb' (filter type)
N=50;// order of filter
n=0:N-1;
wtype="hm";
wp=[0 0];
[h,hm,fr]=wfir(ftype,N,[fc 0],wtype,wp);
subplot(1,2,1);
plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);
plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");
Plot:

54
 π
FIR linear phase lowpass filter of ωc= 4
using Hanning

Program:

clc;clear;close;
wc=%pi/4;//cutoff frequency rad/sec
fc=wc/(2*%pi);//normalized cutoff frequency
ftype='lp'//string : 'lp','hp','bp','sb' (filter type)
N=50;// order of filter
n=0:N-1;
wtype="hn";
wp=[0 0];//
[h,hm,fr]=wfir(ftype,N,[fc 0],wtype,wp);
subplot(1,2,1);plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");
Plot:

55
 Part(b)
π
50th order FIR linear phase bandpass filter with ωcl = 4 and higher cutoff
3π
frequency ωch = 8

Program:

clc;clear;close;
wcl=%pi/4; //lower cutoff frequency rad/sec
wch=3*%pi/4; //higher cutoff frequency rad/sec
fcl=wcl/(2*%pi); //normalized lower cutoff frequency
fch=wch/(2*%pi); //normlized higher cutoff frequency
ftype='bp' //string : 'lp','hp','bp','sb' (filter type)
N=50; //order of filter
n=0:N-1;
wtype="re";
wp=[0 0];//2-vector of window parameters. Kaiser window fpar(1)>0 fpar(2)=0.
//for other types of window this parameters are [0 0]
[h,hm,fr]=wfir(ftype,N,[fcl fch],wtype,wp);
subplot(1,2,1);
plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);
plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");

56
Plot:

π
100th order FIR linear phase bandpass filter with ωcl = 4 and higher cutoff
3π
frequency ωch = 8

Program:

clc;clear;close;
wcl=%pi/4; //lower cutoff frequency rad/sec
wch=3*%pi/4; //higher cutoff frequency rad/sec
fcl=wcl/(2*%pi); //normalized lower cutoff frequency
fch=wch/(2*%pi); //normlized higher cutoff frequency
ftype='bp' //string : 'lp','hp','bp','sb' (filter type)
N=100; //order of filter
n=0:N-1;
wtype="re";

57
wp=[0 0];//2-vector of window parameters. Kaiser window fpar(1)>0 fpar(2)=0.
//for other types of window this parameters are [0 0]
[h,hm,fr]=wfir(ftype,N,[fcl fch],wtype,wp);
subplot(1,2,1);
plot(fr,hm);xlabel("Normalized frequency");ylabel("|H(w)|");
title("Frequency Response of filter");
subplot(1,2,2);
plot2d3(n,h);xlabel("n");ylabel("h(n)");
title("Impulse Response of filter");
Plot:

CONCLUSION :

Marks/Grade Teacher’s Signature

EXPERIMENT 11 Date ___________

58
 IIR filter design

AIM : (a) Design a 7th order digital IIR lowpass Butterworth,

Chebyshev-I and II prototype filter with a cutoff frequency ωc
π
= 4 , passband ripple of 0.1, stopband ripple of 0.1 using
bilinear transformations.

(b) Consider the analog resonator having a resonance

frequency of 3 rad/sec.

𝑠 +0.1
H(s)= 2
𝑠 +0.2𝑠+9.01

Design the digital IIR resonance filter with the resonance

frequency ω𝑟 = π/2. User bilinear transformation for
conversion. Also, plot their poles and zeros.

Part(a)

Program:

/* hz=iir(n,ftype,fdesign,frq,delta)
where,
n = positive number witn integer value, the filter order.
ftype = string specifying the filter type, the possible values are: 'lp' for
low-pass,'hp' for high pass,'bp' for band pass and 'sb' for stop band.
fdesign = string specifying the analog filter design, the possible values are:
'butt', 'cheb1', 'cheb2' and 'ellip'
frq = 2-vector of discrete cut-off frequencies (i.e., 0<frq<.5). For 'lp' and 'hp'
filters only frq(1) is used (in this case, frq can be a scalar). For 'bp' and 'sb'
filters frq(1) is the upper cut-off frequency and frq(2) is the lower cut-off
frequency.
delta = 2-vector of error values for cheb1, cheb2, and ellip filters where only
delta(1) is used for cheb1 case, only delta(2) is used for cheb2 case, and delta(1)
and delta(2) are both used for ellip case. 0<delta(1),delta(2)<1
- for cheb1 filters 1-delta(1)<ripple<1 in passband
- for cheb2 filters 0<ripple<delta(2) in stopband
- for ellip filters 1-delta(1)<ripple<1 in passband and 0<ripple<delta(2) in
stopband

hz = a single input single output discrete transfer function, the low pass filter

*/
59
clc;clear;close;
L=256; //L-point frequency response
Wn=1/4; //Normalized frequency Wc/pi
N=7; // Order of filter
Rp=0.1; //Ripple in passband delta(1)
Rs=0.1; //Ripple in stopband delta(2)
ftype='lp';
hz1=iir(N,ftype,'butt',[Wn 0],[0 0]); //Butterworth filter
disp("Transfer function of Butterworth filter is",hz1);
hz2=iir(N,ftype,'cheb1',[Wn 0],[Rp 0]); //Chebyshev-I filter
disp("Transfer function of Chebyshev-I filter is",hz2);
hz3=iir(N,ftype,'cheb2',[Wn 0],[0 Rs]); //Chebyshev-II filter
disp("Transfer function of Chebyshev-II filter is",hz3);
hz4=iir(N,ftype,'ellip',[Wn 0],[Rp Rs]); // Elliptical filter
disp("Transfer function of Elliptical filter is",hz4);
[hzm1,fr]=frmag(hz1,L);
[hzm2,fr]=frmag(hz2,L);
[hzm3,fr]=frmag(hz3,L);
[hzm4,fr]=frmag(hz4,L);
plot(fr,[hzm1; hzm2; hzm3; hzm4]);
xlabel("Normalized frequency");ylabel("|H(w)|");
title('Discrete IIR filter: low pass');
legend(["Butterworth","Chebyshev-I","Chebyshev-II","Elliptical"]);

60
Result
"Transfer function of Butterworth filter is"

0.0165653 +0.1159571z +0.3478712z² +0.5797853z³ +0.5797853z⁴

+0.3478712z⁵ +0.1159571z⁶ +0.0165653z⁷
------------------------------------------------------------------------
6.874D-35 +0.0076835z -3.278D-18z² +0.1927012z³ +8.084D-18z⁴
+0.919973z⁵
+2.082D-17z⁶ +z⁷

"Transfer function of Chebyshev-I filter is"

0.0041079 +0.0287554z +0.0862661z² +0.1437768z³ +0.1437768z⁴

+0.0862661z⁵ +0.0287554z⁶ +0.0041079z⁷
------------------------------------------------------------------------
-0.2260859 +0.8677288z -1.992739z² +3.1832893z³ -3.7048774z⁴
+3.361089z⁵
-1.9625924z⁶ +z⁷

"Transfer function of Chebyshev-II filter is"

0.2159855 +0.6262675z +1.2209053z² +1.638417z³ +1.638417z⁴

+1.2209053z⁵
+0.6262675z⁶ +0.2159855z⁷
-----------------------------------------------------------------------
0.0465689 +0.2252325z +0.6003597z² +1.1764307z³ +1.4821834z⁴
+1.8753842z⁵ +0.9969913z⁶ +z⁷

"Transfer function of Elliptical filter is"

0.2158778 +0.3198685z +0.7569526z² +0.8609813z³ +0.8609813z⁴

+0.7569526z⁵ +0.3198685z⁶ +0.2158778z⁷
-------------------------------------------------------------
-0.1769321 +0.7105228z -0.7995413z² +2.4147285z³ -1.0840763z⁴
+2.7041794z⁵ -0.4615206z⁶ +z⁷

61
Plot:

Part(b)

Program:

clear;clc;close;
Omegar=3; // Analog resonance freqnecy
Wr=%pi/2; // Digital resonance frequency
// Relationship between frequency variable in two domain
T=(2/Omegar) * tan(Wr/2)
L = 512; // L-point frequency response
s = poly(0,'s');
H = (s+0.1)/(s^2+0.2*s+9.01)
z = poly(0,'z');
Hz = horner(H,(2/T)*((1-z^-1)/(1+z^-1)))
[HW fr] =frmag(Hz(2),Hz(3),L);
62
subplot(1,3,1);plot(fr,20*log10(HW));xlabel("Normalized
frequency");ylabel("|H(w)|");
title("Digital Resonator");
subplot(1,3,2);plzr(H);title("Analog Resonator");
subplot(1,3,3);plzr(Hz);title("Digital Resonator");
disp("Transfer function of Digital Resonator is",Hz);
Result:
"Transfer function of Digital Resonator is"

-0.1558302 +0.0107469z +0.1665771z²

-----------------------------------
0.9355185 +0.0010747z +z²
Plot:

CONCLUSION :

63
Marks/Grade Teacher’s Signature

64