% LQG control of the system G(s) = (1+0.5s)/s^2

alpha = 0.5;
A = [0 0 ; 1 0];
B = [1 ; 0];
C = [alpha 1];
D = 0;
sys = ss(A,B,C,D);

Q=diag([0 1]);
R=0.7e-5;

[k,s,e]=lqr(A,B,Q,R);

G=B;
H=0;
QN = 1;
RN = R*10;
syse = ss(A,[B G],C,[D H]);
[kest,l,p]=kalman(syse,QN,RN);

rlqg = lqgreg(kest,k);
PC = -rlqg*sys;

Ckj = tf([3628 40000],[1 80.28 156.56]),
PCkj = sys*Ckj;
frv= logspace(-2,3,2000);
figure(1)
[re1,im1]=nyquist(PC,2*pi*frv);
plot(re1(:),im1(:),'b','Linewidth',2);hold on
[re2,im2]=nyquist(PCkj,2*pi*frv);
plot(re2(:),im2(:),'r','Linewidth',2)
axis([-5 5 -5 5])
grid
hold on
set(gca,'Fontsize',16)
set(gca,'Fontweight','bold')
print -depsc exkjfig2a.eps

figure(2)
am1=squeeze(freqresp(1/(1+PC),2*pi*frv));
am2=squeeze(freqresp(1/(1+PCkj),2*pi*frv));
semilogx(frv,20*log10(am1),'b','Linewidth',2);hold on
semilogx(frv,20*log10(am2),'r','Linewidth',2)
hold on
grid on
set(gca,'Fontsize',16)
set(gca,'Fontweight','bold')
print -depsc exkjfig2b.eps

figure(3)
am3=squeeze(freqresp(rlqg,2*pi*frv));
am4=squeeze(freqresp(Ckj,2*pi*frv));
semilogx(frv,20*log10(am3),'b','Linewidth',2);hold on
semilogx(frv,20*log10(am4),'r','Linewidth',2)
hold on
grid on
set(gca,'Fontsize',16)
set(gca,'Fontweight','bold')
print -depsc exkjfig2c.eps

zpk(rlqg)