num=400; den=[1,2,400]; G=tf(num,den);
nW1=100*conv([0.005,1],[0.005,1]); 
dW1=conv([0.2,1],[0.2,1]);
nW3=[1,0,0]; dW3=[40000];
[a,b,c,d]=tf2ss(num,den); S=mksys(a,b,c,d); 
W1=[nW1; dW1]; W3=[nW3; 0,0,dW3];
TSS_=augtf(S,W1,[],W3);
[cFH2,ccL]=h2lqg(TSS_); 
[af,bf,cf,df]=branch(cFH2);
F2=zpk(ss(af,bf,cf,df))

delta=-10:1:10; 
f1=figure; f2=figure; ngrid('new');
for i=1:length(delta)
   G.den{1}(2)=delta(i);
   G_o2=G*F2; G_c2=feedback(G_o2,1); 
   figure(f1); step(G_c2), hold on
   figure(f2), nichols(G_o2)
end
figure(f2), axis([-360,0,-40,40]),
figure(f1), set(gca,'Xlim',[0,0.2]),