1-Objectifs du projet 2-Description de l'algorithme d'appariement manuel 3-Description de l'algorithme d'appariement automatique 4-Résultats
1-Objectifs du projet
Dans le cadre de ce quatrième projet du cours de photographie algorithmique, il nous faut implémenter une application
qui permet de créer un panorama d'images à l'aide de prises de vue se chevauchant. Pour effectuer l'appariement d'image, il faudra d'abord déterminer
des points d'intérêts commun entre les images. Pour ce faire, il faudra utiliser une technique de sélection manuelle de points
et une autre qui sélectionnera des points d'intérêts automatiquement.
Pour plus de détails, voyez l'énoncé complet ICI.
2-Description de l'algorithme d'appariement manuel
Étape 1 - Appariement de caractérisque
Étape 1 - Appariement de caractérisque
Tout d'abord, il faut sélectionner manuellement des points d'intérêts entre 2 images. Habituellement, on utilisera l'image du centre comme image de départ et on s'assurera que les points
d'intérêts sélectionnés dans les autres images (de part et d'autre du panorama) sont présents dans cette photo. Cependant, il est également possible de sélectionner des points d'intérêts
que pour les paires d'images adjacentes et de faire le calcul de l'homographie «en cascade». Il est important de bien déterminer les pixels (point x,y) qui correspondent dans chaques images
puisqu'une simple erreur peu engendrer une grande erreur qui en résulterait d'une déformation peu réaliste de l'image vers le panorama.
Il est important de sélectionner certains nombres de points afin d'obtenir de meilleurs résultats comme nous le verrons dans une prochaine étape. Pour ce faire, une fonction à été créer afin
de faire la sélection de points un-à-un sur deux images et de sauvegarder la sélection sur un fichier texte au nom des deux images combinées. Il faut d'abord savoir combien de points on compte
utiliser.
saveAndDisplayFeatures(imgfilename1, imgfilename2, n)
Étape 2 - Calculer l'homographie
Puis, une fois les points d'intérêts déterminées, il faut déterminer la matrice d'homographie d'images «H» qui permettra la déformation d'une image vers l'autre. Pour ce faire, il faut créer un système
de «n» équations linéaires de la forme Ah = b
où h correspondera à un vecteur des 8 inconnus de la matrice H. Comme la matrice de déformation est de type projective, c'est la raison pour
laquelle on se retrouve avec 8 inconnus (car le 9e élément de la matrice pourra être simplifier et être mis à 1 puisqu'il correspond à un facteur d'échelle). Pour résoudre de le système d'équations, il
faudra donc au minimum 4 points d'intérêts.
Cependant, comme il a déjà été mentionné, il faut sélectionner le plus de points possible, le plus précisément possible entre les images et de s'assurer que ceux-ci ne sont pas tous trop près spatialement.
Des points tous sur la même ligne de l'image ou un amassi de points au même endroit n'est pas tolérable et de prendre un minimum de points également. Dans la figure suivante, on voit le résultat que l'on
obtient pour une sélection de quelques points le long d'une clotûre en bordure d'image.
Mauvais choix de points | Déformation non-désirée | Bon choix de points | Déformation désirée |
Pas très satisfaisant si on s'en tient au minimum !
Il faut donc chercher un homographie en surestimant celle-ci à l'aide de plusieurs points qui seront utilisées pour résoudre plusieurs systèmes d'équations et de résoudre ceux-ci en utilisant la méthode des moindres carrée. Dans Matlab, il suffit donc de créer la matrice d'équations linéaires en utilisant les points d'intérêts de l'image source que l'on utilise sous forme homogène. La matrice de sortie correspond aux points de l'image de destination et il faut donc trouver les inconnus qui permettent cette transformation de façon fidèle et juste. On utilisera l'opérateur
mldivide «\»
et il
faut également savoir qu'il est plus simple de travailler sous forme vectorielle pour définir les équations linéaires et d'utiliser celle-ci sous forme d'une matrice par la suite.
Il nous à donc fallu créer la fonction ;
H = computeH(im1_pts,im2_pts)
qui permet de déterminer l'homographie et donc, pour chaque points, on il faut reconstituer le système d'équations ce qui mènera à la matrice suivante:
|
h = |
|
Comme il à été mentionné précédemment, il est plus simple de toujours calculer l'homographie pour un même plan (image de référence) et de déformer les autres images par rapport à celui-ci. Cependant, une meilleure technique serait de calculer l'homographie pour chaque paires d'images adjacentes dans le panorama et d'établir l'homographie en cascade. Pour ce faire, par exemple si on 5 images et que la 3e est celle du centre, normalement, on choisirait des points d'intérêts pour les paires d'images suivantes (1-3, 2-3, 4-3, 5-3) pour ainsi toujours trouver des homographies en fonction de l'image centrale et donc de déformer vers l'image 3. Cependant il est possible d'établir les correspondances suivante (1-2, 2-3 ,4-3, 5-4) et de déterminer les homographies des images non-adjacentes à l'image du centre par le biais des autres homographies. Ainsi, par exemple, pour l'image 1, on calculerait l'homographie en cascade de la manière suivante H12 * H23 = H13.
Étape 3 - Déformation projective
Résultat de déformation d'image |
La déformation peut aisémemt se calculer à l'aide de la fonction Matlab
interp2
qui sera utilisé sur tout les canaux (RGB) de l'image à l'aide de la fonction ;
imwarped = warpImage(allImages, allMatrixH ,X,Y)
De plus les coordonnées d'interpolation seront donc établie à l'aide de l'homographie en utilisant notre fonction ;
[pts_x2, pts_y2] = useH(matrixH, pts_x, pts_y)
Il est à noter, que pour établir la déformation et surtout la nouvelle image, il faut savoir le format (taille) de l'image résultante qui correspond à la mosaïque d'images. Pour déterminer la format du panorama final, la fonction suivante à été codé ;
[X, Y] = findPanoSize(allImages, allMatrixH)
Cette fonction permet, à l'aide du lot d'images et des homographies, de déterminer l'emplacement des quatres coins de l'image résultante et par le fait même, de connaître sa dimension. En effet, en effectuant un «forward-warping» que sur les coins des images et il est donc possible de déterminer les coins minimaux et maximaux en X et Y quicorrespondent aux quatres coins du panorama.
Image destination (référence) | Image source (à déformer) |
Étape 4 - Fusionner les images en une mosaïque
Pour fusionner les images entre elles, il ne suffit pas toujours de les superposer car on remarque trop l'effet de «collage» de photos. Une manière de contrer cet effet, est d'utiliser une moyenne pondérée sur chacuns
des pixels qui se chevauchent. Une fonction à donc été programmé dans le but d'accomplir cette tâche :
[panoramaFinal] = createPano(imgWarped)
En y insérant les images déformées, cette fonction additionne les pixels qui se chevaucheront et détermine le nombre de chevauchement par pixels afin d'obtenir la moyenne pondérée des images. Il suffit de superposer les images par la
suite.
3-Description de l'algorithme d'appariement automatique
Étape 1 - Détection de caractérisques de coins dans une image
Étape 1 - Détection de caractérisques de coins dans une image
Points d'Harris |
Dans un premier temps, pour trouver les coins caractérisques d'une image, on utilisera l'algorithme d'Harris qui est un classique en matière de détection de coins.
Il es possible de faire la détection à plusieurs échelle bien que dans le but de ce projet ce n'est pas nécessaire que nous utiliseront la fonction existante suivante:
function [newy,newx,newv] = harris(imrgb, nbPtsHarris, robust)
Cette fonction a été modifiée légèrement de son originale puisqu'elle comprend maintenant une section qui permet de choisir le nombre de points à retenir pour faire
l'homographie. De plus, elle implémente maintenant la répression adaptive non-maximale qui permet de faire la sélection le plus judicieusement possible, c'est à dire de manière uniforme et pour des points très distincts. Cette technique
est un peu complexe mais fonctionne très bien. On tri d'abord les coins/points selon leurs «forces» et le maximum global sera notre points d'entrée. En formant un certain rayon autour de ce point qu'on élargit peu à peu, on obtient d'autres
points maximaux qui doivent respecter le critère de robustesse (variable d'entrée de la fonction). On trouve donc le rayon minimal en utilisant la fonction pdist2
et ainsi obtenir le nombre de points désirés avec des «forces» suffisament élevés.
Dans l'image à droite, on peut voir les milliers de points en ROUGE qui correspondent aux points de Harris créer par la fonction, puis en BLEUE, les points sélectionnés en utilisant la technique de répression adaptive non-minimale.
Étape 2 - Extraire un descripteur pour chaque point caractéristique
Une fois les séries de points déterminés dans chacunes des images du panorama, il est important dans déterminer un descripteur. Les descripteurs sont des groupements de pixels de 8x8 pris autour des points caractéristique d'Harris qui décrivent
en quelque sorte, l'identité du pixel dans le but de faire l'appariement. Bien qu'il est favorable d'utiliser l'invariance en rotation de ces descripteurs puisque la prise de photos n'est pas toujours parfaitement «droit» et sans rotation, cette technique
n'a pas été utilisé. Cependant dans le but de rendre plus rapide les calculs, on utilise ces descripteurs sur une plus grande fenêtre de dimension 40x40 soit 5x plus grand et légèrement flou. La fonction suivante permet de faire cette étape :
[descripteurs] = getDescript(image, ptsHarris)
Comme indiquer cette fonction applique d'abord un flou sur l'image à l'aide de imfilter
, puis on transforme un domaine spatiale de afin d'obtenir 5 pixels pour un (8x8 - 40x40) à l'aide de meshgrid
.
Finalement, il suffit d'interpoler les points en utilisant interp2
et surtout, de ne pas oublier de normaliser les points avec le gain et le biais (moyenne (mean)
et variance sqrt(var(...))
des descripteurs).
Étape 3 - Apparier ces descripteurs de caractéristiques entre les deux images
À l'aide de la fonction suivante il est possible d'apparier les différents descripteurs provenant des différents lots de points de Harris déterminés auparavant.
[ptsDest, ptsSource] = matchDescript(descript1,descript2,ptsHarris1,ptsHarris2)
Cette fonction utilise la fonction de distance (pdist2)
pour déterminer la distance entre les paires de descripteurs. Pour l'appariement il suffit d'utiliser le plus proche voisin ou d'utiliser la technique de Lowe qui calcule
le rapport entre les 2 plus proches voisins du descripteur. Voir le fichier suivant pour plus de détails.
Étape 4 - Utiliser une méthode robuste (RANSAC) pour calculer l'homographie
Finalement, après avoir apparier les descripteurs caractérisques les uns avec les autres, il faut estimer l'homographie. Cette tâche peut être très longue si on utilise tous les points de Harris, ou même tous les points d'intérêts sélectionnées. Pour
ce faire, il existe l'algorithme de RANSAC (RANdom SAmple Consensus) qui n'est pas très difficile à implémenter. Par itération (plus le nombre sera grand, plus l'algorithme sera robuste), on prend aléatoirement 4 points pour en déterminer l'homographie
en plus de déterminer le nombre de points qui satisfait l'appariement. La tolérance d'erreur des pixels peut être modifié dans la fonction suivante qui effectue l'algorithme de RANSAC et qui retourne la meilleure homographie d'image :
[matrixH] = findRANSAC(ptsDest, ptsSource)
4-Résultats
Dans le but d'alléger le site web, les panoramas ont été simplifiés
Stade du Rouge et Or |
Grand Axe de l'université |
Linges |
Comme on le remarque sur cette dernière photo, lorsque la lumière varie beaucoup, il est difficile pour l'algorithme de bien détecter les points qui sont semblable et aussi le fait qu'il manque un peu de coins dans l'image n'aide pas beaucoup à avoir un beau panorama final.
Dominos |
Comme on le remarque sur la photo précédente, en présence de forts points d'intérêts/coins et le fait que ceux-ci soient en bons nombres, on trouve une correspondance quasi-parfaite.