299 lines
10 KiB
Python
299 lines
10 KiB
Python
# coding: utf-8
|
|
|
|
import sys
|
|
import math
|
|
|
|
|
|
def distanceAuCarre(v1, v2):
|
|
dx = (v1.x() - v2.x())
|
|
dy = (v1.y() - v2.y())
|
|
return dx**2 + dy**2
|
|
|
|
|
|
def conditionDeDelaunay(A, B, C, D):
|
|
dab2 = distanceAuCarre(A, B)
|
|
dad2 = distanceAuCarre(A, D)
|
|
dbc2 = distanceAuCarre(B, C)
|
|
dbd2 = distanceAuCarre(B, D)
|
|
dcd2 = distanceAuCarre(C, D)
|
|
cosa = (dab2 + dad2 - dbd2) / (2 * math.sqrt(dab2) * math.sqrt(dad2))
|
|
cosc = (dbc2 + dcd2 - dbd2) / (2 * math.sqrt(dbc2) * math.sqrt(dcd2))
|
|
sina2 = 1. - cosa * cosa
|
|
if(sina2 > 0.):
|
|
sina = math.sqrt(sina2)
|
|
else:
|
|
sina = 0.
|
|
sinc2 = 1. - cosc * cosc
|
|
if(sinc2 > 0.):
|
|
sinc = math.sqrt(sinc2)
|
|
else:
|
|
sinc = 0.
|
|
return sina * cosc + cosa * sinc >= 0
|
|
|
|
|
|
class Triangulation:
|
|
|
|
class QuadEdge:
|
|
'''Classe définissant un quad-edge : une arête connectée à deux faces et deux sommets'''
|
|
|
|
def __init__(self, f1, f2, v1, v2):
|
|
self.f1 = f1
|
|
self.f2 = f2
|
|
self.v1 = v1
|
|
self.v2 = v2
|
|
|
|
def __eq__(self, e):
|
|
'''Méthode qui teste si deux quad-edges ont même extrémités'''
|
|
if((self.v1 == e.v1) and (self.v2 == e.v2)):
|
|
return True
|
|
elif((self.v1 == e.v2) and (self.v2 == e.v1)):
|
|
return True
|
|
else:
|
|
return False
|
|
|
|
def __str__(self):
|
|
return str(self.v1) + " - " + str(self.v2)
|
|
|
|
def remplaceFace(self, ancienneFace, nouvelleFace):
|
|
'''Méthode qui remplace une face adjacente au quad-edge par une autre'''
|
|
if(self.f1 == ancienneFace):
|
|
self.f1 = nouvelleFace
|
|
else:
|
|
self.f2 = nouvelleFace
|
|
|
|
def lanceAlgorithme(self):
|
|
# Construit la triangulation initiale. Complexité en O(n log(n))
|
|
self.construitTriangulationInitiale()
|
|
|
|
# Bascule les arêtes de la triangulation jusqu'à obtenir la triangulation de Delaunay
|
|
self.basculeAretes()
|
|
|
|
# Efface les sommets à l'infini
|
|
for _ in range(3):
|
|
self.graphe.sommets.pop()
|
|
|
|
def __init__(self, g):
|
|
'''Constructeur de la classe triangulation qui prend en entrée un nuage de points
|
|
et crée une triangulation de Delaunay'''
|
|
self.racine = None
|
|
self.quadEdges = []
|
|
self.graphe = g
|
|
self.lanceAlgorithme()
|
|
|
|
class Triangle:
|
|
def __init__(self, v1, v2, v3):
|
|
self.v1 = v1
|
|
self.v2 = v2
|
|
self.v3 = v3
|
|
self.t12 = None
|
|
self.t13 = None
|
|
self.t23 = None
|
|
self.e12 = None
|
|
self.e13 = None
|
|
self.e23 = None
|
|
|
|
def contient(self, v):
|
|
''' Méthode qui teste si un sommet est à l'intérieur d'un triangle'''
|
|
v1 = self.v1
|
|
v2 = self.v2
|
|
v3 = self.v3
|
|
det = (v1.x() - v3.x()) * (v2.y() - v3.y()) - \
|
|
(v1.y() - v3.y()) * (v2.x() - v3.x())
|
|
c1 = ((v2.y() - v3.y()) * (v.x() - v3.x()) +
|
|
(v3.x() - v2.x()) * (v.y() - v3.y())) * det
|
|
c2 = ((v3.y() - v1.y()) * (v.x() - v3.x()) +
|
|
(v1.x() - v3.x()) * (v.y() - v3.y())) * det
|
|
return (c1 > 0) and (c2 > 0) and (c1 + c2 < det * det)
|
|
|
|
def trouve(self, v):
|
|
'''Méthode récursive qui renvoie le triangle d'une triangulation qui contient un sommet.
|
|
Renvoie le triangle contenant v, None si v n'est pas dans la triangulation'''
|
|
|
|
t12 = self.t12
|
|
t13 = self.t13
|
|
t23 = self.t23
|
|
if(t12 == None):
|
|
return self
|
|
else:
|
|
if(t12.contient(v)):
|
|
return t12.trouve(v)
|
|
elif(t13.contient(v)):
|
|
return t13.trouve(v)
|
|
elif(t23.contient(v)):
|
|
return t23.trouve(v)
|
|
else:
|
|
return None
|
|
|
|
def aPourSommet(self, v):
|
|
'''Méthode qui teste si un sommet est un sommet d'un triangle
|
|
Renvoie vrai si v est un sommet du triangle'''
|
|
return self.v1 == v or self.v2 == v or self.v3 == v
|
|
|
|
def troisiemeSommet(self, a, b):
|
|
''' Renvoie le troisième sommet d'un triangle à partir des deux autres'''
|
|
|
|
if(self.v1 != a and self.v1 != b):
|
|
return self.v1
|
|
elif(self.v2 != a and self.v2 != b):
|
|
return self.v2
|
|
else:
|
|
return self.v3
|
|
|
|
def remplaceSommet(self, ancienSommet, nouveauSommet):
|
|
''' Méthode qui substitue un sommet par un autre dans un triangle '''
|
|
if(self.v1 == ancienSommet):
|
|
self.v1 = nouveauSommet
|
|
elif(self.v2 == ancienSommet):
|
|
self.v2 = nouveauSommet
|
|
elif(self.v3 == ancienSommet):
|
|
self.v3 = nouveauSommet
|
|
else:
|
|
raise Exception("Bad vertex" + ancienSommet)
|
|
|
|
def retrouveArete(self, v1, v2):
|
|
'''Méthode qui renvoie le quad-edge d'un triangle ayant deux sommets pour extrêmités'''
|
|
e = Triangulation.QuadEdge(None, None, v1, v2)
|
|
if(self.e12 == e):
|
|
return self.e12
|
|
elif(self.e13 == e):
|
|
return self.e13
|
|
elif(self.e23 == e):
|
|
return self.e23
|
|
else:
|
|
raise Exception("getEdge")
|
|
|
|
def remplaceArete(self, ancienneArete, nouvelleArete):
|
|
'''Méthode qui substitue le quad-edge d'un triangle par un autre'''
|
|
if(self.e12 == ancienneArete):
|
|
self.e12 = nouvelleArete
|
|
elif(self.e13 == ancienneArete):
|
|
self.e13 = nouvelleArete
|
|
elif(self.e23 == ancienneArete):
|
|
self.e23 = nouvelleArete
|
|
else:
|
|
raise Exception("changeEdge")
|
|
|
|
def construitTriangleInitial(self):
|
|
'''Méthode qui construit une triangulation initiale qui contient tous les sommets'''
|
|
xmin = sys.float_info.max
|
|
xmax = sys.float_info.min
|
|
ymin = sys.float_info.max
|
|
ymax = sys.float_info.min
|
|
|
|
for v in self.graphe.sommets:
|
|
if(v.x() > xmax):
|
|
xmax = v.x()
|
|
if(v.x() < xmin):
|
|
xmin = v.x()
|
|
if(v.y() > ymax):
|
|
ymax = v.y()
|
|
if(v.y() < ymin):
|
|
ymin = v.y()
|
|
|
|
largeur = xmax - xmin
|
|
hauteur = ymax - ymin
|
|
margex = 0.1 * largeur
|
|
margey = 0.1 * hauteur
|
|
|
|
v1 = self.graphe.ajouteSommet(
|
|
xmin - largeur / 2 - margex, ymin - margey)
|
|
v2 = self.graphe.ajouteSommet(
|
|
xmax + largeur / 2 + margex, ymin - margey)
|
|
v3 = self.graphe.ajouteSommet(
|
|
(xmin + xmax) / 2, ymax + hauteur + margey)
|
|
|
|
self.racine = Triangulation.Triangle(v1, v2, v3)
|
|
|
|
e12 = Triangulation.QuadEdge(None, self.racine, v1, v2)
|
|
e13 = Triangulation.QuadEdge(None, self.racine, v1, v3)
|
|
e23 = Triangulation.QuadEdge(None, self.racine, v2, v3)
|
|
self.racine.e12 = e12
|
|
self.racine.e13 = e13
|
|
self.racine.e23 = e23
|
|
|
|
def supprimeAretes(self):
|
|
for s in self.graphe.sommets:
|
|
s.aretes.clear()
|
|
|
|
def construitTriangulationInitiale(self):
|
|
self.quadEdges = []
|
|
self.construitTriangleInitial()
|
|
|
|
for v in self.graphe.sommets:
|
|
t = self.racine.trouve(v)
|
|
if(t != None):
|
|
t.t12 = Triangulation.Triangle(t.v1, t.v2, v)
|
|
t.t13 = Triangulation.Triangle(t.v1, t.v3, v)
|
|
t.t23 = Triangulation.Triangle(t.v2, t.v3, v)
|
|
|
|
e1v = Triangulation.QuadEdge(t.t12, t.t13, t.v1, v)
|
|
e2v = Triangulation.QuadEdge(t.t12, t.t23, t.v2, v)
|
|
e3v = Triangulation.QuadEdge(t.t13, t.t23, t.v3, v)
|
|
self.quadEdges.append(e1v)
|
|
self.quadEdges.append(e2v)
|
|
self.quadEdges.append(e3v)
|
|
|
|
t.t12.e12 = t.e12
|
|
t.e12.remplaceFace(t, t.t12)
|
|
t.t12.e13 = e1v
|
|
t.t12.e23 = e2v
|
|
|
|
t.t13.e12 = t.e13
|
|
t.e13.remplaceFace(t, t.t13)
|
|
t.t13.e13 = e1v
|
|
t.t13.e23 = e3v
|
|
|
|
t.t23.e12 = t.e23
|
|
t.e23.remplaceFace(t, t.t23)
|
|
t.t23.e13 = e2v
|
|
t.t23.e23 = e3v
|
|
|
|
def basculeAretes(self):
|
|
pile = []
|
|
for qe in self.quadEdges:
|
|
pile.append(qe)
|
|
|
|
while(len(pile) > 0):
|
|
e = pile.pop()
|
|
|
|
f1 = e.f1
|
|
f2 = e.f2
|
|
if(f1 != None and f2 != None):
|
|
v1 = e.v1
|
|
v2 = e.v2
|
|
|
|
v3 = f1.troisiemeSommet(v1, v2)
|
|
v4 = f2.troisiemeSommet(v1, v2)
|
|
|
|
if(not conditionDeDelaunay(v3, v1, v4, v2)):
|
|
v1v3 = f1.retrouveArete(v1, v3)
|
|
v2v3 = f1.retrouveArete(v2, v3)
|
|
v1v4 = f2.retrouveArete(v1, v4)
|
|
v2v4 = f2.retrouveArete(v2, v4)
|
|
f1.remplaceArete(v2v3, v1v4)
|
|
f2.remplaceArete(v1v4, v2v3)
|
|
v2v3.remplaceFace(f1, f2)
|
|
v1v4.remplaceFace(f2, f1)
|
|
f1.remplaceSommet(v2, v4)
|
|
f2.remplaceSommet(v1, v3)
|
|
e.v1 = v3
|
|
e.v2 = v4
|
|
pile.append(v1v3)
|
|
pile.append(v2v3)
|
|
pile.append(v1v4)
|
|
pile.append(v2v4)
|
|
|
|
def construitGrapheDeSortie(self, selectionneArete):
|
|
'''Méthode qui construit la triangulation de Delaunay à partir du nuage de points passé au constructeur'''
|
|
|
|
# Sélectionne les arêtes
|
|
graphe = self.graphe
|
|
racine = self.racine
|
|
for qe in self.quadEdges:
|
|
v1 = qe.v1
|
|
v2 = qe.v2
|
|
if(not (racine.aPourSommet(v1) or racine.aPourSommet(v2))):
|
|
v3 = qe.f1.troisiemeSommet(v1, v2)
|
|
v4 = qe.f2.troisiemeSommet(v1, v2)
|
|
selectionneArete(graphe, v1, v2, v3, v4)
|
|
return self.graphe
|