SAMPOORNA INSTITUTE OF TECHNOLOGY AND RESEARCH,

Belekere, Channapatna

DEPARTMENT OF E &C ENGINEERING

V SEMESTER

BEC502

DIGITAL SIGNAL PROCESSING LAB MANUAL

ACADEMIC YEAR: 2022 Scheme

NAME OF THE STUDENT:____________________________

UNIVERSITY SEAT NO :____________________________

BATCH :____________________________

PREPARED BY
LATHAMANI H.N
ASSISTANT PROFESSOR
DEPARTMENT OF ECE
SITAR
 LIST OF EXPERIMENTS

1. Program to generate the following discrete time signals.

a) Unit sample sequence, b)Unit step sequence, c) Exponential sequence, d)Sinusoidal
sequence, e)Random sequence
2. Program to perform the following operations on signals.
a) Signal addition, b) Signal multiplication, c)Scaling, d) Shifting, e)Folding
3. Program to perform convolution of two given sequences (without using built-in function) and
display the signals.
4. Consider a causal system y(n) = 0.9y(n-1)+x(n).

b) Plot |H(ejω)| and ∠ H(ejω)

a) Determine H(z) and sketch its pole zero plot.

c) Determine the impulse response h(n).

5. Computation of N point DFT of a given sequence (without using built-in function) and to plot the
magnitude and phase spectrum.
6. Using the DFT and IDFT, compute the following for any two given sequences
a)Circular convolution
b) Linear convolution
7. Verification of Linearity property, circular time shift property & circular frequency shift property
of DFT.
8. Develop decimation in time radix-2 FFT algorithm without using built-in functions.
9. Design and implementation of digital low pass FIR filter using a window to meet the given
specifications
10. Design and implementation of digital high pass FIR filter using a window to meet the given
specifications
11. Design and implementation of digital IIR Butterworth low pass filter to meet the given
specifications.
12.Design and implementation of digital IIR Butterworth high pass filter to meet the given
specifications
 DSP USING MATLAB

MATLAB:MATLAB is a 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. The MATLAB name
stands for MATrix Laboratory The diagram shows the main features and capabilities of
MATLAB.

MATLAB

GRAPHICS COMPUTATIONS EXTERNAL INTERFACE

2-D Graphics Linear Algebra Interface with C and
3 –D Graphics Signal Processing FORTAN language
Animations Polynomial and
interpolation
Quadrature solution
of ODE’S

TOOL BOX
Signal Processing Image Processing
Statistics Control Systems
Neural Network Communications

At its core ,MATLAB is essentially a set (a “toolbox”) of routines (called “m files”

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 and/or get data back out (i.e. the list of commands is like a function in
most programming languages). Once you save a function, it becomes part of your toolbox
(i.e. it now looks to you as if it were part of the basic toolbox that you started with).
For those with computer programming backgrounds: Note that MATLAB runs as
aninterpretive language (like the old BASIC). That is, it does not need to be compiled. It
simply reads through each line of the function, executes it, and then goes on to the next
line. (In practice, a form of compilation occurs when you first run a function, so that it can
run faster the next time you run it.

MATLAB Windows:
MATLAB works with through three 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

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
filescalled 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 –type in MATLAB is the array. It encompasses several
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 MATLAB. No dimension statements are
required for vectors or arrays .we can find the dimensions of an existing matrix or a vector
with the size and length commands.

BASIC INSTRUCTIONS IN MAT LAB

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
Therefore T = [0 1 2 3 4 5 6 7 8 9 10]

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

3. T= 0:1/pi: 1
Therefore 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. ones (5,2)
The above instruction creates a vector of five rows and two
columns Output = 11

11

11

11

11

6. zeros (2, 4)
Output = 0000
0000
7. a = [ 1 2 3] , b = [4 5 6]
 THEN a.*b = [4 10 18]

i.e. [4*1 5*2 6*3]

8. if C= [2 2 2] and b = [4 5 6]
b.*C results in [8 10 12]

9. plot (t, x)
If x = [6 7 8 9] and t = [1 2 3 4]; This instruction will display a figure window which
indicates the plot of x versus t

10. stem (t,x)

If x(n) = [0 1 2 3];This instruction will display a figure window as shown

11. Subplot: 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
5. Filter:Syntax: y = filter(b,a,X)
Description:y = filter(b,a,X) filters the data in vector X with the filter described by
numerator coefficient vector b and denominator coefficient vector a. If a(1) is not equal to
1, filter normalizes the filter coefficients by a(1). If a(1) equals 0, filter returns an error.

6. Fliplr: Syntax: B = fliplr(A)

Description:B = fliplr(A) returns A with columns flipped in the left-right direction,

that is, about a vertical axis.If A is a row vector, then fliplr(A) returns a vector of the same
length with the order of its elements reversed. If A is a column vector, then fliplr(A)
simply returns A.

7. Conv: Syntax: w = conv(u,v)

Description:w = conv(u,v) convolves vectors u and v. Algebraically, convolution is the

same operation as multiplying the polynomials whose coefficients are the elements of u
and v.

8. Impz:Syntax: [h,t] = impz(b,a,n)

Description:[h,t] = impz(b,a,n) computes the impulse response of the filter with
numerator coefficients b and denominator coefficients a and computes n samples of the
impulse response when n is an integer (t = [0:n-1]'). If n is a vector of integers, impz
computes the impulse response at those integer locations, starting the response
computation from 0 (and t = n or t = [0 n]).If, instead of n, you include the empty
vector for the second argument, the number of samples is computed automatically by
default.

9. 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.

10. xlabel:Syntax: xlabel('string')

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

11. ylabel: Syntax:ylabel('string')

Description:ylabel('string') labels the y-axis of the current axes.

12. Title: Syntax:title('string')

Description:title('string') outputs the string at the top and in the center of the current
axes.

13. grid on: Syntax:grid on

Description:grid on adds major grid lines to the current axes.
1. Program to generate the following discrete time signals.
a) Unit sample sequence, b)Unit step sequence, c) Exponential sequence, d)Sinusoidal
sequence, e)Random sequence
Unit step sequence
n = -10:10; % Generate the unit sample sequence
x = zeros(size(n));
x(n == 0) = 1; % Plot the unit sample sequence
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Unit Sample Sequence');
grid on;

