%--------------------------------------------------------------------------------------
%   Impulse response of
%           y(n)=0.5y(n-1)-0.7y(n-2)+0.75x(n)+0.75x(n-3)
%   Created by Aly El-Osery
%   Modified on : Sept. 1, 2006
%
%--------------------------------------------------------------------------------------

%               0.75 + 0.75 z^(-3)          z^(-3)     0.75 z^3 + 0.75      
%   H(z)= ----------------------------- = ---------  ------------------- = z^(-1) H4(z)
%          1 - 0.5 z^(-1) + 0.7 z^(-2)      z^(-2)     z^2 -0.5 z + 0.7     
%
%
%
%                           0.75 - 1.0714 z^(-1) + 0.535 z^(-2)
%   H1(z)= 1.0714 z^(-1) + -------------------------------------
%                               1 - 0.5 z^(-1) + 0.7 z^(-2)
%
%
%
%                    0.375 z^2 - 0.525 z + 0.75 
%   H4(z)= 0.75 z + ----------------------------
%                        z^2 - 0.5 z + 0.7
%
%
%--------------------------------------------------------------------------------------
clear all
close all
%--------------------------------------------------------------------------------------
% The following is the computation of H(z).
h=dimpulse([0.75 0 0 0.75],[1 -0.5 0.7 0]);
n1=0:length(h)-1;
subplot(3,2,1),stem(n1,h)        
xlabel('n')
ylabel('h(n)')
grid on
%--------------------------------------------------------------------------------------
% The following is the computation of H1(z).
h=dimpulse([0.75 -1.0714 0.535],[1 -0.5 0.7],length(n1));
h(2)=h(2)+1.0714;             % This addes the extra impulse with amplitude 1.0714
subplot(3,2,2),stem(n1,h)        
xlabel('n')
ylabel('h_1(n)')
grid on

%--------------------------------------------------------------------------------------
% The following is the computation of h(n) using the direct method

n=0:length(n1)-1;
h2(n+1)=(-0.00765+0.43372j)*(0.25+0.7984j).^n+(-0.00765-0.43372j)*(0.25-0.7984j).^n;
h2(1)=h2(1)+0.7653;h2(2)=h2(2)+1.0714;
subplot(3,2,3),stem(n,h2,'r')
xlabel('n')
ylabel('h_2(n)')
grid on

%--------------------------------------------------------------------------------------
% The following is the computation of h(n) using state-space representation

a=[0 1 0;0 0 1;0 -0.7 0.5];b=[0;0;1];c=[0.75 -0.7*0.75 0.5*0.75];an=eye(3);
h3(1)=0.75;
h3(2)=c*an*b;
for i=3:length(n1),an=a*an;h3(i)=c*an*b;end
n3=0:length(h3)-1;
subplot(3,2,4),stem(n3,h3,'g')
xlabel('n')
ylabel('h_3(n)')
grid on

%--------------------------------------------------------------------------------------
% The following is the computation of H4(z)

h4=dimpulse([0.375 -0.525 0.75],[1 -0.5 0.7],length(n1));
h4prime(2:length(h4)+1)=h4;     % This shifts the sequence to the right
h4prime(1)=0.75;                % Adds the impulse at n=0
n4=0:length(h4prime)-2;
subplot(3,2,5),stem(n4,h4prime(1:length(n4)),'k')
xlabel('n')
ylabel('h_4(n)')
grid on
