clc;
clear;
load u.mat; %输入序列
load y.mat; %输出序列
figure(4);%第11张图
subplot(2,1,1)
plot(u);
title('脉冲输入序列');
subplot(2,1,2)
plot(y);
title('脉冲输出序列');
figure(5)
subplot(2,1,1)
[a,b] = autocorr(y) ; %a 为各阶的相关系数，b 为滞后阶数
autocorr(y); %原序列的自相关函数图MA(q),观察系数是否在区间(-2T^(1/2),-2T^(1/2))内
subplot(2,1,2);
[c,d] = parcorr(y) ; %c 为各阶的偏自相关系数，d 为滞后阶数
parcorr(y); %画出偏自相关图
z=iddata(y,u,0.01);%将DX转化为matlab接受的格式
test = [];
for p= 1:3 %自回归AR对应PACF
for q = 1:4 %移动平均对应ACF
m = armax(z,'na',p,'nb',q);
AIC = aic(m); %armax(p,q),选择对应FPE最小，AIC值最小的模型
%[H, P, Qstat, CV] = lbqtest(z, [p;q], 0.05) %Ljung-Box Q-statistic lack-of-fit hypothesis test
test= [test;p q AIC];
end
end
for k = 1:size(test,1)
if test(k,3) == min(test(:,3)) %选择AIC值最小的模型
p_test = test(k,1);
q_test = test(k,2);
break;
end
end
m = armax(z,'na',p_test,'nb',q_test,'nk',4);
t1=1:0.01:6;
ltiview('impulse',m,t1);
H=impulse(m,t1);
figure(6);
h=log(H);
plot(h);
title('衰减曲线');