Exponential sequence
n = -10:10;
a = 0.8; % Adjust the value of 'a' to change the decay rate
x = a.^n;
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Exponential Sequence');
grid on;

Sinusoidal sequence
n = 0:99;
A = 2; % Amplitude
w = pi/4; % Angular frequency
phi = pi/3; % Phase shift
x = A*cos(w*n + phi);
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Sinusoidal Sequence');
grid on;

Unit sample sequence

n = -10:10;
x = zeros(size(n));
x(n == 0) = 1;
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Unit Sample Sequence');
grid on;
Random sequence
N = 20; % Specify the length of the random sequence
x = randn(1, N); % Generate a random sequence of N samples
stem(0:N-1, x); % Plot the random sequence
xlabel('Time index (n)');
ylabel('Amplitude');
title('Random Sequence');
grid on;

RESULT:
2. Program to perform the following operations on signals.
a) Signal addition, b) Signal multiplication, c)Scaling, d) Shifting, e)Folding

Signal addition
n = 0:10; % Define the time index range
x1 = 2*sin(0.5*pi*n); % Define the two input signals
x2 = 3*cos(0.2*pi*n);
y = x1 + x2; % Add the two signals
subplot(3,1,1); % Plot the input and output signals
stem(n, x1);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x1');
grid on;

subplot(3,1,2);
stem(n, x2);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x2');
grid on;

