%ECE 480 Lab #3

%Problem #1 a.
n = [0:50]';
x = ones(size(n));
l = 5;
ind = find(n >= (l-1));
x(ind) = 0;
y = x .* (-1).^n;
subplot(2,1,1);
stem(n,x), xlabel('Samples'), ylabel('Amplitude'), title('Problem #1 Input');
subplot(2,1,2);
stem(n,y), xlabel('Samples'), ylabel('Amplitude'), title('Problem #1 Output');
figure;

%Problem #1 b.
w = [0:pi/1024;pi]';
Xspect = freqz(x,1,w);
magX = abs(Xspect);
Yspect = freqz(y,1,w);
magY = abs(Yspect);
subplot(2,1,1);
plot(w,magX), xlabel('Omega'), ylabel('Amplitude'), title('Problem #1 Magnitude Spectrum of x(n)');
grid on;
axis([0,pi,0,4]);
subplot(2,1,2);
plot(w, magY), xlabel('Omega'), ylabel('Amplitude'), title('Problem #1 Magnitude Spectrum of y(n)');
grid on;
axis([0,pi,0,4]);
figure;
phaseX = angle(Xspect);
phaseY = angle(Yspect);
subplot(2,1,1);
plot(w, phaseX), xlabel('Omega'), ylabel('Radians'), title('Problem #1 Phase Spectrum of x(n)');
grid on;
subplot(2,1,2);
plot(w, phaseY), xlabel('Omega'), ylabel('Radians'), title('Problem #1 Phase Spectrum of y(n)');
grid on;
figure;
%Problem #1 d.
w = [0:pi/1024;pi]';
g = conv(y,y);
Gspect = freqz(g,1,w);
magG = abs(Gspect);
subplot(2,1,1);
plot(w,magG), xlabel('Omega'), ylabel('Amplitude'), title('Problem #1 Magnitude Spectrum of g(n)');
grid on;
phaseG = angle(Gspect);
subplot(2,1,2);
plot(w, phaseG), xlabel('Omega'), ylabel('Radians'), title('Problem #1 Phase Spectrum of g(n)');
grid on;
figure;

%Problem #2a.
n = [0:50]';
x1 = ones(size(n));
l = 8;
ind = find(n >= (l));
x1(ind) = 0;
stem(n,x1), xlabel('Samples'), ylabel('Amplitude'), title('Problem #2 Inpulse Response Graph');
figure;
%Problem #2 c.
n = [0:50]';
x1 = ones(size(n));
l = 3;
ind = find(n >= (l-1));
x1(ind) = 0;
a = [1,0,0,0,0,0];
b = [1,1,1,1,1,1];
y1 = filter(b, a, x1);
subplot(2,1,1);
stem(n,x1), xlabel('Samples'), ylabel('Amplitude'), title('Problem #2 Input');
subplot(2,1,2);
stem(n,y1), xlabel('Samples'), ylabel('Amplitude'), title('Problem #2 Output');
figure;
yen = (y1 .^2);
sum(yen)

%Problem #2 d.
w = [0:pi/1024;pi]';
Xspect = freqz(x1,1,w);
magX = abs(Xspect);
Yspect = freqz(y1,1,w);
magY = abs(Yspect);
subplot(2,1,1);
plot(w,magX), xlabel('Omega'), ylabel('Amplitude'), title('Problem #2 Magnitude Spectrum of x(n)');
grid on;
axis([0,pi,0,2]);
subplot(2,1,2);
plot(w, magY), xlabel('Omega'), ylabel('Amplitude'), title('Problem #2 Magnitude Spectrum of y(n)');
grid on;
axis([0,pi,0,12]);
figure;
phaseX = angle(Xspect);
phaseY = angle(Yspect);
subplot(2,1,1);
plot(w, phaseX), xlabel('Omega'), ylabel('Radians'), title('Problem #2 Phase Spectrum of x(n)');
grid on;
axis([0,pi,-1.572,0]);
subplot(2,1,2);
plot(w, phaseY), xlabel('Omega'), ylabel('Radians'), title('Problem #2 Phase Spectrum of y(n)');
grid on;
axis([0,pi,-.5,.5]);
figure;
%group delay
w = [0:pi/1024:pi]';
grd = grpdelay(b,a,w);
plot (w, grd), xlabel('Omega'), ylabel('Samples'), title('Problem #2 Group Delay of Filter');
grid on;
axis([0,pi,-1,3]);
figure;

