
%------graphical illustration----------------------------
disp('Example to show gain in different directions for matrix ')
A=[2 4 ; 0 3]
[UU,SS,VV]=svd(A)
eig(A)
U=[]; Y=[];
uyfig=10; ufig=211; yfig=212;
figure(uyfig); title('test')
subplot(ufig); hold on; grid on 
axis([-2 2 -2 2]/2*1.5); axis('equal');
subplot(yfig); hold on; grid on  
s1=max(max(SS));
axis(1.1*[-s1 s1 -s1 s1]); axis('equal');
for theta=0:pi/100:2*pi
 u1=cos(theta);
 u2=sin(theta);
 u=[u1 u2]';
 U=[U;u];
 subplot(ufig);  
   %title(['Input u= [' num2str(u',3) ']^T,   |u|=1])
 hold on; axis('equal')
 %plot([0 u1],[0 u2],'b');
 plot([0 u1],[0 u2],'b.');
 y=A*u;
 subplot(yfig);  
   %title(['Output y=Au = [' num2str(y',3) ']^T'])
 hold on; axis('equal')
 %plot([0 y(1)],[0 y(2)],'b');
 plot([0 y(1)],[0 y(2)],'b.');
 Y=[Y; y];
 %pause
 subplot(ufig);
 xlabel('u_1'); 
 ylabel('u_2'); 
 %plot([0 u1],[0 u2],'g');
 subplot(yfig);
 xlabel('y_1'); 
 ylabel('y_2'); 
 %plot([0 y(1)],[0 y(2)],'g');

end

[X, lambda]=eig(A);
% plot eigenvectors
scale_Y=s1*1.2;
plot(scale_Y*[X(1,1) -X(1,1)],scale_Y*[X(2,1) -X(2,1)],'r--')
plot(scale_Y*[X(1,2) -X(1,2)],scale_Y*[X(2,2) -X(2,2)],'r--')

%% plot eigenvectors
plot(scale_Y*[UU(1,1) -UU(1,1)],scale_Y*[UU(2,1) -UU(2,1)],'g--')
plot(scale_Y*[UU(1,2) -UU(1,2)],scale_Y*[UU(2,2) -UU(2,2)],'g--')
%
%% plot eigenvectors
plot(scale_Y*[VV(1,1) -VV(1,1)],scale_Y*[VV(2,1) -VV(2,1)],'b--')
plot(scale_Y*[VV(1,2) -VV(1,2)],scale_Y*[VV(2,2) -VV(2,2)],'b--')

xlabel('(red):eigenvectors ; (blue): V ; (green): U    A=U*S*V^T ')

subplot(ufig);
% plot eigenvectors
scale=2
plot(scale*[X(1,1) -X(1,1)],scale*[X(2,1) -X(2,1)],'r--')
plot(scale*[X(1,2) -X(1,2)],scale*[X(2,2) -X(2,2)],'r--')

%% plot eigenvectors
plot(scale*[UU(1,1) -UU(1,1)],scale*[UU(2,1) -UU(2,1)],'g--')
plot(scale*[UU(1,2) -UU(1,2)],scale*[UU(2,2) -UU(2,2)],'g--')
%
%% plot eigenvectors
plot(scale*[VV(1,1) -VV(1,1)],scale*[VV(2,1) -VV(2,1)],'b--')
plot(scale*[VV(1,2) -VV(1,2)],scale*[VV(2,2) -VV(2,2)],'b--')

disp('Press <SPACE> to iterate over inputs |u|=1')

for theta=0:pi/40:2*pi
 u1=cos(theta);
 u2=sin(theta);
 u=[u1 u2]';
 U=[U;u];
 subplot(ufig);  title(['Input u= [' num2str(u',3) ']^T,   |u|= ' num2str(norm(u),3)])
 hold on; axis('equal')
 plot([0 u1],[0 u2],'b');
 plot([0 u1],[0 u2],'b*');
 y=A*u;
 subplot(yfig);  title(['y=Gu = [' num2str(y',3) ']^T,      |y|= ' num2str(norm(y),3)])
 hold on; axis('equal')
 plot([0 y(1)],[0 y(2)],'b');
 plot([0 y(1)],[0 y(2)],'b*');
 Y=[Y; y];
 pause
 subplot(ufig);  
 plot([0 u1],[0 u2],'g');
 subplot(yfig);
 plot([0 y(1)],[0 y(2)],'g');
 
end
  

subplot(ufig); hold off;
subplot(yfig); hold off;