%《MATLAB及其在理工課程中的應用指南》第九章例9－2－2程序
% 由差分方程遞推求輸出
% 西安電子科技大學出版社出版  陳懷琛編著 1999年10月初版，2004年11月第二版
%
clear
a = input('差分方程左端的系數向量a=[a(1),...a(na)]=(書上取 [1,-1,0.9] ) ');
b = input('差分方程右端的系數向量b=[b(1),...b(na)]=(書上取 [2,3] ) ');
u = input('輸入信號序列u=(書上取 cos(0.5*[1:59]) ) ');
na=length(a);nb=length(b);nu=length(u);
s=['起算點前',int2str(na-1),'點ym的值 =[ym(1),..,ym(na-1)]=(書上取 [0,0] )'];
ym=zeros(1,na+nu);ym(1:na-1) = input(s);    % 建立ym序列并賦予初值
um = [zeros(1,na-1),u];					    % 建立um序列并賦予初值
for n=na:na+nu-1							% 這個n以ym的起點為準
   ys = ym(n-1:-1:n-na+1); us = um(n:-1:n-nb+1);	% 生成us及ys
   ym(n) = (b*us'-a(2:na)*ys')/a(1);		% 差分方程遞推求ym
end
% 把ym時間坐標右移na-1位,求出y
y = ym(na:na+nu-1);                         % 給出起點和終點的y數組
stem(y),grid
line([0,60],[0,0])		                    % 畫橫坐標軸
set(gcf,'color','w')                        % 置圖形背景色為白色