Keywords: augmented reality, point pattern matching, registration, texture-less tracking, calibration, projector-camera, augmented map, surface of revolution

Registration is a very important task in Augmented Reality (AR). It provides the spatial alignment between the real environment and virtual objects. Unlike tracking (which relies on previous frame information), wide baseline localization finds the correct solution from a wide search space, so as to overcome the initialization or tracking failure problems. Nowadays, various wide baseline localization methods have been applied successfully. But for objects with no or little texture, there is still no promising method. One possible solution is to rely on the geometric information, which sometimes does not vary as much as texture or colour.

This dissertation focuses on new wide baseline localization methods entirely based on geometric information, and more specifically on points. Point patterns can be identified by using geometric relationships between points. Based on those assumptions, this thesis proposes:

• A novel point pattern matching algorithm (RRDM) to robustly and quickly register two 2D point sets under perspective transformation. It is validated by augmenting different kinds of paper maps of the same city with information retrieved from a Geographic Information System (GIS) with the help of road intersections.

• Another novel and general point pattern matching algorithm (LGC). It registers 2D or 3D point patterns under any known transformation type and supports multi-pattern recognitions. It has a linear behavior with respect to the number of points, which allows for real-time tracking. LGC is applied to multi targets tracking/augmentation, as well as to 3D model registration.

• A practical method for projector-camera system calibration by using random dots calibration patterns: the user only needs to manipulate a small calibration board (B4 size) for 30s to get precise and repeatable calibration results. It can be useful for large scale Spatial Augmented Reality (SAR).

Besides, I also developed a method to estimate the rotation axis of a surface of revolution quickly and precisely on 3D data. It is integrated in a novel framework to reconstruct the surface of revolution on dense SLAM in real-time.

Résumé

La Réalité Augmentée (RA) vise à afficher des informations numériques virtuelles sur des images réelles. Pour la RA, le problème du recalage est d'une importance primordiale, puisqu'il consiste à calculer la position et l'orientation de la caméra filmant la scène réelle afin d'aligner de manière correcte les objets virtuels dans le monde réel. Contrairement au tracking qui résout ce problème en utilisant les informations de l'image précédente, la localisation à grande échelle (wide baseline localization)

calcule la solution en utilisant uniquement les informations présentes dans l'image courante. Il permet ainsi de trouver des solutions initiales au problème de recalage et, n'est pas sujet aux problèmes de perte de tracking qui peuvent survenir lorsque deux images successives sont trop différentes l'une de l'autre. Le problème du recalage en RA est relativement bien étudié dans la littérature, mais les méthodes existantes fonctionnent principalement lorsque la scène augmentée présente des textures. Pourtant, pour les objets peu ou pas texturés, il est possible d'utiliser leurs informations géométriques qui représentent des caractéristiques plus stables que les textures et qui sont indépendantes des changements d'éclairage de la scène.

Cette thèse s'attache donc au problème de recalage basé sur des informations géométriques, et plus précisément sur les points. Nous proposons deux nouvelles méthodes de recalage de points (RRDM et LGC) robustes et rapides. LGC est une amélioration de la méthode RRDM et peut mettre en correspondance des ensembles de motifs de points 2D ou 3D subissant une transformation dont le type est connu. LGC présente un comportement linéaire en fonction du nombre de points, ce qui permet un tracking en temps-réel. La pertinence de LGC a été illustrée en développant une application de calibration de système projecteur-caméra dont les résultats sont comparables avec l'état de l'art tout en présentant des avantages pour l'utilisateur en termes de taille de mire de calibration, i.e. une mire de taille B4 seulement est nécessaire quelle que soit le distance de mise au point du projecteur lorsque l'état de l'art utilise des mires de plus de 1.5 × 2m.

Mots-clés : réalité augmenté, recalage des motifs de points, suivi, peu texturé, calibration, projecteur-caméra, carte augmentée, surface de revolution

Résumé substantiel

Chapitre 1 : Introduction La Réalité Augmentée (RA) vise à afficher des informations numériques virtuelles sur des images réelles. Pour la RA, le problème du recalage est d'une importance primordiale, puisqu'il consiste à calculer la position et l'orientation de la caméra filmant la scène réelle afin d'aligner de manière correcte les objets virtuels dans le monde réel.

Ce problème de recalage en 3D pour être résolu de plusieurs manières différentes en fonction des dispositifs réalisant la capture de la scène réelle. Pour la RA, la capture de la scène se fait le plus souvent via des caméras et dans ce cas, la solution au problème de recalage 3D revient à déterminer la pose (i.e. la position et l'orientation) de la caméra dans le monde réel à partir des images capturées par cette dernière.

Il y a deux grandes familles qui regroupent les différentes méthodes d'estimation de pose : les méthodes dites Narrow Baseline Tracking (NBT) et les méthodes dites Wide Baseline Localization (WBL). Les algorithmes NBT utilisent la pose calculée dans l'image précédente et fait l'hypothèse de petits mouvements entre deux images consécutives afin de calculer la pose de l'image courante. Ces méthodes ne sont donc pas utilisables quand la pose de l'image précédente n'existe pas (e.g. initialisation) ou n'a pas été trouvée. A l'inverse, les méthodes WBL ne connaissent pas de telles restrictions. Plusieurs méthodes de WBL existent, par exemple celles basées sur des marqueurs artificiels ou des textures. Pourtant, ces méthodes présentes également des limitations et ne peuvent pas marcher quand les marqueurs ne sont pas acceptés et en même temps il n'y pas assez de textures.

L'idée centrale de notre thèse provient du problème dit des cartes augmentées, où l'on souhaite rajouter des informations numériques (issues par exemple d'un Système d'Information Géographique -SIG) sur des cartes en papier non préparées (c'est à dire pas de marqueur). L'idée est de pouvoir augmenter plusieurs cartes dont les couleurs sont différentes de la même ville, à partir du même modèle SIG enregistré.

Parce que les couleurs des cartes sont différentes, on ne peut que compter sur des caractéristiques géométriques pour les recaler. Les méthodes existantes basées sur les textures ne peuvent pas résoudre ce problème sans avoir à apprendre au préalable toutes les cartes.

La distribution des intersections des routes d'une ville est, dans le cas général en europe, assez caractéristique de cette dernière. Pour résoudre le problème des cartes augmentées, il est théoriquement possible d'effectuer un recalage en se basant sur cet ensemble de points. De plus, le point est la géométrie la plus basique et peut donc se trouver partout. Par conséquent, le recalage par mise en correspondance d'ensembles de points peut potentiellement servir à la wide baseline localization dans beaucoup de situations.

Chapitre 2 : L'état de l'art Les méthodes wide baseline localization peuvent être séparées en trois : les approches basées sur des marqueurs artificiels, celles basées sur des textures et celles basées sur des géométries.

Le premier système de suivi de marqueur artificiel en RA a été proposé il y a presque 20 ans [Rekimoto 1998]. Ces marqueurs sont souvent planaires, en noir et blanc, et ils contiennent des motifs prédéfinis, distinguables les uns des autres. Par conséquent, ils sont très faciles à détecter dans des images capturées par caméras.

Une fois ces marqueurs détectés, la pose de la caméra peut être calculée de manière simple et efficace. Quoique cette approche soit rapide et robuste, la présence de ces marqueurs un peu partout dans la scène est nécessaire, ce qui est envahissant.

De plus, ces méthodes résistent mal à l'occlusion des marqueurs : si une partie du marqueur est cachée alors la détection échoue. Enfin, les marqueurs artificiels doivent occuper une certaine taille dans l'image ce qui rend leur utilisation difficile en extérieur.

Les approches basées textures utilisent les points caractéristiques que l'on appelle aussi points d'intérêts. Lors d'un prétraitement hors-ligne, les points d'intérêts sont extraits dans un modèle (généralement, une image). Un descripteur de grande dimension est créé à partir des textures locales autour de chaque point d'intérêt.

Ensuite, lors d'une phase en ligne, des points d'intérêts sont détectés dans l'image courante capturée depuis la caméra, puis sont mis en correspondance avec les points d'intérêts extraits de modèle en fonction de la distance entre leurs descripteurs.

Grâce aux correspondances, la pose de caméra peut être trouvée. Bien que ces approches soient implémentées dans beaucoup de SDKs commerciaux, elles ne peuvent résoudre les problèmes similaires à celui des cartes augmentées présenté dans le Chapitre 1.

Les approches basées géométries sont moins matures que les deux autres approches. Peu de méthodes temps-réel existent à cause de la complexité du problème v à résoudre. Ces méthodes essayent d'identifier des géométries dans l'image courante et de les mettre en correspondance avec la géométrie d'un modèle connu que l'on souhaite détecter dans l'image courante. Les géométries les plus utilisées sont les régions, les lignes et les points. Les méthodes basées sur les points ont attiré le plus d'attention car le point est la géométrie la plus répandue et donc la plus facile à détecter. La méthode dite des Random Dot Markers (RDM) [Uchiyama 2011b], permet de recaler des ensembles de points aléatoirement distribués en temps-réel. Cette technique a notamment été utilisée par [Uchiyama 2011c] pour réaliser une carte augmentée à partir d'un SIG et en utilisant les intersections routières d'une carte détectées manuellement sur la carte. Toutefois, cette méthode échoue lorsque l'on utilise des intersections de routes détectées automatiquement par un programme.

Cet échec est dû à la faible robustesse de la méthode quand il y a du bruit dans la détection. Ce bruit va en effet ajouter ou retirer des points à l'ensemble, ainsi qu'introduire des légères différences aux coordonnées (ce que l'on appelle dans la suite le jitter) de chaque points détectés dans l'image courante de la caméra. Ceci va faire échouer la mise en correspondance de cet ensemble de points image avec l'ensemble des intersections issues du SIG. L'un des buts de notre thèse est de résoudre ce problème en proposant un algorithme à la fois efficace et robuste qui permettra de prendre en compte ce bruit issu de la détection automatique des intersections de routes.

Chapitre 3 : Marqueurs consistués des points pseudo aléatoires pour des recalages robustes L'algorithme "Robust Random Dot Markers" (RRDM) est proposé afin de résoudre le recalage de deux ensembles de points 2D sous une transformation perspective de manière robuste et efficace. Par exemple, si on enregistre un ensemble de points modèle 𝑃 dans RRDM, et on lui donne un autre ensemble de points 𝑄 extrait de la scène (i.e. de l'image courante filmée par une caméra), RRDM va calculer une transformation perspective permettant de transformer l'ensemble de points 𝑃 en (approximativement) Q.

Dans une première phase, RRDM sépare les deux ensembles de points originaux 𝑃 et 𝑄 à mettre en correspondance (i.e. à recaler) en sous-ensembles locaux. Ceci a pour objectif d'éviter une complexité trop grande lors de la mise en correspondance due à l'explosion combinatoire. Un descripteur robuste s'appelant ratio de deux surfaces (pour Two Surface Ratio -TSR) est proposé et permet de décrire les positions relatives entre les points dans un sous-ensemble. Le TSR reste stable vi même s'il y a du jitter et des points ajoutés ou supprimés aux ensembles. Une fois une mesure de similarité calculée grâce au descripteur TSR, un processus de vote local est réalisé pour chaque couple de sous-ensembles, où l'un des membres du couple est issu de 𝑃 et l'autre de Q. Ce processus de vote a pour le but de mettre en correspondance des points dans les deux sous-ensembles. Une transformation géométrique locale peut alors être estimé à partir de ces correspondances.

RRDM propose pour la première fois d'utiliser un consensus sur des transformations géométriques locales afin d'éliminer les mauvaises correspondances. Ceci nous permet d'identifier les correspondances correctes même si leur nombre est très faible par rapport aux mauvaises.

Des études menées sur des points générés aléatoirement par ordinateur montrent que RRDM est plus robuste que RDM dans les cas où il y a du jitter, lorsque des points sont ajoutés ou perdus, ainsi que lorsque la distorsion perspective est grande.

En termes de complexité algorithmique, RRDM a un comportement quadratique par rapport au nombre de points. Pour le problème de carte augmentée, RRDM arrive à recaler les intersections de routes issues d'un SIG et celles détectées automatiquement sur les cartes. Il prouve ainsi que le recalage de points peut être suffisamment robuste et efficace pour effectuer un recalage purement géométrique.

Bien que RRDM résolve le recalage de deux ensembles de points de manière robuste, il présente des défauts non négligeables. En effet, RRDM ne peut pas gérer les cas où il existe plusieurs modèles différents (par exemple rechercher parmi plusieurs ensembles représentant les intersections de routes de plusieurs villes). De part son comportement quadratique, il n'est pas assez rapide pour une application en temps-réel lorsque le nombre de points est plus grand que 100.

Chapitre 4 : Consensus géométrique local L'algorithme "Consensus Géométrique Local" (LGC) est une suite de la méthode RRDM et peut mettre en correspondance des ensembles de motifs de points 2D ou 3D subissant une transformation dont le type est connu (par exemple, perspective, similitude, etc.). LGC possède un comportement linéaire en fonction du nombre de points, ce qui permet un tracking en temps-réel, même lorsque les ensembles à mettre en correspondance contiennent un nombre élevé de points.

Tout comme RRDM, LGC découpe tout d'abord les ensembles à mettre en correspondance en sous-ensembles locaux. Rappelons que LGC supporte plusieurs ensembles de points modèles que l'on va essayer de détecter dans les images issues de la caméra. L'algorithme LGC contient trois modules : un générateur d'hypothèses, vii un validateur d'hypothèses et un raffineur. Étant donnée un sous-ensemble local aléatoirement choisi à partir de l'ensemble des points détecté dans l'image de la caméra, le générateur d'hypothèses estime tout d'abord des correspondances présumées et une transformation locale présumée pour chaque sous-ensemble de points modèles. Ces correspondances et cette transformation locale résumées forment les hypothèse. Pour chaque hypothèse, le validateur applique la transformation présumée aux sous-ensembles de points présents dans le voisinage de cette hypothèse.

S'il y a assez de sous-ensembles en consensus avec cette transformation présumée, l'hypothèse est alors considérée comme correcte, c'est-à-dire, elle contient des correspondances correctes. Dans une dernière étape, un raffineur est appliqué avec pour objectif de trouver le plus de correspondances correctes possibles pour estimer une transformation précise.

Des études menées sur des ensembles points générés aléatoirement montrent que

LGC est plus robuste que RDM tout en étant aussi efficace, et qu'il est plus efficace que le Hashing Géométrique classique et que RRDM tout en étant aussi voire plus robuste. Nous montrons également que LGC est capable de détecter et de suivre divers objets avec ou sans texture en temps-réel, en augmentant par exemple des croquis d'ingénierie avec leurs modèles 3D.

Chapitre 5 : Calibration d'un système de projecteur-caméra En combinant une caméra et un projecteur reliés de manière rigide ou non, un système projecteur-caméra (ProCam) permet de projeter des informations numériques directement sur des surfaces ou des objets réels (Spatial Augmented Reality -SAR).

La calibration de système ProCam consiste à déterminer les paramètres intrinsèques du projecteur et de la caméra, ainsi que leurs positions relatives (i.e. paramètres extrinsèques). Ces paramètres sont essentiels pour que les informations numériques soient bien projetées là ou elles le doivent sur les objets réels.

Pour les applications avec une grande zone de projection, le projecteur est souvent mis au point à une distance assez grande. Par conséquent, les projections à courte distance deviennent floues car elles sont faites dans une zone de l'espace pour laquelle le projecteur est hors focus. Les méthodes existantes pour calibrer les systèmes ProCam dans de telles situations ont besoin soit qu'un utilisateur manipule une mire de calibration de traille très importante (e.g. 1.5 × 2m, soit d'utiliser une mire 3D de calibration dont le modèle 3D numérique est connu de manière précise, soit d'attendre quelques minutes pour obtenir les résultats.

Nous proposons une application basée sur l'utilisation de notre algorithme de viii

LGC qui permet de calibrer un système de ProCam de manière pratique et précise et ce quelle que soit la distance de mise au point du projecteur. L'utilisateur ne manipule qu'une mire de taille B4 (soit 250 × 353mm) pendant 30s pour obtenir un résultat aussi stable et précis qu'une méthode représentative de l'état de l'art, alors que cette dernière a besoin d'une mire de calibration 22 fois plus grande que la notre. Des expériences montrent que pour une projection a une distance de 4.5m, l'erreur de re-projection est de 4mm environ.

Chapitre 6 : Reconstruction de surface de révolution La démocratisation des dispositifs permettant à un utilisateur d'acquérir rapidement et simplement un modèle 3D de l'environnement (i.e. scanners 3D, Kinect) a provoqué un regain d'intérêt pour les approches d'extraction automatique de connaissance dans une scène 3D.

La compréhension d'une scène basée sur l'utilisation de données 3D est importante car elle possède de nombreuses applications comme par exemple la préhension robotique, la simplification de modèle et la réalité augmentée. Le but de l'extraction automatique de connaissance dans une scène 3D est dans un premier temps de détecter l'ensemble des formes primitives (par exemple, sphères, cylindres, etc.) qui composent la scène 3D avant d'ensuite tenter d'établir des relations entre elles.

Les formes primitives actuellement détectées automatiquement par les méthodes existantes sont souvent assez limitées, par exemple, planes, sphères, cylindres, etc.

Dans ce chapitre, nous nous intéressons à la détection automatique dans une scène 3D des surfaces de révolution (SoR), puisqu'elles englobent une grande variété de formes primitives. L'estimation de l'axe de rotation est une étape essentielle dans la détection d'une surface de révolution. Les méthodes existantes permettant de détecter les SoR dans une scène 3D ne sont pas satisfaisantes car elles ne sont pas assez précises ou bien ne sont pas suffisamment efficaces. Ce chapitre propose une méthode tirant profit de la symétrie des surfaces de révolution, afin de trouver une formulation simple permettant d'estimer l'axe de rotation de manière précise et efficace. Nous illustrons la pertinence de notre méthode en reconstruisant des surfaces de révolution en temps réel dans des données issues d'un SLAM dense.

Nos résultats sur des troncs de cônes synthétiques montrent que les erreurs d'estimations sur les orientations de l'axe de rotation sont inférieures à 0.5 ∘ et que le processus peut reconstruire différents objets usuels (bouteilles, tasses, mugs, etc.) issus de surfaces de rotations en temps-réel (45fps). Des expériences menées sur des surfaces de révolutions réelles (deux cylindres) montrent que les erreurs d'estimations ix des diamètres sont inférieures à 2mm.

Chapitre 7 : Conclusion

Dans cette thèse nous avons proposé des méthodes pour recaler des motifs de points pseudo aléatoires, puis nous avons illustré la pertinence de nos approches en développant des applications de réalité augmentée et de calibration de système ProCam basées sur nos algorithmes de recalage. Nous avons montré que le recalage des motifs de points est une solution viable pour la localisation à grande échelle dans des applications de réalité augmentée. Les contributions de cette thèse sont :

• Une nouvelle méthode de recalage de points robuste et rapide (RRDM) permettant de recaler deux motifs de points en 2D sous une transformation perspective. Nous avons montré que cette méthode permet d'effectuer le recalage des données issues d'un SIG et de les afficher en réalité augmentée sur des cartes en papier dont les couleurs sont différentes et ce sans préparation préalable.

• LGC, une deuxième méthode de recalage de points plus générale que RRDM.

L'algorithme LGC permet de mettre en correspondance des ensembles de motifs de points 2D ou 3D subissant une transformation dont uniquement le type est connu. LGC présente un comportement linéaire en fonction du nombre de points, ce qui permet un tracking en temps-réel et le rend donc particulièrement intéressant pour des applications de réalité augmentée.

• Une application de LGC à la calibration de systèmes projecteur-caméra dont les résultats sont comparables avec l'état de l'art tout en présentant des avantages pour l'utilisateur en termes de manipulation de mire de calibration qui est 22 fois plus petite que celle nécessaire avec une méthode existante représentative lorsque la projection s'effectue à de grandes distances (≥ 2.5m).

Dans le futur, diverses améliorations sur des applications développées seront envisagées pour les rendre plus faciles à utiliser.

Il est par exemple possible d'optimiser les performances de LGC pour des applications différentes en choisissant de meilleures valeurs de paramètres. Pour ce faire, il sera nécessaire d'étudier attentivement les influences de chaque paramètre sur le comportement global de l'algorithme. De nouvelles applications font aussi partie de nos perspectives de recherche. Nous pensons par exemple à un SLAM monoculaire pour des environnements peu texturés ou bien à une méthode pour recaler des modèles CAO avec ses modèles issus des scans 3D pour lesquels les approches existantes ne sont pas satisfaisantes aujourd'hui. Enfin, afin de proposer une solution plus générale au problème The reality-virtuality continuum [START_REF] Milgram | [END_REF] bead coding [START_REF] Datta | [END_REF]], (b) engineering drawing of an apartment. . 11

2.1

Examples of fiducial markers: (a) ARTag [Fiala 2005], (b) AR-ToolKit [START_REF] Kato | [END_REF]], (c) Circular markers [Naimark 2002

5.5

Calibration Patterns: 𝑃 𝑏 in black are printed on a piece of paper. 𝑃 𝑝 in white are projected. Both form the original pattern 𝑃 𝑜 . . . . . . .

5.6

Point patterns used in Alg. 9. It allows people to perceive some information which is hidden or does not exist in the real-world. One of the most commonly accepted definition of AR is given by Ronald Azuma [Azuma 1997, Billinghurst 2015]: an AR system should have three key properties:

1. It combines real and virtual content

2. It is interactive in real-time 3. It is registered in 3D
Azuma's definition is quite general, not specific to visual information. However, since the vision is usually the most important perception system for humans, most AR research focus on augmenting the real-world environment using visual information. Visual AR is also the focus of this dissertation.

AR is closely related to Virtual Reality (VR), both of which create virtual contents to enhance users' perception. The reality-virtuality continuum [START_REF] Milgram | [END_REF],

Chapter 1. Introduction as shown in Fig. 1.1, clearly illustrates the difference between AR and VR. The continuum extends from Real Environment (pure real-world environment) to Virtual Environment (VR). In VR, people are totally immersed into a virtually created world and no real-world object can be perceived. Mixed Reality (MR) is defined as the mixture between Real Environment and Virtual Environment. AR is a part of MR, where only a part of the Real Environment is replaced by some virtual contents.

In AR, people can interact with both real and virtual objects.

Figure 1.1: The reality-virtuality continuum [START_REF] Milgram | [END_REF] According to Azuma's definition, a basic system of visual AR contains three key components: (1) a 3D registration system, which aligns the virtual environment to the real environment;

(2) a rendering system, which combines the virtual environment and the real environment, taking account of occlusion due to different depths, color changes due to lighting conditions, etc.;

(3) a display, which shows the rendering result to end users. All the three components should run in real-time (normally ≥ 20 frames per second -fps) so that users can notice less artifacts.

Depending on different technologies used to display the augmented contents, AR displays can be classified into three types [Billinghurst 2015]:

1. Video See-Through (VST): VST displays use a camera to capture the real environment. They combine the virtual environment and the real environment numerically and show the final augmented results to users on screens. Users cannot see the real-world directly but only a digital illustration captured by the camera and displayed on screens.

2. Optical See-Through (OST): OST displays use optical systems to combine light from the real environment and light from virtual contents. Users' eyes receive the combined light, so he can see the augmented results.

Registration

From Azuma's definition, another key feature of AR systems is the 3D registration.

It is a common requirement whatever rendering techniques or display devices are used. 3D registration aligns the virtual and the real environments so that virtual objects can be merged into the real-world at their desired positions and orientations.

Although 3D registration exists for some VR systems as well, AR systems require much more accurate registration than VR systems. Since virtual objects are overlaid onto the real environment, even a 0.5 degree of arc mis-registration can easily be noticed by users [Azuma 1997]. 1. Magnetic-based localization relies on the properties of magnetic fields to calculate the pose of the receiver, which is attached to the viewer. with 𝑂 𝑤 being the coordinate system of the real environment. 𝑂 𝑐 is the coordinate system of the camera. Assume that the pinhole camera has zero lens distortion.

Light from a world point x traverses through the camera's optical center and forms an image x 𝑐* on the camera's actual image plane. For the sake of simplicity, the image point x 𝑐 on the camera's equivalent image plane is often used instead. The relationship between x = (𝑥, 𝑦, 𝑧) and x 𝑐 = (𝑥 𝑐 , 𝑦 𝑐 ) in homogeneous coordinates can be expressed as:

Chapter 1. Introduction ⎛ ⎜ ⎜ ⎝ 𝑥 𝑐 𝑦 𝑐 1 ⎞ ⎟ ⎟ ⎠ = K [R|t] ⎛ ⎜ ⎜ ⎜ ⎜ ⎝ 𝑥 𝑦 𝑧 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎠ (1.1)
In (1.1), K is the camera's intrinsic matrix. It is assumed to be invariant as time and position change. [R|t] is the extrinsic matrix, where R is a rotation matrix and t is a vector representing the translation. They describe the transformation from the world coordinate system O 𝑤 to the camera coordinate system O 𝑐 . As a consequence, performing registration for AR systems with a vision-based approach boils down to finding the camera extrinsic parameter matrix [R|t].

In images captured by real cameras, the lens distortion is unavoidable. If x𝑐 = (x 𝑐 , ỹ𝑐 ) is the image point distorted by the lens distortions, the distortion-free image point x 𝑐 = (𝑥 𝑐 , 𝑦 𝑐 ) can be usually derived by solving equations below: [Zhang 2014[START_REF][END_REF]]:

x𝑐 = 𝑥 𝑐 (1 + 𝑘 1 𝑟 2 + 𝑘 2 𝑟 4 + 𝑘 3 𝑟 6 ) + 2𝑝 1 𝑥 𝑐 𝑦 𝑐 + 𝑝 2 [𝑟 2 + 2(𝑥 𝑐 ) 2 ] ỹ𝑐 = 𝑦 𝑐 (1 + 𝑘 1 𝑟 2 + 𝑘 2 𝑟 4 + 𝑘 3 𝑟 6 ) + 𝑝 1 [𝑟 2 + 2(𝑦 𝑐 ) 2 ] + 2𝑝 2 𝑥 𝑐 𝑦 𝑐 (1.2)
where 𝑟 2 = (𝑥 𝑐 ) 2 + (𝑦 𝑐 ) 2 , 𝑘 1 , 𝑘 2 , 𝑘 3 are radial distortion coefficients and 𝑝 1 , 𝑝 2 are tangential distortion coefficients.

Fortunately, since intrinsic matrix and the distortion coefficients can be assumed invariant for a camera as time and position change, cameras in AR systems are often calibrated beforehand. Therefore, localization is to estimate the extrinsic matrix with known the intrinsic matrix K and lens distortion coefficients.

Most of the time, an AR system seamlessly estimates the location of the viewer (i.e. the camera) since the viewer moves with respect to the environment. It consists of estimating the extrinsic matrix [R|t] for each frame in the video sequence captured by the camera. If the viewer's pose in the previous frame is known, it can be used as additional information for estimating the pose of the current frame. The usage of the previous pose to estimate current pose is called narrow baseline tracking (also called recursive tracking [Lepetit 2005] or tracking by tracking [Uchiyama 2012]).

For the first frame of the input video sequence, or if the pose of previous frame is not found (i.e. the tracking failure), the viewer's pose can only be estimated from information in current frame. In such a case, it is called wide baseline localization (also called tracking by detection [Lepetit 2005, Uchiyama 2012]).

1.1. Registration

Narrow baseline tracking

Since the previous camera pose is available, and the current pose is assumed to be changed only little from the previous one, narrow baseline tracking is mostly tackled by optimization, where the previous pose serves as the initial value. Two approaches are commonly used:

1. Bayesian tracking recursively estimates the probability density function of camera pose over time using incoming images captured by the camera and a mathematical process model. The two most important applications are Kalman [Gennery 1992] and Particle [START_REF] Choi | [END_REF]] filters.

2. Local optimization calculates the camera pose [R|t] by minimizing the disparity between the model and the observed data, with the camera pose in previous frame being the initial value. For example, optical flow or template matching minimizes the difference between a template warped by a transformation 𝑇 ([R|t]) and the input image [START_REF] Jurie | [END_REF]]. Visual Servo System (VSS) tackles the dual problem of 2D visual servoing and minimizes the error between the command and the output of a control system [Comport 2006].

Wide baseline localization

Since no presumed position is available, wide baseline localization consists of finding the camera pose by solving (1.1) in most of the time. As the intrinsic matrix K and lens distortion coefficients are known, if there are enough world-image points correspondences, (1.1) can be solved and [R|t] can be uniquely defined. Depending on the nature of world points' spatial arrangement, there are two main approaches:

1. Perspective-n-Point (PnP) is used when world points are not coplanar.

When 𝑛 ≥ 4, [R|t] has a unique solution [Gruen 2002].

2. Pose from 3D plane is used when world points are coplanar. The transformation from the plane containing all world points to the image plane can be represented by a homography 𝐻. It can be uniquely determined when 𝑛 ≥ 4. When both 𝐻 and K are known, [R|t] can be uniquely determined [Lepetit 2005].

Compared to the "local optimization" in Section 1.1.2, these two approaches aim at finding the global optimal solution. Whichever approach used, the most crucial step in wide baseline localization is to find the correct correspondences between Chapter 1. Introduction world points and image points. This is often done by matching points from a previously known model object and its image captured by the camera, i.e. points from the scene object. In most applications, the number of points to match is much greater than 4. Additionally, there is always noise in image points due to imperfect detections. These constraints make this crucial step difficult to solve. Some other work treat the wide baseline matching as pose classification plus narrow baseline tracking. Programs learn roughly the pose of the model under different views offline. During online matching, the learnt pose being most similar to the scene is found via classification. Then, the pose is used as the initial value for narrow baseline tracking [START_REF] Holzer | [END_REF]]. This kind of approaches does not need to know point correspondences explicitly.

Thesis motivation

Although narrow baseline tracking demands less computation efforts and can be relatively easy to solve thanks to information from the previous frame, it suffers from some drawbacks, such as its difficulty to recover from tracking failure and to deal with fast movements (i.e. large movements between two consecutive frames).

Moreover, it cannot estimate the initial pose. Wide baseline localization comes with a more expensive computation cost but it does not suffer from any of above drawbacks.

Finding point correspondences

The key of wide baseline localization is to find correspondences between world points and image points. Currently, there are commonly two kinds of approaches to find these correspondences, i.e. marker based approaches and natural feature based approaches.

Fiducial markers, such as ARToolKit Markers [START_REF] Kato | [END_REF]], contain predefined self-distinguishable patterns and are attached to target objects at specific locations.

During tracking, they are recognized by the AR system in each frame, thus leading to points correspondences. They can achieve frame-rate wide baseline localization, but they require special preparation on each target, and fiducial markers are intrusive for users.

Nature feature methods are proposed for marker-less wide baseline localization, such as SIFT [START_REF] Lowe | [END_REF]], SURF [START_REF] Bay | [END_REF] or BRIEF [Calonder 2012]. These methods rely on local textures to distinguish different interest points to establish 1.2. Thesis motivation point correspondences. Recently, they become a standard approach for wide baseline localization in many cases.

Problem risen from augmented maps

The motivation of this dissertation comes from the idea of Augmented Maps. Traditional paper map offers lots of advantages: users can have natural interactions with them such as adding annotations or drawing trajectories. It is also easier for several users to collaborate on a common paper map which facilitates communication.

However, paper maps remain static: their content cannot be changed according to users' preferences, no computation result (such as flow simulation or shortest path) can be added without adding layers or reprinting the map. These drawbacks can be removed if paper maps were augmented with dynamic digital information. Since most information of cities is now recorded in Geographic Information Systems (GIS), it is a natural way to base augmented maps on GIS. Relying on existing methods, one possible approach is to associate natural features of target paper map with their geo-coordinates manually, then the map can be tracked by traditional natural feature methods. But manual registrations are inconvenient. Especially when working with different types of maps of the same city (cf. Fig. 5.3), manual registrations should be done for each type of map. Therefore, an automatic geo-referencing of paper map with GIS is preferable.

Contents on paper maps and information stored in GIS are very different. For example, Fig. 1.5a-1.5b show two paper maps of different styles of Edinburgh, which are raster images. Fig. 1.5c is an illustration of the road network of the same area retrieved from a GIS, where raw data are vectorial (e.g. lines, points, polygons, etc.). Natural feature methods are not able to establish correspondences between these two sources, since the nature of data is completely different.

In fact, the only common features between a paper map and the data from GIS of the same area are geometric quantities, such as road networks. Uchiyama et al. have successfully registered paper maps with a GIS by using a point pattern matching algorithm (Random Dot Markers -RDM [Uchiyama 2011b]) on road intersections [Uchiyama 2011c]. However, they need to explicitly mark the road intersections manually so that they can be easily detected. In order to augment natural, unprepared maps, it is possible to use algorithms, such as [Callier 2011], to automatically extract road intersections instead of performing a manual detection.

Unfortunately, such algorithms are prone to massive under and over detection as well as approximation errors of positions of the road intersections. As a result, it is impossible for the RDM method to compute a correct matching with these automatically detected road intersections.

Thus, automatic geo-referencing of paper map with GIS still remains unsolved.

More generally speaking, the problem of wide baseline localization when only geometric quantities are available, without any marker or texture, is still not fully solved.

Point pattern matching and its difficulties

