diff pkg/octree/triangulation.go @ 3639:2a4216c81e7b single-beam

Extract the removed triangles from first triangulation, too. Useful to build a artifical DEM for second pass.
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Wed, 12 Jun 2019 12:56:18 +0200
parents db7136854a51
children ec86a7155377
line wrap: on
line diff
--- a/pkg/octree/triangulation.go	Tue Jun 11 18:20:41 2019 +0200
+++ b/pkg/octree/triangulation.go	Wed Jun 12 12:56:18 2019 +0200
@@ -39,7 +39,7 @@
 	return &Triangulation{points, t.convexHull(), t.triangles, t.halfedges}, err
 }
 
-func (t *Triangulation) ConcaveHull(tooLong float64) []int32 {
+func (t *Triangulation) ConcaveHull(tooLong float64) (LineStringZ, map[int32]struct{}) {
 
 	tooLong *= tooLong
 
@@ -60,13 +60,17 @@
 		}
 	}
 
-	removed := map[int32]bool{}
+	removed := map[int32]struct{}{}
 
 	isBorder := func(n int32) bool {
 		n *= 3
 		for i := int32(0); i < 3; i++ {
 			e := n + i
-			if o := t.Halfedges[e]; o < 0 || removed[o/3] {
+			o := t.Halfedges[e]
+			if o < 0 {
+				return true
+			}
+			if _, found := removed[o/3]; found {
 				return true
 			}
 		}
@@ -83,7 +87,7 @@
 		for _, i := range candidates {
 
 			if isBorder(i) {
-				removed[i] = true
+				removed[i] = struct{}{}
 			} else {
 				newCandidates = append(newCandidates, i)
 			}
@@ -119,51 +123,55 @@
 	var rings []*edge
 
 	for i, num := int32(0), int32(len(t.Triangles)/3); i < num; i++ {
-		if removed[i] {
+		if _, found := removed[i]; found {
 			continue
 		}
 		n := i * 3
 		for j := int32(0); j < 3; j++ {
 			e := n + j
-			if f := t.Halfedges[e]; f < 0 || removed[f/3] {
-				a := t.Triangles[e]
-				b := t.Triangles[n+(j+1)%3]
-
-				curr := &edge{a: a, b: b}
+			f := t.Halfedges[e]
+			if f >= 0 {
+				if _, found := removed[f/3]; !found {
+					continue
+				}
+			}
+			a := t.Triangles[e]
+			b := t.Triangles[n+(j+1)%3]
 
-				if old := open[a]; old != nil {
-					delete(open, a)
-					if old.a == a {
-						old.prev = curr
-						curr.next = old
-					} else {
-						old.next = curr
-						curr.prev = old
-					}
+			curr := &edge{a: a, b: b}
 
-					if isClosed(old) {
-						rings = append(rings, old)
-					}
+			if old := open[a]; old != nil {
+				delete(open, a)
+				if old.a == a {
+					old.prev = curr
+					curr.next = old
 				} else {
-					open[a] = curr
+					old.next = curr
+					curr.prev = old
 				}
 
-				if old := open[b]; old != nil {
-					delete(open, b)
-					if old.b == b {
-						old.next = curr
-						curr.prev = old
-					} else {
-						old.prev = curr
-						curr.next = old
-					}
+				if isClosed(old) {
+					rings = append(rings, old)
+				}
+			} else {
+				open[a] = curr
+			}
 
-					if isClosed(old) {
-						rings = append(rings, old)
-					}
+			if old := open[b]; old != nil {
+				delete(open, b)
+				if old.b == b {
+					old.next = curr
+					curr.prev = old
 				} else {
-					open[b] = curr
+					old.prev = curr
+					curr.next = old
 				}
+
+				if isClosed(old) {
+					rings = append(rings, old)
+				}
+			} else {
+				open[b] = curr
 			}
 		}
 	}
@@ -174,19 +182,25 @@
 
 	if len(rings) == 0 {
 		log.Println("warn: no ring found")
-		return nil
+		return nil, removed
 	}
 
 	curr := rings[0]
-	vertices := []int32{curr.a, curr.b}
+
+	polygon := LineStringZ{
+		t.Points[curr.a],
+		t.Points[curr.b],
+	}
 
 	for curr = curr.next; curr != rings[0]; curr = curr.next {
-		vertices = append(vertices, curr.b)
+		polygon = append(polygon, t.Points[curr.b])
 	}
 
-	log.Printf("length of boundary: %d\n", len(vertices))
+	polygon = append(polygon, t.Points[rings[0].a])
 
-	return vertices
+	log.Printf("length of boundary: %d\n", len(polygon))
+
+	return polygon, removed
 }
 
 func (t *Triangulation) TriangleSlices() [][]int32 {