t = 2*pi*rand(1,1000);
rho = rand(1,1000);
xc=3;yc=3;
a=3;b=2;phi=atan(2/3);

x = xc + (rho*a) .* cos(t) * cos(phi) - (rho*b) .* sin(t) * sin(phi);
y = yc + (rho*a) .* cos(t) * sin(phi) + (rho*b) .* sin(t) * cos(phi);

%data = [x'-repmat(mean(x),1000,1) y'-repmat(mean(y),1000,1)];
data = [x' y'];
covdata = cov(data);

[autovetor,score,autovalor] = princomp(data);

%% cov(data) * autovetor(:,1) % A * v1
 % lambda * v1
fprintf(1,'v1 = (%7.4f,%7.4)\n',autovalor(1) * autovetor(1,1),autovalor(1) * autovetor(1,2));

%%cov(data) * autovetor(:,2) % A * v2
% lambda * v2
fprintf(1,'v2 = (%7.4f,%7.4)\n',autovalor(1) * autovetor(2,1),autovalor(1) * autovetor(2,2));

hold;
axis([0 6 0 6]);
plot(x,y,'o');
plot([xc xc+autovalor(1) * autovetor(1,1)],[yc yc+autovalor(1) * autovetor(2,1)],'r-d');
plot([xc xc+autovalor(2) * autovetor(1,2)],[yc yc+autovalor(2) * autovetor(2,2)],'r-d');