Point is the most basic geometric element. In some cases, only points are available, such as a starry sky (cf. Fig. 1.6a) and bead coding used in massive simultaneous assays in drug screening and drug discovery (cf. Fig. 1.6b). In some cases, objects may not have much texture information but it is easy to retrieve lots of points from them, such as engineering drawings (cf. Fig. This inaccuracy of point position is called jitter afterwards.

The difficulty of PPM is that every point is the same. When considering individual points rather than a group of points, any two points can be potentially matched. This gives a huge number of potential matching candidates. The available information to distinguish different points is only the geometrical distribution of points, i.e. the spatial relationship between points. However, due to the presence of under-detection, over-detection and jitter, the only usable information becomes biased. Given the poorness of the quality of the information, robust matching algorithms are generally time consuming. This is the reason why PPM algorithms in the literature are either fast but not robust against point jitter [START_REF] Nakai | [END_REF]], or robust but not fast enough [START_REF] Wolfson | [END_REF]] for tracking purposes in complex AR applications.

Thesis statement and contribution

Thesis statement

Wide baseline point pattern matching can be both fast and robust enough for augmented reality applications in texture-less environment.

The main contributions of this dissertation are: In the past, two surveys on AR have been published in 1997 [Azuma 1997] and 2001 [Azuma 2001]. An analysis of 10 years of ISMAR publications by [START_REF] Zhou | [END_REF] has shown the important position of registration, i.e. more than 20% of the total papers are related to "tracking". Recently, a survey on pose estimation for augmented reality was published [START_REF] Marchand | [END_REF]].

This chapter describes previous work on wide baseline registration methods.

It begins with a short summary on fiducial marker based registration methods in Section 2.1, followed by marker-less registration methods based on local textures in Section 2.2. It ends up drawing attention on marker-less registration methods based on geometric information in Section 2.3, which is not totally solved yet but hardly mentioned in latest survey [START_REF] Marchand | [END_REF]].

Fiducial marker based localization

Fiducial markers are man-made objects, often planar, placed in the field of view of the camera to serve as reference points. Most of these markers contain selfdistinguishable patterns, so point correspondences are easy to establish.

Chapter 2. Related work Such markers were firstly introduced by [Rekimoto 1998] in AR, followed by many other types of designs. The two most representative types of fiducial markers are square and circular markers. Square markers use four corners of a rectangle as points of interest to estimate a perspective transformation. As a single marker is sufficient to compute the camera pose, they are also denoted as "planar" markers [Lepetit 2005]. While ARTag [Fiala 2005] also use matrix codes as in [Rekimoto 1998] for marker identification, other patterns such as alphabet letters and simple geometric shapes are also used in ARToolKit [START_REF] Kato | [END_REF]]. Regarding circular fiducial markers, at least four of them need to be positioned in a coplanar but non collinear way to compute a unique perspective transformation, thus they are referred to as "point" markers [Lepetit 2005]. Contrasting concentric circles (CCC) [Gatrell 1992] is the main trend of this type of markers for which sub-pixel accuracy can be achieved. [Naimark 2002] adds "data rings" between outer and inner rings to identify different markers.

[Uchiyama 2011b] proposed Random Dot Markers (RDM). They use local arrangements of point patterns, known as Locally Likely Arrangement Hashing (LLAH) [START_REF] Nakai | [END_REF]], to establish point correspondences. Their method proves to be robust to partial occlusions. As this method also can be regarded as geometry based registration, it will be detailed in Section 2. There are three basic tasks in the whole procedure above: feature extraction is dedicated to localizing visual salient features; feature description is used to construct the high dimension descriptor from the texture of the local patch around each feature; feature matching is used to establish correspondences according to descriptors.

Feature extraction

Interest points include geometric corners as well as the points in an image where the image intensity changes significantly. Usually, we refer to both of them as corners [START_REF] Rosten | [END_REF]] in a general sense. Curvature maxima can be extracted to detect geometric corners from the contours of objects [START_REF] Awrangjeb | [END_REF]]. Secondorder derivatives or auto-correlation are used to extract general corners, such as Difference of Gaussian (DoG) [START_REF] David G Lowe | Object recognition from local scale-invariant features[END_REF]], Harris corner detector [START_REF] Harris | [END_REF]], Shi-Tomasi detector [Shi 1994]. SURF uses the determinant of the Hessian matrix and relies on integral images for speeding up [START_REF] Bay | [END_REF]]. Pixel intensity comparison is also used to detect corners more efficiently (FAST) [Rosten 2006] Thick line and interest regions can be used as features as well. [START_REF] Grompone | [END_REF] (LSD), [Wang 2009a] (Line Signature), [Zhang 2013] offer some methods for line segment extraction. [Donoser 2006] gives a good method for extracting interest regions (MSER).

Chapter 2. Related work

Regarding the wide baseline localization problem, feature repeatability is an important challenge for the extraction method. A feature is said to be repeatable when the same feature is extracted from different images of the same scene, despite perspective distortions or imaging noise. An evaluation of repeatability of different methods on point feature extraction can be found in [Ballesta 2008], which concludes that Harris corner detector [START_REF] Harris | [END_REF]] works the best under different scales and viewpoints.

Feature description

SIFT proposed by David G. Lowe [START_REF] David G Lowe | Object recognition from local scale-invariant features[END_REF]] is one of the most famous feature descriptor. After having extracted interest points, the 16 × 16 neighborhood of each interest point is uniformly divided into 4 × 4 subregions. For each subregion, the sum of the gradient magnitudes are calculated in 8 directions. The descriptor of the interest point is created by assembling all the sums of the gradient magnitudes in the 8 directions of all the 4 × 4 subregions, resulting in a 128 dimension vector (cf. Following this idea, other kinds of descriptors were then proposed.

SURF [START_REF] Bay | [END_REF]] is built by relying on the distribution of first order Haar wavelet responses in x and y directions. BRIEF [Calonder 2012] is created by comparing two sequences of pixels located near the extracted feature, which produces a binary descriptor. These two methods are used for interest point features. For line segment features, [Bay 2005] uses color histograms as descriptor, [Wang 2009b] adopts a SIFT-like approach to create MSLD descriptor. [START_REF] Forssén | Shape descriptors for maximally stable extremal regions[END_REF]] applies SIFT to create an interest region descriptor.

Geometry based localization

Feature matching

Since descriptors are high dimension vectors, the most common approach is to calculate a distance between different vectors, such as the Euclidean distance [START_REF] David G Lowe | Object recognition from local scale-invariant features[END_REF] or the Hamming distance [Calonder 2012]. For each feature in the scene, one needs to find its corresponding features in the database of model features. This is equivalent to finding the nearest neighbor of the scene feature from all model features in the descriptor space. To this end, either an exhaustive search [Calonder 2012] or approximate kd-trees [START_REF] Jeffrey | Shape indexing using approximate nearest-neighbour search in high-dimensional spaces[END_REF]] are used. To prevent false matchings, a threshold on the distance value [Bay 2005] or on the distance ratio between the first two nearest neighbors [START_REF] Lowe | [END_REF]] is used. To remove the threshold parameter, one can use PROSAC [START_REF] Chum | [END_REF]], a geometric based approach (cf. Section 2.3.3.1), to help find point correspondences.

From another point of view, feature matching can be regarded as a classification problem given feature descriptors. The view set of a model feature under different transformations defines a class. Classes are trained offline and scene features are classified online. Two features belonging to the same class establish a matching.

Random trees [START_REF] Amit | [END_REF]] and random ferns [START_REF] Ozuysal | [END_REF]] have been used and provide good results.

Geometry based localization

Different from texture based methods, geometry based localization does not rely on any photometric information to establish feature correspondences, or to directly find the camera pose. Since no color information is used, this kind of approach can register information in different representations. For example, the same object being painted in different colors can be aligned with its CAD model with no previously defined color; raster paper maps can be registered with information from GIS; an infrared image of a scene can be aligned with its visual appearance, etc. Due to the lack of discriminative information, methods based on geometric information are often more computational expensive and thus only few work have been proposed for real-time applications, e.g. AR. Therefore, the following of this dissertation is not limited to real-time methods in order to give a wider view in this section.

In 2D images, geometric quantities are points, edges, regions (or shapes), etc.

Methods are divided into region based approaches in Section 2.

Region based approaches

A 3D object usually forms a region (or a shape) in a 2D image. Considering the fact that regions are the most complex geometric structures in 2D images, one can obtain more information than edges or points to perform registration. As a consequence, many texture-less registration methods rely on the shape or the contour of regions.

To register planar objects, [START_REF] Goshtasby | [END_REF]] uses regions' mass centers as control points and applies a point pattern matching method to solve the correspondence problem. [Flusser 1994] creates a 4 dimensional affine invariant descriptor from each region. Each descriptor is associated with its homologous mass center, which facilitates the point pattern matching problem. [Martedi 2013] use joints of line segments extracted from the region's contour as control points. Support Vector Shape (SVS) [ Van Nguyen 2013] is a descriptor of shape constructed using Support Vector Machine (SVM) [Cortes 1995], where the resulting support vectors are control points near the contour of the shape. [Campbell 2015] by using control points of SVS boils down planar shape registration to point pattern matching. [START_REF] Holzer | [END_REF]] virtually generates different views of the model regions, then creates Distance Transform (DT) [Rosenfeld 1968] from contours of each virtual view. These distance transforms are then used to train random ferns. The pose estimation problem is solved by online classification. The method requires the target to contain several regions and can achieve near real-time matching. [Rothwell 1992] proposes a method to find four perspective invariant points on a general continuous concave curve. The method is applied by [START_REF] Hagbi | [END_REF]] to calculate a perspective invariant signature for each concave curve in order to distinguish different concave shapes. It can estimate the perspective transformation with only one concave shape. [START_REF] Donoser | [END_REF]] combines the method of [START_REF] Holzer | [END_REF]] and the one of [START_REF] Hagbi | [END_REF] to work with at least one concave shape. They achieve real-time matching.

To register 3D objects, a common approach is adopted by most researchers.

First, 2D views of the 3D model object are virtually generated from a set of possible viewing angles and stored in a database with their corresponding poses in an offline phase. Then during online matching, the best virtual image with respect to the scene image is found by contour matching. At last, the pose corresponding to the matched model image is used as an initial value and refined by a recursive tracking process. [START_REF] Choi | [END_REF]] uses chamfer matching for contour matching and particle filter for pose refinement. [Reinbacher 2010] proposes to use a hierarchical structure to organize the virtual views in order to speed up the matching. [START_REF] Holzer | [END_REF]] treats the best contour matching problem as a classification problem. He uses the distance transforms of all virtual views to train ferns and the scene image is classified later on. As contour matching is a very essential step, some researchers aim at developing more efficient and robust algorithms. [Damen 2011] discretizes the contour to form edgelets (short straight segments) and creates a constellation of edgelets by casting a ray from the one of the edgelets and following reflections of the ray on other edgelets.

The descriptor of this constellation is composed of the parameters (reflecting angles, zig-zag length ratios) of its generating path, which are used for contour matching.

BORDER [Chan 2016] is another descriptor for contour matching. It samples the contour to create a sequence of linelets (small line segments of equal length). Angles between linelets are used to calculate the contour's descriptor, which are matched by kd-trees.

Line segment based approaches

Due to noise, occlusion and broken detected lines, the detected endpoints of line segments are hardly reliable [Li 2016]. This effect coupled with over or under detections, makes it difficult to register line segments. [START_REF] Coiras | [END_REF]] uses extrapolation of lines to create triangles. An exhaustive search on triangle matching is adopted to find the final solution. [START_REF] Guan | [END_REF]] also uses lines or extrapolation of lines but converts line matching to point matching by using their intersections. [Wang 2013] transforms lines to points by using their coordinates in the parametric space and thus converts the line matching to a point matching problem. But it does not work for affine transformations. [START_REF] Long | [END_REF]] models each line segment as a 6D point and employs

Gaussian Mixture Model and Expectation-Maximization to register the high dimension point set. This method is used to register a map with an aerial image. But it cannot converge when the rotational angular difference between two targets is larger than 40 ∘ . [Gros 1998] uses descriptors constructed by neighboring line segments for similarity or affine transformations. It calculates transformations between matched descriptors and finds the clusters of these transformations in parametric space. Line Signature (LS) [Wang 2009a] creates a 11 dimension descriptor for each line segment by using its 5 neighboring segments. It deals with affine transformations and uses a hash table like approach for matching. The BOLD descriptor [Tombari 2013a] is a vector containing two angles between two line segments and the line connecting the middle point of them. BOLD descriptors are matched with the help of kd-trees in order to find correspondences. This method is not affine invariant since angles change under affine transformation.

Wide baseline point pattern matching

Pure Point Pattern Matching (PPM), i.e. registration of different point sets with only their coordinates, is a common problem in many research field such as image registration, astronomy [START_REF][END_REF]] and biochemical tests [START_REF] Datta | [END_REF]]. Major algorithms can be divided into four categories: RANSAC like approaches, geometric invariant approaches, parametric space approaches and optimization approaches.

RANSAC like approaches

RANSAC [START_REF] Fischler | [END_REF]] is a random approach for PPM. First one subset of correspondences is randomly chosen among all possible point correspondences in order to estimate a first "hypothesis" of the transformation. Then other correspondences are used to check whether this hypothesis is valid. The hypothesis that gives rise to the highest number of inliers is chosen as final output. RANSAC is an efficient method usually when the inlier ratio (in the sense of correspondence) is higher than 50% [START_REF] Lowe | [END_REF]]. Some improvements on RANSAC were then proposed to cooperate with cases where inlier ratios are lower. BaySAC and SimSAC [START_REF] Botterill | [END_REF] reduce the probability of correspondences to be selected when they have already led to a false matching in previous hypotheses. NAPSAC [Myatt 2002] is based on the observation that if two correspondences match from two neighboring points in the point set 𝑃 to two neighboring points in the point set 𝑄, then these two correspondences will probably be correct. It chooses the first correspondence uniformly randomly as in RANSAC. But for other correspondences, it favors the selection of those in the neighborhood of the first correspondence. Nevertheless, these methods can hardly work for pure PPM problems, since the number of all combinations of a model point 𝑝 ∈ 𝑃 and a scene point 𝑞 ∈ 𝑄 is very large [START_REF] Denton | [END_REF]].

If biased information on point set or correspondences is available other than only point coordinates, other sampling mechanisms are available. With available point correspondence reliability metrics, PROSAC [START_REF] Chum | [END_REF]] draws the correspondences with higher reliability first during the creation of subsets. Group-SAC [START_REF] Ni | GroupSAC: Efficient consensus in the presence of groupings[END_REF]] is a better choice when source point sets can be classified into different groups, with the help of additional information from images. [Raguram 2013a] has integrated several other methods mentioned above into a unified package (USAC). 20

Geometry based localization

Geometric invariant based approaches

Geometric invariant approaches include all methods which exploit inter-point spatial relationships for matching. Geometric Hashing (GH) [START_REF] Wolfson | [END_REF]] is one of the most famous and classical methods in this category. It uses subsets of points to construct coordinate bases so that other points' coordinates under such bases are invariant to a predefined transformation type. In a first offline step, it chooses every possible basis in the model, calculates other points' coordinates in each basis and stores them into a hash table. During the matching, random bases are iteratively selected in the scene, the coordinates of the remaining points' are expressed in each basis and matched with the model basis in hash table. The pair of bases which receives the highest number of votes gives the final result. GH is robust against extra/missing points and jitter but its major drawback is the huge computational complexity in 𝑂(𝑛 𝑏 ), where 𝑏 -1 is the number of points needed to form a basis [START_REF] Wolfson | [END_REF]]. [Lamdan 1991] and [START_REF] Rigoutsos | A Bayesian approach to model matching with geometric hashing[END_REF]] discussed the influence of noise on GH. [START_REF] Wolfson | [END_REF]] gives a good review on classical GH. Instead of using the whole point set for each basis, [START_REF][END_REF]] only uses the coordinates of two points under the coordinate system defined by two other points to encode geometric relationship into a "codebook", which reduces the memory cost in case of large amount of points. The proposed "codebook" can be seen as descriptors.

[ Van Wamelen 2004] proposes to match the nearest neighbors of two points to accelerate the process. While these two pieces of work are restricted to similarity transformations, others use affine invariant geometric descriptors. [Heyl 2013] improves over [START_REF][END_REF]] by using surface-ratios to deal with affine transformations and uses kd-trees for descriptor matching. Its time complexity is 𝑂(𝑛 4 log 𝑛) with respective to the number of points 𝑛, which is not possible to be used for real-time application. LLAH [START_REF] Nakai | [END_REF]] also uses surface-ratios but only relies on neighbor points to create descriptors. A hash table is used for descriptor matching. This method will be explained in detail in Section 2.3.4.

Parametric space based approaches

Transformations are uniquely defined by transformation parameters. For example, an affine transformation can be determined by 6 parameters of real value. The parameters of a specific transformation type form a high dimension space. The ground truth transformation corresponds to a point in this parametric space. Some methods aim at finding this parametric point in order to solve the point pattern matching problem. Pose clustering [START_REF] Olson | [END_REF]], also known as Generalized Hough

Chapter 2. Related work Transform (GHT) [Ballard 1981], is a representative method in this category. It discretizes the parametric space into bins, enumerates all possible transformations and votes in those bins. The bin receiving the highest number of votes is chosen as the result. However, the enumeration makes these methods computationally intensive. [START_REF] Huttenlocher | [END_REF]] divides the parametric space into subspaces and uses the partial Hausdorff distance to measure the matching score, where a smaller score means a better matching. If the matching score of a subspace is smaller than a threshold, this subspace itself will be divided to smaller subspace for further search.

[ [START_REF] Mount | Efficient algorithms for robust feature matching[END_REF]] proposes to use branch-and-bound breadth-first search for acceleration by introducing lower bounds of subspaces. The parametric space is divided into small cells (i.e. subspaces) and the upper and lower bounds of partial Hausdorff distances are calculated for each small cell. A best value is managed and updated when a smaller upper bound is found. If a cell's lower bound is larger than the best value, the cell is rejected and will not be processed further. This process continues until a satisfactory transformation is found. This method is still not efficient especially when the dimension of the parameter space is high. [START_REF] Datta | [END_REF]] combines geometric invariant approach and parametric space approach. It constructs 4-points affine invariants (4PAI) for each point by using their nearest points. Segment-ratios (which are affine invariants) are created from each 4PAI and are used for 4PAIs matching.

The matched 4PAIs then generate clusters in the transformation parameter space to find the affine transformation.

Optimization based approaches

For optimization approaches, divergence estimators [Jones 2001] of two point sets are defined as objective functions to be minimized. The most classical method of this group is Iterative Closest Points (ICP) [Besl 1992b]. It assigns to each model point its nearest neighbor in the scene point set with the help of an initial guess of the transformation. Then it estimates an updated transformation with these correspondences. These two steps iterate until a threshold is reached. The process minimizes the sum of correspondence distances and guarantees to find a local minimum. But the method is slow for real-time application when the difference between the poses of two point patterns is big and the basin of convergence of the global minimum is quite narrow. [Fitzgibbon 2003] uses the Levenberg-Marquardt algorithm to minimize the same objective function and updates the intermediate point correspondences for derivative computation as well. It improves the basin of convergence to ±60 ∘ in terms of rotation in some cases. Instead of assuming one-to-one

Geometry based localization

correspondence based on the nearest neighbors in each step, robust point matching [Gold 1998] proposes one-to-many relaxations. Each model point is assigned to several scene points with different weights according to the distances between them.

These weights are arranged in a correspondence matrix, which collapses to a binary (one-to-one) correspondence matrix during a deterministic annealing.

Researchers also adopted Gaussian mixture models (GMM). Point sets are in- Mixtures. At beginning, the size of the Gaussian kernel is set large, and thus a rough alignment can be found by minimizing the L2-distance. Then, the size of the Gaussian kernels is gradually decreased (i.e. annealing) to achieve a finer alignment.

terpreted
Point correspondences are not found explicitly. This method can be seen as Kernel Correlation (KC) [Tsin 2004] specialized with a Gaussian kernel. [Jian 2011] also points out that the objective function of their method, the original ICP, EM-ICP and KC are all concrete forms of a general divergence family [Jones 2001]. Similarly, the performance of these methods depends highly on the rotation difference between two point patterns.

Besides linear transformations (also called rigid transformation in most work), researchers get more and more interested in non-linear (i.e. non-rigid) transformations, which can be used for perspective transformation for example. [START_REF] Chui | [END_REF]] adds a regularization term in the objective function which limits the magnitude of second order derivatives of the transformation function. [START_REF] Myronenko | [END_REF]] uses motion coherence theory, which imposes the preservation of topological structure on point sets, to regularize the transformation function. But these modifications did not improve the robustness against rotation.

Chapter 2. Related work 2.3.3.5 Other approaches [Wang 2012[START_REF] Caetano | [END_REF]] treat point pattern matchings under Euclidean transformations as graph matching problems. As their methods base on the invariance of distances between points, it cannot be generalized to affine or perspective transformations. [Shapiro 1992] proposes the spectral context. It calculates the Gaussian distance between every point pair in the same point set and constructs a proximity matrix from these distances. Singular Value Decomposition (SVD) maps each point to an eigenvector of the proximity matrix. In this way, two point sets can be converted to two sets of eigenvectors and the distances between eigenvectors are used for matching. However, this method breaks down if two patterns have different numbers of points. [Tang 2014] uses the same idea but on local point sets to construct descriptors for matching and introduces an "approximate distance order" to deal with noise. Unfortunately, this method takes only similarity transformations into consideration. Evolutionary algorithms [Agrawal 1994, Ezoji 2006, Yin 2012, Zhang 2003] and Particle Swarm Optimizations (PSO) [Yin 2006] have also been used. But they either deal with only similarity transformations, or require certain model points to appear in the scene.

Random Dot Markers

Although there exist numerous approaches for the point pattern matching, algorithms which are efficient enough for real-time wide baseline localization are quite limited. From this point of view, feasible methods are the local geometric invariant based approaches. The most representative method is LLAH [START_REF] Nakai | [END_REF]],

which can be applied to achieve real-time tracking of Random Dots Markers (RDM) [Uchiyama 2011b] (cf. Fig. 2.3). It is the starting point of this dissertation and is the closest related work. Therefore, this section provides a detailed explanation of this method and discusses its drawbacks.

Locally Likely Arrangement Hashing (LLAH)

LLAH was originally designed as a document retrieval system. The original parameters of LLAH are chosen here for the sake of a simple explanation. LLAH is based on area-ratios, which are affine invariants, to create descriptors. Given four coplanar points 𝐴, 𝐵, 𝐶, 𝐷 and 𝑆(•) being the surface of a triangle, the invariant 𝑓 is defined as follows: For each descriptor, the discretized area-ratios of all the combinations of 4 points out of 5 points are used, i.e. 5 ratios. These 5 ratios 𝑓 0 , .., 𝑓 4 are concatenated head-to-tail to form the descriptor, i.e. an integer 𝑑 = 𝑓 0 𝑓 1 𝑓 2 𝑓 3 𝑓 4 .

The matching process is done in two phases: (1) descriptors of model points are calculated offline and stored in a hash table ; (2) during online matching, descriptors of scene points are calculated. Then model points are found in the hash table by matching descriptors. A vote is casted for each model-scene point matching.

RANSAC is finally used to find the homography. Authors claim that the time complexity is linear to the number of points [START_REF] Nakai | [END_REF]].

RDM is developed upon LLAH to track random point patterns in real-time.

Since LLAH cannot match point pattern under large perspective distortions, RDM additionally manages a dynamic hash table that is updated from the previous frame.

As a consequence, the tracking of RDM should be initialized from an approximate front view.

Limitations

The efficiency of the algorithm comes from the fact that it uses 5 surface ratios to construct one descriptor, which makes it discriminant. The algorithm is ro-Chapter 2. Related work bust against small noise thanks to the use of several descriptors per point. RDM has successfully been used to augment paper maps with GIS when road intersections are manually extracted [START_REF] Uchiyama | [END_REF]]. But when we use an automatic extraction of road intersections instead (i.e. more and larger detection errors), such as [Callier 2011], it does not work any more.

To summarize its limitations, LLAH/RDM is quite sensitive when noise becomes bigger or the perspective distortion is large. One of the most common noise is the one on points' positions. Due to optical or numerical noise, the detection of a point may be located at several pixels away from its ground-truth position (i.e. jitter).

This effect is more obvious when points to be detected are not visually salient, such as road intersections on a map (cf. Fig. 2.5). As to perspective distortions, the position of a point under perspective transformation is slightly different from that under affine transformation if the depth of field is small. This difference has a similar effect as the jitter. Hence if an algorithm can manage well jitters, it is able to deal better with perspective distortions too. is added on each coordinate component of all points in 𝑄. Fig. 2.5 shows that the number of descriptors which can be correctly matched dramatically decreases with respect to the augmentation of noise on point coordinates. This leads to a total failure of the algorithm when the noise 𝜎 becomes larger than 4 pixels.

The reason for the algorithm's sensitivity can be explained as follows: when the position of a point changes slightly, the surface ratio in (2.1) changes but may not influence its discretized value 𝑓 . But if the noise on points' positions is larger, it is more likely that 𝑓 changes as well. Since the descriptor is a concatenation of 5 𝑓 s, the change of any 𝑓 𝑖 leads to a different integer, and thus the matching fails. The discrimination of the descriptor makes it quite sensitive to noise on points' positions. It is difficult to improve the algorithm's performance by simple adjustments because of the following constraint: if one makes the descriptor more robust (by reducing 𝐾 for example), it will become less discriminant and thus increases the number of false matchings. The latter makes the vote result less reliable, which will make RANSAC work slowly, and maybe even lead to a matching failure.

In a word, the dilemma between robustness and speed of PPM is not solved by LLAH or any other method. One needs to develop new algorithms to achieve better performance.

Chapter 3

Robust random dot markers (RRDM) We reformulate the problem mathematically as follows: let 𝑃 be a model point set in 2D coordinates; 𝐻 an unknown 2D homography. Let 𝑄 be an observed scene point set, of which some points belong to 𝐻(𝑃 ) but with small differences in their positions, while others are randomly added noise points (i.e. extra points). The problem is to determine 𝐻 given 𝑃 and 𝑄 (cf. Fig. 3.1). We assume that points in 𝑃 are uniformly randomly distributed. The challenges of the problem are large point jitter, extra and missing points (cf. Section 1.2.3). Moreover, we have to estimate 8 parameters of the 2D homography 𝐻, which is more difficult than 2D similarity (with 4 parameters) and 2D affinity (with 6 parameters).

A robust descriptor

We know that (a) pure point pattern matching can only rely on the geometrical pattern itself; and that (b) the presence of extra/missing points and point jitter can significantly alter local geometry of the point pattern. Our task is difficult since we have to take into account both (a) and (b) at the same time. As seen in Section 2.3.4, the LLAH descriptor is an integer composed of five surface ratios computed from seven neighboring points. Any change in those seven neighbors can modify the descriptor, making it very sensitive to noise. In order to solve this problem, we propose a Two-Surface-Ratios (TSR) descriptor. It is constructed only by 4 points (a so-called Quad-Point, cf. Fig 3 .2), so that it is less likely to be influenced by (b).

A robust descriptor

We choose four points since fewer points cannot provide any affine invariant. In this section, we define the TSR descriptor and show how two jittered Quad-Points can be matched under an affine transformation which lays the groundwork for our new algorithm.

Definition

We call any four-coplanar-point configuration a Quad-Point. We do not consider the situation where three points are aligned since it rarely occurs for randomly distributed points. Each Quad-Point is composed of a main point 𝑀 and three associated points 𝐹, 𝑅, 𝐿. Under such a definition, every Quad-Point can be categorized into one of three possible configurations 𝑐. For each configuration, the spatial relationships between these four points, namely 𝑀 (main), 𝐹 (face), 𝑅 (right), 𝐿 (left), are illustrated in Fig. 3.2. Especially, when 𝑐 = 3, 𝐹 is chosen so that the surface of triangle Δ𝑀 𝑅𝐿 is the second smallest among the four triangles. As a consequence, a Quad-Point, denoted as Q = ⟨𝑀 ; 𝐹, 𝑅, 𝐿⟩ (note that the order of the points is important), is uniquely defined when four coplanar points are given and its main point 𝑀 chosen (we therefore say that the Quad-Point belongs to 𝑀 and denote it as Q(𝑀 ) when all four points are known).

Q ≡ Q ′ ⇔ ∃A s.t. 𝑀 ′ = A (𝑀 ), 𝐹 ′ = A (𝐹 ), 𝐿 ′ = A (𝐿), 𝑅 ′ = A (𝑅) (3.1)
All Affinely Equivalent Quad-Points (AEQ) can be grouped together. Since a Quad-Point contains exactly two independent invariants [START_REF] Hartley | [END_REF]], each AEQ group has two degrees of freedom.

We propose to rely on two surface ratios along with a configuration specifier to describe each AEQ group. Given its nature, we name our descriptor Two-Surface-Chapter 3. Robust random dot markers (RRDM) Ratio (TSR). In this context, Q(𝑀 )'s TSR descriptor, t(Q(𝑀 )), is defined as follows:

t(Q(𝑀 )) = (𝑋, 𝑌, 𝑐) (3.2) 𝑤𝑖𝑡ℎ ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ 𝑋 = 𝑆 Δ𝐹 𝑅𝐿 𝑆 𝑎 , 𝑌 = 𝑆 Δ𝑀 𝐹 𝐿 𝑆 𝑎 , for 𝑐 = 1, 2 (3.3) 𝑋 = 𝑆 Δ𝑀 𝑅𝐿 𝑆 𝑎 , 𝑌 = 𝑆 Δ𝑀 𝐹 𝐿 𝑆 𝑎 , for 𝑐 = 3 (3.4)
where 𝑆 𝑎 is the surface formed by the convex hull of Q(𝑀 ).

We prove in the following section that two Quad-Points are affinely equivalent if and only if their TSR descriptors are equal. Note that in the absence of point jitter, when two Quad-Points are matched under an affine transformation, they have the same TSR descriptor.

Affine invariance

We want to demonstrate that two Quad-Points are affinely equivalent if and only if their TSR descriptors are equal, which means:

t(Q) = t(Q ′ ) ⇔ Q ≡ Q ′ (3.5)
As Q and Q ′ have the same affine invariant, the backward demonstration is direct from (3.1) and (3.2). We focus on the forward demonstration.

We first show that for a Quad-Point Q, the position of 𝐹 is uniquely defined by the positions of 𝑀, 𝐿, 𝑅 and 𝑡(Q) = (𝑋, 𝑌, 𝑐).

We use superscript • (𝑖) to represent a point's i-th coordinate. Without losing generality, we translate a Quad-Point so that 𝑀 locates at (0, 0). The 2D plane is divided into seven regions (cf. Fig. 3.3). According to our definition of TSR (Section 3.1.1), 𝐹 cannot lay in since it would be a "left point" rather than a "face point". For the same reason, regions ∼ do not fit. So, 𝐹 can only lay in (c=1), (c=2) or (c=3).

For 𝑐 = 1: With (3.2), we have:

𝑆 𝑎 = 1 2 ⃒ ⃒ ⃒ ⃒ ⃒ 𝑅 (1) 𝐹 (1) 𝑅 (2) 𝐹 (2) ⃒ ⃒ ⃒ ⃒ ⃒ + 1 2 ⃒ ⃒ ⃒ ⃒ ⃒ 𝐹 (1) 𝐿 (1) 𝐹 (2) 𝐿 (2) ⃒ ⃒ ⃒ ⃒ ⃒ 𝑆 𝐹 𝐿𝑅 = 1 2 ⃒ ⃒ ⃒ ⃒ ⃒ 𝑅 (1) -𝐿 (1) 𝐹 (1) -𝐿 (1) 𝑅 (2) -𝐿 (2) 𝐹 (2) -𝐿 (2) ⃒ ⃒ ⃒ ⃒ ⃒ 𝑆 𝑀 𝐹 𝐿 = 1 2 ⃒ ⃒ ⃒ ⃒ ⃒ 𝐹 (1) 𝐿 (1) 𝐹 (2) 𝐿 (2) ⃒ ⃒ ⃒ ⃒ ⃒ (3.6)

A robust descriptor

Z ( 𝐹 (1)
𝐹 ( 2)

) = ( -𝑅 (1) 𝐿 (2) + 𝑅 (2) 𝐿 (1) 0 ) with Z = [ (1 -𝑋)(𝑅 (2) -𝐿 (2) ) -(1 -𝑋)(𝑅 (1) -𝐿 (1) ) -𝑌 𝑅 (2) -(1 -𝑌 )𝐿 (2) 𝑌 𝑅 (1) + (1 -𝑌 )𝐿 (1)
]

(3.7)

This equation gives a unique solution when |Z| ̸ = 0, meaning that the coordinates of 𝐹 are uniquely defined. It can be shown that:

|Z| = -(1 -𝑋) ⃒ ⃒ ⃒ ⃒ ⃒ 𝑅 (1) 𝐿 (1) 𝑅 (2) 𝐿 (2) ⃒ ⃒ ⃒ ⃒ ⃒ = -2(1 -𝑋)𝑆 𝑀 𝐿𝑅 ̸ = 0 (3.8)
A similar approach can be applied for 𝑐 = 2 and 𝑐 = 3 for which |Z| = 2(1 + 𝑌 )𝑆 𝑀 𝐿𝑅 and |Z| = 2𝑋𝑆 𝑀 𝐿𝑅 respectively and where |Z| cannot be zero.

The forward demonstration is as follows: Assume that Q = ⟨𝑀 ; 𝐹, 𝑅, 𝐿⟩ and

