Mercurial > gemma
view pkg/octree/strtree.go @ 2515:6bcaa8bf2603 octree-diff
STR tree: Fixed sorting. Stored num items per nodes. Simpified code.
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Tue, 05 Mar 2019 15:55:14 +0100 |
parents | 1534df518201 |
children | 1ec4c5633eb6 |
line wrap: on
line source
// This is Free Software under GNU Affero General Public License v >= 3.0 // without warranty, see README.md and license for details. // // SPDX-License-Identifier: AGPL-3.0-or-later // License-Filename: LICENSES/AGPL-3.0.txt // // Copyright (C) 2018 by via donau // – Österreichische Wasserstraßen-Gesellschaft mbH // Software engineering by Intevation GmbH // // Author(s): // * Sascha L. Teichmann <sascha.teichmann@intevation.de> package octree import ( "math" "sort" ) const numEntries = 8 type STRTree struct { tri *Triangulation index []int32 triangles [][]int32 bboxes []Box2D } func (s *STRTree) Build(t *Triangulation) { s.tri = t s.triangles = t.TriangleSlices() all := make([]int32, len(s.triangles)) for i := range all { all[i] = int32(i) } s.index = append(s.index, 0) root := s.build(all) s.index[0] = root } func (s *STRTree) build(items []int32) int32 { sort.Slice(items, func(i, j int) bool { ti := s.triangles[items[i]] tj := s.triangles[items[j]] return s.tri.Points[ti[0]].X < s.tri.Points[tj[0]].X }) P := int(math.Ceil(float64(len(items)) / float64(numEntries))) S := int(math.Ceil(math.Sqrt(float64(P)))) sm := S * numEntries slices := make([][]int32, S) rest := items for i := range slices { var n int if l := len(rest); l < sm { n = l } else { n = sm } slices[i] = rest[:n] rest = rest[n:] } leaves := make([]int32, 0, S*S) for _, slice := range slices { sort.Slice(slice, func(i, j int) bool { ti := s.triangles[slice[i]] tj := s.triangles[slice[j]] return s.tri.Points[ti[0]].Y < s.tri.Points[tj[0]].Y }) for len(slice) > 0 { n := numEntries if l := len(slice); l < numEntries { n = l } leaves = append(leaves, s.allocLeaf(slice[:n])) slice = slice[n:] } } return s.buildNodes(leaves) } func (s *STRTree) bbox(idx int32) Box2D { if idx < 0 { // Don't care if leaf or node. idx = -idx - 1 } return s.bboxes[s.index[idx]] } func (s *STRTree) buildNodes(items []int32) int32 { if len(items) <= numEntries { return s.allocNode(items) } sort.Slice(items, func(i, j int) bool { return s.bbox(items[i]).X1 < s.bbox(items[j]).X1 }) P := int(math.Ceil(float64(len(items)) / float64(numEntries))) S := int(math.Ceil(math.Sqrt(float64(P)))) sm := S * numEntries slices := make([][]int32, S) rest := items for i := range slices { var n int if l := len(rest); l < sm { n = l } else { n = sm } slices[i] = rest[:n] rest = rest[n:] } nodes := make([]int32, 0, S*S) for _, slice := range slices { sort.Slice(slice, func(i, j int) bool { return s.bbox(slice[i]).Y1 < s.bbox(slice[j]).Y1 }) for len(slice) > 0 { n := numEntries if l := len(slice); l < numEntries { n = l } nodes = append(nodes, s.allocNode(slice[:n])) slice = slice[n:] } } return s.buildNodes(nodes) } func (s *STRTree) allocNode(items []int32) int32 { pos := len(s.index) s.index = append(s.index, 0, int32(len(items))) s.index = append(s.index, items...) if len(items) > 0 { box := s.bbox(items[0]) for i := 1; i < len(items); i++ { box = box.Union(s.bbox(items[i])) } s.index[pos] = int32(s.allocBBox(box)) } return int32(pos) } func (s *STRTree) allocBBox(box Box2D) int { pos := len(s.bboxes) s.bboxes = append(s.bboxes, box) return pos } func (s *STRTree) allocLeaf(items []int32) int32 { pos := len(s.index) s.index = append(s.index, 0, int32(len(items))) s.index = append(s.index, items...) if len(items) > 0 { vertices := s.tri.Points ti := s.triangles[items[0]] t := Triangle{ vertices[ti[0]], vertices[ti[1]], vertices[ti[2]], } box := t.BBox() for i := 1; i < len(items); i++ { it := items[i] ti := s.triangles[it] t := Triangle{ vertices[ti[0]], vertices[ti[1]], vertices[ti[2]], } box = box.Union(t.BBox()) } s.index[pos] = int32(s.allocBBox(box)) } return int32(-(pos + 1)) }