% FUNCTION [ x, e ] = LS_SVD( A, b )
%
% Least square solution of an overdetermined system Ax = b
% where A  is  a M-by-N matrix with M > N
%       x  is  a M-by-1 vector
%       b  is  a M-by-1 vector
%
% This is done using Singular Value Decomposition of A.
% Also returned: e (scalar) the euclidian norm of Ax-b

function [ x, e ] = ls_svd( A, b )
  
  error( nargchk( 1, 2, nargin ) );
  
  M = size( A, 1 );
  N = size( A, 2 );
  
  [ U, S, V ] = svd( A );
  
  c = U' * b;
  
  y = c( 1:N ) ./ diag( S );
  
  x = V * y;
  
  e = norm( A*x-b );