Q ′ = ⟨𝑀 ′ ; 𝐹 ′ , 𝑅 ′ , 𝐿 ′ ⟩ have the same TSR descriptor. We can calculate an affine transformation A which makes 𝑀 ′ = A (𝑀 ), 𝑅 ′ = A (𝑅), 𝐿 ′ = A (𝐿). Then let 𝐹 * = A (𝐹 ), then Q * = ⟨𝑀 ′ ; 𝐹 * , 𝑅 ′ , 𝐿 ′ ⟩ and Q are affinely equivalent and thus t(Q * ) = t(Q) = t(Q ′ ). Since Q ′ and Q * have

Point jitter and descriptor variance

There are two possible sources for point jitter: (i) distortions from the imaging system and (ii) misdetection of points in the image. The former case can be induced by lens distortions and is both scene and scale independent. In the latter case, point Chapter 3. Robust random dot markers (RRDM) jitter is scale dependent, becoming smaller when the scene gets further away from the camera. For narrative convenience, in the following we only consider point jitter as arising from the second source.

With the presence of point jitter, the shape of Quad-Points can be changed.

The original Quad-Point and its jittered version are generally not affine equivalent, thus they cannot be exactly matched according to (3.1). However, when jitter is not too significant, these two Quad-Points should be "quasi-affine-equivalent". We show that the TSR descriptor can deal with this quasi-affine-equivalence.

We assume that point jitter can be modeled as an additive Gaussian noise to the points' original positions before perspective transformation. Mathematically, each point in the Quad-Point is affected by an additive Gaussian noise 𝛿 ∼ 𝑁 (0, 𝜎 2 ) on two coordinates independently (cf. Fig. 3.4), where 𝜎 is small compared to distances between points forming the Quad-Point.

𝑑𝑋 =𝑋 ′ -𝑋 = 𝑆 𝐹 𝑑𝑆 𝑀 -𝑆 𝑀 𝑑𝑆 𝐹 (𝑆 𝑀 + 𝑆 𝐹 ) 2 ∼ 𝑁 (0, Σ 2 𝑋 ) Σ 2 𝑋 = 𝜎 2 4 (𝑆 𝑀 + 𝑆 𝐹 ) 4 {𝑆 2 𝑀 ( 𝐹 𝑅 2 + 𝐹 𝐿 2 + 𝑅𝐿 2 ) + 𝑆 2 𝐹 ( 𝑀 𝑅 2 + 𝑀 𝐿 2 + 𝑅𝐿 2 ) -𝑆 𝑀 𝑆 𝐹 ( 𝑀 𝑅 2 + 𝑀 𝐿 2 + 𝐹 𝑅 2 + 𝐹 𝐿 2 -2𝑀 𝐹 2 ) } (3.9)
where 𝑆 𝑀 = 𝑆 Δ𝐹 𝑅𝐿 and 𝑆 𝐹 = 𝑆 Δ𝑀 𝑅𝐿 .

So, when a Quad-Point is slightly perturbed by 𝛿 with a small 𝜎, its TSR descriptor will mostly remain unchanged. Thus, we can define two Quad-Points Q

and Q ′ such that they are quasi-affinely-equivalent (Q ≃ Q ′ ) thanks to our TSR descriptor.

Q ≃ Q ′ ⇔ |𝑋 -𝑋 ′ | ≤ 2 Σ𝑋 , |𝑌 -𝑌 ′ | ≤ 2 Σ𝑌 and 𝑐 = 𝑐 ′ (3.10) with 𝑡(Q) = (𝑋, 𝑌, 𝑐) and 𝑡(Q ′ ) = (𝑋 ′ , 𝑌 ′ , 𝑐 ′ ).
Practically, we set Σ𝑋 (resp. Σ𝑌 ) as follows, with Σ 𝑚𝑎𝑥 being a parameter of the method:

Σ𝑋 = { Σ 𝑋 , if Σ 𝑋 ≤ Σ 𝑚𝑎𝑥 (3.11) Σ 𝑚𝑎𝑥 , otherwise (3.12)
For clarity sake, we say that two Quad-Points are quasi-affinely-equivalent when the TSR descriptor of one Quad-Point is inside the matching region of the other.

The matching region of Q is represented by a gray box on the right of Fig. 3.4. In the following, we consider that two Quad-Points are matched if they satisfy (3.10).

Algorithm

Recall that 𝑃 is the model point set, and it contains the coordinates of 𝑚 model points. 𝑄 is the scene point set which contains the coordinates of 𝑛 image points that are detected in the image. Each pair (𝑝, 𝑞), with 𝑝 ∈ 𝑃, 𝑞 ∈ 𝑄 is a tentative correspondence (or in short, a correspondence) that could either be an inlier (i.e. a true correspondence) or an outlier (i.e. a false correspondence).

Both point jitter 𝜎 and point pattern density 𝜌 (i.e. the number of points per unit area) play an important role in Quad-Point matching. Given a fixed jitter value 𝜎, the denser the point pattern, the bigger the descriptor variance, cf. (3.9).

Chapter 3. Robust random dot markers (RRDM)

We define the jitter factor 𝜂 = 𝜎 √ 𝜌 to describe the amount of point jitter (cf.

Section 3.1.3).

The main idea of our approach is that even under a large perspective distortion, the transformation between 𝑃 and 𝑄 can be approximated by an affine transformation in a small region. We use local affine information and consistency between neighboring regions to establish point correspondences before estimating an homography.

Our method can be decomposed into two stages. The first one is in charge of registering every model point into hash tables. This stage is performed once offline before online tracking. During the second stage, we retrieve 𝑄 from each image, and match it with 𝑃 by using a Local Voting and Coherency (LVC) strategy.

Offline pre-registration

The offline pre-registration stage is detailed in Alg. 1. 𝑃 is loaded from a database representing the point pattern. For each point 𝑝 ∈ 𝑃 , we find its k-nearest neighbors 𝑁 𝑘 (𝑝). Each point 𝑝 is associated with a set of Quad-Points QS(𝑝), defined as follows:

QS(𝑝) = {Q(𝑝) = ⟨𝑝; 𝑝 1 , 𝑝 2 , 𝑝 3 ⟩ , ∀𝑝 1 , 𝑝 2 , 𝑝 3 ∈ 𝑁 𝑘 (𝑝)} (3.13)
The order of 𝑝 1 , 𝑝 2 , 𝑝 3 is predefined, so each set of Quad-Points QS(𝑝) contains 𝑘! 3!(𝑘-3)! Quad-Points. We use three bi-dimensional 𝑊 × 𝑊 (𝑊 is an integer parameter of the algorithm) hash tables to register model points so that they can be quickly matched during the online stage with image points according to their descriptors. Each hash 36 3.2. Algorithm table corresponds to one configuration of Quad-Points illustrated in Fig. 3.2. We determine the set of hash bins 𝑆 𝑖𝑑𝑥 = {(𝑖, 𝑗, 𝑐)} into which a Quad-Point Q 𝑚 with TSR descriptor t(Q 𝑚 ) = (𝑋, 𝑌, 𝑐) is registered as follows, with [•] represents a floor function:

[(𝑋 -2 Σ𝑋 )𝑊 ] ≤ 𝑖 ≤ [(𝑋 + 2 Σ𝑋 )𝑊 ] + 1, 𝑖 ∈ 𝑁 [(𝑌 -2 Σ𝑌 )𝑊 ] ≤ 𝑗 ≤ [(𝑌 + 2 Σ𝑌 )𝑊 ] + 1, 𝑗 ∈ 𝑁 (3.14)
For clarity sake, imagine that the square [0, 1]×[0, 1] is uniformly divided into 𝑊 2 subregions by vertical and horizontal lines as illustrated in Fig. 3.5. Each subregion represents a hash bin. We register the Quad-Point with t = (𝑋, 𝑌, 𝑐) into the hash bins which intersect the matching region 𝑅(𝑋, 𝑌, 𝑐). 𝑃 and 𝑄. During the rejection phase, we calculate a local affinity for each correspondence using their neighboring points. We discard correspondences for which affinities are not consistent with those of their neighbors. During the recovery step, we try to find inliers that were potentially wrongly rejected in phase (ii).

Local Voting

Similar to Section 3.2.1, a set of Quad-Points QS(𝑞) is generated for each point 𝑞 ∈ 𝑄. The computational cost is 𝑂(𝑘𝑛 log 𝑛 + 𝑘 3 𝑛).

Then for each Quad-Point Q 𝑠 ∈ QS(𝑞), a set of model Quad-Points QM is retrieved from hash bin 𝐻 𝑖𝑑𝑥 according to its TSR descriptor (𝑋, 𝑌, 𝑐). 𝐸 𝑠 : A 𝑝,𝑞 is the correct transformation between local correspondences 𝑠 𝑝,𝑞

q 2 q 4 q • • • • p 5 p 4 p 2 q 1 q 2 q 4 1 1 4+1 p 1 p 3 q 3 q 5 p 1 p 3 q 3 q 5 4 2 3 5+1 4 3 2 4+1 3 1 V p,q Figure
In other words, 𝐸 𝑠 means that q can be back-projected close enough to its correspondence p in 𝑠 𝑝,𝑞 by A -1 𝑝,𝑞 . Since A -1 𝑝,𝑞 q-p can be regarded as the noise before transformation (i.e. jitter), it follows a Gaussian distribution (cf. Section 3.1.3). As a consequence, 𝜒 2 tests for the sum of squared Gaussian errors can be used to verify 𝐸 𝑠 .

𝐸 𝑠 = ∑

(p,q)∈𝑠𝑝,𝑞

||A -1 𝑝,𝑞 q -p|| 2 < 𝜒 2 0.05,𝑠𝑖𝑧𝑒(𝑠𝑝,𝑞)-6

(3.17)

We use 𝑠𝑖𝑧𝑒(𝑠 𝑝,𝑞 ) -6 as the degree of freedom since the estimation of affinity needs 3 points and each point has 2 degrees of freedom. Tentative correspondences which fail the 𝜒 2 check (3.17) will be rejected.

We then group tentative correspondences that have similar affinities into an affinity group. Using a Singular Value Decomposition (SVD), we extract two scaling factors 𝜆 (1) , 𝜆 (2) along two principal axes and one rigid rotation 𝜃 from any affinity A . We consider two affinities to be similar (A 1 ≈ A 2 ) if they satisfy:

Δ𝜃 ≤ 𝜃 𝑡 , with Δ𝜃 = |𝜃 1 -𝜃 2 | 𝛼 -1 𝑡 ≤ 𝛼 (1) , 𝛼 (2) ≤ 𝛼 𝑡 , with 𝛼 (𝑖) = ⃒ ⃒ ⃒ ⃒ ⃒ 𝜆 (𝑖) 1 𝜆 (𝑖) 2 ⃒ ⃒ ⃒ ⃒ ⃒ (3.18)
Note that based on (3.18), 𝜃 𝑡 and 𝛼 𝑡 are two parameters required by our algorithm and will be discussed in Section 3.3.3.

Chapter 3. Robust random dot markers (RRDM)

We use the Hungarian algorithm for the maximum assignment problem during

Step 1 of Alg. 3. The Hungarian algorithm has a worst-case time complexity of 𝑂(𝑘 3 ). Nevertheless, it should be noted that since there is a large difference in the number of votes between inliers and outliers in 𝑉 𝑝,𝑞 (cf. Section 3.3.2), the algorithm is able to find the optimal solution in a few iterations.

Step 1 of Alg. 3 is the most expensive part of the algorithm with a worst time complexity of 𝑂(𝑚𝑛𝑘 3 ).

In Step 2, neighboring tentative correspondences having similar affinity are grouped together at first to create small affinity groups, which filters out a large amount of "isolated" tentative correspondences that do not share similar correspondences with their neighbors. Then these small groups are fused if they have similar affinity in order to create bigger affinity groups. for all (𝑝, 𝑞)

∈ 𝑠𝑙 do for all (p, q) ∈ 𝑠 𝑝,𝑞 do if |p -𝐻 -1 (q)| ≤ 2𝜎 then 𝑠𝑙 = 𝑠𝑙 ∪ {(p, q)}
until No new element added to 𝑠𝑙 Return H

Choice of parameters

The algorithm presented in Section 3.2 relies on several parameters. In this Section, we detail the meaning of each parameter and propose some default values so that one can easily adapt our algorithm to one's own applications and needs. So, robustness with respect to additional or missing points can be improved by increasing 𝑘, at the expense of impacting the computational efficiency of the algorithm (cf. Section 3.2.2), since most steps have a computational cost in the order of 𝑘 3 .

Parameters 𝜂 and 𝐷 𝑚𝑎𝑥 are used to manage point jitter. When points in the image are likely to have a high incertitude in their positions, 𝜂 should be set to a higher value. When incertitude is so big that (3.12) dominates over (3.11) (resp. for 𝐷 𝑌 ), a larger 𝐷 𝑚𝑎𝑥 should be used. Following those guidelines, one can augment the chances to match a jittered Quad-Point correctly, thus increasing the robustness of the algorithm against large point jitter. However, increasing 𝜂 and 𝐷 𝑚𝑎𝑥 results in enlarging the matching region and in a higher value of the redundancy number 𝑟 (cf. Section 3.2.2) which in turn impacts the efficiency of the algorithm. It is worth noting that the choice of 𝑣 𝑡 has an impact on both the robustness of the algorithm and its speed. Hence, its value should be chosen according to a

Threshold

Results

trade-off between robustness (lower 𝑣 𝑡 values) and computational efficiency (higher 𝑣 𝑡 values).

Influence of point jitter on A and 𝛼 𝑡 , 𝜃 𝑡

Let A 1 and A 2 be the affinities of two inliers estimated during step 1 of Alg. 3. Due to the presence of point jitter and perspective distortion, they are usually different even if we have used inliers for their computation. As a consequence, 𝜃 𝑡 and 𝛼 𝑡 have to be carefully chosen in order to decide whether two affinities are similar or not. In the following, we only study the influence of point jitter on the difference between affinities since perspective distortion in a small region can be regarded as a special case of point jitter.

As one can expect, the smaller the distance between points composing a Quad-Point, the more impact point jitter has on the computation of the affinities. Moreover, the fewer correspondences we use to estimate an affinity, the less precise the result. Both observations lead to the fact that the difference between two affinities A 1 and A 2 becomes more important when 𝑘 gets smaller. In other word, 𝑘 = 3 is probabilistically speaking the worst situation where the difference between A 1 and A 2 might be the greatest.

Results

In this section, we compare the efficiency and the robustness of our technique (denoted as RRDM), an algorithm combining local voting (cf. Section 3.2.2.1) and USAC [Raguram 2013b] (denoted as USAC) and Random Dot Markers (denoted as RDM). An homography is estimated only if at least 9 inliers are found. This value is given by RDM and we apply this rule to all the three methods for a fair comparison. Due to the presence of jitter, an exact estimation is impossible. We use Φ 3 mentioned in [START_REF] Huynh | [END_REF]] to measure if a matching is "precise". With q 0 being the true quaternion of the virtual marker plan and q the estimated one, a matching is "precise" if:

Φ 3 = arccos(‖q 0 • q‖) ≤ 1.5 ∘ (3.19)
In the "Point quantity experiment" (cf. Fig. 3.9), the influence of the amount of points on performance of methods is studied with 𝛽 = 15% and 𝜂 = 3%. Results are averaged from 1000 experiments. We show that even with a moderate point jitter and a small quantity of extra points, RDM gives a very poor estimation.

USAC's estimation is comparable with the one obtained by RRDM but the execution time explodes when 𝑚 increases due to lower inlier ratios. Our method remains consistently precise starting at 𝑚 = 100 even if RDM remains slightly more efficient than RRDM. This characteristic is very important since for a city like Edinburgh, for example, a map of the city center contains about 200 road intersections. RDM and USAC work best at 𝑚 = 100, so we choose this value for later comparison aiming at a more comparative study.

In the "Jitter experiment" (cf. Fig. 3.10), the LS curve is obtained by using the classical least square estimation with known point correspondence and can be considered as an upper bound [START_REF] Hartley | [END_REF]]. Statistics are obtained from 1000 experiments. It can be seen that both RDM and USAC suffer from a very steep performance decrease when 𝜂 increases.

Real markers

We studied the robustness of the three methods against extra points, missing points and partial occlusion with real markers. Two kinds of equipment are used: a 42 inches (930 × 523mm) LCD screen filmed at 1920 × 1080 by an iPad Air video camera and an A4 paper (297 × 210𝑚𝑚) filmed at 1280 × 720 with a Logitech C270 camera, as illustrated in Fig. 3.11. A video sequence of 200 to 300 frames is recorded for each case. A visual assessment is provided for all homographies estimated by these experiments. We have to notice that in these real images, no point jitter occurs, which cannot highlight the robustness of our method. We show that even in situations which are in favor of RDM, our method works better as well.

In the "Extra point experiment", 𝛽 × 𝑚 extra points are added inside markers to

Natural map tracking

We apply RRDM to a map geo-referencing problem in the "Map tracking experiment". Three different unprepared printed maps of downtown Edinburgh containing different textures are used, cf. Fig. 3.15a-3.15c. We manually extracted from map (𝑎) 233 intersections coordinates, which serve as the model point set (i.e. 𝑃 ). For each map, a video sequence containing about 100 frames is recorded with an iPad (with resolution 1280 × 720). Although only one geo-referencing operation is needed for each sequence of video in order to initialize tracking (if no tracking failure appears), we choose to apply it to each frame of the videos to show RRDM's behavior under diverse circumstances. For each frame, an automatic intersection detection technique based on color segmentation (similar to [Callier 2011]) with the help of manually generated filter offline is used to produce 𝑄 on the fly. In this experiment 𝜂 is set to 5%.

To give an idea of the complexity of the problem we tackle, we manually analyzed the detection results on 10 randomly picked frames for each video (cf. Table 3.1).

We consider a detected point to be a real intersection if it does not lay too far

Application: Augmented Maps

In this section, we propose a framework for augmented maps by using RRDM. It During online tracking, the paper map is tracked with a classical method based on SURF interest points and descriptors.

Intersection detection on real maps

Since the map-GIS registration is based on road intersections, the quality of road intersection extraction is crucial. False road intersections can lead to a poor performance of the method, or even to total failure. A usual way to find road intersections on a map is to identify the road network, extract its skeleton before finding lines intersections. The difficulty lies in the treatment of texts on maps, such as road names, which severely reduce the quality of the automatically extracted skeleton.

[ [START_REF] Callier | [END_REF]] relied on lines with homogeneous color in the image to find road pixels The contribution of our new method is mainly about road layer extraction improvement, and is shown in Fig. 3.19. At first, the user picks a single road pixel from the image (a seed point). This is very easy to do with a mouse or on a tactile screen. From this pixel, road color and width can be estimated. Then, a classic level set method is used to find the map in 𝐼 0 to avoid false detections coming from out of the region of interest. At last, a modified level set method is applied to find the road layer inside the map.

We use a 50 × 50 local image portion centered at the seed point to find the road width since road width on a camera image of a map is usually far less than 50 pixels. A Canny edge detector [Canny 1986] is used to extract road edges and a flood fill at the seed point is applied to find the local road layer, which contains 𝑅 road pixels. We then dilate the local road layer by 1 pixel and the result contains 𝑄 road pixels. Since roads can be modeled as thin rectangles, road width 𝑤 can be easily calculated as:

𝑤 = 𝑤 1 + 𝑤 2 2 with 𝑤 1 = 2𝑅/(𝑄 -𝑅) , 𝑤 2 = 2𝑅/(𝑄 -𝑅 + 4) (3.20)
where 𝑤 1 is for the case where a single road is present in the local image portion while 𝑤 2 is for the case where a crossing is present. These are the two most common cases in a local map.

The level set method, commonly known as active contour, is often used in image segmentation. It separates an image into different segments (i.e. subregions) such that colors of all pixels in the same segment are similar. It usually starts with an initial contour Γ 0 and then moves it iteratively in order to minimize 𝐸(Γ) in (3.21)

where Γ is the contour which separates segments 𝐴 1 and 𝐴 2 . Shi's implementation is used as our classic level set method [START_REF] Shi | [END_REF]]:

𝐸(Γ) = ∫ 𝐴 1 |𝐼 0 (𝑥, 𝑦) -C 1 | + ∫ 𝐴 2 |𝐼 0 (𝑥, 𝑦) -C 2 | (3.21)
where 𝐼 0 (𝑥, 𝑦) is the color at (𝑥, 𝑦). C 1 and C 2 are expected colors for the two segmentations respectively.

