Vérifier l'existence d'un point dans un nuage de points

Version imprimable

Voir 40 message(s) de cette discussion en une page

11/05/2015, 12h13
MilWolf

Vérifier l'existence d'un point dans un nuage de points

Bonjour à tous,

je suis actuellement à la recherche d'un algorithme performant pour vérifier l'existence ou non d'un point parmi un nuage de points.

Les contraintes sont:
-le nuage de points peut contenir des millions de points (~3)
-la recherche doit s'exécuter des millions de fois (~20)
-le programme est en Java (pas forcément une contrainte)
-il peut y avoir des erreurs de précision, le nuage de points étant lu depuis un fichier texte et le point cherché étant recalculé

Parmi les algorithmes possibles j'ai vu que l'on peut:
-partitionner l'espace avec un Octree
-utiliser le locality sensitive hashing (cf: wikipedia)

Mais je me demandais si il n'y avait pas un moyen plus simple les algorithmes ci-dessus sont tous fait pour rechercher les voisins les plus proches dans mon cas il s'agit de déterminer si le point appartient ou non au nuage de point sachant que le nuage de points est un extrait des points cherchés par la suite.

Je vous remercie pour votre aide.
Bonne journée.

MilWolf.
11/05/2015, 13h53
dourouc05

:salut:

À lire ton problème, je pense au filtre de Bloom, une structure de données probabiliste pour vérifier rapidement l'appartenance à un ensemble : si la réponse est négative, le point n'appartient pas à l'ensemble ; si elle est positive, il y appartient avec une très grande probabilité (liée à la taille en mémoire de la structure : https://en.wikipedia.org/wiki/Bloom_...alse_positives) ; la complexité est constante pour une recherche et une insertion. Maintenant, je ne suis pas sûr d'avoir compris cette partie de ton message (un peu de ponctuation pourrait aider) :

Citation:

Envoyé par MilWolf

Mais je me demandais si il n'y avait pas un moyen plus simple les algorithmes ci-dessus sont tous fait pour rechercher les voisins les plus proches dans mon cas il s'agit de déterminer si le point appartient ou non au nuage de point sachant que le nuage de points est un extrait des points cherchés par la suite.
11/05/2015, 14h28
MilWolf

Merci de ta réponse.
Je reprends donc la partie qui n'est pas claire:

Le nuage de points a été extrait d'un nuage de points plus grand, dans mon cas il s'agit d'un arbre dans une parcelle de forêt.
Par la suite je cherche pour chaque point de la parcelle si il appartient au sous-ensemble représentant mon arbre.
Comme il s'agit d'un point appartenant ou non à l'ensemble, mais pouvant contenir des erreurs dues à de la simple précision, y-aurait il un algorithme permettant par exemple dans un temps O(1) de connaître la réponse?

Ensuite j'ai regardé le filtre de Bloom, et comme le locality sensitive hashing il fait appel à des fonctions de hashage, et j'avoue ne pas être familier avec ces notions sans compter qu'on obtient qu'une probabilité au final.

J'avais pensé à un algorithme simple mais sûrement pas très rapide:

-Tri de mon nuage de points selon x, y, z.
-Récupération des points minimum et maximum
-Ajout des points dans un octree
-Recherche d'un point dans l'octree

Mais le problème avec cette méthode est qu'il faut quand même à un certain point fait une recherche incrémentielle dans un sous-ensemble.

J'ai pensé à une autre technique:

Création d'une clé calculée à partir de mes points x, y, z, pour ceux qui connaissent la "pairing function (cf: http://en.wikipedia.org/wiki/Pairing_function)"
Il s'agirait là de trouver comme générer une clé unique à partir de mes coordonnées puis de retrouver les coordonnées à partir de la clé.

En attente de vos réponses.
Bon après-midi.

MilWolf.
11/05/2015, 19h44
souviron34

et pourquoi pas une triple recherche dichotomique ??

Tu classes tes points suivant par exemple z croissant, y croissant, x croissant..

Et tu fais une recherche dichotomique..
11/05/2015, 19h58
MilWolf

La recherche dichotomique est une bonne idée.
Toutefois je me demande quelles seraient les performances en comparaison d'une recherche dans une structure Octree.

Pour l'octree il faut prendre en compte sa création.
Pour la recherche dichotomique le temps de tri des éléments.

J'attends de voir si des personnes ont une idée du coût que ça représente.
Pendant ce temps je peux implémenter la recherche dichotomique ce qui sera le plus rapide.

