GCC Code Coverage Report

Directory: ./
File: openvdb/openvdb/tools/MeshToVolume.h
Date: 2022-07-25 17:40:05
Exec Total Coverage
Lines: 953 973 97.9%
Functions: 87 178 48.9%
Branches: 1074 1676 64.1%
Line Branch Exec Source
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3
4 /// @file MeshToVolume.h
5 ///
6 /// @brief Convert polygonal meshes that consist of quads and/or triangles
7 /// into signed or unsigned distance field volumes.
8 ///
9 /// @note The signed distance field conversion requires a closed surface
10 /// but not necessarily a manifold surface. Supports surfaces with
11 /// self intersections and degenerate faces and is independent of
12 /// mesh surface normals / polygon orientation.
13 ///
14 /// @author Mihai Alden
15
16 #ifndef OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
17 #define OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
18
19 #include "openvdb/Platform.h" // for OPENVDB_HAS_CXX11
20 #include "openvdb/Types.h"
21 #include "openvdb/math/FiniteDifference.h" // for GodunovsNormSqrd
22 #include "openvdb/math/Proximity.h" // for closestPointOnTriangleToPoint
23 #include "openvdb/util/NullInterrupter.h"
24 #include "openvdb/util/Util.h"
25 #include "openvdb/thread/Threading.h"
26 #include <openvdb/openvdb.h>
27
28 #include "ChangeBackground.h"
29 #include "Prune.h" // for pruneInactive and pruneLevelSet
30 #include "SignedFloodFill.h" // for signedFloodFillWithValues
31
32 #include <tbb/blocked_range.h>
33 #include <tbb/enumerable_thread_specific.h>
34 #include <tbb/parallel_for.h>
35 #include <tbb/parallel_reduce.h>
36 #include <tbb/partitioner.h>
37 #include <tbb/task_group.h>
38 #include <tbb/task_arena.h>
39
40 #include <algorithm> // for std::sort()
41 #include <cmath> // for std::isfinite(), std::isnan()
42 #include <deque>
43 #include <limits>
44 #include <memory>
45 #include <sstream>
46 #include <type_traits>
47 #include <vector>
48
49 namespace openvdb {
50 OPENVDB_USE_VERSION_NAMESPACE
51 namespace OPENVDB_VERSION_NAME {
52 namespace tools {
53
54
55 ////////////////////////////////////////
56
57
58 /// @brief Mesh to volume conversion flags
59 enum MeshToVolumeFlags {
60
61 /// Switch from the default signed distance field conversion that classifies
62 /// regions as either inside or outside the mesh boundary to a unsigned distance
63 /// field conversion that only computes distance values. This conversion type
64 /// does not require a closed watertight mesh.
65 UNSIGNED_DISTANCE_FIELD = 0x1,
66
67 /// Disable the cleanup step that removes voxels created by self intersecting
68 /// portions of the mesh.
69 DISABLE_INTERSECTING_VOXEL_REMOVAL = 0x2,
70
71 /// Disable the distance renormalization step that smooths out bumps caused
72 /// by self intersecting or overlapping portions of the mesh
73 DISABLE_RENORMALIZATION = 0x4,
74
75 /// Disable the cleanup step that removes active voxels that exceed the
76 /// narrow band limits. (Only relevant for small limits)
77 DISABLE_NARROW_BAND_TRIMMING = 0x8
78 };
79
80
81 /// @brief Convert polygonal meshes that consist of quads and/or triangles into
82 /// signed or unsigned distance field volumes.
83 ///
84 /// @note Requires a closed surface but not necessarily a manifold surface.
85 /// Supports surfaces with self intersections and degenerate faces
86 /// and is independent of mesh surface normals.
87 ///
88 /// @interface MeshDataAdapter
89 /// Expected interface for the MeshDataAdapter class
90 /// @code
91 /// struct MeshDataAdapter {
92 /// size_t polygonCount() const; // Total number of polygons
93 /// size_t pointCount() const; // Total number of points
94 /// size_t vertexCount(size_t n) const; // Vertex count for polygon n
95 ///
96 /// // Return position pos in local grid index space for polygon n and vertex v
97 /// void getIndexSpacePoint(size_t n, size_t v, openvdb::Vec3d& pos) const;
98 /// };
99 /// @endcode
100 ///
101 /// @param mesh mesh data access class that conforms to the MeshDataAdapter
102 /// interface
103 /// @param transform world-to-index-space transform
104 /// @param exteriorBandWidth exterior narrow band width in voxel units
105 /// @param interiorBandWidth interior narrow band width in voxel units
106 /// (set to std::numeric_limits<float>::max() to fill object
107 /// interior with distance values)
108 /// @param flags optional conversion flags defined in @c MeshToVolumeFlags
109 /// @param polygonIndexGrid optional grid output that will contain the closest-polygon
110 /// index for each voxel in the narrow band region
111 template <typename GridType, typename MeshDataAdapter>
112 typename GridType::Ptr
113 meshToVolume(
114 const MeshDataAdapter& mesh,
115 const math::Transform& transform,
116 float exteriorBandWidth = 3.0f,
117 float interiorBandWidth = 3.0f,
118 int flags = 0,
119 typename GridType::template ValueConverter<Int32>::Type * polygonIndexGrid = nullptr);
120
121
122 /// @brief Convert polygonal meshes that consist of quads and/or triangles into
123 /// signed or unsigned distance field volumes.
124 ///
125 /// @param interrupter a callback to interrupt the conversion process that conforms
126 /// to the util::NullInterrupter interface
127 /// @param mesh mesh data access class that conforms to the MeshDataAdapter
128 /// interface
129 /// @param transform world-to-index-space transform
130 /// @param exteriorBandWidth exterior narrow band width in voxel units
131 /// @param interiorBandWidth interior narrow band width in voxel units (set this value to
132 /// std::numeric_limits<float>::max() to fill interior regions
133 /// with distance values)
134 /// @param flags optional conversion flags defined in @c MeshToVolumeFlags
135 /// @param polygonIndexGrid optional grid output that will contain the closest-polygon
136 /// index for each voxel in the active narrow band region
137 template <typename GridType, typename MeshDataAdapter, typename Interrupter>
138 typename GridType::Ptr
139 meshToVolume(
140 Interrupter& interrupter,
141 const MeshDataAdapter& mesh,
142 const math::Transform& transform,
143 float exteriorBandWidth = 3.0f,
144 float interiorBandWidth = 3.0f,
145 int flags = 0,
146 typename GridType::template ValueConverter<Int32>::Type * polygonIndexGrid = nullptr);
147
148
149 ////////////////////////////////////////
150
151
152 /// @brief Contiguous quad and triangle data adapter class
153 ///
154 /// @details PointType and PolygonType must provide element access
155 /// through the square brackets operator.
156 /// @details Points are assumed to be in local grid index space.
157 /// @details The PolygonType tuple can have either three or four components
158 /// this property must be specified in a static member variable
159 /// named @c size, similar to the math::Tuple class.
160 /// @details A four component tuple can represent a quads or a triangle
161 /// if the fourth component set to @c util::INVALID_INDEX
162 template<typename PointType, typename PolygonType>
163 struct QuadAndTriangleDataAdapter {
164
165
3/4
✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 1 times.
 6 QuadAndTriangleDataAdapter(const std::vector<PointType>& points,
166 const std::vector<PolygonType>& polygons)
167
3/4
✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 1 times.
 6 : mPointArray(points.empty() ? nullptr : &points[0])
168 , mPointArraySize(points.size())
169
3/4
✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 1 times.
 6 , mPolygonArray(polygons.empty() ? nullptr : &polygons[0])
170
5/8
✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 4 times.
✓ Branch 5 taken 1 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 5 times.
 11 , mPolygonArraySize(polygons.size())
171 {
172 }
173
174 29 QuadAndTriangleDataAdapter(const PointType * pointArray, size_t pointArraySize,
175 const PolygonType* polygonArray, size_t polygonArraySize)
176 : mPointArray(pointArray)
177 , mPointArraySize(pointArraySize)
178 , mPolygonArray(polygonArray)
179 29 , mPolygonArraySize(polygonArraySize)
180 {
181 }
182
183
2/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 32 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
 35 size_t polygonCount() const { return mPolygonArraySize; }
184 size_t pointCount() const { return mPointArraySize; }
185
186 /// @brief Vertex count for polygon @a n
187 size_t vertexCount(size_t n) const {
188 3283240 return (PolygonType::size == 3 || mPolygonArray[n][3] == util::INVALID_IDX) ? 3 : 4;
189 }
190
191 /// @brief Returns position @a pos in local grid index space
192 /// for polygon @a n and vertex @a v
193 25875896 void getIndexSpacePoint(size_t n, size_t v, Vec3d& pos) const {
194
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 12937948 times.
 25875896 const PointType& p = mPointArray[mPolygonArray[n][int(v)]];
195 25875896 pos[0] = double(p[0]);
196 25875896 pos[1] = double(p[1]);
197 25875896 pos[2] = double(p[2]);
198 25875896 }
199
200 private:
201 PointType const * const mPointArray;
202 size_t const mPointArraySize;
203 PolygonType const * const mPolygonArray;
204 size_t const mPolygonArraySize;
205 }; // struct QuadAndTriangleDataAdapter
206
207
208 ////////////////////////////////////////
209
210
211 // Convenience functions for the mesh to volume converter that wrap stl containers.
212 //
213 // Note the meshToVolume() method declared above is more flexible and better suited
214 // for arbitrary data structures.
215
216
217 /// @brief Convert a triangle mesh to a level set volume.
218 ///
219 /// @return a grid of type @c GridType containing a narrow-band level set
220 /// representation of the input mesh.
221 ///
222 /// @throw TypeError if @c GridType is not scalar or not floating-point
223 ///
224 /// @note Requires a closed surface but not necessarily a manifold surface.
225 /// Supports surfaces with self intersections and degenerate faces
226 /// and is independent of mesh surface normals.
227 ///
228 /// @param xform transform for the output grid
229 /// @param points list of world space point positions
230 /// @param triangles triangle index list
231 /// @param halfWidth half the width of the narrow band, in voxel units
232 template<typename GridType>
233 typename GridType::Ptr
234 meshToLevelSet(
235 const openvdb::math::Transform& xform,
236 const std::vector<Vec3s>& points,
237 const std::vector<Vec3I>& triangles,
238 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
239
240 /// Adds support for a @a interrupter callback used to cancel the conversion.
241 template<typename GridType, typename Interrupter>
242 typename GridType::Ptr
243 meshToLevelSet(
244 Interrupter& interrupter,
245 const openvdb::math::Transform& xform,
246 const std::vector<Vec3s>& points,
247 const std::vector<Vec3I>& triangles,
248 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
249
250
251 /// @brief Convert a quad mesh to a level set volume.
252 ///
253 /// @return a grid of type @c GridType containing a narrow-band level set
254 /// representation of the input mesh.
255 ///
256 /// @throw TypeError if @c GridType is not scalar or not floating-point
257 ///
258 /// @note Requires a closed surface but not necessarily a manifold surface.
259 /// Supports surfaces with self intersections and degenerate faces
260 /// and is independent of mesh surface normals.
261 ///
262 /// @param xform transform for the output grid
263 /// @param points list of world space point positions
264 /// @param quads quad index list
265 /// @param halfWidth half the width of the narrow band, in voxel units
266 template<typename GridType>
267 typename GridType::Ptr
268 meshToLevelSet(
269 const openvdb::math::Transform& xform,
270 const std::vector<Vec3s>& points,
271 const std::vector<Vec4I>& quads,
272 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
273
274 /// Adds support for a @a interrupter callback used to cancel the conversion.
275 template<typename GridType, typename Interrupter>
276 typename GridType::Ptr
277 meshToLevelSet(
278 Interrupter& interrupter,
279 const openvdb::math::Transform& xform,
280 const std::vector<Vec3s>& points,
281 const std::vector<Vec4I>& quads,
282 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
283
284
285 /// @brief Convert a triangle and quad mesh to a level set volume.
286 ///
287 /// @return a grid of type @c GridType containing a narrow-band level set
288 /// representation of the input mesh.
289 ///
290 /// @throw TypeError if @c GridType is not scalar or not floating-point
291 ///
292 /// @note Requires a closed surface but not necessarily a manifold surface.
293 /// Supports surfaces with self intersections and degenerate faces
294 /// and is independent of mesh surface normals.
295 ///
296 /// @param xform transform for the output grid
297 /// @param points list of world space point positions
298 /// @param triangles triangle index list
299 /// @param quads quad index list
300 /// @param halfWidth half the width of the narrow band, in voxel units
301 template<typename GridType>
302 typename GridType::Ptr
303 meshToLevelSet(
304 const openvdb::math::Transform& xform,
305 const std::vector<Vec3s>& points,
306 const std::vector<Vec3I>& triangles,
307 const std::vector<Vec4I>& quads,
308 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
309
310 /// Adds support for a @a interrupter callback used to cancel the conversion.
311 template<typename GridType, typename Interrupter>
312 typename GridType::Ptr
313 meshToLevelSet(
314 Interrupter& interrupter,
315 const openvdb::math::Transform& xform,
316 const std::vector<Vec3s>& points,
317 const std::vector<Vec3I>& triangles,
318 const std::vector<Vec4I>& quads,
319 float halfWidth = float(LEVEL_SET_HALF_WIDTH));
320
321
322 /// @brief Convert a triangle and quad mesh to a signed distance field
323 /// with an asymmetrical narrow band.
324 ///
325 /// @return a grid of type @c GridType containing a narrow-band signed
326 /// distance field representation of the input mesh.
327 ///
328 /// @throw TypeError if @c GridType is not scalar or not floating-point
329 ///
330 /// @note Requires a closed surface but not necessarily a manifold surface.
331 /// Supports surfaces with self intersections and degenerate faces
332 /// and is independent of mesh surface normals.
333 ///
334 /// @param xform transform for the output grid
335 /// @param points list of world space point positions
336 /// @param triangles triangle index list
337 /// @param quads quad index list
338 /// @param exBandWidth the exterior narrow-band width in voxel units
339 /// @param inBandWidth the interior narrow-band width in voxel units
340 template<typename GridType>
341 typename GridType::Ptr
342 meshToSignedDistanceField(
343 const openvdb::math::Transform& xform,
344 const std::vector<Vec3s>& points,
345 const std::vector<Vec3I>& triangles,
346 const std::vector<Vec4I>& quads,
347 float exBandWidth,
348 float inBandWidth);
349
350 /// Adds support for a @a interrupter callback used to cancel the conversion.
351 template<typename GridType, typename Interrupter>
352 typename GridType::Ptr
353 meshToSignedDistanceField(
354 Interrupter& interrupter,
355 const openvdb::math::Transform& xform,
356 const std::vector<Vec3s>& points,
357 const std::vector<Vec3I>& triangles,
358 const std::vector<Vec4I>& quads,
359 float exBandWidth,
360 float inBandWidth);
361
362
363 /// @brief Convert a triangle and quad mesh to an unsigned distance field.
364 ///
365 /// @return a grid of type @c GridType containing a narrow-band unsigned
366 /// distance field representation of the input mesh.
367 ///
368 /// @throw TypeError if @c GridType is not scalar or not floating-point
369 ///
370 /// @note Does not requires a closed surface.
371 ///
372 /// @param xform transform for the output grid
373 /// @param points list of world space point positions
374 /// @param triangles triangle index list
375 /// @param quads quad index list
376 /// @param bandWidth the width of the narrow band, in voxel units
377 template<typename GridType>
378 typename GridType::Ptr
379 meshToUnsignedDistanceField(
380 const openvdb::math::Transform& xform,
381 const std::vector<Vec3s>& points,
382 const std::vector<Vec3I>& triangles,
383 const std::vector<Vec4I>& quads,
384 float bandWidth);
385
386 /// Adds support for a @a interrupter callback used to cancel the conversion.
387 template<typename GridType, typename Interrupter>
388 typename GridType::Ptr
389 meshToUnsignedDistanceField(
390 Interrupter& interrupter,
391 const openvdb::math::Transform& xform,
392 const std::vector<Vec3s>& points,
393 const std::vector<Vec3I>& triangles,
394 const std::vector<Vec4I>& quads,
395 float bandWidth);
396
397
398 ////////////////////////////////////////
399
400
401 /// @brief Return a grid of type @c GridType containing a narrow-band level set
402 /// representation of a box.
403 ///
404 /// @param bbox a bounding box in world units
405 /// @param xform world-to-index-space transform
406 /// @param halfWidth half the width of the narrow band, in voxel units
407 template<typename GridType, typename VecType>
408 typename GridType::Ptr
409 createLevelSetBox(const math::BBox<VecType>& bbox,
410 const openvdb::math::Transform& xform,
411 typename VecType::ValueType halfWidth = LEVEL_SET_HALF_WIDTH);
412
413
414 ////////////////////////////////////////
415
416
417 /// @brief Traces the exterior voxel boundary of closed objects in the input
418 /// volume @a tree. Exterior voxels are marked with a negative sign,
419 /// voxels with a value below @c 0.75 are left unchanged and act as
420 /// the boundary layer.
421 ///
422 /// @note Does not propagate sign information into tile regions.
423 template <typename FloatTreeT>
424 void
425 traceExteriorBoundaries(FloatTreeT& tree);
426
427
428 ////////////////////////////////////////
429
430
431 /// @brief Extracts and stores voxel edge intersection data from a mesh.
432 class MeshToVoxelEdgeData
433 {
434 public:
435
436 //////////
437
438 ///@brief Internal edge data type.
439 struct EdgeData {
440 EdgeData(float dist = 1.0)
441 : mXDist(dist), mYDist(dist), mZDist(dist)
442 , mXPrim(util::INVALID_IDX)
443 , mYPrim(util::INVALID_IDX)
444 , mZPrim(util::INVALID_IDX)
445 {
446 }
447
448 //@{
449 /// Required by several of the tree nodes
450 /// @note These methods don't perform meaningful operations.
451 bool operator< (const EdgeData&) const { return false; }
452 bool operator> (const EdgeData&) const { return false; }
453 template<class T> EdgeData operator+(const T&) const { return *this; }
454 template<class T> EdgeData operator-(const T&) const { return *this; }
455 EdgeData operator-() const { return *this; }
456 //@}
457
458 bool operator==(const EdgeData& rhs) const
459 {
460 return mXPrim == rhs.mXPrim && mYPrim == rhs.mYPrim && mZPrim == rhs.mZPrim;
461 }
462
463 float mXDist, mYDist, mZDist;
464 Index32 mXPrim, mYPrim, mZPrim;
465 };
466
467 using TreeType = tree::Tree4<EdgeData, 5, 4, 3>::Type;
468 using Accessor = tree::ValueAccessor<TreeType>;
469
470
471 //////////
472
473
474 MeshToVoxelEdgeData();
475
476
477 /// @brief Threaded method to extract voxel edge data, the closest
478 /// intersection point and corresponding primitive index,
479 /// from the given mesh.
480 ///
481 /// @param pointList List of points in grid index space, preferably unique
482 /// and shared by different polygons.
483 /// @param polygonList List of triangles and/or quads.
484 void convert(const std::vector<Vec3s>& pointList, const std::vector<Vec4I>& polygonList);
485
486
487 /// @brief Returns intersection points with corresponding primitive
488 /// indices for the given @c ijk voxel.
489 void getEdgeData(Accessor& acc, const Coord& ijk,
490 std::vector<Vec3d>& points, std::vector<Index32>& primitives);
491
492 /// @return An accessor of @c MeshToVoxelEdgeData::Accessor type that
493 /// provides random read access to the internal tree.
494 Accessor getAccessor() { return Accessor(mTree); }
495
496 private:
497 void operator=(const MeshToVoxelEdgeData&) {}
498 TreeType mTree;
499 class GenEdgeData;
500 };
501
502
503 ////////////////////////////////////////////////////////////////////////////////
504 ////////////////////////////////////////////////////////////////////////////////
505
506
507 // Internal utility objects and implementation details
508
509 /// @cond OPENVDB_DOCS_INTERNAL
510
511 namespace mesh_to_volume_internal {
512
513 template<typename PointType>
514 struct TransformPoints {
515
516 15 TransformPoints(const PointType* pointsIn, PointType* pointsOut,
517 const math::Transform& xform)
518
1/4
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 15 times.
✗ Branch 5 not taken.
 15 : mPointsIn(pointsIn), mPointsOut(pointsOut), mXform(&xform)
519 {
520 }
521
522 142 void operator()(const tbb::blocked_range<size_t>& range) const {
523
524 Vec3d pos;
525
526
2/2
✓ Branch 0 taken 142 times.
✓ Branch 1 taken 142 times.
 284 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
527
528 142 const PointType& wsP = mPointsIn[n];
529 142 pos[0] = double(wsP[0]);
530 142 pos[1] = double(wsP[1]);
531 142 pos[2] = double(wsP[2]);
532
533 142 pos = mXform->worldToIndex(pos);
534
535 142 PointType& isP = mPointsOut[n];
536 142 isP[0] = typename PointType::value_type(pos[0]);
537 142 isP[1] = typename PointType::value_type(pos[1]);
538 142 isP[2] = typename PointType::value_type(pos[2]);
539 }
540 142 }
541
542 PointType const * const mPointsIn;
543 PointType * const mPointsOut;
544 math::Transform const * const mXform;
545 }; // TransformPoints
546
547
548 template<typename ValueType>
549 struct Tolerance
550 {
551 static ValueType epsilon() { return ValueType(1e-7); }
552 static ValueType minNarrowBandWidth() { return ValueType(1.0 + 1e-6); }
553 };
554
555
556 ////////////////////////////////////////
557
558
559 template<typename TreeType>
560 class CombineLeafNodes
561 {
562 public:
563
564 using Int32TreeType = typename TreeType::template ValueConverter<Int32>::Type;
565
566 using LeafNodeType = typename TreeType::LeafNodeType;
567 using Int32LeafNodeType = typename Int32TreeType::LeafNodeType;
568
569 34 CombineLeafNodes(TreeType& lhsDistTree, Int32TreeType& lhsIdxTree,
570 LeafNodeType ** rhsDistNodes, Int32LeafNodeType ** rhsIdxNodes)
571 : mDistTree(&lhsDistTree)
572 , mIdxTree(&lhsIdxTree)
573 , mRhsDistNodes(rhsDistNodes)
574
1/4
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 34 times.
✗ Branch 5 not taken.
 34 , mRhsIdxNodes(rhsIdxNodes)
575 {
576 }
577
578 4646 void operator()(const tbb::blocked_range<size_t>& range) const {
579
580 4646 tree::ValueAccessor<TreeType> distAcc(*mDistTree);
581
1/2
✓ Branch 1 taken 2323 times.
✗ Branch 2 not taken.
 4646 tree::ValueAccessor<Int32TreeType> idxAcc(*mIdxTree);
582
583 using DistValueType = typename LeafNodeType::ValueType;
584 using IndexValueType = typename Int32LeafNodeType::ValueType;
585
586
2/2
✓ Branch 0 taken 2671 times.
✓ Branch 1 taken 2323 times.
 9988 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
587
588 5342 const Coord& origin = mRhsDistNodes[n]->origin();
589
590 LeafNodeType* lhsDistNode = distAcc.probeLeaf(origin);
591 Int32LeafNodeType* lhsIdxNode = idxAcc.probeLeaf(origin);
592
593
1/2
✓ Branch 1 taken 2671 times.
✗ Branch 2 not taken.
 5342 DistValueType* lhsDistData = lhsDistNode->buffer().data();
594
1/2
✓ Branch 1 taken 2671 times.
✗ Branch 2 not taken.
 5342 IndexValueType* lhsIdxData = lhsIdxNode->buffer().data();
595
596
1/2
✓ Branch 1 taken 2671 times.
✗ Branch 2 not taken.
 5342 const DistValueType* rhsDistData = mRhsDistNodes[n]->buffer().data();
597
1/2
✓ Branch 1 taken 2671 times.
✗ Branch 2 not taken.
 5342 const IndexValueType* rhsIdxData = mRhsIdxNodes[n]->buffer().data();
598
599
600
2/2
✓ Branch 0 taken 1367552 times.
✓ Branch 1 taken 2671 times.
 2740446 for (Index32 offset = 0; offset < LeafNodeType::SIZE; ++offset) {
601
602
2/2
✓ Branch 0 taken 276792 times.
✓ Branch 1 taken 1090760 times.
 2735104 if (rhsIdxData[offset] != Int32(util::INVALID_IDX)) {
603
604 553584 const DistValueType& lhsValue = lhsDistData[offset];
605 553584 const DistValueType& rhsValue = rhsDistData[offset];
606
607
2/2
✓ Branch 0 taken 169518 times.
✓ Branch 1 taken 107274 times.
 553584 if (rhsValue < lhsValue) {
608 lhsDistNode->setValueOn(offset, rhsValue);
609 lhsIdxNode->setValueOn(offset, rhsIdxData[offset]);
610
2/2
✓ Branch 0 taken 57263 times.
✓ Branch 1 taken 50011 times.
 214548 } else if (math::isExactlyEqual(rhsValue, lhsValue)) {
611 lhsIdxNode->setValueOn(offset,
612
2/2
✓ Branch 0 taken 711 times.
✓ Branch 1 taken 56552 times.
 114526 std::min(lhsIdxData[offset], rhsIdxData[offset]));
613 }
614 }
615 }
616
617
1/2
✓ Branch 0 taken 2671 times.
✗ Branch 1 not taken.
 10684 delete mRhsDistNodes[n];
618
1/2
✓ Branch 0 taken 2671 times.
✗ Branch 1 not taken.
 10684 delete mRhsIdxNodes[n];
619 }
620 }
621
622 private:
623
624 TreeType * const mDistTree;
625 Int32TreeType * const mIdxTree;
626
627 LeafNodeType ** const mRhsDistNodes;
628 Int32LeafNodeType ** const mRhsIdxNodes;
629 }; // class CombineLeafNodes
630
631
632 ////////////////////////////////////////
633
634
635 template<typename TreeType>
636 struct StashOriginAndStoreOffset
637 {
638 using LeafNodeType = typename TreeType::LeafNodeType;
639
640
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✗ Branch 3 not taken.
 36 StashOriginAndStoreOffset(std::vector<LeafNodeType*>& nodes, Coord* coordinates)
641
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✗ Branch 3 not taken.
 36 : mNodes(nodes.empty() ? nullptr : &nodes[0]), mCoordinates(coordinates)
642 {
643 }
644
645 void operator()(const tbb::blocked_range<size_t>& range) const {
646
4/8
✓ Branch 0 taken 907 times.
✓ Branch 1 taken 374 times.
✓ Branch 2 taken 11524 times.
✓ Branch 3 taken 3369 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 16174 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
647 12431 Coord& origin = const_cast<Coord&>(mNodes[n]->origin());
648 12431 mCoordinates[n] = origin;
649 12431 origin[0] = static_cast<int>(n);
650 }
651 }
652
653 LeafNodeType ** const mNodes;
654 Coord * const mCoordinates;
655 };
656
657
658 template<typename TreeType>
659 struct RestoreOrigin
660 {
661 using LeafNodeType = typename TreeType::LeafNodeType;
662
663 36 RestoreOrigin(std::vector<LeafNodeType*>& nodes, const Coord* coordinates)
664
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✗ Branch 3 not taken.
 36 : mNodes(nodes.empty() ? nullptr : &nodes[0]), mCoordinates(coordinates)
665 {
666 }
667
668 void operator()(const tbb::blocked_range<size_t>& range) const {
669
4/8
✓ Branch 0 taken 902 times.
✓ Branch 1 taken 312 times.
✓ Branch 2 taken 11529 times.
✓ Branch 3 taken 2929 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 15672 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
670 12431 Coord& origin = const_cast<Coord&>(mNodes[n]->origin());
671 12431 origin[0] = mCoordinates[n][0];
672 }
673 }
674
675 LeafNodeType ** const mNodes;
676 Coord const * const mCoordinates;
677 };
678
679
680 template<typename TreeType>
681 class ComputeNodeConnectivity
682 {
683 public:
684 using LeafNodeType = typename TreeType::LeafNodeType;
685
686 36 ComputeNodeConnectivity(const TreeType& tree, const Coord* coordinates,
687 size_t* offsets, size_t numNodes, const CoordBBox& bbox)
688 : mTree(&tree)
689 , mCoordinates(coordinates)
690 , mOffsets(offsets)
691 , mNumNodes(numNodes)
692 36 , mBBox(bbox)
693 {
694 }
695
696 ComputeNodeConnectivity(const ComputeNodeConnectivity&) = default;
697
698 // Disallow assignment
699 ComputeNodeConnectivity& operator=(const ComputeNodeConnectivity&) = delete;
700
701 6455 void operator()(const tbb::blocked_range<size_t>& range) const {
702
703 6455 size_t* offsetsNextX = mOffsets;
704 6455 size_t* offsetsPrevX = mOffsets + mNumNodes;
705 6455 size_t* offsetsNextY = mOffsets + mNumNodes * 2;
706 6455 size_t* offsetsPrevY = mOffsets + mNumNodes * 3;
707 6455 size_t* offsetsNextZ = mOffsets + mNumNodes * 4;
708 6455 size_t* offsetsPrevZ = mOffsets + mNumNodes * 5;
709
710 6455 tree::ValueAccessor<const TreeType> acc(*mTree);
711 Coord ijk;
712 const Int32 DIM = static_cast<Int32>(LeafNodeType::DIM);
713
714
2/2
✓ Branch 0 taken 12431 times.
✓ Branch 1 taken 3360 times.
 30465 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
715 24010 const Coord& origin = mCoordinates[n];
716 24010 offsetsNextX[n] = findNeighbourNode(acc, origin, Coord(DIM, 0, 0));
717 24010 offsetsPrevX[n] = findNeighbourNode(acc, origin, Coord(-DIM, 0, 0));
718 24010 offsetsNextY[n] = findNeighbourNode(acc, origin, Coord(0, DIM, 0));
719 24010 offsetsPrevY[n] = findNeighbourNode(acc, origin, Coord(0, -DIM, 0));
720 24010 offsetsNextZ[n] = findNeighbourNode(acc, origin, Coord(0, 0, DIM));
721 24010 offsetsPrevZ[n] = findNeighbourNode(acc, origin, Coord(0, 0, -DIM));
722 }
723 6455 }
724
725 144060 size_t findNeighbourNode(tree::ValueAccessor<const TreeType>& acc,
726 const Coord& start, const Coord& step) const
727 {
728 144060 Coord ijk = start + step;
729 144060 CoordBBox bbox(mBBox);
730
731 while (bbox.isInside(ijk)) {
732 296578 const LeafNodeType* node = acc.probeConstLeaf(ijk);
733
2/2
✓ Branch 0 taken 60132 times.
✓ Branch 1 taken 90519 times.
 296578 if (node) return static_cast<size_t>(node->origin()[0]);
734 ijk += step;
735 }
736
737 return std::numeric_limits<size_t>::max();
738 }
739
740
741 private:
742 TreeType const * const mTree;
743 Coord const * const mCoordinates;
744 size_t * const mOffsets;
745
746 const size_t mNumNodes;
747 const CoordBBox mBBox;
748 }; // class ComputeNodeConnectivity
749
750
751 template<typename TreeType>
752 struct LeafNodeConnectivityTable
753 {
754 enum { INVALID_OFFSET = std::numeric_limits<size_t>::max() };
755
756 using LeafNodeType = typename TreeType::LeafNodeType;
757
758 72 LeafNodeConnectivityTable(TreeType& tree)
759
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 {
760
2/4
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 37 times.
✗ Branch 5 not taken.
 72 mLeafNodes.reserve(tree.leafCount());
761 tree.getNodes(mLeafNodes);
762
763
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 36 times.
 72 if (mLeafNodes.empty()) return;
764
765 70 CoordBBox bbox;
766
1/2
✓ Branch 1 taken 36 times.
✗ Branch 2 not taken.
 70 tree.evalLeafBoundingBox(bbox);
767
768 const tbb::blocked_range<size_t> range(0, mLeafNodes.size());
769
770 // stash the leafnode origin coordinate and temporarily store the
771 // linear offset in the origin.x variable.
772
4/6
✓ Branch 0 taken 36 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 36 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 12431 times.
✓ Branch 7 taken 36 times.
 24080 std::unique_ptr<Coord[]> coordinates{new Coord[mLeafNodes.size()]};
773
2/4
✓ Branch 1 taken 36 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36 times.
✗ Branch 4 not taken.
 70 tbb::parallel_for(range,
774 StashOriginAndStoreOffset<TreeType>(mLeafNodes, coordinates.get()));
775
776 // build the leafnode offset table
777
2/4
✓ Branch 0 taken 36 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 36 times.
✗ Branch 4 not taken.
 70 mOffsets.reset(new size_t[mLeafNodes.size() * 6]);
778
779
780
2/4
✓ Branch 1 taken 36 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36 times.
✗ Branch 4 not taken.
 70 tbb::parallel_for(range, ComputeNodeConnectivity<TreeType>(
781 tree, coordinates.get(), mOffsets.get(), mLeafNodes.size(), bbox));
782
783 // restore the leafnode origin coordinate
784
2/6
✓ Branch 1 taken 36 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 70 tbb::parallel_for(range, RestoreOrigin<TreeType>(mLeafNodes, coordinates.get()));
785 }
786
787 size_t size() const { return mLeafNodes.size(); }
788
789 std::vector<LeafNodeType*>& nodes() { return mLeafNodes; }
790 const std::vector<LeafNodeType*>& nodes() const { return mLeafNodes; }
791
792
793 const size_t* offsetsNextX() const { return mOffsets.get(); }
794 15605 const size_t* offsetsPrevX() const { return mOffsets.get() + mLeafNodes.size(); }
795
796 13735 const size_t* offsetsNextY() const { return mOffsets.get() + mLeafNodes.size() * 2; }
797 15617 const size_t* offsetsPrevY() const { return mOffsets.get() + mLeafNodes.size() * 3; }
798
799 13735 const size_t* offsetsNextZ() const { return mOffsets.get() + mLeafNodes.size() * 4; }
800 13735 const size_t* offsetsPrevZ() const { return mOffsets.get() + mLeafNodes.size() * 5; }
801
802 private:
803 std::vector<LeafNodeType*> mLeafNodes;
804 std::unique_ptr<size_t[]> mOffsets;
805 }; // struct LeafNodeConnectivityTable
806
807
808 template<typename TreeType>
809 class SweepExteriorSign
810 {
811 public:
812
813 enum Axis { X_AXIS = 0, Y_AXIS = 1, Z_AXIS = 2 };
814
815 using ValueType = typename TreeType::ValueType;
816 using LeafNodeType = typename TreeType::LeafNodeType;
817 using ConnectivityTable = LeafNodeConnectivityTable<TreeType>;
818
819
3/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✓ Branch 3 taken 1 times.
 37 SweepExteriorSign(Axis axis, const std::vector<size_t>& startNodeIndices,
820 ConnectivityTable& connectivity)
821 : mStartNodeIndices(startNodeIndices.empty() ? nullptr : &startNodeIndices[0])
822 , mConnectivity(&connectivity)
823
15/24
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 2 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 34 times.
✓ Branch 16 taken 1 times.
✓ Branch 18 taken 35 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 34 times.
✓ Branch 21 taken 1 times.
✓ Branch 23 taken 35 times.
✗ Branch 24 not taken.
✓ Branch 25 taken 34 times.
✓ Branch 26 taken 1 times.
✓ Branch 28 taken 35 times.
✗ Branch 29 not taken.
 145 , mAxis(axis)
824 {
825 }
826
827 10925 void operator()(const tbb::blocked_range<size_t>& range) const {
828
829 constexpr Int32 DIM = static_cast<Int32>(LeafNodeType::DIM);
830
831
2/2
✓ Branch 0 taken 1882 times.
✓ Branch 1 taken 3765 times.
 10925 std::vector<LeafNodeType*>& nodes = mConnectivity->nodes();
832
833 // Z Axis
834 size_t idxA = 0, idxB = 1;
835 Int32 step = 1;
836
837 const size_t* nextOffsets = mConnectivity->offsetsNextZ();
838 const size_t* prevOffsets = mConnectivity->offsetsPrevZ();
839
840
2/2
✓ Branch 0 taken 1882 times.
✓ Branch 1 taken 3765 times.
 10925 if (mAxis == Y_AXIS) {
841
842 idxA = 0;
843 idxB = 2;
844 step = DIM;
845
846 nextOffsets = mConnectivity->offsetsNextY();
847 prevOffsets = mConnectivity->offsetsPrevY();
848
849
2/2
✓ Branch 0 taken 1870 times.
✓ Branch 1 taken 1895 times.
 7310 } else if (mAxis == X_AXIS) {
850
851 idxA = 1;
852 idxB = 2;
853 step = DIM*DIM;
854
855 nextOffsets = mConnectivity->offsetsNextX();
856 prevOffsets = mConnectivity->offsetsPrevX();
857 }
858
859 Coord ijk(0, 0, 0);
860
861 Int32& a = ijk[idxA];
862 Int32& b = ijk[idxB];
863
864
2/2
✓ Branch 0 taken 7227 times.
✓ Branch 1 taken 5647 times.
 25010 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
865
866 14085 size_t startOffset = mStartNodeIndices[n];
867 size_t lastOffset = startOffset;
868
869 Int32 pos(0);
870
871
2/2
✓ Branch 0 taken 57816 times.
✓ Branch 1 taken 7227 times.
 126765 for (a = 0; a < DIM; ++a) {
872
2/2
✓ Branch 0 taken 462528 times.
✓ Branch 1 taken 57816 times.
 1014120 for (b = 0; b < DIM; ++b) {
873
874 901440 pos = static_cast<Int32>(LeafNodeType::coordToOffset(ijk));
875 size_t offset = startOffset;
876
877 // sweep in +axis direction until a boundary voxel is hit.
878
4/4
✓ Branch 0 taken 1485056 times.
✓ Branch 1 taken 117887 times.
✓ Branch 3 taken 1140415 times.
✓ Branch 4 taken 344641 times.
 3119604 while ( offset != ConnectivityTable::INVALID_OFFSET &&
879 traceVoxelLine(*nodes[offset], pos, step) ) {
880
881 lastOffset = offset;
882 2218164 offset = nextOffsets[offset];
883 }
884
885 // find last leafnode in +axis direction
886 offset = lastOffset;
887
2/2
✓ Branch 0 taken 671923 times.
✓ Branch 1 taken 462528 times.
 2206101 while (offset != ConnectivityTable::INVALID_OFFSET) {
888 lastOffset = offset;
889 1304661 offset = nextOffsets[offset];
890 }
891
892 // sweep in -axis direction until a boundary voxel is hit.
893 offset = lastOffset;
894 901440 pos += step * (DIM - 1);
895
4/4
✓ Branch 0 taken 1473001 times.
✓ Branch 1 taken 118266 times.
✓ Branch 3 taken 344262 times.
✓ Branch 4 taken 1128739 times.
 3093432 while ( offset != ConnectivityTable::INVALID_OFFSET &&
896 traceVoxelLine(*nodes[offset], pos, -step)) {
897 2191992 offset = prevOffsets[offset];
898 }
899 }
900 }
901 }
902 10925 }
903
904
905 5758662 bool traceVoxelLine(LeafNodeType& node, Int32 pos, const Int32 step) const {
906
907 5758662 ValueType* data = node.buffer().data();
908
909 bool isOutside = true;
910
911
2/2
✓ Branch 0 taken 23664456 times.
✓ Branch 1 taken 2958057 times.
 51827958 for (Index i = 0; i < LeafNodeType::DIM; ++i) {
912
913
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 23664456 times.
 46069296 assert(pos >= 0);
914 46069296 ValueType& dist = data[pos];
915
916
2/2
✓ Branch 0 taken 6484258 times.
✓ Branch 1 taken 17180198 times.
 46069296 if (dist < ValueType(0.0)) {
917 isOutside = true;
918 } else {
919 // Boundary voxel check. (Voxel that intersects the surface)
920
2/2
✓ Branch 0 taken 5428983 times.
✓ Branch 1 taken 1055275 times.
 12581139 if (!(dist > ValueType(0.75))) isOutside = false;
921
922
2/2
✓ Branch 0 taken 4270949 times.
✓ Branch 1 taken 1158034 times.
 10599018 if (isOutside) dist = ValueType(-dist);
923 }
924
925 46069296 pos += step;
926 }
927
928 5758662 return isOutside;
929 }
930
931
932 private:
933 size_t const * const mStartNodeIndices;
934 ConnectivityTable * const mConnectivity;
935
936 const Axis mAxis;
937 }; // class SweepExteriorSign
938
939
940 template<typename LeafNodeType>
941 inline void
942 seedFill(LeafNodeType& node)
943 {
944 using ValueType = typename LeafNodeType::ValueType;
945 using Queue = std::deque<Index>;
946
947
948 ValueType* data = node.buffer().data();
949
950 // find seed points
951 Queue seedPoints;
952 for (Index pos = 0; pos < LeafNodeType::SIZE; ++pos) {
953 if (data[pos] < 0.0) seedPoints.push_back(pos);
954 }
955
956 if (seedPoints.empty()) return;
957
958 // clear sign information
959 for (Queue::iterator it = seedPoints.begin(); it != seedPoints.end(); ++it) {
960 ValueType& dist = data[*it];
961 dist = -dist;
962 }
963
964 // flood fill
965
966 Coord ijk(0, 0, 0);
967 Index pos(0), nextPos(0);
968
969 while (!seedPoints.empty()) {
970
971 pos = seedPoints.back();
972 seedPoints.pop_back();
973
974 ValueType& dist = data[pos];
975
976 if (!(dist < ValueType(0.0))) {
977
978 dist = -dist; // flip sign
979
980 ijk = LeafNodeType::offsetToLocalCoord(pos);
981
982 if (ijk[0] != 0) { // i - 1, j, k
983 nextPos = pos - LeafNodeType::DIM * LeafNodeType::DIM;
984 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
985 }
986
987 if (ijk[0] != (LeafNodeType::DIM - 1)) { // i + 1, j, k
988 nextPos = pos + LeafNodeType::DIM * LeafNodeType::DIM;
989 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
990 }
991
992 if (ijk[1] != 0) { // i, j - 1, k
993 nextPos = pos - LeafNodeType::DIM;
994 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
995 }
996
997 if (ijk[1] != (LeafNodeType::DIM - 1)) { // i, j + 1, k
998 nextPos = pos + LeafNodeType::DIM;
999 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
1000 }
1001
1002 if (ijk[2] != 0) { // i, j, k - 1
1003 nextPos = pos - 1;
1004 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
1005 }
1006
1007 if (ijk[2] != (LeafNodeType::DIM - 1)) { // i, j, k + 1
1008 nextPos = pos + 1;
1009 if (data[nextPos] > ValueType(0.75)) seedPoints.push_back(nextPos);
1010 }
1011 }
1012 }
1013 } // seedFill()
1014
1015
1016 template<typename LeafNodeType>
1017 inline bool
1018 24054 scanFill(LeafNodeType& node)
1019 {
1020 bool updatedNode = false;
1021
1022 using ValueType = typename LeafNodeType::ValueType;
1023 24054 ValueType* data = node.buffer().data();
1024
1025 Coord ijk(0, 0, 0);
1026
1027 bool updatedSign = true;
1028
2/2
✓ Branch 0 taken 12667 times.
✓ Branch 1 taken 12475 times.
 48300 while (updatedSign) {
1029
1030 updatedSign = false;
1031
1032
2/2
✓ Branch 0 taken 6485504 times.
✓ Branch 1 taken 12667 times.
 12438198 for (Index pos = 0; pos < LeafNodeType::SIZE; ++pos) {
1033
1034 12413952 ValueType& dist = data[pos];
1035
1036
4/4
✓ Branch 0 taken 2004822 times.
✓ Branch 1 taken 4480682 times.
✓ Branch 2 taken 1181335 times.
✓ Branch 3 taken 718237 times.
 12413952 if (!(dist < ValueType(0.0)) && dist > ValueType(0.75)) {
1037
1038 2467920 ijk = LeafNodeType::offsetToLocalCoord(pos);
1039
1040 // i, j, k - 1
1041
4/4
✓ Branch 0 taken 1082630 times.
✓ Branch 1 taken 203955 times.
✓ Branch 2 taken 15054 times.
✓ Branch 3 taken 1067576 times.
 2467920 if (ijk[2] != 0 && data[pos - 1] < ValueType(0.0)) {
1042 updatedSign = true;
1043 15054 dist = ValueType(-dist);
1044
1045 // i, j, k + 1
1046
4/4
✓ Branch 0 taken 1070973 times.
✓ Branch 1 taken 200558 times.
✓ Branch 2 taken 1062640 times.
✓ Branch 3 taken 8333 times.
 2452866 } else if (ijk[2] != (LeafNodeType::DIM - 1) && data[pos + 1] < ValueType(0.0)) {
1047 updatedSign = true;
1048 8333 dist = ValueType(-dist);
1049
1050 // i, j - 1, k
1051
4/4
✓ Branch 0 taken 1060574 times.
✓ Branch 1 taken 202624 times.
✓ Branch 2 taken 1060070 times.
✓ Branch 3 taken 504 times.
 2444533 } else if (ijk[1] != 0 && data[pos - LeafNodeType::DIM] < ValueType(0.0)) {
1052 updatedSign = true;
1053 504 dist = ValueType(-dist);
1054
1055 // i, j + 1, k
1056 } else if (ijk[1] != (LeafNodeType::DIM - 1)
1057
4/4
✓ Branch 0 taken 1059949 times.
✓ Branch 1 taken 202745 times.
✓ Branch 2 taken 1059869 times.
✓ Branch 3 taken 80 times.
 2444029 && data[pos + LeafNodeType::DIM] < ValueType(0.0))
1058 {
1059 updatedSign = true;
1060 80 dist = ValueType(-dist);
1061
1062 // i - 1, j, k
1063 2443949 } else if (ijk[0] != 0
1064
4/4
✓ Branch 0 taken 1058932 times.
✓ Branch 1 taken 203682 times.
✓ Branch 2 taken 1058890 times.
✓ Branch 3 taken 42 times.
 2443949 && data[pos - LeafNodeType::DIM * LeafNodeType::DIM] < ValueType(0.0))
1065 {
1066 updatedSign = true;
1067 42 dist = ValueType(-dist);
1068
1069 // i + 1, j, k
1070 } else if (ijk[0] != (LeafNodeType::DIM - 1)
1071
4/4
✓ Branch 0 taken 1055653 times.
✓ Branch 1 taken 206919 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 1055643 times.
 2443907 && data[pos + LeafNodeType::DIM * LeafNodeType::DIM] < ValueType(0.0))
1072 {
1073 updatedSign = true;
1074 10 dist = ValueType(-dist);
1075 }
1076 }
1077 } // end value loop
1078
1079 24246 updatedNode |= updatedSign;
1080 } // end update loop
1081
1082 24054 return updatedNode;
1083 } // scanFill()
1084
1085
1086 template<typename TreeType>
1087 class SeedFillExteriorSign
1088 {
1089 public:
1090 using ValueType = typename TreeType::ValueType;
1091 using LeafNodeType = typename TreeType::LeafNodeType;
1092
1093
3/4
✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✓ Branch 3 taken 1 times.
 39 SeedFillExteriorSign(std::vector<LeafNodeType*>& nodes, const bool* changedNodeMask)
1094
3/4
✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✓ Branch 3 taken 1 times.
 39 : mNodes(nodes.empty() ? nullptr : &nodes[0])
1095 39 , mChangedNodeMask(changedNodeMask)
1096 {
1097 }
1098
1099 7238 void operator()(const tbb::blocked_range<size_t>& range) const {
1100
2/2
✓ Branch 0 taken 13735 times.
✓ Branch 1 taken 3929 times.
 32552 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1101
2/2
✓ Branch 0 taken 12475 times.
✓ Branch 1 taken 1260 times.
 25314 if (mChangedNodeMask[n]) {
1102 //seedFill(*mNodes[n]);
1103 // Do not update the flag in mChangedNodeMask even if scanFill
1104 // returns false. mChangedNodeMask is queried by neighboring
1105 // accesses in ::SeedPoints which needs to know that this
1106 // node has values propagated on a previous iteration.
1107 24054 scanFill(*mNodes[n]);
1108 }
1109 }
1110 7238 }
1111
1112 LeafNodeType ** const mNodes;
1113 const bool * const mChangedNodeMask;
1114 };
1115
1116
1117 template<typename ValueType>
1118 struct FillArray
1119 {
1120 111 FillArray(ValueType* array, const ValueType v) : mArray(array), mValue(v) { }
1121
1122 void operator()(const tbb::blocked_range<size_t>& range) const {
1123 142 const ValueType v = mValue;
1124
2/2
✓ Branch 0 taken 6389534 times.
✓ Branch 1 taken 142 times.
 6389676 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1125 6389534 mArray[n] = v;
1126 }
1127 }
1128
1129 ValueType * const mArray;
1130 const ValueType mValue;
1131 };
1132
1133
1134 template<typename ValueType>
1135 inline void
1136 111 fillArray(ValueType* array, const ValueType val, const size_t length)
1137 {
1138 111 const auto grainSize = std::max<size_t>(
1139
2/2
✓ Branch 1 taken 77 times.
✓ Branch 2 taken 34 times.
 188 length / tbb::this_task_arena::max_concurrency(), 1024);
1140 const tbb::blocked_range<size_t> range(0, length, grainSize);
1141 111 tbb::parallel_for(range, FillArray<ValueType>(array, val), tbb::simple_partitioner());
1142 111 }
1143
1144
1145 template<typename TreeType>
1146 class SyncVoxelMask
1147 {
1148 public:
1149 using ValueType = typename TreeType::ValueType;
1150 using LeafNodeType = typename TreeType::LeafNodeType;
1151
1152
1/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 2 SyncVoxelMask(std::vector<LeafNodeType*>& nodes,
1153 const bool* changedNodeMask, bool* changedVoxelMask)
1154
1/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 2 : mNodes(nodes.empty() ? nullptr : &nodes[0])
1155 , mChangedNodeMask(changedNodeMask)
1156 2 , mChangedVoxelMask(changedVoxelMask)
1157 {
1158 }
1159
1160 324 void operator()(const tbb::blocked_range<size_t>& range) const {
1161
2/2
✓ Branch 0 taken 1304 times.
✓ Branch 1 taken 324 times.
 1628 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1162
1163
2/2
✓ Branch 0 taken 44 times.
✓ Branch 1 taken 1260 times.
 1304 if (mChangedNodeMask[n]) {
1164 44 bool* mask = &mChangedVoxelMask[n * LeafNodeType::SIZE];
1165
1166 44 ValueType* data = mNodes[n]->buffer().data();
1167
1168
2/2
✓ Branch 0 taken 22528 times.
✓ Branch 1 taken 44 times.
 22572 for (Index pos = 0; pos < LeafNodeType::SIZE; ++pos) {
1169
2/2
✓ Branch 0 taken 3011 times.
✓ Branch 1 taken 19517 times.
 22528 if (mask[pos]) {
1170 3011 data[pos] = ValueType(-data[pos]);
1171 3011 mask[pos] = false;
1172 }
1173 }
1174 }
1175 }
1176 324 }
1177
1178 LeafNodeType ** const mNodes;
1179 bool const * const mChangedNodeMask;
1180 bool * const mChangedVoxelMask;
1181 };
1182
1183
1184 template<typename TreeType>
1185 class SeedPoints
1186 {
1187 public:
1188 using ValueType = typename TreeType::ValueType;
1189 using LeafNodeType = typename TreeType::LeafNodeType;
1190 using ConnectivityTable = LeafNodeConnectivityTable<TreeType>;
1191
1192 39 SeedPoints(ConnectivityTable& connectivity,
1193 bool* changedNodeMask, bool* nodeMask, bool* changedVoxelMask)
1194 : mConnectivity(&connectivity)
1195 , mChangedNodeMask(changedNodeMask)
1196 , mNodeMask(nodeMask)
1197 39 , mChangedVoxelMask(changedVoxelMask)
1198 {
1199 }
1200
1201 6673 void operator()(const tbb::blocked_range<size_t>& range) const {
1202
1203
2/2
✓ Branch 0 taken 13735 times.
✓ Branch 1 taken 3619 times.
 31987 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1204 bool changedValue = false;
1205
1206 25314 changedValue |= processZ(n, /*firstFace=*/true);
1207 25314 changedValue |= processZ(n, /*firstFace=*/false);
1208
1209 25314 changedValue |= processY(n, /*firstFace=*/true);
1210 25314 changedValue |= processY(n, /*firstFace=*/false);
1211
1212 25314 changedValue |= processX(n, /*firstFace=*/true);
1213 25314 changedValue |= processX(n, /*firstFace=*/false);
1214
1215 25314 mNodeMask[n] = changedValue;
1216 }
1217 6673 }
1218
1219
1220 50628 bool processZ(const size_t n, bool firstFace) const
1221 {
1222 const size_t offset =
1223
2/2
✓ Branch 0 taken 13735 times.
✓ Branch 1 taken 13735 times.
 50628 firstFace ? mConnectivity->offsetsPrevZ()[n] : mConnectivity->offsetsNextZ()[n];
1224
4/4
✓ Branch 0 taken 22260 times.
✓ Branch 1 taken 5210 times.
✓ Branch 2 taken 20196 times.
✓ Branch 3 taken 2064 times.
 50628 if (offset != ConnectivityTable::INVALID_OFFSET && mChangedNodeMask[offset]) {
1225
1226 38898 bool* mask = &mChangedVoxelMask[n * LeafNodeType::SIZE];
1227
1228 38898 const ValueType* lhsData = mConnectivity->nodes()[n]->buffer().data();
1229 38898 const ValueType* rhsData = mConnectivity->nodes()[offset]->buffer().data();
1230
1231 const Index lastOffset = LeafNodeType::DIM - 1;
1232
2/2
✓ Branch 0 taken 10098 times.
✓ Branch 1 taken 10098 times.
 38898 const Index lhsOffset =
1233
2/2
✓ Branch 0 taken 10098 times.
✓ Branch 1 taken 10098 times.
 38898 firstFace ? 0 : lastOffset, rhsOffset = firstFace ? lastOffset : 0;
1234
1235 Index tmpPos(0), pos(0);
1236 bool changedValue = false;
1237
1238
2/2
✓ Branch 0 taken 161568 times.
✓ Branch 1 taken 20196 times.
 350082 for (Index x = 0; x < LeafNodeType::DIM; ++x) {
1239 311184 tmpPos = x << (2 * LeafNodeType::LOG2DIM);
1240
2/2
✓ Branch 0 taken 1292544 times.
✓ Branch 1 taken 161568 times.
 2800656 for (Index y = 0; y < LeafNodeType::DIM; ++y) {
1241 2489472 pos = tmpPos + (y << LeafNodeType::LOG2DIM);
1242
1243
2/2
✓ Branch 0 taken 400031 times.
✓ Branch 1 taken 892513 times.
 2489472 if (lhsData[pos + lhsOffset] > ValueType(0.75)) {
1244
2/2
✓ Branch 0 taken 2627 times.
✓ Branch 1 taken 397404 times.
 786492 if (rhsData[pos + rhsOffset] < ValueType(0.0)) {
1245 changedValue = true;
1246 2627 mask[pos + lhsOffset] = true;
1247 }
1248 }
1249 }
1250 }
1251
1252 return changedValue;
1253 }
1254
1255 return false;
1256 }
1257
1258 50628 bool processY(const size_t n, bool firstFace) const
1259 {
1260 const size_t offset =
1261
2/2
✓ Branch 0 taken 13735 times.
✓ Branch 1 taken 13735 times.
 50628 firstFace ? mConnectivity->offsetsPrevY()[n] : mConnectivity->offsetsNextY()[n];
1262
4/4
✓ Branch 0 taken 22190 times.
✓ Branch 1 taken 5280 times.
✓ Branch 2 taken 20126 times.
✓ Branch 3 taken 2064 times.
 50628 if (offset != ConnectivityTable::INVALID_OFFSET && mChangedNodeMask[offset]) {
1263
1264 38758 bool* mask = &mChangedVoxelMask[n * LeafNodeType::SIZE];
1265
1266 38758 const ValueType* lhsData = mConnectivity->nodes()[n]->buffer().data();
1267 38758 const ValueType* rhsData = mConnectivity->nodes()[offset]->buffer().data();
1268
1269 const Index lastOffset = LeafNodeType::DIM * (LeafNodeType::DIM - 1);
1270
2/2
✓ Branch 0 taken 10063 times.
✓ Branch 1 taken 10063 times.
 38758 const Index lhsOffset =
1271
2/2
✓ Branch 0 taken 10063 times.
✓ Branch 1 taken 10063 times.
 38758 firstFace ? 0 : lastOffset, rhsOffset = firstFace ? lastOffset : 0;
1272
1273 Index tmpPos(0), pos(0);
1274 bool changedValue = false;
1275
1276
2/2
✓ Branch 0 taken 161008 times.
✓ Branch 1 taken 20126 times.
 348822 for (Index x = 0; x < LeafNodeType::DIM; ++x) {
1277 310064 tmpPos = x << (2 * LeafNodeType::LOG2DIM);
1278
2/2
✓ Branch 0 taken 1288064 times.
✓ Branch 1 taken 161008 times.
 2790576 for (Index z = 0; z < LeafNodeType::DIM; ++z) {
1279 2480512 pos = tmpPos + z;
1280
1281
2/2
✓ Branch 0 taken 400473 times.
✓ Branch 1 taken 887591 times.
 2480512 if (lhsData[pos + lhsOffset] > ValueType(0.75)) {
1282
2/2
✓ Branch 0 taken 224 times.
✓ Branch 1 taken 400249 times.
 787671 if (rhsData[pos + rhsOffset] < ValueType(0.0)) {
1283 changedValue = true;
1284 224 mask[pos + lhsOffset] = true;
1285 }
1286 }
1287 }
1288 }
1289
1290 return changedValue;
1291 }
1292
1293 return false;
1294 }
1295
1296 50628 bool processX(const size_t n, bool firstFace) const
1297 {
1298 const size_t offset =
1299
2/2
✓ Branch 0 taken 13735 times.
✓ Branch 1 taken 13735 times.
 50628 firstFace ? mConnectivity->offsetsPrevX()[n] : mConnectivity->offsetsNextX()[n];
1300
4/4
✓ Branch 0 taken 22450 times.
✓ Branch 1 taken 5020 times.
✓ Branch 2 taken 20074 times.
✓ Branch 3 taken 2376 times.
 50628 if (offset != ConnectivityTable::INVALID_OFFSET && mChangedNodeMask[offset]) {
1301
1302 38498 bool* mask = &mChangedVoxelMask[n * LeafNodeType::SIZE];
1303
1304 38498 const ValueType* lhsData = mConnectivity->nodes()[n]->buffer().data();
1305 38498 const ValueType* rhsData = mConnectivity->nodes()[offset]->buffer().data();
1306
1307 const Index lastOffset = LeafNodeType::DIM * LeafNodeType::DIM * (LeafNodeType::DIM-1);
1308
2/2
✓ Branch 0 taken 10037 times.
✓ Branch 1 taken 10037 times.
 38498 const Index lhsOffset =
1309
2/2
✓ Branch 0 taken 10037 times.
✓ Branch 1 taken 10037 times.
 38498 firstFace ? 0 : lastOffset, rhsOffset = firstFace ? lastOffset : 0;
1310
1311 Index tmpPos(0), pos(0);
1312 bool changedValue = false;
1313
1314
2/2
✓ Branch 0 taken 160592 times.
✓ Branch 1 taken 20074 times.
 346482 for (Index y = 0; y < LeafNodeType::DIM; ++y) {
1315 307984 tmpPos = y << LeafNodeType::LOG2DIM;
1316
2/2
✓ Branch 0 taken 1284736 times.
✓ Branch 1 taken 160592 times.
 2771856 for (Index z = 0; z < LeafNodeType::DIM; ++z) {
1317 2463872 pos = tmpPos + z;
1318
1319
2/2
✓ Branch 0 taken 405757 times.
✓ Branch 1 taken 878979 times.
 2463872 if (lhsData[pos + lhsOffset] > ValueType(0.75)) {
1320
2/2
✓ Branch 0 taken 224 times.
✓ Branch 1 taken 405533 times.
 797278 if (rhsData[pos + rhsOffset] < ValueType(0.0)) {
1321 changedValue = true;
1322 224 mask[pos + lhsOffset] = true;
1323 }
1324 }
1325 }
1326 }
1327
1328 return changedValue;
1329 }
1330
1331 return false;
1332 }
1333
1334 ConnectivityTable * const mConnectivity;
1335 bool * const mChangedNodeMask;
1336 bool * const mNodeMask;
1337 bool * const mChangedVoxelMask;
1338 };
1339
1340
1341 ////////////////////////////////////////
1342
1343 template<typename TreeType, typename MeshDataAdapter>
1344 747 struct ComputeIntersectingVoxelSign
1345 {
1346 using ValueType = typename TreeType::ValueType;
1347 using LeafNodeType = typename TreeType::LeafNodeType;
1348 using Int32TreeType = typename TreeType::template ValueConverter<Int32>::Type;
1349 using Int32LeafNodeType = typename Int32TreeType::LeafNodeType;
1350
1351 using PointArray = std::unique_ptr<Vec3d[]>;
1352 using MaskArray = std::unique_ptr<bool[]>;
1353 using LocalData = std::pair<PointArray, MaskArray>;
1354 using LocalDataTable = tbb::enumerable_thread_specific<LocalData>;
1355
1356
2/2
✓ Branch 0 taken 34 times.
✓ Branch 1 taken 1 times.
 70 ComputeIntersectingVoxelSign(
1357 std::vector<LeafNodeType*>& distNodes,
1358 const TreeType& distTree,
1359 const Int32TreeType& indexTree,
1360 const MeshDataAdapter& mesh)
1361
2/2
✓ Branch 0 taken 34 times.
✓ Branch 1 taken 1 times.
 70 : mDistNodes(distNodes.empty() ? nullptr : &distNodes[0])
1362 , mDistTree(&distTree)
1363 , mIndexTree(&indexTree)
1364 , mMesh(&mesh)
1365 70 , mLocalDataTable(new LocalDataTable())
1366 {
1367 }
1368
1369
1370 6088 void operator()(const tbb::blocked_range<size_t>& range) const {
1371
1372 6088 tree::ValueAccessor<const TreeType> distAcc(*mDistTree);
1373
1/2
✓ Branch 1 taken 3044 times.
✗ Branch 2 not taken.
 6088 tree::ValueAccessor<const Int32TreeType> idxAcc(*mIndexTree);
1374
1375 ValueType nval;
1376 6088 CoordBBox bbox;
1377 Index xPos(0), yPos(0);
1378 Coord ijk, nijk, nodeMin, nodeMax;
1379 Vec3d cp, xyz, nxyz, dir1, dir2;
1380
1381 LocalData& localData = mLocalDataTable->local();
1382
1383 PointArray& points = localData.first;
1384
5/6
✓ Branch 0 taken 65 times.
✓ Branch 1 taken 2979 times.
✓ Branch 3 taken 65 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 66560 times.
✓ Branch 6 taken 65 times.
 139338 if (!points) points.reset(new Vec3d[LeafNodeType::SIZE * 2]);
1385
1386 MaskArray& mask = localData.second;
1387
3/4
✓ Branch 0 taken 65 times.
✓ Branch 1 taken 2979 times.
✓ Branch 3 taken 65 times.
✗ Branch 4 not taken.
 6088 if (!mask) mask.reset(new bool[LeafNodeType::SIZE]);
1388
1389
1390 typename LeafNodeType::ValueOnCIter it;
1391
1392
2/2
✓ Branch 0 taken 11579 times.
✓ Branch 1 taken 3044 times.
 29246 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1393
1394
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 LeafNodeType& node = *mDistNodes[n];
1395
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 ValueType* data = node.buffer().data();
1396
1397 const Int32LeafNodeType* idxNode = idxAcc.probeConstLeaf(node.origin());
1398
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 const Int32* idxData = idxNode->buffer().data();
1399
1400 23158 nodeMin = node.origin();
1401 23158 nodeMax = nodeMin.offsetBy(LeafNodeType::DIM - 1);
1402
1403 // reset computed voxel mask.
1404 23158 memset(mask.get(), 0, sizeof(bool) * LeafNodeType::SIZE);
1405
1406
2/2
✓ Branch 0 taken 1997311 times.
✓ Branch 1 taken 11579 times.
 4040938 for (it = node.cbeginValueOn(); it; ++it) {
1407 Index pos = it.pos();
1408
1409 3994622 ValueType& dist = data[pos];
1410
4/4
✓ Branch 0 taken 1315320 times.
✓ Branch 1 taken 681991 times.
✓ Branch 2 taken 597083 times.
✓ Branch 3 taken 718237 times.
 3994622 if (dist < 0.0 || dist > 0.75) continue;
1411
1412 1436474 ijk = node.offsetToGlobalCoord(pos);
1413
1414 1436474 xyz[0] = double(ijk[0]);
1415 1436474 xyz[1] = double(ijk[1]);
1416 1436474 xyz[2] = double(ijk[2]);
1417
1418
1419 1436474 bbox.min() = Coord::maxComponent(ijk.offsetBy(-1), nodeMin);
1420 1436474 bbox.max() = Coord::minComponent(ijk.offsetBy(1), nodeMax);
1421
1422 bool flipSign = false;
1423
1424
4/4
✓ Branch 0 taken 1889680 times.
✓ Branch 1 taken 534930 times.
✓ Branch 2 taken 1706373 times.
✓ Branch 3 taken 183307 times.
 4849220 for (nijk[0] = bbox.min()[0]; nijk[0] <= bbox.max()[0] && !flipSign; ++nijk[0]) {
1425 3412746 xPos = (nijk[0] & (LeafNodeType::DIM - 1u)) << (2 * LeafNodeType::LOG2DIM);
1426
4/4
✓ Branch 0 taken 4645561 times.
✓ Branch 1 taken 1528774 times.
✓ Branch 2 taken 177599 times.
✓ Branch 3 taken 4467962 times.
 12348670 for (nijk[1]=bbox.min()[1]; nijk[1] <= bbox.max()[1] && !flipSign; ++nijk[1]) {
1427 8935924 yPos = xPos + ((nijk[1] & (LeafNodeType::DIM-1u)) << LeafNodeType::LOG2DIM);
1428
2/2
✓ Branch 0 taken 12168439 times.
✓ Branch 1 taken 4227185 times.
 32791248 for (nijk[2] = bbox.min()[2]; nijk[2] <= bbox.max()[2]; ++nijk[2]) {
1429 24336878 pos = yPos + (nijk[2] & (LeafNodeType::DIM - 1u));
1430
1431 24336878 const Int32& polyIdx = idxData[pos];
1432
1433
3/4
✓ Branch 0 taken 12168439 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 8917530 times.
✓ Branch 3 taken 3250909 times.
 24336878 if (polyIdx == Int32(util::INVALID_IDX) || !(data[pos] < -0.75))
1434 17835060 continue;
1435
1436
2/2
✓ Branch 0 taken 580971 times.
✓ Branch 1 taken 2669938 times.
 6501818 const Index pointIndex = pos * 2;
1437
1438
2/2
✓ Branch 0 taken 580971 times.
✓ Branch 1 taken 2669938 times.
 6501818 if (!mask[pos]) {
1439
1440 1161942 mask[pos] = true;
1441
1442 1161942 nxyz[0] = double(nijk[0]);
1443 1161942 nxyz[1] = double(nijk[1]);
1444 1161942 nxyz[2] = double(nijk[2]);
1445
1446
1/2
✓ Branch 1 taken 580971 times.
✗ Branch 2 not taken.
 1161942 Vec3d& point = points[pointIndex];
1447
1448
1/2
✓ Branch 1 taken 580971 times.
✗ Branch 2 not taken.
 1161942 point = closestPoint(nxyz, polyIdx);
1449
1450 1161942 Vec3d& direction = points[pointIndex + 1];
1451 1161942 direction = nxyz - point;
1452 1161942 direction.normalize();
1453 }
1454
1455 6501818 dir1 = xyz - points[pointIndex];
1456 6501818 dir1.normalize();
1457
1458
2/2
✓ Branch 0 taken 3010132 times.
✓ Branch 1 taken 240777 times.
 6501818 if (points[pointIndex + 1].dot(dir1) > 0.0) {
1459 flipSign = true;
1460 break;
1461 }
1462 }
1463 }
1464 }
1465
1466
2/2
✓ Branch 0 taken 477460 times.
✓ Branch 1 taken 240777 times.
 1436474 if (flipSign) {
1467 481554 dist = -dist;
1468 } else {
1469
2/2
✓ Branch 0 taken 11977175 times.
✓ Branch 1 taken 450153 times.
 24854656 for (Int32 m = 0; m < 26; ++m) {
1470
2/2
✓ Branch 0 taken 11554596 times.
✓ Branch 1 taken 422579 times.
 23954350 nijk = ijk + util::COORD_OFFSETS[m];
1471
1472
4/6
✓ Branch 1 taken 2305505 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2305505 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 637466 times.
✓ Branch 6 taken 1668039 times.
 4611010 if (!bbox.isInside(nijk) && distAcc.probeValue(nijk, nval) && nval<-0.75) {
1473 1274932 nxyz[0] = double(nijk[0]);
1474 1274932 nxyz[1] = double(nijk[1]);
1475 1274932 nxyz[2] = double(nijk[2]);
1476
1477
2/4
✓ Branch 1 taken 637466 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 637466 times.
✗ Branch 5 not taken.
 1274932 cp = closestPoint(nxyz, idxAcc.getValue(nijk));
1478
1479 1274932 dir1 = xyz - cp;
1480 1274932 dir1.normalize();
1481
1482 1274932 dir2 = nxyz - cp;
1483 1274932 dir2.normalize();
1484
1485
2/2
✓ Branch 0 taken 27307 times.
✓ Branch 1 taken 610159 times.
 1274932 if (dir2.dot(dir1) > 0.0) {
1486 54614 dist = -dist;
1487 54614 break;
1488 }
1489 }
1490 }
1491 }
1492
1493 } // active voxel loop
1494 } // leaf node loop
1495 }
1496
1497 private:
1498
1499 2436874 Vec3d closestPoint(const Vec3d& center, Int32 polyIdx) const
1500 {
1501 Vec3d a, b, c, cp, uvw;
1502
1503 2436874 const size_t polygon = size_t(polyIdx);
1504 2436874 mMesh->getIndexSpacePoint(polygon, 0, a);
1505 2436874 mMesh->getIndexSpacePoint(polygon, 1, b);
1506 2436874 mMesh->getIndexSpacePoint(polygon, 2, c);
1507
1508 2436874 cp = closestPointOnTriangleToPoint(a, c, b, center, uvw);
1509
1510
2/2
✓ Branch 0 taken 1198018 times.
✓ Branch 1 taken 10305 times.
 2416646 if (4 == mMesh->vertexCount(polygon)) {
1511
1512 2396036 mMesh->getIndexSpacePoint(polygon, 3, b);
1513
1514 2396036 c = closestPointOnTriangleToPoint(a, b, c, center, uvw);
1515
1516
2/2
✓ Branch 0 taken 601415 times.
✓ Branch 1 taken 596603 times.
 2396036 if ((center - c).lengthSqr() < (center - cp).lengthSqr()) {
1517 1202830 cp = c;
1518 }
1519 }
1520
1521 2436874 return cp;
1522 }
1523
1524
1525 LeafNodeType ** const mDistNodes;
1526 TreeType const * const mDistTree;
1527 Int32TreeType const * const mIndexTree;
1528 MeshDataAdapter const * const mMesh;
1529
1530 SharedPtr<LocalDataTable> mLocalDataTable;
1531 }; // ComputeIntersectingVoxelSign
1532
1533
1534 ////////////////////////////////////////
1535
1536
1537 template<typename LeafNodeType>
1538 inline void
1539 2399168 maskNodeInternalNeighbours(const Index pos, bool (&mask)[26])
1540 {
1541 using NodeT = LeafNodeType;
1542
1543 2399168 const Coord ijk = NodeT::offsetToLocalCoord(pos);
1544
1545 // Face adjacent neighbours
1546 // i+1, j, k
1547 2399168 mask[0] = ijk[0] != (NodeT::DIM - 1);
1548 // i-1, j, k
1549 2399168 mask[1] = ijk[0] != 0;
1550 // i, j+1, k
1551 2399168 mask[2] = ijk[1] != (NodeT::DIM - 1);
1552 // i, j-1, k
1553 2399168 mask[3] = ijk[1] != 0;
1554 // i, j, k+1
1555 2399168 mask[4] = ijk[2] != (NodeT::DIM - 1);
1556 // i, j, k-1
1557 2399168 mask[5] = ijk[2] != 0;
1558
1559 // Edge adjacent neighbour
1560 // i+1, j, k-1
1561
4/4
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 92015 times.
✓ Branch 3 taken 922290 times.
 2399168 mask[6] = mask[0] && mask[5];
1562 // i-1, j, k-1
1563
4/4
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 101946 times.
✓ Branch 3 taken 986884 times.
 2399168 mask[7] = mask[1] && mask[5];
1564 // i+1, j, k+1
1565
4/4
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 160613 times.
✓ Branch 3 taken 853692 times.
 2399168 mask[8] = mask[0] && mask[4];
1566 // i-1, j, k+1
1567
4/4
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 169909 times.
✓ Branch 3 taken 918921 times.
 2399168 mask[9] = mask[1] && mask[4];
1568 // i+1, j+1, k
1569
4/4
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 162931 times.
✓ Branch 3 taken 851374 times.
 2399168 mask[10] = mask[0] && mask[2];
1570 // i-1, j+1, k
1571
4/4
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 172260 times.
✓ Branch 3 taken 916570 times.
 2399168 mask[11] = mask[1] && mask[2];
1572 // i+1, j-1, k
1573
4/4
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 93692 times.
✓ Branch 3 taken 920613 times.
 2399168 mask[12] = mask[0] && mask[3];
1574 // i-1, j-1, k
1575
4/4
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 103580 times.
✓ Branch 3 taken 985250 times.
 2399168 mask[13] = mask[1] && mask[3];
1576 // i, j-1, k+1
1577
4/4
✓ Branch 0 taken 1087256 times.
✓ Branch 1 taken 112328 times.
✓ Branch 2 taken 169779 times.
✓ Branch 3 taken 917477 times.
 2399168 mask[14] = mask[3] && mask[4];
1578 // i, j-1, k-1
1579
4/4
✓ Branch 0 taken 1087256 times.
✓ Branch 1 taken 112328 times.
✓ Branch 2 taken 101746 times.
✓ Branch 3 taken 985510 times.
 2399168 mask[15] = mask[3] && mask[5];
1580 // i, j+1, k+1
1581
4/4
✓ Branch 0 taken 1008105 times.
✓ Branch 1 taken 191479 times.
✓ Branch 2 taken 159925 times.
✓ Branch 3 taken 848180 times.
 2399168 mask[16] = mask[2] && mask[4];
1582 // i, j+1, k-1
1583
4/4
✓ Branch 0 taken 1008105 times.
✓ Branch 1 taken 191479 times.
✓ Branch 2 taken 91386 times.
✓ Branch 3 taken 916719 times.
 2399168 mask[17] = mask[2] && mask[5];
1584
1585 // Corner adjacent neighbours
1586 // i-1, j-1, k-1
1587
6/6
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 985250 times.
✓ Branch 3 taken 103580 times.
✓ Branch 4 taken 93636 times.
✓ Branch 5 taken 891614 times.
 2399168 mask[18] = mask[1] && mask[3] && mask[5];
1588 // i-1, j-1, k+1
1589
6/6
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 985250 times.
✓ Branch 3 taken 103580 times.
✓ Branch 4 taken 152627 times.
✓ Branch 5 taken 832623 times.
 2399168 mask[19] = mask[1] && mask[3] && mask[4];
1590 // i+1, j-1, k+1
1591
6/6
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 920613 times.
✓ Branch 3 taken 93692 times.
✓ Branch 4 taken 144566 times.
✓ Branch 5 taken 776047 times.
 2399168 mask[20] = mask[0] && mask[3] && mask[4];
1592 // i+1, j-1, k-1
1593
6/6
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 920613 times.
✓ Branch 3 taken 93692 times.
✓ Branch 4 taken 85028 times.
✓ Branch 5 taken 835585 times.
 2399168 mask[21] = mask[0] && mask[3] && mask[5];
1594 // i-1, j+1, k-1
1595
6/6
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 916570 times.
✓ Branch 3 taken 172260 times.
✓ Branch 4 taken 84653 times.
✓ Branch 5 taken 831917 times.
 2399168 mask[22] = mask[1] && mask[2] && mask[5];
1596 // i-1, j+1, k+1
1597
6/6
✓ Branch 0 taken 1088830 times.
✓ Branch 1 taken 110754 times.
✓ Branch 2 taken 916570 times.
✓ Branch 3 taken 172260 times.
✓ Branch 4 taken 144077 times.
✓ Branch 5 taken 772493 times.
 2399168 mask[23] = mask[1] && mask[2] && mask[4];
1598 // i+1, j+1, k+1
1599
6/6
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 851374 times.
✓ Branch 3 taken 162931 times.
✓ Branch 4 taken 135940 times.
✓ Branch 5 taken 715434 times.
 2399168 mask[24] = mask[0] && mask[2] && mask[4];
1600 // i+1, j+1, k-1
1601
6/6
✓ Branch 0 taken 1014305 times.
✓ Branch 1 taken 185279 times.
✓ Branch 2 taken 851374 times.
✓ Branch 3 taken 162931 times.
✓ Branch 4 taken 75977 times.
✓ Branch 5 taken 775397 times.
 2399168 mask[25] = mask[0] && mask[2] && mask[5];
1602 }
1603
1604
1605 template<typename Compare, typename LeafNodeType>
1606 inline bool
1607 2399168 checkNeighbours(const Index pos, const typename LeafNodeType::ValueType * data, bool (&mask)[26])
1608 {
1609 using NodeT = LeafNodeType;
1610
1611 // i, j, k - 1
1612
4/4
✓ Branch 0 taken 1089179 times.
✓ Branch 1 taken 110405 times.
✓ Branch 2 taken 886169 times.
✓ Branch 3 taken 203010 times.
 2399168 if (mask[5] && Compare::check(data[pos - 1])) return true;
1613 // i, j, k + 1
1614
4/4
✓ Branch 0 taken 888777 times.
✓ Branch 1 taken 107797 times.
✓ Branch 2 taken 689706 times.
✓ Branch 3 taken 199071 times.
 1993148 if (mask[4] && Compare::check(data[pos + 1])) return true;
1615 // i, j - 1, k
1616
4/4
✓ Branch 0 taken 736397 times.
✓ Branch 1 taken 61106 times.
✓ Branch 2 taken 552492 times.
✓ Branch 3 taken 183905 times.
 1595006 if (mask[3] && Compare::check(data[pos - NodeT::DIM])) return true;
1617 // i, j + 1, k
1618
4/4
✓ Branch 0 taken 554307 times.
✓ Branch 1 taken 59291 times.
✓ Branch 2 taken 371895 times.
✓ Branch 3 taken 182412 times.
 1227196 if (mask[2] && Compare::check(data[pos + NodeT::DIM])) return true;
1619 // i - 1, j, k
1620
4/4
✓ Branch 0 taken 419361 times.
✓ Branch 1 taken 11825 times.
✓ Branch 2 taken 230501 times.
✓ Branch 3 taken 188860 times.
 862372 if (mask[1] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM])) return true;
