function [rho_estime, N_estime] = estimation(I, S, masque)

    m = zeros(3, length(I));
    masque2 = find(masque);
    m(:, masque2(:)) = S \ I(:, masque2(:));

    rho_estime = sqrt(sum(m.^2, 1));
    N_estime = m ./ rho_estime;

    rho_estime = reshape(rho_estime, size(masque));

end