% Script fig12_30old.m; landing flare for 747 using exponential model-
% following; x=[u w q theta h dt]'; y=[u h]'; u=[de dtc]';
% units ft, sec, crad; xdot = Fx + Gu; y=Hx+Lu; 7/91, 5/28/05
%
F=[-.021 .122 0 -.322 0 1; -.209 -.530 2.21 0 0 -.044;...
.017 -.164 -.412 0 0 .544; 0 0 1 0 0 0; 0 -1 0 2.21 0 0;...
0 0 0 0 0 -.25]; G=[.010 -.064 -.378 0 0 0; 0 0 0 0 0 .25]';
H3=[1 0 0 0 0 0; 0 0 0 0 1 0; 0 0 0 0 0 1]; L3=zeros(3,2);
A=diag([1 .03 1]); B=diag([1 1/3]); k=lqr(F,G,H3'*A*H3,B);
a=1/6.5; hsw=58; hf=7; hdot0=-221*2.5*pi/180; du=15;
H=H3([1 2],:); L=zeros(2); % yc = yc0*exp(-at)
S=[F+a*eye(6) G;H L]; T=[zeros(6,2);eye(2)];
MN=S\T; M=MN([1:6],:); N=MN([7 8],:);
% Calculate initial & final equilibrium states/controls:
F4=F([1:4],[1:4]); G4=[G([1:4],1) F([1:4],6)];
H4=[1 0 0 0;0 -1 0 2.21]; L4=zeros(2);
xu0=[F4 G4;H4 L4]\[0 0 0 0 0 hdot0]';
xuf=[F4 G4;H4 L4]\[0 0 0 0 -du 0]';
dxu=xu0-xuf; tf=log((hsw+hf)/hf)/a; t=tf*[0:.01:1]';
%
% Feedforward/feedback control law: u = N*yc(t) - k*[x - M*yc(t)];
Fcl=F-G*k; Gcl=G*(N+k*M); % xdot = Fcl*x + Gcl*yc
H1=[H;[0 -1 0 2.21 0 0];-k]; L1=[L;[0 0];N+k*M];
yc=exp(-a*t)*[du hsw+hf]; x0=[dxu([1:4]); hsw+hf; dxu(6)];
y=lsim(Fcl,Gcl,H1,L1,yc,t,x0);
xc=yc*M'; uc=yc*N'; y(:,2)=y(:,2)-hf*ones(101,1);
yc(:,2)=yc(:,2)-hf*ones(101,1);
%
figure(1); clf; subplot(211), plot(t,[y(:,2) yc(:,2)]); grid;
text(2,20,'h'); text(5,30,'hc (ft)'); text(7,-18,'-hf');
axis([0 16 -20 60]); xlabel('Time (sec)')
title('747 Landing Flare Using Exponential Model Following')
subplot(212), plot(t,[y(:,1) yc(:,1)]); grid; axis([ 0 16 0 15])
xlabel('Time (sec)'); text(5,8,'du (ft/sec)'); text(3,6,'duc');
%
figure(2); clf; subplot(211), plot(t,[y(:,3) xc*H1(3,:)'])
grid; axis([0 16 -10 0]); xlabel('Time (sec)')
text(3,-3,'hdot'); text(6,-6,'hdotc');
title('747 Landing Flare Using Exponential Model Following')
subplot(212), plot(t,y(:,[4 5]),t,uc); grid; axis([ 0 16 -8 2])
xlabel('Time (sec)'); text(1,-7,'dtc'); text(3,-4,'de')
text(1,.2,'dtcc'); text(.5,-3,'dec')