1621 // i + 1, j, k
1622
4/4
✓ Branch 0 taken 232018 times.
✓ Branch 1 taken 10308 times.
✓ Branch 2 taken 43756 times.
✓ Branch 3 taken 188262 times.
 484652 if (mask[0] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM])) return true;
1623 // i+1, j, k-1
1624
3/4
✓ Branch 0 taken 37305 times.
✓ Branch 1 taken 16759 times.
✓ Branch 2 taken 37305 times.
✗ Branch 3 not taken.
 108128 if (mask[6] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM])) return true;
1625 // i-1, j, k-1
1626
4/4
✓ Branch 0 taken 38212 times.
✓ Branch 1 taken 15852 times.
✓ Branch 2 taken 35514 times.
✓ Branch 3 taken 2698 times.
 108128 if (mask[7] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM - 1])) return true;
1627 // i+1, j, k+1
1628
4/4
✓ Branch 0 taken 32240 times.
✓ Branch 1 taken 19126 times.
✓ Branch 2 taken 29557 times.
✓ Branch 3 taken 2683 times.
 102732 if (mask[8] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM + 1])) return true;
1629 // i-1, j, k+1
1630
4/4
✓ Branch 0 taken 30402 times.
✓ Branch 1 taken 18281 times.
✓ Branch 2 taken 27680 times.
✓ Branch 3 taken 2722 times.
 97366 if (mask[9] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM + 1])) return true;
