view arrows/stats/stats_pca.m @ 61:eff6bddf82e3 tip

Finally implemented perceptual brightness thing.

function pca=stats_pca(stats,I)
	m1=mstep(stats);
	[V,D]=eig(m1.cov);
	J=(1+size(D,1))-I;
	pca.eigvecs=V(:,J);
	pca.eigvals=diag(D);
	pca.eigvals=max(0,pca.eigvals(J));
	pca.mean=m1.mu;
end
function M=mstep(S)
	M.mu=S.amu/S.ag;
	M.cov=msym(S.aco/S.ag - M.mu*M.mu'); 
end

%function M=mstep_z(S)
%	M.mu=zeros(size(S.aco,1),1);
%	M.cov=S.aco/S.ag; 
%end

function A=msym(X), A=(X+X')/2; end