subplot(3,1,3);
stem(n, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Sum of x1 and x2');
grid on;

Signal multiplication

n = 0:10; % Define the time index range

x1 = 2*sin(0.5*pi*n); % Define the two input signals
x2 = 3*cos(0.2*pi*n);
y = x1 .* x2; % Multiply the two signals
subplot(3,1,1); % Plot the input and output signals
stem(n, x1);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x1');
grid on;

subplot(3,1,2);
stem(n, x2);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x2');
grid on;

subplot(3,1,3);
stem(n, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Product of x1 and x2');
grid on;

Signal Scaling
n = 0:10; % Define the time index range
x = sin(0.5*pi*n); % Define the input signal
alpha = 2; % Scaling factor
y = alpha * x; % Scale the signal
subplot(2,1,1); % Plot the input and output signals
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Original Signal');
grid on;

subplot(2,1,2);
stem(n, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Scaled Signal');
grid on;

Time Shifting

n = 0:10; % Define the time index range

x = 2*sin(0.5*pi*n); % Define the input signal
shift = 3; % Define the shift amount (positive for right shift, negative for left shift)
y = circshift(x, shift); % Shift the signal
subplot(2,1,1); % Plot the input and shifted signals
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Original Signal');
grid on;

subplot(2,1,2);
stem(n, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Shifted Signal');
grid on;

Signal Folding

n = -5:5; % Define the time index range

x = [1 2 3 4 5 6 7 8 9 10 11]; % Define the input signal
y = fliplr(x); % Fold the signal
subplot(2,1,1); % Plot the input and folded signals
stem(n, x);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Original Signal');
grid on;

subplot(2,1,2);
stem(n, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Folded Signal');
grid on;

RESULT:
3. Program to perform convolution of two given sequences (without using built-in
function) and display the signals.
x1 = [1 2 3]; % Define the two input sequences
x2 = [4 5 6];
N = length(x1) + length(x2) - 1; % Determine the length of the output sequence
y = zeros(1, N); % Initialize the output sequence
for n = 0:N-1 % Perform convolution

for k = 0:n
y(n+1) = y(n+1) + x1(k+1)*x2(n-k+1);
end
end

subplot(3,1,1); % Plot the input and output sequences

stem(0:length(x1)-1, x1);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x1');
grid on;

subplot(3,1,2);
stem(0:length(x2)-1, x2);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Signal x2');
grid on;

subplot(3,1,3);
stem(0:N-1, y);
xlabel('Time index (n)');
ylabel('Amplitude');
title('Convolution of x1 and x2');
grid on;

RESULT:
4. Consider a causal system y(n) = 0.9y(n-1)+x(n).

b) Plot |H(ejω)| and ∠ H(ejω)

a) Determine H(z) and sketch its pole zero plot.

c) Determine the impulse response h(n).

a) Determining H(z) and Pole-Zero Plot

Taking the Z-transform of the given difference equation:

Y(z) = 0.9z^(-1)Y(z) + X(z)

Rearranging to solve for the transfer function H(z):

H(z) = Y(z)/X(z) = 1 / (1 - 0.9z^(-1))

Pole-Zero Plot: The pole of the system is at z = 0.9. There are no zeros.

b) Plotting |H(ejω)| and ∠ H(ejω)

To plot the magnitude and phase response, we substitute z = e^(jω) into H(z):

H(e^jω) = 1 / (1 - 0.9e^(-jω))

Magnitude Response |H(e^jω)|:

|H(e^jω)| = |1 / (1 - 0.9e^(-jω))|

Phase Response ∠H(e^jω):

∠H(e^jω) = -∠(1 - 0.9e^(-jω))

c) Determining the Impulse Response h(n)

The impulse response h(n) is the inverse Z-transform of H(z). For a system with a single pole at
z = a, the impulse response is given by:

h(n) = a^n * u(n)

MATLAB implementation

% Define the frequency range

omega = linspace(-pi, pi, 1000);

% Calculate the frequency response

H = 1 ./ (1 - 0.9*exp(-1j*omega));
% Plot the magnitude response

figure;

plot(omega, abs(H));

xlabel('Frequency (rad/sample)');

ylabel('Magnitude |H(e^jω)|');

title('Magnitude Response');

% Plot the phase response

figure;

plot(omega, angle(H));

xlabel('Frequency (rad/sample)');

ylabel('Phase ∠H(e^jω)');

title('Phase Response');

RESULT:
5. Computation of N point DFT of a given sequence (without using built-in function) and
to plot the magnitude and phase spectrum.

The Discrete Fourier Transform (DFT) is a mathematical transform that converts a sequence of
N complex numbers into another sequence of N complex numbers. The DFT is defined as:

X(k) = Σ(n=0 to N-1) x(n) * exp(-j*2π*k*n/N)

where:

 X(k) is the k-th DFT coefficient

 x(n) is the n-th sample of the input sequence
 N is the length of the input sequence

MATLAB Implementation:

function X = my_dft(x)

N = length(x);

X = zeros(1, N);

for k = 0:N-1

for n = 0:N-1

X(k+1) = X(k+1) + x(n+1) * exp(-1j*2*pi*k*n/N);

end

end

end
% Example usage:

x = [1 2 3 4];

X = my_dft(x);

% Plotting magnitude and phase spectrum

magnitude = abs(X);

phase = angle(X);

subplot(2,1,1);

stem(0:length(X)-1, magnitude);

xlabel('Frequency Bin');

ylabel('Magnitude');

title('Magnitude Spectrum');

grid on;

subplot(2,1,2);

stem(0:length(X)-1, phase);

xlabel('Frequency Bin');

ylabel('Phase (radians)');

title('Phase Spectrum');

grid on;

RESULT:
6. Using the DFT and IDFT, compute the following for any two given sequences
a)Circular convolution
b) Linear convolution

