1. No category

Lycée Jules Ferry
PTSI1&2
2014-2015
Informatique
TP n°10 : Construction du plan des moindres carrés associé à
une surface nominalement plane
Lors de la phase de contrôle d’une pièce mécanique, le technicien est amené à palper les surfaces de la pièce afin
de pouvoir en reconstruire le profil. Cette opération se fait à l’aide d’une machine à mesurer tridimensionnelle.
Figure 1 – Exemple de machine à mesurer tridimensionnelle
On s’intéresse ici à la vérification au sens de la norme G.P.S de la planéité d’une surface nominalement plane.
Pour cela, on doit construire, à partir du nuage de points, le plan P passant par l’isobarycentre et qui minimise la
somme des carrés des écarts. La tolérance de planéité sera donnée par la distance entre 2 plans parallèles au plan
P et qui englobent le nuage de points. Une illustration du principe est donnée ci-dessous :
Figure 2 – Nuage de points et plan P associé
1
Figure 3 – Zoom et définition de la tolérance de planéité
Modélisation du problème
On va ici mettre en équation le problème de recherche de l’équation du plan P des moindres carrés à partir
d’un nuages de points. En fait le plan va être décrit par un point (isobarycentre des points palpés) et la normale
au plan.
1
1.1
Isobarycentre du nuage de points
Commençons par trouver les coordonnées du centre de gravité G des points palpés. Soient Mi l’ensemble des
points de coordonnées (xi , yi , zi ) dans une base de l’espace (O, ~x0 , ~y0 , ~z0 ). On rappelle que le barycentre est défini
par :
X −−→
GMi = ~0
i
D’où on tire facilement :
−−→
1 X −−→
OG =
OMi
N i
avec N le nombre de points palpés.
1.2
Détermination de la normale au plan des moindres carrés
Soit ~n la normale au plan P. L’idée est de déterminer les composantes (xn , yn , zn ) du vecteur ~n qui minimise
la somme des carrés des écarts ei . On calcule donc la somme :
S=
X
ei 2 =
i
X −−→
2
GMi · ~n
i
Figure 4 – Illustration des écarts ei entre
chaque point Mi et le plan P des moindres
carrés
On réécrit l’équation précédente en interprétant le produit scalaire par un produit matriciel. Si ~n1 = xn1 ~x0 +
yn1 ~y0 + zn1 ~z0 et ~n2 = xn2 ~x0 + yn2 ~y0 + zn2 ~z0 alors :
~n1 · ~n2 = xn1 .xn2 + yn1 .yn2 + zn1 .zn2
Ce produit scalaire peut aussi être vu comme un produit matriciel. En notant [n1 ] la matrice unicolonne des
coordonnées de ~n1 dans la base canonique (~x0 , ~y0 , ~z0 ) :
~n1 · ~n2 = t [n1 ] × [n2 ] =
t



xn2
xn1


 
 yn1  ×  yn2 
zn2
zn1


xn2


= [xn1 , yn1 , zn1 ] ×  yn2 
zn2
= xn1 .xn2 + yn1 .yn2 + zn1 .zn2
Si on s’imposeque la normale au plan P est unitaire (||~n||2 = t [n] × [n] = 1), alors le problème est finalement de
xn
 