At first, a classic level set method is used to extract the map from the image 𝐼 0 . An initial rectangle contour Γ 0 which is just a little smaller than 𝐼 0 is used (cf. 𝐼 0 (𝑥, 𝑦) = arg min

|𝑥 ′ -𝑥|≤2,|𝑦 ′ -𝑦|≤2 |𝐼 0 (𝑥 ′ , 𝑦 ′ ) -𝐶 1 | (3.22)

Results

We recorded three videos of different paper maps of Munich with an iPad-air camera (with resolution 1920 × 1080) to test our method. We use OpenStreetMap data as our GIS database and raster maps found on the Internet and printed (see Fig. 3.21). the end (cf. Frame 167 in Fig. 3.24). As to time performance, although registration takes 500∼1000ms, tracking works at interactive rate, i.e. 20∼30Hz.

Conclusion

In this chapter, we have introduced RRDM, a method to robustly track random dot patterns, where over and under point detection, large perspective distortion and inaccurate detection (jitter) can be taken into account. Unlike traditional feature point methods, RRDM is also able to robustly track 2D texture-less object. Recommendations on parameters are also given in order to facilitate user's implementation.

Synthetic images are used to demonstrate the robustness of RRDM against point jitter and when the marker is composed of a large number of points. In those cases, RRDM clearly outperforms RDM. We also proved through experiments with real videos that RRDM is more robust than RDM in presence of partial occlusion, randomly missing pattern points or under large perspective distortion, and especially with extra points in the pattern. Regarding efficiency, a real-time tracking (30Hz) can be achieved when the number of points in a marker is moderate (about 100 points). We also applied RRDM on registration of unprepared printed maps by using purely geometric coordinates of road intersections as the point pattern. Finally, we 56 3.6. Conclusion 𝑑 is the number of coordinates of each point. As a consequence, 𝑙 𝑖 = 𝜌 -1 𝑖 is a measure of inter-point distance. Without losing generality, we assume that all point sets are normalized beforehand so that the center of mass of each point set is at origin and 𝑙 𝑖 = 𝑙 is now the same for all models. We assume also that the jitter follows an uncorrelated multivariate normal distribution 𝑁 (0, 𝜎I) for all models,

where 𝜎 = 𝑙𝜂 and 𝜂 is the jitter factor (cf. Section 3.1.3). From now on, all point sets are normalized ones unless otherwise stated.

Our algorithm is composed of three modules: hypotheses generator, hypotheses validator, and results refiner (cf. Fig. 4.4). 𝑀 containers 𝐶 𝑖 (𝑖 = 1..𝑀 ) are used to store correspondences coming from each model 𝑃 𝑖 respectively.

The hypotheses generator is a two-stage procedure: offline registration (cf. The time complexity of the method is about 𝑂(𝑚 + 𝑛), where 𝑚 = ∑ 𝑀 𝑖=1 𝑚 𝑖 is the number of total model points, and 𝑛 is the number of scene points. We use the term "about" since it is not formally proved but is only an estimation and validated empirically. It will be explained separately for each module in Section 4.1.2-4.1.4.

Hypotheses generator

The hypotheses generator is basically a geometric hashing module. During offline registration, each 𝑝-patch (i.e. a set of 𝑘 + 1 points) is registered into a hash table along with its model number. During online generation, the input 𝑞-patch is considered as a query set and the corresponding 𝑝-patches that may belong to different models are retrieved.

We chose to use geometric hashing because it can drastically reduce the search space, it works well with moderate point jitter and in the presence of extra/missing points, and it can deal with many different types of transformations. By restricting it to local patches, the number of points in each 𝑝-patch is small enough for efficient matchings.

Geometric hashing uses a local basis. Let 𝑓 (𝑇 𝐿 ) be the degree of freedom of the

X B 1 = X B 2 ⇔ ∃!𝑇 𝐿 s.t. B 2 = 𝑇 𝐿 (B 1 ) (4.1)
where ∃! represents the unique existential quantification operator meaning "there exists one and only one".

Jitter can influence the value of X as well. If the geometric descriptor of basis B is X, and if B ′ which gives rise to a geometric descriptor X ′ is a jittered version of B , then X ′ is a random variable. Since 𝑇 𝐿 is a linear transformation, we can

analytically find Σ = (Σ 1 , Σ 2 , ..., Σ 𝐷 ) so that 𝑋 ′ 𝑖 ∼ 𝑁 (𝑋 𝑖 , Σ 2 𝑖 )
when the variance of jitter 𝜎 is comparably small to inter-point distance 𝑙. Note that for simplicity, Σ is not a covariance matrix but a vector of length 𝐷. As a consequence we say that B and B ′ with descriptors X and X ′ respectively, are matched under jittered conditions if (4.2) is satisfied:

|𝑋 𝑖 -𝑋 ′ 𝑖 | ≤ 2Σ 𝑖 (4.2)
Given 𝑏 points, there are many ways to construct a basis. For a 𝑝-patch, we choose 𝑝 as the first point (i.e. the origin), other points are chosen so that X lies in a bounded region 𝑅 𝑋 = {X : Given the randomly selected input scene point 𝑞, the generation step begins by creating a list of (𝑝 𝑗 -patch, 𝑞-patch) with all model points 𝑝 𝑗 ∈ ∪𝑃 𝑖 (𝑖 = 1..𝑀 ).

|𝑋 𝑖 | ≤ 1, 1 ≤ 𝑖 ≤ 𝐷}.
They are called pre-hypotheses since points correspondences are unknown in these two patches. Each pre-hypothesis has a local voting table for its neighbor points 4.1. General algorithm correspondences. Then the voting schema detailed in Alg. 6 is applied.

Let us take one pre-hypothesis (𝑝-patch, 𝑞-patch) as an example to show the creation of an hypothesis from the voting results. The process is similar to Alg. 3.

First, all votes less than 𝑣 𝑡 (set to 4) are set to 0. Since the generator is a multi-model GH module, its time complexity is 𝑂(𝑚(𝑘 + 1) 2𝑏 ) in the worst case when all 𝑝-patches are hashed in a few bins [START_REF] Wolfson | [END_REF]],

where 𝑚 is the total number of model points, 𝑘 is the number of nearest neighbors used and 𝑏 is the number of points in 𝐵, which only depends on the local transformation type 𝑇 𝐿 .

Hypotheses validator

The core work of the validator is to check whether the local transformation 𝑇 𝐿 of an input hypothesis is an approximation of 𝑇 . If a paired-patch contains point correspondences of the real transformation, it is called an in-paired-patch and its local transformation is a good approximation of 𝑇 . Otherwise, it is called an outpaired-patch.

Local transformations of two neighboring in-paired-patches should be similar (cf. Section 4.2 for details) as they are both approximations of 𝑇 in the same neighborhood. Such a relationship cannot be established for out-paired-patches.

We use this consensus between neighbor paired-patches' transformations to validate hypotheses. Obviously, the more paired-patches share a consensus with 𝑇 𝐿 , the more likely 𝑇 𝐿 is a good approximation of 𝑇 . The measurement of consensus is 𝑇 𝐿 -dependent, it will be addressed when we discuss concrete implementations, cf. Section 4.2.

The idea of consensus between neighboring paired-patches was already seen in RRDM (Chapter 3), however this approach is different. In RRDM, transformations of neighboring paired-patches are first estimated independently by voting before being compared, which consumes a lot of time. Furthermore, when there are lots of noise, voting results will become unreliable and thus correct point correspondences can hardly be found only via voting. In LGC, if the transformation of one pairedpatch is found by voting, we use it as an initial guess to estimate transformations of its neighboring paired-patches. This allows us to improve the possibility and the speed of finding similar transformations between neighboring paired-patches. In addition, by relying on the consensus and on a proper termination criterion 𝑁 𝑙𝑎𝑟𝑔𝑒 ,

LGC can almost guarantee to find a correct estimation of 𝑇 (if it exists), which will be further improved by the refiner module.

Alg. 7 describes the validator module. Recall that a paired-patch contains correspondences between two patches, thus it can be seen as a small correspondence list as well. The validator takes an hypothesis paired-patch (𝑝 0 , 𝑞 0 ) as input, finds its neighboring paired-patches having a similar local transformation. Since these neighboring paired-patches are in geometric consensus with the hypothesis, they can be seen as supporters of it. They are added to a supporterlist. Then the validator module does the same thing to these newly added paired-patches (i.e. supporters).

This process continues until no paired-patch can be added anymore. The resulting supporterlist is a list of correspondences in the same subregion and in geometric consensus. If the local transformation 𝑇 𝐿 of the input hypothesis is considered as a 66 4.1. General algorithm good approximation of 𝑇 , supporterlist is registered into container 𝐶 𝑖 , 𝑖 being the model number of point 𝑝 0 . If 𝐶 𝑖 is not empty, supporterlist and the stored list are merged to give a larger correspondences list sharing the same 𝑇 .

Algorithm 7 Validator

Input: One hypothesis: paired-patch (𝑝 0 , 𝑞 0 ) Output: A list of correspondences: supporterlist supporterlist = {(𝑝 0 , 𝑞 0 )} for each (𝑝, 𝑞) ∈ supporterlist do

𝑇 𝐿 = local transformation of paired-patch (𝑝, 𝑞)
for all (𝑝 ′ , 𝑞 ′ ) ∈ paired-patch (𝑝, 𝑞) do

Find correspondences 𝑠 𝑝 ′ ,𝑞 ′ of (𝑝 ′ -patch, 𝑞 ′ -patch) by using nearest neighbors between 𝑝 ′ -patch and 𝑇 -1

𝐿 (𝑞 ′ -patch). if 𝑠𝑖𝑧𝑒(𝑠 𝑝 ′ ,𝑞 ′ ) > 𝑏 then Estimate 𝑇 ′ 𝐿 using 𝑠 𝑝 ′ ,𝑞 ′ if 𝑠 𝑝 ′ ,𝑞 ′ and 𝑇 ′ 𝐿 satisfy (4.3) and 𝑇 𝐿 ≃ 𝑇 ′ 𝐿 then supporterlist ← supporterlist ∪{(𝑝 ′ , 𝑞 ′ )}
Time complexity of the validator depends on the number of input hypotheses.

In the worst case, 𝑁 𝑚𝑎𝑥 scene points have been fed to the generator and each scene point generates 𝑚 hypotheses, thus generates 𝑚𝑁 𝑚𝑎𝑥 hypotheses. As a consequence, the time complexity is 𝑂(𝑚𝑁 𝑚𝑎𝑥 ) in this worst case.

Result refiner

After loops of "generator-validator", many correspondences are filled into containers.

If the container 𝐶 𝑖 has the most correspondences, then model 𝑃 𝑖 is considered to appear in the scene. The refiner takes the correspondence list in 𝐶 𝑖 as input, noted as 𝑠𝑙. It works on the scene point set 𝑄 and the model point set 𝑃 𝑖 to find transformation 𝑇 . The overall refining process is similar to Section 3.2.2.3, but Delaunay triangulations and a threshold 𝑑 𝑡 are used to improve its performance.

First of all, Delaunay triangulations are generated for both point sets. A point 𝑝's mesh neighbors are defined as the points who directly connect to 𝑝 by Delaunay mesh (cf. Fig. 4.21). As a consequence, each point is connected with its mesh neighbors by the Delaunay mesh. From 𝑠𝑙, a transformation 𝑇 can be estimated.

Then we try to find matching points for all mesh neighbors of points in 𝑠𝑙. When a new matching is found, the new correspondence is added to 𝑠𝑙. With these correspondences, the transformation can be better estimated. This process continues until 𝑠𝑙 stops growing. Note that at the end of Alg. 8, 𝑇 is the transformation between normalized 𝑃 𝑖 and 𝑄. However, since point correspondence between original The number of nearest neighbors 𝑘 defines the size of local patches, the bigger it is, the more robust the algorithm is against extra/missing points. Since our method is based on Geometric Hashing, the efficiency related to 𝑘 in the generator is 𝑂((𝑘 + 1) 𝑏 ), where 𝑏 is the number of points in a basis B (cf. Sec. 4.1.2). 𝜂 = 𝜎/𝑙 is the jitter factor. When it is larger, the algorithm is more robust against acquisition noise, but each 𝑝-patch will be registered into more hash bins as well. The latter impacts the efficiency of the generator at 𝑂(𝜂 𝐷 ), where 𝐷 is the dimension of the geometric descriptor (cf. Sec. 4.1.2). We found that 𝜂 ≈ 10% is the limit for the jittered point pattern matching problem since beyond that value, precise transformation is difficult to obtain even with pre-known point correspondences (cf. Fig 4.10). Fortunately, we found 𝜂 = 5% -7% to be sufficient for many applications. 𝐸 1 = ∃𝑁 𝑙𝑎𝑟𝑔𝑒 consensus correspondences in a local region between 𝑄 and 𝑃 𝑖 Reducing false alerts means reducing the probability of 𝐸 0 knowing 𝐸 1 , i.e. 𝑃 (𝐸 0 |𝐸 1 ). We can show that 𝑃 (𝐸 0 |𝐸 1 ) ≈ 0 with a very small 𝑁 𝑙𝑎𝑟𝑔𝑒 . In a con- servative estimation (see Section 4.1.6 for demonstration) where one has 10 4 models to track and each model has 10 6 points, 𝑃 (𝐸 0 |𝐸 1 ) is smaller than 10 -6 with 𝑁 𝑙𝑎𝑟𝑔𝑒 = 20. However, this result cannot guarantee that the method can achieve an accurate estimation of 𝑇 when a model is present, because the refiner may find false inliers (cf.

.3 𝑁 𝑚𝑎𝑥 termination rules

The second condition on 𝑁 𝑚𝑎𝑥 means that we are almost sure that no model exists for the point set 𝑄. This is a direct inspiration from RANSAC. If 𝑄 has 𝑡 inliers, and if an in-paired-patch has a probability 𝑝 1 to be proposed to the validator and a probability 𝑝 2 to pass the validation phase, then the probability (𝜆) of no result from the validator after trying 𝑁 image points is:

𝜆 = (1 - 𝑡 𝑛 𝑝 1 𝑝 2 ) 𝑁 (4.4)
If we want the probability of no matching results 𝜆 be smaller than a predefined value 𝜆 𝑚𝑎𝑥 , i.e. 𝜆 < 𝜆 𝑚𝑎𝑥 . We have:

𝑁 𝑚𝑎𝑥 = log(𝜆 𝑚𝑎𝑥 ) log(1 -𝑡 𝑛 𝑝 1 𝑝 2 ) (4.5)
From our experiences, the correct correspondences in an in-paired-patch may be quite difficult to find in very noisy cases. So 𝑝 1 can be very small. However, once an in-paired-patch is found, it is easy to find other in-paired-patches in the neighborhood and thus it will probably pass the validation. As a consequence, 𝑝 2 is large. Let us take 𝑝 1 = 0.2 and 𝑝 2 = 0.7. In a scene with 50% inliers, if we want 𝜆 < 𝜆 𝑚𝑎𝑥 = 0.05, we will get 𝑁 𝑚𝑎𝑥 = 41. So no matter how large the point set in the scene, one needs only to do at most 41 simple geometric hashing queries with 𝑘 + 1 points. In all our experiments, 𝑁 𝑚𝑎𝑥 is set to be 45. One can adjust 𝑁 𝑚𝑎𝑥 to a more proper value according to (4.5) to adapt to different outlier ratios.

Local consensus: a guarantee of low false alert

We want to find probability of false alert 𝑃 (𝐸 0 |𝐸 1 ) mentioned in Section 4.1.5.2.

According to Bayes' theory, we have:

𝑃 (𝐸 0 |𝐸 1 ) = 1 1 + 𝑃 (𝐸 1 |𝐸 0 )𝑃 (𝐸 0 ) 𝑃 (𝐸 1 |𝐸 0 )𝑃 (𝐸 0 ) (4.6)
To have a conservative estimation, we want to find the upper bound of 𝑃 (𝐸 0 |𝐸 1 ).

If the upper bound is small enough, it is less likely to have false alerts. 𝑃 (𝐸 1 |𝐸 0 )

and 𝑃 (𝐸 1 |𝐸 0 ) are two events independent of 𝑃 (𝐸 0 ). So 𝑃 (𝐸 0 |𝐸 1 ) will be larger if 𝑃 (𝐸 0 ) is smaller. 𝑃 (𝐸 0 ) can be seen as the probability that model 𝑃 𝑖 appears in the scene. When there are more models to track, 𝑃 (𝐸 0 ) will be smaller. Considering we have 10 4 models to track, which is a huge number, the probability that 𝑃 𝑖 appears in the scene is roughly 10 -4 . In (4.6), 𝑃 (𝐸 1 |𝐸 0 ) is very large, say 𝑃 (𝐸 1 |𝐸 0 ) ≈ 1.

We have to find 𝑃 (𝐸 1 |𝐸 0 ).

𝐸(𝐿) = 𝑙𝑜𝑔 1/𝑝 𝑏 (𝑁 𝑏 (1 -𝑝 𝑏 )) Σ(𝐿) = 𝜋 √ 6𝑙𝑛(1/𝑝 𝑏 ) (4.7)
where the expectation 𝐸(𝐿) is proportional to 𝑙𝑜𝑔𝑁 𝑏 , it varies little with respect to the length of the sequence. The variance Σ(𝐿) is independent of 𝑁 𝑏 .

If we define random field as a field containing many independent and identically distributed (i. For a given transformation 𝑇 , a point correspondence is considered coherent with 𝑇 if the scene point is back-projected within a 3𝜎 radius neighborhood of the model point (cf. Alg. 8). If the scene point has a random position, the probability of the coherence is about (3𝜎/𝑙) 𝑑 . Remember that 𝑙 is the "inter-point" distance defined in Section 4.1.1. We choose the limit of jittered point pattern matching, 𝜎 = 0.1𝑙 (cf. Section 4.1.5). If at least every two correspondences contain one inlier and 𝑑 = 2, which is a conservative estimation, this results in 𝑝 𝑏 ≈ 0.2. Let 𝑁 𝑏 = 10 6 points, we have 𝐸(𝐿) ≈ 8.4 and Σ(𝐿) ≈ 0.8.

If the transformation 𝑇 needs 𝑏 𝑇 points as basis which gives us 𝑏 𝑇 "free" points in 𝑁 𝑙𝑎𝑟𝑔𝑒 , we have:

𝑃 (𝐸 1 |𝐸 0 ) = 𝑃 (𝐿 ≥ 𝑁 𝑙𝑎𝑟𝑔𝑒 -𝑏 𝑇 ) (4.8)
Take as example a 2D homography (i.e. 𝑏 𝑇 = 4) and 𝑁 𝑙𝑎𝑟𝑔𝑒 = 20, thus 𝑁 𝑙𝑎𝑟𝑔𝑒 -𝑏 𝑇 -𝐸(𝐿) ≈ 7.6 ≈ 9.5Σ(𝐿). [START_REF] Schilling | [END_REF] shows that the distribution of 𝐿 is approximately a normal distribution. As a consequence, 𝑃 (𝐸 1 |𝐸 0 ) = 𝑃 (𝐿 ≥ 𝐸(𝐿) + 9.5Σ(𝐿)) ≪ 10 -10 . Thus 𝑃 (𝐸 0 |𝐸 1 ) ≈ 10 -6 and it is a conservative upper 4.2. Specific Implementations bound estimation.

Specific Implementations

In this section, we present the detailed implementation of the algorithm for 2D homography and 3D similarity. 2D homography is used to match coplanar point sets (markers) while 3D similarity matching is very useful in 3D model registration.

Given the general algorithm, one needs three more pieces of information for a detailed implementation: (a) How to construct the basis 𝐵 and compute the geometric descriptor X (b) How to measure the similarity between 𝑇 𝐿 s. (c) How to find 𝑇 𝐿 or 𝑇 with known point correspondences. We detail each of them for both cases.

2D homography

A 2D homography is basically a perspective transformation, and each point has only two coordinates. So 𝑇 belongs to 2D perspective transformation and 𝑇 𝐿 belongs to 2D affine transformations.

(a) The basis 𝐵 is ordered as (𝑝 0 , 𝑝 1 , 𝑝 2 , 𝑝 3 ) so that Δ𝑝 0 𝑝 1 𝑝 2 has the largest surface among all triangles containing 𝑝 0 . 𝑝 1 and 𝑝 2 are chosen so that --→ 𝑝 0 𝑝 1 × --→ 𝑝 0 𝑝 2 > 0. This gives us a local right-hand basis and |𝑋 1 |, |𝑋 2 | ≤ 1 (cf. Fig. 4.7), same as in [START_REF] Lamdan | [END_REF]]. Thus, if 𝑝 (𝑖) stands for the i-th coordinate of point 𝑝, we have:

• • • • p 0 p 1 p 2 p 3 =X 1 (p 1 -p 0 )+X 2 (p 2 -p 0 )+p 0 X 1 X 2 • • • p 0 p 1 X 1 X 2
A = (𝑎 𝑖𝑗 ) = ( 𝑝 (1)
1 -𝑝

(1) 0

𝑝

(1)

-𝑝

(1) 0

𝑝

(2)

1 -𝑝 (2) 0 𝑝 (2) 2 -𝑝 (2) 0 ) -1 ( 𝑋 1 𝑋 2 ) = A ( 𝑝 (1) 3 𝑝 (2) 3 ) ( Σ 1 Σ 2 ) = ((𝑋 1 + 𝑋 2 -1) 2 + 𝑋 2 1 + 𝑋 2 2 + 1)𝜎 ( 𝑎 2 11 + 𝑎 2 12 𝑎 2 21 + 𝑎 2 22
) (4.9) (b) A 2D affinity is a 2 × 3 matrix. After a singular decomposition of the left 2 × 2 square matrix, a rotation 𝜃 and two scale factors 𝛼 (1) , 𝛼 (2) can be retrieved.

By using the result from Section 3.2.2.2, we say that two affine transformations are similar if the difference of their rotation is less than 10 ∘ and the ratio of each scaling is less than 1.3.

(c) With known point correspondences, 𝑇 𝐿 can be found by the least square method while 𝑇 can be found in various ways including OpenCV's findHomography method.

3D similarity

A similarity is a linear transformation itself, so 𝑇 𝐿 = 𝑇 . (a) A basis 𝐵 is ordered as (𝑝 0 , 𝑝 1 , 𝑝 2 ) so that ‖ --→

𝑝 0 𝑝 1 ‖ ≥ | --→ 𝑝 0 𝑝 2 ‖. --→ 𝑝 0 𝑝 2 defines
the positive direction of the second axis of this local basis (cf. Fig. 4.8). Thus the descriptor X (i.e. the coordinates of 𝑝 2 ) can be expressed as:

𝑋 1 = --→ 𝑝 0 𝑝 2 • --→ 𝑝 0 𝑝 1 ‖ --→ 𝑝 0 𝑝 1 ‖ 2 𝑋 2 = √ (‖ --→ 𝑝 0 𝑝 2 ‖/‖ --→ 𝑝 0 𝑝 1 ‖) 2 -𝑋 2 1 𝑋 3 = 0 (4.10)
It is easy to verify that |𝑋 1 |, |𝑋 2 | ≤ 1 and thus variances can be computed as follows: factor. 𝑅 can be represented as a quaternion q. We say two similarities are equal if the difference of rotation is less than 10 ∘ measured by Φ 3 defined in [START_REF] Huynh | [END_REF]] and the ratio of scales is less than 1.3, that is:

Σ 1 = Σ 2 = √ ‖ --→ 𝑝 0 𝑝 1 ‖ 2 + ‖ --→ 𝑝 0 𝑝 2 ‖ 2 + ‖ --→ 𝑝 1 𝑝 2 ‖ 2 ‖ --→ 𝑝 0 𝑝 1 ‖ 2 Σ 3 = 0
• • • • p 0 p 1 p 2 p 3 =X 1 (p 1 -p 0 )+X 2 (p 2 -p 0 )+p 0 X 1 X 2 • • • p 0 p 1 X 1 X 2 p q p' q'
2 arccos(‖q 1 • q 2 ‖) ≤ 10 ∘ , 𝑚𝑎𝑥(| 𝛼 1 𝛼 2 |, | 𝛼 2 𝛼 1 |) ≤ 1.3
𝑅 = 𝑈 𝑆𝑉 𝑇 𝑠 = 𝑛 ∑ 𝑛 𝑖=1 ‖𝑝 𝑖 -p‖ 2 𝑡𝑟(𝐷𝑆) 𝑡 = q -𝑠𝑅p (4.13)
where 𝑈 , 𝑉 and 𝐷 are the results of singular decomposition of covariance matrix Σ (Σ = 𝑈 𝐷𝑉 𝑇 ), and Fig. 4.10 shows results of robustness studies. We increase the jitter factor 𝜂 to study robustness to acquisition noise which can be due to image processing artifacts, etc. Here, we stick to point sets containing 100 points with no extra or missing points. Clearly, RDM can handle some noise, but is outperformed by both RRDM and LGC.

p = 1 𝑛 𝑛 ∑ 𝑖=1 𝑝 𝑖 q = 1 𝑛 𝑛 ∑ 𝑖=1 𝑞 𝑖 Σ = 1 𝑛 ∑ 𝑖=1 (𝑞 𝑖 -q)(𝑝 𝑖 -p) 𝑇 (4.14) 𝑆 = { 𝑑𝑖𝑎𝑔(1, 1, ...,
LGC is the fastest and most robust. LS represents results directly found with pre-known point correspondences, which can be seen as the result obtained using ground-truth point correspondences. Note that 𝜂 = 10% is about the limit of obtaining a precise matching even with ground-truth point correspondences for jittered point pattern matching.

LGC sometimes outperforms LS (cf. 𝜂 = 7% and 8%), because LGC does not keep correspondences when the reprojection distance between the scene and the model points is greater than 3𝜎 (cf. Alg. 8). As a consequence, LGC may use less noisy point correspondences than LS thus leading Then, we present results with extra points (such as outliers that could be detected by image processing techniques) in Fig. 4.11, results with randomly removed points (which simulates underdetections) in Fig. 4.12, results when points in a region are masked (i.e. occlusion) in Fig. 4.13 and results with respect to the perspective angle of the camera in Fig. 4.14. In every case LGC is the fastest and the most robust.

At last, Fig. 4.15 shows the discriminative capabilities of both RDM and LGC (RRDM works for only one model). RDM is a bit faster when no jitter occurs while both techniques fully discriminate between several models. RDM fails to retrieve the correct model as soon as jitter is involved.

LGC can almost always retrieve the correct one whereas its capability to find an accurate transformation remains the same as in the jitter experiment (Fig. 4.10).

As a conclusion to this theoretical study, our proposal outperforms GH, RDM and our previous work RRDM for both robustness and speed, showing a linear behavior in terms of computation time, whatever the conditions. It has also a great discriminative capability even with large jitter. Next we show that LGC can also be applied to 3D transforms.

3D model registration

We show in this section that our method implemented in 3D, allows for registering models undergoing similarity transformations. Most 3D models are composed of 3D point clouds, forming small surfaces which represent their envelopes. Such point clouds are often obtained from 3D scanners, so they are very dense. As a consequence, the geometric distribution of neighbor points does not contain enough discriminative information for LGC.

We use 3D interest point detectors to find geometrical keypoints in such point clouds. The position of each interest point generally depends on the local topology of the 3D model. So, the spatial distribution of these interest points is more or less random on complex models, which is in favor of LGC's requirements. We perform 3D model registrations by matching detected interest points from these models under similarity transformations. Please remember that we only use point coordinates for matching, neither color information nor normal direction is used.

For illustration purpose, we use the Stanford Bunny with Zhong's Intrinsic Shape

Signatures [Zhong 2009] (ISS) since it has a good repeatability [Tombari 2013b[START_REF][END_REF]]. We first create the model point set with ISS key points: for the 35947 𝜂 is set to 3%. Results of registration are listed in Table 4.1. In all cases, rotation errors are less than 0.1 ∘ while scale errors are less than 0.1%. We can conclude that

LGC has the potential to be also used for 3D similarities estimation.

LGC with additional information

So far, we have only considered pure point matching problems where each point is indistinguishable from another. However, points may have their own properties including color or local feature in textures. We show here that such properties can be used in our algorithm.

Some traditional feature point matching methods (e.g. SIFT, SURF or BRIEF) are based on these different properties of points. They first find rough correspondences by using point properties before applying RANSAC-like methods. Since such methods need high inlier ratios, the rough correspondences should contain a lot of inliers. This implies that point properties should be discriminant enough, which results in high dimensional descriptors. We show that using our method instead of RANSAC greatly reduces the descriptors' dimension.

Results on synthetic point sets

In our experiment we use ORB [START_REF] Rublee | [END_REF]]. The original ORB descriptor is composed of a sequence of 32 bits generated from random pixel comparisons. As a consequence, it is very easy to get lower dimension descriptors from the original ones. Image points are no longer fed to the generator randomly but according to their nearest Hamming distances to model points. Let the length of the descriptor in bits be 𝐿 𝑂𝑅𝐵 . In the generator, we do not process (𝑝 𝑗 -patch, 𝑞-patch) which

Hamming distance is larger than 𝐿 𝑂𝑅𝐵 /2 in Alg. 6, thus reducing the number of false hypotheses. The condition in Alg. 8 is also relaxed to find more corresponding points: we add (𝑝 ′ , 𝑞 ′ ) to the solution list only if the Hamming distance between 𝑝 ′ and 𝑞 ′ is less than 𝐿 𝑂𝑅𝐵 /3 and if |𝑝 ′ -𝑇 -1 (𝑞 ′ )| ≤ 3𝑑 𝑡 . We use these values for illustration purpose, they are not optimized. The original 32 bytes ORB are chopped to have lower dimensional descriptors which are used to compare our method with USAC [Raguram 2013a]. The experiment is repeated 5 times for each case. Results are presented in Fig. 4.17 where we 4.4. Applications Similar experiments are repeated twice for a total of 1753 frames and manual visual assessments are provided for all estimated homographies. Both SURF and

�� ��� ��� ��� ��� ���� � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � � �������������������� ������ ��� ��� ����
LGC match at about 15-18ms/frame while RDM matches at around 50ms/frame. LGC outperforms RDM in all cases. Although SURF performs better than LGC in general, it works badly on cross-word and point-pattern since these two models contain less textures. As RDM can be used for "augmenting everything" [Uchiyama 2011a], we claim our method is a better alternative to RDM for this purpose. Mesh neighbors connect points better than nearest neighbors, the latter giving rise to four "important edges" represented in red.

Let u 0 and u be two points which may not belong to point patterns. u is in the neighborhood of u 0 . We assume that 𝑇 can be expanded by the Taylor theorem in the neighborhood of u 0 :

𝑇 (u) = 𝑇 (u 0 ) + 𝜕𝑇 (u 0 ) 𝜕u Δu + 𝜕 2 𝑇 (u 0 + 𝜆Δu) 𝜕u 2 (Δu) 2 (4.17)
with Δu = u-u 0 and 𝜆 ∈ [0, 1]. The first two terms on the right are the linear local transformation 𝑇 𝐿 (u) that we mentioned before, the last term is a small quantity of the same order as (Δu) 2 , which represents the difference between 𝑇 (u) and 𝑇 𝐿 (u).

If the difference between 𝑇 (u) and 𝑇 𝐿 (u) is not very big, the last term can be easily managed by considering that it represents the "jitter". The difference can be expressed as:

‖𝑇 (u) -𝑇 𝐿 (u)‖ ≤ ⃦ ⃦ ⃦ ⃦ 𝜕 2 𝑇 (u 0 + 𝜆Δu) 𝜕u 2 ⃦ ⃦ ⃦ ⃦ 𝑙 2 (4.18) ‖ 𝜕 2 𝑇 (u 0 +𝜆Δu) 𝜕u 2
‖ depends only on transformation 𝑇 , not on inter-point distance 𝑙. Recall that 0.1𝑙 is the limit of algorithm in terms of jitter (cf. Section 4.3). If the inter-point distance 𝑙 is small enough so that 𝑙 ≤ 0.1‖ 𝜕 2 𝑇 (u 0 +𝜆Δu)

𝜕u 2
‖ -1 , then we have: 4.19) The difference between 𝑇 (u) and 𝑇 𝐿 (u) can be smaller than 0.1𝑙. So the algorithm can theoretically deal with any transformation if the distribution of points is dense enough.

‖𝑇 (u) -𝑇 𝐿 (u)‖ ≤ ⃦ ⃦ ⃦ ⃦ 𝜕 2 𝑇 (u 0 + 𝜆Δu) 𝜕u 2 ⃦ ⃦ ⃦ ⃦ 𝑙 2 ≤ 0.1𝑙 (

Repetitive structures

As our method relies on discriminating local geometric structures, repetitive patterns are an issue. They impact both the "generator" and the "validator" since they work on local patterns, but do not affect the "refiner". If the hypothesis (cf. Sec. 4.1.3) is generated inside a repetitive pattern, the algorithm will probably give a false result;

otherwise both "generator" and "validator" have enough discriminative information to estimate a true correspondence. As the generator selects scene points randomly,

Conclusion

In this chapter, we have presented LGC, an algorithm based on Local Geometric Consensus that can be used for several transformation types, both in 2D and 3D.

Compared to RRDM in Chapter 3 and RDM [Uchiyama 2011b], it is more generic and can be applied to any known 2D/3D transformation. It has the ability to recognize and match several different models, which cannot be solved by RRDM.

As to the case of 2D homography, LGC shows a linear time complexity behavior, which is much faster and more robust than RRDM (quadratic time complexity).

LGC is much more robust than RDM [Uchiyama 2011b] in a noisy circumstance. 88 4.

Conclusion

We have also shown two real-time applications developed by using LGC tracking.

It concludes that wide baseline point pattern matching performed by LGC can be both fast and robust enough for texture-less augmented reality applications. In the following chapter, we will show one more application by applying LGC.

Chapter 5

Defocused projector calibration for projector-camera systems Projectors are important display devices for large scale augmented reality applications. Calibration is one of the most essential elements of projector-camera systems (referred to as ProCam systems in the following), namely determining the intrinsic matrix and the distortion coefficients of cameras and projectors, as well as their relative position. A precise calibration will produce a precise mapping between Chapter 5. Defocused projector calibration for projector-camera systems the physical world and camera/projector images, which permits the augmented information to be projected to their correct positions in SAR. However, precisely calibrating projectors with large focus distances implies a trade-off between practicality and accuracy. People either need a huge calibration board or a precise 3D model [START_REF] Resch | [END_REF]].

With the help of LGC introduced in Chapter 4, we present a practical ProCam calibration method to solve this problem. The user only needs a small calibration board to calibrate the system regardless of the focus distance of the projector. Results show that the Root-Mean-Squared reprojection Error (RMSE) for a 450cm

projection distance are only about 4mm, even though it is calibrated using a small B4 (250 × 353mm) calibration board. This work has been presented at International Symposium on Mixed and Augmented Reality (ISMAR) in 2016.

Related work

There exist two families of methods for projector calibration: (i) calibrating each light stripe (cf.

Inverse camera methods

This kind of methods treats projectors as inverse cameras, so that they can be calibrated using Zhang's method [Zhang 2000]. More precisely, a calibration board with a pre-defined pattern is used to calibrate the camera. In the meantime, the projector projects other patterns onto the calibration board. These projected patterns are used to find camera-projector point correspondences in order to calibrate the projector. Several methods exist for finding such correspondences. Structured-light methods [Zhang 2006, Li 2008] employ projections of series of coded gray-bars (cf.

Fig. 5.3a). Even though methods for defocused projector exists [Li 2014], the board needs to be frozen at the same position by the user for several seconds, which is not convenient. Other methods capture fewer or even a single image for each board

Limitations

In [START_REF] Gao | [END_REF]] and [Audet 2009], the calibration board has to be positioned in the focus zone (where projected images are clear), or the pattern becomes too blurred to be recognized. Moreover, Zhang's method [Zhang 2000] requires the pattern to occupy a large part of the field of view of the projector in order to give a better result. Consequently if a wide-angle projector is focusing at about 2m, the user needs to manipulate at least an A0 size (841 × 1189mm) calibration board, which is a burden. The user can also choose to calibrate the projector at a short distance and then change the focus depending on the application, but this will induce a loss of precision, as showed later in Section 5.5. We also experimentally noticed that for Audet et al.'s method [Audet 2009], projectors' brightness greatly influences markers recognition.

Calibration Method

Ouellet et al. [Ouellet 2008] need to take 3 images for each board position while the projector is successively projecting (1) a white image, (2) a regular pattern, (3) a calibration regular pattern. These images are highly related: (1) is subtracted from (2) to find the homography between the projector and the camera (𝐻 𝑐𝑝 ). If the board moves or ambient lighting changes, detection may fail.

(3) needs 𝐻 𝑐𝑝 to interleave projected points between printed ones. If the board moves, points may interfere. Moreover, a few missing points in the regular pattern can lead to detection failures. Thus the whole process is not very robust. At last, this method detects projected dot centers by back projecting them onto the board plane. When the board is out of focus, back projected dots may be heavily deformed, which can make the center detection difficult and introduce large errors.

Contribution

Compared to [Audet 2009] and [Ouellet 2008], our work has the following advantages which makes it much more practical and usable:

• It can work at large distances without having to manipulate unpractical huge calibration boards, unlike [Audet 2009] • It uses random dot patterns, which are robust to pattern interference and insensitive to lighting, unlike both [Audet 2009] and [Ouellet 2008] • It only needs one image for each board position, thus has a faster recovery from detection failure, unlike [Ouellet 2008] • It gives more stable results for intrinsic matrices of cameras and projectors, compared to [Audet 2009]

Calibration Method

Our method relies on simple manipulations of a calibration board whatever the projector's focus distance is: the user holds a small calibration board before the camera and the projector. If the board is still for a while (≈ 1s), an image is automatically captured for calibration. Once a pre-defined number (𝐾 , cf. Section 5.4) of images are acquired, the system can be calibrated using Zhang's method [Zhang 2000] (cf. In the following, we first introduce some basic notations and the generic "inverse camera" method. Then we present our calibration patterns and how to find feature correspondences. Finally, the whole calibration procedure is presented.

Chapter 5. Defocused projector calibration for projector-camera systems

Basic notations and inverse camera method

The ProCam system is described by the camera's intrinsic matrix 𝐾 𝑐 , its distortion coefficients 𝐷 𝑐 , the projector's intrinsic matrix 𝐾 𝑝 , its distortion coefficients 𝐷 𝑝 , the projector rotation in camera frame 𝑅 and the position of the projector's optical center in camera frame 𝑡. More precisely, the intrinsic matrix is defined by [START_REF] Bradski | [END_REF]] as:

K 𝜀 = ⎛ ⎜ ⎜ ⎝ 𝑓 𝜀𝑥 𝛾 𝜀 𝑢 𝜀 0 𝑓 𝜀𝑦 𝑣 𝜀 0 0 1 ⎞ ⎟ ⎟ ⎠ (5.1)
where 𝜀 stands for 𝑐 or 𝑝, 𝑓 𝜀𝑥 and 𝑓 𝜀𝑦 are focal lengths, (𝑢 𝜀 , 𝑣 𝜀 ) is the position of the principal point, 𝛾 𝜀 is the skew.

The "inverse camera" method starts with two sets of correspondences for each view of the calibration board: 𝑐 𝑐𝑏 = {(x, x (𝑐) )} and 𝑐 𝑐𝑝 = {(y (𝑝) , y (𝑐) )}, where x is an interest point (corner/center of markers, etc.) on the calibration board, x (𝑐) is its image in the camera; y (𝑝) is an interest point in the projected image, y (𝑐) is its image viewed by camera. With 𝑐 𝑐𝑏 obtained from 𝐾(𝐾 ≥ 3) views, the camera can be first calibrated using Zhang's method [Zhang 2000]. Lens distortions in x (𝑐) and y (𝑐) can then be removed.

For the projector, the homography from the board to the camera 𝐻 𝑐𝑏 is found with 𝑐 𝑐𝑏 . The location on the calibration board of projected features y (𝑝) can then be expressed as 𝐻 -1 𝑐𝑏 y (𝑐) , thus a new correspondence set 𝑐 𝑝𝑏 = {(y (𝑝) , 𝐻 -1 𝑐𝑏 y (𝑐) )}

Calibration Method

can be established for each view. At last, the projector can be calibrated with 𝑐 𝑝𝑏 obtained from 𝐾 views.

Calibration pattern

We choose to use randomly distributed circular dots as both physical patterns printed on the calibration board and projected patterns, with the following advantages:

• Circular points have no internal structures unlike ARToolKit markers, they are thus less influenced by defocus or blur or lighting changes.

• Different random dots-based markers can be distinguished, real-time tracked and accurately located by LGC even in the case of partial occlusion and over/under-detections, unlike regular patterns. This makes the method robust against pattern overlapping (cf. Chapter 4).

• Points occupy a small surface: the same area can contain more points than other geometries, so more correspondences can be established The first half of these points (𝑃 𝑏 ) are printed on a paper and attached to a rigid board, the second half (𝑃 𝑝 ) is pre-warped and projected onto the board in order to form 𝑃 𝑜 during the calibration procedure (cf. Section 5.2.3). We set a minimum inter-point distance for all 2𝑁 points to have a more homogeneous distribution, which improves camera calibration according to [Shih 1996]. Since our objective

Defocusing error

the homography 𝐻 𝑐𝑏 because of pinhole camera, rectified views 𝐻 -1 𝑐𝑏 (𝐼) are used for physical points 𝑃 ′ 𝑏 . As we use the pinhole projector model, the relationship between the projector image and the camera image is the homography 𝐻 𝑐𝑝 induced by the calibration board [START_REF] Hartley | [END_REF]]. Therefore, rectified views 𝐻 -1 𝑐𝑝 (𝐼) are used for projected points 𝑃 ′ 𝑝 . This latter approach introduces the defocusing error, which is discussed in the following section. To detect ellipse centers in rectified views, the OpenCV's MSER detector [START_REF] Nistér | Linear Time Maximally Stable Extremal Regions[END_REF]] is used with its default parameters. The initialization step takes ≈ 1s, the manipulation sequence is performed in real-time, the improvement step can process 2 -3 frames per second.

Defocusing error

As mentioned in Section 5.2.3, the center of a circle is not perspective invariant.

This can be dealt with view rectification. But if the transformed circle is defocused, the view rectification cannot compensate the error introduced by the defocus effect. In this section, we study the error in center localization by defocus.

The origin of the defocusing error

An optical schema of a single point projection is shown in Fig. 5.7. Since the depth of field of pinhole cameras/projectors is infinite, a projected dot is never defocused despite of the distance between the calibration board and the projector (cf. the green cone in Fig. 5.7). But when the light follows a lens model inside the projector (cf.

the blue cone in Fig. 5.7), there is only one plane of focus. Therefore, a projected dot is defocused if the calibration board is not exactly on the plane of focus. Therefore, in order to estimate the defocusing error, the lens model should be considered.

In Fig. on the camera's image plane. To detect the centers of s 𝑏 , homography 𝐻 -1

𝑐𝑝 is used to rectify the image 𝐻 𝑐𝑏 (s 𝑏 ). This process implies that a pinhole projector model is used (cf. Section 5.2.3). The rectified spot is 𝐻 -1 𝑐𝑝 𝐻 𝑐𝑏 (s 𝑏 ) = 𝐻 𝑝𝑏 (s 𝑏 ). 𝐻 𝑝𝑏 is the homography from positions on the calibration board 𝑝 𝑏 to those on the projector image plane 𝑝 𝑜 induced by a pinhole projector model. As a consequence, if we think that the projector works inversely, then the rectified 𝐻 𝑝𝑏 (s 𝑏 ) can be regarded as an image of s 𝑏 "seen" by the pinhole inverse-projector. In denoted as 𝐸 𝑟 , and we use (𝐸 𝑜 ,𝐻 -1 𝑝𝑏 (𝐸 𝑟 )) as one correspondence for the projector calibration.

However, if there is no defocus in the optical system (i.e. with the pinhole projector model), the light from s 𝑜 should form a green cone and produces an image s ′ 𝑏 = 𝐻 -1 𝑝𝑏 (s 𝑜 ) on the calibration board 𝑝 𝑏 . If the center of s 𝑜 is denoted as 𝐸 𝑜 , then the correspondence (𝐸 𝑜 ,𝐻 𝑝𝑏 (𝐸 𝑜 )) should be used theoretically without any defocus. Unfortunately, s ′ 𝑏 is physically not observable, so 𝐻 -1 𝑝𝑏 (𝐸 𝑟 ) is used instead of 𝐻 -1 𝑝𝑏 (𝐸 𝑜 ) which introduces a defocusing error 𝐻 -1 𝑝𝑏 (𝐸 𝑟 )-𝐻 -1 𝑝𝑏 (𝐸 𝑜 ). The projection of this defocusing error on the projector's object plane 𝑝 𝑜 is then 𝑒 = 𝐸 𝑟 -𝐸 𝑜 .

Chapter 5. Defocused projector calibration for projector-camera systems

An estimation of the defocusing error

To derive the error 𝑒, a coordinate system is constructed so that the angle between 𝑝 𝑏 and the projector's 𝑧-axis is 𝜙. Let 𝑂 be the origin of the coordinate system, as well as the projector's optical center. The 𝑧-axis intersects 𝑝 𝑏 at point (0, 0, 𝑏), meaning that 𝑏 represents the distance between 𝑝 𝑏 and the projector. Therefore, the equation of the board plane 𝑝 𝑏 is:

