rand('state',0)
zos=0.32891; tt=4.2354; zeta=0.51452;
u0=0.2; du=-0.5; pos=-2.36; 
[y,t]=step(tf(zos,[tt^2 2*tt*zeta 1]));
y = du * y + pos+ 0.005*zos*du*rand(size(y));

plot(t,y)
ylabel('y')
xlabel('t')
set(1, 'PaperUnits', 'centimeters')
set(1, 'PaperPosition', [3 10 11 8.25])
print -deps 'id21exm2a.eps'

y0 = y(1);
yn = (y-y0)/du;
plot(t,yn)

[y1,i1]=max(yn);
t1=t(i1);

yy=yn(i1:end);
[y2,j]=min(yy);
i2=i1+j-1;
t2=t(i2);

K = yn(end);
M=(y1-y2)/y1;
zeta = abs( log(M)/sqrt(pi^2+log(M)^2));
omega0 = pi/((t2-t1)*sqrt(1-zeta^2));
T=1/omega0;

K = fix(K*100+0.5)/100
zeta = fix(zeta*100+0.5)/100
T = fix(T*100+0.5)/100

g1=tf(K, [T^2 2*T*zeta 1]);
[yo,to]=step(g1,t(end));
plot(t,yn,'k',to,yo,'k--')
text(t1,y1,'+', 'HorizontalAlignment', 'center')
text(t2,y2,'+', 'HorizontalAlignment', 'center')
ylabel('y')
xlabel('t')
legend('y', 'y^n',4)
print -deps 'id21exm2b.eps'
