function [x] = ART(A,y,x0,b,lambda)

    x = x0;
    for i=1:size(A,1)
        ai = A(i,:)-y;
        normSqrd = ai * ai';
        if(normSqrd > 0)
%            x = x + (lambda * (b(i,1) - ai*x) / normSqrd * ai)';
            x = x + (((b(i,1) - ai*x) / normSqrd) * ai)';
        end
    end
    