commit ca6da606526bb044fc21c2a39671921645ee68c4 Author: Laureηt Date: Sat Jun 10 21:05:32 2023 +0200 init diff --git a/tp1/creation_cercle_et_donnees_bruitees.p b/tp1/creation_cercle_et_donnees_bruitees.p new file mode 100755 index 0000000..f2bb35a Binary files /dev/null and b/tp1/creation_cercle_et_donnees_bruitees.p differ diff --git a/tp1/donnees.m b/tp1/donnees.m new file mode 100755 index 0000000..9acbdd0 --- /dev/null +++ b/tp1/donnees.m @@ -0,0 +1,38 @@ +clear; +close all; + +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +% Fenetre d'affichage : +figure('Name','Points situes au voisinage d''un cercle', ... + 'Position',[0.4*L,0.05*H,0.6*L,0.7*H]); +axis equal; +hold on; +set(gca,'FontSize',20); +hx = xlabel('$x$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$y$','FontSize',30); +set(hy,'Interpreter','Latex'); + +% Bornes d'affichage des donnees centrees en (0,0) : +taille = 20; +bornes = [-taille taille -taille taille]; + +% Creation du cercle et des donnees bruitees : +n = 50; +sigma = 0.5; +[x_cercle,y_cercle,x_donnees_bruitees,y_donnees_bruitees,theta_donnees_bruitees] ... + = creation_cercle_et_donnees_bruitees(taille,n,sigma); + +% Affichage du cercle : +plot(x_cercle([1:end 1]),y_cercle([1:end 1]),'r','LineWidth',3); + +% Affichage des donnees bruitees : +plot(x_donnees_bruitees,y_donnees_bruitees,'k+','MarkerSize',10,'LineWidth',2); +axis(bornes); +lg = legend(' Cercle', ... + ' Donnees bruitees', ... + 'Location','Best'); +grid on; diff --git a/tp1/donnees_occultees.m b/tp1/donnees_occultees.m new file mode 100755 index 0000000..82f1bf1 --- /dev/null +++ b/tp1/donnees_occultees.m @@ -0,0 +1,52 @@ +clear; +close all; + +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +% Fenetre d'affichage : +figure('Name','Points situes au voisinage d''un cercle', ... + 'Position',[0.4*L,0.05*H,0.6*L,0.7*H]); +axis equal; +hold on; +set(gca,'FontSize',20); +hx = xlabel('$x$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$y$','FontSize',30); +set(hy,'Interpreter','Latex'); + +% Bornes d'affichage des donnees centrees en (0,0) : +taille = 20; +bornes = [-taille taille -taille taille]; + +% Creation du cercle et des donnees bruitees : +n = 50; +sigma = 0.5; +[x_cercle,y_cercle,x_donnees_bruitees,y_donnees_bruitees,theta_donnees_bruitees] ... + = creation_cercle_et_donnees_bruitees(taille, n, sigma); + +theta = rand(1, 2)*2*pi; +theta = sort(theta); + +while theta(2)-theta(1) < pi/2 + + theta = rand(1, 2)*2*pi; + theta = sort(theta); + +end + +x_donnees_bruitees = x_donnees_bruitees( theta(1) <= theta_donnees_bruitees & theta_donnees_bruitees <= theta(2) ); +y_donnees_bruitees = y_donnees_bruitees( theta(1) <= theta_donnees_bruitees & theta_donnees_bruitees <= theta(2) ); +theta_donnees_bruitees = theta_donnees_bruitees( theta(1) <= theta_donnees_bruitees & theta_donnees_bruitees <= theta(2) ); + +% Affichage du cercle : +plot(x_cercle([1:end 1]),y_cercle([1:end 1]),'r','LineWidth',3); + +% Affichage des donnees bruitees : +plot(x_donnees_bruitees,y_donnees_bruitees,'k+','MarkerSize',10,'LineWidth',2); +axis(bornes); +lg = legend(' Cercle', ... + ' Donnees bruitees', ... + 'Location','Best'); +grid on; diff --git a/tp1/exercice_1.m b/tp1/exercice_1.m new file mode 100755 index 0000000..a198875 --- /dev/null +++ b/tp1/exercice_1.m @@ -0,0 +1,36 @@ +donnees; + +n_tests = 500; + +% Estimation de la position du centre : +[C_estime,R_moyen] = estimation_1(x_donnees_bruitees,y_donnees_bruitees,n_tests); + +% Affichage du cercle estime : +n_points_cercle = 100; +theta_cercle = 2*pi/n_points_cercle:2*pi/n_points_cercle:2*pi; +x_cercle_estime = C_estime(1)+R_moyen*cos(theta_cercle); +y_cercle_estime = C_estime(2)+R_moyen*sin(theta_cercle); +plot(x_cercle_estime([1:end 1]),y_cercle_estime([1:end 1]),'b','LineWidth',3); +lg = legend(' Cercle initial', ... + ' Donnees bruitees', ... + ' Cercle estime', ... + 'Location','Best'); + +function [C_estime, R_moyen] = estimation_1(x_donnees_bruitees, y_donnees_bruitees, n_tests) + + G = mean( [ x_donnees_bruitees.' y_donnees_bruitees.' ] ); + R_moyen = mean( sqrt((x_donnees_bruitees.'-G(1)).^2 + (y_donnees_bruitees.'-G(2)).^2) ); + + x = repmat(x_donnees_bruitees.', 1, n_tests); + y = repmat(y_donnees_bruitees.', 1, n_tests); + + x_rand = repmat((rand(1, n_tests)*R_moyen - R_moyen/2) + G(1), length(x_donnees_bruitees), 1); + y_rand = repmat((rand(1, n_tests)*R_moyen - R_moyen/2) + G(2), length(y_donnees_bruitees), 1); + + dist = (sqrt((x-x_rand).^2 + (y-y_rand).^2) - R_moyen).^2; + + [min_val, min_index] = min( sum( dist ) ); + + C_estime = [ x_rand(1, min_index) y_rand(1, min_index) ]; + +end \ No newline at end of file diff --git a/tp1/exercice_2.m b/tp1/exercice_2.m new file mode 100755 index 0000000..cd31d0a --- /dev/null +++ b/tp1/exercice_2.m @@ -0,0 +1,42 @@ +donnees; + +n_tests = 1000; + +% Estimation du rayon et de la position du centre : +[C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + +% Affichage du cercle estime : +n_points_cercle = 100; +theta_cercle = 2*pi/n_points_cercle:2*pi/n_points_cercle:2*pi; +x_cercle_estime = C_estime(1)+R_estime*cos(theta_cercle); +y_cercle_estime = C_estime(2)+R_estime*sin(theta_cercle); +plot(x_cercle_estime([1:end 1]),y_cercle_estime([1:end 1]),'b','LineWidth',3); +lg = legend(' Cercle initial', ... + ' Donnees bruitees', ... + ' Cercle estime', ... + 'Location','Best'); + + function [C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + + G = mean( [ x_donnees_bruitees.' y_donnees_bruitees.' ] ); + R_moyen = mean( sqrt((x_donnees_bruitees.'-G(1)).^2 + (y_donnees_bruitees.'-G(2)).^2) ); + + x = repmat(x_donnees_bruitees.', 1, n_tests, n_tests); + y = repmat(y_donnees_bruitees.', 1, n_tests, n_tests); + + x_rand = repmat((rand(1, n_tests)*2 - 1) + G(1), length(x_donnees_bruitees), 1, n_tests); + y_rand = repmat((rand(1, n_tests)*2 - 1) + G(2), length(y_donnees_bruitees), 1, n_tests); + + r_rand = rand(1, 1, n_tests) - 1/2 + R_moyen; + R = repmat(r_rand, length(x_donnees_bruitees), n_tests, 1); + + dist = (sqrt((x-x_rand).^2 + (y-y_rand).^2) - R).^2; + + somme = reshape(sum(dist), [n_tests, n_tests]); + + [min_val, idx] = min(somme(:)); + [row, col] = ind2sub(size(somme), idx); + + C_estime = [ x_rand(1, row ) y_rand(1, row ) ]; + R_estime = r_rand(col); +end \ No newline at end of file diff --git a/tp1/exercice_3.m b/tp1/exercice_3.m new file mode 100755 index 0000000..0d9e649 --- /dev/null +++ b/tp1/exercice_3.m @@ -0,0 +1,42 @@ +donnees_occultees; + +n_tests = 1000; + +% Estimation du rayon et de la position du centre : +[C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + +% Affichage du cercle estime : +n_points_cercle = 100; +theta_cercle = 2*pi/n_points_cercle:2*pi/n_points_cercle:2*pi; +x_cercle_estime = C_estime(1)+R_estime*cos(theta_cercle); +y_cercle_estime = C_estime(2)+R_estime*sin(theta_cercle); +plot(x_cercle_estime([1:end 1]),y_cercle_estime([1:end 1]),'b','LineWidth',3); +lg = legend(' Cercle initial', ... + ' Donnees bruitees', ... + ' Cercle estime', ... + 'Location','Best'); + + function [C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + + G = mean( [ x_donnees_bruitees.' y_donnees_bruitees.' ] ); + R_moyen = mean( sqrt((x_donnees_bruitees.'-G(1)).^2 + (y_donnees_bruitees.'-G(2)).^2) ); + + x = repmat(x_donnees_bruitees.', 1, n_tests, n_tests); + y = repmat(y_donnees_bruitees.', 1, n_tests, n_tests); + + x_rand = repmat((rand(1, n_tests)*3*R_moyen - 1.5*R_moyen) + G(1), length(x_donnees_bruitees), 1, n_tests); + y_rand = repmat((rand(1, n_tests)*3*R_moyen - 1.5*R_moyen) + G(2), length(y_donnees_bruitees), 1, n_tests); + + r_rand = rand(1, 1, n_tests)*3*R_moyen - 1.5*R_moyen; + R = repmat(r_rand, length(x_donnees_bruitees), n_tests, 1); + + dist = (sqrt((x-x_rand).^2 + (y-y_rand).^2) - R).^2; + + somme = reshape(sum(dist), [n_tests, n_tests]); + + [min_val, idx] = min(somme(:)); + [row, col] = ind2sub(size(somme), idx); + + C_estime = [ x_rand(1, row ) y_rand(1, row ) ]; + R_estime = r_rand(col); +end \ No newline at end of file diff --git a/tp1/exercice_4.m b/tp1/exercice_4.m new file mode 100755 index 0000000..6bbad55 --- /dev/null +++ b/tp1/exercice_4.m @@ -0,0 +1,42 @@ +donnees_occultees; + +n_tests = 1000; + +% Estimation du rayon et de la position du centre : +[C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + +% Affichage du cercle estime : +n_points_cercle = 100; +theta_cercle = 2*pi/n_points_cercle:2*pi/n_points_cercle:2*pi; +x_cercle_estime = C_estime(1)+R_estime*cos(theta_cercle); +y_cercle_estime = C_estime(2)+R_estime*sin(theta_cercle); +plot(x_cercle_estime([1:end 1]),y_cercle_estime([1:end 1]),'b','LineWidth',3); +lg = legend(' Cercle initial', ... + ' Donnees bruitees', ... + ' Cercle estime', ... + 'Location','Best'); + + function [C_estime,R_estime] = estimation_2(x_donnees_bruitees,y_donnees_bruitees,n_tests); + + G = mean( [ x_donnees_bruitees.' y_donnees_bruitees.' ] ); + R_moyen = mean( sqrt((x_donnees_bruitees.'-G(1)).^2 + (y_donnees_bruitees.'-G(2)).^2) ); + + x = repmat(x_donnees_bruitees.', 1, n_tests, n_tests); + y = repmat(y_donnees_bruitees.', 1, n_tests, n_tests); + + x_rand = repmat((randn(1, n_tests)*R_moyen - R_moyen/2) + G(1), length(x_donnees_bruitees), 1, n_tests); + y_rand = repmat((randn(1, n_tests)*R_moyen - R_moyen/2) + G(2), length(y_donnees_bruitees), 1, n_tests); + + r_rand = randn(1, 1, n_tests)*R_moyen - R_moyen/2; + R = repmat(r_rand, length(x_donnees_bruitees), n_tests, 1); + + dist = (sqrt((x-x_rand).^2 + (y-y_rand).^2) - R).^2; + + somme = reshape(sum(dist), [n_tests, n_tests]); + + [min_val, idx] = min(somme(:)); + [row, col] = ind2sub(size(somme), idx); + + C_estime = [ x_rand(1, row ) y_rand(1, row ) ]; + R_estime = r_rand(col); +end \ No newline at end of file diff --git a/tp2/donnees.m b/tp2/donnees.m new file mode 100755 index 0000000..c1363ec --- /dev/null +++ b/tp2/donnees.m @@ -0,0 +1,70 @@ +clear; +close all; + +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +% Fenetre d'affichage : +figure('Name','Points situes au voisinage d''une droite', ... + 'Position',[0.4*L,0,0.6*L,0.8*H]); +axis equal; +hold on; +set(gca,'FontSize',20); +hx = xlabel('$x$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$y$','FontSize',30); +set(hy,'Interpreter','Latex'); + +% Bornes d'affichage des donnees centrees en (0,0) : +taille = 20; +bornes = [-taille taille -taille taille]; + +% Parametres de la droite : +theta_D = 2*pi*(rand-0.5); +cos_theta_D = cos(theta_D); +sin_theta_D = sin(theta_D); +marge = 5; +rho_D = (taille-marge)*rand; + +% Affichage de la droite : +droite_horizontale = abs(cos_theta_D)-taille & ... + x_donnees_bruitees-taille & ... + y_donnees_bruitees=y_min & y<=y_max + x_extremites(i) = x; + y_extremites(i) = y; + i = i+1; + end + + % Cas ou D_k coupe le bord droit : + x = x_max; + y = (rho_k-x*cos_theta_k)/sin_theta_k; + if y>=y_min & y<=y_max + x_extremites(i) = x; + y_extremites(i) = y; + i = i+1; + end + + % Cas ou D_k coupe le bord superieur : + y = y_min; + x = (rho_k-y*sin_theta_k)/cos_theta_k; + if x>=x_min & x<=x_max + x_extremites(i) = x; + y_extremites(i) = y; + i = i+1; + end + + % Cas ou D_k coupe le bord inferieur : + y = y_max; + x = (rho_k-y*sin_theta_k)/cos_theta_k; + if x>=x_min & x<=x_max + x_extremites(i) = x; + y_extremites(i) = y; + end + + % Affichage de la droite D_k : + plot(x_extremites,y_extremites,couleur,'LineWidth',2); +end diff --git a/tp3/cercle_3_points.m b/tp3/cercle_3_points.m new file mode 100755 index 0000000..19464a7 --- /dev/null +++ b/tp3/cercle_3_points.m @@ -0,0 +1,13 @@ +function [C,R] = cercle_3_points(x,y) + +x1 = x(1); +y1 = y(1); +x2 = x(2); +y2 = y(2); +x3 = x(3); +y3 = y(3); + +A = 2*[-x1+x2 -y1+y2 ; -x1+x3 -y1+y3]; +B = [-x1^2+x2^2-y1^2+y2^2 ; -x1^2+x3^2-y1^2+y3^2]; +C = A\B; +R = sqrt((x1-C(1))^2+(y1-C(2))^2); diff --git a/tp3/creation_cercle_et_donnees_bruitees.p b/tp3/creation_cercle_et_donnees_bruitees.p new file mode 100755 index 0000000..5066819 Binary files /dev/null and b/tp3/creation_cercle_et_donnees_bruitees.p differ diff --git a/tp3/donnees_aberrantes.m b/tp3/donnees_aberrantes.m new file mode 100755 index 0000000..26320f8 --- /dev/null +++ b/tp3/donnees_aberrantes.m @@ -0,0 +1,44 @@ +clear; +close all; + +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +% Fenetre d'affichage : +figure('Name','Points situes au voisinage d''un cercle', ... + 'Position',[0.4*L,0.05*H,0.6*L,0.7*H]); +axis equal; +hold on; +set(gca,'FontSize',20); +hx = xlabel('$x$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$y$','FontSize',30); +set(hy,'Interpreter','Latex'); + +% Limites de l'affichage des donnees centrees en (0,0) : +taille = 20; +limites_affichage = [-taille taille -taille taille]; + +% Creation du cercle et des donnees bruitees : +n = 50; +sigma = 0.5; +[x_cercle,y_cercle,theta_cercle,x_donnees_bruitees,y_donnees_bruitees,theta_donnees_bruitees] ... + = creation_cercle_et_donnees_bruitees(taille,n,sigma); + +% Affichage du cercle : +plot(x_cercle([1:end 1]),y_cercle([1:end 1]),'r','LineWidth',3); + +% Donnees aberrantes : +proportion = 0.25; +n_donnees_aberrantes = floor(proportion*n); +x_donnees_bruitees = [x_donnees_bruitees taille*(2*rand(1,n_donnees_aberrantes)-1)]; +y_donnees_bruitees = [y_donnees_bruitees taille*(2*rand(1,n_donnees_aberrantes)-1)]; + +% Affichage des donnees bruitees : +plot(x_donnees_bruitees,y_donnees_bruitees,'k+','MarkerSize',10,'LineWidth',2); +axis(limites_affichage); +lg = legend(' Cercle', ... + ' Donnees (bruitees + aberrantes)', ... + 'Location','Best'); +grid on; diff --git a/tp3/exercice_0.m b/tp3/exercice_0.m new file mode 100755 index 0000000..26ed58e --- /dev/null +++ b/tp3/exercice_0.m @@ -0,0 +1,86 @@ +clear; +close all; +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +% Parametres : +alpha = pi/180*12; % Seuil sur l'orientation du gradient (en radians) +cos_alpha = cos(alpha); +seuil_norme = 2/sin(alpha); % Seuil sur la norme du gradient (decoule d'une etude theorique) + +% Chargement des donnees : +% load Donnees/parquet; +%load Donnees/parquet; +load Donnees/bateau; +[n_lignes,n_colonnes] = size(I); +limites_affichages = [1 n_colonnes 1 n_lignes]; + +figure('Name','Detection des alignements','Position',[0.2*L,0,0.8*L,H]); + +% Affichage de l'image : +subplot(2,2,1); +imagesc(I); +axis equal; +axis off; +colormap gray; +title('Image originale','FontSize',30); +axis(limites_affichages); + +% Fenetre d'affichage des pixels de contour : +subplot(2,2,3); +imagesc(I); +axis equal; +axis off; +colormap gray; +hold on; + +% Gradient du niveau de gris (x vers la droite, y vers le bas) : +I = double(I); +[G_x,G_y] = gradient(I); +G_norme = sqrt(G_x.^2+G_y.^2); + +% Selection des pixels de contour : +contour = G_norme>seuil_norme; + +% Pas de pixel de contour sur le bord de l'image => Traitement simplifie ! +contour([1 n_lignes],:) = 0; +contour(:,[1 n_colonnes]) = 0; + +% Tri des pixels de contour : +indices_contour = find(contour); +[~,indices] = sort(G_norme(indices_contour),'descend'); +indices_contour = indices_contour(indices); + +% Affichage d'une petite fleche sur les pixels de contour : +[i,j] = ind2sub(size(I),indices_contour); +quiver(j,i,G_x(indices_contour),G_y(indices_contour),'r'); % Attention : x = j et y = i +axis equal; +title('Pixels de contour','FontSize',30); +axis(limites_affichages); +drawnow; + +% Affichage des ensembles E (la fonction label2rgb donne a chaque ensemble E une couleur differente) : +subplot(2,2,2); +imagesc(I_resultat); +axis(limites_affichages); +axis equal; +axis off; +title('Ensembles candidats','FontSize',30); + +% Affichage de l'esquisse : +subplot(2,2,4); +imagesc(120*ones(size(I)),[0,255]); +axis equal; +axis off; +title('Alignements','FontSize',30); +hold on; +axis(limites_affichages); + +% Boucle sur les ensembles E : +for k = 1:size(extremites,3) + % Affichage du segment : + line(extremites(1,:,k),extremites(2,:,k),'Color','w','LineWidth',2); +end + +save exercice_0; diff --git a/tp3/exercice_0.mat b/tp3/exercice_0.mat new file mode 100755 index 0000000..d9a96a3 Binary files /dev/null and b/tp3/exercice_0.mat differ diff --git a/tp3/exercice_1.m b/tp3/exercice_1.m new file mode 100755 index 0000000..33a8bc6 --- /dev/null +++ b/tp3/exercice_1.m @@ -0,0 +1,81 @@ +clear; +close all; +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +load exercice_0; + +% Estimation du point de fuite : +[rho_F,theta_F] = estimation_F(rho,theta); + +% Coordonnees cartesiennes du point de fuite : +x_F = rho_F*cos(theta_F); +y_F = rho_F*sin(theta_F); + +figure('Name','Estimation du point de fuite','Position',[0.7*L,0,0.3*L,H]); + +% Affichage des points de coordonnees (rho,theta) : +subplot(2,1,1); +plot(theta,rho,'k+','MarkerSize',10,'LineWidth',2); +axis([-pi pi -1.2*rho_F 1.2*rho_F]); +set(gca,'FontSize',20); +hx = xlabel('$\theta$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$\rho$','FontSize',30); +set(hy,'Interpreter','Latex'); +grid; +hold on; + +% Affichage de la sinusoide correspondant au point de fuite : +pas = 0.01; +theta_affichage = -pi:pas:pi; +rho_affichage = rho_F*cos(theta_affichage-theta_F); +plot(theta_affichage,rho_affichage,'b-','LineWidth',3); + +% Affichage des points de coordonnees (rho,theta) : +plot(theta,rho,'r+','MarkerSize',10,'LineWidth',2); +title('Sinusoide estimee'); + +% Affichage de l'image : +subplot(2,1,2); +imagesc(I); +set(gca,'FontSize',20); +axis ij equal off; +colormap gray; +hold on; + +% Limites des affichages : +marge = round(min(n_lignes,n_colonnes)/10); +x_min = min(1,x_F)-marge; +x_max = max(n_colonnes,x_F)+marge; +y_min = min(1,y_F)-marge; +y_max = max(n_lignes,y_F)+marge; +limites_affichages = [x_min x_max y_min y_max]; +axis(limites_affichages); + +% Affichage d'une selection des droites formant le premier faisceau : +taille_selection = min(length(rho),10); +longueurs_segments_au_carre = (extremites(1,1,:)-extremites(1,2,:)).^2 ... + +(extremites(2,1,:)-extremites(2,2,:)).^2; +[~,indices_tries] = sort(longueurs_segments_au_carre,'descend'); +selection = indices_tries(1:taille_selection); +affichage_faisceau(rho(selection),theta(selection),limites_affichages,'r'); + +% Affichage du point de fuite : +plot(x_F,y_F,'bx','MarkerSize',20,'LineWidth',5); +title('Point de fuite'); + +function [rho_F, theta_F] = estimation_F(rho, theta) + + B = rho; + A = [cos(theta.') ; sin(theta.')].'; + + X = (A.'*A)\(A.'*B); + + rho_F = sqrt( X(1)^2 + X(2)^2 ); + theta_F = atan2(X(2), X(1)); + +end + + diff --git a/tp3/exercice_2.m b/tp3/exercice_2.m new file mode 100755 index 0000000..1b388cd --- /dev/null +++ b/tp3/exercice_2.m @@ -0,0 +1,201 @@ +clear; +close all; +taille_ecran = get(0,'ScreenSize'); +L = taille_ecran(3); +H = taille_ecran(4); + +load exercice_0; + +% Parametres de l'algorithme RANSAC : +n_donnees = length(rho); +S1 = 5; +S2 = 0.3; +k_max = floor(nchoosek(n_donnees,2)/n_donnees); +parametres = [S1 S2 k_max]; + +% Estimation du premier point de fuite : +[rho_F_1,theta_F_1] = RANSAC_2(rho,theta,parametres); + +% Coordonnees cartesiennes du premier point de fuite : +x_F_1 = rho_F_1*cos(theta_F_1); +y_F_1 = rho_F_1*sin(theta_F_1); + +% Droites conformes au premier point de fuite : +conformes_1 = abs(rho-rho_F_1*cos(theta-theta_F_1))<=S1; +rho_conformes_1 = rho(conformes_1); +theta_conformes_1 = theta(conformes_1); + +% Droites restantes : +theta = theta(~conformes_1); +rho = rho(~conformes_1); + +% Estimation du deuxieme point de fuite : +[rho_F_2,theta_F_2] = RANSAC_2(rho,theta,parametres); + +% Coordonnees cartesiennes du deuxieme point de fuite : +x_F_2 = rho_F_2*cos(theta_F_2); +y_F_2 = rho_F_2*sin(theta_F_2); + +% Droites conformes au deuxieme point de fuite : +conformes_2 = abs(rho-rho_F_2*cos(theta-theta_F_2))<=S1; +rho_conformes_2 = rho(conformes_2); +theta_conformes_2 = theta(conformes_2); + +figure('Name','Estimation de la ligne de fuite','Position',[0.3*L,0,0.7*L,H]); + +% Affichage des points de coordonnees (rho,theta) : +subplot(2,2,2); +plot(theta,rho,'k+','MarkerSize',10,'LineWidth',2); +rho_F_max = max(rho_F_1,rho_F_2); +axis([-pi pi -1.2*rho_F_max 1.2*rho_F_max]); +set(gca,'FontSize',20); +hx = xlabel('$\theta$','FontSize',30); +set(hx,'Interpreter','Latex'); +hy = ylabel('$\rho$','FontSize',30); +set(hy,'Interpreter','Latex'); +grid; +hold on; + +% Affichage de la sinusoide correspondant au premier point de fuite : +pas = 0.01; +theta_affichage = -pi:pas:pi; +rho_affichage_1 = rho_F_1*cos(theta_affichage-theta_F_1); +plot(theta_affichage,rho_affichage_1,'b-','LineWidth',3); + +% Affichage des points conformes a la premiere sinusoide : +plot(theta_conformes_1,rho_conformes_1,'r+','MarkerSize',10,'LineWidth',2); + +% Affichage de la deuxieme sinusoide : +rho_affichage_2 = rho_F_2*cos(theta_affichage-theta_F_2); +plot(theta_affichage,rho_affichage_2,'m-','LineWidth',3); + +% Affichage des points conformes a la deuxieme sinusoide : +plot(theta_conformes_2,rho_conformes_2,'g+','MarkerSize',10,'LineWidth',2); +title('Paire de sinusoides estimees'); + +% Affichage de l'image : +subplot(2,2,3); +imagesc(I); +set(gca,'FontSize',20); +axis ij equal off; +colormap gray; +hold on; + +% Limites des affichages : +marge = round(min(n_lignes,n_colonnes)/10); +x_min = min(1,x_F_1); +x_min = min(x_min,x_F_2)-marge; +x_max = max(n_colonnes,x_F_1); +x_max = max(x_max,x_F_2)+marge; +y_min = min(1,y_F_1); +y_min = min(y_min,y_F_2)-marge; +y_max = max(n_lignes,y_F_1); +y_max = max(y_max,y_F_2)+marge; +limites_affichages = [x_min x_max y_min y_max]; +axis(limites_affichages); + +% Affichage d'une selection de droites formant le premier faisceau : +taille_selection_1 = min(length(rho_conformes_1),10); +longueurs_segments_au_carre_1 = (extremites(1,1,conformes_1)-extremites(1,2,conformes_1)).^2 ... + +(extremites(2,1,conformes_1)-extremites(2,2,conformes_1)).^2; +[~,indices_tries_1] = sort(longueurs_segments_au_carre_1,'descend'); +selection_1 = indices_tries_1(1:taille_selection_1); +affichage_faisceau(rho_conformes_1(selection_1),theta_conformes_1(selection_1),limites_affichages,'r'); + +% Affichage du premier point de fuite : +plot(x_F_1,y_F_1,'bx','MarkerSize',20,'LineWidth',5); + +% Affichage d'une selection de droites formant le deuxieme faisceau : +taille_selection_2 = min(length(rho_conformes_2),10); +longueurs_segments_au_carre_2 = (extremites(1,1,conformes_2)-extremites(1,2,conformes_2)).^2 ... + +(extremites(2,1,conformes_2)-extremites(2,2,conformes_2)).^2; +[~,indices_tries_2] = sort(longueurs_segments_au_carre_2,'descend'); +selection_2 = indices_tries_2(1:taille_selection_2); +affichage_faisceau(rho_conformes_2(selection_2),theta_conformes_2(selection_2),limites_affichages,'g'); + +% Affichage du deuxieme point de fuite : +plot(x_F_2,y_F_2,'mx','MarkerSize',20,'LineWidth',5); +title('Paire de points de fuite'); + +% Affichage de l'image : +subplot(2,2,4); +imagesc(I); +set(gca,'FontSize',20); +axis ij equal off; +colormap gray; +hold on; +axis(limites_affichages); + +% Affichage du premier point de fuite : +plot(x_F_1,y_F_1,'bx','MarkerSize',20,'LineWidth',5); + +% Affichage du deuxieme point de fuite : +plot(x_F_2,y_F_2,'mx','MarkerSize',20,'LineWidth',5); + +% Affichage de la ligne de fuite : +line([x_F_1 x_F_2],[y_F_1 y_F_2],'Color','c','LineWidth',2); +title('Ligne de fuite'); + +% Affichage de l'image originale : +subplot(2,2,1); +imagesc(I); +set(gca,'FontSize',20); +axis ij equal off; +colormap gray; +hold on; +title('Image originale'); +axis(limites_affichages); + + +function [rho_F_1, theta_F_1] = RANSAC_2(rho, theta, parametres) + + n = length(rho); + erreur_F_1 = +inf; + + for j = 1:parametres(3) + + points_model = randperm(n, 2); + points_candidat = []; + + for i = setdiff(1:n, points_model) + + [rho_candidat, theta_candidat, erreur_candidat] = estimation_F_2(rho([points_model i]), theta([points_model i])); + + if erreur_candidat < parametres(1) + + points_candidat = [ points_candidat i ]; + + end + + end + + if length(points_candidat)/n > parametres(2) + + [rho_new, theta_new, erreur_new] = estimation_F_2(rho([points_model points_candidat]), theta([points_model points_candidat])); + + if erreur_new < erreur_F_1 + + rho_F_1 = rho_new; + theta_F_1 = theta_new; + erreur_F_1 = erreur_new; + + end + + end + + end +end + +function [rho_F, theta_F, erreur_F] = estimation_F_2(rho, theta) + + B = rho; + A = [cos(theta.') ; sin(theta.')].'; + + X = (A.'*A)\(A.'*B); + + rho_F = sqrt( X(1)^2 + X(2)^2 ); + theta_F = atan2(X(2), X(1)); + + erreur_F = mean( abs( rho - rho_F*cos(theta - theta_F) ) ); + +end \ No newline at end of file diff --git a/tp3/exercice_3.m b/tp3/exercice_3.m new file mode 100755 index 0000000..d6f6a92 --- /dev/null +++ b/tp3/exercice_3.m @@ -0,0 +1,103 @@ +donnees_aberrantes; + +n_tests = 1000; + +% Parametres de l'algorithme RANSAC : +n_donnees = length(x_donnees_bruitees); +S1 = 2; +S2 = 0.5; +k_max = floor(nchoosek(n_donnees,3)/n_donnees); +parametres = [S1 S2 k_max]; + +% Estimation du rayon et de la position du centre : +[C_estime,R_estime] = RANSAC_3(x_donnees_bruitees,y_donnees_bruitees,parametres,n_tests); + +% Affichage du cercle estime : +x_cercle_estime = C_estime(1)+R_estime*cos(theta_cercle); +y_cercle_estime = C_estime(2)+R_estime*sin(theta_cercle); +plot(x_cercle_estime([1:end 1]),y_cercle_estime([1:end 1]),'b','LineWidth',3); + +% Affichage des points conformes au modele : +conformes = abs(sqrt((x_donnees_bruitees-C_estime(1)).^2+ ... + (y_donnees_bruitees-C_estime(2)).^2)-R_estime)<=S1; +plot(x_donnees_bruitees(conformes), ... + y_donnees_bruitees(conformes),'b+','MarkerSize',10,'LineWidth',2); +lg = legend(' Cercle initial', ... + ' Donnees (bruitees + aberrantes)', ... + ' Cercle estime', ... + ' Donnees conformes', ... + 'Location','Best'); + +function [C_estime,R_estime] = RANSAC_3(x_donnees_bruitees,y_donnees_bruitees,parametres,n_tests) + + n = length(x_donnees_bruitees); + erreur = +inf; + + for j = 1:parametres(3) + + j/parametres(3)*100 + + points_model = randperm(n, 3); + [C_model, R_model] = cercle_3_points(x_donnees_bruitees(points_model), y_donnees_bruitees(points_model)); + + points_candidat = []; + + for i = setdiff(1:n, points_model) + + [C_candidat, R_candidat] = estimation(x_donnees_bruitees([points_model i]), y_donnees_bruitees([points_model i]), n_tests); + + % ma fonction d'erreur ne doit pas être la bonne mais ça marche + % a peu près bien ~~ + erreur_candidat = mean( abs( sqrt((C_model(1)-x_donnees_bruitees([points_model i])).^2 + (C_model(2)-y_donnees_bruitees([points_model i])).^2 )) - R_model ); + + if erreur_candidat < parametres(1) + + points_candidat = [ points_candidat i ]; + + end + + end + + if length(points_candidat)/n > parametres(2) + + [C_new, R_new] = estimation(x_donnees_bruitees([points_model points_candidat]), y_donnees_bruitees([points_model points_candidat]), n_tests); + + erreur_new = mean( abs( sqrt((C_model(1)-x_donnees_bruitees([points_model points_candidat])).^2 + (C_model(2)-y_donnees_bruitees([points_model points_candidat])).^2 )) - R_model ); + + if erreur_new < erreur + + C_estime = C_new; + R_estime = R_new; + erreur = erreur_new; + + end + + end + + end +end + +function [C_estime,R_estime] = estimation(x_donnees_bruitees,y_donnees_bruitees,n_tests) + + G = mean( [ x_donnees_bruitees.' y_donnees_bruitees.' ] ); + R_moyen = mean( sqrt((x_donnees_bruitees.'-G(1)).^2 + (y_donnees_bruitees.'-G(2)).^2) ); + + x = repmat(x_donnees_bruitees.', 1, n_tests, n_tests); + y = repmat(y_donnees_bruitees.', 1, n_tests, n_tests); + + x_rand = repmat((randn(1, n_tests)*R_moyen - R_moyen/2) + G(1), length(x_donnees_bruitees), 1, n_tests); + y_rand = repmat((randn(1, n_tests)*R_moyen - R_moyen/2) + G(2), length(y_donnees_bruitees), 1, n_tests); + + r_rand = randn(1, 1, n_tests)*R_moyen - R_moyen/2; + R = repmat(r_rand, length(x_donnees_bruitees), n_tests, 1); + + dist = (sqrt((x-x_rand).^2 + (y-y_rand).^2) - R).^2; + + somme = reshape(sum(dist), [n_tests, n_tests]); + + [min_val, idx] = min(somme(:)); + [row, col] = ind2sub(size(somme), idx); + + C_estime = [ x_rand(1, row ) y_rand(1, row ) ]; + R_estime = r_rand(col); + end \ No newline at end of file