% [lambda,V] = varimax(lambda[,sgn,tol]) Varimax rotation
%
% Performs varimax rotation (Kaiser 1958) on the column vectors contained
% in lambda. We follow Kaiser's notation.
%
% In:
%   lambda: DxL matrix (L<D) of column vectors (possibly the loadings from
%      a factor analysis).
%   sgn: preferred majoritary sign for the vector components (-1:
%      negative, 1: positive, 0: don't care). Default 0.
%   tol: minimum relative increase of the objective function to keep on
%      iterating (default 1e-5).
% Out:
%   lambda: DxL matrix with the rotated loadings.
%   V: training curve for the objective function.
%
% Any non-mandatory argument can be given the value [] to force it to take
% its default value.

% Copyright (c) 1997,2003 by Miguel A. Carreira-Perpinan

function [lambda,V] = varimax(lambda,sgn,tol)

% Argument defaults
if ~exist('sgn','var') | isempty(sgn) sgn = 0; end;
if ~exist('tol','var') | isempty(tol) tol = 1e-5; end;

[D,L] = size(lambda);

h = sqrt(sum(lambda.^2,2)) + eps;			% Communalities
temp = lambda./repmat(h,1,L);
V = [sum(temp(:).^4)-sum(sum(temp.^2,1).^2)/D];		% Objective function

V_old = V*(1-2*tol);

while abs(V(end)-V_old) > tol*abs(V(end))
  
  V_old = V(end);
  
  for i=1:L-1
    for j=i+1:L
      
      % Optimal angle to rotate columns i, j
      x = lambda(:,i)./h;
      y = lambda(:,j)./h;
      u = x.*x - y.*y;
      v = 2*x.*y;
      t = atan2( 2*(D*u'*v-sum(u)*sum(v)), D*(u'*u-v'*v)-sum(u)^2+sum(v)^2 )/4;
      
      % Anticlockwise rotation of angle t (t+pi is valid, too)
      temp = lambda(:,[i j])*[cos(t) -sin(t); sin(t) cos(t)];
      lambda(:,[i j]) = temp;
      
    end
  end
  
  % New value of the objective function
  h = sqrt(sum(lambda.^2,2)) + eps;			% Communalities
  temp = lambda./repmat(h,1,L);
  V = [V sum(temp(:).^4)-sum(sum(temp.^2,1).^2)/D];	% Objective function
  
end

% Sign inversion so that each column vector of lambda has mainly
% components of the same sign
if sgn ~= 0
  temp = find(sgn*sum(lambda,1) < 0);
  lambda(:,temp) = -lambda(:,temp);
end