Computing Circular and Linear Convolution Using DFT and IDFT

Circular Convolution

Method:

1. Zero-pad the sequences: Ensure both sequences have the same length, N. If they are
not, zero-pad the shorter sequence.
2. Compute the DFT of both sequences: Use the DFT formula or a built-in FFT function.
3. Multiply the DFT coefficients: Multiply the corresponding DFT coefficients of the two
sequences.
4. Compute the IDFT of the product: Use the IDFT formula or a built-in IFFT function
to obtain the circular convolution result.

Linear Convolution

Method:

1. Zero-pad the sequences: Pad the sequences with zeros to a length equal to the sum of
their lengths minus one.
2. Compute the DFT of both zero-padded sequences.
3. Multiply the DFT coefficients: Multiply the corresponding DFT coefficients of the two
sequences.
4. Compute the IDFT of the product: Use the IDFT to obtain the linear convolution
result.

MATLAB CODE

function [y_circular, y_linear] = circular_linear_convolution(x1, x2)

 % Circular Convolution
X1 = fft(x1);
X2 = fft(x2);
Y_circular = ifft(X1 .* X2);
y_circular = Y_circular(:); % Ensure column vector
% Linear Convolution using Overlap-Save Method
N = length(x1) + length(x2) - 1;
x1_padded = [x1 zeros(1, N-length(x1))];
x2_padded = [x2 zeros(1, N-length(x2))];
% Divide the input sequences into blocks of length N
L = floor(N/2);
num_blocks = ceil(length(x1_padded)/L);
y_linear = zeros(1, N);
for i = 1:num_blocks
x1_block = x1_padded((i-1)*L+1:min(i*L, end));
x2_block = x2_padded((i-1)*L+1:min(i*L, end));
% Circular convolution of blocks
Y_block = ifft(fft(x1_block) .* fft(x2_block));
y_linear((i-1)*L+1:min(i*L+L-1, N)) = y_linear((i-1)*L+1:min(i*L+L-1, N)) +
Y_block(:);
end
end

RESULT:

7. Verification of Linearity property, circular time shift property & circular frequency
shift property of DFT.
% Define a sequence
x = [1 2 3 4];

% DFT of x
X = fft(x);

% Linearity Verification
y = 2*x + 3*[5 6 7 8];
Y = fft(y);
linearity_check = abs(Y - (2*X + 3*fft([5 6 7 8]))) < 1e-10;
% Circular Time Shift Verification
x_shifted = circshift(x, 1);
X_shifted = fft(x_shifted);
time_shift_check = abs(X_shifted - exp(-1i*2*pi*(0:3)/4)*X) < 1e-10;

% Circular Frequency Shift Verification

X_shifted_freq = circshift(X, 1);
x_shifted_time = ifft(X_shifted_freq);
freq_shift_check = abs(x_shifted_time - exp(1i*2*pi*(0:3)/4)*x) < 1e-10;

% Display results
disp('Linearity Verification:');
disp(linearity_check);
disp('Circular Time Shift Verification:');
disp(time_shift_check);
disp('Circular Frequency Shift Verification:');
disp(freq_shift_check);

RESULT:

8. Develop decimation in time radix-2 FFT algorithm without using built-in functions.

function X = dit_fft(x)
N = length(x);
if N == 1
X = x;
return;
end

% Bit-reversal permutation
x = bitrevorder(x);

% Divide and conquer

x_even = x(1:2:N-1);
x_odd = x(2:2:N);
X_even = dit_fft(x_even);
X_odd = dit_fft(x_odd);

% Combine the results using twiddle factors

W = exp(-2*pi*1i/N);
k = 0:(N/2-1);
W_k = W.^k;
X = [X_even + W_k.*X_odd, X_even - W_k.*X_odd];
End

RESULT:

9. Design and implementation of digital low pass FIR filter using a window to meet the given
specifications.
function y = lowpass_fir(x, fp, fs, Fs, As, Rp)
% Calculate normalized frequencies
wp = 2*pi*fp/Fs;
ws = 2*pi*fs/Fs;

% Determine filter order

delta_w = ws - wp;
N = ceil((As-7.95)/(14.36*delta_w));

% Design ideal low-pass filter impulse response

n = 0:N-1;
hd = sin(wp*n)/(pi*n);
hd(n==0) = wp/pi;

% Apply Hamming window

w = hamming(N);
h = hd .* w';

% Filter the input signal

y = filter(h, 1, x);
end

RESULT: