disp('Running MATLAB script file lmtzcf53.m') ;
%
%    Makes Figure 5.3, for the Locantore, Marron, Tripoli,
%    Zhang, Cohen paper.
%
%    Copied from cs1.m  (part 4)


seed = 20943343 ;
rand('seed',seed) ;
randn('seed',seed) ;

n = 60 ;
xdat = -2 + .6 * randn(n,1) ;
ydat = -2 - (xdat + 2) + 0.2 * randn(n,1) ;
xdat = [xdat; 3] ;
ydat = [ydat; 3] ;

plot(xdat,ydat,'o') ;
  axis([-4,4,-4,4]) ;
  axis('square') ;
  set(gca,'XTickLabels',[]) ;
  set(gca,'YTickLabels',[]) ;

vmean = rmean([xdat, ydat]) ;
r = 1.5 ;
theta = 2 * pi * (0:(1/100):1) ;
cx = vmean(1) + r * cos(theta) ;
cy = vmean(2) + r * sin(theta) ;
hold on ;
  plot(cx,cy,'--') ;
hold off ;

cxdat = xdat - vmean(1) ;
cydat = ydat - vmean(2) ;
vlength = sqrt(cxdat.^2 + cydat.^2) ;
sxdat = vmean(1) + r * cxdat ./ vlength ;
sydat = vmean(2) + r * cydat ./ vlength ;

hold on ;
  plot(sxdat,sydat,'+') ;
hold off ;


mcov = cov([sxdat, sydat]) ;
[meigv,eigval] = eig(mcov) ; 
eigval = diag(eigval) ;
[temp,visort] = sort(eigval) ;
visort = flipud(visort) ;
meigv = meigv(:,visort) ;
eigval = eigval(visort) ;


hold on ;
  plot([vmean(1); (vmean(1) + 2 * eigval(1) * meigv(1,1))], ...
       [vmean(2); (vmean(2) + 2 * eigval(1) * meigv(2,1))],'-') ;
    vhand = get(gca,'Children') ;
    set(vhand(1),'LineWidth',1.5) ;
  plot([vmean(1); (vmean(1) + 2 * eigval(2) * meigv(1,2))], ...
       [vmean(2); (vmean(2) + 2 * eigval(2) * meigv(2,2))],'-') ;
    vhand = get(gca,'Children') ;
    set(vhand(1),'LineWidth',0.5) ;

hold off ;


  orient tall ;
  set(gcf,'PaperPosition',[2.5, 4.5, 4.5, 4.5]) ;
        %  2.0 inches on left and right, 4.0 interior, 8.5 total
        %  4.5 from bottom, 4.5 interior, 11 total
print -deps \matlab\steve\ps\lmtzcf53.eps ;




