% IR: The estimated impulse response (IR(1) corresponds to g(0))
%
% [IR,R,CL] = cra(Z,M,NA) gives access to
% M: The number of lags for which the functions are computed (def 20)
% NA: The order of the whitening filter.(Def 10). With NA=0, no prewhite-
% ning is performed. Then the covariance functions of the original
% data are obtained.
% R: The covariance/correlation information
% R(:,1) contains the lag indices
% R(:,2) contains the covariance function of y (poss. prewhitened)
% R(:,3) contains the covariance function of u (poss. prewhitened)
% R(:,4) contains the correlation function between (poss prewhitened)
% u and y (positive lags correponds to an influence from u to y)
% CL is the 99 % significance level for the impulse response
% These functions are also plotted by the cra. The plots can be
% redisplayed by cra(R);
%
% L. Ljung 10-2-90,1-25-92
% Copyright (c) 1990-92 by the MathWorks, Inc.
% All Rights Reserved.
[Ncap,nz]=size(z);redisp=0;
if nz==4, if z((Ncap+1)/2,1)==0, redisp=1;end,end
if nz>2 & ~redisp,error('This routine is for two data records only. To check correlations between more signals, use cra for two signals at a time.'),end
if redisp, r=z; n=0;end
clf reset
if ~redisp
if nargin<3, n=[];end
if nargin<2, M=[];end
if isempty(n),n=10;end
if isempty(M),M=20;end
if nz==1, error('For a single signal, use covf instead!'),end
if n>0
a=th2poly(ar(z(:,2),n,'fb0'));
z(:,1)=filter(a,1,z(:,1));
z(:,2)=filter(a,1,z(:,2));
end
M=M-1;
Rcap=covf(z,M+1);
r(:,1)=[-M:1:M]';
r(M+1:2*M+1,2:3)=Rcap([1 4],:)';
r(1:M,2:3)=Rcap([1 4],M+1:-1:2)';
scir=Rcap(4,1); sccf=sqrt(Rcap(1,1)*Rcap(4,1));
r(M+1:2*M+1,4)=Rcap(2,:)'/sccf;
r(1:M,4)=Rcap(3,M+1:-1:2)'/sccf;
ir=r(M+1:2*M+1,4)*sccf/scir;
end
subplot(221)
plot(r(:,1),r(:,2)),
if n>0,title('Covf for filtered y'),else title('Covf for y'),end
subplot(222)
plot(r(:,1),r(:,3)),
if n>0,title('Covf for prewhitened u'),else title('Covf for u'),end
subplot(223)
plot(r(:,1),r(:,4))
if n>0,title('Correlation from u to y (prewh)'),else title('Correlation from u to y'),end
