function [ W, U, E ] = ObliProj( V, Wperp, opt, tol, ind )
% ObliProj constructs an oblique projection matrix onto the span of the columns of V along
% the span of the columns of Wperp
%
% for opt=[1,1] computes the projector as E=V*W, where W=pinv(U'*V)*U'
%          with U=V-orth_proj_{Wperp} V
%
% for opt=[2,1] computes the projector as E=V*W, where W=pinv(U'*U)*U'
%          with U=V-OrthProj_{Wperp} V
%
% for opt=[3,1] computes the projector as E=V*W, where W=pinv(U)'
%          with U=V-OrthProj_{Wperp} V
%
% for opt=[1,2] computes the projector as E=V*W, where W=pinv(Q'*V)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=DREOr(U)
%
% for opt=[1,3] computes the projector as E=V*W, where W=pinv(Q'*V)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=qr(U)
%
% for opt=[1,4] computes the projector as E=V*W, where W=pinv(Q'*V)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=orth(U)
%
% for opt=[2,2] computes the projector as E=V*W, where W=pinv(Q'*U)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=DREOr(U)
%
% for opt=[2,3] computes the projector as E=V*W', where W'=pinv(Q'*U)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=qr(U)
%
% for opt=[2,4] computes the projector as E=V*W', where W'=pinv(Q'*U)*Q'
%          with  U=V-OrthProj_{Wperp}  and Q=orth(U)
%
% Usage:    [ W, U, E ] = ObliProj( V, Wperp, opt, tol, ind );
%           [ W, U, E ] = ObliProj( V, Wperp );
%
% Inputs:
%   V       matrix the columns of which span then space to project onto
%   Wperp   matrix the columns of which span then space to project along
%   opt     array to chose the method to calculate the projector (see above) default 
%           opt=[3,1]
%   tol     if opt=[1,2], or [2,2]  tol is the torance set for considering linearly 
%           independent columns [default tol= 1.0000e-7]                 
%   ind     (optional) the indices of vectors to start the orthogonalization (see DREOr)
%   []      can be used for ind, tol and opt
%
% Outputs:
%   W       matrix producing E=U*W;
%   U       matrix producing E=U*W;
%   E       Projector onto span of columns of V onto columns of Wperp   
%
% see orth DREOrp DREOr qr

% See   http://www.ncrg.aston.ac.uk/Projects/HNLApprox/
%       http://www.nonlinear-approx.info/

% setting of default tolerance parameter
if (nargin<5) |(isempty(ind)==1); ind=[];     end
if (nargin<4) |(isempty(ind)==1); tol=1e-7;   end
if (nargin<3) |(isempty(ind)==1); opt=[3,1];  end
%
%
tol_ort=1e-7;%it's needed for opt_ort=2 only see ORTH_PROJ to change opt_ort
opt_ort=2;
Wperp=NormDict(Wperp);
%
[pp,Q]=OrthProj(Wperp,opt_ort,tol_ort);
%
U=V-pp*V;
%
opt1=opt(1);
opt2=opt(2);
if opt2==1
     switch opt1;
     case 1
     W=pinv(U'*V)*U';
     case 2 
     W=pinv(U'*U)*U';
     case 3
     W=pinv(U);
     otherwise
     disp('wrong option compute W')
     end
else
   if opt1==1 & opt2~=1
    switch opt2
     case 2 
     or=DREOr(U,tol,ind);
     case 3 
     [or,rr]=qr(U,0);
     case 4
     or=orth(U);
     otherwise
     disp('wrong option to orthogonalize vectors U')
    end
    W=pinv(or'*V)*or';
  else 
     if opt1==2 & opt2~=1
      switch opt2
       case 2
       or=DREOr(U,tol,ind);
       case 3
       [or,rr]=qr(U,0);
       case 4
       or=orth(U);
       otherwise
       disp('wrong option to orthogonalize vectors U')
      end 
      W=pinv(or'*U)*or';
     end
   end 
end
E=V*W;
%
%Laura REBOLLO-NEIRA 2008
%
