function w = gauss_legendre( L, nodes)
%
% w = gauss_legendre( L, nodes)
%
% Compute generalized Gauss-Legendre quadrature weights for H_L, that is 
% the space spanned by all Legendre polynomials with l = 0 ... L for the 
% given nodes. Note that at least (L+1) nodes are required although we do
% allow more nodes.
%
% @return w  quadrature weights
% 
% @param L   maximum degree of Legendre polynomials
% @param nodes  quadrature node locations (in [-1,1]

  assert( numel(nodes) >= (L+1));
  assert( (min(nodes) >= -1.0) && max(nodes) <= 1.0);
  
  % compute kernel matrix
  K = zeros( numel(nodes), L+1);
  for( i = 1 : numel(nodes))
    for( l = 0 : L)
      p_lm = legendre( l, nodes(i), 'norm')';
      K(i,l+1) = p_lm(1);
    end
  end
    
  % invert matrix
  S = pinv(K);
  
  % the scaled first row of the kernel matrix are the sought weights
  w = sqrt(2) * S(1,:);
  
end