1631 // i+1, j+1, k
1632
4/4
✓ Branch 0 taken 28266 times.
✓ Branch 1 taken 17695 times.
✓ Branch 2 taken 26475 times.
✓ Branch 3 taken 1791 times.
 91922 if (mask[10] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM + NodeT::DIM])) return true;
1633 // i-1, j+1, k
1634
4/4
✓ Branch 0 taken 27121 times.
✓ Branch 1 taken 17049 times.
✓ Branch 2 taken 26170 times.
✓ Branch 3 taken 951 times.
 88340 if (mask[11] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM + NodeT::DIM])) return true;
1635 // i+1, j-1, k
1636
4/4
✓ Branch 0 taken 26847 times.
✓ Branch 1 taken 16372 times.
✓ Branch 2 taken 25044 times.
✓ Branch 3 taken 1803 times.
 86438 if (mask[12] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM - NodeT::DIM])) return true;
1637 // i-1, j-1, k
1638
4/4
✓ Branch 0 taken 26055 times.
✓ Branch 1 taken 15361 times.
✓ Branch 2 taken 25179 times.
✓ Branch 3 taken 876 times.
 82832 if (mask[13] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM - NodeT::DIM])) return true;
1639 // i, j-1, k+1
1640
4/4
✓ Branch 0 taken 23803 times.
✓ Branch 1 taken 16737 times.
✓ Branch 2 taken 22510 times.
✓ Branch 3 taken 1293 times.
 81080 if (mask[14] && Compare::check(data[pos - NodeT::DIM + 1])) return true;
