TP-probleme-inverse-3D/TP3/MVS.m

61 lines
2 KiB
Mathematica
Raw Permalink Normal View History

2023-06-25 14:38:01 +00:00
function erreur_reprojection = MVS(I, k_ref, masque_ref, K, R, t, valeurs_z)
% Indices des pixels du masque de ref
ind_masque = transpose(find(masque_ref > 0));
[i_masque, j_masque] = ind2sub(size(masque_ref), ind_masque);
% Passage du repère pixels au repère image :
u_masque = j_masque;
v_masque = i_masque;
p_tilde = [ u_masque ; v_masque ; ones(1, length(ind_masque)) ];
erreur_reprojection = zeros(size(p_tilde, 2), size(valeurs_z, 2));
% construction de la liste des indices des images sur lesquelles on va reprojeter
valeurs_k = [1:k_ref-1, k_ref+1:size(t, 2)];
% on va tester plusieurs valeurs de déprojection
for j = 1:size(valeurs_z, 2)
z = valeurs_z(j);
% calcul de w à partir de p_tilde
Q = K \ p_tilde * z;
% on déprojete le point dans le repère monde
Qm = R(:,:,k_ref)' * (Q - t(:, k_ref));
% pour chacune des autres image (non référence)
for i = 1:size(valeurs_k, 2)
k = valeurs_k(i);
% on reprojète le point (repère monde) dans le repère de l'image k
Qk = R(:,:,k) * Qm + t(:,k);
wk = Qk ./ Qk(3,:);
% on transforme Q en pixels
pk_tilde = K * wk;
% on arrondi les coords des pixels
pk_tilde = round(pk_tilde);
for o = 1:size(pk_tilde, 2)
i_pixel_reproj = pk_tilde(2, o); % v
j_pixel_reproj = pk_tilde(1, o); % u
i_pixel_ref = p_tilde(2, o); % v
j_pixel_ref = p_tilde(1, o); % u
invisible = i_pixel_reproj <= 0 || i_pixel_reproj > 540 || j_pixel_reproj <= 0 || j_pixel_reproj > 540;
if invisible
erreur_reprojection(o, k) = Inf;
else
erreur = I(i_pixel_ref, j_pixel_ref, k_ref) - I(i_pixel_reproj, j_pixel_reproj, k);
erreur = erreur^2;
erreur_reprojection(o, j) = erreur;
end
end
end
end
end