function [x_k,rho,eta] = tsvd(U,s,V,b,k) %TSVD Truncated SVD regularization. % % [x_k,rho,eta] = 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), ... ] . % % The solution and residual norms are returned in eta and rho.  % Per Christian Hansen, IMM, 12/21/97.  % Initialization. [n,p] = size(V); lk = length(k); if (min(k)<0 | max(k)>p)   error('Illegal truncation parameter k') end x_k = zeros(n,lk); eta = zeros(lk,1); rho = zeros(lk,1); beta = U(:,1:p)'*b; xi = beta./s;  % Treat each k separately. for j=1:lk   i = k(j);   if (i>0)     x_k(:,j) = V(:,1:i)*xi(1:i);     eta(j) = norm(xi(1:i));     rho(j) = norm(beta(i+1:p));   end end  if (nargout > 1 & size(U,1) > p)   rho = sqrt(rho.^2 + norm(b - U(:,1:p)*beta)^2); end