𝑝 𝑏 : 𝑧 -𝑥 cot 𝜙 -𝑏 = 0 (5.
𝐶 𝑏 : ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ 𝑥(𝑡, 𝛼) = 𝑥𝑜 𝑧𝑜 𝑡 + 𝑑 2 [ 1 -(1 -2𝑟 𝑖 𝑑 ) 𝑡 𝑧 𝑖 ] cos 𝛼 𝑦(𝑡, 𝛼) = 𝑦𝑜 𝑧𝑜 𝑡 + 𝑑 2 [ 1 -(1 -2𝑟 𝑖 𝑑 ) 𝑡 𝑧 𝑖 ] sin 𝛼 𝑧(𝑡, 𝛼) = 𝑡 (5.5)
where 𝑡 > 0, 𝛼 ∈ [0, 2𝜋) are two parameters for 𝐶 𝑏 .

By using the fact that s 𝑏 is the intersection of 𝐶 𝑏 and 𝑝 𝑏 ; both s 𝑏 and s 𝑟 lie on 𝐶 𝑟 ; 𝐶 𝑟 is an elliptic cone passing through 𝑂; and s 𝑟 is the intersection between 𝐶 𝑟 and 𝑝 𝑜 . We can find the equation of s 𝑟 with parameter 𝛼 ∈ [0, 2𝜋):

s 𝑟 : ⎧ ⎨ ⎩ 𝑥(𝛼) = 𝑥 𝑜 + 𝑧𝑜𝑑𝑐𝑜𝑠𝛼 2𝑏+𝑑𝑐𝑜𝑡𝜙𝑐𝑜𝑠𝛼 ( 1 -𝑥𝑜 𝑧𝑜 cot 𝜙 -𝑏 𝑧 𝑖 -2𝑏𝑟𝑜 𝑧𝑜𝑑 ) , 𝑦(𝛼) = 𝑦 𝑜 + 𝑧𝑜𝑑𝑠𝑖𝑛𝛼 2𝑏+𝑑𝑐𝑜𝑡𝜙𝑐𝑜𝑠𝛼 ( 1 -𝑥𝑜 𝑧𝑜 cot 𝜙 -𝑏 𝑧 𝑖 -2𝑏𝑟𝑜

𝑧𝑜𝑑

) .

( 5.6) We can prove that s 𝑟 is an ellipse centered at (𝑥 𝑟 , 𝑦 𝑟 ), with:

⎧ ⎨ ⎩ 𝑥 𝑟 = 𝑥 𝑜 -𝑧𝑜𝑑 2 cot 𝜙 4𝑏 2 -𝑑 2 cot 2 𝜙 ( 1 -𝑥𝑜 𝑧𝑜 cot 𝜙 -𝑏 𝑧 𝑖 -2𝑏𝑟𝑜 𝑧𝑜𝑑 ) 𝑦 𝑟 = 𝑦 𝑜 (5.7) Thus the error 𝑒 is 𝑒 = - 𝑧 𝑜 𝑑 2 cot 𝜙 4𝑏 2 -𝑑 2 cot 2 𝜙 ( 1 + cot 𝜙 𝑥 𝑜 𝑧 0 - 𝑏 𝑧 𝑖 - 2𝑏𝑟 𝑜 𝑧 𝑜 𝑑 ) (5.8)
where -𝑧 𝑜 is the object distance since 𝑧 𝑜 < 0; 𝑧 𝑖 is the focus distance; 𝑑 is the lens diameter; 𝜙 is the angle between the 𝑧-axis and the board (cf. Fig. 5.7), and usually lies in [45 ∘ , 135 ∘ ] for calibration; 𝑝 𝑏 intersects the 𝑧-axis at (0, 0, 𝑏); 𝑥 𝑜 measures the distance between s 𝑜 's center and the principle point in the projector

Calibration results

plane along the 𝑥-axis.

We can see from (5.8) that with a pinhole projector model (i.e. 𝑑 = 0), the error will be 0. This corresponds to our intuitive knowledge. Moreover, there is no defocusing error when the board is vertical (i.e. 𝜙 = 90 ∘ ). Because in this case, the defocused spot is enlarged uniformly so that 𝐸 𝑜 and 𝐸 𝑟 coincide.

In our experiment, the projector's width 𝑊 = 1920 px, its height 𝐻 = 1080

𝑒 𝑚𝑎𝑥 ≈ 0.52𝑓 𝑑 2 𝑏 2 ≈ 0.29𝑝𝑥
(5.9)

Note that (5.9) can also be used to estimate the smallest admissible board distance for different projectors 𝑏 𝑚𝑖𝑛 = 𝑑 √ 𝑓 2𝑒𝑚𝑎𝑥 . When the lens diameter 𝑑 of a projector is large, the board should be placed further from the projector. This indicates the drawback of our method for projectors with lenses of large sizes.

Calibration results

Our experimental system consists of a LogiTech C270 webcam (working resolution 640 × 480), an Epson EB-1771W projector (resolution 1920 × 1080, lens diameter 8.35mm), a laptop (Intel i7-4510U@2.00GHz CPU, 8GB of RAM) and calibration boards. We use our method to calibrate the system with the projector focusing at various distances, and compare the results with [Audet 2009]. This method is chosen as a reference since it presents competitive results compared to other ones, the source code is available online and it is easy to use.

[ Audet 2009] proposes to use a B4 paper (250×353mm) attached to a rigid board as calibration pattern. In order to ease comparison, we use a calibration board of the same size with 2𝑁 = 200 points. We found that with this size, using a 2mm radius for points with 16mm minimum inter-point distance works well. The radius of projected points is 𝑟 𝑜 = 6px.

For the projector we use, a B4 board covers almost all the projection view only at about 50cm from the projector. For other focus distances, three different sizes of pattern (A2, A1, 2 × 𝐴0) are printed. Marker centers are used in [Audet 2009] as they report to have a better result. For our method, we only use the B4 calibration board to demonstrate that it is much easier to use (cf. Fig. the effective focal length varies according to focus change. In Fig. 5.11, our principal point positions seem more stable as well, but the difference is not significant in t-test (𝑡(5) = 1.35 for SD of 𝑢 and 𝑡(5) = 1.92 for SD of 𝑣, 𝑝 > 0.05). For the camera, our method gives stable results 𝑓 𝑐 = 812.4 ± 0.9 against Audet's 𝑓 𝑐 = 808.1 ± 4.2, (𝑓 𝑐𝑥 , 𝑓 𝑐𝑦 are averaged together).

To make a short summary, our method gives smaller RMSE on calibration data, and more stable intrinsic estimates, despite using a 22.7 times smaller calibration board.

Augmentation evaluation

In the SAR community, people do not care so much about the true value of intrinsic matrices, nor about the RMSE of calibration data, since real augmentation will hardly lie on these calibration points. They care more about how precise some information can be projected in the focus zone of the projector. In this section, we use calibration-independent data to evaluate this effect.

When relying on [Audet 2009] to augment information at large distances, the user may choose to calibrate either at a short focus distance with a small board (B4 5.5. Augmentation evaluation size) with a loss of precision, or to calibrate at a larger distance with more burden (due to the use of potentially huge calibration boards). As mentioned before, our method offers the advantage to be able to calibrate the ProCam system with only a B4 size calibration board, whatever the focus distance of the projector. In order to draw a fair comparison between our method and Audet's one, three different results are compared here: our method and Audet's calibrating at Correct Focus Distances with different size of boards (i.e. focus distances of the projector during calibration and evaluation are the same, denoted in the figures as Audet-CFD), and Audet's method calibrating with a focus set at 50cm with a B4 board (denoted as .

We use two asymmetric circle patterns (cf. Fig. 5.12) to investigate projection errors in the focus zone of the projector. 𝑄 𝑏 contains the coordinates of all the printed black points while 𝑄 𝑤 contains the coordinates of the projected white ones.

According to the distance between the ProCam system and the calibration board, different sizes of these patterns are generated. 𝑄 𝑏 are printed at exact positions and with known sizes so that we know their true physical positions on the paper. 𝑄 𝑤 contains the ground-truth positions at which white points should be projected. They are drawn in an image 𝐼(𝑄 𝑝 ) to be warped and projected later. After projection, the real positions of the projected circles are measured by the camera and compared to their ground truth values to compute projection errors. A detailed version of this evaluation algorithm is presented in Alg. 10. Its inputs are intrinsic parameters of the camera and the projector, their relative position, as well as two evaluation patterns 𝑄 𝑏 and 𝑄 𝑤 . Its output is a list of reprojection errors { ⃗ 𝑒 𝑟 } of 𝑄 𝑤 . Before showing the results, we need to address an issue: how precise can this evaluation algorithm be? The precision depends on the intrinsics of both the camera and the projector . Let us take an example with the camera (the projector would follow the same reasoning): assuming a board parallel to the camera's image plane is positioned at a distance 𝑧 from the camera's origin (cf. Fig. 5.13). A segment of length 𝐿 is measured by the camera as being 𝑙 pixels long, so we have 𝑙 = 𝑓 𝑐 𝑧 -1 𝐿.

Furthermore, if the angle between the board and the camera's image plane is 𝜙, 𝑙 = 𝑓 𝑐 𝑧 -1 𝐿 cos(𝜙).

We thus see that the smaller the focal length, the smaller 𝑙. Assuming 𝑙 = 0.5 pixel (used for illustrative purpose as an approximation of the MSER detector precision), the minimum on board difference we can measure is:

𝑑 𝑚𝑖𝑛 = 𝑧 2𝑚𝑖𝑛(𝑓 𝑐 , 𝑓 𝑝 ) cos 𝜙

∑

e 1 /𝑛 1 -∑ e 2 /𝑛 2 is less than 𝑑 𝑚𝑖𝑛 , we consider it to be numerical noise.

This effect is represented in Fig. 5.14, 5.17, 5.18 as a gray region encompassing our curve, which means that other RMSE inside this region shows no difference from our result due to the measurement limit. Otherwise, an independent significant test (t-test with different variations) is performed against the result of our method: 5.11) We choose 𝑝 = 0.05, which means the probability that 𝐻 0 is wrongly rejected

𝐻 0 : ∑ e 1 /𝑛 1 = ∑ e 2 /𝑛 2 (
H pre P p P' b & P' p L l f c z O c x y Figure 5
.13: Projection of a segment: 𝑂 𝑐 is the camera's origin. A world segment 𝐿 is on a board parallel to the camera's image plane at 𝑧, its image 𝑙 is on the image plane. We have 𝑙/𝐿 = 𝑓 𝑐 /𝑧. is less than 5%. In our system, 𝑓 𝑐 < 𝑓 𝑝 , so we have 𝑚𝑖𝑛(𝑓 𝑐 , 𝑓 𝑝 ) = 𝑓 𝑐 ≈ 810 (cf. Section 5.4).

Focus distance

We first evaluate RMSE at different focus distances (cf. Fig. 5.14). The circle pattern is placed at the center of the focus zone to have the least lens distortions.

Although only calibrated at 50cm, Audet-50cm works well for short focus distances (< 250cm), but the RMSE grow rapidly for large focus distances. Except for 250cm, there is no significant performance difference between our method and Audet-CFD although we have to remind the results of Audet-CFD are obtained with huge cumbersome calibration boards (i.e. 1189 × 1682mm). For our method, the RMSE for a 450cm focus/projection distance are only about 4mm. To have a clear visual difference, we choose 450cm to show the reprojection difference (cf. Fig. 5.15). However, this is not contradictory, since the error difference reaches the limits of the measurement system.

Error distribution

Perspective and depth

At last, we study the influence of perspective angle and depth in the projector's focus zone. The projector is always focusing at 250cm. Fig. 5.17 shows that all methods perform worse under large oblique angles than under small ones. Audet-50cm is significantly worse in many cases while our method and Audet-CFD give almost the same result, yet again at the expense of the size of the calibration board.

Conclusion

This chapter presented a practical method to calibrate a projector-camera system, by using the LGC method introduced in Chapter 4. The method can calibrate a ProCam system despite the focus distance of the projector. It is user-friendly and interactive in real-time. We evaluated the calibration results with calibrationindependent data. The reprojection error of the system is competitive with regard to the state-of-art, but the method only demands a small calibration board and simple manipulations even at large focus distances (more than 2.5m). Results also show that the calibrated intrinsic matrix of the projector is more stable than other state-of-the-art methods [Audet 2009]. The drawback of our method is that it potentially has larger system error for projector with large lenses according to (5.9).

Nevertheless this can be solved by calibrating the system at a slightly larger distance. Since defocus projections are tolerated by our method, the calibration board would still be smaller than state-of-the-art methods [Audet 2009]. In this chapter, we aim at detecting and reconstructing surfaces of revolution (SoR) in a cluttered scene in real-time. SoR can represent the majority of primitives and are quite common in man-made objects due to their ease of production. Since the most crucial step in detecting a SoR is to estimate its rotation axis, we propose a fast and accurate method, which boils down the estimation of the rotation axis to a one dimension search. The main contributions of this chapter are summarized as follows:

• a fast and accurate estimation of SoR axes from 3D data

• a framework for detecting, localizing and reconstructing SoR in a cluttered scene in real-time.

Related work

For semantic 3D understanding, structural modeling using geometric primitives has been investigated in the literature. To reconstruct SoR in real-time, we propose a new iterative method to improve the efficiency while still ensuring a good accuracy. It is inspired by the method of Mara et al. [START_REF][END_REF]] that uses parallel planes to cut a SoR. However, while [START_REF][END_REF]] uses a 2-dimensional brute force search, our method finds the axis by solving a 1-dimensional objective function. The same cutting operation (cf. Section 6.3.2) is repeated almost a hundred times in [START_REF][END_REF]] to reach a solution, but is only required at most 5 times in our method.

Surface of revolution axis estimation

The procedure of reconstructing a SoR in our method is summarized as follows:

1. Estimate an initial rotation axis 𝑙 0 using Pottmann's method [Pottmann 1999] (denoted as P-method afterwards)

2. Improve the rotation axis estimation (cf. A fast and accurate iterative step 2 is our main contribution. Usually, the direction of a rotation axis 𝑙 can be determined by two Euler angles and many existing methods use two dimensional search, such as Mara et al. [START_REF][END_REF]]. We

show that the search space can be reduced to one dimension, and that the objective function is easy to solve. Both of them make the algorithm more efficient. Steps 1 and 3 can be seen as pre-process and post-process respectively.

In this section, we first explain the idea of step 2 and then mathematically prove it. The implementation of step 2 and step 3 are then presented in Section 6.3.

Basic idea

We first show that, given an initial guess 𝑙 0 obtained from the P-method, the search space of the direction of the rotation axis is unidimensional. In Fig. 6.1, we use a cone R to represent a general SoR for the sake of simplicity. The solid line 𝑙 𝑔 is the ground truth of its rotation axis. Line 𝑙 0 represents the initial guess. A plane 𝑃 𝑠 containing 𝑙 𝑔 can be defined by the normal ⃗ 𝑙 𝑔 × ⃗ 𝑙 0 , where the unit vector ⃗ • represents the line's direction. By rotating ⃗ 𝑙 0 inside 𝑃 𝑠 , it is possible to find ⃗ 𝑙 𝑔 . Therefore, the search space is reduced to one dimension. Note that we use a 2D sketch of 𝑃 𝑠 in Fig. 6.2 to represent Fig. 6.1 because all useful information is contained in 𝑃 𝑠 .

We now show that 𝑃 𝑠 can be determined without knowing 𝑙 𝑔 . We first create 𝑛 parallel planes 𝑃 𝑗 (𝑗 = 1..𝑛) with normal ⃗ 𝑙 0 to cut R, which results in a set of 2D curves 𝐶 𝑗 . Since 𝐶 𝑗 are symmetric curves with respect to 𝑃 𝑠 , their centers 𝑂 𝑗 should lie on 𝑃 𝑠 . We then project the centers of all 𝐶 𝑗 on an arbitrary plane, e.g. Next, we derive the objective function of our method. Let us assume that ⃗ 𝑙 is rotated from ⃗ 𝑙 0 inside 𝑃 𝑠 by angle 𝜃 (cf. Fig. 6.1-6.2). The set of parallel planes 𝑃 𝑗 is also rotated so that their normal is ⃗ 𝑙. A line 𝑣 is fitted to all points 𝑂 𝑗 . Intuitively, if ⃗ 𝑙 = ⃗ 𝑙 𝑔 , all 𝐶 𝑗 should be circles because of symmetry, and their centers 𝑂 𝑗 should lie on 𝑙 𝑔 . It means that ⃗ 𝑙 = ⃗ 𝑙 𝑔 ⇒ ⃗ 𝑣 = ⃗ 𝑙 𝑔 = ⃗ 𝑙. In the next section, we show that the converse is true if R is not a sphere. That is to say ⃗ 𝑣 = ⃗ 𝑙 ⇔ ⃗ 𝑙 = ⃗ 𝑙 𝑔 for non-spheres.

P 1 C 1 O 1 C 2 O 2 l l 0 C 3 O 3 P s l g v θ φ O p 2 O p 3 l c
In Fig. 6.2, unit vectors inside 𝑃 𝑠 can be represented by 2D vectors. By abusing notations, we still use ⃗ • to refer to these 2D vectors, so their cross products are real numbers. Finally, the objective function is:

𝜑(𝜃) = sin -1 ( ⃗ 𝑙 × ⃗ 𝑣) with 𝜃 = sin -1 ( ⃗ 𝑙 0 × ⃗ 𝑙) (6.1)
Both 𝜑(𝜃) and 𝜃 are signed angles. Our objective is to find 𝜃 such that 𝜑(𝜃) = 0.

Once ⃗ 𝑙 𝑔 is found, we can use the average of 𝑂 𝑗 as a point on the axis. This gives

𝑙 𝑔 = ( ⃗ 𝑙 𝑔 , 1 𝑛 𝑛 ∑ 𝑗=1 𝑂 𝑗 ) (6.2)

Approximately linear objective function

In this section, we prove that the objective function 𝜑 (6.1) is approximately linear with respect to 𝜃. This is the core of our method since it guarantees its fast convergence.

The unknown 𝑙 𝑔 is chosen to be the 𝑧-axis in Fig. 6.2 for an easier deduction.

The profile of the SoR in this plane can be expressed as 𝑟 = ±𝑓 (𝑧).

𝜙 = sin -1 ( ⃗ 𝑙 𝑔 × ⃗ 𝑙 0 ) (6.3)
Since 𝜃 is the angle from ⃗ 𝑙 0 to ⃗ 𝑙, so when 𝜃 = -𝜙, we have ⃗ 𝑙 = ⃗ 𝑙 𝑔 and thus the objective function 𝜑(𝜃) = 𝜑(-𝜙) = 0. To derive the formula of 𝜑(𝜃), two intermediate angles are used:

𝛼 = sin -1 ( ⃗ 𝑙 𝑔 × ⃗ 𝑙) 𝛽 = sin -1 ( ⃗ 𝑙 𝑔 × ⃗ 𝑣) (6.4)

Surface of revolution axis estimation

It is easy to derive that:

𝜃 = 𝛼 -𝜙 𝜑 = 𝛽 -𝛼 (6.5)
By assuming that 𝛼 is small, the profile near 𝑧 𝑗 can be approximated by its first order Taylor expansion:

𝑟 = ±𝑓 (𝑧) ≈ ± ( 𝑓 (𝑧 𝑗 ) + 𝑓 ′ (𝑧 𝑗 )(𝑧 -𝑧 𝑗 ) ) (6.6)
The equation of plane 𝑃 𝑗 can be expressed as:

𝑟 = tan( 𝜋 2 + 𝛼)(𝑧 -𝑧 𝑗 ) (6.7)
The intersection points (𝑧 + 𝑗 , 𝑟 + 𝑗 ) can be calculated by solving (6.7) and using 𝑟 = +𝑓 (𝑧) from (6.6). Similarly, (𝑧 - 𝑗 , 𝑟 - 𝑗 ) can be calculated from (6.7) and using 𝑟 = -𝑓 (𝑧) from (6.6). The center point 𝑂 𝑗 (𝑧 * 𝑗 , 𝑟 * 𝑗 ) lies in the middle of the [(𝑧 - 𝑗 , 𝑟 - 𝑗 ), (𝑧 + 𝑗 , 𝑟 + 𝑗 )] line segment and thus can be derived as:

𝑧 * 𝑗 = 1 2 (𝑧 - 𝑗 + 𝑧 + 𝑗 ) = 𝑧 𝑗 + 𝑓 (𝑧 𝑗 )𝑓 ′ (𝑧 𝑗 ) cot 2 𝛼 -𝑓 ′2 (𝑧 𝑗 ) 𝑟 * 𝑗 = 1 2 (𝑟 - 𝑗 + 𝑟 + 𝑗 ) = - cot 𝛼𝑓 (𝑧 𝑗 )𝑓 ′ (𝑧 𝑗 ) cot 2 𝛼 -𝑓 ′2 (𝑧 𝑗 ) (6.8)
Since 𝑣 is an interpolation of all 𝑂 𝑗 , its slope 𝑘 in the plane can be expressed as follows:

𝑘 = tan 𝛽 = - 𝑛 ∑ 1 (𝑧 * 𝑗 ) 2 -𝑛 ( 𝑛 ∑ 1 𝑧 * 𝑗 𝑛 ) 2 𝑛 ∑ 1 𝑧 * 𝑗 𝑟 * 𝑗 -𝑛 ( 𝑛 ∑ 1 𝑧 * 𝑗 𝑛 ) ( 𝑛 ∑ 1 𝑟 * 𝑗 𝑛 ) tan 𝛼 (6.9)
Since 𝛼 is small, the equation above can be simplified by the Taylor expansion and by using the approximation tan 𝛼 ≈ 𝛼 along with (6.8):

𝑘 ≈ - 𝜉 𝑗 = 0. As a consequence, (6.10) can be simplified as: 6.11) where

𝑛 𝑛 ∑ 1 𝑧 𝑗 𝑓 (𝑧 𝑗 )𝑓 ′ (𝑧 𝑗 ) - ( 𝑛 ∑ 1 𝑧 𝑗 ) ( 𝑛 ∑ 1 𝑓 (𝑧 𝑗 )𝑓 ′ (𝑧 𝑗 ) ) 𝑛 𝑛 ∑ 1 𝑧 2 𝑗 - ( 𝑛 ∑ 1 𝑧 𝑗 ) 2 𝛼 (6.
𝛽 ≈ - 𝑛 ∑ 1 𝜉 𝑗 𝑓 (𝜉 𝑗 + z)𝑓 ′ (𝜉 𝑗 + z) 𝑛 ∑ 1 𝜉 2 𝑗 𝛼 = -(𝐴 -1)𝛼 (
𝐴 = 𝑛 ∑ 1 𝜉 𝑗 𝑓 (𝜉 𝑗 + z)𝑓 ′ (𝜉 𝑗 + z) 𝑛 ∑ 1 𝜉 2 𝑗 + 1 (6.12)
Thus, by using (6.5) and ( 6.11), the objective function 𝜑(𝜃) can be written as:

𝜑(𝜃) = 𝛽 -𝛼 = -𝐴𝛼 = -𝐴(𝜃 + 𝜙) (6.13)
When ⃗ 𝑙 is rotated to the ground-truth ⃗ 𝑙 𝑔 , then 𝜃 = -𝜙 and thus 𝜑(𝜃) = 0, which corresponds to our intuitive knowledge. In (6.13), 𝜙 is an unknown variable depending on 𝑙 0 . Although 𝐴 depends on 𝜉 𝑗 (cf. (6.12)), by letting Δ𝜉 = 𝜉𝑛-𝜉 1

𝑛-1 , we can show that when 𝑛 → ∞:

𝐴 = 𝑛 ∑ 1 𝜉 𝑗 𝑓 (𝜉 𝑗 + z)𝑓 ′ (𝜉 𝑗 + z)Δ𝜉 𝑛 ∑ 1 𝜉 2 𝑗 Δ𝜉 + 1 ≈ 𝜉𝑚𝑎𝑥 ∫ 𝜉 𝑚𝑖𝑛 𝜉 𝑗 𝑓 (𝜉 𝑗 + z)𝑓 ′ (𝜉 𝑗 + z)𝑑𝜉 𝜉𝑚𝑎𝑥 ∫ 𝜉 𝑚𝑖𝑛 𝜉 2 𝑗 𝑑𝜉 + 1 (6.14)
It shows that 𝐴 is a discrete approximation of an integral. Thus although 𝐴 depends on 𝜉 𝑗 , it does not vary much when 𝑛 is big. In other words, the objective function 𝜑(𝜃) is approximately linear.

We now show that the solution 𝜑(𝜃) = 0 is unique except for spheres and that it corresponds to the rotation axis. From (6.13), we can see that 𝐴 = 0 ⇒ 𝜑(𝜃) ≡ 0.

By solving 𝐴 = 0 from (6.14), we find the profile function is 𝑓

(𝑧) = ± √ (𝐶 -𝑧 2 ),
where 𝐶 is an arbitrary constant. This equation represents a sphere of radius √ 𝐶.

𝜑 ≡ 0 means a sphere can have its rotation axis in any direction, which is obvious.

But for other surfaces, since 𝐴 ̸ = 0, the solution of 𝜑(𝜃) = 0 is uniquely 𝜃 = -𝜙 by (6.13). It means that by finding the unique solution of 𝜑(𝜃) = 0, one can find its 6.3. Implementation details 𝜏 , a point is considered inside a plane if its distance to the plane is less than 0.5𝜏 .

If a SDF representation (with voxel size 𝜏 ) is used, a ray casting method can be used to find points inside a specified plane. Circles are fitted to points in the same plane. We define the geometric error as the distance between a point and the fitted circle in the same plane. A point is considered as an inlier if its geometric error is less than 𝜏 .

Since scanned 3D data contain lots of defects, such as incomplete scanning and noise, some planes may contain data with better quality than others. To compensate those limitations, a circle is considered as valid only when its inlier ratio of circle fitting is more than 80%, and each valid circle 𝐶 𝑗 is associated with a weight 𝑤 𝑗 according to its fitting result. Additionally, planes near the two ends of a SoR are more easily influenced by inaccurate axis direction estimation, so we only use empirically the middle 70% region of the surface for axis estimation. The weight 𝑤 𝑗 for a valid circle is calculated by: 𝑤 𝑗 = Ω 𝑗 ē𝑗 (6.17) where Ω 𝑗 is the angular span of data points on the fitted circle as shown in Fig. 6.3, ē𝑗 is the average geometric error of all inliers in the current plane.

If findCircleCenters finds less than 2 valid circles, it means that the target surface is not a true SoR or the data is too noisy to process. The algorithm fails in this case.

123 Chapter 6. Surface of revolution reconstruction from 3D data 6.3.3 Determining 𝑃 𝑠 and solving 𝜑(𝜃)

The basic idea of calculating 𝑃 𝑠 has been explained in Section 6.2.1. Assuming that the incertitude of center detection is 𝛿, the center projections 𝑂 𝑝 𝑗 in Fig. 6.1 will lie in a thick line of width 2𝛿. By using Principle Component Analysis (PCA) on these center projections, the major eigenvector can be regarded as a good approximation of ⃗ 𝑙 𝑐 (cf. Fig. 6.1) while the minor eigenvalue 𝜆 𝑚𝑖𝑛 is the sum of squared distances between detected circle centers and 𝑃 𝑠 . Thus, 𝛿 can be approximated by: 𝛿 ≈ √ 𝜆 𝑚𝑖𝑛 /𝑛 𝑣 (6.18) where 𝑛 𝑣 is the number of valid circles.

Although 𝜑(𝜃) has a simple form, since the scanned 3D data is noisy and 𝜑(𝜃) is not exactly linear (cf. Section 6.2.2), we still have to rely on an iterative approach to solve 𝜑(𝜃) = 0. Recall Newton's method:

𝜃 𝑖+1 = 𝜃 𝑖 - 𝜑(𝜃 𝑖 ) 𝜑 ′ (𝜃 𝑖 ) = 𝜃 𝑖 - 𝜑(𝜃 𝑖 ) 𝐴 𝑖 (6.19)
For the i-th step, 𝜃 𝑖 is known and 𝜑(𝜃 𝑖 ) is evaluated. Owing to its approximate linearity, 𝜑 ′ (𝜃 𝑖 ) = 𝐴 𝑖 can be approximated by fitting previous 𝑖 + 1 data by using calculateSlope. For the first step 𝑖 = 0, 𝜃 0 = 0, we solve the primitive function (6.12) to find an approximation of 𝐴 0 : 6.20) where 𝐷 is the distance between the two most distant planes (𝑧 𝑚𝑎𝑥 -𝑧 𝑚𝑖𝑛 ). When 𝛼 is small, 𝑃 𝑗 is almost perpendicular to 𝑙 𝑔 , so we use the the radius of fitted circle 𝑅 𝑗 to approximate 𝑓 (𝑧 𝑗 ) (cf. Fig. 6.2) in the second approximation.

𝐴 0 ≈ 3 𝑓 2 (𝑧 𝑚𝑎𝑥 ) + 𝑓 2 (𝑧 𝑚𝑖𝑛 ) -2 𝑛𝑣 𝑛𝑣 ∑ 1 𝑓 2 (𝑧 𝑗 ) (𝑧 𝑚𝑎𝑥 -𝑧 𝑚𝑖𝑛 ) 2 + 1 ≈ 3 𝑅 2 𝑚𝑎𝑥 + 𝑅 2 𝑚𝑖𝑛 -2 𝑛𝑣 𝑛𝑣 ∑ 1 𝑅 2 𝑗 𝐷 2 + 1 (
The iterative process can be terminated when |𝜑(𝜃)| is small enough, e.g. ≤ 10 -3 𝑟𝑎𝑑. However, since circle center detection cannot be more precise than 𝛿 according to (6.18), the estimation of the rotation axis cannot be more precise than tan -1 (𝛿/𝐷). Since 𝛿 ≪ 𝐷, 𝛿/𝐷 can be used instead of tan -1 (𝛿/𝐷). Therefore, the iterative process is terminated when |𝜑(𝜃 𝑖 )| < 𝜀 = max(𝛿/𝐷, 10 -3 ). Practically, we found that this method usually converges in 2-3 steps. However, in the special case where a point cloud is used, 𝜑(𝜃 𝑖 ) may not be continuous and |𝜑(𝜃 𝑖 )| < 𝜀 may 6.4. Real-time SoR reconstruction never be achieved. This leads the process to oscillate around 𝜑(𝜃 𝑖 ) = 0. In order to avoid this, we set the maximum iteration number 𝑖 𝑚𝑎𝑥 = 5. If 𝑖 𝑚𝑎𝑥 is reached but |𝜑(𝜃 𝑖𝑚𝑎𝑥 )| ≥ 𝜀, the 𝑙 𝑖 giving the smallest |𝜑| will be chosen as 𝑙 𝑒 .

Real-time SoR reconstruction

We propose a workflow to reconstruct SoR from dense SLAM data. It uses real-time segmentation result from the work of [Tateno 2015], where surfel-based representation is used. A surfel is a representation of a small surface patch at position ⃗ 𝑝(𝑥, 𝑦, 𝑧)

in 3D space, with its normal direction ⃗ 𝑛, its confidence radius 𝑟 and one segment label. Surfels having the same segment label make up a segment. Usually, each segment represents a single object in the scene. We should mention that, although we use surfel-based representation in this section, our algorithm can be easily adapted to deal with any other SLAM data structure.

Although this workflow can deal with SoR, planes and spheres, we only focus on SoR to highlight our contribution. A method to distinguish among different types of surfaces is presented in Section 6.4.1. It allows us to find segments which represent SoR. The detailed algorithm to reconstruct SoR using [Tateno 2015] as our SLAM data provider at each frame is presented in Section 6.4.2.

Segment classification

We use a two-step method to classify each segment. A surfel at position ⃗ 𝑝 having a normal ⃗ 𝑛 is denoted by (⃗ 𝑝, ⃗ 𝑛). Let ⃗ 𝑟 = (⃗ 𝑝 × ⃗ 𝑛, ⃗ 𝑛) be a 6 dimension vector. For all 𝑁 points belonging to one segment, the following two matrices can be calculated: Although different frames contain different sets of visible surfels, their rotation axis should be similar. Thus, 𝑐 is incremented by 1 if the angle between two estimations of the rotation axis is smaller than a threshold 𝛾 𝑡 (set to 10 ∘ in our experiment).

When 𝑐 = 𝑐 𝑡 (set to 5 in our experiment), the data is considered to be stable. In later frames, we continue to check the precision of this reconstruction. If the reconstruction is precise, the reconstructed model should fit well all visible surfels 𝑆 𝑣 in new frames. Otherwise, the rotation axis of the SoR should be re-estimated and the SoR should be reconstructed by using the newest global model 𝑆 𝑚 . Details on checking the reconstruction precision are shown in Alg. 13, where 𝑒𝑛𝑢𝑚 indicates the number of surfels whose modeling error is larger than 𝑒 𝑡 . Empirically, 𝑒 𝑡 is set to 1.5𝑒, with 𝑒 being the SoR's global reconstruction error (cf. Section 6.3.1). When 𝑒𝑛𝑢𝑚 is larger than a threshold 𝑁 𝑡 , the reconstruction should be improved. Since consecutive frames have large overlapping area, the same surfel may be frequently recalculated. We found that 𝑁 𝑡 = 30%|𝑆 𝑚 | works well in our experiment.

The calcModelingError process estimates the modeling error 𝑒 𝑖 on surfel 𝑠 𝑖 . 𝑒 𝑖 is computed as the distance between 𝑠 𝑖 and its approximate nearest point on the reconstructed surface. The approximate nearest point is constrained so that itself and 𝑠 𝑖 are lying on the same plane which is perpendicular to the rotation axis. Improve axis estimation and reconstruct the SoR from 𝑆 𝑚 Set 𝑒𝑛𝑢𝑚 = 0

Results

In this section, we first evaluate the method for estimating the rotation axis, before using real data to evaluate reconstruction results and the efficiency of our method.

Synthetic study

We use synthetic truncated cone to investigate the feasibility of the method (cf. To simulate noise, 𝛿 𝑟 ∼ 𝑁 (0, 0.5𝜏 ) is added to the three coordinates of each point of the synthetic truncated cone, and 𝛿 𝑛 ∼ 𝑁 (0, 𝜂) is added to each normal component (normalized afterwards). At last, the synthetic cone is rotated and moved by a random rigid transformation. We compare the results of our method with the one from [Pottmann 1999] (denoted as P-method), since we are interested in efficient methods for real-time application. 6.5.2 (cf. Fig. 6.10). We can see that these 𝜑(𝜃)s are almost linear, when 𝜃 is near the ground truth. Fig. 6.7 shows the rotation axis estimation error cos -1 ( ⃗ 𝑙 𝑔 , ⃗ 𝑙) on synthetic data.

Our method gives better results, especially for large noise and opening angles. Moreover, the error of our method is less than 0.5 ∘ in most cases.

Real data

For real-time reconstruction, 𝑛 1 = 𝑛 2 = 10 is used. Two datasets are used for evaluating our online reconstruction algorithm. We use downsampled images to resolution 160 × 120 to let the whole workflow run in real-time on the CPU as indicated in [Tateno 2015].

We use two cylinders and a cone with known size in Dataset 1. The cylinders' radii and the cone's opening angle are used to measure the precision of the algorithm.

Results are shown in Fig. 6.8. It can be seen that the accuracy of cylinders' radii is about 1-2mm, while the accuracy of the cone's opening angle has a big variation at the beginning but gradually converges near the ground truth in about 80 frames. Dataset 2 is an online dataset provided by Tateno 1 . It contains diverse daily necessities lying on a desktop. We draw the time consumption for each frame in Fig. 6.9. The solid line represents the reconstruction time, with an average of 1.87ms per frame, while the dashed line represents the SLAM+segmentation time, with an average of 20.40ms per frame. Therefore, on average our method runs at 45fps. These results are produced on a laptop equipped with Intel Core i7-4510U CPU @2.00GHz and 8GB RAM.

The final reconstructed objects along with the segmentation map is presented in Fig. 6.10. The reconstructed objects correctly match the original data and are visually similar to the original objects in the RGB image.

Conclusion

We have developed a fast and accurate method to estimate the rotation axis of revolution surfaces as well as a real-time online geometric reconstruction workflow.

Experiments show that the proposed method and workflow work well for real world objects, achieving about 1-2mm accuracy for radius estimation and < 1 ∘ for cone's opening angle estimation.

With these reconstructed geometries, point clouds can be simplified and other interesting work may be performed. For example, SLAM wide-baseline re-localization problems could be solved by 3D rotation axis matching. Other future work may also include improving segmentation results by using reconstructed geometries; establishing spatial relationships between different geometries, etc. 7.2. Projector-camera system calibration 7.2 Projector-camera system calibration

The projector-camera system is used in Spatial Augmented Reality (SAR). In large scale SAR, the projector's is usually set at large distance and hence projections at short distance are defocused. This makes the projector-camera system difficult to calibrate under such situations. People either need to use a huge calibration board, acquire a precise 3D model, or wait several minutes for the calibration using traditional methods. An application for calibrating a projector-camera system is developed to solve this problem. It relies on two random point patterns and LGC for calibration. The application is practical and calibrates accurately whatever the projector's focus distance. The user only needs to manipulate a B4 size (250 × 353mm) rigid board during the calibration. The superiority of the method relies on the fact that (1) In a defocused situation, dots can be more easily detected than targets requiring edges information (such as chessboards, ARToolKit markers, etc.).

(2) Random dots can be robustly tracked by LGC.

Contribution: The proposed application is practical and as precise as one representative state of art. It produces calibration results that are stable and accurate with a small calibration board whatever the projector's focus distance. For example, the reprojection error is only about 4mm for a 4.5m projection distance, which is comparable to the state of art method using a 22 times larger calibration board.

Scene understanding in 3D data

Chapter 6 presented a piece of work related to scene understanding from 3D data.

It proposed a novel method to estimate the rotation axis of a Surface of Revolution (SoR) efficiently. It also designed a workflow to integrate real-time SoR reconstruction into dense SLAM. By taking advantage of the symmetry of SoR, it proves that the rotation axis of SoR can be found using a one dimension search. An approximate objective function is also proposed to obtain fast convergences.

Contribution: A fast method to estimate the rotation axis of SoR is proposed, which makes the real-time reconstruction of SoR possible. The work can provide basic elements for high level structure modeling of the scene.

Perspectives

Although some interesting applications have been developed, there are still many improvements that can be considered to make them more useful. Regarding the Chapter 7. Conclusion

projector-camera calibration application, a fully automated method would be of interest in the future and could be achieved by using a non planar calibration target. Experiments on projectors with extreme focal lengths (e.g. very short/long focal lengths) are also interesting research topics. For the augmented map application, a road layer extraction which does not need any user interaction can be beneficial in order to create a fully automated application. Quick map localization in GIS which contains a much wider area than the map is also an interesting topic (cf. [START_REF] Minster | [END_REF]). As to real-time reconstruction of SoR, more experiments can be performed to study the sensitivity of the method with respect to acquisition noise.

More comparisons against other existing but slower methods are needed to study in detail whether our real-time solution has an impact on the precision of reconstructions. Real-time reconstruction of more generic primitives, such as superquadric primitives, is also an interest of research in the future.

As to LGC, two directions can be exploited, i.e. improvements of the algorithm and development of new applications.

Influence on the performance of some parameters is not fully studied, such as 𝑘 (i.e. the number of nearest neighbors). More comprehensive studies on these parameters are required to find their optimal values in different cases, which would help the algorithm to achieve better performance. Similarly, parameter settings on other kinds of transformation can also be studied. The algorithm can be improved to better deal with partially regular point patterns, for example, by selecting carefully points that are fed to the generator. Similar to the "colored LGC" in Section 4.3.3, other ways of combining both geometric information and "color" information for registration can also be exploited in order to gain performance when points have their own characteristics. The problem of massive extra or missing points is another challenge to be faced as well. It often happens when two point patterns to be matched are detected from the same target in different scales. For example, lots of details may be lost in an image with a scale of lower resolution, and thus it produces massive missing points. The current LGC cannot solve this problem.

Since LGC is a point pattern matching algorithm, which is an elementary task in many problems, it may have many other potential applications. Current monocular SLAM techniques rely on texture-based methods to establish point correspondences, and such methods fail if they work in texture-less environments, such as construction sites. By using LGC as the point matching module, texture-less monocular SLAM may be doable. The registration between a scanned 3D model (i.e. a dense point cloud) and its corresponding CAD model (e.g. triangle meshes with sparse points) is difficult. It could also be solved by LGC if 3D geometric salient 7.4. Perspectives points can be found from both models, for example, by using Support Vector Shape (SVS) [ Van Nguyen 2013]. Similarly, a solution to re-localizations of dense SLAM in case of tracking failure may also be proposed based on LGC. Other transformation types can further be considered. By applying mapping functions such as equisolid angle, one may be able to match a planar object with its image captured by a fisheye camera.

Besides LGC and point pattern matchings, other related research directions are interesting as well in the future. Since LGC needs reliable and repeatable point patterns as inputs, studies on stable and repeatable point detectors will be very beneficial. As to texture-less tracking, LGC could be integrated with Virtual Visual Servoing (VVS) [Comport 2006]. New methods on texture-less wide baseline localization can also be developed by using other types of geometries, such as regions and line segments, together with points. Or, a combination of different texture-less registration methods may be interesting in order to propose a general solution.

  1999]), (c) with OST head-mounted display ([Maimone 2013]), (d) with VST hand-held display ([Ridden 2013]), (e) OST display free of head-wore or hand-held device ([Hilliges 2012]), (f) with SAR ([Jones 2014]). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.3 AR displays in different position ([Bimber 2005]) . . . . . . . . . . . 4 1.4 Pin-hole camera model . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.5 Edinburgh: (a)-(b) paper maps of different styles. (c) road network extracted from GIS. . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.6 Some examples of point patterns (a) star matching [Lang 2010], (b)

  not work), lower band for screen-based markers. . . . . . . . . . . . . 47 3.12 Extra point experiment, 𝛽 × 𝑚 extra points are added inside the marker. Solid bars represent correct reprojections, empty bars represent incorrect ones. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 3.13 Tilted view experiment. Solid bars represent correct reprojections, empty bars represent incorrect ones. . . . . . . . . . . . . . . . . . . 49 3.14 Missing point experiment. Solid bars represent correct reprojections, empty bars represent incorrect ones. . . . . . . . . . . . . . . . . . . 50 3.15 Augmenting natural maps of Edinburgh. Maps (a), (b) and (c) are three different maps registered with a single model and a 3D model of the castle is rendered onto them. . . . . . . . . . . . . . . . . . . . 50 3.16 Schema of the application. 𝐼 0 is the first image used for initialization. 𝐻 0 is the homography between the GIS and the map contained in 𝐼 0 52 3.17 Workflow of road layer extraction . . . . . . . . . . . . . . . . . . . . 53 3.18 Classic level set method for map extraction. (a) The initial contour Γ 0 . (b) The final contour. (c) The result of two segments represented in different uniform colors 𝐴 1 (white) and 𝐴 2 (black). . . . . . . . . 54 xviii List of Figures 3.19 (b)-(d) Results of various level sets. The segment 𝐴 1 represents road layer in white, the segment 𝐴 2 represents the background in black. (a) Original map portion. (b) Classic level set method. (c) The presence of road labels disturbs contours. (d) Our modified level set method. . 3.20 Left: an example of rectangular stripes. Right: Initial contours are created from the edges between white stripes and black stripes. . . .