1641 // i, j-1, k-1
1642
4/4
✓ Branch 0 taken 26105 times.
✓ Branch 1 taken 13142 times.
✓ Branch 2 taken 25499 times.
✓ Branch 3 taken 606 times.
 78494 if (mask[15] && Compare::check(data[pos - NodeT::DIM - 1])) return true;
1643 // i, j+1, k+1
1644
4/4
✓ Branch 0 taken 21512 times.
✓ Branch 1 taken 17129 times.
✓ Branch 2 taken 20959 times.
✓ Branch 3 taken 553 times.
 77282 if (mask[16] && Compare::check(data[pos + NodeT::DIM + 1])) return true;
1645 // i, j+1, k-1
1646
4/4
✓ Branch 0 taken 24294 times.
✓ Branch 1 taken 13794 times.
✓ Branch 2 taken 23123 times.
✓ Branch 3 taken 1171 times.
 76176 if (mask[17] && Compare::check(data[pos + NodeT::DIM - 1])) return true;
1647 // i-1, j-1, k-1
1648
4/4
✓ Branch 0 taken 19162 times.
✓ Branch 1 taken 17755 times.
✓ Branch 2 taken 18295 times.
✓ Branch 3 taken 867 times.
 73834 if (mask[18] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM - NodeT::DIM - 1])) return true;
1649 // i-1, j-1, k+1
1650
4/4
✓ Branch 0 taken 15343 times.
✓ Branch 1 taken 20707 times.
✓ Branch 2 taken 14728 times.
✓ Branch 3 taken 615 times.
 72100 if (mask[19] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM - NodeT::DIM + 1])) return true;
1651 // i+1, j-1, k+1
1652
4/4
✓ Branch 0 taken 14152 times.
✓ Branch 1 taken 21283 times.
✓ Branch 2 taken 13515 times.
✓ Branch 3 taken 637 times.
 70870 if (mask[20] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM - NodeT::DIM + 1])) return true;
1653 // i+1, j-1, k-1
1654
4/4
✓ Branch 0 taken 16539 times.
✓ Branch 1 taken 18259 times.
✓ Branch 2 taken 15579 times.
✓ Branch 3 taken 960 times.
 69596 if (mask[21] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM - NodeT::DIM - 1])) return true;
1655 // i-1, j+1, k-1
1656
4/4
✓ Branch 0 taken 15931 times.
✓ Branch 1 taken 17907 times.
✓ Branch 2 taken 15248 times.
✓ Branch 3 taken 683 times.
 67676 if (mask[22] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM + NodeT::DIM - 1])) return true;
1657 // i-1, j+1, k+1
1658
4/4
✓ Branch 0 taken 12503 times.
✓ Branch 1 taken 20652 times.
✓ Branch 2 taken 11548 times.
✓ Branch 3 taken 955 times.
 66310 if (mask[23] && Compare::check(data[pos - NodeT::DIM * NodeT::DIM + NodeT::DIM + 1])) return true;
1659 // i+1, j+1, k+1
1660
4/4
✓ Branch 0 taken 11218 times.
✓ Branch 1 taken 20982 times.
✓ Branch 2 taken 10304 times.
✓ Branch 3 taken 914 times.
 64400 if (mask[24] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM + NodeT::DIM + 1])) return true;
1661 // i+1, j+1, k-1
1662
4/4
✓ Branch 0 taken 13163 times.
✓ Branch 1 taken 18123 times.
✓ Branch 2 taken 897 times.
✓ Branch 3 taken 12266 times.
 62572 if (mask[25] && Compare::check(data[pos + NodeT::DIM * NodeT::DIM + NodeT::DIM - 1])) return true;
1663
1664 return false;
1665 }
1666
1667
1668 template<typename Compare, typename AccessorType>
1669 inline bool
1670 60778 checkNeighbours(const Coord& ijk, AccessorType& acc, bool (&mask)[26])
1671 {
1672
2/2
✓ Branch 0 taken 177107 times.
✓ Branch 1 taken 127 times.
 354468 for (Int32 m = 0; m < 26; ++m) {
1673
4/4
✓ Branch 0 taken 67828 times.
✓ Branch 1 taken 109279 times.
✓ Branch 3 taken 30262 times.
✓ Branch 4 taken 37566 times.
 354214 if (!mask[m] && Compare::check(acc.getValue(ijk + util::COORD_OFFSETS[m]))) {
1674 60524 return true;
1675 }
1676 }
1677
1678 return false;
1679 }
1680
1681
1682 template<typename TreeType>
1683 struct ValidateIntersectingVoxels
1684 {
1685 using ValueType = typename TreeType::ValueType;
1686 using LeafNodeType = typename TreeType::LeafNodeType;
1687
1688 struct IsNegative { static bool check(const ValueType v) { return v < ValueType(0.0); } };
1689
1690 35 ValidateIntersectingVoxels(TreeType& tree, std::vector<LeafNodeType*>& nodes)
1691 : mTree(&tree)
1692
3/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 31 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
 35 , mNodes(nodes.empty() ? nullptr : &nodes[0])
1693 {
1694 }
1695
1696 6154 void operator()(const tbb::blocked_range<size_t>& range) const
1697 {
1698 6154 tree::ValueAccessor<const TreeType> acc(*mTree);
1699 bool neighbourMask[26];
1700
1701
2/2
✓ Branch 0 taken 11579 times.
✓ Branch 1 taken 3077 times.
 29312 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1702
1703
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 LeafNodeType& node = *mNodes[n];
1704
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 ValueType* data = node.buffer().data();
1705
1706 typename LeafNodeType::ValueOnCIter it;
1707
2/2
✓ Branch 0 taken 1997311 times.
✓ Branch 1 taken 11579 times.
 4040938 for (it = node.cbeginValueOn(); it; ++it) {
1708
1709 const Index pos = it.pos();
1710
1711 3994622 ValueType& dist = data[pos];
1712
4/4
✓ Branch 0 taken 1199584 times.
✓ Branch 1 taken 797727 times.
✓ Branch 2 taken 597083 times.
✓ Branch 3 taken 602501 times.
 3994622 if (dist < 0.0 || dist > 0.75) continue;
1713
1714 // Mask node internal neighbours
1715 1205002 maskNodeInternalNeighbours<LeafNodeType>(pos, neighbourMask);
1716
1717 const bool hasNegativeNeighbour =
1718
3/4
✓ Branch 0 taken 20660 times.
✓ Branch 1 taken 581841 times.
✓ Branch 2 taken 20660 times.
✗ Branch 3 not taken.
 1246322 checkNeighbours<IsNegative, LeafNodeType>(pos, data, neighbourMask) ||
1719
1/4
✓ Branch 2 taken 20660 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 1246322 checkNeighbours<IsNegative>(node.offsetToGlobalCoord(pos), acc, neighbourMask);
1720
1721 if (!hasNegativeNeighbour) {
1722 // push over boundary voxel distance
1723 ✗ dist = ValueType(0.75) + Tolerance<ValueType>::epsilon();
1724 }
1725 }
1726 }
1727 }
1728
1729 TreeType * const mTree;
1730 LeafNodeType ** const mNodes;
1731 }; // ValidateIntersectingVoxels
1732
1733
1734 template<typename TreeType>
1735 struct RemoveSelfIntersectingSurface
1736 {
1737 using ValueType = typename TreeType::ValueType;
1738 using LeafNodeType = typename TreeType::LeafNodeType;
1739 using Int32TreeType = typename TreeType::template ValueConverter<Int32>::Type;
1740
1741 struct Comp { static bool check(const ValueType v) { return !(v > ValueType(0.75)); } };
1742
1743 35 RemoveSelfIntersectingSurface(std::vector<LeafNodeType*>& nodes,
1744 TreeType& distTree, Int32TreeType& indexTree)
1745
3/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 31 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
 35 : mNodes(nodes.empty() ? nullptr : &nodes[0])
1746 , mDistTree(&distTree)
1747 35 , mIndexTree(&indexTree)
1748 {
1749 }
1750
1751 6278 void operator()(const tbb::blocked_range<size_t>& range) const
1752 {
1753 6278 tree::ValueAccessor<const TreeType> distAcc(*mDistTree);
1754
1/2
✓ Branch 1 taken 3139 times.
✗ Branch 2 not taken.
 6278 tree::ValueAccessor<Int32TreeType> idxAcc(*mIndexTree);
1755 bool neighbourMask[26];
1756
1757
2/2
✓ Branch 0 taken 11579 times.
✓ Branch 1 taken 3139 times.
 29436 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1758
1759
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 LeafNodeType& distNode = *mNodes[n];
1760
1/2
✓ Branch 1 taken 11579 times.
✗ Branch 2 not taken.
 23158 ValueType* data = distNode.buffer().data();
1761
1762 typename Int32TreeType::LeafNodeType* idxNode =
1763 idxAcc.probeLeaf(distNode.origin());
1764
1765 typename LeafNodeType::ValueOnCIter it;
1766
2/2
✓ Branch 0 taken 1997311 times.
✓ Branch 1 taken 11579 times.
 4040938 for (it = distNode.cbeginValueOn(); it; ++it) {
1767
1768 const Index pos = it.pos();
1769
1770
2/2
✓ Branch 0 taken 1400228 times.
✓ Branch 1 taken 597083 times.
 3994622 if (!(data[pos] > 0.75)) continue;
1771
1772 // Mask node internal neighbours
1773 1194166 maskNodeInternalNeighbours<LeafNodeType>(pos, neighbourMask);
1774
1775 const bool hasBoundaryNeighbour =
1776
4/4
✓ Branch 0 taken 9729 times.
✓ Branch 1 taken 587354 times.
✓ Branch 2 taken 9602 times.
✓ Branch 3 taken 127 times.
 1213624 checkNeighbours<Comp, LeafNodeType>(pos, data, neighbourMask) ||
1777
1/4
✓ Branch 2 taken 9729 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 1213370 checkNeighbours<Comp>(distNode.offsetToGlobalCoord(pos),distAcc,neighbourMask);
1778
1779 if (!hasBoundaryNeighbour) {
1780 254 distNode.setValueOff(pos);
1781 254 idxNode->setValueOff(pos);
1782 }
1783 }
1784 }
1785 }
1786
1787 LeafNodeType * * const mNodes;
1788 TreeType * const mDistTree;
1789 Int32TreeType * const mIndexTree;
1790 }; // RemoveSelfIntersectingSurface
1791
1792
1793 ////////////////////////////////////////
1794
1795
1796 template<typename NodeType>
1797 struct ReleaseChildNodes
1798 {
1799 ReleaseChildNodes(NodeType ** nodes) : mNodes(nodes) {}
1800
1801 void operator()(const tbb::blocked_range<size_t>& range) const {
1802
1803 using NodeMaskType = typename NodeType::NodeMaskType;
1804
1805 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1806 const_cast<NodeMaskType&>(mNodes[n]->getChildMask()).setOff();
1807 }
1808 }
1809
1810 NodeType ** const mNodes;
1811 };
1812
1813
1814 template<typename TreeType>
1815 inline void
1816 releaseLeafNodes(TreeType& tree)
1817 {
1818 using RootNodeType = typename TreeType::RootNodeType;
1819 using NodeChainType = typename RootNodeType::NodeChainType;
1820 using InternalNodeType = typename NodeChainType::template Get<1>;
1821
1822 std::vector<InternalNodeType*> nodes;
1823 tree.getNodes(nodes);
1824
1825 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
1826 ReleaseChildNodes<InternalNodeType>(nodes.empty() ? nullptr : &nodes[0]));
1827 }
1828
1829
1830 template<typename TreeType>
1831 struct StealUniqueLeafNodes
1832 {
1833 using LeafNodeType = typename TreeType::LeafNodeType;
1834
1835 136 StealUniqueLeafNodes(TreeType& lhsTree, TreeType& rhsTree,
1836 std::vector<LeafNodeType*>& overlappingNodes)
1837 : mLhsTree(&lhsTree)
1838 , mRhsTree(&rhsTree)
1839 136 , mNodes(&overlappingNodes)
1840 {
1841 }
1842
1843
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 272 void operator()() const {
1844
1845 std::vector<LeafNodeType*> rhsLeafNodes;
1846
1847
2/4
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 136 times.
✗ Branch 5 not taken.
 272 rhsLeafNodes.reserve(mRhsTree->leafCount());
1848 //mRhsTree->getNodes(rhsLeafNodes);
1849 //releaseLeafNodes(*mRhsTree);
1850
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 272 mRhsTree->stealNodes(rhsLeafNodes);
1851
1852
1/2
✓ Branch 1 taken 136 times.
✗ Branch 2 not taken.
 272 tree::ValueAccessor<TreeType> acc(*mLhsTree);
1853
1854
2/2
✓ Branch 0 taken 28500 times.
✓ Branch 1 taken 136 times.
 57272 for (size_t n = 0, N = rhsLeafNodes.size(); n < N; ++n) {
1855
2/2
✓ Branch 1 taken 23158 times.
✓ Branch 2 taken 5342 times.
 57000 if (!acc.probeLeaf(rhsLeafNodes[n]->origin())) {
1856
1/2
✓ Branch 1 taken 23158 times.
✗ Branch 2 not taken.
 46316 acc.addLeaf(rhsLeafNodes[n]);
1857 } else {
1858
1/2
✓ Branch 1 taken 5342 times.
✗ Branch 2 not taken.
 10684 mNodes->push_back(rhsLeafNodes[n]);
1859 }
1860 }
1861 }
1862
1863 private:
1864 TreeType * const mLhsTree;
1865 TreeType * const mRhsTree;
1866 std::vector<LeafNodeType*> * const mNodes;
1867 };
1868
1869
1870 template<typename DistTreeType, typename IndexTreeType>
1871 inline void
1872
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 combineData(DistTreeType& lhsDist, IndexTreeType& lhsIdx,
1873 DistTreeType& rhsDist, IndexTreeType& rhsIdx)
1874 {
1875 using DistLeafNodeType = typename DistTreeType::LeafNodeType;
1876 using IndexLeafNodeType = typename IndexTreeType::LeafNodeType;
1877
1878 std::vector<DistLeafNodeType*> overlappingDistNodes;
1879 std::vector<IndexLeafNodeType*> overlappingIdxNodes;
1880
1881 // Steal unique leafnodes
1882 tbb::task_group tasks;
1883
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 tasks.run(StealUniqueLeafNodes<DistTreeType>(lhsDist, rhsDist, overlappingDistNodes));
1884
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 tasks.run(StealUniqueLeafNodes<IndexTreeType>(lhsIdx, rhsIdx, overlappingIdxNodes));
1885
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 tasks.wait();
1886
1887 // Combine overlapping leaf nodes
1888
3/4
✓ Branch 0 taken 34 times.
✓ Branch 1 taken 34 times.
✓ Branch 2 taken 34 times.
✗ Branch 3 not taken.
 136 if (!overlappingDistNodes.empty() && !overlappingIdxNodes.empty()) {
1889
1/2
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
 68 tbb::parallel_for(tbb::blocked_range<size_t>(0, overlappingDistNodes.size()),
1890 CombineLeafNodes<DistTreeType>(lhsDist, lhsIdx,
1891 &overlappingDistNodes[0], &overlappingIdxNodes[0]));
1892 }
1893 }
1894
1895 /// @brief TBB body object to voxelize a mesh of triangles and/or quads into a collection
1896 /// of VDB grids, namely a squared distance grid, a closest primitive grid and an
1897 /// intersecting voxels grid (masks the mesh intersecting voxels)
1898 /// @note Only the leaf nodes that intersect the mesh are allocated, and only voxels in
1899 /// a narrow band (of two to three voxels in proximity to the mesh's surface) are activated.
1900 /// They are populated with distance values and primitive indices.
1901 template<typename TreeType>
1902 struct VoxelizationData {
1903
1904 using Ptr = std::unique_ptr<VoxelizationData>;
1905 using ValueType = typename TreeType::ValueType;
1906
1907 using Int32TreeType = typename TreeType::template ValueConverter<Int32>::Type;
1908 using UCharTreeType = typename TreeType::template ValueConverter<unsigned char>::Type;
1909
1910 using FloatTreeAcc = tree::ValueAccessor<TreeType>;
1911 using Int32TreeAcc = tree::ValueAccessor<Int32TreeType>;
1912 using UCharTreeAcc = tree::ValueAccessor<UCharTreeType>;
1913
1914
1915 136 VoxelizationData()
1916
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 : distTree(std::numeric_limits<ValueType>::max())
1917 , distAcc(distTree)
1918
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 , indexTree(Int32(util::INVALID_IDX))
1919 , indexAcc(indexTree)
1920 , primIdTree(MaxPrimId)
1921 , primIdAcc(primIdTree)
1922
1/2
✓ Branch 4 taken 68 times.
✗ Branch 5 not taken.
 408 , mPrimCount(0)
1923 {
1924 }
1925
1926 TreeType distTree;
1927 FloatTreeAcc distAcc;
1928
1929 Int32TreeType indexTree;
1930 Int32TreeAcc indexAcc;
1931
1932 UCharTreeType primIdTree;
1933 UCharTreeAcc primIdAcc;
1934
1935 69556 unsigned char getNewPrimId() {
1936
1937 /// @warning Don't use parallel methods here!
1938 /// The primIdTree is used as a "scratch" pad to mark visits for a given polygon
1939 /// into voxels which it may contribute to. The tree is kept as lightweight as
1940 /// possible and is reset when a maximum count or size is reached. A previous
1941 /// bug here occurred due to the calling of tree methods with multi-threaded
1942 /// implementations, resulting in nested parallelization and re-use of the TLS
1943 /// from the initial task. This consequently resulted in non deterministic values
1944 /// of mPrimCount on the return of the initial task, and could potentially end up
1945 /// with a mPrimCount equal to that of the MaxPrimId. This is used as the background
1946 /// value of the scratch tree.
1947 /// @see jira.aswf.io/browse/OVDB-117, PR #564
1948 /// @todo Consider profiling this operator with tree.clear() and Investigate the
1949 /// chosen value of MaxPrimId
1950
1951
3/4
✓ Branch 0 taken 34455 times.
✓ Branch 1 taken 323 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34455 times.
 138466 if (mPrimCount == MaxPrimId || primIdTree.leafCount() > 1000) {
1952 646 mPrimCount = 0;
1953 646 primIdTree.root().clear();
1954 646 primIdTree.clearAllAccessors();
1955
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 323 times.
 646 assert(mPrimCount == 0);
1956 }
1957
1958 69556 return mPrimCount++;
1959 }
1960
1961 private:
1962
1963 enum { MaxPrimId = 100 };
1964
1965 unsigned char mPrimCount;
1966 };
1967
1968
1969 template<typename TreeType, typename MeshDataAdapter, typename Interrupter = util::NullInterrupter>
1970 class VoxelizePolygons
1971 {
1972 public:
1973
1974 using VoxelizationDataType = VoxelizationData<TreeType>;
1975 using DataTable = tbb::enumerable_thread_specific<typename VoxelizationDataType::Ptr>;
1976
1977 35 VoxelizePolygons(DataTable& dataTable,
1978 const MeshDataAdapter& mesh,
1979 Interrupter* interrupter = nullptr)
1980 : mDataTable(&dataTable)
1981 , mMesh(&mesh)
1982 35 , mInterrupter(interrupter)
1983 {
1984 }
1985
1986 1498 void operator()(const tbb::blocked_range<size_t>& range) const {
1987
1988 1498 typename VoxelizationDataType::Ptr& dataPtr = mDataTable->local();
1989
2/2
✓ Branch 0 taken 68 times.
✓ Branch 1 taken 681 times.
 1634 if (!dataPtr) dataPtr.reset(new VoxelizationDataType());
1990
1991 Triangle prim;
1992
1993
2/2
✓ Branch 0 taken 1886 times.
✓ Branch 1 taken 749 times.
 5270 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
1994
1995 if (this->wasInterrupted()) {
1996 ✗ thread::cancelGroupExecution();
1997 break;
1998 }
1999
2000
2/2
✓ Branch 0 taken 1830 times.
✓ Branch 1 taken 24 times.
 3708 const size_t numVerts = mMesh->vertexCount(n);
2001
2002 // rasterize triangles and quads.
2003
1/2
✓ Branch 0 taken 1854 times.
✗ Branch 1 not taken.
 3708 if (numVerts == 3 || numVerts == 4) {
2004
2005 3772 prim.index = Int32(n);
2006
2007 3772 mMesh->getIndexSpacePoint(n, 0, prim.a);
2008 3772 mMesh->getIndexSpacePoint(n, 1, prim.b);
2009 3772 mMesh->getIndexSpacePoint(n, 2, prim.c);
2010
2011 3772 evalTriangle(prim, *dataPtr);
2012
2013
2/2
✓ Branch 0 taken 1830 times.
✓ Branch 1 taken 24 times.
 3708 if (numVerts == 4) {
2014 3660 mMesh->getIndexSpacePoint(n, 3, prim.b);
2015 3660 evalTriangle(prim, *dataPtr);
2016 }
2017 }
2018 }
2019 }
2020
2021 private:
2022
2023
4/16
✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 32 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 1854 times.
✗ Branch 11 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 1854 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
 1886 bool wasInterrupted() const { return mInterrupter && mInterrupter->wasInterrupted(); }
2024
2025 struct Triangle { Vec3d a, b, c; Int32 index; };
2026
2027 struct SubTask
2028 {
2029 enum { POLYGON_LIMIT = 1000 };
2030
2031 41416 SubTask(const Triangle& prim, DataTable& dataTable,
2032 int subdivisionCount, size_t polygonCount,
2033 Interrupter* interrupter = nullptr)
2034 : mLocalDataTable(&dataTable)
2035 , mPrim(prim)
2036 , mSubdivisionCount(subdivisionCount)
2037 , mPolygonCount(polygonCount)
2038 41416 , mInterrupter(interrupter)
2039 {
2040 }
2041
2042 82832 void operator()() const
2043 {
2044
4/4
✓ Branch 0 taken 37672 times.
✓ Branch 1 taken 3744 times.
✓ Branch 2 taken 27648 times.
✓ Branch 3 taken 10024 times.
 82832 if (mSubdivisionCount <= 0 || mPolygonCount >= POLYGON_LIMIT) {
2045
2046 62784 typename VoxelizationDataType::Ptr& dataPtr = mLocalDataTable->local();
2047
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 31392 times.
 62784 if (!dataPtr) dataPtr.reset(new VoxelizationDataType());
2048
2049 62784 voxelizeTriangle(mPrim, *dataPtr, mInterrupter);
2050
2051
2/4
✓ Branch 0 taken 10024 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 10024 times.
✗ Branch 4 not taken.
 20048 } else if (!(mInterrupter && mInterrupter->wasInterrupted())) {
2052 20048 spawnTasks(mPrim, *mLocalDataTable, mSubdivisionCount, mPolygonCount, mInterrupter);
2053 }
2054 }
2055
2056 DataTable * const mLocalDataTable;
2057 Triangle const mPrim;
2058 int const mSubdivisionCount;
2059 size_t const mPolygonCount;
2060 Interrupter * const mInterrupter;
2061 }; // struct SubTask
2062
2063
2/2
✓ Branch 0 taken 1360 times.
✓ Branch 1 taken 2356 times.
 7432 inline static int evalSubdivisionCount(const Triangle& prim)
