function [B] =REF (A,b )
% This function computes the row echelon form of 
% the augemented matrix [A b] by using Gauss elimination 
% with pivoting by row (partial pivoting)
% 
% Input: A -> matrix m x n with real coefficients
%        b -> column vector m x 1 with real coefficients
%
% Output B -> row echelon matrix (size m x (n+1))

if length(b)~=size(A,1)
    error('Please make sure that b and A are compatible')
end


A =[A b];

[m,n]=size(A);

for kk=1:m-1

    [a,idx]=max(abs(A(kk:end,kk)));
    
    if idx~=kk   % PIVOTING STEP
        q = A(kk, kk:end);
        A(kk,kk:end)=A(idx+kk-1,kk:end);
        A(idx+kk-1,kk:end)=q;
    end
    
     for jj=kk+1:m % ELIMINATION STEP
         A(jj,kk:end)=A(jj,kk:end)-A(kk,kk:end)*A(jj,kk)/A(kk,kk);
     end
end
    
B=A;        
    
    

end

