Mercurial > gemma

comparison pkg/octree/triangulation.go @ 2484:4fa92d468164 octree-diff

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Use fogleman's triangulation.
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Fri, 01 Mar 2019 17:42:46 +0100
parents 620038ade708
children 6248a4bc10c7
comparison
equal deleted inserted replaced
2483:620038ade708 2484:4fa92d468164
25 ) 25 )
26 26
27 type Triangulation struct { 27 type Triangulation struct {
28 Points []Vertex 28 Points []Vertex
29 ConvexHull []Vertex 29 ConvexHull []Vertex
30 Triangles []int 30 Triangles []int32
31 Halfedges []int 31 Halfedges []int32
32 } 32 }
33 33
34 // Triangulate returns a Delaunay triangulation of the provided points. 34 // Triangulate returns a Delaunay triangulation of the provided points.
35 func Triangulate(points []Vertex) (*Triangulation, error) { 35 func Triangulate(points []Vertex) (*Triangulation, error) {
36 t := newTriangulator(points) 36 t := newTriangulator(points)
37 err := t.triangulate() 37 err := t.triangulate()
38 return &Triangulation{points, t.convexHull(), t.triangles, t.halfedges}, err 38 return &Triangulation{points, t.convexHull(), t.triangles, t.halfedges}, err
39 }
40
41 func (t *Triangulation) TriangleSlices() [][]int32 {
42 sl := make([][]int32, len(t.Triangles)/3)
43 var j int
44 for i := range sl {
45 sl[i] = t.Triangles[j : j+3]
46 j += 3
47 }
48 return sl
49 }
50
51 func (t *Triangulation) Tin() *Tin {
52
53 min := Vertex{math.MaxFloat64, math.MaxFloat64, math.MaxFloat64}
54 max := Vertex{-math.MaxFloat64, -math.MaxFloat64, -math.MaxFloat64}
55
56 vertices := t.Points
57
58 for _, v := range vertices {
59 min.Minimize(v)
60 max.Maximize(v)
61 }
62
63 return &Tin{
64 Vertices: vertices,
65 Triangles: t.TriangleSlices(),
66 Min: min,
67 Max: max,
68 }
39 } 69 }
40 70
41 func (t *Triangulation) area() float64 { 71 func (t *Triangulation) area() float64 {
42 var result float64 72 var result float64
43 points := t.Points 73 points := t.Points
58 // verify halfedges 88 // verify halfedges
59 for i1, i2 := range t.Halfedges { 89 for i1, i2 := range t.Halfedges {
60 if i1 != -1 && t.Halfedges[i1] != i2 { 90 if i1 != -1 && t.Halfedges[i1] != i2 {
61 return fmt.Errorf("invalid halfedge connection") 91 return fmt.Errorf("invalid halfedge connection")
62 } 92 }
63 if i2 != -1 && t.Halfedges[i2] != i1 { 93 if i2 != -1 && t.Halfedges[i2] != int32(i1) {
64 return fmt.Errorf("invalid halfedge connection") 94 return fmt.Errorf("invalid halfedge connection")
65 } 95 }
66 } 96 }
67 97
68 // verify convex hull area vs sum of triangle areas 98 // verify convex hull area vs sum of triangle areas