  function xBasic = LowRankLS(W,Y,b)
% W,Y are m-by-r with full column rank
% b is m-by-1
% xbasic minimizes norm(WY'x - b,2) and has at most r nonzero components

  [m,r] = size(W);
% First solve min norm(Wz - b)

% Get the thin qr factorization...
  [Q,R] = qr(W,0);
  z = R\(Q'*b);

% Now determine a suitable x so Y'*x = z.

% QR with column pivoting...
  [U,S,E] = qr(Y');

% We want z = Y'*xBasic = USE'*xBasic
  ztilde = U'*z;

% If 
  xtilde = S(:,1:r)\ztilde;
  xBasic = E*[xtilde;zeros(m-r,1)];
  
% -3 if dont use qr(*,0)
% -8 if you use svd
% -17 if you for WY'