clear all
clc
l=[50];
m=[1 5 50];
G=tf(l,m);
subplot(4,1,1)
step(G)
T=5000;
dt=0.01
t=0:dt:T-dt;
x=randn(size(t));
subplot(4,1,2)
plot(t,x)
y=lsim(G,x,t);
y=y';
subplot(4,1,3)
lsim(G,x,t)

X=fft(x);
Y=fft(y);

PSDxy=(2*dt^2*conj(X).*Y)/T;
PSDxx=(2*dt^2*conj(X).*X)/T;

Gxy=PSDxy./PSDxx;
Gxy=decimate(Gxy,10);
subplot(4,1,4)
plot(Gxy)