Mercurial > gemma
view cmd/octree2contour/loader.go @ 687:be90ab542aa7 octree
octree: contouring: Do the math to calculate the intersection points of the triangles and the planes.
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Thu, 20 Sep 2018 11:32:03 +0200 |
parents | a8672ba9ebad |
children | 9db4ae29ded9 |
line wrap: on
line source
package main import ( "bufio" "encoding/binary" "io" "log" "math" "os" "github.com/golang/snappy" ) type octree struct { epsg uint32 vertices []vertex triangles [][]int32 index []int32 min vertex max vertex } func (v *vertex) read(r io.Reader) error { var buf [8]byte b := buf[:] if _, err := io.ReadFull(r, b); err != nil { return nil } v.x = math.Float64frombits(binary.LittleEndian.Uint64(b)) if _, err := io.ReadFull(r, b); err != nil { return nil } v.y = math.Float64frombits(binary.LittleEndian.Uint64(b)) if _, err := io.ReadFull(r, b); err != nil { return nil } v.z = math.Float64frombits(binary.LittleEndian.Uint64(b)) return nil } func loadOctreeReader(r *bufio.Reader) (*octree, error) { tree := new(octree) if err := binary.Read(r, binary.LittleEndian, &tree.epsg); err != nil { return nil, err } log.Printf("EPSG: %d\n", tree.epsg) if err := tree.min.read(r); err != nil { return nil, err } if err := tree.max.read(r); err != nil { return nil, err } log.Printf("BBOX: [[%f, %f, %f], [%f, %f, %f]]\n", tree.min.x, tree.min.y, tree.min.z, tree.max.x, tree.max.y, tree.max.z) var numVertices uint32 if err := binary.Read(r, binary.LittleEndian, &numVertices); err != nil { return nil, err } log.Printf("vertices: %d\n", numVertices) vertices := make([]vertex, numVertices) tree.vertices = vertices for i := range vertices { if err := vertices[i].read(r); err != nil { return nil, err } } var numTriangles uint32 if err := binary.Read(r, binary.LittleEndian, &numTriangles); err != nil { return nil, err } log.Printf("triangles: %d\n", numTriangles) indices := make([]int32, 3*numTriangles) triangles := make([][]int32, numTriangles) tree.triangles = triangles var last int32 for i := range triangles { tri := indices[:3] indices = indices[3:] triangles[i] = tri for j := range tri { v, err := binary.ReadVarint(r) if err != nil { return nil, err } value := int32(v) + last tri[j] = value last = value } } var numNodes uint32 if err := binary.Read(r, binary.LittleEndian, &numNodes); err != nil { return nil, err } log.Printf("num nodes: %d\n", numNodes) tree.index = make([]int32, numNodes) entries := tree.index[1:] last = 0 for i := range entries { v, err := binary.ReadVarint(r) if err != nil { return nil, err } value := int32(v) + last entries[i] = value last = value } return tree, nil } func (ot *octree) horizontal(h float64, fn func(*triangle)) { type frame struct { pos int32 min float64 max float64 } if h < ot.min.z || ot.max.z < h { return } stack := []frame{{1, ot.min.z, ot.max.z}} dupes := map[int32]struct{}{} for len(stack) > 0 { top := stack[len(stack)-1] stack = stack[:len(stack)-1] pos := top.pos if pos == 0 { continue } min, max := top.min, top.max if pos > 0 { // node if mid := (max-min)*0.5 + min; h >= mid { pos += 4 // nodes with z-bit set min = mid } else { max = mid } stack = append(stack, frame{ot.index[pos+0], min, max}, frame{ot.index[pos+1], min, max}, frame{ot.index[pos+2], min, max}, frame{ot.index[pos+3], min, max}) } else { // leaf pos = -pos - 1 n := ot.index[pos] //log.Printf("%d %d %d\n", pos, n, len(ot.index)) indices := ot.index[pos+1 : pos+1+n] for _, idx := range indices { if _, found := dupes[idx]; found { continue } tri := ot.triangles[idx] t := triangle{ ot.vertices[tri[0]], ot.vertices[tri[1]], ot.vertices[tri[2]], } if !(math.Min(t[0].z, math.Min(t[1].z, t[2].z)) > h || math.Max(t[0].z, math.Max(t[1].z, t[2].z)) < h) { dupes[idx] = struct{}{} fn(&t) } } } } } func loadOctree(fname string) (*octree, error) { f, err := os.Open(fname) if err != nil { return nil, err } defer f.Close() return loadOctreeReader(bufio.NewReader(snappy.NewReader(f))) }