GCC Code Coverage Report

Directory: ./
File: openvdb/openvdb/tools/VolumeToSpheres.h
Date: 2022-07-25 17:40:05
Exec Total Coverage
Lines: 260 268 97.0%
Functions: 17 26 65.4%
Branches: 227 354 64.1%
Line Branch Exec Source
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3
4 /// @file tools/VolumeToSpheres.h
5 ///
6 /// @brief Fill a closed level set or fog volume with adaptively-sized spheres.
7
8 #ifndef OPENVDB_TOOLS_VOLUME_TO_SPHERES_HAS_BEEN_INCLUDED
9 #define OPENVDB_TOOLS_VOLUME_TO_SPHERES_HAS_BEEN_INCLUDED
10
11 #include <openvdb/tree/LeafManager.h>
12 #include <openvdb/math/Math.h>
13 #include "Morphology.h" // for erodeActiveValues()
14 #include "PointScatter.h"
15 #include "LevelSetRebuild.h"
16 #include "LevelSetUtil.h"
17 #include "VolumeToMesh.h"
18 #include <openvdb/openvdb.h>
19
20 #include <tbb/blocked_range.h>
21 #include <tbb/parallel_for.h>
22 #include <tbb/parallel_reduce.h>
23
24 #include <algorithm> // for std::min(), std::max()
25 #include <cmath> // for std::sqrt()
26 #include <limits> // for std::numeric_limits
27 #include <memory>
28 #include <random>
29 #include <utility> // for std::pair
30 #include <vector>
31
32
33 namespace openvdb {
34 OPENVDB_USE_VERSION_NAMESPACE
35 namespace OPENVDB_VERSION_NAME {
36 namespace tools {
37
38 /// @brief Fill a closed level set or fog volume with adaptively-sized spheres.
39 ///
40 /// @param grid a scalar grid that defines the surface to be filled with spheres
41 /// @param spheres an output array of 4-tuples representing the fitted spheres<BR>
42 /// The first three components of each tuple specify the sphere center,
43 /// and the fourth specifies the radius.
44 /// The spheres are ordered by radius, from largest to smallest.
45 /// @param sphereCount lower and upper bounds on the number of spheres to be generated<BR>
46 /// The actual number will be somewhere within the bounds.
47 /// @param overlapping toggle to allow spheres to overlap/intersect
48 /// @param minRadius the smallest allowable sphere size, in voxel units<BR>
49 /// @param maxRadius the largest allowable sphere size, in voxel units
50 /// @param isovalue the voxel value that determines the surface of the volume<BR>
51 /// The default value of zero works for signed distance fields,
52 /// while fog volumes require a larger positive value
53 /// (0.5 is a good initial guess).
54 /// @param instanceCount the number of interior points to consider for the sphere placement<BR>
55 /// Increasing this count increases the chances of finding optimal
56 /// sphere sizes.
57 /// @param interrupter pointer to an object adhering to the util::NullInterrupter interface
58 ///
59 /// @note The minimum sphere count takes precedence over the minimum radius.
60 template<typename GridT, typename InterrupterT = util::NullInterrupter>
61 void
62 fillWithSpheres(
63 const GridT& grid,
64 std::vector<openvdb::Vec4s>& spheres,
65 const Vec2i& sphereCount = Vec2i(1, 50),
66 bool overlapping = false,
67 float minRadius = 1.0,
68 float maxRadius = std::numeric_limits<float>::max(),
69 float isovalue = 0.0,
70 int instanceCount = 10000,
71 InterrupterT* interrupter = nullptr);
72
73
74 ////////////////////////////////////////
75
76
77 /// @brief Accelerated closest surface point queries for narrow band level sets
78 /// @details Supports queries that originate at arbitrary world-space locations,
79 /// is not confined to the narrow band region of the input volume geometry.
80 template<typename GridT>
81 class ClosestSurfacePoint
82 {
83 public:
84 using Ptr = std::unique_ptr<ClosestSurfacePoint>;
85 using TreeT = typename GridT::TreeType;
86 using BoolTreeT = typename TreeT::template ValueConverter<bool>::Type;
87 using Index32TreeT = typename TreeT::template ValueConverter<Index32>::Type;
88 using Int16TreeT = typename TreeT::template ValueConverter<Int16>::Type;
89
90 /// @brief Extract surface points and construct a spatial acceleration structure.
91 ///
92 /// @return a null pointer if the initialization fails for any reason,
93 /// otherwise a unique pointer to a newly-allocated ClosestSurfacePoint object.
94 ///
95 /// @param grid a scalar level set or fog volume
96 /// @param isovalue the voxel value that determines the surface of the volume
97 /// The default value of zero works for signed distance fields,
98 /// while fog volumes require a larger positive value
99 /// (0.5 is a good initial guess).
100 /// @param interrupter pointer to an object adhering to the util::NullInterrupter interface.
101 template<typename InterrupterT = util::NullInterrupter>
102 static Ptr create(const GridT& grid, float isovalue = 0.0,
103 InterrupterT* interrupter = nullptr);
104
105 /// @brief Compute the distance from each input point to its closest surface point.
106 /// @param points input list of points in world space
107 /// @param distances output list of closest surface point distances
108 bool search(const std::vector<Vec3R>& points, std::vector<float>& distances);
109
110 /// @brief Overwrite each input point with its closest surface point.
111 /// @param points input/output list of points in world space
112 /// @param distances output list of closest surface point distances
113 bool searchAndReplace(std::vector<Vec3R>& points, std::vector<float>& distances);
114
115 /// @brief Tree accessor
116 ✗ const Index32TreeT& indexTree() const { return *mIdxTreePt; }
117 /// @brief Tree accessor
118 ✗ const Int16TreeT& signTree() const { return *mSignTreePt; }
119
120 private:
121 using Index32LeafT = typename Index32TreeT::LeafNodeType;
122 using IndexRange = std::pair<size_t, size_t>;
123
124 std::vector<Vec4R> mLeafBoundingSpheres, mNodeBoundingSpheres;
125 std::vector<IndexRange> mLeafRanges;
126 std::vector<const Index32LeafT*> mLeafNodes;
127 PointList mSurfacePointList;
128 size_t mPointListSize = 0, mMaxNodeLeafs = 0;
129 typename Index32TreeT::Ptr mIdxTreePt;
130 typename Int16TreeT::Ptr mSignTreePt;
131
132
1/4
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 12 ClosestSurfacePoint() = default;
133 template<typename InterrupterT = util::NullInterrupter>
134 bool initialize(const GridT&, float isovalue, InterrupterT*);
135 bool search(std::vector<Vec3R>&, std::vector<float>&, bool transformPoints);
136 };
137
138
139 ////////////////////////////////////////
140
141
142 // Internal utility methods
143
144 /// @cond OPENVDB_DOCS_INTERNAL
145
146 namespace v2s_internal {
147
148 struct PointAccessor
149 {
150 PointAccessor(std::vector<Vec3R>& points)
151 10 : mPoints(points)
152 {
153 }
154
155 void add(const Vec3R &pos)
156 {
157 65000 mPoints.push_back(pos);
158 }
159 private:
160 std::vector<Vec3R>& mPoints;
161 };
162
163
164 template<typename Index32LeafT>
165 class LeafOp
166 {
167 public:
168
169 LeafOp(std::vector<Vec4R>& leafBoundingSpheres,
170 const std::vector<const Index32LeafT*>& leafNodes,
171 const math::Transform& transform,
172 const PointList& surfacePointList);
173
174 void run(bool threaded = true);
175
176
177 void operator()(const tbb::blocked_range<size_t>&) const;
178
179 private:
180 std::vector<Vec4R>& mLeafBoundingSpheres;
181 const std::vector<const Index32LeafT*>& mLeafNodes;
182 const math::Transform& mTransform;
183 const PointList& mSurfacePointList;
184 };
185
186 template<typename Index32LeafT>
187 12 LeafOp<Index32LeafT>::LeafOp(
188 std::vector<Vec4R>& leafBoundingSpheres,
189 const std::vector<const Index32LeafT*>& leafNodes,
190 const math::Transform& transform,
191 const PointList& surfacePointList)
192 : mLeafBoundingSpheres(leafBoundingSpheres)
193 , mLeafNodes(leafNodes)
194 , mTransform(transform)
195 12 , mSurfacePointList(surfacePointList)
196 {
197 12 }
198
199 template<typename Index32LeafT>
200 void
201 12 LeafOp<Index32LeafT>::run(bool threaded)
202 {
203
1/2
✓ Branch 0 taken 12 times.
✗ Branch 1 not taken.
 12 if (threaded) {
204 12 tbb::parallel_for(tbb::blocked_range<size_t>(0, mLeafNodes.size()), *this);
205 } else {
206 ✗ (*this)(tbb::blocked_range<size_t>(0, mLeafNodes.size()));
207 }
208 12 }
209
210 template<typename Index32LeafT>
211 void
212 776 LeafOp<Index32LeafT>::operator()(const tbb::blocked_range<size_t>& range) const
213 {
214 typename Index32LeafT::ValueOnCIter iter;
215 Vec3s avg;
216
217
2/2
✓ Branch 0 taken 1354 times.
✓ Branch 1 taken 776 times.
 2130 for (size_t n = range.begin(); n != range.end(); ++n) {
218 1354 avg[0] = 0.0;
219 1354 avg[1] = 0.0;
220 1354 avg[2] = 0.0;
221
222 int count = 0;
223
2/2
✓ Branch 1 taken 86624 times.
✓ Branch 2 taken 1354 times.
 89332 for (iter = mLeafNodes[n]->cbeginValueOn(); iter; ++iter) {
224 86624 avg += mSurfacePointList[iter.getValue()];
225 86624 ++count;
226 }
227
2/2
✓ Branch 0 taken 1346 times.
✓ Branch 1 taken 8 times.
 1354 if (count > 1) avg *= float(1.0 / double(count));
228
229 float maxDist = 0.0;
230
2/2
✓ Branch 1 taken 86624 times.
✓ Branch 2 taken 1354 times.
 89332 for (iter = mLeafNodes[n]->cbeginValueOn(); iter; ++iter) {
231
2/2
✓ Branch 1 taken 3726 times.
✓ Branch 2 taken 82898 times.
 86624 float tmpDist = (mSurfacePointList[iter.getValue()] - avg).lengthSqr();
232
2/2
✓ Branch 0 taken 3726 times.
✓ Branch 1 taken 82898 times.
 86624 if (tmpDist > maxDist) maxDist = tmpDist;
233 }
234
235 1354 Vec4R& sphere = mLeafBoundingSpheres[n];
236 1354 sphere[0] = avg[0];
237 1354 sphere[1] = avg[1];
238 1354 sphere[2] = avg[2];
239 1354 sphere[3] = maxDist * 2.0; // padded radius
240 }
241 776 }
242
243
244 class NodeOp
245 {
246 public:
247 using IndexRange = std::pair<size_t, size_t>;
248
249 NodeOp(std::vector<Vec4R>& nodeBoundingSpheres,
250 const std::vector<IndexRange>& leafRanges,
251 const std::vector<Vec4R>& leafBoundingSpheres);
252
253 inline void run(bool threaded = true);
254
255 inline void operator()(const tbb::blocked_range<size_t>&) const;
256
257 private:
258 std::vector<Vec4R>& mNodeBoundingSpheres;
259 const std::vector<IndexRange>& mLeafRanges;
260 const std::vector<Vec4R>& mLeafBoundingSpheres;
261 };
262
263 inline
264 NodeOp::NodeOp(std::vector<Vec4R>& nodeBoundingSpheres,
265 const std::vector<IndexRange>& leafRanges,
266 12 const std::vector<Vec4R>& leafBoundingSpheres)
267 : mNodeBoundingSpheres(nodeBoundingSpheres)
268 , mLeafRanges(leafRanges)
269 12 , mLeafBoundingSpheres(leafBoundingSpheres)
270 {
271 }
272
273 inline void
274 12 NodeOp::run(bool threaded)
275 {
276
1/2
✓ Branch 0 taken 12 times.
✗ Branch 1 not taken.
 12 if (threaded) {
277 12 tbb::parallel_for(tbb::blocked_range<size_t>(0, mLeafRanges.size()), *this);
278 } else {
279 ✗ (*this)(tbb::blocked_range<size_t>(0, mLeafRanges.size()));
280 }
281 12 }
282
283 inline void
284 41 NodeOp::operator()(const tbb::blocked_range<size_t>& range) const
285 {
286 Vec3d avg, pos;
287
288
2/2
✓ Branch 0 taken 41 times.
✓ Branch 1 taken 41 times.
 82 for (size_t n = range.begin(); n != range.end(); ++n) {
289
290 41 avg[0] = 0.0;
291 41 avg[1] = 0.0;
292 41 avg[2] = 0.0;
293
294 41 int count = int(mLeafRanges[n].second) - int(mLeafRanges[n].first);
295
296
2/2
✓ Branch 0 taken 1354 times.
✓ Branch 1 taken 41 times.
 1395 for (size_t i = mLeafRanges[n].first; i < mLeafRanges[n].second; ++i) {
297 1354 avg[0] += mLeafBoundingSpheres[i][0];
298 1354 avg[1] += mLeafBoundingSpheres[i][1];
299 1354 avg[2] += mLeafBoundingSpheres[i][2];
300 }
301
302
2/2
✓ Branch 0 taken 33 times.
✓ Branch 1 taken 8 times.
 41 if (count > 1) avg *= float(1.0 / double(count));
303
304
305 double maxDist = 0.0;
306
307
2/2
✓ Branch 0 taken 1354 times.
✓ Branch 1 taken 41 times.
 1395 for (size_t i = mLeafRanges[n].first; i < mLeafRanges[n].second; ++i) {
308
2/2
✓ Branch 0 taken 87 times.
✓ Branch 1 taken 1267 times.
 1354 pos[0] = mLeafBoundingSpheres[i][0];
309 pos[1] = mLeafBoundingSpheres[i][1];
310 pos[2] = mLeafBoundingSpheres[i][2];
311 const auto radiusSqr = mLeafBoundingSpheres[i][3];
312
313 1354 double tmpDist = (pos - avg).lengthSqr() + radiusSqr;
314
2/2
✓ Branch 0 taken 87 times.
✓ Branch 1 taken 1267 times.
 1354 if (tmpDist > maxDist) maxDist = tmpDist;
315 }
316
317 41 Vec4R& sphere = mNodeBoundingSpheres[n];
318
319 41 sphere[0] = avg[0];
320 41 sphere[1] = avg[1];
321 41 sphere[2] = avg[2];
322 41 sphere[3] = maxDist * 2.0; // padded radius
323 }
324 41 }
325
326
327 ////////////////////////////////////////
328
329
330 template<typename Index32LeafT>
331 class ClosestPointDist
332 {
333 public:
334 using IndexRange = std::pair<size_t, size_t>;
335
336 ClosestPointDist(
337 std::vector<Vec3R>& instancePoints,
338 std::vector<float>& instanceDistances,
339 const PointList& surfacePointList,
340 const std::vector<const Index32LeafT*>& leafNodes,
341 const std::vector<IndexRange>& leafRanges,
342 const std::vector<Vec4R>& leafBoundingSpheres,
343 const std::vector<Vec4R>& nodeBoundingSpheres,
344 size_t maxNodeLeafs,
345 bool transformPoints = false);
346
347
348 void run(bool threaded = true);
349
350
351 void operator()(const tbb::blocked_range<size_t>&) const;
352
353 private:
354
355 void evalLeaf(size_t index, const Index32LeafT& leaf) const;
356 void evalNode(size_t pointIndex, size_t nodeIndex) const;
357
358
359 std::vector<Vec3R>& mInstancePoints;
360 std::vector<float>& mInstanceDistances;
361
362 const PointList& mSurfacePointList;
363
364 const std::vector<const Index32LeafT*>& mLeafNodes;
365 const std::vector<IndexRange>& mLeafRanges;
366 const std::vector<Vec4R>& mLeafBoundingSpheres;
367 const std::vector<Vec4R>& mNodeBoundingSpheres;
368
369 std::vector<float> mLeafDistances, mNodeDistances;
370
371 const bool mTransformPoints;
372 size_t mClosestPointIndex;
373 };// ClosestPointDist
374
375
376 template<typename Index32LeafT>
377 14 ClosestPointDist<Index32LeafT>::ClosestPointDist(
378 std::vector<Vec3R>& instancePoints,
379 std::vector<float>& instanceDistances,
380 const PointList& surfacePointList,
381 const std::vector<const Index32LeafT*>& leafNodes,
382 const std::vector<IndexRange>& leafRanges,
383 const std::vector<Vec4R>& leafBoundingSpheres,
384 const std::vector<Vec4R>& nodeBoundingSpheres,
385 size_t maxNodeLeafs,
386 bool transformPoints)
387 : mInstancePoints(instancePoints)
388 , mInstanceDistances(instanceDistances)
389 , mSurfacePointList(surfacePointList)
390 , mLeafNodes(leafNodes)
391 , mLeafRanges(leafRanges)
392 , mLeafBoundingSpheres(leafBoundingSpheres)
393 , mNodeBoundingSpheres(nodeBoundingSpheres)
394 28 , mLeafDistances(maxNodeLeafs, 0.0)
395 28 , mNodeDistances(leafRanges.size(), 0.0)
396 , mTransformPoints(transformPoints)
397
1/4
✓ Branch 2 taken 14 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 14 , mClosestPointIndex(0)
398 {
399 14 }
400
401
402 template<typename Index32LeafT>
403 void
404 14 ClosestPointDist<Index32LeafT>::run(bool threaded)
405 {
406
1/2
✓ Branch 0 taken 14 times.
✗ Branch 1 not taken.
 14 if (threaded) {
407 14 tbb::parallel_for(tbb::blocked_range<size_t>(0, mInstancePoints.size()), *this);
408 } else {
409 ✗ (*this)(tbb::blocked_range<size_t>(0, mInstancePoints.size()));
410 }
411 14 }
412
413 template<typename Index32LeafT>
414 void
415 365473 ClosestPointDist<Index32LeafT>::evalLeaf(size_t index, const Index32LeafT& leaf) const
416 {
417 typename Index32LeafT::ValueOnCIter iter;
418 365473 const Vec3s center = mInstancePoints[index];
419 size_t& closestPointIndex = const_cast<size_t&>(mClosestPointIndex);
420
421
2/2
✓ Branch 0 taken 29458060 times.
✓ Branch 1 taken 365473 times.
 30189006 for (iter = leaf.cbeginValueOn(); iter; ++iter) {
422
423
2/2
✓ Branch 1 taken 1491138 times.
✓ Branch 2 taken 27966922 times.
 29458060 const Vec3s& point = mSurfacePointList[iter.getValue()];
424 float tmpDist = (point - center).lengthSqr();
425
426
2/2
✓ Branch 0 taken 1491138 times.
✓ Branch 1 taken 27966922 times.
 29458060 if (tmpDist < mInstanceDistances[index]) {
427 1491138 mInstanceDistances[index] = tmpDist;
428 1491138 closestPointIndex = iter.getValue();
429 }
430 }
431 365473 }
432
433
434 template<typename Index32LeafT>
435 void
436 129462 ClosestPointDist<Index32LeafT>::evalNode(size_t pointIndex, size_t nodeIndex) const
437 {
438
2/2
✓ Branch 0 taken 124462 times.
✓ Branch 1 taken 5000 times.
 129462 if (nodeIndex >= mLeafRanges.size()) return;
439
440 124462 const Vec3R& pos = mInstancePoints[pointIndex];
441 124462 float minDist = mInstanceDistances[pointIndex];
442 size_t minDistIdx = 0;
443 Vec3R center;
444 bool updatedDist = false;
445
446 4681982 for (size_t i = mLeafRanges[nodeIndex].first, n = 0;
447
2/2
✓ Branch 0 taken 4557520 times.
✓ Branch 1 taken 124462 times.
 4681982 i < mLeafRanges[nodeIndex].second; ++i, ++n)
448 {
449 float& distToLeaf = const_cast<float&>(mLeafDistances[n]);
450
451
2/2
✓ Branch 0 taken 4395390 times.
✓ Branch 1 taken 162130 times.
 4557520 center[0] = mLeafBoundingSpheres[i][0];
452 center[1] = mLeafBoundingSpheres[i][1];
453 center[2] = mLeafBoundingSpheres[i][2];
454 const auto radiusSqr = mLeafBoundingSpheres[i][3];
455
456
2/2
✓ Branch 0 taken 4395390 times.
✓ Branch 1 taken 162130 times.
 4557520 distToLeaf = float(std::max(0.0, (pos - center).lengthSqr() - radiusSqr));
457
458
2/2
✓ Branch 0 taken 454835 times.
✓ Branch 1 taken 4102685 times.
 4557520 if (distToLeaf < minDist) {
459 minDist = distToLeaf;
460 minDistIdx = i;
461 updatedDist = true;
462 }
463 }
464
465
2/2
✓ Branch 0 taken 99817 times.
✓ Branch 1 taken 24645 times.
 124462 if (!updatedDist) return;
466
467 99817 evalLeaf(pointIndex, *mLeafNodes[minDistIdx]);
468
469 3866781 for (size_t i = mLeafRanges[nodeIndex].first, n = 0;
470
2/2
✓ Branch 0 taken 3766964 times.
✓ Branch 1 taken 99817 times.
 3866781 i < mLeafRanges[nodeIndex].second; ++i, ++n)
471 {
472
4/4
✓ Branch 0 taken 359240 times.
✓ Branch 1 taken 3407724 times.
✓ Branch 2 taken 265656 times.
✓ Branch 3 taken 93584 times.
 3766964 if (mLeafDistances[n] < mInstanceDistances[pointIndex] && i != minDistIdx) {
473 265656 evalLeaf(pointIndex, *mLeafNodes[i]);
474 }
475 }
476 }
477
478
479 template<typename Index32LeafT>
480 void
481 1369 ClosestPointDist<Index32LeafT>::operator()(const tbb::blocked_range<size_t>& range) const
482 {
483 Vec3R center;
484
2/2
✓ Branch 0 taken 65938 times.
✓ Branch 1 taken 1369 times.
 67307 for (size_t n = range.begin(); n != range.end(); ++n) {
485
486 65938 const Vec3R& pos = mInstancePoints[n];
487 65938 float minDist = mInstanceDistances[n];
488 size_t minDistIdx = 0;
489
490
2/2
✓ Branch 0 taken 231064 times.
✓ Branch 1 taken 65938 times.
 297002 for (size_t i = 0, I = mNodeDistances.size(); i < I; ++i) {
491 float& distToNode = const_cast<float&>(mNodeDistances[i]);
492
493
2/2
✓ Branch 0 taken 118733 times.
✓ Branch 1 taken 112331 times.
 231064 center[0] = mNodeBoundingSpheres[i][0];
494 center[1] = mNodeBoundingSpheres[i][1];
495 center[2] = mNodeBoundingSpheres[i][2];
496 const auto radiusSqr = mNodeBoundingSpheres[i][3];
497
498
2/2
✓ Branch 0 taken 118733 times.
✓ Branch 1 taken 112331 times.
 231064 distToNode = float(std::max(0.0, (pos - center).lengthSqr() - radiusSqr));
499
500
2/2
✓ Branch 0 taken 127797 times.
✓ Branch 1 taken 103267 times.
 231064 if (distToNode < minDist) {
501 minDist = distToNode;
502 minDistIdx = i;
503 }
504 }
505
506 65938 evalNode(n, minDistIdx);
507
508
2/2
✓ Branch 0 taken 231064 times.
✓ Branch 1 taken 65938 times.
 297002 for (size_t i = 0, I = mNodeDistances.size(); i < I; ++i) {
509
4/4
✓ Branch 0 taken 124462 times.
✓ Branch 1 taken 106602 times.
✓ Branch 2 taken 63524 times.
✓ Branch 3 taken 60938 times.
 231064 if (mNodeDistances[i] < mInstanceDistances[n] && i != minDistIdx) {
510 63524 evalNode(n, i);
511 }
512 }
513
514
2/2
✓ Branch 0 taken 18 times.
✓ Branch 1 taken 65920 times.
 65938 mInstanceDistances[n] = std::sqrt(mInstanceDistances[n]);
515
516
2/2
✓ Branch 0 taken 18 times.
✓ Branch 1 taken 65920 times.
 65938 if (mTransformPoints) mInstancePoints[n] = mSurfacePointList[mClosestPointIndex];
517 }
518 1369 }
519
520
521 class UpdatePoints
522 {
523 public:
524 UpdatePoints(
525 const Vec4s& sphere,
526 const std::vector<Vec3R>& points,
527 std::vector<float>& distances,
528 std::vector<unsigned char>& mask,
529 bool overlapping);
530
531 43 float radius() const { return mRadius; }
532 43 int index() const { return mIndex; }
533
534 inline void run(bool threaded = true);
535
536
537 UpdatePoints(UpdatePoints&, tbb::split);
538 inline void operator()(const tbb::blocked_range<size_t>& range);
539 void join(const UpdatePoints& rhs)
540 {
541
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 161 times.
 162 if (rhs.mRadius > mRadius) {
542 1 mRadius = rhs.mRadius;
543 1 mIndex = rhs.mIndex;
544 }
545 }
546
547 private:
548 const Vec4s& mSphere;
549 const std::vector<Vec3R>& mPoints;
550 std::vector<float>& mDistances;
551 std::vector<unsigned char>& mMask;
552 bool mOverlapping;
553 float mRadius;
554 int mIndex;
555 };
556
557 inline
558 UpdatePoints::UpdatePoints(
559 const Vec4s& sphere,
560 const std::vector<Vec3R>& points,
561 std::vector<float>& distances,
562 std::vector<unsigned char>& mask,
563 43 bool overlapping)
564 : mSphere(sphere)
565 , mPoints(points)
566 , mDistances(distances)
567 , mMask(mask)
568 , mOverlapping(overlapping)
569 , mRadius(0.0)
570 43 , mIndex(0)
571 {
572 }
573
574 inline
575 162 UpdatePoints::UpdatePoints(UpdatePoints& rhs, tbb::split)
576 162 : mSphere(rhs.mSphere)
577 162 , mPoints(rhs.mPoints)
578 162 , mDistances(rhs.mDistances)
579 162 , mMask(rhs.mMask)
580 162 , mOverlapping(rhs.mOverlapping)
581 162 , mRadius(rhs.mRadius)
582 162 , mIndex(rhs.mIndex)
583 {
584 }
585
586 inline void
587 43 UpdatePoints::run(bool threaded)
588 {
589
1/2
✓ Branch 0 taken 43 times.
✗ Branch 1 not taken.
 43 if (threaded) {
590 43 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, mPoints.size()), *this);
591 } else {
592 ✗ (*this)(tbb::blocked_range<size_t>(0, mPoints.size()));
593 }
594 43 }
595
596 inline void
597 7064 UpdatePoints::operator()(const tbb::blocked_range<size_t>& range)
598 {
599 Vec3s pos;
600
2/2
✓ Branch 0 taken 242820 times.
✓ Branch 1 taken 7064 times.
 249884 for (size_t n = range.begin(); n != range.end(); ++n) {
601
2/2
✓ Branch 0 taken 172342 times.
✓ Branch 1 taken 70478 times.
 242820 if (mMask[n]) continue;
602
603
2/2
✓ Branch 0 taken 59569 times.
✓ Branch 1 taken 10909 times.
 70478 pos.x() = float(mPoints[n].x()) - mSphere[0];
604
2/2
✓ Branch 0 taken 59569 times.
✓ Branch 1 taken 10909 times.
 70478 pos.y() = float(mPoints[n].y()) - mSphere[1];
605
2/2
✓ Branch 0 taken 59569 times.
✓ Branch 1 taken 10909 times.
 70478 pos.z() = float(mPoints[n].z()) - mSphere[2];
606
607 float dist = pos.length();
608
609
2/2
✓ Branch 0 taken 59569 times.
✓ Branch 1 taken 10909 times.
 70478 if (dist < mSphere[3]) {
610 59569 mMask[n] = 1;
611 59569 continue;
612 }
613
614
2/2
✓ Branch 0 taken 5713 times.
✓ Branch 1 taken 5196 times.
 10909 if (!mOverlapping) {
615
2/2
✓ Branch 0 taken 1570 times.
✓ Branch 1 taken 4143 times.
 7283 mDistances[n] = std::min(mDistances[n], (dist - mSphere[3]));
616 }
617
618
2/2
✓ Branch 0 taken 287 times.
✓ Branch 1 taken 10622 times.
 10909 if (mDistances[n] > mRadius) {
619 287 mRadius = mDistances[n];
620 287 mIndex = int(n);
621 }
622 }
623 7064 }
624
625
626 } // namespace v2s_internal
627
628 /// @endcond
629
630 ////////////////////////////////////////
631
632
633 template<typename GridT, typename InterrupterT>
634 void
635
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
 20 fillWithSpheres(
636 const GridT& grid,
637 std::vector<openvdb::Vec4s>& spheres,
638 const Vec2i& sphereCount,
639 bool overlapping,
640 float minRadius,
641 float maxRadius,
642 float isovalue,
643 int instanceCount,
644 InterrupterT* interrupter)
645 {
646 spheres.clear();
647
648
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 if (grid.empty()) return;
649
650 const int
651 minSphereCount = sphereCount[0],
652 20 maxSphereCount = sphereCount[1];
653
2/4
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
 20 if ((minSphereCount > maxSphereCount) || (maxSphereCount < 1)) {
654 ✗ OPENVDB_LOG_WARN("fillWithSpheres: minimum sphere count ("
655 << minSphereCount << ") exceeds maximum count (" << maxSphereCount << ")");
656 ✗ return;
657 }
658 20 spheres.reserve(maxSphereCount);
659
660 20 auto gridPtr = grid.copy(); // shallow copy
661
662
3/4
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 7 times.
✓ Branch 4 taken 3 times.
 20 if (gridPtr->getGridClass() == GRID_LEVEL_SET) {
663 // Clamp the isovalue to the level set's background value minus epsilon.
664 // (In a valid narrow-band level set, all voxels, including background voxels,
665 // have values less than or equal to the background value, so an isovalue
666 // greater than or equal to the background value would produce a mask with
667 // effectively infinite extent.)
668 14 isovalue = std::min(isovalue,
669
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 6 times.
 16 static_cast<float>(gridPtr->background() - math::Tolerance<float>::value()));
670
3/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 1 times.
 6 } else if (gridPtr->getGridClass() == GRID_FOG_VOLUME) {
671 // Clamp the isovalue of a fog volume between epsilon and one,
672 // again to avoid a mask with infinite extent. (Recall that
673 // fog volume voxel values vary from zero outside to one inside.)
674
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 8 isovalue = math::Clamp(isovalue, math::Tolerance<float>::value(), 1.f);
675 }
676
677 // ClosestSurfacePoint is inaccurate for small grids.
678 // Resample the input grid if it is too small.
679
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 auto numVoxels = gridPtr->activeVoxelCount();
680
2/2
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
 20 if (numVoxels < 10000) {
681
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 const auto scale = 1.0 / math::Cbrt(2.0 * 10000.0 / double(numVoxels));
682
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 auto scaledXform = gridPtr->transform().copy();
683
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 scaledXform->preScale(scale);
684
685
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 auto newGridPtr = levelSetRebuild(*gridPtr, isovalue,
686 LEVEL_SET_HALF_WIDTH, LEVEL_SET_HALF_WIDTH, scaledXform.get(), interrupter);
687
688
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 const auto newNumVoxels = newGridPtr->activeVoxelCount();
689
2/2
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 1 times.
 10 if (newNumVoxels > numVoxels) {
690 OPENVDB_LOG_DEBUG_RUNTIME("fillWithSpheres: resampled input grid from "
691 << numVoxels << " voxel" << (numVoxels == 1 ? "" : "s")
692 << " to " << newNumVoxels << " voxel" << (newNumVoxels == 1 ? "" : "s"));
693 gridPtr = newGridPtr;
694 numVoxels = newNumVoxels;
695 }
696 }
697
698 const bool addNarrowBandPoints = (numVoxels < 10000);
699 20 int instances = std::max(instanceCount, maxSphereCount);
700
701 using TreeT = typename GridT::TreeType;
702 using BoolTreeT = typename TreeT::template ValueConverter<bool>::Type;
703 using Int16TreeT = typename TreeT::template ValueConverter<Int16>::Type;
704
705 using RandGen = std::mersenne_twister_engine<uint32_t, 32, 351, 175, 19,
706 0xccab8ee7, 11, 0xffffffff, 7, 0x31b6ab00, 15, 0xffe50000, 17, 1812433253>; // mt11213b
707 RandGen mtRand(/*seed=*/0);
708
709 const TreeT& tree = gridPtr->tree();
710
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 math::Transform transform = gridPtr->transform();
711
712 std::vector<Vec3R> instancePoints;
713 {
714 // Compute a mask of the voxels enclosed by the isosurface.
715
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 typename Grid<BoolTreeT>::Ptr interiorMaskPtr;
716
3/4
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 7 times.
✓ Branch 4 taken 3 times.
 20 if (gridPtr->getGridClass() == GRID_LEVEL_SET) {
717
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 14 interiorMaskPtr = sdfInteriorMask(*gridPtr, isovalue);
718 } else {
719 // For non-level-set grids, the interior mask comprises the active voxels.
720
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 6 interiorMaskPtr = typename Grid<BoolTreeT>::Ptr(Grid<BoolTreeT>::create(false));
721
3/8
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 3 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
 12 interiorMaskPtr->setTransform(transform.copy());
722 interiorMaskPtr->tree().topologyUnion(tree);
723 }
724
725
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 20 if (interrupter && interrupter->wasInterrupted()) return;
726
727 // If the interior mask is small and eroding it results in an empty grid,
728 // use the uneroded mask instead. (But if the minimum sphere count is zero,
729 // then eroding away the mask is acceptable.)
730
2/2
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
 20 if (!addNarrowBandPoints || (minSphereCount <= 0)) {
731
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 tools::erodeActiveValues(interiorMaskPtr->tree(), /*iterations=*/1, tools::NN_FACE, tools::IGNORE_TILES);
732
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 tools::pruneInactive(interiorMaskPtr->tree());
733 } else {
734 auto& maskTree = interiorMaskPtr->tree();
735
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 auto copyOfTree = StaticPtrCast<BoolTreeT>(maskTree.copy());
736
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 tools::erodeActiveValues(maskTree, /*iterations=*/1, tools::NN_FACE, tools::IGNORE_TILES);
737
1/2
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
 10 tools::pruneInactive(maskTree);
738
3/6
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 12 if (maskTree.empty()) { interiorMaskPtr->setTree(copyOfTree); }
739 }
740
741 // Scatter candidate sphere centroids (instancePoints)
742
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 instancePoints.reserve(instances);
743 v2s_internal::PointAccessor ptnAcc(instancePoints);
744
745
2/2
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
 20 const auto scatterCount = Index64(addNarrowBandPoints ? (instances / 2) : instances);
746
747 UniformPointScatter<v2s_internal::PointAccessor, RandGen, InterrupterT> scatter(
748 ptnAcc, scatterCount, mtRand, 1.0, interrupter);
749
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 scatter(*interiorMaskPtr);
750 }
751
752
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 20 if (interrupter && interrupter->wasInterrupted()) return;
753
754
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 40 auto csp = ClosestSurfacePoint<GridT>::create(*gridPtr, isovalue, interrupter);
755
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
 20 if (!csp) return;
756
757 // Add extra instance points in the interior narrow band.
758
2/2
✓ Branch 0 taken 6 times.
✓ Branch 1 taken 4 times.
 20 if (instancePoints.size() < size_t(instances)) {
759 const Int16TreeT& signTree = csp->signTree();
760
2/2
✓ Branch 0 taken 136 times.
✓ Branch 1 taken 6 times.
 284 for (auto leafIt = signTree.cbeginLeaf(); leafIt; ++leafIt) {
761
2/2
✓ Branch 0 taken 8264 times.
✓ Branch 1 taken 136 times.
 16800 for (auto it = leafIt->cbeginValueOn(); it; ++it) {
762
1/2
✓ Branch 1 taken 8264 times.
✗ Branch 2 not taken.
 16528 const int flags = int(it.getValue());
763
2/2
✓ Branch 0 taken 920 times.
✓ Branch 1 taken 7344 times.
 16528 if (!(volume_to_mesh_internal::EDGES & flags)
764 && (volume_to_mesh_internal::INSIDE & flags))
765 {
766
1/2
✓ Branch 1 taken 920 times.
✗ Branch 2 not taken.
 3680 instancePoints.push_back(transform.indexToWorld(it.getCoord()));
767 }
768
1/2
✓ Branch 0 taken 8264 times.
✗ Branch 1 not taken.
 16528 if (instancePoints.size() == size_t(instances)) break;
769 }
770
1/2
✓ Branch 0 taken 136 times.
✗ Branch 1 not taken.
 272 if (instancePoints.size() == size_t(instances)) break;
771 }
772 }
773
774
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 20 if (interrupter && interrupter->wasInterrupted()) return;
775
776 // Assign a radius to each candidate sphere. The radius is the world-space
777 // distance from the sphere's center to the closest surface point.
778 std::vector<float> instanceRadius;
779
2/4
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
 20 if (!csp->search(instancePoints, instanceRadius)) return;
780
781 20 float largestRadius = 0.0;
782 int largestRadiusIdx = 0;
783
2/2
✓ Branch 0 taken 65920 times.
✓ Branch 1 taken 10 times.
 131860 for (size_t n = 0, N = instancePoints.size(); n < N; ++n) {
784
2/2
✓ Branch 0 taken 292 times.
✓ Branch 1 taken 65628 times.
 131840 if (instanceRadius[n] > largestRadius) {
785 584 largestRadius = instanceRadius[n];
786 584 largestRadiusIdx = int(n);
787 }
788 }
789
790
1/4
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 20 std::vector<unsigned char> instanceMask(instancePoints.size(), 0);
791
792
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 minRadius = float(minRadius * transform.voxelSize()[0]);
793
3/4
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 9 times.
 20 maxRadius = float(maxRadius * transform.voxelSize()[0]);
794
795
4/4
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 9 times.
✓ Branch 2 taken 52 times.
✓ Branch 3 taken 1 times.
 108 for (size_t s = 0, S = std::min(size_t(maxSphereCount), instancePoints.size()); s < S; ++s) {
796
797
1/8
✗ Branch 0 not taken.
✓ Branch 1 taken 52 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
 104 if (interrupter && interrupter->wasInterrupted()) return;
798
799 104 largestRadius = std::min(maxRadius, largestRadius);
800
801
3/4
✓ Branch 0 taken 9 times.
✓ Branch 1 taken 43 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 9 times.
 104 if ((int(s) >= minSphereCount) && (largestRadius < minRadius)) break;
802
803
1/2
✓ Branch 1 taken 43 times.
✗ Branch 2 not taken.
 86 const Vec4s sphere(
804 float(instancePoints[largestRadiusIdx].x()),
805 float(instancePoints[largestRadiusIdx].y()),
806
1/2
✓ Branch 1 taken 43 times.
✗ Branch 2 not taken.
 86 float(instancePoints[largestRadiusIdx].z()),
807 largestRadius);
808
809
1/2
✓ Branch 1 taken 43 times.
✗ Branch 2 not taken.
 86 spheres.push_back(sphere);
810 86 instanceMask[largestRadiusIdx] = 1;
811
812 v2s_internal::UpdatePoints op(
813 sphere, instancePoints, instanceRadius, instanceMask, overlapping);
814
1/2
✓ Branch 1 taken 43 times.
✗ Branch 2 not taken.
 86 op.run();
815
816 86 largestRadius = op.radius();
817 largestRadiusIdx = op.index();
818 }
819 } // fillWithSpheres
820
821
822 ////////////////////////////////////////
823
824
825 template<typename GridT>
826 template<typename InterrupterT>
827 typename ClosestSurfacePoint<GridT>::Ptr
828 12 ClosestSurfacePoint<GridT>::create(const GridT& grid, float isovalue, InterrupterT* interrupter)
829 {
830 12 auto csp = Ptr{new ClosestSurfacePoint};
831
2/4
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 12 times.
 12 if (!csp->initialize(grid, isovalue, interrupter)) csp.reset();
832 12 return csp;
833 }
834
835
836 template<typename GridT>
837 template<typename InterrupterT>
838 bool
839 12 ClosestSurfacePoint<GridT>::initialize(
840 const GridT& grid, float isovalue, InterrupterT* interrupter)
841 {
842 using Index32LeafManagerT = tree::LeafManager<Index32TreeT>;
843 using ValueT = typename GridT::ValueType;
844
845 const TreeT& tree = grid.tree();
846 const math::Transform& transform = grid.transform();
847
848 { // Extract surface point cloud
849
850 24 BoolTreeT mask(false);
851
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 volume_to_mesh_internal::identifySurfaceIntersectingVoxels(mask, tree, ValueT(isovalue));
852
853
2/4
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 12 times.
✗ Branch 6 not taken.
 12 mSignTreePt.reset(new Int16TreeT(0));
854
3/6
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 12 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 12 times.
✗ Branch 9 not taken.
 12 mIdxTreePt.reset(new Index32TreeT(std::numeric_limits<Index32>::max()));
855
856
857
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 volume_to_mesh_internal::computeAuxiliaryData(
858 *mSignTreePt, *mIdxTreePt, mask, tree, ValueT(isovalue));
859
860
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 12 if (interrupter && interrupter->wasInterrupted()) return false;
861
862 // count unique points
863
864 using Int16LeafNodeType = typename Int16TreeT::LeafNodeType;
865 using Index32LeafNodeType = typename Index32TreeT::LeafNodeType;
866
867 std::vector<Int16LeafNodeType*> signFlagsLeafNodes;
868 mSignTreePt->getNodes(signFlagsLeafNodes);
869
870 const tbb::blocked_range<size_t> auxiliaryLeafNodeRange(0, signFlagsLeafNodes.size());
871
872
2/4
✓ Branch 0 taken 12 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 12 times.
✗ Branch 4 not taken.
 12 std::unique_ptr<Index32[]> leafNodeOffsets(new Index32[signFlagsLeafNodes.size()]);
873
874
1/4
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 12 tbb::parallel_for(auxiliaryLeafNodeRange,
875 volume_to_mesh_internal::LeafNodePointCount<Int16LeafNodeType::LOG2DIM>
876 (signFlagsLeafNodes, leafNodeOffsets));
877
878 {
879 Index32 pointCount = 0;
880
2/2
✓ Branch 0 taken 1354 times.
✓ Branch 1 taken 12 times.
 1366 for (size_t n = 0, N = signFlagsLeafNodes.size(); n < N; ++n) {
881 1354 const Index32 tmp = leafNodeOffsets[n];
882 1354 leafNodeOffsets[n] = pointCount;
883 1354 pointCount += tmp;
884 }
885
886 12 mPointListSize = size_t(pointCount);
887
3/4
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 86624 times.
✓ Branch 4 taken 12 times.
 86636 mSurfacePointList.reset(new Vec3s[mPointListSize]);
888 }
889
890
891 std::vector<Index32LeafNodeType*> pointIndexLeafNodes;
892 mIdxTreePt->getNodes(pointIndexLeafNodes);
893
894
4/8
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 10 times.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 24 tbb::parallel_for(auxiliaryLeafNodeRange, volume_to_mesh_internal::ComputePoints<TreeT>(
895 mSurfacePointList.get(), tree, pointIndexLeafNodes,
896 signFlagsLeafNodes, leafNodeOffsets, transform, ValueT(isovalue)));
897 }
898
899
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 12 if (interrupter && interrupter->wasInterrupted()) return false;
900
901 24 Index32LeafManagerT idxLeafs(*mIdxTreePt);
902
903 using Index32RootNodeT = typename Index32TreeT::RootNodeType;
904 using Index32NodeChainT = typename Index32RootNodeT::NodeChainType;
905 static_assert(Index32NodeChainT::Size > 1,
906 "expected tree depth greater than one");
907 using Index32InternalNodeT = typename Index32NodeChainT::template Get<1>;
908
909 typename Index32TreeT::NodeCIter nIt = mIdxTreePt->cbeginNode();
910 nIt.setMinDepth(Index32TreeT::NodeCIter::LEAF_DEPTH - 1);
911 12 nIt.setMaxDepth(Index32TreeT::NodeCIter::LEAF_DEPTH - 1);
912
913 std::vector<const Index32InternalNodeT*> internalNodes;
914
915 12 const Index32InternalNodeT* node = nullptr;
916
2/2
✓ Branch 0 taken 41 times.
✓ Branch 1 taken 12 times.
 53 for (; nIt; ++nIt) {
917 nIt.getNode(node);
918
2/4
✓ Branch 0 taken 41 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 41 times.
✗ Branch 4 not taken.
 41 if (node) internalNodes.push_back(node);
919 }
920
921
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
 12 std::vector<IndexRange>().swap(mLeafRanges);
922
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 mLeafRanges.resize(internalNodes.size());
923
924
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
 12 std::vector<const Index32LeafT*>().swap(mLeafNodes);
925
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 mLeafNodes.reserve(idxLeafs.leafCount());
926
927 typename Index32InternalNodeT::ChildOnCIter leafIt;
928 12 mMaxNodeLeafs = 0;
929
2/2
✓ Branch 0 taken 41 times.
✓ Branch 1 taken 12 times.
 53 for (size_t n = 0, N = internalNodes.size(); n < N; ++n) {
930
931 41 mLeafRanges[n].first = mLeafNodes.size();
932
933 41 size_t leafCount = 0;
934
2/2
✓ Branch 1 taken 1354 times.
✓ Branch 2 taken 41 times.
 1436 for (leafIt = internalNodes[n]->cbeginChildOn(); leafIt; ++leafIt) {
935
1/4
✓ Branch 1 taken 1354 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 1354 mLeafNodes.push_back(&(*leafIt));
936 1354 ++leafCount;
937 }
938
939
2/2
✓ Branch 0 taken 6 times.
✓ Branch 1 taken 35 times.
 47 mMaxNodeLeafs = std::max(leafCount, mMaxNodeLeafs);
940
941 41 mLeafRanges[n].second = mLeafNodes.size();
942 }
943
944
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
 12 std::vector<Vec4R>().swap(mLeafBoundingSpheres);
945
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 mLeafBoundingSpheres.resize(mLeafNodes.size());
946
947 12 v2s_internal::LeafOp<Index32LeafT> leafBS(
948 mLeafBoundingSpheres, mLeafNodes, transform, mSurfacePointList);
949
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 leafBS.run();
950
951
952
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 12 times.
 12 std::vector<Vec4R>().swap(mNodeBoundingSpheres);
953
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 mNodeBoundingSpheres.resize(internalNodes.size());
954
955 v2s_internal::NodeOp nodeBS(mNodeBoundingSpheres, mLeafRanges, mLeafBoundingSpheres);
956
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 nodeBS.run();
957 return true;
958 } // ClosestSurfacePoint::initialize
959
960
961 template<typename GridT>
962 bool
963
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
 28 ClosestSurfacePoint<GridT>::search(std::vector<Vec3R>& points,
964 std::vector<float>& distances, bool transformPoints)
965 {
966 distances.clear();
967 28 distances.resize(points.size(), std::numeric_limits<float>::infinity());
968
969 28 v2s_internal::ClosestPointDist<Index32LeafT> cpd(points, distances, mSurfacePointList,
970 28 mLeafNodes, mLeafRanges, mLeafBoundingSpheres, mNodeBoundingSpheres,
971 mMaxNodeLeafs, transformPoints);
972
973
1/2
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
 28 cpd.run();
974
975 28 return true;
976 }
977
978
979 template<typename GridT>
980 bool
981 20 ClosestSurfacePoint<GridT>::search(const std::vector<Vec3R>& points, std::vector<float>& distances)
982 {
983 20 return search(const_cast<std::vector<Vec3R>& >(points), distances, false);
984 }
985
986
987 template<typename GridT>
988 bool
989 4 ClosestSurfacePoint<GridT>::searchAndReplace(std::vector<Vec3R>& points,
990 std::vector<float>& distances)
991 {
992 4 return search(points, distances, true);
993 }
994
995
996 ////////////////////////////////////////
997
998
999 // Explicit Template Instantiation
1000
1001 #ifdef OPENVDB_USE_EXPLICIT_INSTANTIATION
1002
1003 #ifdef OPENVDB_INSTANTIATE_VOLUMETOSPHERES
1004 #include <openvdb/util/ExplicitInstantiation.h>
1005 #endif
1006
1007 OPENVDB_INSTANTIATE_CLASS ClosestSurfacePoint<FloatGrid>;
1008 OPENVDB_INSTANTIATE_CLASS ClosestSurfacePoint<DoubleGrid>;
1009
1010 #define _FUNCTION(TreeT) \
1011 void fillWithSpheres(const Grid<TreeT>&, std::vector<openvdb::Vec4s>&, const Vec2i&, \
1012 bool, float, float, float, int, util::NullInterrupter*)
1013 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
1014 #undef _FUNCTION
1015
1016 #endif // OPENVDB_USE_EXPLICIT_INSTANTIATION
1017
1018
1019 } // namespace tools
1020 } // namespace OPENVDB_VERSION_NAME
1021 } // namespace openvdb
1022
1023 #endif // OPENVDB_TOOLS_VOLUME_TO_MESH_HAS_BEEN_INCLUDED
1024