function x_k = tsvd(U,s,V,b,k)
%TSVD Truncated SVD regularization.
%
% x_k = tsvd(U,s,V,b,k)
%
% Computes the truncated SVD solution
%    x_k = V(:,1:k)*inv(diag(s(1:k)))*U(:,1:k)'*b .
% If k is a vector, then x_k is a matrix such that
%    x_k = [ x_k(1), x_k(2), ... ] .

% Per Christian Hansen, UNI-C, 11/18/91.

% Initialization.
[n,p] = size(V); lk = length(k);
if (min(k)<1 | max(k)>p)
  error('Illegal truncation parameter k')
end
x_k = zeros(n,lk); xi = (U(:,1:p)'*b)./s;

% Treat each k separately.
for j=1:lk
  i = k(j);
  x_k(:,j) = V(:,1:i)*xi(1:i);
end
