**********************
**********************
*LABORATORY SESSION 5*
**********************
**********************

%FIR
%task-a
hl = [0.5 0.5];
hl4 = conv(conv(conv(hl,hl),hl),hl); %=[1/16 1/4 3/8 1/4 1/16] 
%verify lowpass filtering of hl4  http://www.fourier-series.com under digital convolution tab

%task-b
hh = [0.5 -0.5];
hh4 = conv(conv(conv(hl,hl),hl),hl); %=[1/16 -1/4 3/8 -1/4 1/16]
%verify highpass filtering effect of hh4  http://www.fourier-series.com under digital convolution tab


%IIR
%task-a
% y(n) = cy(n-1) + x(n)
% Derive the impulse response (only first few values) using y(-1) = 0
y(0) = cy(-1)  + delta(0) = 0 + 1 = 1 
y(1) = cy(0) + delta(1) = c + 0 = c
y(2) = cy(1) + delta(2) = c^2 + 0 = c^2
y(3) = cy(2) + delta(3) = c^3 + 0 = c^3 
y(4) = cy(3) + delta(4) = c^4 + 0 = c^4
y(5) = cy(4) + delta(5) = c^5 + 0 = c^5

% derive the unit step response (only first few values) using y(-1) = 0
y(0) = cy(-1)  + u(0) = 0 + 1 = 1 
y(1) = cy(0) + u(1) = c + 1 = c + 1
y(2) = cy(1) + u(2) = c(c+1) + 1 = c^2 + c + 1
y(3) = cy(2) + u(3) = c(c^2 + c + 1) + 1 = c^3 + c^2 + c + 1 
y(4) = cy(3) + u(4) = c(c^3 + c^2 + c + 1) + 1 = c^4 + c^3 + c^2 + c + 1 
y(5) = cy(4) + u(5) = c(c^4 + c^3 + c^2 + c + 1) + 1 = c^5 + c^4 + c^3 + c^2 + c + 1 

% Draw the impulse response and unit step response for c = 0.5 (only first few values) using y(-1) = 0
% for impulse response we can use Matlab function impz
b = [1]; % coefficient of x's in system
a = [1 -0.5]; % coefficient of y's in system
impz(b,a)

% unit step response
u = ones(1,10); % snippet from unit step for n>=0
b = [1]; % coefficient of x's in system
a = [1 -0.5]; % coefficient of y's in system

resp = filter(b,a,u);
plot((0:length(u)-1), resp, '*'); % first part sets x axis to start from 0 (default by Matlab is 1)
axis([0 length(u)-1 0 2]); % set x and y coordinates

% For what values of c is IIR system stable
% (using some signal and different values of c try to figure when the IIR system is stable)
sig = [1 zeros(1,10)];
c = 0.5;
b = [1]; % coefficient of x's in system
a = [1 -c]; % coefficient of y's in system
plot(filter(b,a,sig))

c = 2;
b = [1]; % coefficient of x's in system
a = [1 -c]; % coefficient of y's in system
plot(filter(b,a,sig))

c = -2;
b = [1]; % coefficient of x's in system
a = [1 -c]; % coefficient of y's in system
plot(filter(b,a,sig))

IIR system is stable for following |c| < 1


%task-b
% y(n)  + y(n-1) - 6y(n-2) = x(n)
% derive first few values of system response for the following input signal x(n) = 8u(n), where u(n) is unit step, using initial conditions: y(-1) = 1, y(-2) = -1

b = [1]; % coefficient of x's in system
a = [1 1 -6]; % coefficient of y's in system
x = 8*(ones(1,10)) % snippet from the input signal for n>=0
resp = filter(b, a, x, filtic(b,a,[1 -1])); % filtic(b,a,[1 -1]) sets the initial conditions
stem((0:length(x)-1), resp);

% initial conditions: y(-1) = y(-2) = 0
resp = filter(b, a, x, filtic(b,a,[0 0])); % filtic(b,a,[0 0]) sets the initial conditions, when initial conditions are all zero, we can also use following resp = filter(b, a, x);
stem((0:length(x)-1), resp);

%cross-correlation
[x fs] = wavread('Generated_Rain_Scene_8000.wav');
c = wavread('Cuckoo_mono_8000_detajl.wav');
a = xcorr2(x(:,1),c); % using firs channel from record
plot(a);
