NATIONAL INSTITUTE OF TECHNOLOGY CALICUT
Department of Electronics and Communication Engineering
Monsoon 2024
EC6405E Statistical Signal Processing
Assignment - 2
SUBMITTED BY
SYAMANTAK SARKAR (M230922EC)
M.Tech – Signal Processing
Part I and ii
Theory:
Code:
clc;
clear;
close all;
%Generation of random signal
N = 100;
v = randn(1, N);
x = zeros(1, N);
for n=1:N
if n==1
x(n) = 0.447 *[v(3)+v(2)+v(1)];
elseif n==2
x(n)= 0.447 *[v(4)+v(3)+v(2)+v(1)];
elseif n== 3:N
x(n) = 0.447 *[v(n+2)+v(n+1)+v(n)+v(n-1)+v(n-2)] ;
elseif n==N-1
x(n)= 0.447 *[v(n+1)+v(n)+v(n-1)+v(n-2)];
else
x(n)= 0.447 *[v(n)+v(n-1)+v(n-2)];
end
end
figure;
subplot (2,1,1);
stem(1:N, v, 'LineWidth', 2);
xlabel('n');
ylabel('v(n)');
title(' random noise v(n)');
grid on;
subplot (2,1,2);
stem(1:N, x, 'LineWidth', 2);
xlabel('n');
ylabel('x(n)');
title('Generated random signal x(n)');
% Prediction using winer 2 tap filter
x_hat = zeros(N);
A = [1 0.8; 0.8 1];
B = [0.8 ;0.6];
W = inv(A )*B;
w0 = W(1);
w1 = W(2);
for i = 1:1:N
if i==1
x_hat(i) = x(i);
elseif i==2
x_hat(i) = w0*x(i-1);
else
x_hat(i) = w0*x(i-1) + w1*x(i-2);
end
end
figure;
subplot(3,1,1);
stem(1:N,x);
title('Generated Signal');
subplot(3,1,2);
stem(1:N,x_hat);
title('x_h (Estimated Signal)');
e = zeros(N);
for i = 1:1:N
e(i) = x(i) - x_hat(i);
end
subplot(3,1,3);
stem(1:N,e);
title('error signal');
% PCM quantization error calculation(4 bit case)
Mx = max(x_hat,[],"all");
del_2 = Mx/2.^3;
quan_dpcm = (del_2.*del_2) /12;
SNR_dpcm = 1/quan_dpcm;
% DPCM quantization error calculation (4 bit case)
Me = max(e,[],"all");
del_1 = Me/2.^3;
quan_pcm = (del_1.*del_1) /12;
SNR_pcm = 1/quan_pcm;
fprintf('Signal to Noise Ratio for PCM\n');
disp(SNR_pcm);
fprintf('Signal to Noise Ratio for DPCM\n');
disp(SNR_dpcm);
Results:
Signal to Noise Ratio for PCM
397.2769
Signal to Noise Ratio for DPCM
513.441
Inference:
• The autocorrelation function taken here is 1-0.2|k|
• As we can see that the SNR for DPCM is much higher than SNR for PCM.
• As there is no noise while estimation the signal so the SNR value on an average is higher
compared to the case where there is noise present.
Part iii
Theory:
Code:
% Define k values
k = 0:10;
N=2;
% Calculate autocorrelation function
r_x = autocorr_func(k);
% Define the matrix W
Rx_mat = autocorrelationMatrix(N);
% Define the vector v
rx_vect = autocorrelationvector(N);
% Perform matrix multiplication
s=inv(Rx_mat);
y=transpose(rx_vect);
w=s*y;
% Display the result
disp('w matrix(L=2)');
disp(w);
%error power
% L=3 tap calculations
P_2n = r_x(1)- [w(1,1)*r_x(2)] - [w(2,1)*r_x(3)];
r_3=[conj(r_x(4))] -[w(1,1)*conj(r_x(3))] -[w(2,1)*conj(r_x(2))];
T_3 =(-r_3/P_2n);
v=T_3;
w_21 = w(1,1)+(T_3*w(2,1));
z=w(2,1);
c=w(1,1);
v=T_3;
w_22=w(2,1)+(T_3*w(1,1));
w_23 = T_3;
P_3 = [1-((T_3)*(T_3))]*P_2n;
P_3n = r_x(1)- [w_21*r_x(2)] - [w_22*r_x(3)]-[w_23*r_x(4)];
fprintf('P2 = ');
disp(P_2n);
fprintf('gamma_3 = ');
disp(r_3);
fprintf('T_3 = ');
disp(v);
disp('w matrix(L=3)');
disp(w_21);
disp(w_22);
disp(w_23);
% L=4 tap calculations
P_3n = r_x(1)- [w_21*r_x(2)] - [w_22*r_x(3)]-[w_23*r_x(4)];
P_3 = [1-((T_3)*(T_3))]*P_2n
r_4=[conj(r_x(5))] -[w_21*conj(r_x(4))] -[w_22*conj(r_x(3))] -
[w_23*conj(r_x(2))];
T_4 =(-r_4/P_3);
w_31 = w_21+(T_4*w_23);
w_32=w_22+(T_4*w_22);
w_33 = w_23+(T_4*w_21);
w_34 = T_4;
fprintf('P3 = ');
disp(P_3);
fprintf('gamma_4 = ');
disp(r_4);
fprintf('T_4 = ');
disp(T_4);
disp('w matrix(L=4)');
disp(w_31);
disp(w_32);
disp(w_33);
disp(w_34);
% L=5 tap calculations
P_4n = r_x(1)- [w_31*r_x(2)] - [w_32*r_x(3)]-[w_33*r_x(4)] -[w_34*r_x(5)];
P_4 = [1-((T_4)*(T_4))]*P_3;
r_5=[conj(r_x(6))] -[w_34*conj(r_x(5))] -[w_33*conj(r_x(4))] -
[w_32*conj(r_x(3))] - [w_31*conj(r_x(3))];
T_5 =(-r_5/P_4);
w_41 = w_31+(T_5*w_34);
w_42=w_32+(T_5*w_33);
w_43 = w_33+(T_5*w_32);
w_44 = w_34+(T_5*w_31);
w_45 = T_5;
fprintf('P4 = ');
disp(P_4);
fprintf('gamma_5 = ');
disp(r_5);
fprintf('T_5 = ');
disp(T_5);
disp('w matrix(L=5)');
disp(w_41);
disp(w_42);
disp(w_43);
disp(w_44);
disp(w_45);
% Prediction using winer 5 tap filter
N = 100;
v = randn(1, N);
x = zeros(1, N);
for n=1:N
if n==1
x(n) = 0.447 *[v(3)+v(2)+v(1)];
elseif n==2
x(n)= 0.447 *[v(4)+v(3)+v(2)+v(1)];
elseif n== 3:N
x(n) = 0.447 *[v(n+2)+v(n+1)+v(n)+v(n-1)+v(n-2)] ;
elseif n==N-1
x(n)= 0.447 *[v(n+1)+v(n)+v(n-1)+v(n-2)];
else
x(n)= 0.447 *[v(n)+v(n-1)+v(n-2)];
end
end
x_hat = zeros(N);
w0 = w_41;
w1 = w_42;
w2 = w_43;
w3 = w_44;
w4 = w_45;
for i = 1:1:N
if i==1
x_hat(i) = x(i);
elseif i==2
x_hat(i) = w0*x(i-1);
elseif i==3
x_hat(i) = w0*x(i-1) + w1*x(i-2);
elseif i==4
x_hat(i) = w0*x(i-1) + w1*x(i-2)+ w2*x(i-3);
elseif i==5
x_hat(i) = w0*x(i-1) + w1*x(i-2) + w2*x(i-3) + w3*x(i-4);
else
x_hat(i) = w0*x(i-1) + w1*x(i-2) + w2*x(i-3) + w3*x(i-4) + w4*x(i-
5);
end
end
figure;
subplot(3,1,1);
stem(1:N,x);
title('Generated Signal');
subplot(3,1,2);
stem(1:N,x_hat);
title('x_h (Predicted Signal)');
e = zeros(N);
for i = 1:1:N
e(i) = x(i) - x_hat(i);
end
subplot(3,1,3);
stem(1:N,e);
title('error signal');
% PCM quantization error calculation(4 bit case)
Mx_max = max(x_hat,[],"all");
Mx_min = min(x_hat,[],"all");
del_2 = (Mx_max - Mx_min)/2.^4;
quan_dpcm = (del_2.*del_2) /12;
SNR_dpcm = 1/quan_dpcm;
% DPCM quantization error calculation (4 bit case)
Me_max = max(e,[],"all");
Me_min = min(e,[],"all");
del_1 = (Me_max - Me_min)/2.^4;
quan_pcm = (del_1.*del_1) /12;
SNR_pcm = 1/quan_pcm;
fprintf('Signal to Noise Ratio for PCM\n');
disp(SNR_pcm);
fprintf('Signal to Noise Ratio for DPCM\n');
disp(SNR_dpcm);
function Rx = autocorrelationMatrix(N)
Rx = zeros(N);
for i = 1:N
for j = i:N
Rx(i, j) = autocorr_func(i - j);
Rx(j, i) = Rx(i, j);
end
end
end
%rx vector creation
function rx = autocorrelationvector(N)
rx = zeros(1:N);
for i = 1:N
rx(i) = autocorr_func(i);
end
end
%function wie = wiener_func(i,j,T_3)
%wie = w(i,i) + T_3 * conj(w(i,j));
function r_x = autocorr_func(k)
r_x = 1 - 0.2*abs(k);
end
Results:
>> exp2_2
w matrix(L=2)
0.8889
-0.1111
P2 = 0.3556
gamma_3 = -0.0444
T_3 = 0.1250
w matrix(L=3)
0.8750
-1.0686e-15
0.1250
P_3 =
0.3500
P3 = 0.3500
gamma_4 = -0.2500
T_4 = 0.7143
w matrix(L=4)
0.9643
-1.8319e-15
0.7500
0.7143
P4 = 0.1714
gamma_5 = -1.0214
T_5 = 5.9583
w matrix(L=5)
5.2202
4.4687
0.7500
6.4598
5.9583
Signal to Noise Ratio for PCM
433.95
Signal to Noise Ratio for DPCM
550.23
Inference:
• As we can see the coefficients for the 5 tap weiner estimer is much higher.
• As we are increasing the no of taps from L=2 to L=5 at each of the step the filter coefficeents are
changing. So there is need for recalculation for all the filter coefficients at each step.
• The SNR values are higher than L = 2 case in both PCM and DPCM. This is due to the fact that it is
a memory less correlation function and the the prediction will be more and more accurate as
the number of tape(i.e. the dependencies on the previous inputs) are increased.
Part iv
Theory:
Code:
For White Noise
%White Noise
clc;
clear;
close all;
%Generation of random signal
N = 100;
v = randn(1, N);
x = zeros(1, N);
y = zeros(1, N);
for n=1:N
if n==1
x(n) = 0.447 *[v(3)+v(2)+v(1)];
elseif n==2
x(n)= 0.447 *[v(4)+v(3)+v(2)+v(1)];
elseif n== 3:N
x(n) = 0.447 *[v(n+2)+v(n+1)+v(n)+v(n-1)+v(n-2)] ;
elseif n==N-1
x(n)= 0.447 *[v(n+1)+v(n)+v(n-1)+v(n-2)];
else
x(n)= 0.447 *[v(n)+v(n-1)+v(n-2)];
end
end
figure;
subplot (2,1,1);
stem(1:N, x, 'LineWidth', 2);
xlabel('n');
ylabel('x(n)');
title('Generated random signal x(n)')
% white Noise
u = normrnd(0,sqrt(0.25),N);
for i = 1:1:N
y(i) = x(i)+ u(i);
end
subplot(2,1,2);
stem(1:N,y);
title('Y (x+white noise)');
% Prediction using winer 2 tap filter
x_hat = zeros(N);
A = [1.25 0.8; 0.8 1.25];
B = [0.8 ;0.6];
W = inv(A )*B;
w0 = W(1);
w1 = W(2);
for i = 1:1:N
if i==1
x_hat(i) = y(i);
elseif i==2
x_hat(i) = w0*y(i-1);
else
x_hat(i) = w0*y(i-1) + w1*y(i-2);
end
end
figure;
subplot(3,1,1);
stem(1:N,x);
title('Generated Signal');
subplot(3,1,2);
stem(1:N,x_hat);
title('x_h (predicted signal)');
e = zeros(N);
for i = 1:1:N
e(i) = x(i) - x_hat(i);
end
subplot(3,1,3);
stem(1:N,e);
title('error signal');
% PCM quantization error calculation(4 bit case)
Mx = max(x_hat,[],"all");
del_2 = Mx/2.^3;
quan_dpcm = (del_2.*del_2) /12;
SNR_dpcm = 1/quan_dpcm;
% DPCM quantization error calculation (4 bit case)
Me = max(e,[],"all");
del_1 = Me/2.^3;
quan_pcm = (del_1.*del_1) /12;
SNR_pcm = 1/quan_pcm;
fprintf('Signal to Noise Ratio for PCM\n');
disp(SNR_pcm);
fprintf('Signal to Noise Ratio for DPCM\n');
disp(SNR_dpcm);
For Colored Noise
%coloured Noise
clc;
clear;
close all;
%Generation of random signal
N = 100;
v = randn(1, N);
x = zeros(1, N);
y = zeros(1, N);
for n=1:N
if n==1
x(n) = 0.447 *[v(3)+v(2)+v(1)];
elseif n==2
x(n)= 0.447 *[v(4)+v(3)+v(2)+v(1)];
elseif n== 3:N
x(n) = 0.447 *[v(n+2)+v(n+1)+v(n)+v(n-1)+v(n-2)] ;
elseif n==N-1
x(n)= 0.447 *[v(n+1)+v(n)+v(n-1)+v(n-2)];
else
x(n)= 0.447 *[v(n)+v(n-1)+v(n-2)];
end
end
figure;
subplot (2,1,1);
stem(1:N, x, 'LineWidth', 2);
xlabel('n');
ylabel('x(n)');
title('Generated random signal x(n)')
u = zeros(N);% coloured Noise
u(1) = sqrt(0.25/2)*v(1);
for i = 2:1:N
u(i) = sqrt(0.25/2)*v(i-1) + sqrt(0.25/2)*v(i);
end
for i = 1:1:N
y(i) = x(i)+ u(i);
end
subplot(2,1,2);
stem(1:N,y);
title('Y (x+coloured noise)');
% Prediction using winer 2 tap filter
x_hat = zeros(N);
A = [1.25 0.925; 0.8 0.925];
B = [0.8 ;0.6];
W = inv(A )*B;
w0 = W(1);
w1 = W(2);
for i = 1:1:N
if i==1
x_hat(i) = y(i);
elseif i==2
x_hat(i) = w0*y(i-1);
else
x_hat(i) = w0*y(i-1) + w1*y(i-2);
end
end
figure;
subplot(3,1,1);
stem(1:N,x);
title('Generated Signal');
subplot(3,1,2);
stem(1:N,x_hat);
title('x_h (predicted signal)');
e = zeros(N);
for i = 1:1:N
e(i) = x(i) - x_hat(i);
end
subplot(3,1,3);
stem(1:N,e);
title('error signal');
% PCM quantization error calculation(4 bit case)
Mx = max(x_hat,[],"all");
del_2 = Mx/2.^3;
quan_dpcm = (del_2.*del_2) /12;
SNR_dpcm = 1/quan_dpcm;
% DPCM quantization error calculation (4 bit case)
Me = max(e,[],"all");
del_1 = Me/2.^3;
quan_pcm = (del_1.*del_1) /12;
SNR_pcm = 1/quan_pcm;
fprintf('Signal to Noise Ratio for PCM\n');
disp(SNR_pcm);
fprintf('Signal to Noise Ratio for DPCM\n');
disp(SNR_dpcm);
Results:
For White Noise
Signal to Noise Ratio for PCM
230.5587
Signal to Noise Ratio for DPCM
338.3980
For Colored Noise
Signal to Noise Ratio for PCM
196.3035
Signal to Noise Ratio for DPCM
249.7851
Inference:
• The SNR value is less in case when we are predicting from noisey signals.
• The SNR value is lesser when the signal is coloured one when compared to the white noise one.
• As noise is introduced the SNR value for prediction also decreases.
A comparison table is given below for different DPCM and PCM SNR values.
Prediction using L = Prediction using L = Prediction using L = Prediction using L =
2 Weiner filter 5 Weiner filter 2 Weiner filter 2 Weiner filter
(No noise) (No noise) (White noise) (Colored noise)
DPCM 513.441 550.23 338.3980 249.7851
PCM 397.2769 433.95 230.5587 196.3035