Merci pour vos réponses.

MilWolf.
11/05/2015, 20h10
souviron34

ben le tri : entre O(N) and O(N logN)

dichotomie : O(log(N))

C'est très rapide et efficace, la dichotomie... (et sans structures complexes)
11/05/2015, 20h20
MilWolf

Et bien ce serait une bonne nouvelle dans ce cas.

Bon j'ai implémenté le tri sur un de mes fichiers contenant ~320000 points et ça met à peu près 1/10 de seconde.
Reste à voir le comportement de la recherche dichotomique et surtout sur des datasets beaucoup plus grands.

J'attends quand même d'autres avis.
Merci d'avance.

Bonne soirée.
MilWolf.
13/05/2015, 12h02
dinobogan

Une info n'est pas très claire pour moi : le point que tu recherches dans ton ensemble, a-t-il une marge d'erreur ?
Si tu as un point X(a ; b ; c) à chercher dans l'ensemble, t'attends-tu à trouver exactement X ? Si non, quelle est la marge d'erreur sur les 3 axes ?

Autre question : quelle est l'amplitude sur les 3 axes ?
La pairing function n'est pas très économe en espace mémoire. Pour un cas similaire au problème résolu par la pairing function, j'ai choisi une autre approche : je prend l'amplitude maximale des clés. Chaque clé commence à zéro. Ayant l'amplitude, je connais le nombre de bit maximum pour représenter chaque clé. Enfin, je colle les clés les unes aux autres en remplissant pour chaque clé avec des zéros afin d'arriver au bon nombre de bit pour chaque clé. La nouvelle valeur contenant toutes les clés est donc compactée au maximum et faire le calcul dans les deux sens est extrêmement rapide puisque ce ne sont que des décalages de bits.

Bonjour dinobogan,

Le point recherché dans mon ensemble peut soit être à la valeur exacte, soit contenir une marge d'erreur en x et/ou y et/ou z dû à un arrondi suite à l'export.
Ces coordonnées étant exportées 3 chiffres après la virgule et arrondi et/ou tronqué (je me renseigne encore là dessus) tandis que les coordonnées d'origine contiennent 8 chiffres après la virgule.

J'ai implémenté la recherche dichotomique (binary search) et j'obtiens des performances correctes.
Je recherche pour ~65 000 0000 de points si ils existent dans dans nuage de ~1 150 000 points.
Sur mon pc la recherche s'effectue en 200s.

