num=1; den=[1,0,0]; G=tf(num,den); 
wn=[10,50,100]; t=0:.01:1; 
f1=figure;  f2=figure; ngrid('new');
[a,b,c,d]=tf2ss(num,den); 
p1=-0.1; a_shift=a-p1*eye(size(a));
S_shift=mksys(a_shift,b,c,d);
for i=1:length(wn)
   G1=std_tf(1,wn(i),2); G2=1-G1; 
   w1=[G2.den{1};G2.num{1}]; w1(2,3)=0.001;
   TSS_shift=augtf(S_shift,w1,[1e-5;1],[]);
   [gg,ss_Fopt_shift]=hinfopt(TSS_shift);
   [af_shift,bf,cf,df]=branch(ss_Fopt_shift); 
   af=af_shift+p1*eye(size(af_shift));
   Gc=ss(af,bf,cf,df); 
   G_o=G*Gc; G_c=feedback(G_o,1);
   figure(f2); nichols(G_o);  pause
   figure(f1); step(G_c,t); hold on;
end
figure(f2), axis([-360,0,-40,40])

zpk(Gc)