Mercurial > gemma

comparison cmd/octree2contour/loader.go @ 678:7bb961d750b6 octree

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
octree: traverse horizontally over tree to find out which triangles to process.
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Wed, 19 Sep 2018 17:34:19 +0200
parents 120a82bd9953
children c79c7be29a7a
comparison
equal deleted inserted replaced
675:1cb565d244cf 678:7bb961d750b6
50 50
51 if err := binary.Read(r, binary.LittleEndian, &tree.epsg); err != nil { 51 if err := binary.Read(r, binary.LittleEndian, &tree.epsg); err != nil {
52 return nil, err 52 return nil, err
53 } 53 }
54 54
55 log.Printf("EPSG: %d\n", tree.epsg)
56
55 if err := tree.min.read(r); err != nil { 57 if err := tree.min.read(r); err != nil {
56 return nil, err 58 return nil, err
57 } 59 }
58 60
59 if err := tree.max.read(r); err != nil { 61 if err := tree.max.read(r); err != nil {
60 return nil, err 62 return nil, err
61 } 63 }
64
65 log.Printf("BBOX: [[%f, %f, %f], [%f, %f, %f]]\n",
66 tree.min.x, tree.min.y, tree.min.z,
67 tree.max.x, tree.max.y, tree.max.z)
62 68
63 var numVertices uint32 69 var numVertices uint32
64 if err := binary.Read(r, binary.LittleEndian, &numVertices); err != nil { 70 if err := binary.Read(r, binary.LittleEndian, &numVertices); err != nil {
65 return nil, err 71 return nil, err
66 } 72 }
126 } 132 }
127 133
128 return tree, nil 134 return tree, nil
129 } 135 }
130 136
137 func (ot *octree) horizontal(h float64, fn func([]int32)) {
138
139 type frame struct {
140 pos int32
141 min float64
142 max float64
143 }
144
145 if h < ot.min.z || h > ot.max.z {
146 return
147 }
148
149 stack := []frame{{1, ot.min.z, ot.max.z}}
150
151 for len(stack) > 0 {
152 top := stack[len(stack)-1]
153 stack = stack[:len(stack)-1]
154
155 pos, min, max := top.pos, top.min, top.max
156 if pos == 0 {
157 continue
158 }
159
160 if pos > 0 { // node
161 if mid := (max-min)*0.5 + min; h <= mid {
162 stack = append(stack,
163 frame{ot.index[pos+0], min, mid},
164 frame{ot.index[pos+1], min, mid},
165 frame{ot.index[pos+4], min, mid},
166 frame{ot.index[pos+5], min, mid})
167 } else {
168 stack = append(stack,
169 frame{ot.index[pos+2], mid, max},
170 frame{ot.index[pos+3], mid, max},
171 frame{ot.index[pos+6], mid, max},
172 frame{ot.index[pos+7], mid, max})
173 }
174 } else { // leaf
175 pos = -pos - 1
176 n := ot.index[pos]
177 //log.Printf("%d %d %d\n", pos, n, len(ot.index))
178 indices := ot.index[pos+1 : pos+1+n]
179
180 for idx := range indices {
181 tri := ot.triangles[idx]
182 v0 := ot.vertices[tri[0]]
183 v1 := ot.vertices[tri[1]]
184 v2 := ot.vertices[tri[2]]
185
186 if !(math.Min(v0.z, math.Min(v1.z, v2.z)) > h ||
187 math.Max(v0.z, math.Max(v1.z, v2.z)) < h) {
188 fn(tri)
189 }
190 }
191 }
192 }
193 }
194
131 func loadOctree(fname string) (*octree, error) { 195 func loadOctree(fname string) (*octree, error) {
132 196
133 f, err := os.Open(fname) 197 f, err := os.Open(fname)
134 if err != nil { 198 if err != nil {
135 return nil, err 199 return nil, err