Mercurial > gemma

changeset 4654:3eda5a7215ab stree-experiment

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Generate STRTrees instaead of Octree during sounding result imports.
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Mon, 14 Oct 2019 13:12:38 +0200
parents 8efc6b3289f3
children a2f16987911b
files pkg/imports/sr.go pkg/octree/areas.go pkg/octree/strtree.go pkg/octree/tree.go
diffstat 4 files changed, 102 insertions(+), 105 deletions(-) [+]
line wrap: on
line diff
--- a/pkg/imports/sr.go	Mon Oct 14 12:26:00 2019 +0200
+++ b/pkg/imports/sr.go	Mon Oct 14 13:12:38 2019 +0200
@@ -24,7 +24,6 @@
 	"errors"
 	"fmt"
 	"io"
-	"log"
 	"math"
 	"os"
 	"path"
@@ -156,7 +155,7 @@
 
 	insertOctreeSQL = `
 UPDATE waterway.sounding_results SET
-  octree_checksum = $2, octree_index = $3
+  mesh_checksum = $2, mesh_index = $3
 WHERE id = $1`
 
 	repairBoundarySQL = `
@@ -498,12 +497,9 @@
 
 		start = time.Now()
 
-		var tree *octree.Tree
-		{
-			builder := octree.NewBuilder(tri.Tin())
-			builder.Build(removed)
-			tree = builder.Tree()
-		}
+		tin := tri.Tin()
+		var virtual octree.STRTree
+		virtual.BuildWithout(tin, removed)
 
 		feedback.Info("Building took %v", time.Since(start))
 
@@ -518,9 +514,13 @@
 
 		generated := make(octree.LineStringZ, 0, numPoints+clippingPolygon.NumVertices(0))
 
-		tree.GenerateRandomVertices(numPoints, func(vertices []octree.Vertex) {
-			generated = append(generated, vertices...)
-		})
+		octree.GenerateRandomVertices(
+			numPoints,
+			tin.Min, tin.Max,
+			virtual.Value,
+			func(vertices []octree.Vertex) {
+				generated = append(generated, vertices...)
+			})
 
 		feedback.Info("Generating %d points took %v.", len(generated), time.Since(start))
 
@@ -532,7 +532,7 @@
 				return
 			}
 			dupes[key] = struct{}{}
-			if z, ok := tree.Value(x, y); ok {
+			if z, ok := virtual.Value(x, y); ok {
 				generated = append(generated, octree.Vertex{X: x, Y: y, Z: z})
 			}
 		})
@@ -548,7 +548,6 @@
 		feedback.Info("Second triangulation took %v.", time.Since(start))
 		feedback.Info("Number triangles: %d.", len(tri.Triangles)/3)
 		feedback.Info("Clipping triangles from new mesh.")
-
 	}
 
 	start = time.Now()
@@ -567,23 +566,9 @@
 	feedback.Info("Number of triangles to clip %d.", len(removed))
 
 	start = time.Now()
-	s := octree.STRTree{Entries: 16}
-	s.BuildWithout(tin, removed)
-	if _, err2 := s.Bytes(); err2 != nil {
-		log.Printf("serializing STRTree failed: %v\n", err2)
-	}
-	log.Printf("Building strtree took: %v.\n", time.Since(start))
-
-	// return nil, UnchangedError("nothing to do")
+	final := octree.STRTree{Entries: 16}
+	final.BuildWithout(tin, removed)
 
-	feedback.Info("Build final octree index")
-
-	builder := octree.NewBuilder(tin)
-	builder.Build(removed)
-	octreeIndex, err := builder.Bytes()
-	if err != nil {
-		return nil, err
-	}
 	feedback.Info("Building octree took %v.", time.Since(start))
 	feedback.Info("Store octree.")
 
@@ -623,15 +608,20 @@
 		return nil, err
 	}
 
+	index, err := final.Bytes()
+	if err != nil {
+		return nil, err
+	}
+
 	h := sha1.New()
-	h.Write(octreeIndex)
+	h.Write(index)
 	checksum := hex.EncodeToString(h.Sum(nil))
-	_, err = tx.ExecContext(ctx, insertOctreeSQL, id, checksum, octreeIndex)
+	_, err = tx.ExecContext(ctx, insertOctreeSQL, id, checksum, index)
 	if err != nil {
 		return nil, err
 	}
 	feedback.Info("Storing octree index took %s.", time.Since(start))
-	err = generateIsos(ctx, tx, feedback, builder.Tree(), id)
+	err = generateIsos(ctx, tx, feedback, &final, id)
 	if err != nil {
 		return nil, err
 	}
@@ -838,7 +828,7 @@
 	ctx context.Context,
 	tx *sql.Tx,
 	feedback Feedback,
-	tree *octree.Tree,
+	tree *octree.STRTree,
 	id int64,
 ) error {
 
@@ -847,16 +837,18 @@
 		"morphology_classbreaks",
 	)
 
+	minZ, maxZ := tree.Min().Z, tree.Max().Z
+
 	if err != nil {
 		feedback.Warn("Loading class breaks failed: %v", err)
 		feedback.Info("Using default class breaks")
 		heights = nil
-		h := contourStepWidth * math.Ceil(tree.Min.Z/contourStepWidth)
-		for ; h <= tree.Max.Z; h += contourStepWidth {
+		h := contourStepWidth * math.Ceil(minZ/contourStepWidth)
+		for ; h <= maxZ; h += contourStepWidth {
 			heights = append(heights, h)
 		}
 	} else {
-		heights = octree.ExtrapolateClassBreaks(heights, tree.Min.Z, tree.Max.Z)
+		heights = octree.ExtrapolateClassBreaks(heights, minZ, maxZ)
 		// We set steps for InBetweenClassBreaks to 1, so it
 		// becomes a null operation.  The extra class breaks
 		// were considered unexpected and confusing by the
@@ -881,7 +873,7 @@
 	ctx context.Context,
 	tx *sql.Tx,
 	feedback Feedback,
-	tree *octree.Tree,
+	tree *octree.STRTree,
 	heights []float64,
 	id int64,
 ) error {
@@ -892,11 +884,11 @@
 			time.Since(total))
 	}()
 
-	areas := octree.TraceAreas(heights, isoCellSize, tree.Min, tree.Max, tree.Value)
+	areas := octree.TraceAreas(heights, isoCellSize, tree.Min(), tree.Max(), tree.Value)
 
 	return storeAreas(
 		ctx, tx, feedback,
-		areas, tree.EPSG, heights, id)
+		areas, tree.EPSG(), heights, id)
 }
 
 func storeAreas(
--- a/pkg/octree/areas.go	Mon Oct 14 12:26:00 2019 +0200
+++ b/pkg/octree/areas.go	Mon Oct 14 13:12:38 2019 +0200
@@ -26,6 +26,73 @@
 	"gemma.intevation.de/gemma/pkg/wkb"
 )
 
+func GenerateRandomVertices(
+	n int,
+	min, max Vertex,
+	eval func(float64, float64) (float64, bool),
+	callback func([]Vertex),
+) {
+	var wg sync.WaitGroup
+
+	jobs := make(chan int)
+	out := make(chan []Vertex)
+	done := make(chan struct{})
+
+	cpus := runtime.NumCPU()
+
+	free := make(chan []Vertex, cpus)
+
+	for i := 0; i < cpus; i++ {
+		wg.Add(1)
+		go func() {
+			defer wg.Done()
+			xRange := common.Random(min.X, max.X)
+			yRange := common.Random(min.Y, max.Y)
+
+			for size := range jobs {
+				var vertices []Vertex
+				select {
+				case vertices = <-free:
+				default:
+					vertices = make([]Vertex, 0, 1000)
+				}
+				for len(vertices) < size {
+					x, y := xRange(), yRange()
+					if z, ok := eval(x, y); ok {
+						vertices = append(vertices, Vertex{X: x, Y: y, Z: z})
+					}
+				}
+				out <- vertices
+			}
+		}()
+	}
+
+	go func() {
+		defer close(done)
+		for vertices := range out {
+			callback(vertices)
+			select {
+			case free <- vertices[:0]:
+			default:
+			}
+		}
+	}()
+
+	for remain := n; remain > 0; {
+		if remain > 1000 {
+			jobs <- 1000
+			remain -= 1000
+		} else {
+			jobs <- remain
+			remain = 0
+		}
+	}
+	close(jobs)
+	wg.Wait()
+	close(out)
+	<-done
+}
+
 func TraceAreas(
 	heights []float64,
 	cellSize float64,
--- a/pkg/octree/strtree.go	Mon Oct 14 12:26:00 2019 +0200
+++ b/pkg/octree/strtree.go	Mon Oct 14 13:12:38 2019 +0200
@@ -32,6 +32,10 @@
 	bboxes  []Box2D
 }
 
+func (s *STRTree) Min() Vertex  { return s.tin.Min }
+func (s *STRTree) Max() Vertex  { return s.tin.Max }
+func (s *STRTree) EPSG() uint32 { return s.tin.EPSG }
+
 func (s *STRTree) Value(x, y float64) (float64, bool) {
 
 	stack := []int32{s.index[0]}
--- a/pkg/octree/tree.go	Mon Oct 14 12:26:00 2019 +0200
+++ b/pkg/octree/tree.go	Mon Oct 14 13:12:38 2019 +0200
@@ -15,10 +15,6 @@
 
 import (
 	"math"
-	"runtime"
-	"sync"
-
-	"gemma.intevation.de/gemma/pkg/common"
 )
 
 // Tree is an Octree holding triangles.
@@ -364,65 +360,3 @@
 
 	return result
 }
-
-func (ot *Tree) GenerateRandomVertices(n int, callback func([]Vertex)) {
-	var wg sync.WaitGroup
-
-	jobs := make(chan int)
-	out := make(chan []Vertex)
-	done := make(chan struct{})
-
-	cpus := runtime.NumCPU()
-
-	free := make(chan []Vertex, cpus)
-
-	for i := 0; i < cpus; i++ {
-		wg.Add(1)
-		go func() {
-			defer wg.Done()
-			xRange := common.Random(ot.Min.X, ot.Max.X)
-			yRange := common.Random(ot.Min.Y, ot.Max.Y)
-
-			for size := range jobs {
-				var vertices []Vertex
-				select {
-				case vertices = <-free:
-				default:
-					vertices = make([]Vertex, 0, 1000)
-				}
-				for len(vertices) < size {
-					x, y := xRange(), yRange()
-					if z, ok := ot.Value(x, y); ok {
-						vertices = append(vertices, Vertex{X: x, Y: y, Z: z})
-					}
-				}
-				out <- vertices
-			}
-		}()
-	}
-
-	go func() {
-		defer close(done)
-		for vertices := range out {
-			callback(vertices)
-			select {
-			case free <- vertices[:0]:
-			default:
-			}
-		}
-	}()
-
-	for remain := n; remain > 0; {
-		if remain > 1000 {
-			jobs <- 1000
-			remain -= 1000
-		} else {
-			jobs <- remain
-			remain = 0
-		}
-	}
-	close(jobs)
-	wg.Wait()
-	close(out)
-	<-done
-}