2064 {
2065
2/2
✓ Branch 0 taken 1360 times.
✓ Branch 1 taken 2356 times.
 7432 const double ax = prim.a[0], bx = prim.b[0], cx = prim.c[0];
2066
4/4
✓ Branch 0 taken 1339 times.
✓ Branch 1 taken 2377 times.
✓ Branch 2 taken 821 times.
✓ Branch 3 taken 2895 times.
 14864 const double dx = std::max(ax, std::max(bx, cx)) - std::min(ax, std::min(bx, cx));
2067
2068
2/2
✓ Branch 0 taken 821 times.
✓ Branch 1 taken 2895 times.
 7432 const double ay = prim.a[1], by = prim.b[1], cy = prim.c[1];
2069
4/4
✓ Branch 0 taken 1886 times.
✓ Branch 1 taken 1830 times.
✓ Branch 2 taken 1093 times.
✓ Branch 3 taken 2623 times.
 14864 const double dy = std::max(ay, std::max(by, cy)) - std::min(ay, std::min(by, cy));
2070
2071
2/2
✓ Branch 0 taken 1093 times.
✓ Branch 1 taken 2623 times.
 7432 const double az = prim.a[2], bz = prim.b[2], cz = prim.c[2];
2072
4/4
✓ Branch 0 taken 1595 times.
✓ Branch 1 taken 2121 times.
✓ Branch 2 taken 1768 times.
✓ Branch 3 taken 1948 times.
 14864 const double dz = std::max(az, std::max(bz, cz)) - std::min(az, std::min(bz, cz));
2073
2074 7432 return int(std::max(dx, std::max(dy, dz)) / double(TreeType::LeafNodeType::DIM * 2));
2075 }
2076
2077 7432 void evalTriangle(const Triangle& prim, VoxelizationDataType& data) const
2078 {
2079 7432 const size_t polygonCount = mMesh->polygonCount();
2080 const int subdivisionCount =
2081
1/2
✓ Branch 0 taken 3716 times.
✗ Branch 1 not taken.
 7432 polygonCount < SubTask::POLYGON_LIMIT ? evalSubdivisionCount(prim) : 0;
2082
2083
2/2
✓ Branch 0 taken 3386 times.
✓ Branch 1 taken 330 times.
 7432 if (subdivisionCount <= 0) {
2084 6772 voxelizeTriangle(prim, data, mInterrupter);
2085 } else {
2086 660 spawnTasks(prim, *mDataTable, subdivisionCount, polygonCount, mInterrupter);
2087 }
2088 }
2089
2090 20708 static void spawnTasks(
2091 const Triangle& mainPrim,
2092 DataTable& dataTable,
2093 int subdivisionCount,
2094 size_t polygonCount,
2095 Interrupter* const interrupter)
2096 {
2097 20708 subdivisionCount -= 1;
2098 20708 polygonCount *= 4;
2099
2100 tbb::task_group tasks;
2101
2102 20708 const Vec3d ac = (mainPrim.a + mainPrim.c) * 0.5;
2103 20708 const Vec3d bc = (mainPrim.b + mainPrim.c) * 0.5;
2104 20708 const Vec3d ab = (mainPrim.a + mainPrim.b) * 0.5;
2105
2106 Triangle prim;
2107 20708 prim.index = mainPrim.index;
2108
2109 20708 prim.a = mainPrim.a;
2110 20708 prim.b = ab;
2111 20708 prim.c = ac;
2112
1/2
✓ Branch 1 taken 10354 times.
✗ Branch 2 not taken.
 20708 tasks.run(SubTask(prim, dataTable, subdivisionCount, polygonCount, interrupter));
2113
2114 20708 prim.a = ab;
2115 20708 prim.b = bc;
2116 20708 prim.c = ac;
2117
1/2
✓ Branch 1 taken 10354 times.
✗ Branch 2 not taken.
 20708 tasks.run(SubTask(prim, dataTable, subdivisionCount, polygonCount, interrupter));
2118
2119 20708 prim.a = ab;
2120 20708 prim.b = mainPrim.b;
2121 20708 prim.c = bc;
2122
1/2
✓ Branch 1 taken 10354 times.
✗ Branch 2 not taken.
 20708 tasks.run(SubTask(prim, dataTable, subdivisionCount, polygonCount, interrupter));
2123
2124 20708 prim.a = ac;
2125 20708 prim.b = bc;
2126 20708 prim.c = mainPrim.c;
2127
1/2
✓ Branch 1 taken 10354 times.
✗ Branch 2 not taken.
 20708 tasks.run(SubTask(prim, dataTable, subdivisionCount, polygonCount, interrupter));
2128
2129
1/2
✓ Branch 1 taken 10354 times.
✗ Branch 2 not taken.
 20708 tasks.wait();
2130 }
2131
2132 69556 static void voxelizeTriangle(const Triangle& prim, VoxelizationDataType& data, Interrupter* const interrupter)
2133 {
2134 std::deque<Coord> coordList;
2135 Coord ijk, nijk;
2136
2137 69556 ijk = Coord::floor(prim.a);
2138
1/2
✓ Branch 1 taken 34778 times.
✗ Branch 2 not taken.
 69556 coordList.push_back(ijk);
2139
2140 // The first point may not be quite in bounds, and rely
2141 // on one of the neighbours to have the first valid seed,
2142 // so we cannot early-exit here.
2143
1/2
✓ Branch 1 taken 34778 times.
✗ Branch 2 not taken.
 69556 updateDistance(ijk, prim, data);
2144
2145
1/2
✓ Branch 1 taken 34778 times.
✗ Branch 2 not taken.
 69556 unsigned char primId = data.getNewPrimId();
2146
1/2
✓ Branch 1 taken 34778 times.
✗ Branch 2 not taken.
 69556 data.primIdAcc.setValueOnly(ijk, primId);
2147
2148
2/2
✓ Branch 0 taken 34778 times.
✓ Branch 1 taken 34778 times.
 139112 while (!coordList.empty()) {
2149
3/6
✓ Branch 0 taken 34778 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 34778 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 34778 times.
 69556 if (interrupter && interrupter->wasInterrupted()) {
2150 ✗ thread::cancelGroupExecution();
2151 break;
2152 }
2153
3/4
✗ Branch 0 not taken.
✓ Branch 1 taken 1467770 times.
✓ Branch 2 taken 34778 times.
✓ Branch 3 taken 1432992 times.
 2935540 for (Int32 pass = 0; pass < 1048576 && !coordList.empty(); ++pass) {
2154
2/2
✓ Branch 0 taken 1427219 times.
✓ Branch 1 taken 5773 times.
 2865984 ijk = coordList.back();
2155 coordList.pop_back();
2156
2157
2/2
✓ Branch 0 taken 37257792 times.
✓ Branch 1 taken 1432992 times.
 77381568 for (Int32 i = 0; i < 26; ++i) {
2158
1/2
✓ Branch 1 taken 37257792 times.
✗ Branch 2 not taken.
 74515584 nijk = ijk + util::COORD_OFFSETS[i];
2159
3/4
✓ Branch 1 taken 37257792 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 7696414 times.
✓ Branch 4 taken 29561378 times.
 74515584 if (primId != data.primIdAcc.getValue(nijk)) {
2160
1/2
✓ Branch 1 taken 7696414 times.
✗ Branch 2 not taken.
 15392828 data.primIdAcc.setValueOnly(nijk, primId);
2161
4/6
✓ Branch 1 taken 7696414 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1398214 times.
✓ Branch 4 taken 6298200 times.
✓ Branch 6 taken 1398214 times.
✗ Branch 7 not taken.
 15392828 if(updateDistance(nijk, prim, data)) coordList.push_back(nijk);
2162 }
2163 }
2164 }
2165 }
2166 }
2167
2168 15462384 static bool updateDistance(const Coord& ijk, const Triangle& prim, VoxelizationDataType& data)
2169 {
2170 15462384 Vec3d uvw, voxelCenter(ijk[0], ijk[1], ijk[2]);
2171
2172 using ValueType = typename TreeType::ValueType;
2173
2174 15462384 const ValueType dist = ValueType((voxelCenter -
2175 15462384 closestPointOnTriangleToPoint(prim.a, prim.c, prim.b, voxelCenter, uvw)).lengthSqr());
2176
2177 // Either the points may be NAN, or they could be far enough from
2178 // the origin that computing distance fails.
2179
1/2
✓ Branch 0 taken 7731192 times.
✗ Branch 1 not taken.
 15462384 if (std::isnan(dist))
2180 return false;
2181
2182 15462384 const ValueType oldDist = data.distAcc.getValue(ijk);
2183
2184
2/2
✓ Branch 0 taken 4060747 times.
✓ Branch 1 taken 3670445 times.
 15462384 if (dist < oldDist) {
2185 8121494 data.distAcc.setValue(ijk, dist);
2186 8121494 data.indexAcc.setValue(ijk, prim.index);
2187
2/2
✓ Branch 0 taken 404873 times.
✓ Branch 1 taken 3265572 times.
 7340890 } else if (math::isExactlyEqual(dist, oldDist)) {
2188 // makes reduction deterministic when different polygons
2189 // produce the same distance value.
2190
2/2
✓ Branch 1 taken 44073 times.
✓ Branch 2 taken 360800 times.
 897892 data.indexAcc.setValueOnly(ijk, std::min(prim.index, data.indexAcc.getValue(ijk)));
2191 }
2192
2193 15462384 return !(dist > 0.75); // true if the primitive intersects the voxel.
2194 }
2195
2196 DataTable * const mDataTable;
2197 MeshDataAdapter const * const mMesh;
2198 Interrupter * const mInterrupter;
2199 }; // VoxelizePolygons
2200
2201
2202 ////////////////////////////////////////
2203
2204
2205 template<typename TreeType>
2206 struct DiffLeafNodeMask
2207 {
2208 using AccessorType = typename tree::ValueAccessor<TreeType>;
2209 using LeafNodeType = typename TreeType::LeafNodeType;
2210
2211 using BoolTreeType = typename TreeType::template ValueConverter<bool>::Type;
2212 using BoolLeafNodeType = typename BoolTreeType::LeafNodeType;
2213
2214 40 DiffLeafNodeMask(const TreeType& rhsTree,
2215 std::vector<BoolLeafNodeType*>& lhsNodes)
2216
1/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 40 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 40 : mRhsTree(&rhsTree), mLhsNodes(lhsNodes.empty() ? nullptr : &lhsNodes[0])
2217 {
2218 }
2219
2220 6668 void operator()(const tbb::blocked_range<size_t>& range) const {
2221
2222 6668 tree::ValueAccessor<const TreeType> acc(*mRhsTree);
2223
2224
2/2
✓ Branch 0 taken 19574 times.
✓ Branch 1 taken 3334 times.
 45816 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2225
2226 39148 BoolLeafNodeType* lhsNode = mLhsNodes[n];
2227 const LeafNodeType* rhsNode = acc.probeConstLeaf(lhsNode->origin());
2228
2229
2/2
✓ Branch 0 taken 15435 times.
✓ Branch 1 taken 4139 times.
 39148 if (rhsNode) lhsNode->topologyDifference(*rhsNode, false);
2230 }
2231 }
2232
2233 private:
2234 TreeType const * const mRhsTree;
2235 BoolLeafNodeType ** const mLhsNodes;
2236 };
2237
2238
2239 template<typename LeafNodeTypeA, typename LeafNodeTypeB>
2240 struct UnionValueMasks
2241 {
2242 40 UnionValueMasks(std::vector<LeafNodeTypeA*>& nodesA, std::vector<LeafNodeTypeB*>& nodesB)
2243
2/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 34 times.
✓ Branch 3 taken 6 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 40 : mNodesA(nodesA.empty() ? nullptr : &nodesA[0])
2244
3/16
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 34 times.
✓ Branch 6 taken 6 times.
✓ Branch 8 taken 40 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
 74 , mNodesB(nodesB.empty() ? nullptr : &nodesB[0])
2245 {
2246 }
2247
2248 6960 void operator()(const tbb::blocked_range<size_t>& range) const {
2249
2/2
✓ Branch 0 taken 13669 times.
✓ Branch 1 taken 3480 times.
 34298 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2250 27338 mNodesA[n]->topologyUnion(*mNodesB[n]);
2251 }
2252 }
2253
2254 private:
2255 LeafNodeTypeA ** const mNodesA;
2256 LeafNodeTypeB ** const mNodesB;
2257 };
2258
2259
2260 template<typename TreeType>
2261
1/8
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 29 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 124 struct ConstructVoxelMask
2262 {
2263 using LeafNodeType = typename TreeType::LeafNodeType;
2264
2265 using BoolTreeType = typename TreeType::template ValueConverter<bool>::Type;
2266 using BoolLeafNodeType = typename BoolTreeType::LeafNodeType;
2267
2268 29 ConstructVoxelMask(BoolTreeType& maskTree, const TreeType& tree,
2269 std::vector<LeafNodeType*>& nodes)
2270 : mTree(&tree)
2271
1/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 29 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 29 , mNodes(nodes.empty() ? nullptr : &nodes[0])
2272 58 , mLocalMaskTree(false)
2273
2/16
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 29 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 29 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
 58 , mMaskTree(&maskTree)
2274 {
2275 }
2276
2277 95 ConstructVoxelMask(ConstructVoxelMask& rhs, tbb::split)
2278 95 : mTree(rhs.mTree)
2279 95 , mNodes(rhs.mNodes)
2280 190 , mLocalMaskTree(false)
2281 95 , mMaskTree(&mLocalMaskTree)
2282 {
2283 }
2284
2285 5298 void operator()(const tbb::blocked_range<size_t>& range)
2286 {
2287 using Iterator = typename LeafNodeType::ValueOnCIter;
2288
2289 5298 tree::ValueAccessor<const TreeType> acc(*mTree);
2290
1/2
✓ Branch 1 taken 2649 times.
✗ Branch 2 not taken.
 5298 tree::ValueAccessor<BoolTreeType> maskAcc(*mMaskTree);
2291
2292 Coord ijk, nijk, localCorod;
2293 Index pos, npos;
2294
2295
2/2
✓ Branch 0 taken 11139 times.
✓ Branch 1 taken 2649 times.
 27576 for (size_t n = range.begin(); n != range.end(); ++n) {
2296
2297 22278 LeafNodeType& node = *mNodes[n];
2298
2299 22278 CoordBBox bbox = node.getNodeBoundingBox();
2300 bbox.expand(-1);
2301
2302
1/2
✓ Branch 1 taken 11139 times.
✗ Branch 2 not taken.
 22278 BoolLeafNodeType& maskNode = *maskAcc.touchLeaf(node.origin());
2303
2304
2/2
✓ Branch 0 taken 1929892 times.
✓ Branch 1 taken 11139 times.
 3882062 for (Iterator it = node.cbeginValueOn(); it; ++it) {
2305
1/2
✓ Branch 1 taken 1929892 times.
✗ Branch 2 not taken.
 3859784 ijk = it.getCoord();
2306 pos = it.pos();
2307
2308 3859784 localCorod = LeafNodeType::offsetToLocalCoord(pos);
2309
2310
2/2
✓ Branch 0 taken 1667878 times.
✓ Branch 1 taken 262014 times.
 3859784 if (localCorod[2] < int(LeafNodeType::DIM - 1)) {
2311 3335756 npos = pos + 1;
2312
2/2
✓ Branch 0 taken 123407 times.
✓ Branch 1 taken 1544471 times.
 3335756 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2313 } else {
2314 524028 nijk = ijk.offsetBy(0, 0, 1);
2315
2/2
✓ Branch 1 taken 77762 times.
✓ Branch 2 taken 184252 times.
 524028 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2316 }
2317
2318
2/2
✓ Branch 0 taken 1740828 times.
✓ Branch 1 taken 189064 times.
 3859784 if (localCorod[2] > 0) {
2319 3481656 npos = pos - 1;
2320
2/2
✓ Branch 0 taken 196357 times.
✓ Branch 1 taken 1544471 times.
 3481656 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2321 } else {
2322 378128 nijk = ijk.offsetBy(0, 0, -1);
2323
2/2
✓ Branch 1 taken 4812 times.
✓ Branch 2 taken 184252 times.
 378128 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2324 }
2325
2326
2/2
✓ Branch 0 taken 1665582 times.
✓ Branch 1 taken 264310 times.
 3859784 if (localCorod[1] < int(LeafNodeType::DIM - 1)) {
2327 3331164 npos = pos + LeafNodeType::DIM;
2328
2/2
✓ Branch 0 taken 124424 times.
✓ Branch 1 taken 1541158 times.
 3331164 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2329 } else {
2330 528620 nijk = ijk.offsetBy(0, 1, 0);
2331
2/2
✓ Branch 1 taken 78522 times.
✓ Branch 2 taken 185788 times.
 528620 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2332 }
2333
2334
2/2
✓ Branch 0 taken 1741999 times.
✓ Branch 1 taken 187893 times.
 3859784 if (localCorod[1] > 0) {
2335 3483998 npos = pos - LeafNodeType::DIM;
2336
2/2
✓ Branch 0 taken 200841 times.
✓ Branch 1 taken 1541158 times.
 3483998 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2337 } else {
2338 375786 nijk = ijk.offsetBy(0, -1, 0);
2339
2/2
✓ Branch 1 taken 2105 times.
✓ Branch 2 taken 185788 times.
 375786 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2340 }
2341
2342
2/2
✓ Branch 0 taken 1676206 times.
✓ Branch 1 taken 253686 times.
 3859784 if (localCorod[0] < int(LeafNodeType::DIM - 1)) {
2343 3352412 npos = pos + LeafNodeType::DIM * LeafNodeType::DIM;
2344
2/2
✓ Branch 0 taken 136461 times.
✓ Branch 1 taken 1539745 times.
 3352412 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2345 } else {
2346 507372 nijk = ijk.offsetBy(1, 0, 0);
2347
2/2
✓ Branch 1 taken 78079 times.
✓ Branch 2 taken 175607 times.
 507372 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2348 }
2349
2350
2/2
✓ Branch 0 taken 1751884 times.
✓ Branch 1 taken 178008 times.
 3859784 if (localCorod[0] > 0) {
2351 3503768 npos = pos - LeafNodeType::DIM * LeafNodeType::DIM;
2352
2/2
✓ Branch 0 taken 212139 times.
✓ Branch 1 taken 1539745 times.
 3503768 if (!node.isValueOn(npos)) maskNode.setValueOn(npos);
2353 } else {
2354 356016 nijk = ijk.offsetBy(-1, 0, 0);
2355
2/2
✓ Branch 1 taken 2401 times.
✓ Branch 2 taken 175607 times.
 356016 if (!acc.isValueOn(nijk)) maskAcc.setValueOn(nijk);
2356 }
2357 }
2358 }
2359 }
2360
2361 95 void join(ConstructVoxelMask& rhs) { mMaskTree->merge(*rhs.mMaskTree); }
2362
2363 private:
2364 TreeType const * const mTree;
2365 LeafNodeType ** const mNodes;
2366
2367 BoolTreeType mLocalMaskTree;
2368 BoolTreeType * const mMaskTree;
2369 };
2370
2371
2372 /// @note The interior and exterior widths should be in world space units and squared.
2373 template<typename TreeType, typename MeshDataAdapter>
2374 struct ExpandNarrowband
2375 {
2376 using ValueType = typename TreeType::ValueType;
2377 using LeafNodeType = typename TreeType::LeafNodeType;
2378 using NodeMaskType = typename LeafNodeType::NodeMaskType;
2379 using Int32TreeType = typename TreeType::template ValueConverter<Int32>::Type;
2380 using Int32LeafNodeType = typename Int32TreeType::LeafNodeType;
2381 using BoolTreeType = typename TreeType::template ValueConverter<bool>::Type;
2382 using BoolLeafNodeType = typename BoolTreeType::LeafNodeType;
2383
2384 struct Fragment
2385 {
2386 Int32 idx, x, y, z;
2387 ValueType dist;
2388
2389 Fragment() : idx(0), x(0), y(0), z(0), dist(0.0) {}
2390
2391 3375974 Fragment(Int32 idx_, Int32 x_, Int32 y_, Int32 z_, ValueType dist_)
2392 3375974 : idx(idx_), x(x_), y(y_), z(z_), dist(dist_)
2393 {
2394 }
2395
2396
18/96
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✗ Branch 47 not taken.
✗ Branch 48 not taken.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✗ Branch 51 not taken.
✓ Branch 52 taken 150920 times.
✓ Branch 53 taken 3357450 times.
✓ Branch 54 taken 435329 times.
✓ Branch 55 taken 7210603 times.
✓ Branch 56 taken 446701 times.
✓ Branch 57 taken 7210603 times.
✓ Branch 58 taken 10800 times.
✓ Branch 59 taken 282171 times.
✓ Branch 60 taken 916 times.
✓ Branch 61 taken 9884 times.
✓ Branch 62 taken 7287 times.
✓ Branch 63 taken 2597 times.
✓ Branch 64 taken 7341 times.
✓ Branch 65 taken 274830 times.
✓ Branch 66 taken 2903 times.
✓ Branch 67 taken 271927 times.
✗ Branch 68 not taken.
✗ Branch 69 not taken.
✗ Branch 70 not taken.
✗ Branch 71 not taken.
✗ Branch 72 not taken.
✗ Branch 73 not taken.
✗ Branch 74 not taken.
✗ Branch 75 not taken.
✗ Branch 76 not taken.
✗ Branch 77 not taken.
✗ Branch 78 not taken.
✗ Branch 79 not taken.
✗ Branch 80 not taken.
✗ Branch 81 not taken.
✗ Branch 82 not taken.
✗ Branch 83 not taken.
✗ Branch 84 not taken.
✗ Branch 85 not taken.
✗ Branch 86 not taken.
✗ Branch 87 not taken.
✗ Branch 88 not taken.
✗ Branch 89 not taken.
✗ Branch 90 not taken.
✗ Branch 91 not taken.
✓ Branch 92 taken 714 times.
✓ Branch 93 taken 265448 times.
✗ Branch 94 not taken.
✗ Branch 95 not taken.
 19370739 bool operator<(const Fragment& rhs) const { return idx < rhs.idx; }
2397 }; // struct Fragment
2398
2399 ////////////////////
2400
2401 40 ExpandNarrowband(
2402 std::vector<BoolLeafNodeType*>& maskNodes,
2403 BoolTreeType& maskTree,
2404 TreeType& distTree,
2405 Int32TreeType& indexTree,
2406 const MeshDataAdapter& mesh,
2407 ValueType exteriorBandWidth,
2408 ValueType interiorBandWidth,
2409 ValueType voxelSize)
2410
1/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 40 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 40 : mMaskNodes(maskNodes.empty() ? nullptr : &maskNodes[0])
2411 , mMaskTree(&maskTree)
2412 , mDistTree(&distTree)
2413 , mIndexTree(&indexTree)
2414 , mMesh(&mesh)
2415
1/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 40 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 40 , mNewMaskTree(false)
2416 , mDistNodes()
2417 , mUpdatedDistNodes()
2418 , mIndexNodes()
2419 , mUpdatedIndexNodes()
2420 , mExteriorBandWidth(exteriorBandWidth)
2421 , mInteriorBandWidth(interiorBandWidth)
2422
1/8
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 6 taken 40 times.
✗ Branch 7 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 40 , mVoxelSize(voxelSize)
2423 {
2424 }
2425
2426 296 ExpandNarrowband(const ExpandNarrowband& rhs, tbb::split)
2427 296 : mMaskNodes(rhs.mMaskNodes)
2428 296 , mMaskTree(rhs.mMaskTree)
2429 296 , mDistTree(rhs.mDistTree)
2430 296 , mIndexTree(rhs.mIndexTree)
2431 296 , mMesh(rhs.mMesh)
2432 592 , mNewMaskTree(false)
2433 , mDistNodes()
2434 , mUpdatedDistNodes()
2435 , mIndexNodes()
2436 , mUpdatedIndexNodes()
2437 296 , mExteriorBandWidth(rhs.mExteriorBandWidth)
2438 296 , mInteriorBandWidth(rhs.mInteriorBandWidth)
2439 296 , mVoxelSize(rhs.mVoxelSize)
2440 {
2441 }
2442
2443 296 void join(ExpandNarrowband& rhs)
2444 {
2445 296 mDistNodes.insert(mDistNodes.end(), rhs.mDistNodes.begin(), rhs.mDistNodes.end());
2446 296 mIndexNodes.insert(mIndexNodes.end(), rhs.mIndexNodes.begin(), rhs.mIndexNodes.end());
2447
2448 296 mUpdatedDistNodes.insert(mUpdatedDistNodes.end(),
2449 rhs.mUpdatedDistNodes.begin(), rhs.mUpdatedDistNodes.end());
2450
2451 296 mUpdatedIndexNodes.insert(mUpdatedIndexNodes.end(),
2452 rhs.mUpdatedIndexNodes.begin(), rhs.mUpdatedIndexNodes.end());
2453
2454 296 mNewMaskTree.merge(rhs.mNewMaskTree);
2455 }
2456
2457 7838 void operator()(const tbb::blocked_range<size_t>& range)
2458 {
2459 7838 tree::ValueAccessor<BoolTreeType> newMaskAcc(mNewMaskTree);
2460
1/2
✓ Branch 1 taken 3919 times.
✗ Branch 2 not taken.
 7838 tree::ValueAccessor<TreeType> distAcc(*mDistTree);
2461
1/2
✓ Branch 1 taken 3919 times.
✗ Branch 2 not taken.
 7838 tree::ValueAccessor<Int32TreeType> indexAcc(*mIndexTree);
2462
2463 std::vector<Fragment> fragments;
2464
1/2
✓ Branch 1 taken 3919 times.
✗ Branch 2 not taken.
 7838 fragments.reserve(256);
2465
2466 7838 std::unique_ptr<LeafNodeType> newDistNodePt;
2467 15676 std::unique_ptr<Int32LeafNodeType> newIndexNodePt;
2468
2469
2/2
✓ Branch 0 taken 19574 times.
✓ Branch 1 taken 3919 times.
 46986 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2470
2471 39148 BoolLeafNodeType& maskNode = *mMaskNodes[n];
2472
2/2
✓ Branch 0 taken 1764 times.
✓ Branch 1 taken 17810 times.
 39148 if (maskNode.isEmpty()) continue;
2473
2474 // Setup local caches
2475
2476 const Coord& origin = maskNode.origin();
2477
2478 35620 LeafNodeType * distNodePt = distAcc.probeLeaf(origin);
2479 35620 Int32LeafNodeType * indexNodePt = indexAcc.probeLeaf(origin);
2480
2481
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 17810 times.
 35620 assert(!distNodePt == !indexNodePt);
2482
2483 bool usingNewNodes = false;
2484
2485
3/4
✓ Branch 0 taken 4139 times.
✓ Branch 1 taken 13671 times.
✓ Branch 2 taken 4139 times.
✗ Branch 3 not taken.
 35620 if (!distNodePt && !indexNodePt) {
2486
2487
2/4
✓ Branch 1 taken 4139 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4139 times.
✗ Branch 4 not taken.
 8278 const ValueType backgroundDist = distAcc.getValue(origin);
2488
2489
2/4
✓ Branch 0 taken 4139 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4139 times.
✗ Branch 3 not taken.
 8278 if (!newDistNodePt.get() && !newIndexNodePt.get()) {
2490
2/4
✓ Branch 1 taken 4139 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4139 times.
✗ Branch 5 not taken.
 8278 newDistNodePt.reset(new LeafNodeType(origin, backgroundDist));
2491
3/6
✓ Branch 1 taken 4139 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4139 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 4139 times.
✗ Branch 8 not taken.
 8278 newIndexNodePt.reset(new Int32LeafNodeType(origin, indexAcc.getValue(origin)));
2492 } else {
2493
2494 ✗ if ((backgroundDist < ValueType(0.0)) !=
2495 ✗ (newDistNodePt->getValue(0) < ValueType(0.0))) {
2496 newDistNodePt->buffer().fill(backgroundDist);
2497 }
2498
2499 newDistNodePt->setOrigin(origin);
2500 newIndexNodePt->setOrigin(origin);
2501 }
2502
2503 8278 distNodePt = newDistNodePt.get();
2504 8278 indexNodePt = newIndexNodePt.get();
2505
2506 usingNewNodes = true;
2507 }
2508
2509
2510 // Gather neighbour information
2511
2512 35620 CoordBBox bbox(Coord::max(), Coord::min());
2513
2/2
✓ Branch 0 taken 1324184 times.
✓ Branch 1 taken 17810 times.
 2683988 for (typename BoolLeafNodeType::ValueOnIter it = maskNode.beginValueOn(); it; ++it) {
2514
1/2
✓ Branch 1 taken 1324184 times.
✗ Branch 2 not taken.
 2648368 bbox.expand(it.getCoord());
2515 }
2516
2517 bbox.expand(1);
2518
2519
1/2
✓ Branch 1 taken 17810 times.
✗ Branch 2 not taken.
 35620 gatherFragments(fragments, bbox, distAcc, indexAcc);
2520
2521
2522 // Compute first voxel layer
2523
2524 35620 bbox = maskNode.getNodeBoundingBox();
2525 NodeMaskType mask;
2526 35620 bool updatedLeafNodes = false;
2527
2528
2/2
✓ Branch 0 taken 1324184 times.
✓ Branch 1 taken 17810 times.
 2683988 for (typename BoolLeafNodeType::ValueOnIter it = maskNode.beginValueOn(); it; ++it) {
2529
2530
1/2
✓ Branch 1 taken 1324184 times.
✗ Branch 2 not taken.
 2648368 const Coord ijk = it.getCoord();
2531
2532
3/4
✓ Branch 1 taken 1324184 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 160241 times.
✓ Branch 4 taken 1163943 times.
 2648368 if (updateVoxel(ijk, 5, fragments, *distNodePt, *indexNodePt, &updatedLeafNodes)) {
2533
2534
2/2
✓ Branch 0 taken 160241 times.
✓ Branch 1 taken 961446 times.
 2243374 for (Int32 i = 0; i < 6; ++i) {
2535
2/2
✓ Branch 0 taken 938455 times.
✓ Branch 1 taken 22991 times.
 1922892 const Coord nijk = ijk + util::COORD_OFFSETS[i];
2536 if (bbox.isInside(nijk)) {
2537 1648844 mask.setOn(BoolLeafNodeType::coordToOffset(nijk));
2538 } else {
2539 newMaskAcc.setValueOn(nijk);
2540 }
2541 }
2542
2543
2/2
✓ Branch 0 taken 3204820 times.
✓ Branch 1 taken 160241 times.
 6730122 for (Int32 i = 6; i < 26; ++i) {
2544
2/2
✓ Branch 0 taken 3020892 times.
✓ Branch 1 taken 183928 times.
 6409640 const Coord nijk = ijk + util::COORD_OFFSETS[i];
2545 if (bbox.isInside(nijk)) {
2546 4441332 mask.setOn(BoolLeafNodeType::coordToOffset(nijk));
2547 }
2548 }
2549 }
2550 }
2551
2552
2/2
✓ Branch 0 taken 17808 times.
✓ Branch 1 taken 2 times.
 35620 if (updatedLeafNodes) {
2553
2554 // Compute second voxel layer
2555 35616 mask -= indexNodePt->getValueMask();
2556
2557
2/2
✓ Branch 1 taken 219730 times.
✓ Branch 2 taken 17808 times.
 510692 for (typename NodeMaskType::OnIterator it = mask.beginOn(); it; ++it) {
2558
2559 const Index pos = it.pos();
2560 439460 const Coord ijk = maskNode.origin() + LeafNodeType::offsetToLocalCoord(pos);
2561
2562
3/4
✓ Branch 1 taken 219730 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 188080 times.
✓ Branch 4 taken 31650 times.
 439460 if (updateVoxel(ijk, 6, fragments, *distNodePt, *indexNodePt)) {
2563
2/2
✓ Branch 0 taken 1128480 times.
✓ Branch 1 taken 188080 times.
 2633120 for (Int32 i = 0; i < 6; ++i) {
2564
1/2
✓ Branch 1 taken 1128480 times.
✗ Branch 2 not taken.
 2256960 newMaskAcc.setValueOn(ijk + util::COORD_OFFSETS[i]);
2565 }
2566 }
2567 }
2568
2569 // Export new distance values
2570
2/2
✓ Branch 0 taken 4139 times.
✓ Branch 1 taken 13669 times.
 35616 if (usingNewNodes) {
2571 newDistNodePt->topologyUnion(*newIndexNodePt);
2572
1/2
✓ Branch 1 taken 4139 times.
✗ Branch 2 not taken.
 8278 mDistNodes.push_back(newDistNodePt.release());
2573
1/2
✓ Branch 1 taken 4139 times.
✗ Branch 2 not taken.
 8278 mIndexNodes.push_back(newIndexNodePt.release());
2574 } else {
2575
1/2
✓ Branch 1 taken 13669 times.
✗ Branch 2 not taken.
 27338 mUpdatedDistNodes.push_back(distNodePt);
2576
1/2
✓ Branch 1 taken 13669 times.
✗ Branch 2 not taken.
 27338 mUpdatedIndexNodes.push_back(indexNodePt);
2577 }
2578 }
2579 } // end leafnode loop
2580 }
2581
2582 //////////
2583
2584 BoolTreeType& newMaskTree() { return mNewMaskTree; }
2585
2586 std::vector<LeafNodeType*>& newDistNodes() { return mDistNodes; }
2587 std::vector<LeafNodeType*>& updatedDistNodes() { return mUpdatedDistNodes; }
2588
2589 std::vector<Int32LeafNodeType*>& newIndexNodes() { return mIndexNodes; }
2590 std::vector<Int32LeafNodeType*>& updatedIndexNodes() { return mUpdatedIndexNodes; }
2591
2592 private:
2593
2594 /// @note The output fragment list is ordered by the primitive index
2595 void
2596
2/2
✓ Branch 0 taken 14059 times.
✓ Branch 1 taken 3751 times.
 35620 gatherFragments(std::vector<Fragment>& fragments, const CoordBBox& bbox,
2597 tree::ValueAccessor<TreeType>& distAcc, tree::ValueAccessor<Int32TreeType>& indexAcc)
2598 {
2599 fragments.clear();
2600 const Coord nodeMin = bbox.min() & ~(LeafNodeType::DIM - 1);
2601 const Coord nodeMax = bbox.max() & ~(LeafNodeType::DIM - 1);
2602
2603 35620 CoordBBox region;
2604 Coord ijk;
2605
2606
2/2
✓ Branch 0 taken 44480 times.
✓ Branch 1 taken 17810 times.
 124580 for (ijk[0] = nodeMin[0]; ijk[0] <= nodeMax[0]; ijk[0] += LeafNodeType::DIM) {
2607
2/2
✓ Branch 0 taken 109083 times.
✓ Branch 1 taken 44480 times.
 307126 for (ijk[1] = nodeMin[1]; ijk[1] <= nodeMax[1]; ijk[1] += LeafNodeType::DIM) {
2608
2/2
✓ Branch 0 taken 258005 times.
✓ Branch 1 taken 109083 times.
 734176 for (ijk[2] = nodeMin[2]; ijk[2] <= nodeMax[2]; ijk[2] += LeafNodeType::DIM) {
2609
2/2
✓ Branch 1 taken 178144 times.
✓ Branch 2 taken 79861 times.
 516010 if (LeafNodeType* distleaf = distAcc.probeLeaf(ijk)) {
2610 356288 region.min() = Coord::maxComponent(bbox.min(), ijk);
2611 356288 region.max() = Coord::minComponent(bbox.max(),
2612 ijk.offsetBy(LeafNodeType::DIM - 1));
2613 356288 gatherFragments(fragments, region, *distleaf, *indexAcc.probeLeaf(ijk));
2614 }
2615 }
2616 }
2617 }
2618
2619 35620 std::sort(fragments.begin(), fragments.end());
2620 }
2621
2622 void
2623 356288 gatherFragments(std::vector<Fragment>& fragments, const CoordBBox& bbox,
2624 const LeafNodeType& distLeaf, const Int32LeafNodeType& idxLeaf) const
2625 {
2626 const typename LeafNodeType::NodeMaskType& mask = distLeaf.getValueMask();
2627 356288 const ValueType* distData = distLeaf.buffer().data();
2628 356288 const Int32* idxData = idxLeaf.buffer().data();
2629
2630
2/2
✓ Branch 0 taken 608496 times.
✓ Branch 1 taken 178144 times.
 1573280 for (int x = bbox.min()[0]; x <= bbox.max()[0]; ++x) {
2631 1216992 const Index xPos = (x & (LeafNodeType::DIM - 1u)) << (2 * LeafNodeType::LOG2DIM);
2632
2/2
✓ Branch 0 taken 2049922 times.
✓ Branch 1 taken 608496 times.
 5316836 for (int y = bbox.min()[1]; y <= bbox.max()[1]; ++y) {
2633 4099844 const Index yPos = xPos + ((y & (LeafNodeType::DIM - 1u)) << LeafNodeType::LOG2DIM);
2634
2/2
✓ Branch 0 taken 6905748 times.
✓ Branch 1 taken 2049922 times.
 17911340 for (int z = bbox.min()[2]; z <= bbox.max()[2]; ++z) {
2635 13811496 const Index pos = yPos + (z & (LeafNodeType::DIM - 1u));
2636
2/2
✓ Branch 1 taken 3375974 times.
✓ Branch 2 taken 3529774 times.
 13811496 if (mask.isOn(pos)) {
2637 6751948 fragments.push_back(Fragment(idxData[pos],x,y,z, std::abs(distData[pos])));
2638 }
2639 }
2640 }
2641 }
2642 }
2643
2644 /// @note This method expects the fragment list to be ordered by the primitive index
2645 /// to avoid redundant distance computations.
2646 ValueType
2647 3087828 computeDistance(const Coord& ijk, const Int32 manhattanLimit,
2648 const std::vector<Fragment>& fragments, Int32& closestPrimIdx) const
2649 {
2650 3087828 Vec3d a, b, c, uvw, voxelCenter(ijk[0], ijk[1], ijk[2]);
2651 double primDist, tmpDist, dist = std::numeric_limits<double>::max();
2652 3087828 Int32 lastIdx = Int32(util::INVALID_IDX);
2653
2654
2/2
✓ Branch 0 taken 355178889 times.
✓ Branch 1 taken 1543914 times.
 713445606 for (size_t n = 0, N = fragments.size(); n < N; ++n) {
2655
2656 const Fragment& fragment = fragments[n];
2657
2/2
✓ Branch 0 taken 324681202 times.
✓ Branch 1 taken 30497687 times.
 710357778 if (lastIdx == fragment.idx) continue;
2658
2659
2/2
✓ Branch 0 taken 28424624 times.
✓ Branch 1 taken 2073063 times.
 60995374 const Int32 dx = std::abs(fragment.x - ijk[0]);
2660 60995374 const Int32 dy = std::abs(fragment.y - ijk[1]);
2661 60995374 const Int32 dz = std::abs(fragment.z - ijk[2]);
2662
2663 60995374 const Int32 manhattan = dx + dy + dz;
2664
2/2
✓ Branch 0 taken 28424624 times.
✓ Branch 1 taken 2073063 times.
 60995374 if (manhattan > manhattanLimit) continue;
2665
2666 lastIdx = fragment.idx;
2667
2668 4146126 const size_t polygon = size_t(lastIdx);
2669
2670 4146126 mMesh->getIndexSpacePoint(polygon, 0, a);
2671 4146126 mMesh->getIndexSpacePoint(polygon, 1, b);
2672 4146126 mMesh->getIndexSpacePoint(polygon, 2, c);
2673
2674 4146126 primDist = (voxelCenter -
2675 4146126 closestPointOnTriangleToPoint(a, c, b, voxelCenter, uvw)).lengthSqr();
2676
2677 // Split quad into a second triangle
2678
2/2
✓ Branch 0 taken 1857942 times.
✓ Branch 1 taken 215121 times.
 4146126 if (4 == mMesh->vertexCount(polygon)) {
2679
2680 3715884 mMesh->getIndexSpacePoint(polygon, 3, b);
2681
2682 3715884 tmpDist = (voxelCenter - closestPointOnTriangleToPoint(
2683 a, b, c, voxelCenter, uvw)).lengthSqr();
2684
2685
2/2
✓ Branch 0 taken 950225 times.
✓ Branch 1 taken 907717 times.
 3715884 if (tmpDist < primDist) primDist = tmpDist;
2686 }
2687
2688
2/2
✓ Branch 0 taken 1732047 times.
✓ Branch 1 taken 341016 times.
 4146126 if (primDist < dist) {
2689 dist = primDist;
2690 3464094 closestPrimIdx = lastIdx;
2691 }
2692 }
2693
2694 3087828 return ValueType(std::sqrt(dist)) * mVoxelSize;
2695 }
2696
2697 /// @note Returns true if the current voxel was updated and neighbouring
2698 /// voxels need to be evaluated.
2699 bool
2700 3087828 updateVoxel(const Coord& ijk, const Int32 manhattanLimit,
2701 const std::vector<Fragment>& fragments,
2702 LeafNodeType& distLeaf, Int32LeafNodeType& idxLeaf, bool* updatedLeafNodes = nullptr)
2703 {
2704 3087828 Int32 closestPrimIdx = 0;
2705 3087828 const ValueType distance = computeDistance(ijk, manhattanLimit, fragments, closestPrimIdx);
2706
2707 const Index pos = LeafNodeType::coordToOffset(ijk);
2708 3087828 const bool inside = distLeaf.getValue(pos) < ValueType(0.0);
2709
2710 bool activateNeighbourVoxels = false;
2711
2712
4/4
✓ Branch 0 taken 887838 times.
✓ Branch 1 taken 656076 times.
✓ Branch 2 taken 880166 times.
✓ Branch 3 taken 7672 times.
 3087828 if (!inside && distance < mExteriorBandWidth) {
2713
2/2
✓ Branch 0 taken 742649 times.
✓ Branch 1 taken 137517 times.
 1760332 if (updatedLeafNodes) *updatedLeafNodes = true;
2714 1760332 activateNeighbourVoxels = (distance + mVoxelSize) < mExteriorBandWidth;
2715 1760332 distLeaf.setValueOnly(pos, distance);
2716 idxLeaf.setValueOn(pos, closestPrimIdx);
2717
4/4
✓ Branch 0 taken 656076 times.
✓ Branch 1 taken 7672 times.
✓ Branch 2 taken 654110 times.
✓ Branch 3 taken 1966 times.
 1327496 } else if (inside && distance < mInteriorBandWidth) {
2718
2/2
✓ Branch 0 taken 578847 times.
✓ Branch 1 taken 75263 times.
 1308220 if (updatedLeafNodes) *updatedLeafNodes = true;
2719 1308220 activateNeighbourVoxels = (distance + mVoxelSize) < mInteriorBandWidth;
2720 1308220 distLeaf.setValueOnly(pos, -distance);
2721 idxLeaf.setValueOn(pos, closestPrimIdx);
2722 }
2723
2724 3087828 return activateNeighbourVoxels;
2725 }
2726
2727 //////////
2728
2729 BoolLeafNodeType ** const mMaskNodes;
2730 BoolTreeType * const mMaskTree;
2731 TreeType * const mDistTree;
2732 Int32TreeType * const mIndexTree;
2733
2734 MeshDataAdapter const * const mMesh;
2735
2736 BoolTreeType mNewMaskTree;
2737
2738 std::vector<LeafNodeType*> mDistNodes, mUpdatedDistNodes;
2739 std::vector<Int32LeafNodeType*> mIndexNodes, mUpdatedIndexNodes;
2740
2741 const ValueType mExteriorBandWidth, mInteriorBandWidth, mVoxelSize;
2742 }; // struct ExpandNarrowband
2743
2744
2745 template<typename TreeType>
2746 struct AddNodes {
2747 using LeafNodeType = typename TreeType::LeafNodeType;
2748
2749 80 AddNodes(TreeType& tree, std::vector<LeafNodeType*>& nodes)
2750 80 : mTree(&tree) , mNodes(&nodes)
2751 {
2752 }
2753
2754 160 void operator()() const {
2755 160 tree::ValueAccessor<TreeType> acc(*mTree);
2756 160 std::vector<LeafNodeType*>& nodes = *mNodes;
2757
2/2
✓ Branch 0 taken 8278 times.
✓ Branch 1 taken 80 times.
 16716 for (size_t n = 0, N = nodes.size(); n < N; ++n) {
2758
1/2
✓ Branch 1 taken 8278 times.
✗ Branch 2 not taken.
 16556 acc.addLeaf(nodes[n]);
2759 }
2760 }
2761
2762 TreeType * const mTree;
2763 std::vector<LeafNodeType*> * const mNodes;
2764 }; // AddNodes
2765
2766
2767 template<typename TreeType, typename Int32TreeType, typename BoolTreeType, typename MeshDataAdapter>
2768 inline void
2769
1/2
✓ Branch 0 taken 40 times.
✗ Branch 1 not taken.
 80 expandNarrowband(
2770 TreeType& distTree,
2771 Int32TreeType& indexTree,
2772 BoolTreeType& maskTree,
2773 std::vector<typename BoolTreeType::LeafNodeType*>& maskNodes,
2774 const MeshDataAdapter& mesh,
2775 typename TreeType::ValueType exteriorBandWidth,
2776 typename TreeType::ValueType interiorBandWidth,
2777 typename TreeType::ValueType voxelSize)
2778 {
2779 80 ExpandNarrowband<TreeType, MeshDataAdapter> expandOp(maskNodes, maskTree,
2780 distTree, indexTree, mesh, exteriorBandWidth, interiorBandWidth, voxelSize);
2781
2782
3/4
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✓ Branch 4 taken 6 times.
 80 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, maskNodes.size()), expandOp);
