De Ensiwiki.

Recalage de données tomographiques RX 3D

Labo	Gipsa-lab
Equipe	sigma-phy
Encadrants	Francois.Cokelaer@gipsa-lab.grenoble-inp.fr jocelyn.chanussot@gipsa-lab.grenoble-inp.fr, Francois.Cokelaer@gipsa-lab.grenoble-inp.fr

Etudiant

Allain, Benjamin

Introduction

Le projet DELPIX

Ce TER a été réalisé au sein du laboratoire Gipsa-lab, dans le cadre du projet DELPIX, réunissant 6 entreprises et 2 laboratoires. DELPIX signifie Détecteur, Equipement et Logiciels pour contrôle de Production avec Imagerie par tomographie en rayons X. La tomographie à rayons X (connue sous le non de « scanner ») permet l’acquisition de volumes 3D. Ce projet doit fournir les outils pour l’utilisation de cette technologie dans l’industrie dans le but de détecter les défauts de fabrication sur les chaines de production.

Rôle du recalage

Le processus de détection de défaut sur un objet se divise en plusieurs étapes successives :

acquisition : acquisition des projections de l’objet selon plusieurs directions.
reconstruction : reconstruction du volume (une grille 3D de voxels) à partir des projections
Débruitage : suppression du bruit (d’origine physique et numérique)
Recalage : repositionnement géométrique de l’objet (image 3D) pour qu’il concorde avec le modèle
Détection de défauts : comparaison entre le modèle et l’objet

Pourquoi introduit-on cette étape de recalage ?
Lorsque l'objet est scanné, sa position dans l'espace ne correspond pas exactement à celle qu'avait le modèle lorsqu'il a été scanné. Il faut donc "recaler" l'image de l'objet sur celle du modèle. Sans ce recalage, la comparaison avec le modèle risque d'être biaisée, car elle portera sur des zones qui ne correspondent pas.

Sujet

