Re: Schnelle Triangulation und konvexe Hülle in 2D
Verfasst: 12.01.2014 22:11
Ja danke.
Sollte ich so oder so optional von außerhalb anwenden, teils ist der Effekt auch gewünscht.

Das deutsche PureBasic-Forum
https://www.purebasic.fr/german/
Code: Alles auswählen
CompilerIf Not #PB_Compiler_Thread
CompilerError "Threadsicherer Modus muss aktiviert sein."
CompilerEndIf
#WIDTH = 1400
#HEIGHT = 800
;-- START OF TRIANGULATION STRUCTURES AND FUNCTIONS
DeclareModule Triangulation
EnableExplicit
#PRECISION = #PB_Float ;oder #PB_Double
#MAKE_DELAUNAY = #True
#IGNORE_SECOND_SORTING = #False
#MAX_FLIPPING = 0
#DEBUG = #True
CompilerIf #PRECISION = #PB_Float
Macro prec
f
EndMacro
#EPSILON = 0.00001
CompilerElseIf #PRECISION = #PB_Double
Macro prec
d
EndMacro
#EPSILON = 0.00000001
CompilerElse
CompilerError "Precision not supported."
CompilerEndIf
Structure Point2D
x.prec
y.prec
id.i
EndStructure
Structure Vector2D Extends Point2D
EndStructure
; Eine Kante merkt sich seine Endpunkte und die maximal zwei Dreicke, zu denen sie gehört
Structure Edge2D
*p.Point2DPC[2]
*t.Triangle2D[2]
mark.i ; Nur für makeDelauney. #True, wenn die Kante in der ToFlip-Liste ist.
flipped.i ; Nur für makeDelauney. Zähler wie häufig eine Kante geflipt wurde.
EndStructure
; Ein Punkt merkt sich natürlich seine Koordinaten und die Kanten und Dreiecken, zu denen er gehört.
Structure Point2DPC Extends Point2D
List *edges.Edge2D()
EndStructure
Structure CircumCircle
center.Point2D
;radius.prec
radiusSq.prec
EndStructure
; Ein Dreieck merkt sich seine drei Eckpunkte und seine drei Kanten
Structure Triangle2D
*p.Point2DPC[3]
*e.Edge2D[3]
cc.CircumCircle
EndStructure
Structure Point2DDiff
*p.Point2D
diff.prec
EndStructure
Structure ConvexHullPoint2D
*p.Point2DPC
used.i
EndStructure
Structure Triangulation
n.i
time.i
Array points.Point2DPC(1)
List edges.Edge2D()
List triangles.Triangle2D()
List convexHull.ConvexHullPoint2D()
nextTriangle.Triangle2D[2]
EndStructure
Prototype StepCallback(*pc.Triangulation, *userData)
Declare clearTriangulation(*pc.Triangulation)
Declare getCircumCircle(*a.Point2D, *b.Point2D, *c.Point2D, *cc.CircumCircle)
Declare.i getPoint(*pc.Triangulation, x.prec, y.prec)
Declare setPoint(*pc.Triangulation, i.i, x.prec, y.prec)
Declare.i isRightHand(*a.Point2D, *b.Point2D, *c.Point2D)
Declare.i isEdgeDelaunay(*edge.Edge2D)
Declare.i isTriangleDelaunay(*triangle.Triangle2D)
Declare triangulate(*pc.Triangulation, seed.i = 0, *cb.StepCallback = 0, *userData = 0)
EndDeclareModule
Module Triangulation
Procedure.s strT(*t.Triangle2D)
ProcedureReturn "T" + *t\p[0]\id + "-" + *t\p[1]\id + "-" + *t\p[2]\id
EndProcedure
Procedure clearTriangulation(*pc.Triangulation) ;correct
Protected i.i
ClearList(*pc\triangles())
ClearList(*pc\edges())
ClearList(*pc\convexHull())
For i = 0 To *pc\n - 1
ClearList(*pc\points(i)\edges())
Next
For i = 0 To 5
*pc\nextTriangle[i % 2]\p[i / 2] = 0
Next
EndProcedure
Procedure.i getPoint(*pc.Triangulation, x.prec, y.prec) ;correct
Protected i.i, diffSq.prec = Infinity(), actualDiffSq.prec
Protected nearest.i = -1
With *pc
For i = 0 To \n - 1
actualDiffSq = Pow(\points(i)\x - x, 2) + Pow(\points(i)\y - y, 2)
If (actualDiffSq < diffSq)
nearest = i
diffSq = actualDiffSq
EndIf
Next
EndWith
ProcedureReturn nearest
EndProcedure
Procedure setPoint(*pc.Triangulation, i.i, x.prec, y.prec) ;correct
If (i >= 0 And i < *pc\n)
*pc\points(i)\x = x
*pc\points(i)\y = y
EndIf
EndProcedure
Procedure getCircumCircle(*a.Point2D, *b.Point2D, *c.Point2D, *cc.CircumCircle) ;correct
Protected i.i
For i = 0 To 1
Protected m.Point2D
m\x = 0.5 * (*b\x - *c\x)
m\y = 0.5 * (*b\y - *c\y)
Protected v1.Vector2D
v1\x = *b\x - *a\x
v1\y = *b\y - *a\y
Protected v2.Vector2D
v2\x = *c\x - *a\x
v2\y = *c\y - *a\y
Protected lambda.prec = NaN()
Protected t.prec = v2\y * v1\x - v1\y * v2\x
If (Abs(t) > #EPSILON)
lambda = (v2\x * m\x + v2\y * m\y) / t
Break
EndIf
Swap *a, *b
Next
*cc\center\x = lambda * v1\y + 0.5 * (*a\x + *b\x)
*cc\center\y = - lambda * v1\x + 0.5 * (*a\y + *b\y)
*cc\radiusSq = Pow(*cc\center\x - *a\x, 2) + Pow(*cc\center\y - *a\y, 2)
If *cc\radiusSq >= 1.0 / (#EPSILON * #EPSILON)
Debug *cc\radiusSq
EndIf
ProcedureReturn Bool(*cc\radiusSq < 1.0 / (#EPSILON * #EPSILON))
EndProcedure
Procedure.i isRightHand(*a.Point2D, *b.Point2D, *c.Point2D) ;correct
Protected ab.Vector2D, bc.Vector2D
ab\x = *b\x - *a\x
ab\y = *b\y - *a\y
bc\x = *c\x - *b\x
bc\y = *c\y - *b\y
Protected z.prec = ab\x * bc\y - ab\y * bc\x
ProcedureReturn Bool(z < -#EPSILON)
EndProcedure
Procedure orderRightHand(*a.Point2D, *p_b.Integer, *p_c.Integer) ;correct
If (Not isRightHand(*a, *p_b\i, *p_c\i))
Swap *p_b\i, *p_c\i
EndIf
EndProcedure
Procedure.i addEdge(*pc.Triangulation, *a.Point2DPC, *b.Point2DPC) ;correct
;Zwei identische Punkte ergeben keine Linie
If (*a = *b)
ProcedureReturn #False
EndIf
;Existiert bereits eine Linie mit diesen beiden Punkten?
ForEach *a\edges()
;If (*a\edges()\p[0] = *b Or *a\edges()\p[1] = *b) ;changed
If (*a\edges()\p[0] = *a And *a\edges()\p[1] = *b) XOr (*a\edges()\p[0] = *b And *a\edges()\p[1] = *a)
ProcedureReturn *a\edges()
EndIf
Next
;Füge die Linie hinzu und speichere ihre Referenz in den Punkten
If AddElement(*pc\edges())
*pc\edges()\p[0] = *a
*pc\edges()\p[1] = *b
*pc\edges()\t[0] = 0
*pc\edges()\t[1] = 0
AddElement(*a\edges())
*a\edges() = @*pc\edges()
AddElement(*b\edges())
*b\edges() = @*pc\edges()
ProcedureReturn @*pc\edges()
EndIf
ProcedureReturn #False
EndProcedure
Procedure.i addTriangle(*pc.Triangulation, *a.Point2DPC, *b.Point2DPC, *c.Point2DPC, force.i = #False) ;correct
;Zwei identische Punkte ergeben kein Dreieck
If (*a = *b Or *b = *c Or *a = *c)
ProcedureReturn #False
EndIf
;Make sure the triangle's points are ordered in right hand order
orderRightHand(*a, @*b, @*c)
Protected cc.CircumCircle
If Not getCircumCircle(*a, *b, *c, cc)
If Not force
ProcedureReturn #False
EndIf
EndIf
Protected Dim *l.Edge2D(2)
*l(0) = addEdge(*pc, *a, *b)
*l(1) = addEdge(*pc, *b, *c)
*l(2) = addEdge(*pc, *c, *a)
;Prüfe, ob die drei Linien zufällig schon ein Dreieck bilden.
;Falls ja, dann gib einfach das zurück.
Protected j.i, k.i, *triangle.Triangle2D, used.i
For j = 0 To 1
used = 0
*triangle = *l(0)\t[j]
If (*triangle)
For k = 0 To 1
used + Bool(*l(1)\t[k] = *triangle)
used + Bool(*l(2)\t[k] = *triangle)
Next
If (used = 2)
Debug "Doppeltes Dreieck?"
ProcedureReturn *triangle
EndIf
EndIf
Next
If (AddElement(*pc\triangles()))
With *pc\triangles()
For j = 0 To 2
\e[j] = *l(j)
\e[j]\t[Bool(\e[j]\t[0])] = @*pc\triangles()
Next
\p[0] = *a
\p[1] = *b
\p[2] = *c
\cc = cc
EndWith
ProcedureReturn @*pc\triangles()
EndIf
ProcedureReturn #False
EndProcedure
Procedure.i isEdgeDelaunay(*edge.Edge2D) ;correct
Protected i.i, j.i, *p.Point2D
With *edge
If (\t[1])
For i = 0 To 1
For j = 0 To 2 ;Iterate over points of other triangle
*p = \t[1 - i]\p[j]
If (*p <> \p[0] And *p <> \p[1]) ;Is the actual point not belonging to the actual edge?
; Is this point from the other triangle within the circum circle of the actual triangle?
If (Pow(*p\x - \t[i]\cc\center\x, 2) + Pow(*p\y - \t[i]\cc\center\y, 2) - \t[i]\cc\radiusSq < -#EPSILON)
ProcedureReturn #False
EndIf
EndIf
Next
Next
EndIf
EndWith
ProcedureReturn #True
EndProcedure
Procedure.i isTriangleDelaunay(*triangle.Triangle2D) ;correct
Protected cc.CircumCircle, i.i, j.i, k.i, *edge.Edge2D, *p.Point2D
For k = 0 To 2
If (Not isEdgeDelaunay(*triangle\e[k]))
ProcedureReturn #False
EndIf
Next
ProcedureReturn #True
EndProcedure
Procedure makeDelaunay(*pc.Triangulation, *cb.StepCallback = 0, *userData = 0) ;correct
Protected NewList *ndEdges.Edge2D()
Protected Dim cc.CircumCircle(1), *p.Point2D
Protected Dim p_i.i(1) ;p_i(x) : Index to point of triangle x which not belongs to the actual edge
Protected i.i, j.i, doFlip.i
Protected *actualEdge.Edge2D, *newEdge.Edge2D
; Iteriere über alle Kanten
ForEach *pc\edges()
; Gehört die Kante zu zwei Dreiecken, prüfe, ob sie Delauney ist
If (*pc\edges()\t[1] And (Not isEdgeDelaunay(*pc\edges())))
; Wenn sie es nicht ist, füge sie zu der Prüfliste hinzu
If (AddElement(*ndEdges()))
*ndEdges() = @*pc\edges()
EndIf
; Und markiere sie
*pc\edges()\mark = #True
Else
*pc\edges()\mark = #False
EndIf
*pc\edges()\flipped = 0
Next
If *cb
*cb(*pc, *userData)
EndIf
While FirstElement(*ndEdges())
*actualEdge = *ndEdges()
DeleteElement(*ndEdges())
With *actualEdge
\mark = #False
doFlip = #False
For i = 0 To 1
For j = 0 To 2 ;Iterate over points of other triangle
*p = \t[1 - i]\p[j]
If (*p <> \p[0] And *p <> \p[1]) ;Is the actual point not belonging to the actual edge?
p_i(1 - i) = j
; Is this point from the other triangle within the circum circle of the actual triangle?
If Pow(*p\x - \t[i]\cc\center\x, 2) + Pow(*p\y - \t[i]\cc\center\y, 2) < \t[i]\cc\radiusSq - #EPSILON * #EPSILON
doFlip = #True
EndIf
EndIf
Next
Next
If (doFlip)
CopyStructure(\t[0], *pc\nextTriangle[0], Triangle2D)
CopyStructure(\t[1], *pc\nextTriangle[1], Triangle2D)
If *cb
*cb(*pc, *userData)
EndIf
CompilerIf #DEBUG
Protected error.i = #False
;Be sure p_i is correct
For i = 0 To 1
If (\t[i]\p[p_i(i)] = \p[0] Or \t[i]\p[p_i(i)] = \p[1])
Debug "p_i(" + i + ") ist nicht korrekt! (Points) [" + \t[i] + "]"
error = #True
EndIf
If (\t[i]\e[(p_i(i) + 1) % 3] <> *actualEdge)
Debug "p_i(" + i + ") ist nicht korrekt! (Edges) [" + \t[i] + "]"
error = #True
EndIf
If (Not isRightHand(\t[i]\p[0], \t[i]\p[1], \t[i]\p[2]))
Debug "Triangle " + strT(\t[i]) + " is not in the right order! [" + \t[i] + "]"
error = #True
EndIf
Next
If (error)
Debug "Error on " + *actualEdge\p[0]\id + "-" + *actualEdge\p[1]\id
ProcedureReturn #False
EndIf
Debug "flip on " + *actualEdge\p[0]\id + "-" + *actualEdge\p[1]\id + " With [" + strT(\t[0]) + "] And [" + strT(\t[1]) + "]"
CompilerEndIf
;Delete Edge from Points-Array
For i = 0 To 1
ForEach \p[i]\edges()
If (\p[i]\edges() = *actualEdge)
DeleteElement(\p[i]\edges())
Break
EndIf
Next
Next
;Swap points to create new triangles.
;This loop runs without error if the triangle's points are ordered in right hand order
;and the edge's order correlates to the points.
For i = 0 To 1
;Give actual edge the new coordinates
\p[i] = \t[i]\p[p_i(i)]
;Make the edge known to the point
AddElement(\p[i]\edges())
\p[i]\edges() = *actualEdge
;Change third point of triangle i to first point of triangle 1 - i
\t[i]\p[(p_i(i) + 2) % 3] = \t[1 - i]\p[p_i(1 - i)]
;Change second edge of triangle i to third edge of triangle 1 - i
\t[i]\e[(p_i(i) + 1) % 3] = \t[1 - i]\e[(p_i(1 - i) + 2) % 3]
;Correct neighbours of new second edge of triangle i
Protected *e.Edge2D
*e = \t[i]\e[(p_i(i) + 1) % 3]
For j = 0 To 1
If (*e\t[j] = \t[1 - i])
*e\t[j] = \t[i]
EndIf
Next
Next
For i = 0 To 1
;Change third edge of triangle i to actual edge
\t[i]\e[(p_i(i) + 2) % 3] = *actualEdge
Next
CompilerIf #DEBUG
If (addEdge(*pc, \t[0]\p[p_i(0)], \t[1]\p[p_i(1)]) <> *actualEdge)
Debug "actualEdge problem!"
EndIf
For i = 0 To 1
If (Not isRightHand(\t[i]\p[0], \t[i]\p[1], \t[i]\p[2]))
Debug "Triangle " + i + " is no more in the right order!"
EndIf
Next
CompilerEndIf
\flipped + 1
;Edge is now Delauney. Add adjacent Edges to List.
LastElement(*ndEdges())
For i = 0 To 1
For j = 0 To 1
*newEdge = \t[i]\e[(p_i(i) + j) % 3]
If ((Not *newEdge\mark) And *newEdge\t[1] And (#MAX_FLIPPING = 0 Or *newEdge\flipped < #MAX_FLIPPING))
If (AddElement(*ndEdges()))
*ndEdges() = *newEdge
*newEdge\mark = #True
EndIf
EndIf
Next
Next
; Update CircumCirles
For i = 0 To 1
If Not getCircumCircle(\t[i]\p[0], \t[i]\p[1], \t[i]\p[2], \t[i]\cc)
Debug "oh oh"
EndIf
Next
EndIf
EndWith
Wend
*pc\nextTriangle[0]\p[0] = 0
*pc\nextTriangle[0]\p[1] = 0
*pc\nextTriangle[0]\p[2] = 0
*pc\nextTriangle[1]\p[0] = 0
*pc\nextTriangle[1]\p[1] = 0
*pc\nextTriangle[1]\p[2] = 0
If *cb
*cb(*pc, *userData)
EndIf
EndProcedure
Procedure triangulate(*pc.Triangulation, seed.i = 0, *cb.StepCallback = 0, *userData = 0) ;correct
Protected Dim sortedPoints.Point2DDiff(*pc\n - 1)
Protected i.i
With *pc
\time = ElapsedMilliseconds()
If (\n < 3)
ProcedureReturn #False
EndIf
;1. Select a seed point x_0 from x_i.
If (seed < 0 Or seed >= \n)
ProcedureReturn #False
EndIf
Protected *x0.Point2DPC = @\points(seed)
clearTriangulation(*pc)
;TODO Hier würde eine lineare Suche reichen, da in Schrift 5 sowieso wieder sortiert wird
;2. Sort according to |x_i - x_0|^2.
For i = 0 To \n - 1
sortedPoints(i)\p = @\points(i)
sortedPoints(i)\diff = Pow(*x0\x - \points(i)\x, 2) + Pow(*x0\y - \points(i)\y, 2)
Next
SortStructuredArray(sortedPoints(), #PB_Sort_Ascending, OffsetOf(Point2DDiff\diff), #PRECISION)
;3. Find the point x_j closest to x_0.
Protected *xj.Point2DPC = sortedPoints(1)\p
;4. Find the point x_k that creates the smallest circum-circle
; with x_0 and x_j and record the center of the circum-circle C.
Protected bestIndex.i
Protected C.CircumCircle, bestC.CircumCircle
bestC\radiusSq = Infinity()
For i = 2 To \n - 1
getCircumCircle(*x0, *xj, sortedPoints(i)\p, @C)
If (C\radiusSq < bestC\radiusSq)
bestIndex = i
bestC = C
EndIf
; Wir suchen nur solange bis der (sortedPoints(i)-x_0)^2 < (bestC\radius * 2) ^ 2 ist
If sortedPoints(i)\diff > bestC\radiusSq * 4
Break
EndIf
Next
;5. Resort the remaining points according to |x_i - C|^2 to
; give points s_i.
If (bestIndex <> 2)
Swap sortedPoints(2)\p, sortedPoints(bestIndex)\p
EndIf
CompilerIf Not #IGNORE_SECOND_SORTING
For i = 3 To \n - 1
sortedPoints(i)\diff = Pow(bestC\center\x - sortedPoints(i)\p\x, 2) + Pow(bestC\center\y - sortedPoints(i)\p\y, 2)
Next
SortStructuredArray(sortedPoints(), #PB_Sort_Ascending, OffsetOf(Point2DDiff\diff), #PRECISION, 3, \n - 1)
CompilerEndIf
;6. Order point x_0, x_j, x_k to give a right handed system.
; This is the initial seed convex hull.
Protected *xk.Point2DPC = sortedPoints(2)\p
orderRightHand(*x0, @*xk, @*xj)
;7. Sequentially add the points s_i to the propagating 2D convex
; hull that is seeded with the triangle formed from x_0, x_j, x_k.
; As a new point is added the facets of the 2D-hull that are visible
; to it form new triangles.
addTriangle(*pc, *x0, *xk, *xj)
ClearList(\convexHull())
AddElement(\convexHull()) : \convexHull()\p = *x0
AddElement(\convexHull()) : \convexHull()\p = *xk
AddElement(\convexHull()) : \convexHull()\p = *xj
If *cb
*cb(*pc, *userData)
EndIf
Protected *p.Point2D, *last.ConvexHullPoint2D, alreadyAdded.i, convexHullSize.i, j.i
Protected *cEdges, *cTriangles, *cConvexHull
;Protected *bestLast.ConvexHullPoint2D, *best.ConvexHullPoint2D, bestRadius.prec = Infinity(), cc.CircumCircle
For i = 3 To \n - 1
*p = sortedPoints(i)\p
Debug "Punkt " + *p\id
alreadyAdded = #False
convexHullSize = ListSize(\convexHull())
FirstElement(\convexHull())
*last = @\convexHull()
NextElement(\convexHull())
For j = 0 To convexHullSize - 1
\nextTriangle[0]\p[0] = *last\p
\nextTriangle[0]\p[1] = *p
\nextTriangle[0]\p[2] = \convexHull()\p
If *cb
*cConvexHull = @\convexHull()
*cEdges = @\edges()
*cTriangles = @\triangles()
*cb(*pc, *userData)
ChangeCurrentElement(\convexHull(), *cConvexHull)
ChangeCurrentElement(\edges(), *cEdges)
ChangeCurrentElement(\triangles(), *cTriangles)
EndIf
If (isRightHand(*last\p, *p, \convexHull()\p))
If addTriangle(*pc, *last\p, *p, \convexHull()\p)
Debug " + ADD " + *last\p\id + "-" + *p\id + "-" + \convexHull()\p\id
\convexHull()\used + 1
*last\used + 1
If (Not alreadyAdded)
InsertElement(\convexHull())
\convexHull()\p = *p
NextElement(\convexHull())
convexHullSize + 1
alreadyAdded = #True
EndIf
Else
Debug " - not"
EndIf
Else
Debug " - not " + *last\p\id + "-" + *p\id + "-" + \convexHull()\p\id
EndIf
*last = @\convexHull()
If (Not NextElement(\convexHull()))
FirstElement(\convexHull())
EndIf
Next
;TODO vielleicht muss das schon viel früher ausgeführt werden.
ForEach \convexHull()
If (\convexHull()\used >= 2)
DeleteElement(\convexHull(), 1)
Else
\convexHull()\used = 0
EndIf
Next
Next
\nextTriangle[0]\p[0] = 0
\nextTriangle[0]\p[1] = 0
\nextTriangle[0]\p[2] = 0
;8. A non-overlapping triangulation of the set of points is created.
; (This is an extremely fast method for creating an non-overlapping
; triangualtion of a 2D point set).
;9: Adjacent pairs of triangles of this triangulation must be 'flipped'
; to create a Delaunay triangulation from the initial non-overlapping
; triangulation.
CompilerIf #MAKE_DELAUNAY
makeDelaunay(*pc, *cb, *userData)
CompilerEndIf
\time = ElapsedMilliseconds() - \time
EndWith
EndProcedure
EndModule
;-- END OF TRIANGULATION FUNCTIONS
Structure Window
id.i
width.i
height.i
title.s
canvasId.i
*pc.Triangulation::Triangulation
clicked.i
leftDown.i
rightDown.i
mutex.i
thread.i
newThread.i
EndStructure
#MAX_INTEGER = 1 << (SizeOf(Integer) * 8 - 1) - 1
#MIN_INTEGER = ~#MAX_INTEGER
;Erstellt zufällige Punkte
Procedure CreateRandomizedPointCloud(*pc.Triangulation::Triangulation, n.i, minX.d = 0.0, minY.d = 0.0, maxX.d = 1.0, maxY.d = 1.0)
Protected i.i, w.i
With *pc
Triangulation::clearTriangulation(*pc)
\n = n
ReDim \points(\n - 1)
w = Round(Sqr(\n), #PB_Round_Up)
For i = 0 To \n - 1
\points(i)\id = i
\points(i)\x = (Random(#MAX_INTEGER) * (maxX - minX)) / #MAX_INTEGER + minX
\points(i)\y = (Random(#MAX_INTEGER) * (maxY - minY)) / #MAX_INTEGER + minY
\points(i)\x = Mod(i * (maxX - minX) * (maxY - minY) / \n, maxX - minX)
\points(i)\y = Mod((i * (maxX - minX) * (maxY - minY) / \n) / (maxX - minX), maxY - minY)
;\points(i)\x = (0.5 + Mod(i, w)) * (maxX - minX) / w
;\points(i)\y = (0.5 + Int(i / w)) * (maxY - minY) / w
Next
EndWith
EndProcedure
#VECTORDRAWING = #True
Procedure DrawPoints(*main.Window)
Protected i.i, width.i, height.i
Static font.i
With *main\pc
CompilerIf #VECTORDRAWING
If Not IsFont(font)
font = LoadFont(#PB_Any, "Courier New", 16, #PB_Font_Bold)
EndIf
If StartVectorDrawing(CanvasVectorOutput(*main\canvasId))
VectorFont(FontID(font), 16)
CompilerElse
If StartDrawing(CanvasOutput(*main\canvasId))
DrawingMode(#PB_2DDrawing_Default)
CompilerEndIf
width = GadgetWidth(*main\canvasId)
height = GadgetHeight(*main\canvasId)
CompilerIf #VECTORDRAWING
VectorSourceColor(RGBA(255, 255, 255, 255))
FillVectorOutput()
CompilerElse
Box(0, 0, GadgetWidth(*main\canvasId), GadgetHeight(*main\canvasId), $ffffff)
DrawingMode(#PB_2DDrawing_Outlined)
CompilerEndIf
CompilerIf #VECTORDRAWING
If LastElement(\convexHull())
MovePathCursor(\convexHull()\p\x, \convexHull()\p\y)
VectorSourceColor(RGBA($ff, 0, 0, $3f))
ForEach \convexHull()
AddPathLine(\convexHull()\p\x, \convexHull()\p\y)
Next
DashPath(5, 10, #PB_Path_RoundCorner)
EndIf
CompilerElse
Protected ConvexHullSize.i = ListSize(\convexHull())
If (ConvexHullSize > 2)
Protected *last.Triangulation::ConvexHullPoint2D = 0
LastElement(\convexHull())
*last = @\convexHull()
i = 0
ForEach \convexHull()
LineXY(*last\p\x, *last\p\y, \convexHull()\p\x, \convexHull()\p\y, $0000ff)
DrawText((*last\p\x + \convexHull()\p\x) / 2, (*last\p\y + \convexHull()\p\y) / 2, Str(i), 0)
i + 1
*last = @\convexHull()
Next
EndIf
CompilerEndIf
Protected cc.Triangulation::CircumCircle, color.i
ForEach \triangles()
If (Triangulation::isTriangleDelaunay(@\triangles()))
color = RGBA(0, $ff, 0, $5f)
Else
color = RGBA($ff, 0, 0, $5f)
EndIf
Triangulation::getCircumCircle(\triangles()\p[0], \triangles()\p[1], \triangles()\p[2], @cc)
CompilerIf #VECTORDRAWING
VectorSourceColor(color)
AddPathCircle(cc\center\x, cc\center\y, Sqr(cc\radiusSq))
StrokePath(1)
CompilerElse
; Ohne das hier, kann Kreise malen ganz schön langsam sein
If cc\center\x >= 0 And cc\center\x < width And cc\center\y >= 0 And cc\center\y < height
Circle(cc\center\x, cc\center\y, Sqr(cc\radiusSq), color)
EndIf
CompilerEndIf
Next
ForEach \edges()
If \edges()\mark
color = RGBA(0, 0, $ff, $ff)
Else
color = RGBA($7f, $7f, $7f, $ff)
EndIf
CompilerIf #VECTORDRAWING
MovePathCursor(\edges()\p[0]\x, \edges()\p[0]\y)
VectorSourceColor(color)
AddPathLine(\edges()\p[1]\x, \edges()\p[1]\y)
StrokePath(1 + \edges()\mark * 2)
CompilerElse
LineXY(\edges()\p[0]\x, \edges()\p[0]\y, \edges()\p[1]\x, \edges()\p[1]\y, color)
CompilerEndIf
Next
CompilerIf Not #VECTORDRAWING
DrawingMode(#PB_2DDrawing_Transparent)
CompilerEndIf
For i = 0 To \n - 1
Protected s.s = Str(\points(i)\id)
CompilerIf #VECTORDRAWING
If \points(i) = \nextTriangle\p[1]
VectorSourceColor(RGBA($7f, $7f, $ff, $ff))
AddPathCircle(\points(i)\x, \points(i)\y, 10)
FillPath()
color = RGBA($ff, $3f, $3f, $ff)
Else
color = RGBA(0, 0, 0, $ff)
EndIf
VectorSourceColor(color)
AddPathCircle(\points(i)\x, \points(i)\y, 8)
FillPath()
MovePathCursor(\points(i)\x - VectorTextWidth(s) / 2 + 1, \points(i)\y - VectorTextHeight(s) / 2 + 1)
VectorSourceColor(RGBA($ff, $ff, $ff, $ff))
DrawVectorText(s)
CompilerElse
Circle(\points(i)\x, \points(i)\y, 1, $000000)
DrawText(\points(i)\x, \points(i)\y, s, $ff0000)
CompilerEndIf
Next
CompilerIf #VECTORDRAWING
For i = 0 To 1
If \nextTriangle[i]\p[0] And \nextTriangle[i]\p[1] And \nextTriangle[i]\p[2]
If Triangulation::isRightHand(\nextTriangle[i]\p[0], \nextTriangle[i]\p[1], \nextTriangle[i]\p[2])
color = RGBA($7f, $ff, $7f, $7f)
Else
color = RGBA($ff, $7f, $7f, $7f)
EndIf
MovePathCursor(\nextTriangle[i]\p[0]\x, \nextTriangle[i]\p[0]\y)
AddPathLine(\nextTriangle[i]\p[1]\x, \nextTriangle[i]\p[1]\y)
AddPathLine(\nextTriangle[i]\p[2]\x, \nextTriangle[i]\p[2]\y)
ClosePath()
VectorSourceColor(color)
FillPath()
EndIf
Next
CompilerEndIf
CompilerIf Not #VECTORDRAWING
DrawingMode(#PB_2DDrawing_Default)
DrawText(0, height - 20, " Points: " + Str(\n) + " " +
"Edges: " + Str(ListSize(\edges())) + " " +
"Triangles: " + Str(ListSize(\triangles())) + " " +
"Time: " + Str(\time) + " ms ", $0000ff, $ffffff)
CompilerEndIf
CompilerIf #VECTORDRAWING
StopVectorDrawing()
CompilerElse
StopDrawing()
CompilerEndIf
EndIf
EndWith
ProcedureReturn #True
EndProcedure
Procedure callback(*pc, *main.Window)
PostEvent(#PB_Event_FirstCustomValue, *main\id, 0)
WaitSemaphore(*main\mutex)
EndProcedure
Procedure triangulate(*main.Window)
Triangulation::triangulate(*main\pc, 0, @callback(), *main)
EndProcedure
Procedure CanvasEvent()
Protected x.i, y.i, gadgetId.i, *main.Window, seed.i
gadgetId = EventGadget()
*main = GetGadgetData(gadgetId)
With *main
x = GetGadgetAttribute(\canvasId, #PB_Canvas_MouseX)
y = GetGadgetAttribute(\canvasId, #PB_Canvas_MouseY)
Select (EventType())
Case #PB_EventType_RightClick
If Not \newThread
If IsThread(\thread)
SignalSemaphore(\mutex)
EndIf
EndIf
\newThread = #False
Case #PB_EventType_LeftButtonDown
\leftDown = #True
\clicked = Triangulation::getPoint(\pc, x, y)
Case #PB_EventType_Input
If GetGadgetAttribute(\canvasId, #PB_Canvas_Input) = 'f'
If IsThread(\thread)
SignalSemaphore(\mutex)
EndIf
EndIf
Case #PB_EventType_RightButtonDown
\rightDown = #True
If Not IsThread(\thread)
\thread = CreateThread(@triangulate(), *main)
\newThread = #True
EndIf
;Triangulation::triangulate(\pc, Triangulation::getPoint(\pc, x, y), @callback(), *main)
;DrawPoints(*main)
Case #PB_EventType_MouseMove
If (\leftDown)
Triangulation::setPoint(\pc, \clicked, x, y)
DrawPoints(*main)
EndIf
If Not IsThread(\thread)
If (\rightDown)
If (\clicked >= 0)
seed = \clicked
Else
seed = Triangulation::getPoint(\pc, x, y)
EndIf
Triangulation::triangulate(\pc, seed)
EndIf
;DrawPoints(*main)
EndIf
Case #PB_EventType_LeftButtonUp
\leftDown = #False
\clicked = -1
Case #PB_EventType_RightButtonUp
\rightDown = #False
\clicked = -1
EndSelect
EndWith
EndProcedure
Procedure Redraw()
Protected *main.Window = GetWindowData(EventWindow())
If Event() = #PB_Event_FirstCustomValue
DrawPoints(*main)
EndIf
EndProcedure
Procedure CreateMainWindow(*main.Window)
With *main
\id = OpenWindow(#PB_Any, 0, 0, \width, \height, \title, #PB_Window_MinimizeGadget | #PB_Window_SystemMenu | #PB_Window_ScreenCentered)
If (Not \id) : ProcedureReturn #False : EndIf
SetWindowData(\id, *main)
\canvasId = CanvasGadget(#PB_Any, 0, 0, \width, \height, #PB_Canvas_Keyboard)
\clicked = -1
SetGadgetData(\canvasId, *main)
If (Not \canvasId)
CloseWindow(\id)
ProcedureReturn #False
EndIf
BindGadgetEvent(\canvasId, @CanvasEvent())
BindEvent(#PB_Event_FirstCustomValue, @Redraw(), \id)
ProcedureReturn \id
EndWith
EndProcedure
Define.Triangulation::Triangulation pc
CreateRandomizedPointCloud(@pc, 53, 0, 0, #WIDTH - 1, #HEIGHT - 1)
Define.Window main
main\width = #WIDTH
main\height = #HEIGHT
main\title = "Triangulation Test"
main\pc = @pc
main\mutex = CreateSemaphore(0)
CreateMainWindow(@main)
DrawPoints(@main)
Repeat : Until WaitWindowEvent() = #PB_Event_CloseWindow