2783
2784
2/4
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 40 times.
✗ Branch 5 not taken.
 80 tbb::parallel_for(tbb::blocked_range<size_t>(0, expandOp.updatedIndexNodes().size()),
2785 UnionValueMasks<typename TreeType::LeafNodeType, typename Int32TreeType::LeafNodeType>(
2786 expandOp.updatedDistNodes(), expandOp.updatedIndexNodes()));
2787
2788 tbb::task_group tasks;
2789
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 tasks.run(AddNodes<TreeType>(distTree, expandOp.newDistNodes()));
2790
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 tasks.run(AddNodes<Int32TreeType>(indexTree, expandOp.newIndexNodes()));
2791
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 tasks.wait();
2792
2793
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 maskTree.clear();
2794
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 maskTree.merge(expandOp.newMaskTree());
2795 }
2796
2797
2798 ////////////////////////////////////////
2799
2800
2801 // Transform values (sqrt, world space scaling and sign flip if sdf)
2802 template<typename TreeType>
2803 struct TransformValues
2804 {
2805 using LeafNodeType = typename TreeType::LeafNodeType;
2806 using ValueType = typename TreeType::ValueType;
2807
2808 34 TransformValues(std::vector<LeafNodeType*>& nodes,
2809 ValueType voxelSize, bool unsignedDist)
2810 : mNodes(&nodes[0])
2811 , mVoxelSize(voxelSize)
2812
2/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 31 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
 34 , mUnsigned(unsignedDist)
2813 {
2814 }
2815
2816 6350 void operator()(const tbb::blocked_range<size_t>& range) const {
2817
2818 typename LeafNodeType::ValueOnIter iter;
2819
2820 6350 const bool udf = mUnsigned;
2821 6350 const ValueType w[2] = { -mVoxelSize, mVoxelSize };
2822
2823
2/2
✓ Branch 0 taken 11579 times.
✓ Branch 1 taken 3175 times.
 29508 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2824
2825
2/2
✓ Branch 1 taken 1997184 times.
✓ Branch 2 taken 11579 times.
 4040684 for (iter = mNodes[n]->beginValueOn(); iter; ++iter) {
2826 3994368 ValueType& val = const_cast<ValueType&>(iter.getValue());
2827
3/4
✓ Branch 0 taken 1997184 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1199457 times.
✓ Branch 3 taken 797727 times.
 6393282 val = w[udf || (val < ValueType(0.0))] * std::sqrt(std::abs(val));
2828 }
2829 }
2830 }
2831
2832 private:
2833 LeafNodeType * * const mNodes;
2834 const ValueType mVoxelSize;
2835 const bool mUnsigned;
2836 };
2837
2838
2839 // Inactivate values outside the (exBandWidth, inBandWidth) range.
2840 template<typename TreeType>
2841 struct InactivateValues
2842 {
2843 using LeafNodeType = typename TreeType::LeafNodeType;
2844 using ValueType = typename TreeType::ValueType;
2845
2846 32 InactivateValues(std::vector<LeafNodeType*>& nodes,
2847 ValueType exBandWidth, ValueType inBandWidth)
2848
2/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 29 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
 32 : mNodes(nodes.empty() ? nullptr : &nodes[0])
2849 , mExBandWidth(exBandWidth)
2850
2/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 29 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
 32 , mInBandWidth(inBandWidth)
2851 {
2852 }
2853
2854 6204 void operator()(const tbb::blocked_range<size_t>& range) const {
2855
2856 typename LeafNodeType::ValueOnIter iter;
2857 6204 const ValueType exVal = mExBandWidth;
2858 6204 const ValueType inVal = -mInBandWidth;
2859
2860
2/2
✓ Branch 0 taken 14090 times.
✓ Branch 1 taken 3102 times.
 34384 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2861
2862
2/2
✓ Branch 1 taken 3008851 times.
✓ Branch 2 taken 14090 times.
 6074062 for (iter = mNodes[n]->beginValueOn(); iter; ++iter) {
2863
2864 6017702 ValueType& val = const_cast<ValueType&>(iter.getValue());
2865
2866 6017702 const bool inside = val < ValueType(0.0);
2867
2868
4/4
✓ Branch 0 taken 1659482 times.
✓ Branch 1 taken 1349369 times.
✓ Branch 2 taken 3345 times.
✓ Branch 3 taken 1656137 times.
 6017702 if (inside && !(val > inVal)) {
2869 6690 val = inVal;
2870 6690 iter.setValueOff();
2871
4/4
✓ Branch 0 taken 1349369 times.
✓ Branch 1 taken 1656137 times.
✓ Branch 2 taken 6727 times.
✓ Branch 3 taken 1342642 times.
 6011012 } else if (!inside && !(val < exVal)) {
2872 13454 val = exVal;
2873 13454 iter.setValueOff();
2874 }
2875 }
2876 }
2877 }
2878
2879 private:
2880 LeafNodeType * * const mNodes;
2881 const ValueType mExBandWidth, mInBandWidth;
2882 };
2883
2884
2885 template<typename TreeType>
2886 struct OffsetValues
2887 {
2888 using LeafNodeType = typename TreeType::LeafNodeType;
2889 using ValueType = typename TreeType::ValueType;
2890
2891 60 OffsetValues(std::vector<LeafNodeType*>& nodes, ValueType offset)
2892
8/32
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 27 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 27 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 27 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 27 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✓ Branch 30 taken 3 times.
✗ Branch 31 not taken.
✓ Branch 33 taken 3 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 3 times.
✗ Branch 36 not taken.
✓ Branch 38 taken 3 times.
✗ Branch 39 not taken.
 120 : mNodes(nodes.empty() ? nullptr : &nodes[0]), mOffset(offset)
2893 {
2894 }
2895
2896 13030 void operator()(const tbb::blocked_range<size_t>& range) const {
2897
2898 13030 const ValueType offset = mOffset;
2899
2900
2/2
✓ Branch 0 taken 31036 times.
✓ Branch 1 taken 6515 times.
 75102 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2901
2902 62072 typename LeafNodeType::ValueOnIter iter = mNodes[n]->beginValueOn();
2903
2904
2/2
✓ Branch 0 taken 6994568 times.
✓ Branch 1 taken 31036 times.
 14051208 for (; iter; ++iter) {
2905 13989136 ValueType& val = const_cast<ValueType&>(iter.getValue());
2906 13989136 val += offset;
2907 }
2908 }
2909 }
2910
2911 private:
2912 LeafNodeType * * const mNodes;
2913 const ValueType mOffset;
2914 };
2915
2916
2917 template<typename TreeType>
2918 struct Renormalize
2919 {
2920 using LeafNodeType = typename TreeType::LeafNodeType;
2921 using ValueType = typename TreeType::ValueType;
2922
2923 30 Renormalize(const TreeType& tree, const std::vector<LeafNodeType*>& nodes,
2924 ValueType* buffer, ValueType voxelSize)
2925 : mTree(&tree)
2926
2/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 27 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
 30 , mNodes(nodes.empty() ? nullptr : &nodes[0])
2927 , mBuffer(buffer)
2928
2/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 4 taken 27 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
 30 , mVoxelSize(voxelSize)
2929 {
2930 }
2931
2932 6220 void operator()(const tbb::blocked_range<size_t>& range) const
2933 {
2934 using Vec3Type = math::Vec3<ValueType>;
2935
2936 6220 tree::ValueAccessor<const TreeType> acc(*mTree);
2937
2938 Coord ijk;
2939 Vec3Type up, down;
2940
2941 6220 const ValueType dx = mVoxelSize, invDx = ValueType(1.0) / mVoxelSize;
2942
2943
2/2
✓ Branch 0 taken 15518 times.
✓ Branch 1 taken 3110 times.
 37256 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2944
2945 31036 ValueType* bufferData = &mBuffer[n * LeafNodeType::SIZE];
2946
2947 31036 typename LeafNodeType::ValueOnCIter iter = mNodes[n]->cbeginValueOn();
2948
2/2
✓ Branch 0 taken 3497284 times.
✓ Branch 1 taken 15518 times.
 7025604 for (; iter; ++iter) {
2949
2950 6994568 const ValueType phi0 = *iter;
2951
2952
1/2
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
 6994568 ijk = iter.getCoord();
2953
2954
2/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3497284 times.
✗ Branch 5 not taken.
 6994568 up[0] = acc.getValue(ijk.offsetBy(1, 0, 0)) - phi0;
2955
2/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3497284 times.
✗ Branch 5 not taken.
 6994568 up[1] = acc.getValue(ijk.offsetBy(0, 1, 0)) - phi0;
2956
2/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3497284 times.
✗ Branch 5 not taken.
 6994568 up[2] = acc.getValue(ijk.offsetBy(0, 0, 1)) - phi0;
2957
2958
2/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3497284 times.
✗ Branch 5 not taken.
 6994568 down[0] = phi0 - acc.getValue(ijk.offsetBy(-1, 0, 0));
2959
2/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3497284 times.
✗ Branch 5 not taken.
 6994568 down[1] = phi0 - acc.getValue(ijk.offsetBy(0, -1, 0));
2960
1/4
✓ Branch 1 taken 3497284 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 6994568 down[2] = phi0 - acc.getValue(ijk.offsetBy(0, 0, -1));
2961
2962 6994568 const ValueType normSqGradPhi = math::GodunovsNormSqrd(phi0 > 0.0, down, up);
2963
2964 6994568 const ValueType diff = math::Sqrt(normSqGradPhi) * invDx - ValueType(1.0);
2965 6994568 const ValueType S = phi0 / (math::Sqrt(math::Pow2(phi0) + normSqGradPhi));
2966
2967 6994568 bufferData[iter.pos()] = phi0 - dx * S * diff;
2968 }
2969 }
2970 }
2971
2972 private:
2973 TreeType const * const mTree;
2974 LeafNodeType const * const * const mNodes;
2975 ValueType * const mBuffer;
2976
2977 const ValueType mVoxelSize;
2978 };
2979
2980
2981 template<typename TreeType>
2982 struct MinCombine
2983 {
2984 using LeafNodeType = typename TreeType::LeafNodeType;
2985 using ValueType = typename TreeType::ValueType;
2986
2987 30 MinCombine(std::vector<LeafNodeType*>& nodes, const ValueType* buffer)
2988
4/16
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 27 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 27 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 3 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 3 times.
✗ Branch 19 not taken.
 60 : mNodes(nodes.empty() ? nullptr : &nodes[0]), mBuffer(buffer)
2989 {
2990 }
2991
2992 6200 void operator()(const tbb::blocked_range<size_t>& range) const {
2993
2994
2/2
✓ Branch 0 taken 15518 times.
✓ Branch 1 taken 3100 times.
 37236 for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
2995
2996 31036 const ValueType* bufferData = &mBuffer[n * LeafNodeType::SIZE];
2997
2998 31036 typename LeafNodeType::ValueOnIter iter = mNodes[n]->beginValueOn();
2999
3000
2/2
✓ Branch 0 taken 3497284 times.
✓ Branch 1 taken 15518 times.
 7025604 for (; iter; ++iter) {
3001 6994568 ValueType& val = const_cast<ValueType&>(iter.getValue());
3002
2/2
✓ Branch 0 taken 198274 times.
✓ Branch 1 taken 3299010 times.
 7391116 val = std::min(val, bufferData[iter.pos()]);
3003 }
3004 }
3005 }
3006
3007 private:
3008 LeafNodeType * * const mNodes;
3009 ValueType const * const mBuffer;
3010 };
3011
3012
3013 } // mesh_to_volume_internal namespace
3014
3015 /// @endcond
3016
3017
3018 ////////////////////////////////////////
3019
3020 // Utility method implementation
3021
3022
3023 template <typename FloatTreeT>
3024 void
3025 72 traceExteriorBoundaries(FloatTreeT& tree)
3026 {
3027 using ConnectivityTable = mesh_to_volume_internal::LeafNodeConnectivityTable<FloatTreeT>;
3028
3029 // Build a node connectivity table where each leaf node has an offset into a
3030 // linearized list of nodes, and each leaf stores its six axis aligned neighbor
3031 // offsets
3032 144 ConnectivityTable nodeConnectivity(tree);
3033
3034 std::vector<size_t> zStartNodes, yStartNodes, xStartNodes;
3035
3036 // Store all nodes which do not have negative neighbors i.e. the nodes furthest
3037 // in -X, -Y, -Z. We sweep from lowest coordinate positions +axis and then
3038 // from the furthest positive coordinate positions -axis
3039
2/2
✓ Branch 0 taken 12431 times.
✓ Branch 1 taken 37 times.
 24082 for (size_t n = 0; n < nodeConnectivity.size(); ++n) {
3040
2/2
✓ Branch 0 taken 2438 times.
✓ Branch 1 taken 9993 times.
 24010 if (ConnectivityTable::INVALID_OFFSET == nodeConnectivity.offsetsPrevX()[n]) {
3041
1/2
✓ Branch 1 taken 2438 times.
✗ Branch 2 not taken.
 4805 xStartNodes.push_back(n);
3042 }
3043
3044
2/2
✓ Branch 0 taken 2412 times.
✓ Branch 1 taken 10019 times.
 24010 if (ConnectivityTable::INVALID_OFFSET == nodeConnectivity.offsetsPrevY()[n]) {
3045
1/2
✓ Branch 1 taken 2412 times.
✗ Branch 2 not taken.
 4675 yStartNodes.push_back(n);
3046 }
3047
3048
2/2
✓ Branch 0 taken 2377 times.
✓ Branch 1 taken 10054 times.
 24010 if (ConnectivityTable::INVALID_OFFSET == nodeConnectivity.offsetsPrevZ()[n]) {
3049
1/2
✓ Branch 1 taken 2377 times.
✗ Branch 2 not taken.
 4605 zStartNodes.push_back(n);
3050 }
3051 }
3052
3053 using SweepingOp = mesh_to_volume_internal::SweepExteriorSign<FloatTreeT>;
3054
3055 // Sweep the exterior value signs (make them negative) up until the voxel intersection
3056 // with the isosurface. Do this in both lowest -> + and largest -> - directions
3057
3058
3/4
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36 times.
✓ Branch 4 taken 1 times.
 72 tbb::parallel_for(tbb::blocked_range<size_t>(0, zStartNodes.size()),
3059 SweepingOp(SweepingOp::Z_AXIS, zStartNodes, nodeConnectivity));
3060
3061
3/4
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 36 times.
✓ Branch 4 taken 1 times.
 72 tbb::parallel_for(tbb::blocked_range<size_t>(0, yStartNodes.size()),
3062 SweepingOp(SweepingOp::Y_AXIS, yStartNodes, nodeConnectivity));
3063
3064
2/6
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 37 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 72 tbb::parallel_for(tbb::blocked_range<size_t>(0, xStartNodes.size()),
3065 SweepingOp(SweepingOp::X_AXIS, xStartNodes, nodeConnectivity));
3066
3067 const size_t numLeafNodes = nodeConnectivity.size();
3068 72 const size_t numVoxels = numLeafNodes * FloatTreeT::LeafNodeType::SIZE;
3069
3070
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 std::unique_ptr<bool[]> changedNodeMaskA{new bool[numLeafNodes]};
3071
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 std::unique_ptr<bool[]> changedNodeMaskB{new bool[numLeafNodes]};
3072
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 std::unique_ptr<bool[]> changedVoxelMask{new bool[numVoxels]};
3073
3074
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 mesh_to_volume_internal::fillArray(changedNodeMaskA.get(), true, numLeafNodes);
3075
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 mesh_to_volume_internal::fillArray(changedNodeMaskB.get(), false, numLeafNodes);
3076
1/2
✓ Branch 1 taken 37 times.
✗ Branch 2 not taken.
 72 mesh_to_volume_internal::fillArray(changedVoxelMask.get(), false, numVoxels);
3077
3078 const tbb::blocked_range<size_t> nodeRange(0, numLeafNodes);
3079
3080 bool nodesUpdated = false;
3081
2/2
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 37 times.
 74 do {
3082 // Perform per leaf node localized propagation of signs by looping over
3083 // all voxels and checking to see if any of their neighbors (within the
3084 // same leaf) are negative
3085
2/4
✓ Branch 1 taken 39 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 39 times.
✗ Branch 5 not taken.
 74 tbb::parallel_for(nodeRange, mesh_to_volume_internal::SeedFillExteriorSign<FloatTreeT>(
3086 nodeConnectivity.nodes(), changedNodeMaskA.get()));
3087
3088 // For each leaf, check its axis aligned neighbors and propagate any changes
3089 // which occurred previously (in SeedFillExteriorSign OR in SyncVoxelMask) to
3090 // the leaf faces. Note that this operation stores the propagated face results
3091 // in a separate buffer (changedVoxelMask) to avoid writing to nodes being read
3092 // from other threads. Additionally mark any leaf nodes which will absorb any
3093 // changes from its neighbors in changedNodeMaskB
3094
1/2
✓ Branch 1 taken 39 times.
✗ Branch 2 not taken.
 74 tbb::parallel_for(nodeRange, mesh_to_volume_internal::SeedPoints<FloatTreeT>(
3095 nodeConnectivity, changedNodeMaskA.get(), changedNodeMaskB.get(),
3096 changedVoxelMask.get()));
3097
3098 // Only nodes where a value was influenced by an adjacent node need to be
3099 // processed on the next pass.
3100 changedNodeMaskA.swap(changedNodeMaskB);
3101
3102 nodesUpdated = false;
3103
2/2
✓ Branch 0 taken 12832 times.
✓ Branch 1 taken 37 times.
 24483 for (size_t n = 0; n < numLeafNodes; ++n) {
3104 24411 nodesUpdated |= changedNodeMaskA[n];
3105
2/2
✓ Branch 0 taken 12830 times.
✓ Branch 1 taken 2 times.
 24411 if (nodesUpdated) break;
3106 }
3107
3108 // Use the voxel mask updates in ::SeedPoints to actually assign the new values
3109 // across leaf node faces
3110
2/2
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 37 times.
 74 if (nodesUpdated) {
3111
1/4
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 2 tbb::parallel_for(nodeRange, mesh_to_volume_internal::SyncVoxelMask<FloatTreeT>(
3112 nodeConnectivity.nodes(), changedNodeMaskA.get(), changedVoxelMask.get()));
3113 }
3114 } while (nodesUpdated);
3115
3116 72 } // void traceExteriorBoundaries()
3117
3118
3119 ////////////////////////////////////////
3120
3121
3122 template <typename GridType, typename MeshDataAdapter, typename Interrupter>
3123 typename GridType::Ptr
3124 70 meshToVolume(
3125 Interrupter& interrupter,
3126 const MeshDataAdapter& mesh,
3127 const math::Transform& transform,
3128 float exteriorBandWidth,
3129 float interiorBandWidth,
3130 int flags,
3131 typename GridType::template ValueConverter<Int32>::Type * polygonIndexGrid)
3132 {
3133 using GridTypePtr = typename GridType::Ptr;
3134 using TreeType = typename GridType::TreeType;
3135 using LeafNodeType = typename TreeType::LeafNodeType;
3136 using ValueType = typename GridType::ValueType;
3137
3138 using Int32GridType = typename GridType::template ValueConverter<Int32>::Type;
3139 using Int32TreeType = typename Int32GridType::TreeType;
3140
3141 using BoolTreeType = typename TreeType::template ValueConverter<bool>::Type;
3142
3143 //////////
3144
3145 // Setup
3146
3147
3/6
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 35 times.
✗ Branch 8 not taken.
 140 GridTypePtr distGrid(new GridType(std::numeric_limits<ValueType>::max()));
3148
2/6
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 70 distGrid->setTransform(transform.copy());
3149
3150 70 ValueType exteriorWidth = ValueType(exteriorBandWidth);
3151
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 ValueType interiorWidth = ValueType(interiorBandWidth);
3152
3153 // Note: inf interior width is all right, this value makes the converter fill
3154 // interior regions with distance values.
3155
2/4
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 35 times.
 70 if (!std::isfinite(exteriorWidth) || std::isnan(interiorWidth)) {
3156 ✗ std::stringstream msg;
3157 ✗ msg << "Illegal narrow band width: exterior = " << exteriorWidth
3158 ✗ << ", interior = " << interiorWidth;
3159 OPENVDB_LOG_DEBUG(msg.str());
3160 return distGrid;
3161 }
3162
3163
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 const ValueType voxelSize = ValueType(transform.voxelSize()[0]);
3164
3165
2/4
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 35 times.
✗ Branch 3 not taken.
 70 if (!std::isfinite(voxelSize) || math::isZero(voxelSize)) {
3166 ✗ std::stringstream msg;
3167 msg << "Illegal transform, voxel size = " << voxelSize;
3168 OPENVDB_LOG_DEBUG(msg.str());
3169 return distGrid;
3170 }
3171
3172 // Convert narrow band width from voxel units to world space units.
3173 70 exteriorWidth *= voxelSize;
3174 // Avoid the unit conversion if the interior band width is set to
3175 // inf or std::numeric_limits<float>::max().
3176
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 if (interiorWidth < std::numeric_limits<ValueType>::max()) {
3177 70 interiorWidth *= voxelSize;
3178 }
3179
3180 70 const bool computeSignedDistanceField = (flags & UNSIGNED_DISTANCE_FIELD) == 0;
3181 70 const bool removeIntersectingVoxels = (flags & DISABLE_INTERSECTING_VOXEL_REMOVAL) == 0;
3182 const bool renormalizeValues = (flags & DISABLE_RENORMALIZATION) == 0;
3183 70 const bool trimNarrowBand = (flags & DISABLE_NARROW_BAND_TRIMMING) == 0;
3184
3185 Int32GridType* indexGrid = nullptr;
3186
3187 70 typename Int32GridType::Ptr temporaryIndexGrid;
3188
3189
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 if (polygonIndexGrid) {
3190 indexGrid = polygonIndexGrid;
3191 } else {
3192
3/6
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 35 times.
✗ Branch 8 not taken.
 70 temporaryIndexGrid.reset(new Int32GridType(Int32(util::INVALID_IDX)));
3193 indexGrid = temporaryIndexGrid.get();
3194 }
3195
3196
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 indexGrid->newTree();
3197
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
 70 indexGrid->setTransform(transform.copy());
3198
3199
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 if (computeSignedDistanceField) {
3200
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 distGrid->setGridClass(GRID_LEVEL_SET);
3201 } else {
3202 ✗ distGrid->setGridClass(GRID_UNKNOWN);
3203 ✗ interiorWidth = ValueType(0.0);
3204 }
3205
3206 TreeType& distTree = distGrid->tree();
3207 Int32TreeType& indexTree = indexGrid->tree();
3208
3209
3210 //////////
3211
3212 // Voxelize mesh
3213
3214 {
3215 using VoxelizationDataType = mesh_to_volume_internal::VoxelizationData<TreeType>;
3216 using DataTable = tbb::enumerable_thread_specific<typename VoxelizationDataType::Ptr>;
3217
3218
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 140 DataTable data;
3219 using Voxelizer =
3220 mesh_to_volume_internal::VoxelizePolygons<TreeType, MeshDataAdapter, Interrupter>;
3221
3222 const tbb::blocked_range<size_t> polygonRange(0, mesh.polygonCount());
3223
3224
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 tbb::parallel_for(polygonRange, Voxelizer(data, mesh, &interrupter));
3225
3226
2/2
✓ Branch 0 taken 68 times.
✓ Branch 1 taken 35 times.
 206 for (typename DataTable::iterator i = data.begin(); i != data.end(); ++i) {
3227 VoxelizationDataType& dataItem = **i;
3228 136 mesh_to_volume_internal::combineData(
3229
1/2
✓ Branch 1 taken 68 times.
✗ Branch 2 not taken.
 136 distTree, indexTree, dataItem.distTree, dataItem.indexTree);
3230 }
3231 }
3232
3233 // The progress estimates are based on the observed average time for a few different
3234 // test cases and is only intended to provide some rough progression feedback to the user.
3235
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 35 times.
✗ Branch 4 not taken.
 70 if (interrupter.wasInterrupted(30)) return distGrid;
3236
3237
3238 //////////
3239
3240 // Classify interior and exterior regions
3241
3242
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 if (computeSignedDistanceField) {
3243
3244 // Determines the inside/outside state for the narrow band of voxels.
3245
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 traceExteriorBoundaries(distTree);
3246
3247 std::vector<LeafNodeType*> nodes;
3248
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
 70 nodes.reserve(distTree.leafCount());
3249 distTree.getNodes(nodes);
3250
3251 const tbb::blocked_range<size_t> nodeRange(0, nodes.size());
3252
3253 using SignOp =
3254 mesh_to_volume_internal::ComputeIntersectingVoxelSign<TreeType, MeshDataAdapter>;
3255
3256
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 35 times.
✗ Branch 5 not taken.
 70 tbb::parallel_for(nodeRange, SignOp(nodes, distTree, indexTree, mesh));
3257
3258
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 35 times.
✗ Branch 4 not taken.
 70 if (interrupter.wasInterrupted(45)) return distGrid;
3259
3260 // Remove voxels created by self intersecting portions of the mesh.
3261
1/2
✓ Branch 0 taken 35 times.
✗ Branch 1 not taken.
 70 if (removeIntersectingVoxels) {
3262
3263
3/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✓ Branch 4 taken 1 times.
 70 tbb::parallel_for(nodeRange,
3264 mesh_to_volume_internal::ValidateIntersectingVoxels<TreeType>(distTree, nodes));
3265
3266
1/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 70 tbb::parallel_for(nodeRange,
3267 mesh_to_volume_internal::RemoveSelfIntersectingSurface<TreeType>(
3268 nodes, distTree, indexTree));
3269
3270
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 tools::pruneInactive(distTree, /*threading=*/true);
3271
1/2
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
 70 tools::pruneInactive(indexTree, /*threading=*/true);
3272 }
3273 }
3274
3275
2/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 35 times.
✗ Branch 4 not taken.
 70 if (interrupter.wasInterrupted(50)) return distGrid;
3276
3277
3/4
✓ Branch 1 taken 35 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 34 times.
 70 if (distTree.activeVoxelCount() == 0) {
3278
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 2 distTree.clear();
3279
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 2 distTree.root().setBackground(exteriorWidth, /*updateChildNodes=*/false);
3280 return distGrid;
3281 }
3282
3283 // Transform values (world space scaling etc.).
3284 {
3285 std::vector<LeafNodeType*> nodes;
3286
2/4
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 34 times.
✗ Branch 5 not taken.
 68 nodes.reserve(distTree.leafCount());
3287 distTree.getNodes(nodes);
3288
3289
2/6
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 68 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3290 mesh_to_volume_internal::TransformValues<TreeType>(
3291
1/2
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
 68 nodes, voxelSize, !computeSignedDistanceField));
3292 }
3293
3294 // Propagate sign information into tile regions.
3295
1/2
✓ Branch 0 taken 34 times.
✗ Branch 1 not taken.
 68 if (computeSignedDistanceField) {
3296
1/2
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
 68 distTree.root().setBackground(exteriorWidth, /*updateChildNodes=*/false);
3297
1/4
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 68 tools::signedFloodFillWithValues(distTree, exteriorWidth, -interiorWidth);
3298 } else {
3299 ✗ tools::changeBackground(distTree, exteriorWidth);
3300 }
3301
3302
2/4
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✗ Branch 4 not taken.
 68 if (interrupter.wasInterrupted(54)) return distGrid;