trouver [n] =  yn  telle que la somme S soit minimale. On écrit :
zn
t
min
X t t
[GMi ] × [n] ×
[n]×[n]=1 i
2
t
[GMi ] × [n]
Ou :
t
min
t
[n] ×
[n]×[n]=1
X
[GMi ] × t [GMi ] × [n]
i
P En notant [A] la matrice 3 × 3 telle que [A] = i [GMi ] × t [GMi ] , alors on démontre que la normale [n]
recherchée est solution de l’équation dite "aux valeurs propres" suivante :
[A] × [n] = λ[n]
où le nombre λ est appelé valeur propre associée au vecteur propre [n].
On suppose ici que pour notre matrice [A], de taille 3×3, il existe 3 couples (λi , [ni ]) vérifiant l’équation précédente
(les conditions d’existence et le calcul de ces "couples propres" dépassent largement le cadre du TP, et seront vus
prochainement en cours de Mathématiques...).
Le vecteur propre qui nous intéresse ici et qui correspond donc à la normale au plan des moindres carrés sera celui
associé à la plus petite valeur propre.
2
2.1
Programmation
Tracé des points palpés et de l’isobarycentre dans l’espace
Exercice 1
Tracé des points palpés
Rappels (commandes vues au TP8) : Pour travailler sur le fichier texte avec les coordonnées des points palpés,
il faut se placer dans le bon dossier.Pour cela on utilise un module nommé os avec une commande chdir qui
s’emploie de comme suit :
import os
os.chdir("nom␣du␣dossier")
où nom du dossier est le nom du répertoire de travail.
Vous pouvez faire apparaître le contenu du dossier avec la commande os.listdir(), et connaître le dossier de
travail actif par os.getcwd().
Les coordonnées des points palpés se trouvent dans un fichier texte nommé plan. Chaque ligne est un point
différent, et pour chaque ligne les composantes sont séparées par des virgules. Il faudra tout d’abord aller lire
sur le fichier externe les cordonnées de chaque point et construire les listes listeX, listeY, listeZ contenant
respectivement les composantes sur (~x0 , ~y0 , ~z0 ) des points palpés par la MMT.
Pour cela on utilisera les commandes suivantes :
#Ouverture du fichier contenant les coordonnees des points en lecture
#(ne pas oublier l'extension du fichier a savoir .txt)
plan1=open("nom␣du␣fichier","r")
for ligne in plan1:
point=ligne.split(',')
...
Attention,pensez au type des éléments de la liste point.
Il faudra ensuite fermer le fichier par la commande plan1.close().
L’affichage des points dans une fenêtre graphique 3D se fait en appelant des commandes de la bibliothèque
mpl_toolkits.mplot3d. Le tracé utilisera donc les commandes suivantes :
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure()
ax=Axes3D(fig)
ax.set_xlabel("x")
ax.set_ylabel("y")
ax.set_zlabel("z")
ax.scatter(X,Y,Z)
Exercice 2
Calcul et tracé de l’isobarycentre du nuage de points
On cherche maintenant à calculer les coordonnées de l’isobarycentre des points palpés. Construire une fonction
3
prenant les listes listeX, listeY, et listeZ et renvoyant les coordonnées XG, YG, ZG de l’isobarycentre. Vous
pourrez placer le centre de gravité des points palpés par la commande suivante :
ax.scatter(XG,YG,ZG,color='r',marker='^')
2.2
Calcul de la normale au plan des moindres carrés P
Dans cette section, on fera appel à la bibliothèque d’algèbre linéaire de Python contenue dans le package
numpy. Pensez à importer cette bibliothèque.
Exercice 3
Construction de la matrice [A], et résolution du problème aux valeurs propres
La matrice [A] sera construite comme le produit [A] = [P ] × t [P ] (voir relation page 2).
[P ] est alors une matrice 3 × N où N est le nombre de points palpés. Pour obtenir cette matrice positionnant
chaque point palpé par rapport à l’isobarycentre, vous pourrez construire une liste PListe telle que :
PListe =
hh
XGM1
i h
XGM2 · · · XGMN , YGM1
i h
YGM2 · · · YGMN , ZGM1
ZGM2
· · · ZGMN
ii
Afin d’utiliser les outils d’algèbre linéaire de la librairie Numpy, il est nécessaire de transformer la liste PListe en
une matrice de type array. Pour cela vous utiliserez la commande :
P=np.array(PListe)
Le calcul de la matrice [A] est maintenant réalisable par :
A=np.dot(P,P.transpose())
Enfin, la solution du problème aux valeurs propres est obtenue directement grâce à la commande suivante :
W=np.linalg.eig(A)
Cette commande place la solution dans un tableau numpy dont le premier argument contient les valeurs propres,
tandis que le second retourne une matrice dont chaque colonne constitue un vecteur propre (ces vecteurs étant
bien entendu classés d’une manière respective par rapport aux valeurs propres).
Afficher cette solution W pour visualiser sa forme.
Écrire une fonction qui extrait de W les composantes du vecteur propre associé à la plus petite des valeurs propres.
Exercice 4
Tracé du plan des moindres carrés dans la fenêtre graphique
— Écrire une fonction qui renvoie la plus petite valeur et la plus grande valeur d’une liste.
— Ćonstruire 2 listes Xpl et Ypl qui serviront à réaliser un maillage du plan des moindres carrés. Ces listes
devront contenir un ensemble de coordonnées respectivement sur ~x et ~y régulièrement réparties allant de la
plus petite composantes à la plus grande composante des points palpés (on pourra prendre par exemple 30
points).
— Construire le maillage par la commande : Xplan,Yplan=np.meshgrid(Xpl,Ypl)
— Construire une matrice (voir cours 8) des coordonnées Zplan sur ~z0 connaissant les coordonnées de la normale
au plan, et les listes Xplan et Yplan.
— Tracer le plan dans la fenêtre graphique 3D par le code suivant :
ax.plot_wireframe(Xplan,Yplan,Zplan,rstride = 2 , cstride = 2,color='g')
Exercice 5
Détermination de la tolérance de planéité de la surface nominalement plane
— Déterminer le plus grand écart (mesure algébrique ) entre les points palpés et le plan des moindres carrés
— Déterminer le plus petit écart (mesure algébrique ) entre les points palpés et le plan des moindres carrés
— En déduire la plus petite tolérance de forme possible du dessin d’ensemble de la pièce pour que celle-ci soit
considérée comme bonne.
4