function [AT,BT,T]=varimax2(X,Y);

% varimax rotatie:


% als start I, dan Identiteits matrix eerste rotatie,
% als start R, dan random orthogonoale matrix.
% Varimax rotatie voor 2 deel matrices.

start='I';

delta=eps;
e=1;
count=0;

[n,k]=size(X);
[p,k]=size(Y);
V=[X;Y];

Mx=eye(n)-(1/n)*ones(n);
My=eye(p)-(1/p)*ones(p);

if start=='I' Ts=eye(k); end;
if start=='R' Ts=orth(randn(k,k)-.5); end;

B=X*Ts;
D=Y*Ts;

while e>delta,

	C=B.*B;
	E=D.*D;

	MC=Mx*C;
	ME=My*E;

	Q=[(2/n)*(B.*MC);(2/p)*(D.*ME)];

	QVVQ=Q'*V*V'*Q;
        Z=QVVQ^0.5;

	Vi=V*V'*Q*inv(Z);
	B=Vi(1:n,:);
	D=Vi(n+1:n+p,:);

	e=trace((Vi-V)'*(Vi-V));
        V=Vi;
	
	count=count+1;
        if count==100, e=0; end

	end;


C=B.*B;
E=D.*D;
T=[X;Y]'*Q*inv(Z);
max=trace((1/n)*(C'*Mx*C))+trace((1/p)*(E'*My*E));
count;
e;
AT=V(1:n,:);
BT=V(n+1:n+p,:);