3303
3304
3305 //////////
3306
3307 // Expand the narrow band region
3308
3309 68 const ValueType minBandWidth = voxelSize * ValueType(2.0);
3310
3311
3/4
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 29 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 5 times.
 68 if (interiorWidth > minBandWidth || exteriorWidth > minBandWidth) {
3312
3313 // Create the initial voxel mask.
3314
1/2
✓ Branch 2 taken 29 times.
✗ Branch 3 not taken.
 116 BoolTreeType maskTree(false);
3315
3316 {
3317 std::vector<LeafNodeType*> nodes;
3318
2/4
✓ Branch 1 taken 29 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 29 times.
✗ Branch 5 not taken.
 58 nodes.reserve(distTree.leafCount());
3319 distTree.getNodes(nodes);
3320
3321 mesh_to_volume_internal::ConstructVoxelMask<TreeType> op(maskTree, distTree, nodes);
3322
1/2
✓ Branch 1 taken 29 times.
✗ Branch 2 not taken.
 58 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, nodes.size()), op);
3323 }
3324
3325 // Progress estimation
3326 unsigned maxIterations = std::numeric_limits<unsigned>::max();
3327
3328 float progress = 54.0f, step = 0.0f;
3329 58 double estimated =
3330
1/2
✓ Branch 0 taken 29 times.
✗ Branch 1 not taken.
 58 2.0 * std::ceil((std::max(interiorWidth, exteriorWidth) - minBandWidth) / voxelSize);
3331
3332
1/2
✓ Branch 0 taken 29 times.
✗ Branch 1 not taken.
 58 if (estimated < double(maxIterations)) {
3333 58 maxIterations = unsigned(estimated);
3334 58 step = 40.0f / float(maxIterations);
3335 }
3336
3337 std::vector<typename BoolTreeType::LeafNodeType*> maskNodes;
3338
3339 unsigned count = 0;
3340 72 while (true) {
3341
3342
2/6
✓ Branch 1 taken 65 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 65 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 130 if (interrupter.wasInterrupted(int(progress))) return distGrid;
3343
3344 130 const size_t maskNodeCount = maskTree.leafCount();
3345
2/2
✓ Branch 0 taken 40 times.
✓ Branch 1 taken 25 times.
 130 if (maskNodeCount == 0) break;
3346
3347 maskNodes.clear();
3348
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 maskNodes.reserve(maskNodeCount);
3349 maskTree.getNodes(maskNodes);
3350
3351 const tbb::blocked_range<size_t> range(0, maskNodes.size());
3352
3353
1/4
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 80 tbb::parallel_for(range,
3354 mesh_to_volume_internal::DiffLeafNodeMask<TreeType>(distTree, maskNodes));
3355
3356
1/2
✓ Branch 1 taken 40 times.
✗ Branch 2 not taken.
 80 mesh_to_volume_internal::expandNarrowband(distTree, indexTree, maskTree, maskNodes,
3357 mesh, exteriorWidth, interiorWidth, voxelSize);
3358
3359
2/2
✓ Branch 0 taken 36 times.
✓ Branch 1 taken 4 times.
 80 if ((++count) >= maxIterations) break;
3360 72 progress += step;
3361 }
3362 }
3363
3364
2/4
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✗ Branch 4 not taken.
 68 if (interrupter.wasInterrupted(94)) return distGrid;
3365
3366
2/4
✓ Branch 0 taken 34 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 34 times.
✗ Branch 4 not taken.
 68 if (!polygonIndexGrid) indexGrid->clear();
3367
3368
3369 /////////
3370
3371 // Renormalize distances to smooth out bumps caused by self intersecting
3372 // and overlapping portions of the mesh and renormalize the level set.
3373
3374
2/2
✓ Branch 0 taken 30 times.
✓ Branch 1 taken 4 times.
 68 if (computeSignedDistanceField && renormalizeValues) {
3375
3376 std::vector<LeafNodeType*> nodes;
3377
2/4
✓ Branch 1 taken 30 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 30 times.
✗ Branch 5 not taken.
 60 nodes.reserve(distTree.leafCount());
3378 distTree.getNodes(nodes);
3379
3380
2/4
✓ Branch 0 taken 30 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 30 times.
✗ Branch 4 not taken.
 60 std::unique_ptr<ValueType[]> buffer{new ValueType[LeafNodeType::SIZE * nodes.size()]};
3381
3382 60 const ValueType offset = ValueType(0.8 * voxelSize);
3383
3384
3/6
✓ Branch 0 taken 30 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 30 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 30 times.
✗ Branch 6 not taken.
 120 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3385 mesh_to_volume_internal::OffsetValues<TreeType>(nodes, -offset));
3386
3387
2/4
✓ Branch 1 taken 30 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 30 times.
✗ Branch 4 not taken.
 60 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3388 mesh_to_volume_internal::Renormalize<TreeType>(
3389 distTree, nodes, buffer.get(), voxelSize));
3390
3391
1/2
✓ Branch 1 taken 30 times.
✗ Branch 2 not taken.
 60 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3392 mesh_to_volume_internal::MinCombine<TreeType>(nodes, buffer.get()));
3393
3394
3/8
✓ Branch 0 taken 30 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 30 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 30 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
 120 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3395 mesh_to_volume_internal::OffsetValues<TreeType>(
3396 nodes, offset - mesh_to_volume_internal::Tolerance<ValueType>::epsilon()));
3397 }
3398
3399
2/4
✓ Branch 1 taken 34 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34 times.
✗ Branch 4 not taken.
 68 if (interrupter.wasInterrupted(99)) return distGrid;
3400
3401
3402 /////////
3403
3404 // Remove active voxels that exceed the narrow band limits
3405
3406
3/4
✓ Branch 0 taken 34 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32 times.
✓ Branch 3 taken 2 times.
 136 if (trimNarrowBand && std::min(interiorWidth, exteriorWidth) < voxelSize * ValueType(4.0)) {
3407
3408 std::vector<LeafNodeType*> nodes;
3409
2/4
✓ Branch 1 taken 32 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 32 times.
✗ Branch 5 not taken.
 64 nodes.reserve(distTree.leafCount());
3410 distTree.getNodes(nodes);
3411
3412
3/6
✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 32 times.
✗ Branch 6 not taken.
 128 tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
3413 mesh_to_volume_internal::InactivateValues<TreeType>(
3414 nodes, exteriorWidth, computeSignedDistanceField ? interiorWidth : exteriorWidth));
3415
3416 64 tools::pruneLevelSet(
3417
3/8
✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 32 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 32 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
 64 distTree, exteriorWidth, computeSignedDistanceField ? -interiorWidth : -exteriorWidth);