Pour déterminer la fonction de recalage (fonction qui replace l'objet à dans la même position que le modèle), on commence par détecter les positions des coins de l'objet. Chaque coin détecté est associé à un coin du modèle (par une méthode encore indéterminée). On dispose donc d'un ensemble d'associations de deux points (couples de points de l'espaces). Le but de ce TER est de trouver une transformation rigide (un déplacement) qui soit cohérente avec ces associations, sachant qu'il y a des associations aberrantes (les "outliers").

La méthode imposée est le couplage de la méthode de K.P. Horn (article [1]) avec la méthode des moindres carrés medians (article [2] et [3]), puis une optimisation standard.

Cette méthode doit être implémentée en C, de manière optimisée. Le code ainsi produit sera ultérieurement intégré dans le logiciel industriel de détection de défauts.

Eléments de pré-requis

On appelle 'déplacement' de l'espace toute application de l'espace dans lui-même qui conserve les distances et l'orientation (des bases). La notion mathématique de déplacement correspond à la notion intuitive : un changement de place. Ainsi la fonction de recalage recherchée est un déplacement.

Tout déplacement se décompose de manière unique comme la composée d'une translation par une rotation vectorielle : $D = t \circ{} r$

Une rotation peut se représenter de diverses manières : un couple axe/angle, une matrice de rotation, un quaternion unitaire ou encore avec le triplet des angles d'Euler.

Soit D un déplacement et (P,Q) un couple de points de l'espace. On appelle résidu quadratique du couple (P,Q) la quantité suivante, qui correspond au carré de l'erreur de recalage pour le couple (P,Q) : ${r_{(P,Q)}(D)}^2 = {distance^2(D(P),Q)}$

On notera $n$ le nombre total de couples (associations) dont on dispose.

Travail réalisé

Méthode de Horn

La méthode de recalage de K.P. Horn est une méthode directe qui donne la similitude directe qui minimise la somme des résidus quadratiques de l'ensemble des couples. Elle minimise donc l'erreur quadratique moyenne de recalage.

Un déplacement est un cas particulier de similitude directe. Comme nos données sont issues d'un déplacement, la solution optimale de méthode de Horn est bien un déplacement.

Expression de la solution : La translation optimale correspond au recalage du barycentre des points de l'objet sur le barycentre des points du modèle. La rotation optimale est donnée s'obtient par le calcul d'un vecteur propre associé à la plus grande valeur propre de la matrice N, où N est une matrice symétrique réelle de taille 4x4, définie à partir des coordonnées des points dans leurs repères barycentriques. Le vecteur propre obtenu, après normalisation, est le quaternion unitaire représentant la rotation optimale. (Remarque : l'expression de N n'a pas sa place ici sans la démonstration)

Cette méthode a été implémentée puis testée.

Inconvénient majeur :

La méthode de Horn permet d’obtenir un résultat correct malgré la présence de points aberrants (les « outliers »). Cependant, la détections de défauts nécessite une grande précision du recalage. On souhaite donc minimiser uniquement les résidus des « inliers » (points non aberrant). Il nous faut donc identifier les inliers. Nous utilisons pour cela la méthode des moindres carrés medians (Least Median Squares (LMedS) en anglais), qui a pour but de d’éliminer les outliers.

Méthode des moindres carrés medians (LMedS)

Principe

Soit $ε$ le taux d'ouliers dans notre ensemble de couples.
Hypothèse : On suppose qu'on a plus de données correctes que de données aberrantes ( $ε < 1 / 2$ ).

Le principe est d'utiliser la méthode de Horn sur plusieurs sous-parties de taille $p$ de l'ensemble des couples. On choisit assez de sous-parties pour être "sûr" qu'il y en a au moins une qui ne contienne que des inliers. Chaque sous-partie $k$ nous donne un déplacement $D k$ (qui est optimal au sens de l'erreur quadratique moyenne sur la sous-partie k). La méthode va exploiter la cohérence des inliers, l'incohérence des outliers et la dominance des inliers pour trouver, à partir des $D k$ et de l'ensemble de couples, une sous-partie $k o p t$ ne contenant que des inliers. De plus, à partir des résidus quadratiques associées à $D_{k_{opt}}$ , pour chaque couple, on va pouvoir décider si c'est un outlier ou un inlier. On sera alors en possession de la liste exhaustive des inliers.

Probabilité de réussite

Les sous-parties sont choisies aléatoirement. Pour $m$ sous-parties de taille $p$ la probabilité de réussite de la méthode est : $P_{reusite} = 1 - {(1 - {(1-\epsilon)}^p)}^m$ En inversant la formule, on obtient le nombre de sous-partie à tirer en fonction de la probabilité désirée, du taux d'ouliers et de la taille des sous-parties.

Fonctionnement

Supposons maintenant qu'au moins une des sous-parties tirées est "bonne" (ie elle ne contient que des inliers). Alors la sous-partie $k o p t$ qui minimise $median_{1\leq i\leq n}\; {r_i{(D_k)^2}}$ ne contient que des inliers.

Intéressons nous à la sous-partie $k$

Si elle ne contient que des inliers, alors $D k$ est optimale, donc les résidus quadratiques des inliers sont quasiments nuls, tandis que les résidus des outliers sont élevés. Comme les inliers sont prédominants, la médiane des résidus quadratiques est le résidu d'un inlier, donc cette médiane est quasiment nulle.
Si elle contient un outlier (ou plus), alors $D k$ est incorrecte, donc les résidus quadratiques des inliers sont élevés, ainsi que ceux des outliers (car les outliers ne sont pas cohérents, même entre eux). La médiane des résidus quadratiques est donc élevée.

Ainsi en minimisant la médiane des résidus quadratiques, on sélectionne une sous-partie constituée uniquement d'inliers.

Sélection des inliers

On a vu précédemment que les résidus quadratiques associés à $D_{k_{opt}}$ sont quasiment nuls pour les inliers, et grands pour les outliers. On peut donc identifier les inliers. Le critère est (cf [2] et [3]) :

${r_i}^2(D_{k_{opt}}) \leq 2.5 {\sigma ^2}_{robuste}$

où ${\sigma ^2}_{robuste} = 1.4826 * (1 + \frac{5}{n-p}) \sqrt{median_i\; r_i^2}$

On dispose donc ici de la liste des inliers.

Optimisation des moindres carrés

Maintenant qu'on a éliminé les outliers, il nous reste à optimiser les résidus des inliers seulement. La solution la plus simple est d'utiliser la méthode de Horn sur l'ensemble des inliers. Une autre solution, en cours de développement par François Cokelaer, est d'utiliser une optimisation des moindres carrés non spécifique au recalage, avec la bibliothèque MINPACK.

Simulation et résultats

Génération de données

L'étape de détection de coins et d'association des coins n'est pas encore implémentée, ni bien définie. Pour tester notre méthode, il faut donc simuler des données. On génère un nuage de $n$ points, on lui applique un déplacement connu D. Puis on ajoute du bruit (additif, blanc, gaussien), et on introduit des outliers par interversions du second terme entre certaines paires de couples. On applique ensuite notre méthode, qui nous renvoie le déplacement estimé, et l'erreur quadratique moyenne.

Résultats

L'estimation fonctionne, même en présence de bruit additif raisonnable.
Les outliers sont quasiment tous détectés. Ceux qui ne le sont pas sont issues d'échanges entre points très proches. Ce ne sont donc pas vraiment des outliers.

Conclusion

L’utilisation du LMedS avec la méthode de Horn, suivi d’une optimisation non linéeaire (méthode de Horn sur les inliers) a été comprise, implémentée et testée avec succès. Il reste à mettre en place l’optimisation non linéeaire ”aveugle” pour voir si elle améliore la solution de la méthode de Horn sur les inliers. Certaines optimisations de code peuvent encore être faites dans la méthod de Horn, et le LMedS peut être parallelisée très facilement. Enfin il faudra vérifier l’efficacité de la méthode et régler ses paramètres avec des données réelles.

Références

Berthold K. P. Horn. Closed-form solution of absolute orientation using unit quaternions. 1986.
Z. Zhang G. Xu. Epipolar Geometry in stereo, motion and object recognition (p102 -105). 1996.
Annick M.Leroy Peter J. Rousseeuw. Robuste regression and outlier detection. 1987.

Documents additionnels

Transparents de la soutenance Media:TER0910_soutenance_AllainBenjamin.pdf
Rapport Media:TER0910_rapport_AllainbBenjamin.pdf

Allain Benjamin Recalage de données tomographiques RX 3D Résultats