《DSP using MATLAB》Problem 4.23

代码：

%% ------------------------------------------------------------------------
%%            Output Info about this m-file
fprintf('
***********************************************************
');
fprintf('        <DSP using MATLAB> Problem 4.23 

');

banner();
%% ------------------------------------------------------------------------


% ----------------------------------------------------
%               1        H1(z)
% ----------------------------------------------------

b = [1, -0.866, -1.2278, 0.7938];  nb = [-1:2];
a = [1, -0.7, -0.81, 0.567]*0.7;   na = [0:3];            %  

[p, np, r, nr] = deconv_m(b, nb, a, na)

b1 = [-0.1660, -0.4178, 0.2268]
[R, p, C] = residuez(b1, a)

Mp = (abs(p))'
Ap = (angle(p))'/pi

%% ------------------------------------------------------
%%   START a    determine Y(z) and sketch    
%% ------------------------------------------------------
figure('NumberTitle', 'off', 'Name', 'P4.23 Y(z) its pole-zero plot')
set(gcf,'Color','white'); 
zplane(b,a);
title('pole-zero plot'); grid on;

%% ----------------------------------------------
%%    END
%% ----------------------------------------------

% ------------------------------------
%                  h(n)  
% ------------------------------------

[delta, n] = impseq(0, 0, 19); 
h_check = filter(b, a, delta)                                             % check sequence

[delta1, n1] = sigshift(delta, n, -1);
h_answer1 = 1.4286*delta1;
h_answer2 = 0.2323*(-0.9).^n.*stepseq(0,0,19) - 1.1257*(0.9).^n.*stepseq(0,0,19) + 0.6563*(0.7).^n.*stepseq(0,0,19);         % answer sequence

[h_answer, n2] = sigadd(h_answer1, n1, h_answer2, n)

　　运行结果：

        Y(z)进行部分分式展开，

        其留数及对应极点为：

零极点图如下

        输出y(n)的前20个元素（这里变量名称h_check、h_answer用的不合适。）