function [Q,R] = mgs(A)
[m,n] = size(A);
Q = zeros(m,n); R = zeros(n,n);
V = A;
for ii = 1:n
    R(ii,ii) = norm(V(:,ii));
    Q(:,ii) = V(:,ii)/R(ii,ii);
    for jj = ii+1:n
        R(ii,jj) = Q(:,ii)'*V(:,jj);
        V(:,jj) = V(:,jj)- R(ii,jj)*Q(:,ii);
    end
end
end