4. 22

 22 Impact of regular patterns. (a) chessboard, 165/165 repetitive feature points. (b) office design, 26/109 repetitive feature points. Blue points are detected while red points are projected model points. . . . . . . . 5.1 Projector's light stripe calibration. Left: the setup of the calibration system. Middle: calibration board. Right: light stripe projected onto the calibration board. [Yamauchi 2008] . . . . . . . . . . . . . . . . . 5.2 An example of setup of two-views based method. The car and its corresponding 3D numeric model is the calibration object. Structuredlight patterns are projected onto the car [Resch 2015]. . . . . . . . . 5.3 Different calibration pattern used in inverse camera methods. (a) Regular dot pattern is used and point correspondences between projector and camera are found by projecting structuredlights [Li 2014]. (b) A physical and a projected regular dot pattern is used [Ouellet 2008]. (c) A physical and a projected chessboard is used [Gao 2008]. (d) Two matrix of ARToolKit markers are used [Audet 2009]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.4 Calibration: only a small board is manipulated whatever the focus distance of the projector. . . . . . . . . . . . . . . . . . . . . . . . . .

3.

  Figure 1.2: Augmentation shown (a) with OST retina display ([Magic Leap 2015]), (b) with VST Head-mounted display ([Kato 1999]), (c) with OST head-mounted display ([Maimone 2013]), (d) with VST hand-held display ([Ridden 2013]), (e) OST display free of head-wore or hand-held device ([Hilliges 2012]), (f) with SAR ([Jones 2014]).

Figure 1

 1 Figure 1.3: AR displays in different position ([Bimber 2005])

Figure 1 . 4 :

 14 Figure 1.4: Pin-hole camera model

Figure 1 . 5 :

 15 Figure 1.5: Edinburgh: (a)-(b) paper maps of different styles. (c) road network extracted from GIS.

Figure 1 . 6 :

 16 Figure 1.6: Some examples of point patterns (a) star matching [Lang 2010], (b) bead coding [Datta 2013], (b) engineering drawing of an apartment.

Figure

  Figure 2.1: Examples of fiducial markers: (a) ARTag [Fiala 2005], (b) ARToolKit [Kato 1999], (c) Circular markers [Naimark 2002], (d) RDM [Uchiyama 2011b].

Fig. 2

 2 Fig.2.2). The vector is then normalized to achieve illumination invariance.

Figure 2

 2 Figure 2.2: SIFT descriptor. Left: Image gradients of 16 × 16 neighborhood of interest points. Right: Sums of the gradient magnitudes in 8 directions of 4 × 4 subregions [Lowe 1999].

  as a statistical sample drawing from GMM. This representation naturally takes jitter into account and converts a hard discrete optimization to a relatively easier continuous optimization. [Granger 2002] (EM-ICP) uses model point set to construct Gaussian Mixtures, where each model point serves as the center of a Gaussian kernel. It treats the scene points as samples drawn from model Gaussian Mixture and treats the correspondence matrix as a hidden variable. It applies the Expectation-Maximization (EM) algorithm: the transformation is fixed and correspondences are estimated from expectation during the E-step; then correspondences are fixed and transformation is updated by maximizing a likelihood function during the M-step. Although EM-ICP is more robust than the original ICP, it is about two times slower for the same number of points. [Jian 2011] models both point sets as Gaussian Mixtures and minimizes the L2-distance between the two Gaussian

Figure 2

 2 Figure 2.3: Random dot markers [Uchiyama 2011b]

Figure 2 . 4 :

 24 Figure 2.4: An intersection on a paper map is an area (represented by the red circle for example). It is difficult to tell the precise location of the intersection.

Figure 2 . 5 :

 25 Figure 2.5: The number of correctly matched LLAH descriptors for each point with noise on coordinates. 21 indicates a full match (each point is associated with 21 descriptors).

Figure 3

 3 Figure 3.1: 𝑃 is the model point set in 2D coordinates, 𝑄 is the observed scene point set, they are related by a 2D homography 𝐻.

Figure 3

 3 Figure 3.2: Three Quad-Point configurations.

Figure 3

 3 Figure 3.3: Possible regions for 𝐹 : only feasible regions are , and .

  the same TSR descriptor and have three points in common, they are the same Quad-Point, therefore Q ′ = Q * ≡ Q.

Figure 3 . 4 :

 34 Figure 3.4: Jittered Quad-Point and its TSR descriptor. Left : Original Quad-Point Q (solid blue line) and jittered Quad-Point Q ′ (red dotted line) from a front view. Right : Original (𝑋, 𝑌, 𝑐) and jittered (𝑋 ′ , 𝑌 ′ , 𝑐) descriptors in descriptor space. Gray region represents the matching region 𝑅(𝑋, 𝑌, 𝑐).

Figure 3 . 5 :

 35 Figure 3.5: Quantization of descriptor value for configuration 𝑐. The subregion with red stripes represents 𝐻 𝑖𝑑𝑥 (𝑋, 𝑌, 𝑐) while subregions with green stripes represent 𝑆 𝑖𝑑𝑥 (𝑋, 𝑌, 𝑐).

3. 3 . 1

 31 Robustness of descriptor and 𝑘, 𝜂, 𝐷 𝑚𝑎𝑥 Random addition or removal of points in the close vicinity of a point can deeply modify its k-nearest neighboring points. We call the points not originally present in the neighborhood of a point false neighbors. Such false neighbors result in the creation of new Quad-Points in the image which necessarily lead to a reduced number of correct matchings between model points and scene points. For illustration purposes, let us consider 𝑘 = 7 which leads to the creation of 𝐶 3 7 = 35 Quad-Points for each point in the pattern. In the presence of 25% additional or missing points (uniformly distributed wrt. the original point pattern)in the image, we can assume that each point will have 2 false neighbors in all its 𝑘 = 7 nearest neighbors. Thus the number of remaining original Quad-Points reduces to 𝐶 3 5 = 10, meaning that about 10 correct matching are left for each point.

Figure 3 . 7 :

 37 Figure 3.7: Vote probability distribution functions (𝜙 𝑖 ,𝜙 𝑜 )

  Experimental results are collected from 1000 runs with 200 points uniformly distributed and 𝑘 = 3. Cumulative distributions for 𝛼 (1) , 𝛼 (2) and Δ𝜃 are illustrated in Fig. 3.8 with 𝜂 = 5%. They are computed from inliers according to (3.18). With 𝛼 𝑡 = 1.3, 𝜃 𝑡 = 10 ∘ , 𝑃 (| ln 𝛼 (1) | ≤ ln 𝛼 𝑡 ) = 91%, 𝑃 (| ln 𝛼 (2) | ≤ ln 𝛼 𝑡 ) = 86.8% and 𝑃 (Δ𝜃 ≤ 𝜃 𝑡 ) = 84.7%, we show that about 𝑃 (| ln 𝛼 (1) | ≤ ln 𝛼 𝑡 ) × 𝑃 (| ln 𝛼 (2) | ≤ ln 𝛼 𝑡 ) × 𝑃 (Δ𝜃 ≤ 𝜃 𝑡 ) = 66.9% of inliers affinity pairs are judged similar in the least favourable case.

For

  Figure 3.8: Cumulative distribution for 𝑃 (| ln 𝛼 (1) | ≤ ln 𝛼 𝑡 ), 𝑃 (| ln 𝛼 (2) | ≤ ln 𝛼 𝑡 ) (top) and 𝑃 (Δ𝜃 ≤ 𝜃 𝑡 ) (bottom)

Figure 3 . 9 :Figure 3

 393 Figure 3.9: Point quantity experiment: robustness and speed with respect to 𝑚. 𝛽 = 15%, 𝜂 = 3%.

Figure 3

 3 Figure 3.11: Video snapshots: upper band for paper markers (the rightest one did not work), lower band for screen-based markers.

(Figure 3

 3 Figure 3.13: Tilted view experiment. Solid bars represent correct reprojections, empty bars represent incorrect ones.

Figure 3

 3 Figure 3.14: Missing point experiment. Solid bars represent correct reprojections, empty bars represent incorrect ones.

  also includes a new method to better extract road network from images of paper maps without using manually generated filters. The work has been presented at the conference on Machine Vision Applications (MVA) in 2015.The framework has three modules: road intersection detection, map-GIS registration and tracking. Reference intersections coordinates are calculated from GIS and cleaned beforehand since there are small roads which will never be present on a city map. During the initialization, three tasks are accomplished on the very first frame 𝐼 0 : (1) Map intersections are detected from real maps by our road intersection detection module (cf. Section 3.5.1); (2) map intersections and reference intersections are matched by the map-GIS registration module using RRDM; (3) SURF key points from the paper map in 𝐼 0 are registered in the tracking module for later use.

Figure 3

 3 Figure 3.16: Schema of the application. 𝐼 0 is the first image used for initialization.𝐻 0 is the homography between the GIS and the map contained in 𝐼 0

  Figure 3.17: Workflow of road layer extraction

Fig. 3 .

 3 Fig. 3.18a). 𝐶 1 and 𝐶 2 are the averages of pixels' values in segments 𝐴 1 and 𝐴 2 respectively. They change dynamically as the contour moves. Fig.3.18c shows the segmentation result, with 𝐴 1 in white represents the map segment and 𝐴 2 in black represents the background .Once we have got the map segment in 𝐼 0 , a modified level set method is applied on this segment since the classic level set does not work well for road extraction.

  Figure 3.18: Classic level set method for map extraction. (a) The initial contour Γ 0 . (b) The final contour. (c) The result of two segments represented in different uniform colors 𝐴 1 (white) and 𝐴 2 (black).

Figure 3 .

 3 Figure 3.21: Data sources. Left: GIS road network (green) and road intersections (red). Right: One raster map.

Figure 3 .

 3 Figure 3.23: Initialization results for two different maps.

Figure 3 .

 3 Figure 3.24: Frame 11, 94 and 167 (from left to right) of tracking results in video 3. Top row presents original augmented frames while bottom row shows zooms to illustrate the precision of the matching.

Chapter 4 .

 4 Local geometric consensus (LGC) many transformations such as similarity transformations in 2D/3D, and affine and perspective transformations in 2D. This work has been presented at International Symposium on Mixed and Augmented Reality (ISMAR) in 2015. We reformulate the matching problem as follows: consider having 𝑀 model point sets (𝑃 1 , 𝑃 2 , • • • , 𝑃 𝑀 ), each point set 𝑃 𝑖 containing 𝑚 𝑖 uniformly randomly distributed points in 𝑑 = (2, 3) dimensions. Let 𝑇 : R 𝑑 → R 𝑑 be a geometric transformation on 𝑑 dimensions of which we only know its type (similarity, homography, etc.). Let 𝑄 be a set of observed scene points, some of which belong to 𝑇 (𝑃 ′ 𝑖 ) where 𝑃 ′ 𝑖 is a subset of 𝑃 𝑖 (𝑃 ′ 𝑖 ⊆ 𝑃 𝑖 ) in which each point may have been slightly transformed (i.e. jitter). Other points of 𝑄 are extra points, i.e. acquisition noise. Our objective is to find 𝑃 𝑖 amongst the 𝑀 known model point sets and to determine the transformation 𝑇 between 𝑃 𝑖 and 𝑄. 4.1 General algorithm Our algorithm takes as inputs several model point sets 𝑃 1 , 𝑃 2 , • • • , 𝑃 𝑀 , a query scene point set 𝑄 and the type of transformation. The output is the matched model number 𝑖 (𝑖 ∈ 1..𝑀 ) and 𝑇 between 𝑃 𝑖 and 𝑄, with 𝑇 belonging to the given transformation type (cf. Fig. 4.1).

Figure 4 . 1 :

 41 Figure 4.1: Input and outputs of our algorithm. 𝜂, 𝑘, 𝑁 𝑙𝑎𝑟𝑔𝑒 and 𝑁 𝑚𝑎𝑥 are the parameters of the algorithm.

Fig. 4 .

 4 Fig.4.3) and online generation. During the offline stage, every model patch is registered in the generator. During the online stage, a scene point 𝑞 is randomly chosen and the 𝑞-patch is fed to the hypotheses generator to build several pairedpatches. These generated paired-patches are called "hypotheses" in the following.The validator takes each hypothesis as input and verifies them. If a hypothesis is judged as containing correct correspondences, the validator produces a list of correspondences and adds them to one of the containers 𝐶 𝑖 . The "generation-validation" process is looped until either one of the two following conditions is satisfied: (1) one container has more than 𝑁 𝑙𝑎𝑟𝑔𝑒 correspondences or (2) 𝑁 𝑚𝑎𝑥 scene points have been sent to the generator. We explain these conditions and their values in Section 4.1.5. Finally, the result refiner finds and improves the final result with the largest correspondences set among all containors 𝐶 1 ..𝐶𝑀 (cf. Fig. 4.4). Note that only the hypotheses generator needs to work on all models whereas the hypotheses validator and the result refiner only use the scene point set 𝑄 and the model 𝑃 𝑖 which contains the model points of the input hypothesis.

Figure 4 . 3 :Figure 4 . 4 :

 4344 Figure 4.3: General schema for offline registration: 𝑝-patches are created for each point in all models and are registered into the generator.

Figure 4 . 5 :

 45 Figure 4.5: An illustration before the refiner module: model points (in red) are projected into the scene. Scene points are in blue, correspondences in 𝑠𝑙 are in black. The convex hull of scene points in the correspondence list C is represented by a green polygon, 𝑜 being its center. Region is near C and thus the projected model point and the scene point almost overlap. Region is far from C and thus the distance between the projected model point and the scene point is larger. A correspondence will be established for point 𝑞 ′ outside of C if a model point is projected within 𝑑 𝑡 = 2 𝑜𝑞 ′ 𝑜𝑞 𝐶 𝜎

  Fig. 4.6).

Figure 4 . 6 :

 46 Figure 4.6: Condition 𝑁 𝑙𝑎𝑟𝑔𝑒 confirms the presence of model but refiner finds outliers so the estimation of 𝑇 is not accurate. It can be seen in the yellow ellipse where points are correctly detected but not matched. Correspondences found in validator are inside the green ellipse and are inliers.

  i.d) random experiments organized in space. The problem of long runs describes the probability that all experiments in a region of the random field have the same experiment results. 𝐸 1 |𝐸 0 describes the probability that the most of experiments (instead of all experiments) in a region of the random field have the same experiment results. Both problems have a similar nature, thus (4.7) can be used for a rough estimation.

Figure 4 . 7 :

 47 Figure 4.7: Local basis 𝐵 (𝑝 0 , 𝑝 1 , 𝑝 2 , 𝑝 3 ) for 2D homography. 𝑝 0 is the origin. x𝑖 are axes of the local coordinate system.

Chapter 4 .

 4 Local geometric consensus (LGC)

  similarity can be decomposed into a rotation matrix 𝑅 and a scale 𝛼

Figure 4 . 8 :

 48 Figure 4.8: Local basis 𝐵 (𝑝 0 , 𝑝 1 , 𝑝 2 ) for 3D similarity. 𝑝 0 is the origin. x𝑖 are axes of the local coordinate system, with x𝑖 ⊥ x𝑗 , (𝑖 ̸ = 𝑗). 𝑋 𝑖 are the coordinates of 𝑝 2 under this local coordinate system, with 𝑋 3 ≡ 0.

  We use the method from[Umeyama 1991] to estimate the similarity 𝑇 with known point correspondences. Knowing 𝑛 point correspondences (𝑝 𝑖 , 𝑞 𝑖 ), 𝑝 𝑖 ∈ 𝑃 , 𝑞 𝑖 ∈ 𝑄, 𝑖 = 1..𝑛, the rotation 𝑅, the scale 𝑠 and the translation 𝑡 minimize the sum of squared error ∑ 𝑛 𝑖=1 ‖𝑞 𝑖 -(𝑐𝑅𝑝 𝑖 -𝑡)‖ 2 are:

Chapter 4 .Figure 4 .

 44 Figure 4.10: Jitter experiment (model 𝑃 contains 100 points). 𝛽=0%, missing=0%, occlusion=0%.

  Figure 4.11: Extra points experiment (model 𝑃 contains 100 points). 𝜂=3%, missing=0%, occlusion=0%.

Figure 4 .

 4 Figure 4.15: Discriminative power. Top: with different number of models. Each model contains 100 points. Bottom: with different jitter factor. 50 models are used with 100 points in each model. Full boxes stand for success, empty boxes for good model found with wrong or imprecise transform. (3.19) is used to judge if a transform is precise.

Figure 4 .

 4 Figure 4.16: Left: the model point set (in blue) detected from the original bunny. Right: an example of scene point set (in blue) detected from a transformed instance.

Figure 4 .

 4 Figure 4.17: Graf series and results. Top: graf figures. Bottom left: result for image pair (left, middle). Bottom right: result for image pair (left, right).

Figure 4 .

 4 Figure 4.19: Tracking results: Solid bars stand for correct matching, empty bars for bad matching (i.e. the algorithm matches the scene to a wrong model, or produces an imprecise matching by visual assessments).

Fig. 4 .

 4 Fig. 4.19 shows a summary of the results for each model during tracking. LGC out-

4. 4 . 2

 42 Figure 4.21: Delaunay mesh neighbors (left) and nearest neighbors (right) of the same point set. Neighboring points (black dots) are connected by a black edge.

Figure 4 .

 4 Figure 4.22: Impact of regular patterns. (a) chessboard, 165/165 repetitive feature points. (b) office design, 26/109 repetitive feature points. Blue points are detected while red points are projected model points.

Contents 5 . 1

 51 Related work . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 5.1.1 Two-views based methods . . . . . . . . . . . . . . . . . . . . 92 5.1.2 Inverse camera methods . . . . . . . . . . . . . . . . . . . . . 93 5.1.3 Limitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 5.1.4 Contribution . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 5.2 Calibration Method . . . . . . . . . . . . . . . . . . . . . . . . 95 5.2.1 Basic notations and inverse camera method . . . . . . . . . . 96 5.2.2 Calibration pattern . . . . . . . . . . . . . . . . . . . . . . . . 97 5.2.3 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 5.3 Defocusing error . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.3.1 The origin of the defocusing error . . . . . . . . . . . . . . . . 100 5.3.2 An estimation of the defocusing error . . . . . . . . . . . . . 102 5.4 Calibration results . . . . . . . . . . . . . . . . . . . . . . . . . 103 5.5 Augmentation evaluation . . . . . . . . . . . . . . . . . . . . . 106 5.5.1 Focus distance . . . . . . . . . . . . . . . . . . . . . . . . . . 109 5.5.2 Error distribution . . . . . . . . . . . . . . . . . . . . . . . . 109 5.5.3 Perspective and depth . . . . . . . . . . . . . . . . . . . . . . 110 5.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112

Fig. 5 .

 5 2)[Yamauchi 2008] or each pixel[Luo 2014] of the projector, and (ii) calibrating the projector as a whole using the pinhole model. The former is very time consuming since each light stripe or pixel needs to be calibrated individually. The latter mainly features two categories of techniques: two-views based methods and inverse camera methods.

Figure 5 .

 5 Figure 5.1: Projector's light stripe calibration. Left: the setup of the calibration system. Middle: calibration board. Right: light stripe projected onto the calibration board. [Yamauchi 2008]

Figure 5 . 2 :

 52 Figure 5.2: An example of setup of two-views based method. The car and its corresponding 3D numeric model is the calibration object. Structured-light patterns are projected onto the car [Resch 2015].

Chapter 5 .Figure 5 . 3 :

 553 Fig. 5.3d) [Audet 2009]. Our method follows a similar approach.

Fig. 5 .

 5 Fig. 5.4).

Figure 5 . 4 :

 54 Figure 5.4: Calibration: only a small board is manipulated whatever the focus distance of the projector.

Figure 5 . 5 :

 55 Figure 5.5: Calibration Patterns: 𝑃 𝑏 in black are printed on a piece of paper. 𝑃 𝑝 in white are projected. Both form the original pattern 𝑃 𝑜 .

Chapter 5 .Figure 5 . 6 :

 556 Figure 5.6: Point patterns used in Alg. 9. From left to right: 𝑃 𝑝 , projected pattern 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 ), 𝑃 ′ 𝑏 and 𝑃 ′ 𝑝 in camera view, board pattern 𝑃 𝑏 . 𝐻 𝑐𝑝 = 𝐻 𝑐𝑏 if 𝑃 ′ 𝑏 and 𝑃 ′ 𝑝 are

Figure 5 . 7 :

 57 Figure 5.7: Projection geometry for defocusing error: 𝑂 is the projector optical center, 𝑝 𝑜 is the projector plane, 𝑝 𝑖 is the plane of focus, 𝑝 𝑏 is the board plane. The blue cone represents the light that follows a lens projector model. s 𝑜 (green) is projected by the lens and forms a defocused spot s 𝑏 on the board. Both the green cone and red cone represent the light that follows a pinhole projector model. If s 𝑜 is projected by a pinhole projector (the green cone), it should form a clear spot s ′ 𝑏 on 𝑝 𝑏 . s 𝑟 (red) is the defocused spot "seen" by an inverse pinhole projector (the red cone).

Figure 5 . 9 :

 59 Figure 5.9: Average and maximum reprojection errors (RMSE): A red point on the curve indicates a significant difference between Audet's method and ours, at 𝑝 < 0.05 using a Student's t-test.

Figure 5 .Figure 5 .

 55 Figure 5.10: Focal length results: our method gives a significantly more stable estimation at 𝑝 < 0.05 using a Student's t-test. It shows the trend of focal length variation.

Figure 5 .

 5 Figure 5.12: Evaluation pattern: (printed) black points are used for localization, (projected) white points are used to measure projection errors.

  from the list of reprojection errors {⃗ 𝑒 𝑟 } in each evaluation image. Considering two series of RMSE e 1 and e 2 , if the difference of the average value

Fig. 5 .

 5 Fig. 5.16 shows error distribution in a front view at 250cm. The circle board is placed at four different places to cover the whole projector view. Both our method and Audet-CFD have small errors while Audet-50cm gives large errors especially for small 𝑦 values. We can notice that Audet-CFD generally has smaller errors than ours in Fig.5.16 while our method shows smaller RMSE at 250cm in Fig.5.14.

Figure 5 .

 5 Figure 5.14: RMSE with different focus distances. Color dots mean that RMSE significantly differ from our result: 𝐻 0 of (5.11) is rejected.

Figure 5 .

 5 Figure 5.15: Reprojection error in rectified views (focusing at 450cm). From left to right: ours, Audet-CFD, Audet-50cm. Red lines are drawn to show points alignment. Circles' diameter is 48mm.

Fig. 5 .Figure 5 .Figure 5 .Figure 5 .

 5555 Fig. 5.18 shows that when further away from the projector, Audet-50cm has a poor performance.

  Alg. 11) (a) Use parallel planes to cut the SoR along 𝑙 0 (b) For each plane: find the center of the intersecting curve (cf. Fig. 6.1) (c) Fit a line to all centers (d) Rotate the cutting planes and goto step 2a if the line is not perpendicular to the planes 3. Reconstruct the SoR by finding its profile (cf. Section 6.3.1)

Figure 6 . 1 :

 61 Figure 6.1: Symmetry of SoR and planes

Figure 6 .

 6 Figure 6.2: Sketch of 𝑃 𝑠 : All useful quantities are inside 𝑃 𝑠

Figure 6 .

 6 Figure 6.3: Angular span of data points on a fitted circle: black points represent data points in plane 𝑃 𝑗 , 𝐶 𝑗 is the fitted circle, Ω 𝑗 is the angular span of data points.

  else 𝑐 = 𝑐 -1

Fig. 6 .

 6 Fig. 6.5). The cone is created with small radius 𝑟 = 20, height ℎ = 10 and grid size 𝜏 = 0.2. Since we only calculate the error on the rotation axis estimation, 𝑛 1 = 10 is used.

Figure 6 . 5 :

 65 Figure 6.5: A synthetic truncated cone: with small radius 𝑟, height ℎ, angular span Ω, and the opening angle Θ.

Figure 6 . 6 :

 66 Figure 6.6: 𝜑(𝜃) for different objects: synthetic cone (right top), Can (right middle) and Cup (right bottom). They are approximately linear.

Figure 6 . 7 :

 67 Figure 6.7: Synthetic result: errors on the direction of the rotation axis. Top: errors with respect to opening angles Ω, with 𝜂 = 0.03; Bottom: errors with respect to the amount of noise 𝜂, with Ω = 120 ∘ .

Figure 6

 6 Figure 6.8: Dataset 1: Measuring objects of known dimension. Dashed lines represent ground truth values. Solid lines represent measurements from reconstruction results.

Figure 6

 6 Figure 6.9: Dataset 2: Time consumption.

