function  [Q, R, P, r] = varpro_qr (X, defeps);%VARPRO_QR      Special QR decomposition%% [Q, R, P, r] = varpro_qr (X, defeps);% makes Q*X*P == R ( == [R1(r,r), R2(r,n-r); 0(m-r,r), 0(m-r,n-r)]).%% defeps: limit for rank-deficiency%% r: rank of R (with respect to defeps)  [m, n] = size(X);  [Q, R, P] = qr(X);    normr = norm(R);  r = 1;  while (r <= min(m,n)),    if (norm(R(r,r:n)) < defeps*normr), break; end;    r = r + 1;  end  r = r - 1;   Q = Q';end % varpro_qr