%Problem #3 a.
a = -0.6;
w = [0:pi/128:pi]';
x = 1./(1-a*exp(-j*w));
magX = abs(x);
phaseX = angle(x);
subplot(2,1,1);
plot(w, magX), xlabel('Omega'), ylabel('Amplitude'), title('Problem #3 Magnitude Response');
axis([0,pi,.5,2.5]);
grid on;
subplot(2,1,2);
plot(w, phaseX), xlabel('Omega'), ylabel('Radians'), title('Problem #3 Phase Spectrum of x(n)');
grid on;
axis([0,pi,0,.8]);
figure;

a = [1, .6];
b = [1];
x = freqz(b, a, w);
magX = abs(x);
phaseX = angle(x);
subplot(2,1,1);
plot(w, magX), xlabel('Omega'), ylabel('Amplitude'), title('Problem #3 Magnitude Response');
grid on;
axis([0,pi,.5,2.5]);
subplot(2,1,2);
plot(w, phaseX), xlabel('Omega'), ylabel('Radians'), title('Problem #3 Phase Spectrum of x(n)');
grid on;
axis([0,pi,0,.8]);
figure;

%Problem #4
w = [0:pi/128:pi]';
b1 = [1, -1*exp(-j*-.25*pi)];
b2 = [1, -1*exp(-j*.25*pi)];
b = conv(b1, b2);
a1 = [1, -.8*exp(-j*-.25*pi)];
a2 = [1, -.8*exp(-j*.25*pi)];
a = conv(a1, a2);
xspect = freqz(b,a,w);
magx = abs(xspect);
plot(w,magx), xlabel('Omega'), ylabel('Amplitude'), title('Problem #4');
grid on;
axis([0,pi,0,1.3]);
figure;

%Problem #5
n1 = [.287, -.447];
d1 = [1, -1.297, .695];
n2 = [-2.142, 1.145];
d2 = [1, -1.069, .37];
d = conv(d1, d2);
n = conv(n1, d2) + conv(d1, n2);
n3 = [1.856, -.63];
d3 = [1, -.997, .257];
dfin = conv(d, d3);
nfin = conv(n, d3) + conv(d, n3);
w = [0:pi/128:pi]';
xspect = freqz(nfin,dfin,w);
plot(w,xspect), xlabel('Omega'), ylabel('Amplitude'), title('Problem #5, Frequency Response');
%group delay
w = [0:pi/1024:pi]';
grd = grpdelay(nfin,dfin,w);
plot (w, grd), xlabel('Omega'), ylabel('Samples'), title('Problem #5 Group Delay of Filter');
grid on;
axis([0,pi,-2,10]);
figure;

%Problem #6
n1 = [1,-1];
n2 = conv(n1, n1);
n4 = conv(n2, n2);
n5 = [.8, 1];
n6 = conv(n4, n5);
nfin = n6 * .02426;
d1 = [1, 1.0416, .4019];
d2 = [1, .5661, .7657];
d3 = [1, .8];
d4 = conv(d1, d2);
dfin = conv(d4, d3);
w = [0:pi/128:pi]';
x = freqz(nfin, dfin, w);
magX = abs(x);
phaseX = angle(x);
subplot(2,1,1);
plot(w, magX), xlabel('Omega'), ylabel('Amplitude'), title('Problem #6 Magnitude Response');
grid on;
axis([0,pi,0,1.25]);
subplot(2,1,2);
plot(w, phaseX), xlabel('Omega'), ylabel('Radians'), title('Problem #6 Phase Spectrum of x(n)');
grid on;
axis([0,pi,-4,4]);
figure;
%group delay
w = [0:pi/128:pi]';
grd = grpdelay(nfin,dfin,w);
plot (w, grd), xlabel('Omega'), ylabel('Samples'), title('Problem #6 Group Delay of Filter');
grid on;
figure;
%DB Problem
iend = length(x);
magXdb = zeros(size(x));
for i = 1:iend
    magXdb(i) = 20*log(magX(i));
end
plot(w,magXdb), xlabel('Omega in Radians per Second'), ylabel('Decibels'), title('Magnitude Response in Decibels');

%Problem #7
%n = [0:255];