function [f,P] = pisarenko(Ryy)
%
% Compute the the modal frequencies using Pisarenko's method,
% then find the amplitudes
%
% function [f,P] = pisarenko(Ryy)
%
% Ryy = autocorrelation function from observed data
%
% f = vector of frequencies
% P = vector of amplitudes

% Copyright 1999 by Todd K. Moon

% Debugged by G. Lathoud - lathoud@idiap.ch
  
[n,n] = size(Ryy);
[v,u] = eig(Ryy);
sigma2 = u(n,n);
vm = v(:,n);               % get eigenvector from smallest eigenvalue
r = roots(vm);             % find roots of the polynomial
f = angle(r)/(2*pi); % get angle and convert to Hz/sample
% Now compute the amplitudes

% Build up the coefficient matrix
for k=1:n-1
   A(k,:) = exp(2*pi*j*k .* f).';
end
r = Ryy(2:n,1);
P = A\r;