1) Program to generate the following discrete time signals in octave a) unit sample sequence
B) unit step response C) Exponential Sequence D) Sinusoidal Sequence E) Random
sequence
clc;
clear all;
close all;
% Time index
n = -10:10;
% a) Unit Sample Sequence (delta[n])
unit_sample = (n == 0);
subplot(3,2,1);
stem(n, unit_sample, 'filled');
title('a) Unit Sample Sequence \delta[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
% b) Unit Step Sequence (u[n])
unit_step = (n >= 0);
subplot(3,2,2);
stem(n, unit_step, 'filled');
title('b) Unit Step Sequence u[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
% c) Exponential Sequence (a^n u[n]), a = 0.8
a = 0.8;
exp_seq = (a .^ n) .* (n >= 0);
subplot(3,2,3);
stem(n, exp_seq, 'filled');
title('c) Exponential Sequence (0.8^n u[n])');
xlabel('n'); ylabel('Amplitude');
grid on;
% d) Sinusoidal Sequence (sin(0.2πn))
sin_seq = sin(0.2 * pi * n);
subplot(3,2,4);
stem(n, sin_seq, 'filled');
title('d) Sinusoidal Sequence sin(0.2\pin)');
xlabel('n'); ylabel('Amplitude');
grid on;
% e) Random Sequence (between 0 and 1)
rand_seq = rand(1, length(n));
subplot(3,2,5);
stem(n, rand_seq, 'filled');
title('e) Random Sequence');
xlabel('n'); ylabel('Amplitude');
grid on;
Output
2) program to perform following operations on signals. a) Signal Addition B) Signal
Multiplication C) Scaling D) Shifting E) Folding.
clc;
clear all;
close all;
% Define two signals x1[n] and x2[n]
n = -5:5;
x1 = sin(0.2 * pi * n); % Signal 1: Sinusoidal
x2 = (n >= 0); % Signal 2: Unit Step
% a) Signal Addition: y[n] = x1[n] + x2[n]
add = x1 + x2;
subplot(3,2,1);
stem(n, add, 'filled');
title('a) Signal Addition: x1[n] + x2[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
% b) Signal Multiplication: y[n] = x1[n] * x2[n]
mult = x1 .* x2;
subplot(3,2,2);
stem(n, mult, 'filled');
title('b) Signal Multiplication: x1[n] * x2[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
% c) Scaling: y[n] = 2 * x1[n]
scaled = 2 * x1;
subplot(3,2,3);
stem(n, scaled, 'filled');
title('c) Signal Scaling: 2 * x1[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
% d) Shifting: y[n] = x1[n - 2] (Right shift by 2)
shift = sin(0.2 * pi * (n - 2));
subplot(3,2,4);
stem(n, shift, 'filled');
title('d) Right Shifting: x1[n - 2]');
xlabel('n'); ylabel('Amplitude');
grid on;
% e) Folding: y[n] = x1[-n]
folded = sin(0.2 * pi * (-n));
subplot(3,2,5);
stem(n, folded, 'filled');
title('e) Signal Folding: x1[-n]');
xlabel('n'); ylabel('Amplitude');
grid on;
Output
3) Consider Program to perform convolution of two given sequences (without using built-in
function) and display the signals.
clc;
clear all;
close all;
% Input sequences
x = [1 2 3]; % First signal x[n]
h = [4 5 6]; % Second signal h[n]
% Lengths of sequences
Lx = length(x);
Lh = length(h);
% Total length of convolution result
Ly = Lx + Lh - 1;
% Zero-padding the sequences
x = [x, zeros(1, Ly - Lx)];
h = [h, zeros(1, Ly - Lh)];
% Manual convolution
y = zeros(1, Ly);
for n = 1:Ly
for k = 1:n
y(n) = y(n) + x(k) * h(n - k + 1);
end
end
% Time indices
nx = 0:(length(x)-1);
nh = 0:(length(h)-1);
ny = 0:(length(y)-1);
% Plot input and output signals
subplot(3,1,1);
stem(nx, x, 'filled');
title('Input Sequence x[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
subplot(3,1,2);
stem(nh, h, 'filled');
title('Impulse Response h[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
subplot(3,1,3);
stem(ny, y, 'filled');
title('Convolved Output y[n] = x[n] * h[n]');
xlabel('n'); ylabel('Amplitude');
grid on;
Output
4)
plot. b) Plot |H(ejω)| and ∠ H(ejω) c) Determine the impulse response h(n).
Consider a causal system y(n) = 0.9y(n-1)+x(n). a) Determine H(z) and sketch its pole zero
clc;
clear all;
close all;
% System: y[n] = 0.9*y[n-1] + x[n]
% => H(z) = 1 / (1 - 0.9*z^(-1))
b = [1]; % Numerator
a = [1 -0.9]; % Denominator
% a) Pole-Zero Plot using roots
z = roots(b); % Zeros
p = roots(a); % Poles
figure;
plot(real(z), imag(z), 'ob', 'MarkerSize', 10, 'LineWidth', 2); hold on;
plot(real(p), imag(p), 'xr', 'MarkerSize', 10, 'LineWidth', 2);
theta = linspace(0, 2*pi, 200);
plot(cos(theta), sin(theta), 'k--'); % Unit circle
xlabel('Real Part'); ylabel('Imaginary Part');
title('a) Pole-Zero Plot of H(z)');
legend('Zeros','Poles','Unit Circle');
axis equal; grid on;
% b) Frequency Response H(e^jw)
N = 512;
w = linspace(0, pi, N);
H = freqz(b, a, w); % Frequency response
figure;
subplot(2,1,1);
plot(w/pi, abs(H), 'b', 'LineWidth', 2);
title('b) Magnitude Response |H(e^{j\omega})|');
xlabel('\omega/\pi'); ylabel('Magnitude');
grid on;
subplot(2,1,2);
plot(w/pi, angle(H), 'r', 'LineWidth', 2);
title('Phase Response ∠H(e^{j\omega})');
xlabel('\omega/\pi'); ylabel('Phase (radians)');
grid on;
% c) Impulse Response h[n]
n = 0:30;
x = [1, zeros(1, length(n)-1)]; % Unit impulse
h = filter(b, a, x); % System response
figure;
stem(n, h, 'filled');
title('c) Impulse Response h[n]');
xlabel('n'); ylabel('h[n]');
grid on;
Output
5) Computation of N point DFT of a given sequence (without using built-in function) and to
plot the magnitude and phase spectrum.
clc;
clear all;
close all;
% Input sequence
x = [1 2 3 4]; % You can change this
N = length(x); % N-point DFT
% Manual DFT computation
X = zeros(1, N);
for k = 0:N-1
for n = 0:N-1
X(k+1) = X(k+1) + x(n+1) * exp(-j * 2 * pi * k * n / N);
end
end
% Frequency index
k = 0:N-1;
% Plot Magnitude Spectrum
subplot(2,1,1);
stem(k, abs(X), 'filled');
title('Magnitude Spectrum |X[k]|');
xlabel('Frequency Index k'); ylabel('|X[k]|');
grid on;
% Plot Phase Spectrum
subplot(2,1,2);
title('Phase Spectrum ∠X[k]');
stem(k, angle(X), 'filled');
xlabel('Frequency Index k'); ylabel('Phase (radians)');
grid on;
Output
6) Using the DFT and IDFT, compute the following for any two given sequences a) Circular
convolution b) Linear convolution.
clc;
clear all;
close all;
% Input sequences
x1 = [1 2 3]; % First sequence
x2 = [4 5 6]; % Second sequence
% ---- a) Circular Convolution ----
% Use length N = max(length(x1), length(x2)) or you can use a fixed N
N = max(length(x1), length(x2));
x1_circ = [x1, zeros(1, N - length(x1))];
x2_circ = [x2, zeros(1, N - length(x2))];
% DFT of both sequences
X1 = zeros(1,N);
X2 = zeros(1,N);
for k = 0:N-1
for n = 0:N-1
X1(k+1) = X1(k+1) + x1_circ(n+1) * exp(-j*2*pi*k*n/N);
X2(k+1) = X2(k+1) + x2_circ(n+1) * exp(-j*2*pi*k*n/N);
end
end
% Point-wise multiplication in frequency domain
Y_circ = X1 .* X2;
% Inverse DFT (manual)
y_circ = zeros(1, N);
for n = 0:N-1
for k = 0:N-1
y_circ(n+1) = y_circ(n+1) + Y_circ(k+1) * exp(j*2*pi*k*n/N);
end
y_circ(n+1) = y_circ(n+1) / N; % Normalization
end
% ---- b) Linear Convolution ----
L = length(x1) + length(x2) - 1; % Linear convolution length
x1_lin = [x1, zeros(1, L - length(x1))];
x2_lin = [x2, zeros(1, L - length(x2))];
% DFT of zero-padded signals
X1L = zeros(1, L);
X2L = zeros(1, L);
for k = 0:L-1
for n = 0:L-1
X1L(k+1) = X1L(k+1) + x1_lin(n+1) * exp(-j*2*pi*k*n/L);
X2L(k+1) = X2L(k+1) + x2_lin(n+1) * exp(-j*2*pi*k*n/L);
end
end
% Frequency-domain multiplication
Y_lin = X1L .* X2L;
% IDFT for linear convolution
y_lin = zeros(1, L);
for n = 0:L-1
for k = 0:L-1
y_lin(n+1) = y_lin(n+1) + Y_lin(k+1) * exp(j*2*pi*k*n/L);
end
y_lin(n+1) = y_lin(n+1) / L;
end
% ---- Plot Results ----
n1 = 0:N-1;
n2 = 0:L-1;
figure;
subplot(2,1,1);
stem(n1, real(y_circ), 'filled');
title('a) Circular Convolution via DFT');
xlabel('n'); ylabel('Amplitude');
grid on;
subplot(2,1,2);
stem(n2, real(y_lin), 'filled');
title('b) Linear Convolution via DFT');
xlabel('n'); ylabel('Amplitude');
grid on;
Output
7) Verification of Linearity property, circular time shift property & circular frequency shift
property of DFT.
clc;
clear all;
close all;
% Define length and input sequences
N = 8;
n = 0:N-1;
% Define two sequences x1[n] and x2[n]
x1 = [1 2 3 4 0 0 0 0];
x2 = [4 3 2 1 0 0 0 0];
a = 2; b = 3; % Scalars for linearity
% ===== 1. Linearity Property: DFT{a*x1 + b*x2} = a*DFT{x1} + b*DFT{x2} =====
x_comb = a*x1 + b*x2;
% Manual DFT function
function X = myDFT(x, N)
X = zeros(1,N);
for k = 0:N-1
for n = 0:N-1
X(k+1) += x(n+1) * exp(-j*2*pi*k*n/N);
end
end
end
X1 = myDFT(x1, N);
X2 = myDFT(x2, N);
X_comb = myDFT(x_comb, N);
X_verify = a*X1 + b*X2;
fprintf('Linearity Property Verified: %d\n', all(abs(X_comb - X_verify) < 1e-10));
% ===== 2. Circular Time-Shift Property =====
k_shift = 2; % Circular shift by 2 samples
x_shifted = circshift(x1, [0 k_shift]); % Circular shift to the right
X_shifted = myDFT(x_shifted, N);
X1_phase = X1 .* exp(-j*2*pi*(0:N-1)*k_shift/N);
fprintf('Circular Time Shift Property Verified: %d\n', all(abs(X_shifted - X1_phase) < 1e-10));
% ===== 3. Circular Frequency-Shift Property =====
k0 = 2;
% Frequency shift: x[n] * exp(j*2*pi*k0*n/N)
x_freq_shifted = x1 .* exp(j*2*pi*k0*n/N);
X_freq_shifted = myDFT(x_freq_shifted, N);
% Expected result: Circular shift of X1 by k0
X1_shifted = circshift(X1, [0, k0]);
fprintf('Circular Frequency Shift Property Verified: %d\n', all(abs(X_freq_shifted - X1_shifted) < 1e-
10));
% ===== Plot for Visualization =====
figure;
subplot(3,1,1);
stem(0:N-1, abs(X_comb), 'b', 'filled'); hold on;
stem(0:N-1, abs(X_verify), 'r--');
title('1. Linearity Property |DFT|');
legend('|DFT{a*x1+b*x2}|', '|a*DFT{x1}+b*DFT{x2}|');
grid on;
subplot(3,1,2);
stem(0:N-1, abs(X_shifted), 'b', 'filled'); hold on;
stem(0:N-1, abs(X1_phase), 'r--');
title('2. Circular Time Shift |DFT|');
legend('DFT{x[n-k]}', 'DFT{x[n]} * exp(-j2pikn/N)');
grid on;
subplot(3,1,3);
stem(0:N-1, abs(X_freq_shifted), 'b', 'filled'); hold on;
stem(0:N-1, abs(X1_shifted), 'r--');
title('3. Circular Frequency Shift |DFT|');
legend('DFT{x[n]*e^{j2\pi k_0 n/N}}', 'Circular Shift of DFT{x[n]}');
grid on;
Output
Linearity Property Verified: 1
Circular Time Shift Property Verified: 1
Circular Frequency Shift Property Verified: 1
8) Develop decimation in time radix-2 FFT algorithm without using built-in functions.
clc;
clear all;
close all;
% --- Recursive Radix-2 DIT FFT Function ---
function X = myFFT(x)
N = length(x);
if N == 1
X = x;
else
% Divide: even and odd indexed elements
x_even = x(1:2:end);
x_odd = x(2:2:end);
% Conquer: recursively apply FFT
X_even = myFFT(x_even);
X_odd = myFFT(x_odd);
% Combine
X = zeros(1, N);
for k = 0:(N/2 - 1)
W = exp(-j*2*pi*k/N);
X(k+1) = X_even(k+1) + W * X_odd(k+1);
X(k+1 + N/2) = X_even(k+1) - W * X_odd(k+1);
end
end
end
% --- Input Sequence ---
x = [1 2 3 4 0 0 0 0]; % 8-point input sequence (must be power of 2)
N = length(x);
% --- Compute FFT ---
X = myFFT(x);
% --- Plot Results ---
n = 0:N-1;
figure;
subplot(2,1,1);
stem(n, abs(X), 'filled');
title('Magnitude Spectrum |X[k]|');
xlabel('k'); ylabel('|X[k]|');
grid on;
subplot(2,1,2);
stem(n, angle(X), 'filled');
title('Phase Spectrum ∠X[k]');
xlabel('k'); ylabel('Phase (radians)');
grid on;
Output
9) Design and implementation of digital low pass FIR filter using a window to meet the given
specifications.
Problem Statement (Example Specification):
Design a Low Pass FIR Filter with:
Sampling frequency fs=2000f_s = 2000fs=2000 Hz
Cutoff frequency fc=500f_c = 500fc=500 Hz
Filter order N=40N = 40N=40
Use Hamming Window
clc;
clear;
close all;
% Filter specifications
fs = 2000; % Sampling frequency in Hz
fc = 500; % Cutoff frequency in Hz
N = 40; % Filter order
wc = 2 * pi * fc / fs; % Normalized cutoff frequency (radians)
% Time vector for N+1 coefficients (FIR filter has N+1 taps)
n = 0:N;
% Ideal impulse response of low-pass filter (sinc function)
hd = (sin(wc*(n - N/2))) ./ (pi*(n - N/2));
hd(N/2+1) = wc/pi; % Handle NaN at center (divide by zero)
% Apply Hamming window
w = hamming(N+1)'; % Transpose to match dimensions
h = hd .* w;
% Plot impulse response
figure;
stem(n, h, 'filled');
title('Impulse Response of FIR Low Pass Filter');
xlabel('n'); ylabel('h[n]');
grid on;
% Frequency response
[H, f] = freqz(h, 1, 1024, fs); % Frequency response
% Plot magnitude response
figure;
plot(f, abs(H), 'b', 'LineWidth', 2);
title('Magnitude Response of FIR Low Pass Filter');
xlabel('Frequency (Hz)');
ylabel('|H(f)|');
grid on;
% Plot phase response
figure;
plot(f, angle(H), 'r', 'LineWidth', 2);
title('Phase Response of FIR Low Pass Filter');
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
grid on;
% Signal test (Optional)
% Generate a test signal with 2 frequencies: one below and one above cutoff
t = 0:1/fs:1;
x = sin(2*pi*200*t) + sin(2*pi*800*t); % 200 Hz (pass), 800 Hz (stop)
y = filter(h, 1, x);
% Plot input and filtered signals
figure;
subplot(2,1,1);
plot(t, x); title('Input Signal'); xlabel('Time (s)'); ylabel('Amplitude');
subplot(2,1,2);
plot(t, y); title('Filtered Signal (Low-Passed)'); xlabel('Time (s)'); ylabel('Amplitude');
Output
10) Design and implementation of digital high pass FIR filter using a window to meet the given
specifications.
Example Specifications:
Sampling frequency fs=2000f_s = 2000fs=2000 Hz
Cutoff frequency fc=500f_c = 500fc=500 Hz
Filter order N=40N = 40N=40
Use Hamming Window
clc;
clear;
close all;
% Filter specifications
fs = 2000; % Sampling frequency in Hz
fc = 500; % Cutoff frequency in Hz
N = 40; % Filter order
wc = 2 * pi * fc / fs; % Normalized cutoff frequency (radians)
% Time vector
n = 0:N;
% Ideal impulse response for low pass filter
hd_lp = (sin(wc*(n - N/2))) ./ (pi*(n - N/2));
hd_lp(N/2+1) = wc/pi; % Handle divide by zero at center
% Convert to high pass by spectral inversion
hd_hp = -hd_lp;
hd_hp(N/2+1) = 1 - wc/pi;
% Apply Hamming window
w = hamming(N+1)'; % Row vector window
h = hd_hp .* w;
% Plot impulse response
figure;
stem(n, h, 'filled');
title('Impulse Response of FIR High Pass Filter');
xlabel('n'); ylabel('h[n]');
grid on;
% Frequency response
[H, f] = freqz(h, 1, 1024, fs);
% Plot magnitude response
figure;
plot(f, abs(H), 'b', 'LineWidth', 2);
title('Magnitude Response of FIR High Pass Filter');
xlabel('Frequency (Hz)');
ylabel('|H(f)|');
grid on;
% Plot phase response
figure;
plot(f, angle(H), 'r', 'LineWidth', 2);
title('Phase Response of FIR High Pass Filter');
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
grid on;
% Test signal: mixture of 200 Hz (stop) and 800 Hz (pass)
t = 0:1/fs:1;
x = sin(2*pi*200*t) + sin(2*pi*800*t);
y = filter(h, 1, x);
% Plot input and filtered signals
figure;
subplot(2,1,1);
plot(t, x); title('Input Signal'); xlabel('Time (s)'); ylabel('Amplitude');
subplot(2,1,2);
plot(t, y); title('Filtered Signal (High-Passed)'); xlabel('Time (s)'); ylabel('Amplitude');
Output
11) Design and implementation of digital IIR Butterworth low pass filter to meet the given
specifications.
Example Specifications:
Sampling frequency fs=2000f_s = 2000fs=2000 Hz
Cutoff frequency fc=500f_c = 500fc=500 Hz
Filter order N=4N = 4N=4
Type: Butterworth Low Pass
clc;
clear;
close all;
pkg load signal; % Load signal processing package
% Filter Specifications
fs = 2000; % Sampling frequency in Hz
fc = 500; % Cutoff frequency in Hz
N = 4; % Filter order
% Normalized cutoff frequency (0 to 1)
Wn = fc / (fs/2);
% Design Butterworth Low Pass Filter
[b, a] = butter(N, Wn, 'low');
% Frequency Response
[H, f] = freqz(b, a, 1024, fs);
% Plot Magnitude Response
figure;
plot(f, abs(H), 'b', 'LineWidth', 2);
title('Butterworth Low Pass Filter - Magnitude Response');
xlabel('Frequency (Hz)');
ylabel('|H(f)|');
grid on;
% Plot Phase Response
figure;
plot(f, angle(H), 'r', 'LineWidth', 2);
title('Butterworth Low Pass Filter - Phase Response');
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
grid on;
% Impulse Response
figure;
impz(b, a);
title('Impulse Response');
xlabel('Sample Index'); ylabel('Amplitude');
grid on;
% Step Response (manual)
step_input = ones(1, 100);
step_response = filter(b, a, step_input);
figure;
plot(step_response, 'LineWidth', 2);
title('Step Response of Butterworth Low Pass Filter');
xlabel('Sample Index');
ylabel('Amplitude');
grid on;
% Test Signal (Combination of 200 Hz and 800 Hz)
t = 0:1/fs:1;
x = sin(2*pi*200*t) + sin(2*pi*800*t); % 200 Hz should pass, 800 Hz should be attenuated
y = filter(b, a, x); % Apply the filter
% Plot Input and Filtered Output
figure;
subplot(2,1,1);
plot(t, x);
title('Input Signal (200 Hz + 800 Hz)');
xlabel('Time (s)');
ylabel('Amplitude');
subplot(2,1,2);
plot(t, y);
title('Filtered Signal (Low Passed)');
xlabel('Time (s)');
ylabel('Amplitude');
Output
12) Design and implementation of digital IIR Butterworth high pass filter to meet the given
specifications.
Example Specifications:
Sampling frequency fs=2000f_s = 2000fs=2000 Hz
Cutoff frequency fc=500f_c = 500fc=500 Hz
Filter order N=4N = 4N=4
Type: Butterworth High Pass
clc;
clear;
close all;
pkg load signal; % Load the Signal Processing Package
% Filter specifications
fs = 2000; % Sampling frequency in Hz
fc = 500; % Cutoff frequency in Hz
N = 4; % Filter order
% Normalize the cutoff frequency (between 0 and 1)
Wn = fc / (fs/2);
% Design Butterworth High Pass Filter
[b, a] = butter(N, Wn, 'high');
% Frequency Response
[H, f] = freqz(b, a, 1024, fs);
% Plot Magnitude Response
figure;
plot(f, abs(H), 'b', 'LineWidth', 2);
title('Butterworth High Pass Filter - Magnitude Response');
xlabel('Frequency (Hz)');
ylabel('|H(f)|');
grid on;
% Plot Phase Response
figure;
plot(f, angle(H), 'r', 'LineWidth', 2);
title('Butterworth High Pass Filter - Phase Response');
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
grid on;
% Impulse Response
figure;
impz(b, a);
title('Impulse Response of High Pass Filter');
xlabel('Sample Index');
ylabel('Amplitude');
grid on;
% Step Response (manual method using filter on step input)
step_input = ones(1, 100);
step_response = filter(b, a, step_input);
figure;
plot(step_response, 'LineWidth', 2);
title('Step Response of Butterworth High Pass Filter');
xlabel('Sample Index');
ylabel('Amplitude');
grid on;
% Test Signal: 200 Hz (should be blocked), 800 Hz (should pass)
t = 0:1/fs:1;
x = sin(2*pi*200*t) + sin(2*pi*800*t);
y = filter(b, a, x);
% Plot Input and Filtered Output
figure;
subplot(2,1,1);
plot(t, x);
title('Input Signal (200 Hz + 800 Hz)');
xlabel('Time (s)');
ylabel('Amplitude');
subplot(2,1,2);
plot(t, y);
title('Filtered Output (High Passed)');
xlabel('Time (s)');
ylabel('Amplitude');
Output
*****END*****