G=tf(10, [1,10,35,50,24]);
nyquist(G); axis('square'); axis([-0.6,0.6,-0.4,0.4])
[Kc,pp,wg,wp]=margin(G); [Kc,wg], Tc=2*pi/wg;

[Gc1,Kp1]=ziegler(1,[Kc,Tc,10]); Kp1

[Gc2,Kp2,Ti2]=ziegler(2,[Kc,Tc,10]); [Kp2,Ti2]

[Gc3,Kp3,Ti3,Td3]=ziegler(3,[Kc,Tc,10]); [Kp3,Ti3,Td3]

G_c1=feedback(G*Gc1,1); step(G_c1); hold on
G_c2=feedback(G*Gc2,1); step(G_c2); 
G_c3=feedback(G*Gc3,1); step(G_c3);