3418 }
3419
3420 return distGrid;
3421 }
3422
3423
3424 template <typename GridType, typename MeshDataAdapter>
3425 typename GridType::Ptr
3426 20 meshToVolume(
3427 const MeshDataAdapter& mesh,
3428 const math::Transform& transform,
3429 float exteriorBandWidth,
3430 float interiorBandWidth,
3431 int flags,
3432 typename GridType::template ValueConverter<Int32>::Type * polygonIndexGrid)
3433 {
3434 20 util::NullInterrupter nullInterrupter;
3435 return meshToVolume<GridType>(nullInterrupter, mesh, transform,
3436 20 exteriorBandWidth, interiorBandWidth, flags, polygonIndexGrid);
3437 }
3438
3439
3440 ////////////////////////////////////////
3441
3442
3443 //{
3444 /// @cond OPENVDB_DOCS_INTERNAL
3445
3446 /// @internal This overload is enabled only for grids with a scalar, floating-point ValueType.
3447 template<typename GridType, typename Interrupter>
3448 inline typename std::enable_if<std::is_floating_point<typename GridType::ValueType>::value,
3449 typename GridType::Ptr>::type
3450
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 15 times.
 30 doMeshConversion(
3451 Interrupter& interrupter,
3452 const openvdb::math::Transform& xform,
3453 const std::vector<Vec3s>& points,
3454 const std::vector<Vec3I>& triangles,
3455 const std::vector<Vec4I>& quads,
3456 float exBandWidth,
3457 float inBandWidth,
3458 bool unsignedDistanceField = false)
3459 {
3460
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 15 times.
 30 if (points.empty()) {
3461 ✗ return typename GridType::Ptr(new GridType(typename GridType::ValueType(exBandWidth)));
3462 }
3463
3464 const size_t numPoints = points.size();
3465
3/4
✓ Branch 0 taken 15 times.
✗ Branch 1 not taken.
✓ Branch 4 taken 142 times.
✓ Branch 5 taken 15 times.
 314 std::unique_ptr<Vec3s[]> indexSpacePoints{new Vec3s[numPoints]};
3466
3467 // transform points to local grid index space
3468
1/4
✓ Branch 1 taken 15 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 30 tbb::parallel_for(tbb::blocked_range<size_t>(0, numPoints),
3469 mesh_to_volume_internal::TransformPoints<Vec3s>(
3470 &points[0], indexSpacePoints.get(), xform));
3471
3472
1/2
✓ Branch 0 taken 15 times.
✗ Branch 1 not taken.
 30 const int conversionFlags = unsignedDistanceField ? UNSIGNED_DISTANCE_FIELD : 0;
3473
3474
2/2
✓ Branch 0 taken 3 times.
✓ Branch 1 taken 12 times.
 30 if (quads.empty()) {
3475
3476 QuadAndTriangleDataAdapter<Vec3s, Vec3I>
3477 mesh(indexSpacePoints.get(), numPoints, &triangles[0], triangles.size());
3478
3479 return meshToVolume<GridType>(
3480
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 6 interrupter, mesh, xform, exBandWidth, inBandWidth, conversionFlags);
3481
3482
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
 24 } else if (triangles.empty()) {
3483
3484 QuadAndTriangleDataAdapter<Vec3s, Vec4I>
3485 mesh(indexSpacePoints.get(), numPoints, &quads[0], quads.size());
3486
3487 return meshToVolume<GridType>(
3488
1/2
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
 20 interrupter, mesh, xform, exBandWidth, inBandWidth, conversionFlags);
3489 }
3490
3491 // pack primitives
3492
3493 4 const size_t numPrimitives = triangles.size() + quads.size();
3494
4/6
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 48 times.
✓ Branch 7 taken 2 times.
 100 std::unique_ptr<Vec4I[]> prims{new Vec4I[numPrimitives]};
3495
3496
2/2
✓ Branch 0 taken 24 times.
✓ Branch 1 taken 2 times.
 52 for (size_t n = 0, N = triangles.size(); n < N; ++n) {
3497 const Vec3I& triangle = triangles[n];
3498 Vec4I& prim = prims[n];
3499 48 prim[0] = triangle[0];
3500 48 prim[1] = triangle[1];
3501 48 prim[2] = triangle[2];
3502 48 prim[3] = util::INVALID_IDX;
3503 }
3504
3505 const size_t offset = triangles.size();
3506
2/2
✓ Branch 0 taken 24 times.
✓ Branch 1 taken 2 times.
 52 for (size_t n = 0, N = quads.size(); n < N; ++n) {
3507 48 prims[offset + n] = quads[n];
3508 }
3509
3510 QuadAndTriangleDataAdapter<Vec3s, Vec4I>
3511 mesh(indexSpacePoints.get(), numPoints, prims.get(), numPrimitives);
3512
3513 return meshToVolume<GridType>(interrupter, mesh, xform,
3514
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
 4 exBandWidth, inBandWidth, conversionFlags);
3515 }
3516
3517
3518 /// @internal This overload is enabled only for grids that do not have a scalar,
3519 /// floating-point ValueType.
3520 template<typename GridType, typename Interrupter>
3521 inline typename std::enable_if<!std::is_floating_point<typename GridType::ValueType>::value,
3522 typename GridType::Ptr>::type
3523 doMeshConversion(
3524 Interrupter&,
3525 const math::Transform& /*xform*/,
3526 const std::vector<Vec3s>& /*points*/,
3527 const std::vector<Vec3I>& /*triangles*/,
3528 const std::vector<Vec4I>& /*quads*/,
3529 float /*exBandWidth*/,
3530 float /*inBandWidth*/,
3531 bool /*unsignedDistanceField*/ = false)
3532 {
3533 OPENVDB_THROW(TypeError,
3534 "mesh to volume conversion is supported only for scalar floating-point grids");
3535 }
3536
3537 /// @endcond
3538 //}
3539
3540
3541 ////////////////////////////////////////
3542
3543
3544 template<typename GridType>
3545 typename GridType::Ptr
3546 meshToLevelSet(
3547 const openvdb::math::Transform& xform,
3548 const std::vector<Vec3s>& points,
3549 const std::vector<Vec3I>& triangles,
3550 float halfWidth)
3551 {
3552 util::NullInterrupter nullInterrupter;
3553 return meshToLevelSet<GridType>(nullInterrupter, xform, points, triangles, halfWidth);
3554 }
3555
3556
3557 template<typename GridType, typename Interrupter>
3558 typename GridType::Ptr
3559 ✗ meshToLevelSet(
3560 Interrupter& interrupter,
3561 const openvdb::math::Transform& xform,
3562 const std::vector<Vec3s>& points,
3563 const std::vector<Vec3I>& triangles,
3564 float halfWidth)
3565 {
3566 ✗ std::vector<Vec4I> quads(0);
3567 return doMeshConversion<GridType>(interrupter, xform, points, triangles, quads,
3568 ✗ halfWidth, halfWidth);
3569 }
3570
3571
3572 template<typename GridType>
3573 typename GridType::Ptr
3574 1 meshToLevelSet(
3575 const openvdb::math::Transform& xform,
3576 const std::vector<Vec3s>& points,
3577 const std::vector<Vec4I>& quads,
3578 float halfWidth)
3579 {
3580 1 util::NullInterrupter nullInterrupter;
3581 1 return meshToLevelSet<GridType>(nullInterrupter, xform, points, quads, halfWidth);
3582 }
3583
3584
3585 template<typename GridType, typename Interrupter>
3586 typename GridType::Ptr
3587 2 meshToLevelSet(
3588 Interrupter& interrupter,
3589 const openvdb::math::Transform& xform,
3590 const std::vector<Vec3s>& points,
3591 const std::vector<Vec4I>& quads,
3592 float halfWidth)
3593 {
3594 2 std::vector<Vec3I> triangles(0);
3595 return doMeshConversion<GridType>(interrupter, xform, points, triangles, quads,
3596
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 4 halfWidth, halfWidth);
3597 }
3598
3599
3600 template<typename GridType>
3601 typename GridType::Ptr
3602 meshToLevelSet(
3603 const openvdb::math::Transform& xform,
3604 const std::vector<Vec3s>& points,
3605 const std::vector<Vec3I>& triangles,
3606 const std::vector<Vec4I>& quads,
3607 float halfWidth)
3608 {
3609 util::NullInterrupter nullInterrupter;
3610 return meshToLevelSet<GridType>(
3611 nullInterrupter, xform, points, triangles, quads, halfWidth);
3612 }
3613
3614
3615 template<typename GridType, typename Interrupter>
3616 typename GridType::Ptr
3617 28 meshToLevelSet(
3618 Interrupter& interrupter,
3619 const openvdb::math::Transform& xform,
3620 const std::vector<Vec3s>& points,
3621 const std::vector<Vec3I>& triangles,
3622 const std::vector<Vec4I>& quads,
3623 float halfWidth)
3624 {
3625 return doMeshConversion<GridType>(interrupter, xform, points, triangles, quads,
3626 28 halfWidth, halfWidth);
3627 }
3628
3629
3630 template<typename GridType>
3631 typename GridType::Ptr
3632 meshToSignedDistanceField(
3633 const openvdb::math::Transform& xform,
3634 const std::vector<Vec3s>& points,
3635 const std::vector<Vec3I>& triangles,
3636 const std::vector<Vec4I>& quads,
3637 float exBandWidth,
3638 float inBandWidth)
3639 {
3640 util::NullInterrupter nullInterrupter;
3641 return meshToSignedDistanceField<GridType>(
3642 nullInterrupter, xform, points, triangles, quads, exBandWidth, inBandWidth);
3643 }
3644
3645
3646 template<typename GridType, typename Interrupter>
3647 typename GridType::Ptr
3648 ✗ meshToSignedDistanceField(
3649 Interrupter& interrupter,
3650 const openvdb::math::Transform& xform,
3651 const std::vector<Vec3s>& points,
3652 const std::vector<Vec3I>& triangles,
3653 const std::vector<Vec4I>& quads,
3654 float exBandWidth,
3655 float inBandWidth)
3656 {
3657 return doMeshConversion<GridType>(interrupter, xform, points, triangles,
3658 ✗ quads, exBandWidth, inBandWidth);
3659 }
3660
3661
3662 template<typename GridType>
3663 typename GridType::Ptr
3664 meshToUnsignedDistanceField(
3665 const openvdb::math::Transform& xform,
3666 const std::vector<Vec3s>& points,
3667 const std::vector<Vec3I>& triangles,
3668 const std::vector<Vec4I>& quads,
3669 float bandWidth)
3670 {
3671 util::NullInterrupter nullInterrupter;
3672 return meshToUnsignedDistanceField<GridType>(
3673 nullInterrupter, xform, points, triangles, quads, bandWidth);
3674 }
3675
3676
3677 template<typename GridType, typename Interrupter>
3678 typename GridType::Ptr
3679 ✗ meshToUnsignedDistanceField(
3680 Interrupter& interrupter,
3681 const openvdb::math::Transform& xform,
3682 const std::vector<Vec3s>& points,
3683 const std::vector<Vec3I>& triangles,
3684 const std::vector<Vec4I>& quads,
3685 float bandWidth)
3686 {
3687 return doMeshConversion<GridType>(interrupter, xform, points, triangles, quads,
3688 ✗ bandWidth, bandWidth, true);
3689 }
3690
3691
3692 ////////////////////////////////////////////////////////////////////////////////
3693
3694
3695 // Required by several of the tree nodes
3696 inline std::ostream&
3697 operator<<(std::ostream& ostr, const MeshToVoxelEdgeData::EdgeData& rhs)
3698 {
3699 ostr << "{[ " << rhs.mXPrim << ", " << rhs.mXDist << "]";
3700 ostr << " [ " << rhs.mYPrim << ", " << rhs.mYDist << "]";
3701 ostr << " [ " << rhs.mZPrim << ", " << rhs.mZDist << "]}";
3702 return ostr;
3703 }
3704
3705 // Required by math::Abs
3706 inline MeshToVoxelEdgeData::EdgeData
3707 Abs(const MeshToVoxelEdgeData::EdgeData& x)
3708 {
3709 return x;
3710 }
3711
3712
3713 ////////////////////////////////////////
3714
3715
3716 class MeshToVoxelEdgeData::GenEdgeData
3717 {
3718 public:
3719
3720 GenEdgeData(
3721 const std::vector<Vec3s>& pointList,
3722 const std::vector<Vec4I>& polygonList);
3723
3724 void run(bool threaded = true);
3725
3726 GenEdgeData(GenEdgeData& rhs, tbb::split);
3727 inline void operator() (const tbb::blocked_range<size_t> &range);
3728 inline void join(GenEdgeData& rhs);
3729
3730 inline TreeType& tree() { return mTree; }
3731
3732 private:
3733 void operator=(const GenEdgeData&) {}
3734
3735 struct Primitive { Vec3d a, b, c, d; Int32 index; };
3736
3737 template<bool IsQuad>
3738 inline void voxelize(const Primitive&);
3739
3740 template<bool IsQuad>
3741 inline bool evalPrimitive(const Coord&, const Primitive&);
3742
3743 inline bool rayTriangleIntersection( const Vec3d& origin, const Vec3d& dir,
3744 const Vec3d& a, const Vec3d& b, const Vec3d& c, double& t);
3745
3746
3747 TreeType mTree;
3748 Accessor mAccessor;
3749
3750 const std::vector<Vec3s>& mPointList;
3751 const std::vector<Vec4I>& mPolygonList;
3752
3753 // Used internally for acceleration
3754 using IntTreeT = TreeType::ValueConverter<Int32>::Type;
3755 IntTreeT mLastPrimTree;
3756 tree::ValueAccessor<IntTreeT> mLastPrimAccessor;
3757 }; // class MeshToVoxelEdgeData::GenEdgeData
3758
3759
3760 inline
3761 MeshToVoxelEdgeData::GenEdgeData::GenEdgeData(
3762 const std::vector<Vec3s>& pointList,
3763 const std::vector<Vec4I>& polygonList)
3764 : mTree(EdgeData())
3765 , mAccessor(mTree)
3766 , mPointList(pointList)
3767 , mPolygonList(polygonList)
3768 , mLastPrimTree(Int32(util::INVALID_IDX))
3769 , mLastPrimAccessor(mLastPrimTree)
3770 {
3771 }
3772
3773
3774 inline
3775 MeshToVoxelEdgeData::GenEdgeData::GenEdgeData(GenEdgeData& rhs, tbb::split)
3776 : mTree(EdgeData())
3777 , mAccessor(mTree)
3778 , mPointList(rhs.mPointList)
3779 , mPolygonList(rhs.mPolygonList)
3780 , mLastPrimTree(Int32(util::INVALID_IDX))
3781 , mLastPrimAccessor(mLastPrimTree)
3782 {
3783 }
3784
3785
3786 inline void
3787 MeshToVoxelEdgeData::GenEdgeData::run(bool threaded)
3788 {
3789 if (threaded) {
3790 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, mPolygonList.size()), *this);
3791 } else {
3792 (*this)(tbb::blocked_range<size_t>(0, mPolygonList.size()));
3793 }
3794 }
3795
3796
3797 inline void
3798 MeshToVoxelEdgeData::GenEdgeData::join(GenEdgeData& rhs)
3799 {
3800 using RootNodeType = TreeType::RootNodeType;
3801 using NodeChainType = RootNodeType::NodeChainType;
3802 static_assert(NodeChainType::Size > 1, "expected tree height > 1");
3803 using InternalNodeType = typename NodeChainType::template Get<1>;
3804
3805 Coord ijk;
3806 Index offset;
3807
3808 rhs.mTree.clearAllAccessors();
3809
3810 TreeType::LeafIter leafIt = rhs.mTree.beginLeaf();
3811 for ( ; leafIt; ++leafIt) {
3812 ijk = leafIt->origin();
3813
3814 TreeType::LeafNodeType* lhsLeafPt = mTree.probeLeaf(ijk);
3815
3816 if (!lhsLeafPt) {
3817
3818 mAccessor.addLeaf(rhs.mAccessor.probeLeaf(ijk));
3819 InternalNodeType* node = rhs.mAccessor.getNode<InternalNodeType>();
3820 node->stealNode<TreeType::LeafNodeType>(ijk, EdgeData(), false);
3821 rhs.mAccessor.clear();
3822
3823 } else {
3824
3825 TreeType::LeafNodeType::ValueOnCIter it = leafIt->cbeginValueOn();
3826 for ( ; it; ++it) {
3827
3828 offset = it.pos();
3829 const EdgeData& rhsValue = it.getValue();
3830
3831 if (!lhsLeafPt->isValueOn(offset)) {
3832 lhsLeafPt->setValueOn(offset, rhsValue);
3833 } else {
3834
3835 EdgeData& lhsValue = const_cast<EdgeData&>(lhsLeafPt->getValue(offset));
3836
3837 if (rhsValue.mXDist < lhsValue.mXDist) {
3838 lhsValue.mXDist = rhsValue.mXDist;
3839 lhsValue.mXPrim = rhsValue.mXPrim;
3840 }
3841
3842 if (rhsValue.mYDist < lhsValue.mYDist) {
3843 lhsValue.mYDist = rhsValue.mYDist;
3844 lhsValue.mYPrim = rhsValue.mYPrim;
3845 }
3846
3847 if (rhsValue.mZDist < lhsValue.mZDist) {
3848 lhsValue.mZDist = rhsValue.mZDist;
3849 lhsValue.mZPrim = rhsValue.mZPrim;
3850 }
3851
3852 }
3853 } // end value iteration
3854 }
3855 } // end leaf iteration
3856 }
3857
3858
3859 inline void
3860 MeshToVoxelEdgeData::GenEdgeData::operator()(const tbb::blocked_range<size_t> &range)
3861 {
3862 Primitive prim;
3863
3864 for (size_t n = range.begin(); n < range.end(); ++n) {
3865
3866 const Vec4I& verts = mPolygonList[n];
3867
3868 prim.index = Int32(n);
3869 prim.a = Vec3d(mPointList[verts[0]]);
3870 prim.b = Vec3d(mPointList[verts[1]]);
3871 prim.c = Vec3d(mPointList[verts[2]]);
3872
3873 if (util::INVALID_IDX != verts[3]) {
3874 prim.d = Vec3d(mPointList[verts[3]]);
3875 voxelize<true>(prim);
3876 } else {
3877 voxelize<false>(prim);
3878 }
3879 }
3880 }
3881
3882
3883 template<bool IsQuad>
3884 inline void
3885 MeshToVoxelEdgeData::GenEdgeData::voxelize(const Primitive& prim)
3886 {
3887 std::deque<Coord> coordList;
3888 Coord ijk, nijk;
3889
3890 ijk = Coord::floor(prim.a);
3891 coordList.push_back(ijk);
3892
3893 evalPrimitive<IsQuad>(ijk, prim);
3894
3895 while (!coordList.empty()) {
3896
3897 ijk = coordList.back();
3898 coordList.pop_back();
3899
3900 for (Int32 i = 0; i < 26; ++i) {
3901 nijk = ijk + util::COORD_OFFSETS[i];
3902
3903 if (prim.index != mLastPrimAccessor.getValue(nijk)) {
3904 mLastPrimAccessor.setValue(nijk, prim.index);
3905 if(evalPrimitive<IsQuad>(nijk, prim)) coordList.push_back(nijk);
3906 }
3907 }
3908 }
3909 }
3910
3911
3912 template<bool IsQuad>
3913 inline bool
3914 MeshToVoxelEdgeData::GenEdgeData::evalPrimitive(const Coord& ijk, const Primitive& prim)
3915 {
3916 Vec3d uvw, org(ijk[0], ijk[1], ijk[2]);
3917 bool intersecting = false;
3918 double t;
3919
3920 EdgeData edgeData;
3921 mAccessor.probeValue(ijk, edgeData);
3922
3923 // Evaluate first triangle
3924 double dist = (org -
3925 closestPointOnTriangleToPoint(prim.a, prim.c, prim.b, org, uvw)).lengthSqr();
3926
3927 if (rayTriangleIntersection(org, Vec3d(1.0, 0.0, 0.0), prim.a, prim.c, prim.b, t)) {
3928 if (t < edgeData.mXDist) {
3929 edgeData.mXDist = float(t);
3930 edgeData.mXPrim = prim.index;
3931 intersecting = true;
3932 }
3933 }
3934
3935 if (rayTriangleIntersection(org, Vec3d(0.0, 1.0, 0.0), prim.a, prim.c, prim.b, t)) {
3936 if (t < edgeData.mYDist) {
3937 edgeData.mYDist = float(t);
3938 edgeData.mYPrim = prim.index;
3939 intersecting = true;
3940 }
3941 }
3942
3943 if (rayTriangleIntersection(org, Vec3d(0.0, 0.0, 1.0), prim.a, prim.c, prim.b, t)) {
3944 if (t < edgeData.mZDist) {
3945 edgeData.mZDist = float(t);
3946 edgeData.mZPrim = prim.index;
3947 intersecting = true;
3948 }
3949 }
3950
3951 if (IsQuad) {
3952 // Split quad into a second triangle and calculate distance.
3953 double secondDist = (org -
3954 closestPointOnTriangleToPoint(prim.a, prim.d, prim.c, org, uvw)).lengthSqr();
3955
3956 if (secondDist < dist) dist = secondDist;
3957
3958 if (rayTriangleIntersection(org, Vec3d(1.0, 0.0, 0.0), prim.a, prim.d, prim.c, t)) {
3959 if (t < edgeData.mXDist) {
3960 edgeData.mXDist = float(t);
3961 edgeData.mXPrim = prim.index;
3962 intersecting = true;
3963 }
3964 }
3965
3966 if (rayTriangleIntersection(org, Vec3d(0.0, 1.0, 0.0), prim.a, prim.d, prim.c, t)) {
3967 if (t < edgeData.mYDist) {
3968 edgeData.mYDist = float(t);
3969 edgeData.mYPrim = prim.index;
3970 intersecting = true;
3971 }
3972 }
3973
3974 if (rayTriangleIntersection(org, Vec3d(0.0, 0.0, 1.0), prim.a, prim.d, prim.c, t)) {
3975 if (t < edgeData.mZDist) {
3976 edgeData.mZDist = float(t);
3977 edgeData.mZPrim = prim.index;
3978 intersecting = true;
3979 }
3980 }
3981 }
3982
3983 if (intersecting) mAccessor.setValue(ijk, edgeData);
3984
3985 return (dist < 0.86602540378443861);
3986 }
3987
3988
3989 inline bool
3990 MeshToVoxelEdgeData::GenEdgeData::rayTriangleIntersection(
3991 const Vec3d& origin, const Vec3d& dir,
3992 const Vec3d& a, const Vec3d& b, const Vec3d& c,
3993 double& t)
3994 {
3995 // Check if ray is parallel with triangle
3996
3997 Vec3d e1 = b - a;
3998 Vec3d e2 = c - a;
3999 Vec3d s1 = dir.cross(e2);
4000
4001 double divisor = s1.dot(e1);
4002 if (!(std::abs(divisor) > 0.0)) return false;
4003
4004 // Compute barycentric coordinates
4005
4006 double inv_divisor = 1.0 / divisor;
4007 Vec3d d = origin - a;
4008 double b1 = d.dot(s1) * inv_divisor;
4009
4010 if (b1 < 0.0 || b1 > 1.0) return false;
4011
4012 Vec3d s2 = d.cross(e1);
4013 double b2 = dir.dot(s2) * inv_divisor;
4014
4015 if (b2 < 0.0 || (b1 + b2) > 1.0) return false;
4016
4017 // Compute distance to intersection point
4018
4019 t = e2.dot(s2) * inv_divisor;
4020 return (t < 0.0) ? false : true;
4021 }
4022
4023
4024 ////////////////////////////////////////
4025
4026
4027 inline
4028 MeshToVoxelEdgeData::MeshToVoxelEdgeData()
4029 : mTree(EdgeData())
4030 {
4031 }
4032
4033
4034 inline void
4035 MeshToVoxelEdgeData::convert(
4036 const std::vector<Vec3s>& pointList,
4037 const std::vector<Vec4I>& polygonList)
4038 {
4039 GenEdgeData converter(pointList, polygonList);
4040 converter.run();
4041
4042 mTree.clear();
4043 mTree.merge(converter.tree());
4044 }
4045
4046
4047 inline void
4048 MeshToVoxelEdgeData::getEdgeData(
4049 Accessor& acc,
4050 const Coord& ijk,
4051 std::vector<Vec3d>& points,
4052 std::vector<Index32>& primitives)
4053 {
4054 EdgeData data;
4055 Vec3d point;
4056
4057 Coord coord = ijk;
4058
4059 if (acc.probeValue(coord, data)) {
4060
4061 if (data.mXPrim != util::INVALID_IDX) {
4062 point[0] = double(coord[0]) + data.mXDist;
4063 point[1] = double(coord[1]);
4064 point[2] = double(coord[2]);
4065
4066 points.push_back(point);
4067 primitives.push_back(data.mXPrim);
4068 }
4069
4070 if (data.mYPrim != util::INVALID_IDX) {
4071 point[0] = double(coord[0]);
4072 point[1] = double(coord[1]) + data.mYDist;
4073 point[2] = double(coord[2]);
4074
4075 points.push_back(point);
4076 primitives.push_back(data.mYPrim);
4077 }
4078
4079 if (data.mZPrim != util::INVALID_IDX) {
4080 point[0] = double(coord[0]);
4081 point[1] = double(coord[1]);
4082 point[2] = double(coord[2]) + data.mZDist;
4083
4084 points.push_back(point);
4085 primitives.push_back(data.mZPrim);
4086 }
4087
4088 }
4089
4090 coord[0] += 1;
4091
4092 if (acc.probeValue(coord, data)) {
4093
4094 if (data.mYPrim != util::INVALID_IDX) {
4095 point[0] = double(coord[0]);
4096 point[1] = double(coord[1]) + data.mYDist;
4097 point[2] = double(coord[2]);
4098
4099 points.push_back(point);
4100 primitives.push_back(data.mYPrim);
4101 }
4102
4103 if (data.mZPrim != util::INVALID_IDX) {
4104 point[0] = double(coord[0]);
4105 point[1] = double(coord[1]);
4106 point[2] = double(coord[2]) + data.mZDist;
4107
4108 points.push_back(point);
4109 primitives.push_back(data.mZPrim);
4110 }
4111 }
4112
4113 coord[2] += 1;
4114
4115 if (acc.probeValue(coord, data)) {
4116 if (data.mYPrim != util::INVALID_IDX) {
4117 point[0] = double(coord[0]);
4118 point[1] = double(coord[1]) + data.mYDist;
4119 point[2] = double(coord[2]);
4120
4121 points.push_back(point);
4122 primitives.push_back(data.mYPrim);
4123 }
4124 }
4125
4126 coord[0] -= 1;
4127
4128 if (acc.probeValue(coord, data)) {
4129
4130 if (data.mXPrim != util::INVALID_IDX) {
4131 point[0] = double(coord[0]) + data.mXDist;
4132 point[1] = double(coord[1]);
4133 point[2] = double(coord[2]);
4134
4135 points.push_back(point);
4136 primitives.push_back(data.mXPrim);
4137 }
4138
4139 if (data.mYPrim != util::INVALID_IDX) {
4140 point[0] = double(coord[0]);
4141 point[1] = double(coord[1]) + data.mYDist;
4142 point[2] = double(coord[2]);
4143
4144 points.push_back(point);
4145 primitives.push_back(data.mYPrim);
4146 }
4147 }
4148
4149
4150 coord[1] += 1;
4151
4152 if (acc.probeValue(coord, data)) {
4153
4154 if (data.mXPrim != util::INVALID_IDX) {
4155 point[0] = double(coord[0]) + data.mXDist;
4156 point[1] = double(coord[1]);
4157 point[2] = double(coord[2]);
4158
4159 points.push_back(point);
4160 primitives.push_back(data.mXPrim);
4161 }
4162 }
4163
4164 coord[2] -= 1;
4165
4166 if (acc.probeValue(coord, data)) {
4167
4168 if (data.mXPrim != util::INVALID_IDX) {
4169 point[0] = double(coord[0]) + data.mXDist;
4170 point[1] = double(coord[1]);
4171 point[2] = double(coord[2]);
4172
4173 points.push_back(point);
4174 primitives.push_back(data.mXPrim);
4175 }
4176
4177 if (data.mZPrim != util::INVALID_IDX) {
4178 point[0] = double(coord[0]);
4179 point[1] = double(coord[1]);
4180 point[2] = double(coord[2]) + data.mZDist;
4181
4182 points.push_back(point);
4183 primitives.push_back(data.mZPrim);
4184 }
4185 }
4186
4187 coord[0] += 1;
4188
4189 if (acc.probeValue(coord, data)) {
4190
4191 if (data.mZPrim != util::INVALID_IDX) {
4192 point[0] = double(coord[0]);
4193 point[1] = double(coord[1]);
4194 point[2] = double(coord[2]) + data.mZDist;
4195
4196 points.push_back(point);
4197 primitives.push_back(data.mZPrim);
4198 }
4199 }
4200 }
4201
4202
4203 template<typename GridType, typename VecType>
4204 typename GridType::Ptr
4205 14 createLevelSetBox(const math::BBox<VecType>& bbox,
4206 const openvdb::math::Transform& xform,
4207 typename VecType::ValueType halfWidth)
4208 {
4209 14 const Vec3s pmin = Vec3s(xform.worldToIndex(bbox.min()));
4210 14 const Vec3s pmax = Vec3s(xform.worldToIndex(bbox.max()));
4211
4212 126 Vec3s points[8];
4213 14 points[0] = Vec3s(pmin[0], pmin[1], pmin[2]);
4214 14 points[1] = Vec3s(pmin[0], pmin[1], pmax[2]);
4215 14 points[2] = Vec3s(pmax[0], pmin[1], pmax[2]);
4216 14 points[3] = Vec3s(pmax[0], pmin[1], pmin[2]);
4217 14 points[4] = Vec3s(pmin[0], pmax[1], pmin[2]);
4218 14 points[5] = Vec3s(pmin[0], pmax[1], pmax[2]);
4219 14 points[6] = Vec3s(pmax[0], pmax[1], pmax[2]);
4220 14 points[7] = Vec3s(pmax[0], pmax[1], pmin[2]);
4221
4222 98 Vec4I faces[6];
4223 14 faces[0] = Vec4I(0, 1, 2, 3); // bottom
4224 14 faces[1] = Vec4I(7, 6, 5, 4); // top
4225 14 faces[2] = Vec4I(4, 5, 1, 0); // front
4226 14 faces[3] = Vec4I(6, 7, 3, 2); // back
4227 14 faces[4] = Vec4I(0, 3, 7, 4); // left
4228 14 faces[5] = Vec4I(1, 5, 6, 2); // right
4229
4230 QuadAndTriangleDataAdapter<Vec3s, Vec4I> mesh(points, 8, faces, 6);
4231
4232 14 return meshToVolume<GridType>(mesh, xform, static_cast<float>(halfWidth), static_cast<float>(halfWidth));
4233 }
4234
4235
4236 ////////////////////////////////////////
4237
4238
4239 // Explicit Template Instantiation
4240
4241 #ifdef OPENVDB_USE_EXPLICIT_INSTANTIATION
4242
4243 #ifdef OPENVDB_INSTANTIATE_MESHTOVOLUME
4244 #include <openvdb/util/ExplicitInstantiation.h>
4245 #endif
4246
4247 #define _FUNCTION(TreeT) \
4248 Grid<TreeT>::Ptr meshToVolume<Grid<TreeT>>(util::NullInterrupter&, \
4249 const QuadAndTriangleDataAdapter<Vec3s, Vec3I>&, const openvdb::math::Transform&, \
4250 float, float, int, Grid<TreeT>::ValueConverter<Int32>::Type*)
4251 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4252 #undef _FUNCTION
4253
4254 #define _FUNCTION(TreeT) \
4255 Grid<TreeT>::Ptr meshToVolume<Grid<TreeT>>(util::NullInterrupter&, \
4256 const QuadAndTriangleDataAdapter<Vec3s, Vec4I>&, const openvdb::math::Transform&, \
4257 float, float, int, Grid<TreeT>::ValueConverter<Int32>::Type*)
4258 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4259 #undef _FUNCTION
4260
4261 #define _FUNCTION(TreeT) \
4262 Grid<TreeT>::Ptr meshToLevelSet<Grid<TreeT>>(util::NullInterrupter&, \
4263 const openvdb::math::Transform&, const std::vector<Vec3s>&, const std::vector<Vec3I>&, \
4264 float)
4265 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4266 #undef _FUNCTION
4267
4268 #define _FUNCTION(TreeT) \
4269 Grid<TreeT>::Ptr meshToLevelSet<Grid<TreeT>>(util::NullInterrupter&, \
4270 const openvdb::math::Transform&, const std::vector<Vec3s>&, const std::vector<Vec4I>&, \
4271 float)
4272 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4273 #undef _FUNCTION
4274
4275 #define _FUNCTION(TreeT) \
4276 Grid<TreeT>::Ptr meshToLevelSet<Grid<TreeT>>(util::NullInterrupter&, \
4277 const openvdb::math::Transform&, const std::vector<Vec3s>&, \
4278 const std::vector<Vec3I>&, const std::vector<Vec4I>&, float)
4279 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4280 #undef _FUNCTION
4281
4282 #define _FUNCTION(TreeT) \
4283 Grid<TreeT>::Ptr meshToSignedDistanceField<Grid<TreeT>>(util::NullInterrupter&, \
4284 const openvdb::math::Transform&, const std::vector<Vec3s>&, \
4285 const std::vector<Vec3I>&, const std::vector<Vec4I>&, float, float)
4286 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4287 #undef _FUNCTION
4288
4289 #define _FUNCTION(TreeT) \
4290 Grid<TreeT>::Ptr meshToUnsignedDistanceField<Grid<TreeT>>(util::NullInterrupter&, \
4291 const openvdb::math::Transform&, const std::vector<Vec3s>&, \
4292 const std::vector<Vec3I>&, const std::vector<Vec4I>&, float)
4293 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4294 #undef _FUNCTION
4295
4296 #define _FUNCTION(TreeT) \
4297 Grid<TreeT>::Ptr createLevelSetBox<Grid<TreeT>>(const math::BBox<Vec3s>&, \
4298 const openvdb::math::Transform&, float)
4299 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4300 #undef _FUNCTION
4301
4302 #define _FUNCTION(TreeT) \
4303 Grid<TreeT>::Ptr createLevelSetBox<Grid<TreeT>>(const math::BBox<Vec3d>&, \
4304 const openvdb::math::Transform&, double)
4305 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4306 #undef _FUNCTION
4307
4308 #define _FUNCTION(TreeT) \
4309 void traceExteriorBoundaries(TreeT&)
4310 OPENVDB_REAL_TREE_INSTANTIATE(_FUNCTION)
4311 #undef _FUNCTION
4312
4313 #endif // OPENVDB_USE_EXPLICIT_INSTANTIATION
4314
4315
4316 } // namespace tools
4317 } // namespace OPENVDB_VERSION_NAME
4318 } // namespace openvdb
4319
4320 #endif // OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
4321