Figure 6 .

 6 Figure 6.10: Dataset 2. Top: Original segmented point cloud (left) and our reconstruction result (right). Bottom: visual comparison between RGB image (left), original point cloud (middle) and reconstruction result (right) for Cup and Can objects.

  matching . . . . . . . . . . . . . . . . . . . 133 7.1.2 LGC: a general solution for PPM . . . . . . . . . . . . . . . . 134 Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135

			Contents
	de localisation à grande échelle, nous souhaitons proposer de nouvelles manières d'utiliser l'ensemble des informations de textures et géométriques, ainsi qu'étudier d'autres méthodes basées sur des géométries plus complexes que des points, e.g. les 7.1.1 Projector-camera system calibration . . . . . . . . . . . . . . . . . . 135 2D point pattern 7.2 List of Figures
	lignes, les régions, etc. 7.3 Scene understanding in 3D data . . . . . . . . . . . . . . . . . . . . . 135
	7.4	1.1
	Bibliography	139
	7 Conclusion
		x xvi xv xi

7.1

Point Pattern Matching algorithms . . . . . . . . . . . . . . . . . . .

  𝑃 is the model point set in 2D coordinates, 𝑄 is the observed scene point set, they are related by a 2D homography 𝐻. Quantization of descriptor value for configuration 𝑐. The subregion with red stripes represents 𝐻 𝑖𝑑𝑥 (𝑋, 𝑌, 𝑐) while subregions with green stripes represent 𝑆 𝑖𝑑𝑥 (𝑋, 𝑌, 𝑐). . . . . . . . . . . . . . . . . . . . . . . 36 3.6 An increment of vote for table of 𝑉 𝑝,𝑞 . Left : Neighbor correspondences are established from two matched Quad-Points. Right : Cor-

			List of Figures
	3.4	Jittered Quad-Point and its TSR descriptor. Left : Original Quad-
			. 34
	3.5
			.	26
		2.5	The number of correctly matched LLAH descriptors for each point
			with noise on coordinates. 21 indicates a full match (each point is
			associated with 21 descriptors). . . . . . . . . . . . . . . . . . . . . . 27
		3.1

], (d) RDM [Uchiyama 2011b]. . . . . . . . . . . . . . . . . . . . . . . . . . 14 2.2 SIFT descriptor. Left: Image gradients of 16×16 neighborhood of interest points. Right: Sums of the gradient magnitudes in 8 directions of 4 × 4 subregions [Lowe 1999]. . . . . . . . . . . . . . . . . . . . . . 16 2.3 Random dot markers [Uchiyama 2011b] . . . . . . . . . . . . . . . . 25 2.4 An intersection on a paper map is an area (represented by the red circle for example). It is difficult to tell the precise location of the intersection. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.2 Three Quad-Point configurations. . . . . . . . . . . . . . . . . . . . . 31 3.3 Possible regions for 𝐹 : only feasible regions are , and . . . . . . 33 Point Q (solid blue line) and jittered Quad-Point Q ′ (red dotted line) from a front view. Right : Original (𝑋, 𝑌, 𝑐) and jittered (𝑋 ′ , 𝑌 ′ , 𝑐) descriptors in descriptor space. Gray region represents the matching region 𝑅(𝑋, 𝑌, 𝑐). . . . . . . . . . . . . . . . . . . . . . . . . . . . . responding cells are incremented by one vote in 𝑉 𝑝,𝑞 . . . . . . . . . . 38 3.7 Vote probability distribution functions (𝜙 𝑖 ,𝜙 𝑜 ) . . . . . . . . . . . . 42 3.8 Cumulative distribution for 𝑃 (| ln 𝛼 (1) | ≤ ln 𝛼 𝑡 ), 𝑃 (| ln 𝛼 (2) | ≤ ln 𝛼 𝑡 ) (top) and 𝑃 (Δ𝜃 ≤ 𝜃 𝑡 ) (bottom) . . . . . . . . . . . . . . . . . . . . . 44 3.9 Point quantity experiment: robustness and speed with respect to 𝑚. 𝛽 = 15%, 𝜂 = 3%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.10 Jitter experiment: 𝑚 = 100. . . . . . . . . . . . . . . . . . . . . . . . 46 3.11 Video snapshots: upper band for paper markers (the rightest one did

  3.21 Data sources. Left: GIS road network (green) and road intersections (red). Right: One raster map. . . . . . . . . . . . . . . . . . . . . . . Condition 𝑁 𝑙𝑎𝑟𝑔𝑒 confirms the presence of model but refiner finds outliers so the estimation of 𝑇 is not accurate. It can be seen in the yellow ellipse where points are correctly detected but not matched. Local basis 𝐵 (𝑝 0 , 𝑝 1 , 𝑝 2 , 𝑝 3 ) for 2D homography. 𝑝 0 is the origin. x𝑖 are axes of the local coordinate system. . . . . . . . . . . . . . . . . . Local basis 𝐵 (𝑝 0 , 𝑝 1 , 𝑝 2 ) for 3D similarity. 𝑝 0 is the origin. x𝑖 are axes of the local coordinate system, with x𝑖 ⊥ x𝑗 , (𝑖 ̸ = 𝑗). 𝑋 𝑖 are the coordinates of 𝑝 2 under this local coordinate system, with 𝑋 3 ≡ 0. . 75

	List of Figures	List of Figures
	4.6	4.19 Tracking results: Solid bars stand for correct matching, empty bars
		for bad matching (i.e. the algorithm matches the scene to a wrong
		model, or produces an imprecise matching by visual assessments). . .
		Correspondences found in validator are inside the green ellipse and 4.20 Augmenting CAD drawings of a ragum, an apartment and a kart (see
		are inliers. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . supplementary video). . . . . . . . . . . . . . . . . . . . . . . . . . . 70
	4.7
		4.21 Delaunay mesh neighbors (left) and nearest neighbors (right) of the	73
	4.8	same point set. Neighboring points (black dots) are connected by
		3.22 Intersection detection results (blue points). Top: filter-based method a black edge. Mesh neighbors connect points better than nearest
		used in Section 3.4.3. Bottom: Modified level set (new method). . . . neighbors, the latter giving rise to four "important edges" represented
	4.9	Speed experiment.
				. .	83
	4.18 Models used in this experiment.	From left to right:	(Top)
				84
		xix xx

3.23 Initialization results for two different maps. . . . . . . . . . . . . . . 3.24 Frame 11, 94 and 167 (from left to right) of tracking results in video 3. Top row presents original augmented frames while bottom row shows zooms to illustrate the precision of the matching. . . . . . . . 4.1 Input and outputs of our algorithm. 𝜂, 𝑘, 𝑁 𝑙𝑎𝑟𝑔𝑒 and 𝑁 𝑚𝑎𝑥 are the parameters of the algorithm. . . . . . . . . . . . . . . . . . . . . . . .

4.2

An example of a paired-patch (𝑝, 𝑞) with 𝑘 = 6. 𝑝 and its 6 nearest neighbors are in red. 𝑞 and its 6 nearest neighbors are in blue. Black lines with arrows represent correspondences of the paired-patch. Note that not all points have a correspondence due to missing/extra points, and wrong correspondences could also occur.

. . . . . . . . . . . . .

4.3

General schema for offline registration: 𝑝-patches are created for each point in all models and are registered into the generator. . . . . . . .

4.4

General schema for online matching: Generator (Section 4.1.2) is a geometric hashing module. Validator (Section 4.1.3) validates the input hypothesis and creates a list of correspondences from it. Refiner (Section 4.1.4) computes the final result. . . . . . . . . . . . . . . . . 4.5 An illustration before the refiner module: model points (in red) are projected into the scene. Scene points are in blue, correspondences in 𝑠𝑙 are in black. The convex hull of scene points in the correspondence list C is represented by a green polygon, 𝑜 being its center. Region is near C and thus the projected model point and the scene point almost overlap. Region is far from C and thus the distance between the projected model point and the scene point is larger. A correspondence will be established for point 𝑞 ′ outside of C if a model point is projected within 𝑑 𝑡 = 2 𝑜𝑞 ′ 𝑜𝑞 𝐶 𝜎 . . . . . . . . . . . . . . . . . . . . . . . Top: Ideal conditions without jitter nor extra/missing points. Reprojections for all methods are 100% precise. Middle and Bottom: 𝛽=15%, 𝜂=3%, missing=0%, occlusion=0%. . 77 4.10 Jitter experiment (model 𝑃 contains 100 points). 𝛽=0%, miss-ing=0%, occlusion=0%. . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.11 Extra points experiment (model 𝑃 contains 100 points). 𝜂=3%, miss-ing=0%, occlusion=0%. . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.12 Random missing points experiment (model 𝑃 contains 100 points). 𝛽=0%, 𝜂=3%, occlusion=0%. . . . . . . . . . . . . . . . . . . . . . . 79 4.13 Occlusion experiment (model 𝑃 contains 100 points). 𝛽=0%, 𝜂=3%, missing=0%, perspective=30 ∘ . . . . . . . . . . . . . . . . . . . . . . . 79 4.14 Perspective experiment (model 𝑃 contains 100 points). 𝛽=0%, 𝜂=3%, missing=0%, occlusion=0%. . . . . . . . . . . . . . . . . . . . 79 4.15 Discriminative power. Top: with different number of models. Each model contains 100 points. Bottom: with different jitter factor. 50 models are used with 100 points in each model. Full boxes stand for success, empty boxes for good model found with wrong or imprecise transform. (3.19) is used to judge if a transform is precise. . . . . . . 81 4.16 Left: the model point set (in blue) detected from the original bunny. Right: an example of scene point set (in blue) detected from a transformed instance. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.17 Graf series and results. Top: graf figures. Bottom left: result for image pair (left, middle). Bottom right: result for image pair (left, right). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . appartment-plan, city-map, person, advertisement; (Bottom) postcard, cross-word, graffiti, point-pattern. . . . . . . . . . . . . . . . . in red. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

  Projection geometry for defocusing error: 𝑂 is the projector optical center, 𝑝 𝑜 is the projector plane, 𝑝 𝑖 is the plane of focus, 𝑝 𝑏 is the board plane. The blue cone represents the light that follows a lens projector model. s 𝑜 (green) is projected by the lens and forms a defocused spot s 𝑏 on the board. Both the green cone and red cone represent the light that follows a pinhole projector model. If s 𝑜 is projected by a pinhole projector (the green cone), it should form a clear spot s ′ 𝑏 on 𝑝 𝑏 . s 𝑟 (red) is the defocused spot "seen" by an inverse pinhole projector (the red cone). . . . . . . . . . . . . . . . . . . . . 101 .18 RMSE with various depths (focusing at 250cm). Color dots mean RMSE significantly differ from our result: 𝐻 0 of (5.11) is rejected. . Sketch of 𝑃 𝑠 : All useful quantities are inside 𝑃 𝑠 . . . . . . . . . . . . 𝑃 𝑗 , 𝐶 𝑗 is the fitted circle, Ω 𝑗 is the angular span of data points. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . cone: with small radius 𝑟, height ℎ, angular span Ω, and the opening angle Θ. . . . . . . . . . . . . . . . . . . . . Registration results . . . . . . . . . . . . . . . . . . . . . . . . . . . . Different surface types according to 𝑓 𝑡 and 𝑓 𝑟 . . . . . . . . . . . . . Registration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.1.1 Vision-based localization . . . . . . . . . . . . . . . . . . . . .

		List of Figures
	5.7
	6.1	Symmetry of SoR and planes . . . . . . . . . . . . . . . . . . . . . . Introduction
	6.2 3.1	Qualities of detections in map tracking experiment (Mean value ±
	6.3	Angular span of data points on a fitted circle: black points represent
	5.8 Contents Comparison of a 2 × 𝐴0 calibration board used by [Audet 2009] and 3.2 Results of the "Map tracking experiment". Good matching rates is data points in plane 6.4 Workflow of real-time surface reconstruction for SoR. . . . . . . . . . 6.5 4.1 6.1 A synthetic truncated 6.6 shown. Detection time and matching time per frame are reported in the form of "Mean value ± Standard Deviation". . . . . . . . . . . . 1.1
		xxi

From left to right: 𝑃 𝑝 , projected pattern 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 ), 𝑃 ′ 𝑏 and 𝑃 ′ 𝑝 in camera view, board pattern 𝑃 𝑏 . 𝐻 𝑐𝑝 = 𝐻 𝑐𝑏 if 𝑃 ′ 𝑏 and 𝑃 ′ 𝑝 are well aligned on board. . . . . . . . . . . . . . . a 𝐵4 calibration board used by our method (bottom-right) for large focus distances (≥ 250cm). . . . . . . . . . . . . . . . . . . . . . . . . 104 5.9 Average and maximum reprojection errors (RMSE): A red point on the curve indicates a significant difference between Audet's method and ours, at 𝑝 < 0.05 using a Student's t-test. . . . . . . . . . . . . . 105 5.10 Focal length results: our method gives a significantly more stable estimation at 𝑝 < 0.05 using a Student's t-test. It shows the trend of focal length variation. . . . . . . . . . . . . . . . . . . . . . . . . . . 105 5.11 Principal point position results: there is no significant difference between our method and Audet's one, despite the size difference between the calibration patterns. . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.12 Evaluation pattern: (printed) black points are used for localization, (projected) white points are used to measure projection errors. . . . 107 5.13 Projection of a segment: 𝑂 𝑐 is the camera's origin. A world segment 𝐿 is on a board parallel to the camera's image plane at 𝑧, its image 𝑙 is on the image plane. We have 𝑙/𝐿 = 𝑓 𝑐 /𝑧. . . . . . . . . . . . . . 109 5.14 RMSE with different focus distances. Color dots mean that RMSE significantly differ from our result: 𝐻 0 of (5.11) is rejected. . . . . . . 110 5.15 Reprojection error in rectified views (focusing at 450cm). From left to right: ours, Audet-CFD, Audet-50cm. Red lines are drawn to show points alignment. Circles' diameter is 48mm. . . . . . . . . . . . . . 110 5.16 Error distribution in projector's view (focusing at 250cm). Reprojection errors (instead of RMSE) are used in the t-test: Color circles indicate significant differences from our result. . . . . . . . . . . . . . 111 5.17 RMSE with various rotations (focusing at 250cm). Color dots mean RMSE significantly differ from our result: 𝐻 0 of (5.11) is rejected. . 111 xxii List of Figures 5𝜑(𝜃) for different objects: synthetic cone (right top), Can (right middle) and Cup (right bottom). They are approximately linear. . . . . 6.7 Synthetic result: errors on the direction of the rotation axis. Top: errors with respect to opening angles Ω, with 𝜂 = 0.03; Bottom: errors with respect to the amount of noise 𝜂, with Ω = 120 ∘ . . . . . 6.8 Dataset 1: Measuring objects of known dimension. Dashed lines represent ground truth values. Solid lines represent measurements from reconstruction results. . . . . . . . . . . . . . . . . . . . . . . . 6.9 Dataset 2: Time consumption. . . . . . . . . . . . . . . . . . . . . . 6.10 Dataset 2. Top: Original segmented point cloud (left) and our reconstruction result (right). Bottom: visual comparison between RGB image (left), original point cloud (middle) and reconstruction result (right) for Cup and Can objects. . . . . . . . . . . . . . . . . . . . . xxiii List of Tables standard deviation). . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.1.2 Narrow baseline tracking . . . . . . . . . . . . . . . . . . . . 7 1.1.3 Wide baseline localization . . . . . . . . . . . . . . . . . . . . 7 1.2 Thesis motivation . . . . . . . . . . . . . . . . . . . . . . . . . 8 1.2.1 Finding point correspondences . . . . . . . . . . . . . . . . . 8 1.2.2 Problem risen from augmented maps . . . . . . . . . . . . . . 9 1.2.3 Point pattern matching and its difficulties . . . . . . . . . . . 10 1.3 Thesis statement and contribution . . . . . . . . . . . . . . . 11 1.4 Thesis outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Augmented Reality (AR) is a technology which enhances the real-world environment in real-time by adding some extra information, usually generated by computers.

  Choi 2006]. LCS aims at finding the transformation 𝑇 which maximizes the cardinality of the subset 𝑄 𝑠 from 𝑄 and keeps the Hausdorff distance

	1.6c). At last and most important,
	points can be a very good connection between "images" in different modes, such
	as road intersections from paper maps and those retrieved from GIS. Therefore, a
	fast and robust Point Pattern Matching (PPM) method can be a good approach for
	geometric (or texture-less) wide baseline localization for AR. Although sometimes
	point may have their own characteristics such as colors, sizes, etc. Pure geometric
	points with only coordinates are considered here because of their generality.
	Point pattern matching, or point set registration, is the process of finding a
	spatial transformation that aligns two point sets. Let us consider two point sets
	𝑃 and 𝑄, where 𝑃 is the model point set and 𝑄, called the scene point set, is an
	observation of 𝑃 . PPM aims at determining the pair-wise point correspondences
	between 𝑃 and 𝑄, as well as finding the transformation 𝑇 which brings such cor-
	respondences satisfactorily close. PPM has another description known as Largest
	Common Set (LCS) [

  𝑑𝑖𝑠𝑡(𝑇 (𝑃 ), 𝑄 𝑠 ) under a specific small threshold. Due to acquisition noise, points in 𝑄 may not originally exist in 𝑃 (i.e. over-detection or extra points) and points in 𝑃 may not appear in 𝑄 (i.e. under-detection or missing points). Moreover, even though there is neither over-detection nor under-detection and the ground truth transformation 𝑇 is given, the model point set after transformation 𝑇 (𝑃 ) may not yet lie exactly on their corresponding points in 𝑄 because of noisy point positions.

  been presented at the symposium on Virtual Reality Software and Technology (VRST) in 2014.

	In this chapter, we focus on overcoming the limitations of RDM mentioned in
	Section 2.3.4.2. By considering PPM between one model point set and one scene

Contents 3.1 A robust descriptor . . . . . . . . . . . . . . . . . . . . . . . . 30 3.1.1 Definition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.1.2 Affine invariance . . . . . . . . . . . . . . . . . . . . . . . . . 32 3.1.3 Point jitter and descriptor variance . . . . . . . . . . . . . . . 33 3.2 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 3.2.1 Offline pre-registration . . . . . . . . . . . . . . . . . . . . . . 36 3.2.2 Local Voting and Coherency . . . . . . . . . . . . . . . . . . 37 3.3 Choice of parameters . . . . . . . . . . . . . . . . . . . . . . . 41 3.3.1 Robustness of descriptor and 𝑘, 𝜂, 𝐷 𝑚𝑎𝑥 . . . . . . . . . . . . 41 3.3.2 Threshold 𝑣 𝑡 . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 3.3.3 Influence of point jitter on A and 𝛼 𝑡 , 𝜃 𝑡 . . . . . . . . . . . . 43 3.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.4.1 Synthetic images . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.4.2 Real markers . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 3.4.3 Natural map tracking . . . . . . . . . . . . . . . . . . . . . . 48 3.5 Application: Augmented Maps . . . . . . . . . . . . . . . . . 51 3.5.1 Intersection detection on real maps . . . . . . . . . . . . . . . 51 3.5.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 3.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 point set, we propose a robust matching algorithm which is able to register automatically extracted road intersection from paper maps with GIS data. This work Chapter 3. Robust random dot markers (RRDM) has

  Algorithm 1 Offline model points pre-registration for 𝑝 ∈ 𝑃 do

Extract QS(𝑝) according to

(3.13)

for

Q 𝑚 ∈ QS(𝑝) do Calculate t(Q 𝑚 ) from (3.2) Register Q 𝑚 in all hash bins in 𝑆 𝑖𝑑𝑥 according to (3.14)

3.2.2 Local Voting and Coherency

This stage can be decomposed into three phases: (i) pre-alignment by local voting, (ii) rejection by local coherency and (iii) recovery. During the local voting phase, we establish the pre-alignment by computing a score for each correspondence between

  These matching results are recorded into voting tables. 𝐻 𝑖𝑑𝑥 is calculated as follows, where[•] From (3.10), (3.14) and (3.15), we are sure that all model Quad-Points which can be matched with Q 𝑠 are contained in QM. each tentative correspondence (𝑝, 𝑞), we use a 𝑘 × 𝑘 matrix to record their kneighbor's corresponding relations, denoted by 𝑉 𝑝,𝑞 . These relations are established from two matched Quad-Points as presented in Fig.3.6. The detailed process is described in Alg. 2. As each matrix only represents the correspondences between two subsets, we call it local voting.

	Chapter 3. Robust random dot markers (RRDM)
	For p 2
	p 4
	p
	p 5
	q 1
	represents a floor function:
	(3.15)

𝐻 𝑖𝑑𝑥 = ([𝑊 𝑋], [𝑊 𝑌 ], 𝑐)

  represented by the sum of the votes in its voting table. With this measure, a PROSAC-like algorithm can be implemented. However, this approach does not suit our problem as we will show in the experiments of Section 3.4.𝑠 = ⟨𝑞; 𝐹 𝑠 , 𝑅 𝑠 , 𝐿 𝑠 ⟩ ∈ QS(𝑞) do Calculate t(Q 𝑠 ) from (3.2) and 𝐻 𝑖𝑑𝑥 (t(Q 𝑠 )) from (3.15) Retrieve QM at 𝐻 𝑖𝑑𝑥 from hash tableLet us take one tentative correspondence (𝑝, 𝑞) as an example to explain how to find its affinity. (cf.Step 1 of Alg. 3). We first find its local correspondences 𝑠 𝑝,𝑞 by solving a maximum assignment problem. When the number of correspondences 𝑠𝑖𝑧𝑒(𝑠 𝑝,𝑞 ) is larger than 4, 𝑠 𝑝,𝑞 should contain correspondences from at least two different Quad-Points. This cross validation between Quad-Points helps in reducing false matchings. Thus, we only consider such kind of tentative correspondence and estimate an affinity A 𝑝,𝑞 from 𝑠 𝑝,𝑞 . Moreover, we assume a statement 𝐸 𝑠 :

	Algorithm 2 Local voting procedure for 𝑞 ∈ 𝑄 do Extract QS(𝑞) according to (3.13) 3.2.2.2 Rejection This phase aims at rejecting outliers. It begins by clearing votes that are lower or equal to a threshold 𝑣 𝑡 (which value is discussed in Section 3.3.2): for Q 3.2. Algorithm 𝑉 𝑝,𝑞 (p, q) = 0, if 𝑉 𝑝,𝑞 (p, q) ≤ 𝑣 𝑡 (3.16)

3.6: An increment of vote for table of 𝑉 𝑝,𝑞 . Left : Neighbor correspondences are established from two matched Quad-Points. Right : Corresponding cells are incremented by one vote in 𝑉 𝑝,𝑞 .

Since each scene Quad-Point Q 𝑠 has to be tested against each Quad-Point in QM, the time complexity is 𝑂(𝑘 3 𝑛𝑟), where the redundancy number 𝑟 is the average size of QM across all hash entries. After local voting, the "quality" (cf.

[START_REF] Chum | [END_REF]

) of each tentative correspondence can be for all Q 𝑚 = ⟨𝑝; 𝐹 𝑚 , 𝑅 𝑚 , 𝐿 𝑚 ⟩ ∈ QM do 𝑉 𝑝,𝑞 (𝐹 𝑚 , 𝐹 𝑠 )++, 𝑉 𝑝,𝑞 (𝑅 𝑚 , 𝑅 𝑠 )++, 𝑉 𝑝,𝑞 (𝐿 𝑚 , 𝐿 𝑠 )++

  After this procedure, every correspondence in the same affinity group 𝑔 can be described by the same homography 𝐻 𝑔 . Since outliers are randomly matched, it is almost impossible to estimate the same homography for several outliers. Therefore, we consider that the largest affinity group contains inliers and call it the solution list (denoted as 𝑠𝑙 in Alg. 4). More correspondences will be added to this solution Choice of parameters homography between 𝑃 and 𝑄. The detailed procedure of the recovery phase is shown in Alg. 4.

	Algorithm 4 Recovery
	𝑠𝑙 = largest 𝑔 ∈ 𝐺
	Algorithm 3 Rejection
	Step 1. Calculate Affinity for all (𝑝, 𝑞) ∈{All tentative correspondences} do
	list in the following section.
	3.2.2.3 Recovery

𝑠 = Solution of maximum assignment problem for 𝑉 𝑝,𝑞 𝑠 𝑝,𝑞 = 𝑠 -{(p, q)} , ∀𝑉 𝑝,𝑞 (p, q) = 0 Rejection: mark (𝑝, 𝑞) as outlier if 𝑠𝑖𝑧𝑒(𝑠 𝑝,𝑞 ) ≤ 4 Estimate A 𝑝,𝑞 using 𝑠 𝑝,𝑞 Rejection: mark (𝑝, 𝑞) as outlier if 𝐸 𝑠 of (3.17) fails

Step 2. Create a set of affinity groups (𝐺) Group tentative correspondences having similar affinity (3.18) into the same affinity group 𝑔. All 𝑔s make up a set of affinity groups 𝐺. (𝑔 ∈ 𝐺)

Step 3. Estimate homography for each affinity group (𝑔) for all 𝑔 ∈ 𝐺 and 𝑠𝑖𝑧𝑒(𝑔) > 4 do

Estimate 𝐻 𝑔 using RANSAC Discard (𝑝, 𝑞) ∈ 𝑔, if |𝑝 -𝐻 -1 𝑔 (𝑞)| > 2𝜎

We refine the solution list by gradually adding false outliers rejected during Step 2. This increases the number of inliers and results in a better estimation of the 40 3.3. repeat Estimate homography 𝐻 using 𝑠𝑙 with least-square method.

  𝑣 𝑡 Before applying Alg. 3, correspondences that have received less than or equal to 𝑣 𝑡 votes in 𝑉 𝑝,𝑞 have been rejected. The idea behind this thresholding mechanism is to keep as many inliers as possible while rejecting potential outliers. A brief study of the effectiveness of this threshold is now presented.Let 𝑍 𝑖 (𝑣) and 𝑍 𝑜 (𝑣) be two random variables representing the votes respectively received by an inlier and an outlier in voting tables. We find their probability distribution functions 𝑍 𝑖 ∼ 𝜑 𝑖 , 𝑍 𝑜 ∼ 𝜑 𝑜 experimentally by analyzing 200 randomly distributed points 1000 times with 𝑘 = 7, 𝐷 𝑚𝑎𝑥 = 0.2 and 𝜂 = 5%

Using these parameters values, an important difference in the number of votes received by inliers and outliers is presented in Fig.

3

.7. 𝜙 𝑖 (1) ≈ 0.03% and 𝜙 𝑜 (1) ≈ 17.3% indicate that 0.03% inliers and 17.3% outliers receive 1 vote. The best choice is 𝑣 𝑡 = 8, since 𝑣 𝑡 ≥ 9 will not reject more outliers (in proportion) than inliers.

  Synthetic images are constructed as a virtual marker containing 𝑚 randomly distributed points in a 1280 × 720 rectangle. It serves as model point set 𝑃 . This virtual marker is then contaminated by point jitter which is simulated by additive Gaussian noise, and by 𝛽 × 𝑚 randomly generated noisy points. This contaminated virtual marker is placed before a virtual camera so that the angle between the virtual 3.4. Results marker plan and the focal plane of the virtual camera is 30 ∘ . 𝑄 is the camera image of the point set in the contaminated virtual marker. Remember that Gaussian noise is controlled by 𝜂 (cf. Section 3.2).

Table 3

 3 (remember that 𝑃 was constructed based on map (𝑎)) thus preventing some inliers to be reprojected close enough to their corresponding detected points. This prevents the result list from a significant expansion during the recovery phase (cf. Section 3.2.2.3), as a consequence, correct reprojections are confined in a subregion of the map (see top left area in Fig.3.15c).

	Map ID	Point pattern size	Noisy points	Inliers found
	(𝑎)	236.7 ± 14.7	52.7 ± 5.0	102.0 ± 28.1
	(𝑏)	279.6 ± 16.9	114.9 ± 11.9	76.8 ± 19.8
	(𝑐)	237.1 ± 10.3	51.7 ± 7.4	56.9 ± 5.5
	Tracking results of RRDM for the videos of maps (𝑎), (𝑏) and (𝑐) are presented
	in Table 3.2. As both the total number of inliers and the ground truth of the
	homography are unknown, it is difficult to estimate whether a reprojection is correct.
	We consider a reprojection to be good if the basement of the 3D model of Edinburgh

.1: Qualities of detections in map tracking experiment (Mean value ± standard deviation).

to its uniform road color (white) and high image quality (originally 1171 × 795 pixels). Roads in map (𝑏) are white and yellow, those colors are similar to some background elements, leading to a much noisier detection result. Map (𝑐) is also challenging since its original quality is low (553 × 461 pixels). Furthermore, scaling is not homogeneous between maps (𝑐) and (𝑎) castle's 3D model is inside the castle area on the map which is mainly bounded by Princes St., Kings Stables Rd, Johnston Terrace and The Mount on the map. Note

Table 3 .

 3 2: Results of the "Map tracking experiment". Good matching rates is shown.Detection time and matching time per frame are reported in the form of "Mean value

	± Standard Deviation".
	Map ID	Good	Detection	Matching
		matching (%)	time (ms)	time (ms)
	(𝑎)	100.0%	73.0 ± 11.5 167.0 ± 14.2
	(𝑏)	90.2%	210.6 ± 23.8 223.6 ± 18.7
	(𝑐)	96.0%	81.2 ± 6.1	182.4 ± 15.4
	that neither RDM nor USAC work for any of the videos.

  General algorithm . . . . . . . . . . . . . . . . . . . . . . . . . 60 4.1.1 Definitions and a brief description of the algorithm . . . . . . 61

	Chapter 3. Robust random dot markers (RRDM)
	proposed a workflow for augmenting unprepared paper maps with a GIS database. Chapter 4
	Although RRDM is robust and has almost solved one of the problems mentioned
	at the beginning of this dissertation (cf. Section 1.2), i.e. augmenting different paper Local geometric consensus (LGC)
	maps without any preparation using GIS, it has three major drawbacks: (1) It can
	register the scene point set 𝑄 against only one model point set 𝑃 but it cannot judge
	the correctness of a registration. Therefore, recognition of different point patterns
	cannot be performed. (2) It has quadratic time complexity in relation to the number
	of points to be tracked, which makes the algorithm slow when pattern contains more Contents than 100 points. (3) It can only deal with 2D perspective transformations due to the TSR descriptor. In next chapter, we propose another algorithm which can overcome all these drawbacks. 4.1 . 88
	4.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
	In this chapter, we consider the PPM problem between several model point sets
	and one scene point set. We present a method which can quickly and robustly
	57

4.1.2 Hypotheses generator . . . . . . . . . . . . . . . . . . . . . . 62 4.1.3 Hypotheses validator . . . . . . . . . . . . . . . . . . . . . . . 66 4.1.4 Result refiner . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 4.1.5 Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 4.1.6 Local consensus: a guarantee of low false alert . . . . . . . . 71 4.2 Specific Implementations . . . . . . . . . . . . . . . . . . . . . 73 4.2.1 2D homography . . . . . . . . . . . . . . . . . . . . . . . . . . 73 4.2.2 3D similarity . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 4.3 Results on synthetic point sets . . . . . . . . . . . . . . . . . 76 4.3.1 Speed and robustness studies . . . . . . . . . . . . . . . . . . 76 4.3.2 3D model registration . . . . . . . . . . . . . . . . . . . . . . 80 4.3.3 LGC with additional information . . . . . . . . . . . . . . . . 82 4.4 Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 4.4.1 Tracking ordinary planar objects . . . . . . . . . . . . . . . . 84 4.4.2 Augmenting engineering drawings . . . . . . . . . . . . . . . 85 4.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.5.1 Neighbors . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.5.2 Transformation 𝑇 . . . . . . . . . . . . . . . . . . . . . . . . 86 4.5.3 Repetitive structures . . . . . . . . . . . . . . . . . . . . . . match 2D and 3D point patterns based on their sole spatial distribution, but it can also handle other cues if available. This method can be easily adapted to

  𝑚 𝑖 be the number of points in model 𝑃 𝑖 and 𝑉 𝑖 be the space enclosed by 𝑃 𝑖 's convex hull. Then 𝜌 𝑖 = 𝑑 √ 𝑚 𝑖 /𝑉 𝑖 is the point density of 𝑃 𝑖 , where dimension

		•	p2
	X2	•
		• p1
	p0 •	X1
	p
		p'	q	q'

Figure 4.2: An example of a paired-patch (𝑝, 𝑞) with 𝑘 = 6. 𝑝 and its 6 nearest neighbors are in red. 𝑞 and its 6 nearest neighbors are in blue. Black lines with arrows represent correspondences of the paired-patch. Note that not all points have a correspondence due to missing/extra points, and wrong correspondences could also occur. Let

  Σ𝑖 is used instead of Σ 𝑖 in Alg. 5, with Σ𝑖 = 𝑚𝑎𝑥(Σ 𝑖 , 0.05). Algorithm 5 Offline registration for all 𝑝-patches in model 𝑃 𝑖 for B ∈ {All combination of 𝑏 points in 𝑝-patch containing 𝑝} do Compute geometric descriptor X and its variance Σ Register (model number 𝑖, B ) into all bins that intersect [X -2 Σ, X + 2 Σ]

This region is uniformly partitioned into 100 𝐷 bins (meaning that [-1,1] in each dimension is equally divided in 100 segments) to serve as hash entries. To avoid inefficient record in the hash

table,

  Then point correspondences 𝑠 𝑝,𝑞 are estimated by solving the maximum assignment problem. To reduce the presence of bad hypotheses, the local transformation 𝑇 𝐿 is only estimated when the number of correspondences 𝑠𝑖𝑧𝑒(𝑠 𝑝,𝑞 ) is larger than 𝑏, i.e. the number of points in basis 𝐵. 𝑠 𝑝,𝑞 and 𝑇 𝐿 are considered valid if they satisfy (4.3), where 𝑓 (𝑇 𝐿 ) is the degree of freedom of 𝑇 𝐿 . If so, 𝑠 𝑝,𝑞 is used to create the paired-patch(𝑝, 𝑞), thus giving rise to

	one hypothesis.
	∑	||𝑇 -1 𝐿 q -p|| 2 < 𝜒 2 0.05,𝑠𝑖𝑧𝑒(𝑠𝑝,𝑞)-𝑓 (𝑇 𝐿 )	(4.3)
	(q,p)∈𝑠𝑝,𝑞

Algorithm 6 Online hypotheses generation with 𝑞-patch Input: A randomly selected scene point 𝑞 ∈ 𝑄 Output: A list of hypotheses: hypo Create (𝑝-patch, 𝑞-patch) pre-hypotheses, 𝑝 ∈ ∪𝑃 𝑖 , 𝑖 = 1..𝑀 for all basis 𝐵 * ∈ {𝑞-patch} do Find geometric descriptor X for all (model number 𝑖, 𝐵) in the hash bin which covers X do Let 𝑝 = the origin point of 𝐵 Cast one vote for each correspondence between 𝐵 * and 𝐵 in the local voting table of (𝑝-patch, 𝑞-patch) for all pre-hypotheses (𝑝-patch, 𝑞-patch) do Set votes less than 𝑣 𝑡 to 0 𝑠 𝑝,𝑞 = Solution of maximum assignment problem of the voting results (0 vote excluded) if 𝑠𝑖𝑧𝑒(𝑠 𝑝,𝑞 ) > 𝑏 then Estimate 𝑇 𝐿 using 𝑠 𝑝,𝑞 if 𝑠 𝑝,𝑞 and 𝑇 𝐿 satisfy (4.3) then hypo ← paired-patch (𝑝, 𝑞)

  very few false alerts. Given a model 𝑃 𝑖 and the scene 𝑄, let us consider two random events 𝐸 0 and 𝐸 1 : 𝐸 0 = 𝑄 is a random set. (i.e. 𝑄 does not come from 𝑃 𝑖 )

	4.1.5.2 𝑁 𝑙𝑎𝑟𝑔𝑒 termination rule

The first termination condition on 𝑁 𝑙𝑎𝑟𝑔𝑒 means that we are almost sure to have found the model so there is no need to test more hypotheses. A good choice should

Chapter 4. Local geometric consensus (LGC)

give

  Before solving 𝑃 (𝐸 1 |𝐸 0 ), let us consider a simpler problem (the problem of long runs): In a 𝑁 𝑏 random binary string, 1 appears at each bit with probability 𝑝 𝑏 . The length of the largest sequence containing only 1s is a random variable 𝐿. Then, what is the expectation and the variance of 𝐿? Schilling et al.[START_REF] Schilling | [END_REF] show that the result can be approximated by the following equations:

Table 4 .

 4

	1: Registration results
	Instance ID	1	2	3	4	5
	Key point size	485	485	479	488	491
	Inliers found	459	455	447	457	439
	Rotation difference ( ∘ )	0.10	0.06	0.00	0.10	0.10
	Scale difference (%)	0.09	0.04	0.05	0.05	0.09
	Bunny points, we obtain 480 keypoints. Then, random similarity transformations
	are applied on the bunny model to create 5 instances for matching (cf. Fig. 4.16).

  Algorithm 9 Algorithm 𝐿𝐺𝐶.𝑟𝑒𝑔𝑖𝑠𝑡𝑒𝑟(𝑃 𝑏 ) 𝐿𝐺𝐶.𝑟𝑒𝑔𝑖𝑠𝑡𝑒𝑟(𝑃 𝑝 ) set 𝐻 𝑝𝑟𝑒 such that 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 ) exactly fits the projector's resolution 𝑘 = 0 while 𝑘 < 𝐾 do Projector: Draw circles at 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 ) and project them Camera: Get new image 𝐼 [𝑃 ′ 𝑏 , 𝑃 𝑏 , 𝐻 𝑐𝑏 ] = 𝐿𝐺𝐶.𝑚𝑎𝑡𝑐ℎ(𝐼) 𝑝 , 𝑃 𝑝 , 𝐻 𝑐𝑝 ] = 𝐿𝐺𝐶.𝑚𝑎𝑡𝑐ℎ( Ĩ), Ĩ is 𝐼 with inverted color. if 𝑃 ′ 𝑝 found then 𝑠𝑡𝑒𝑎𝑑𝑦 = ∀𝑝 ∈ 𝑃 ′ 𝑏 has moved less than 𝑡 𝑠 pixels since last image 𝑎𝑙𝑖𝑔𝑛𝑒𝑑 = ∀𝑝 ∈ 𝑃 𝑏 , ‖𝐻 𝑐𝑏 (𝑝) -𝐻 𝑐𝑝 (𝑝)‖ < 𝑡 𝑎 pixels if 𝑠𝑡𝑒𝑎𝑑𝑦 && 𝑎𝑙𝑖𝑔𝑛𝑒𝑑 then if last over 1𝑠 then Capture 𝐼, record correspondences (𝑃 ′ 𝑏 , 𝑃 𝑏 ), (𝑃 ′ 𝑝 , 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 )) and homographies 𝐻 𝑐𝑏 , 𝐻 𝑐𝑝 𝑘 + + else 𝐻 𝑝𝑟𝑒 = 𝐻 𝑝𝑟𝑒 𝐻 -1 𝑐𝑝 𝐻 𝑐𝑏 Improve detection of 𝑃 ′ 𝑏 in rectified views 𝐼 𝑏 = 𝐻 -1 𝑐𝑏 𝐼 Improve detection of 𝑃 ′ 𝑝 in rectified views 𝐼 𝑝 = 𝐻 -1 𝑐𝑝 𝐼

	if 𝑃 ′ 𝑏 found then [𝑃 ′
	Calibrate the system (cf. Section 5.2.1)

  5.7, the projector lens is located at 𝑂. A circle s 𝑜 is located on the projector's plane 𝑝 𝑜 . It is projected and forms a clear image s 𝑖 on the plane of focus 𝑝 𝑖 . All light coming from s 𝑜 and being transmitted by the lens form a blue cone 𝐶 𝑏 , which follows a standard lens projection model. 𝐶 𝑏 intersects the board plane 𝑝 𝑏 and forms the light spot s 𝑏 on the board.During the calibration, s 𝑏 is captured by the camera and forms an image 𝐻 𝑐𝑏 (s 𝑏 )

  4)The projector's plane 𝑝 𝑜 is at 𝑧 = 𝑧 𝑜 and the plane of focus 𝑝 𝑖 is at 𝑧 = 𝑧 𝑖 . 𝑑 is the diameter of the projector lens. The circle s 𝑜 has it center at 𝐸 𝑜 (𝑥 𝑜 , 𝑦 𝑜 , 𝑧 𝑜 ) with 𝑟 𝑜 , the circle s 𝑖 has a radius of 𝑟 𝑖 . Then the equation of 𝐶 𝑏 is:

  px, the effective focal length (i.e. object distance) 𝑓 = -𝑧 𝑜 ≈ 2000px, 𝑏 = 50cm, 𝑟 𝑜 = 6px. The maximum error 𝑒 𝑚𝑎𝑥 is obtained when 𝑣 → ∞, 𝜙 = 𝜙 𝑚𝑎𝑥 = 135 ∘ and 𝑥 𝑜 = √ 𝑊 2 /4 + 𝐻 2 . It can be approximated by:

  𝐾 𝑐 , 𝐷 𝑐 , 𝐾 𝑝 , 𝐷 𝑝 , 𝑅, 𝑡 to be evaluated, 𝑄 𝑏 and 𝑄 𝑤 Output: List of reprojection errors { ⃗ 𝑒 𝑟 } 𝐼 ← 𝑢𝑛𝑑𝑖𝑠𝑡𝑜𝑟𝑡𝐶𝑎𝑚𝑒𝑟𝑎𝑉 𝑖𝑒𝑤(𝐾 𝑐 , 𝐷 𝑐 ) 𝑄 ′ 𝑏 ← 𝑑𝑒𝑡𝑒𝑐𝑡𝐵𝑙𝑎𝑐𝑘𝐷𝑜𝑡𝑠(𝐼)* 𝐻 𝑐𝑏 ← 𝑐𝑜𝑚𝑝𝑢𝑡𝑒𝐻𝑜𝑚𝑜𝑔𝑟𝑎𝑝ℎ𝑦(𝑄 ′ 𝑏 , 𝑄 𝑏 ) with 𝑄 ′ 𝑏 = 𝐻 𝑐𝑏 (𝑄 𝑏 ) 𝑅 𝑐𝑏 , 𝑡 𝑐𝑏 ← 𝑐𝑜𝑚𝑝𝑢𝑡𝑒𝑃 𝑎𝑡𝑡𝑒𝑟𝑛𝑃 𝑜𝑠𝑖𝑡𝑖𝑜𝑛(𝐾 𝑐 , 𝐻 𝑐𝑏 ) 𝑅 𝑝𝑏 , 𝑡 𝑝𝑏 ← 𝑐𝑜𝑜𝑟𝑑𝑖𝑛𝑎𝑡𝑒𝑇 𝑟𝑎𝑛𝑠𝑓 𝑜𝑟𝑚𝑎𝑡𝑖𝑜𝑛(𝑅, 𝑡, 𝑅 𝑐𝑏 , 𝑡 𝑐𝑏 ) 𝐻 𝑝𝑏 ← 𝑐𝑜𝑛𝑠𝑡𝑟𝑢𝑐𝑡𝐻𝑜𝑚𝑜𝑔𝑟𝑎𝑝ℎ𝑦(𝑅 𝑝𝑏 , 𝑡 𝑝𝑏 , 𝐾 𝑝 ) 𝐼 𝑑𝑝 ← 𝑝𝑟𝑜𝑗𝑒𝑐𝑡𝐷𝑖𝑠𝑡𝑜𝑟𝑡𝑒𝑑𝐼𝑚𝑎𝑔𝑒(𝐼(𝑄 𝑝 ), 𝐷 𝑝 , 𝐻 𝑝𝑏 ) 𝑄 ′ 𝑤 ← 𝑑𝑒𝑡𝑒𝑐𝑡𝑊 ℎ𝑖𝑡𝑒𝐷𝑜𝑡𝑠(𝐼 𝑑𝑝 )* { ⃗ 𝑒 𝑟 } ← 𝐻 -1 𝑐𝑏 (𝑞 ′ 𝑤 ) -𝑞 𝑤 , for all (𝑞 ′ 𝑤 , 𝑞 𝑤 ) ∈ (𝑄 ′ 𝑤 , 𝑄 𝑤 )

	Chapter 5. Defocused projector calibration for projector-camera systems
	Algorithm 10 Evaluation
	Input: (*)Perspective and lens' distortions are removed from images. This allows us to
	find correct circle centers.

  Related work . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 6.2 Surface of revolution axis estimation . . . . . . . . . . . . . . 116

	Chapter 6. Surface of revolution reconstruction from 3D data
	Chapter 6 detecting pre-known models. Bottom-up approaches are sometimes called geomet-
	ric reconstruction methods. They detect primitive shapes in a scene and establish
	Surface of revolution relationships between these shapes whenever possible.
	Primitive detection and localization is a crucial step in bottom-up approaches:
	it represents the first step of those methods and determines their application reconstruction from 3D data range. However, existing approaches can only deal with a limited type of prim-
	itive shapes such as: planes only [Nguyen 2015], cylinders only [Qiu 2014] or
	thin-structures [Song 2015]. More complex primitives such as spheres, cones or
	tori [Schnabel 2007] can be detected but those techniques still remain limited to
	pre-defined parametric shapes. It is important for applications in reverse engineer-Contents ing or in the construction field to deal with various types of geometric primitive
	shapes. Also, on-line scene understanding is desirable as dense SLAM (i.e. SLAM 6.1 This chapter presents the work realized during a three-month internship in from RGB-D data) runs in real-time from 3D data [Tateno 2016].
	Kyushu University in Japan, and has been presented at International Conference
	on 3D Vision (3DV) in 2016. It tackles scene understanding from 3D data, which
	is completely different from the previous chapters. Scene understanding from 3D
	data is an important research topic in computer vision and robotic perception be-
	cause it can be applied to various types of systems such as robotic grasping, model

