function [lambda1, lambda2, phi] = vcv_to_eig(uu, vv, uv)
% function [lambda1, lambda2, phi] = vcv_to_eig(uu, vv, uv)
%
%  Calculates eigenvalues and angle from the x-axis of
%  the principal eigenvector for each of a set of
%  2-D variance-covariance matrices.
%  For example, suppose we have time series of (u,v)
%  at each of several points.  Then
%     uu    =  vector of mean-squared u at each point;
%     vv    =                         v
%     uv    =            mean product of u and v.
%  Typically, u and v will be zero-mean perturbation quantities.
%
%  Outputs are
%     lambda1 = vector of largest eigenvalues, one for
%                       each time series
%     lambda2 = vector of second eigenvalues
%     phi     = vector of angles from the x-axis to the
%                       eigenvector corresponding to the
%                       largest eigenvalue (lambda1)
%
%  Typical application would be to plot the standard deviation
%  ellipses for a set of current vector time series.  Input
%  arguments to m_vec_ell would then include:
%
%        'semi-major', sqrt(lambda1), 'semi-minor', sqrt(lambda2), ...
%        'angles', phi
%
%  2000/03/05 EF
%

if nargin ~= 3,
   help vcv_to_eig
   error('vcv_to_eig requires exactly 3 input arguments')
end


pm = sqrt(uu.*uu + vv.*vv - 2 * uu.*vv + 4 * uv.*uv);

lambda1 = 0.5 * (uu + vv + pm);
lambda2 = 0.5 * (uu + vv - pm);
lambda2 = max(lambda2, 0);  % Otherwise it can be negative due
                            % to truncation error.
phi = atan2(lambda1 - uu, uv);