Voici l'implémentation:
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 public int nearestPoint(Point3f point, float maxDistance){ int low = 0; int high = points.size()-1; while (low <= high) { int mid = (low + high) >>> 1; Comparable<Point3f> midVal = points.get(mid); int cmp = midVal.compareTo(point); if (cmp < 0){ low = mid + 1; }else if (cmp > 0){ high = mid - 1; }else{ return mid; // key found } } if((low < points.size()) && (point.distanceTo(points.get(low)) < maxDistance)){ return low; }else if(high >= 0 && (point.distanceTo(points.get(high)) < maxDistance)){ return high; } return -(low + 1); // key not found }
avec "points" étant le nuage de points dans lequel je recherche et "poin"t étant le point parmi les 65 000 000 à tester.

Pensez-vous que cette implémentation est correcte?

Bonne journée!

13/05/2015, 14h15
dinobogan

Il y aurait probablement une préparation des données à faire au début, permettant de gagner encore du temps mais pour cela il faudrait un peu plus d'infos sur les données.
Quelles sont les amplitudes (valeurs min et max) sur chaque axe ?

Voici le but de ma question : les calculs flottants sont beaucoup plus lent que des calculs entiers. si l'amplitude est assez petite, tu vas pouvoir convertir tes flottants en entier en multipliant par 1000 (car 3 chiffres significatifs après la virgule). Tu conserves le même algorithme mais il faut travailler avec des entiers.

Second point, les objets. Tu travailles avec des Point3f. Il y a beaucoup de constructions/destructions d'objets. Utilises directement des entiers en paramètres, idem pour le stockage en RAM : 3 tableaux d'entiers pour les 3 axes. Le temps de calcul s'en trouvera fortement diminué.

Quand à ton implémentation par recherche dichotomique puis du plus proche voisin, elle me semble correcte.
13/05/2015, 14h35
MilWolf

Pour l'amplitude, actuellement le nuage de points a une amplitude de:
x: 18
y: 30
z: 40

Les points cherchés eux ont une amplitude de:
x:100
y:100
z:50

Les coordonnées sont en mètres.

Par contre le nuage de points ne sera pas toujours le même donc l'amplitude pourra changer.
Quand aux chiffres après la virgule ce n'est pas encore fixé.

Les coordonnées peuvent très bien se trouver dans un espace local (aux alentours d'une origine à 0) que dans un espace global (origine: (-50000, 30000, 40) par exemple).
Ce qui fait que je peux me retrouver avec une coordonnée en x de -50000,035 (soit -50000035 en entier) donc en fonction du nombre de chiffres significatifs je peux risquer de dépasser la limite d'un entier.

Mais ça dépendra vraiment du contexte, je dois encore me renseigner pour savoir comment les coordonnées seront exportées.

Merci pour ton aide ;)
13/05/2015, 15h03
dinobogan

L'amplitude ne tient pas compte du zéro de départ. Tu peux tout à fait faire un décalage par une valeur fixe pour toutes tes valeurs. Pour un ensemble de points, tu prends la valeur minimum sur chaque axe et pour tous les points de ton ensemble, tu soustrais son min à chaque axe et tu te retrouves avec la même amplitude mais un point d'origine à zéro.
Exemple en une dimension : j'ai un ensemble de valeurs qui va de -10 milliard à 10 milliards. Son amplitude est donc 20 mais n'est pas représentable avec des entiers 32 bits. Maintenant si je retire -10 milliards (donc j'ajoute +10 milliards) à toutes les valeurs, l'amplitude est toujours 20 mais les valeurs iront de 0 à 20. La représentation par des entiers 32 bits devient possible.

L'amplitude d'un entier 32bits est 4 milliards. On peut donc mesurer des distances de 4 milliards de millimètres (donc 3 chiffres après la virgule), soit environ 4 fois la longueur nord-sud de la France.
Pour un entier 64bits, on passe à 18 milliards de milliards. Si on prend l'extrême à 8 chiffres significatifs après la virgule, l'amplitude représente environ 4600 fois le tour de la Terre à l'équateur.

La phase de préparation des données dans un algorithme n'est jamais à négliger car cela peut apporter de très grandes simplifications de calculs mais aussi supprimer des calculs redondants au runtime.
13/05/2015, 15h15
MilWolf

Oui je n'ai pas réfléchi sur le coup, je vais optimiser de ce côté là.
Merci pour ces précisions.
14/05/2015, 16h27
souviron34

je ne sais pas comment fonctionne ta fonction Compare, ni comment sont tes points, mais moi pour optimiser je diviserais en 3 :

1er test par exemple sur Z (tu risques d'avoir moins d'altitudes que de x,y)

SI ce test réussi alors tu testes les 2 autres successivement

C'est ce que je disais plus haut.. Je pense que 3 dichotomies seraient beaucoup plus rapides que 1 seule globale (surtout pour un point ne figurant pas dans l'ensemble) : quand tu as identifié l'intervalle en Z, au lieu d'avoir log(N) sur 3 comparaisons, tu as log(N) sur une, + log(Nx) sur une autre et log (Ny) sur la 3ieme..

N = nbe total de points
Nx = nbe de points de x différents entre les 2 points Z (ou au point Z)
Ny = nbe de points de y différents entre les 2 points x (ou au point x)

C'est particulièrement valable si tu as qque chose comme un modèle 3D. Si réellement TOUS les points ont TOUS des coordonnées différentes alors ça ne gagnera rien.. Mais sinon tu auras un bon facteur d'accéleration,,,,
14/05/2015, 16h43
MilWolf

Bonjour souviron34,

en fait ma classe Point3f implémente l'interface Java "Comparable" et la comparaison ne se fait que sur une coordonnée sauf si les deux sont identiques, dans ce cas on teste la deuxième coordonnée, puis la troisième.

A++ !

MilWolf.
15/05/2015, 17h18
Graffito

Avec une fonctiion associant à x, y ou z un hashcode de n bits (n dépendant de la mémoire disponible), on peut remplir une matrice à 3 dimensions qui pointe sur la liste ordonée des points ayant le même hashcode et on prévoira une recherche dichotomique sur cette liste.
15/05/2015, 23h13
MilWolf

Bonsoir Graffito,

Peux tu développer? Quel est l'avantage en performance? Qu'est ce que cette matrice? Est-ce que c'est une technique connue?
De ce que j'ai compris ça permet de faire une recherche dichotomique sur une liste plus petite mais je ne vois pas l'avantage entre la recherche d'une cordonnée sur un axe ou de son hashcode associé.

Merci pour ton aide.
Bonne soirée.

MilWolf.

Le hashcode permet de "convertir" une coordonée float en entier codé sur n bits. On s'affranchit ainsi de la gestion des valeurs min et max sur l'axe correspondant.
La premiére étape d'une recherche (trouver une petite liste de points) est un accés à un tableau, donc extrémement rapide.
Pour une fonction de hash sur 8 bits et 20 millions de points, la recherche dichotomique se fait sur une liste contenant en moyenne 1.20 points (20 000 000 / 256 * 256 * 256).

Avec une fonction de hash sur 8 bits, le code (en C#) ressemblerait à ceci .
Fonctions à coder

Code:

1
2
3
4
bool ReadPoints(List<Point3D> TabPoints) ; // Lecture du nuage de points
int GetHashCode(float v) ; // retourne le hashcode
void SortPoints3D(List<Point3D> TheSmallPointsList) ; // trie une liste de points
bool IsPointInList(List<Point3D> TheSmallPointsList, Point3D) ; // recherche dichotomique

Initialisations

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
 
// données d'entrée : nuage de points à lire 
List<Point3D> TabPoints = new List<Point3D>(); 
ReadPoints(TabPoints) ;
 
//Initialisation de la matrice
int [,,] TabHash = new int[256,256,256] ;
for (int i=0;i<256;i++)  for (int j=0;j<256;j++)  for (int k=0;k<256;k++) TabHash[i,j,k]=null ;
 
// remplissage de la matrice
for (int n=0;n<TabPoints.Count;n++)
{
   int ii=GetHashCode(TabPoinst[n].x) ;
   int jj=GetHashCode(TabPoinst[n].y) ;
   int kk=GetHashCode(TabPoinst[n].z) ;
   // Retrouver la liste de points correspondant à l'entrée de la matrice,
   // Si non trouvée, la créer
   List<Point3D> TheSmallPointsList=TabHash[ii,jj,kk] ;
   if (TabHash[ii,jj,kk]==null) { TheSmallPointsList = new List<Point3D> () ; TabHash[ii,jj,kk]=TheSmallPointsList ; }
   TheSmallPointsList.Add((TabPoinst[n]) ;
}
 
// tri des petites listes pour future recherche dichotomique 
for (int i=0;i<256;i++)  for (int j=0;j<256;j++)  for (int k=0;k<256;k++)
  if (TabHash[i,j,k]!=null) SortPoints3D(TabHash[i,j,k]) ;

Code de recherche :

Code:

1
2
3
4
5
6
7
bool PointExists(Point3D ThePoint)
{
   int ii=GetHashCode(ThePoint) ;
   int jj=GetHashCode(ThePoint) ;
   int kk=GetHashCode(ThePoint) ;
   return IsPointInList(TabHash[ii,jj,kk],ThePoint) ;
}

Un algo (en C) de pour un hashcode 8 bits (la conversion du float en tableau de bytes sépend du langage de programmation.

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
//CRC-8 - based on the CRC8 formulas by Dallas/Maxim
//code released under the terms of the GNU GPL 3.0 license
byte CRC8(const byte *data, byte len) 
{
  byte crc = 0x00;
  while (len--) 
  {
    byte extract = *data++;
    for (byte tempI = 8; tempI; tempI--) 
    {
      byte sum = (crc ^ extract) & 0x01;
      crc >>= 1;
      if (sum) crc ^= 0x8C;
      extract >>= 1;
    }
  }
  return crc;
}

16/05/2015, 12h34
MilWolf

Bonjour Graffito,

merci pour tous tes efforts, c'est très complet et je suis sûr que ce code servira à beaucoup sur le forum.
Toutefois qu'en est-il de mes coordonnées à valeurs approximatives, l'erreur de précision que j'ai mentionné dans le post?
SI je regarde ton code, le hashcode est obtenu à partir de la valeur en flottant et si le hashcode n'existe pas le point n'existe pas.

Alors certes je pourrai adapter ton code pour faire une première recherche sur les points ayant exactement les mêmes coordonnées que recherchées, et si la coordonnée n'est pas retrouvée faire une recherche sur toute la liste. Mais les points ayant exactement les mêmes coordonnées ne représente qu'un peu moins d'1/10ème du nombre de valeurs.

Y-aurait-t'il quelque chose que je n'aurai pas compris?

Passe une bonne journée.

MilWolf.

Voir 40 message(s) de cette discussion en une page