6.2.1 Basic idea . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 6.2.2 Approximately linear objective function . . . . . . . . . . . . 118 6.3 Implementation details . . . . . . . . . . . . . . . . . . . . . . 121 6.3.1 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121 6.3.2 Plane cutting and circle fitting . . . . . . . . . . . . . . . . . 122 6.3.3 Determining 𝑃 𝑠 and solving 𝜑(𝜃) . . . . . . . . . . . . . . . . 124 6.4 Real-time SoR reconstruction . . . . . . . . . . . . . . . . . . 125 6.4.1 Segment classification . . . . . . . . . . . . . . . . . . . . . . 125 6.4.2 Workflow . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 6.5 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 6.5.1 Synthetic study . . . . . . . . . . . . . . . . . . . . . . . . . . 128 6.5.2 Real data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 6.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 simplification, augmented reality, etc. Existing methods on this topic can be classified into two main families: top-down and bottom-up approaches [Nguyen 2015].

Top-down approaches are often referred to as object detection, since they aim at

  Qiu et al. [Qiu 2014] use parallel cylinder de-Once the position of its rotation axis is determined, the profile can be easily calculated. Existing approaches are dedicated to the estimation of the rotation axis of a SoR, and can be classified as direct, iterative and brute force methods. For direct methods, Pottmann et al.[Pottmann 1999] use 3D line geometry to classify surfaces and find the rotation axis of a SoR by solving a 6D eigenvalue problem. For iterative methods,[Lou 2009],[Willis 2003] and[Pavlakos 2015] use ICP-like method[Besl 1992a]. They first estimate the profile function of a SoR by using a presumed rotation axis, thus creating an initial model of the SoR. Then, they try to align this model with the point cloud, which results in a better estimation of the axis for the next iteration. Brute force methods use the fact that the intersection between the SoR and a plane perpendicular to its rotation axis is a circle. Han et al.[START_REF] Han | [END_REF]] use planes containing the normal of a surface point to cut a SoR and evaluate the intersecting curve. The plane which gives the least curvature variation is used to calculate the rotation axis. Mara et al.[START_REF][END_REF]] use a set of parallel planes to cut a SoR. Then circles are fitted on each plane and the variation of circle centers across all planes is calculated. This set of planes is then rotated in a 2-dimensional search space so that the set which gives the least variation of circle centers is used to calculate the rotation axis.Although direct methods are quite fast, they cannot find accurate results if input data contains even moderate noise (cf. Section 6.5). Iterative methods use the whole point cloud, fit the SoR profile at each step, and are thus very time consuming. Brute force methods are time consuming as well.

	6.1. Related work
	All these methods can detect a limited set of parametric shapes, and are not designed
	to work in real-time.
	A surface of revolution (SoR), or rotational object, is formed by rotating a 2D
	curve (i.e. the profile) in 3D space with respect to a 3D line (i.e. the rotation
	axis), thus it is axis-symmetric and is a more general representation of most geo-
	metric primitives.
	tection to reconstruct an industrial pipeline system. Song et al. [Song 2015] use
	beams and planes to reconstruct thin structural systems, such as chairs. Thanh
	et al. [Nguyen 2015] use planes to find semantic structures in a scene. In order to
	model the structure of more complex objects, Schnabel et al. [Schnabel 2007] pro-
	pose a RANSAC-like method which can detect cones, cylinders, spheres and tori to

reconstruct more types of primitive shapes. Taguchi et al.

[START_REF] Taguchi | [END_REF]

] improved this method by introducing distance field functions for efficient computation. Drost et al.

[START_REF] Drost | [END_REF]

propose to use Hough Transform to detect cylinders and spheres.

  Parallel planes intersect with the 𝑧-axis at 𝑧 𝑗 . The 2D projection of 𝐶 𝑗 s on 𝑃 𝑠 can be represented by segments between (𝑧 + 𝑗 , 𝑟 + 𝑗 ) and (𝑧 - 𝑗 , 𝑟 - 𝑗 ). 𝜙 is defined as the offset angle of initial guess ⃗ 𝑙 0 from ⃗ 𝑙 𝑔 :

Table 6 .

 6 𝐿 is a 3 × 3 matrix and 𝑀 is a 6 × 6 matrix. Let {𝜆 𝐿,𝑖 : 𝑖 = 1, 2, 3} be 1: Different surface types according to 𝑓 𝑡 and 𝑓 𝑟 6.4. Real-time SoR reconstruction threshold 𝑐 𝑡 , we check the stability by estimating the variation of 𝑆 𝑣 's rotation axis.

		𝐿 =	1 𝑁	𝑁 𝑗=1 ∑	⃗ 𝑛 𝑗 ⃗ 𝑛 𝑇 𝑗
		𝑀 =	𝑁 1	𝑗=1 ∑ 𝑁	𝑗 ⃗ 𝑟 𝑗 ⃗ 𝑟 𝑇	(6.21)
	𝑓 𝑡	2				1
	𝑓 𝑟	-		> 2	1	0
	Type	Plane	Sphere	SoR	Complex

where

  The rotation axis of 𝑆 𝑚 is estimated and the SoR is then reconstructed for the first time by applying our methods (cf. Sections 6.2-6.3) on the global model 𝑆 𝑚 , with 𝜏 (cf. Section 6.3.2) being the average of diameters (i.e. 2×radius) of all surfels in 𝑆 𝑚 .More details are shown in Alg. 12.Algorithm 12 Stability Estimation and first reconstructionInput: 𝑆 𝑣 , 𝑆 𝑚 at current frame Output: (Possibly) first reconstruction Calculate ⃗ 𝑙 from 𝑆 𝑣 using P-method[Pottmann 1999] if𝑐 = 0 then Set ⃗ 𝑙 𝑟 = ⃗ 𝑙 else if cos -1 ( ⃗ 𝑙 • ⃗ 𝑙 𝑟 ) < 𝛾 𝑡 then 𝑐 = 𝑐 + 1 if 𝑐 ≥ 𝑐 𝑡 thenEstimate initial ⃗ 𝑙 𝑚 0 from 𝑆 𝑚 with P-method Improve axis estimation & reconstruct from 𝑆 𝑚

  Chapter 6. Surface of revolution reconstruction from 3D data Algorithm 13 Check reconstruction precision and improve reconstruction Input: 𝑆 𝑣 , 𝑆 𝑚 at current frame, previous reconstruction Output: (Possibly) new reconstruction for each 𝑠 𝑖 ∈ 𝑆 𝑣 do 𝑒 𝑖 ← 𝑐𝑎𝑙𝑐𝑀 𝑜𝑑𝑒𝑙𝑖𝑛𝑔𝐸𝑟𝑟𝑜𝑟(𝑠 𝑖 , 𝑆 𝑚 ) 𝑒𝑛𝑢𝑚 = 𝑒𝑛𝑢𝑚 + 1 if 𝑒 𝑖 > 𝑒 𝑡 if 𝑒𝑛𝑢𝑚 > 𝑁 𝑡 then

A novel point pattern matching algorithm (RRDM) to robustly and quickly

(4.16)

http://campar.in.tum.de/Chair/ProjectInSeg

sets can be derived from 𝑠𝑙, it is very easy to calculate the transformation between them. The reason for using mesh neighbors instead of nearest neighbors is explained in Section 4. 5

until 𝑠𝑙 stops growing Return: 𝑇 In Alg. 8, 𝑑 𝑡 is a threshold which describes how far could a projected scene point 𝑇 -1 (𝑞 ′ ) be away from its model corresponding point. This value can vary according to the position of 𝑞 ′ . Let C be the convex hull of scene points in 𝑠𝑙. If 𝑞 ′ is inside C , we can take roughly 𝑑 𝑡 = 2𝜎 with 𝜎 being the variance of the additive Gaussian noise. When 𝑞 ′ is outside C , the farther it is from C , the larger 𝑑 𝑡 will be (cf. Fig. 4.5). Assuming 𝑜 is the center of C and 𝐶 intersects 𝑜𝑞 ′ at 𝑞 𝐶 , we then take 𝑑 𝑡 = 2 𝑜𝑞 ′ 𝑜𝑞 𝐶 𝜎.

In the worst case, each repeat loop only adds one correspondence into 𝑠𝑙, which means that repeat loop runs at most 𝑛 times (i.e. the number of points in 𝑄). Thus the time complexity in the worst case is about 𝑂(𝑛). We use the term "about" since we consider the inner loop is fast enough and thus not influence much the refiner module.

Parameters

All parameters used in LGC will be discussed here so that they can be more easily adjusted according to specific use.

Results on synthetic point sets

In this section, we present results obtained under different conditions. Section

Speed and robustness studies

The performance of our method is compared with other methods using synthetic data.

In the graphs, RDM stands for Uchiyama's Random Dot Markers [Uchiyama 2011b] (the original "tracking by detection" implementation is used), RRDM for the Chapter 3, LGC for our current proposal and GH for traditional Geometric Hashing [START_REF] Wolfson | [END_REF]]. The default parameters of these methods are used.

The 2D homography is set to a 30 ∘ perspective transform unless stated otherwise. We first investigate the performance with different number of points in the model point set 𝑃 . Fig. 4.9 shows the result in an ideal condition, that is without jitter Corners in real textured images are assumed to be randomly distributed. As demonstrated in [Uchiyama 2011a], these corners can be used for tracking. We compare our method with RDM and SURF in this section on real planar targets. However, our purpose is not to compete with textural-based trackers in these situations. SURF serves as a reference. We aim at illustrating the flexibility of our approach to deal with traditional targets and showing its robustness on real targets.

We use the OpenCV "goodFeaturesToTrack" method [Shi 1994] as interest point detector in this experiment which performs at about 10ms/frame.

The max-Corners and minDistance parameters are selected in favor of RDM according to [Uchiyama 2011a] (qualityLevel =0.13 and blockSize=12 for stable point detection). All remaining parameters are set to their default values. As to SURF, a FLANN based matcher in OpenCV is used to build a database of descriptors for tracking.

We use a Logitech C270 camera with resolution set to 640×480.

Discussion

In this section, we discuss some aspects of the algorithm.

Neighbors

We use nearest neighbors to create local patches but use mesh neighbors in the refiner. We tried to use the same type of neighbors in the algorithm but neither performs better. Although mesh neighbors are more robust against perspective

Transformation 𝑇

Theoretically, the algorithm works with any transformation 𝑇 , but some conditions on point sets have to be satisfied.

Chapter 5. Defocused projector calibration for projector-camera systems is to calibrate the projector with a far focus by using a small calibration board, the board does not have to be close to the plane of focus of the projector, so the projected points will be defocused and become larger. The minimum inter-point distance requirement makes the defocused projected points less likely to overlap with printed points.

Algorithm

We use LGC to find point correspondences, the two basic operations of this algorithm are defined as follows:

𝐿𝐺𝐶.𝑟𝑒𝑔𝑖𝑠𝑡𝑒𝑟(𝑃 )

(5.2)

[ P , 𝑃 𝑖 , 𝐻] = 𝐿𝐺𝐶.𝑚𝑎𝑡𝑐ℎ(𝐼)

2) registers a model point pattern 𝑃 into LGC, several models may be registered; (5.3) matches detected black points in image 𝐼 with all registered models. If a model 𝑃 𝑖 is found in the detected pattern, (5.3) returns 𝑃 𝑖 and its corresponding detected point set P . An homography 𝐻 : P = 𝐻(𝑃 𝑖 ) is also returned.

The algorithm contains three parts: initialization, manipulation and improvements (cf. Alg. 9). For initialization, both 𝑃 𝑏 and 𝑃 𝑝 points are registered into

LGC as models. An initial value 𝐻 𝑝𝑟𝑒 is defined so that 𝐻 𝑝𝑟𝑒 (𝑃 𝑝 ) exactly fits the projector resolution.

During manipulation, the algorithm matches 𝑃 𝑏 (i.e. the calibration board) with LGC. Once the board is found, it starts to detect and match 𝑃 𝑝 (i.e. the projected pattern) on the calibration board in the same way. If 𝑃 𝑝 is found as well, an homography 𝐻 𝑝𝑟𝑒 is computed whose role is to warp the projected pattern to the right position for the next iteration (cf. Fig. 5.6). Note that 𝐻 𝑝𝑟𝑒 only warps points' positions, while the points to be projected are always circular. When both patterns are well aligned (meaning that they form well 𝑃 𝑜 ) and the board is steady for ≈ 1s, an image 𝐼 is captured, and both 𝐻 𝑏𝑐 and 𝐻 𝑝𝑐 are recorded. We choose 𝑡 𝑠 = 3.0 pixels and 𝑡 𝑎 = 2.0 pixels in Alg. 9 empirically.

Once images of 𝐾 views are obtained, an improvement step is executed. It is used to deal with the fact that the center of a circle is not perspective invariant:

after applying a perspective transformation, the transformed center of the circle is not the center of the transformed circle. The approach is similar to [START_REF] Datta | [END_REF]]:

each captured image 𝐼 is rectified so that points are nearly circle. More precisely, since the transformation from the calibration board to the camera image plane is Chapter 5. Defocused projector calibration for projector-camera systems calibration, 𝐾 = 10 images are captured, following [Audet 2009]. Both methods call

OpenCV's calibrateCamera method for final calibration, with the same set of options and parameters. To prevent the projector's optical properties from changing with the temperature, the projector is pre-warmed 30 minutes before the experiments.

We first show the results of calibration reprojection root mean square errors (RMSE) in Fig. 5.9. RMSE is the most commonly used error measurement for camera and projector calibrations. Compared to [Audet 2009], our average RMSE is smaller but our maximum RMSE is larger. This is reasonable since the more correspondences a method uses, the more likely extreme values appear.

Both methods rely on Zhang's method [Zhang 2000] to find the projector's intrinsic matrix and distortions coefficients by minimizing the sum of reprojection errors of all correspondences. However, as pointed out by [González 2005, Shih 1996], minimizing reprojection errors on calibration data cannot guarantee that the estimated internal parameters are the best ones. Ideally, when no modification is applied to the ProCam system, separated consecutive calibrations should lead to repeatability, i.e. stable intrinsic matrix and lens distortions for both devices. Focal lengths estimated by our method are very stable, with 𝑡(5) = 2.69 for standard deviation (SD) of 𝑓 𝑝𝑥 and 𝑡(5) = 2.49 for SD of 𝑓 𝑝𝑦 , 𝑝 < 0.05.

Algorithm

The outline of rotation axis estimation is shown in Alg. 11. It takes a 3D representation (e.g. a point cloud, a SDF) of a SoR and an initial axis estimation 𝑙 0 obtained from [Pottmann 1999] as inputs and returns the rotation axis 𝑙 𝑒 .

Algorithm 11 Rotation axis estimation

Input: 3D data of SoR, initial rotation axis 𝑙 0 estimated using [Pottmann 1999] Output: SoR rotation axis 𝑙 𝑒 Parameter: number of cutting planes 𝑛 1 , 𝑖 𝑚𝑎𝑥 𝑖 = 0, 𝜃 0 = 0, 𝑣𝑎𝑙𝑢𝑒𝑃 𝑎𝑖𝑟 = ∅ while 𝑖 < 𝑖 𝑚𝑎𝑥 do if 𝑖 ̸ = 0 then

Estimate 𝐴 0 with (O, r) from (6.20) else

Surface of revolution reconstruction from 3D data

The process of cutting the SoR by parallel planes and finding the centers of intersecting curves is frequently used in our method. Since the angle between the presumed axis 𝑙 𝑖 and the ground-truth 𝑙 𝑔 is usually smaller than 10 ∘ (cf. Fig. 6.7), we use circle fitting to find the centers 𝑂 𝑗 . We call this process findCircleCenters.

O, w, R, ē ← 𝑓 𝑖𝑛𝑑𝐶𝑖𝑟𝑐𝑙𝑒𝐶𝑒𝑛𝑡𝑒𝑟𝑠(𝑙, 𝑛) (6.15) This process takes the presumed axis 𝑙 and the number of planes 𝑛 as inputs and returns 𝑛 𝑣 valid circles as output, with circle centers O = {𝑂 1 , 𝑂 2 , ..., 𝑂 𝑛𝑣 }, circle radii R = {𝑅 1 , 𝑅 2 , ..., 𝑅 𝑛𝑣 }, weights w = {𝑤 1 , 𝑤 2 , ..., 𝑤 𝑛𝑣 } and average geometric errors ē = {ē 1 , ē2 , ..., ē𝑛𝑣 }. 𝑛 𝑣 is the number of valid circles, which may be smaller than 𝑛. A detailed explanation is provided in Section 6.3.2. Different 𝑛 can be used as input in different cases. In Alg. 11, we use 𝑛 1 in order to distinguish 𝑛 2 that is used for reconstruction (see later in this Section).

In Alg. 11, findSymmetricPlane is the process of finding 𝑃 𝑠 . It takes all valid circle centers O as input and gives 𝑃 𝑠 and 𝛿 as outputs. 𝛿 represents the incertitude of circle centers detection, which is explained in Section 6.3.3. findCenter-LineDirection is the process of finding the direction of 𝑣 (cf. Fig. 6.2). It takes circle centers O and their weights w as inputs. It fits a line to points O by using different weights w and returns ⃗ 𝑣 as output. calculateSlope is the process of finding the slope 𝐴 of 𝜑(𝜃). It takes as input the previous 𝑖 + 1 value pairs 𝑣𝑎𝑙𝑢𝑒𝑃 𝑎𝑖𝑟 = {(𝜃 0 , 𝜑 0 ), (𝜃 1 , 𝜑 1 ), .., (𝜃 𝑖 , 𝜑 𝑖 )}, fits a line to them and returns the slope 𝐴 of this fitted line.

After finding ⃗ 𝑙 𝑒 , the whole object can be reconstructed. 𝑓 𝑖𝑛𝑑𝐶𝑖𝑟𝑐𝑙𝑒𝐶𝑒𝑛𝑡𝑒𝑟𝑠 is used as well, but with inputs ⃗ 𝑙 𝑒 and 𝑛 2 . 𝑛 2 can be larger than 𝑛 1 in Alg. 11 for a finer reconstruction result. A surface's intrinsic coordinate system is constructed, with origin 𝑂 𝑒 being the average of 𝑂 𝑗 and the direction of z-axis ẑ𝑒 = ⃗ 𝑙 𝑒 . Within this coordinate system, the profile function 𝑓 (𝑧 𝑒 ) is fitted by splines, thus finishing the SoR reconstruction.

The global reconstruction error 𝑒, is defined as the average of geometrical error ē (cf. Section 6.3.2). 𝑒 = 𝑎𝑣𝑒𝑟𝑎𝑔𝑒(ē) the eigenvalues of 𝐿 and 𝑓 𝑡 be the number of small eigenvalues in {𝜆 𝐿,𝑖 }. According to [Lou 2009], the segment can be translated in 𝑓 𝑡 directions without changing the distance function values nearby. We can say that the segment has 𝑓 𝑡 translation free directions. Let {𝜆 𝑀,𝑖 : 𝑖 = 1, 2, 3} be the eigenvalues of problem 𝑀⃗ 𝑥 = 𝜆𝐷⃗ 𝑥,

where 𝐷 = 𝑑𝑖𝑎𝑔(1, 1, 1, 0, 0, 0), and 𝑓 𝑟 be the number of small eigenvalues in {𝜆 𝑀,𝑖 }.

According to [Pottmann 1999], the segment has 𝑓 𝑟 rotation free directions. With both 𝑓 𝑡 and 𝑓 𝑟 , surfaces can be classified into 4 groups as shown in Table 6.1. Practically, an eigenvalue is considered to be small when one of the following conditions are true:

𝜆 𝐿,𝑖 < 0.05 𝜆 𝑀,𝑖 < 0.5𝜏 2 (6.22)

Workflow

Let us take one segment of SoR as an example to illustrate the workflow. All other segments follow the same procedure.

Two sets of surfels of the SoR, 𝑆 𝑚 and 𝑆 𝑣 , can be obtained from [Tateno 2015] at each frame. 𝑆 𝑚 contains all surfels in the model obtained from the SLAM procedure.

𝑆 𝑣 is the visible subset of 𝑆 𝑚 at current frame. When the SoR has only appeared in front of the camera for a few frames, or when it is occluded or far from the camera, the acquired data may be not stable. With more information coming from new frames, data become more stable and the reconstruction more reliable. We use a confidence label 𝑐 to represent the stability of data points acquired from the SoR.

The workflow to reconstruct the SoR is shown in Fig. LGC est une amélioration de la méthode RRDM et peut mettre en correspondance des ensembles de motifs de points 2D ou 3D subissant une transformation dont le type est connu.

LGC présente un comportement linéaire en fonction du nombre de points, ce qui permet un tracking en temps-réel. La pertinence de LGC a été illustrée en développant une application de calibration de système projecteur-caméra dont les résultats sont comparables avec l'état de l'art tout en présentant des avantages pour l'utilisateur en termes de taille de mire de calibration..

Mots-clés

réalité augmentée, recalage des motifs de points, suivi, peu-texturé, calibration, projecteur-caméra, carte augmentée, surface de révolution

Abstract

Registration is a very important task in Augmented Reality (AR). It provides the spatial alignment between the real environment and virtual objects. Unlike tracking (which relies on previous frame information), wide baseline localization finds the correct solution from a wide search space, so as to overcome the initialization or tracking failure problems. Nowadays, various wide baseline localization methods have been applied successfully. But for objects with no or little texture, there is still no promising method. One possible solution is to rely on the geometric information, which sometimes does not vary as much as texture or color.