OpenVDB  11.0.0
TopologyToLevelSet.h
Go to the documentation of this file.
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3 //
4 /// @file TopologyToLevelSet.h
5 ///
6 /// @brief This tool generates a narrow-band signed distance field / level set
7 /// from the interface between active and inactive voxels in a vdb grid.
8 ///
9 /// @par Example:
10 /// Combine with @c tools::PointsToVolume for fast point cloud to level set conversion.
11 
12 #ifndef OPENVDB_TOOLS_TOPOLOGY_TO_LEVELSET_HAS_BEEN_INCLUDED
13 #define OPENVDB_TOOLS_TOPOLOGY_TO_LEVELSET_HAS_BEEN_INCLUDED
14 
15 #include "LevelSetFilter.h"
16 #include "Morphology.h" // for erodeActiveValues and dilateActiveValues
17 
18 #include <openvdb/Grid.h>
19 #include <openvdb/Types.h>
20 #include <openvdb/math/FiniteDifference.h> // for math::BiasedGradientScheme
22 #include <openvdb/openvdb.h>
24 
25 #include <tbb/task_group.h>
26 
27 #include <algorithm> // for std::min(), std::max()
28 #include <vector>
29 
30 
31 namespace openvdb {
33 namespace OPENVDB_VERSION_NAME {
34 namespace tools {
35 
36 
37 /// @brief Compute the narrow-band signed distance to the interface between
38 /// active and inactive voxels in the input grid.
39 ///
40 /// @return A shared pointer to a new sdf / level set grid of type @c float
41 ///
42 /// @param grid Input grid of arbitrary type whose active voxels are used
43 /// in constructing the level set.
44 /// @param halfWidth Half the width of the narrow band in voxel units.
45 /// @param closingSteps Number of morphological closing steps used to fill gaps
46 /// in the active voxel region.
47 /// @param dilation Number of voxels to expand the active voxel region.
48 /// @param smoothingSteps Number of smoothing interations.
49 template<typename GridT>
50 typename GridT::template ValueConverter<float>::Type::Ptr
51 topologyToLevelSet(const GridT& grid,
52  int halfWidth = 3,
53  int closingSteps = 1,
54  int dilation = 0,
55  int smoothingSteps = 0);
56 
57 
58 /// @brief Compute the narrow-band signed distance to the interface between
59 /// active and inactive voxels in the input grid.
60 ///
61 /// @return A shared pointer to a new sdf / level set grid of type @c float
62 ///
63 /// @param grid Input grid of arbitrary type whose active voxels are used
64 /// in constructing the level set.
65 /// @param halfWidth Half the width of the narrow band in voxel units.
66 /// @param closingSteps Number of morphological closing steps used to fill gaps
67 /// in the active voxel region.
68 /// @param dilation Number of voxels to expand the active voxel region.
69 /// @param smoothingSteps Number of smoothing interations.
70 /// @param interrupt Optional object adhering to the util::NullInterrupter interface.
71 template<typename GridT, typename InterrupterT>
72 typename GridT::template ValueConverter<float>::Type::Ptr
73 topologyToLevelSet(const GridT& grid,
74  int halfWidth = 3,
75  int closingSteps = 1,
76  int dilation = 0,
77  int smoothingSteps = 0,
78  InterrupterT* interrupt = nullptr);
79 
80 
81 ////////////////////////////////////////
82 
83 /// @cond OPENVDB_DOCS_INTERNAL
84 
85 namespace ttls_internal {
86 
87 
88 template<typename TreeType>
89 struct OffsetAndMinComp
90 {
91  using LeafNodeType = typename TreeType::LeafNodeType;
92  using ValueType = typename TreeType::ValueType;
93 
94  OffsetAndMinComp(std::vector<LeafNodeType*>& lhsNodes,
95  const TreeType& rhsTree,
96  ValueType offset)
97  : mLhsNodes(lhsNodes)
98  , mRhsTree(rhsTree)
99  , mOffset(offset) {}
100 
101  void operator()(const tbb::blocked_range<size_t>& range) const
102  {
103  using Iterator = typename LeafNodeType::ValueOnIter;
104 
105  tree::ValueAccessor<const TreeType> rhsAcc(mRhsTree);
106  const ValueType offset = mOffset;
107 
108  for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
109 
110  LeafNodeType& lhsNode = *mLhsNodes[n];
111  const LeafNodeType* rhsNodePt = rhsAcc.probeConstLeaf(lhsNode.origin());
112  if (!rhsNodePt) continue;
113 
114  auto* const data = lhsNode.buffer().data();
115  for (Iterator it = lhsNode.beginValueOn(); it; ++it) {
116  ValueType& val = data[it.pos()];
117  val = std::min(val, offset + rhsNodePt->getValue(it.pos()));
118  }
119  }
120  }
121 
122 private:
123  std::vector<LeafNodeType*>& mLhsNodes;
124  const TreeType& mRhsTree;
125  const ValueType mOffset;
126 }; // struct OffsetAndMinComp
127 
128 
129 template<typename GridType, typename InterrupterType>
130 inline void
131 normalizeLevelSet(GridType& grid, const int halfWidthInVoxels, InterrupterType* interrupt = nullptr)
132 {
133  LevelSetFilter<GridType, GridType, InterrupterType> filter(grid, interrupt);
134  filter.setSpatialScheme(math::FIRST_BIAS);
135  filter.setNormCount(halfWidthInVoxels);
136  filter.normalize();
137  filter.prune();
138 }
139 
140 
141 template<typename GridType, typename InterrupterType>
142 inline void
143 smoothLevelSet(GridType& grid, int iterations, int halfBandWidthInVoxels,
144  InterrupterType* interrupt = nullptr)
145 {
146  using ValueType = typename GridType::ValueType;
147  using TreeType = typename GridType::TreeType;
148  using LeafNodeType = typename TreeType::LeafNodeType;
149 
150  GridType filterGrid(grid);
151 
152  LevelSetFilter<GridType, GridType, InterrupterType> filter(filterGrid, interrupt);
153  filter.setSpatialScheme(math::FIRST_BIAS);
154 
155  for (int n = 0; n < iterations; ++n) {
156  if (interrupt && interrupt->wasInterrupted()) break;
157  filter.mean(1);
158  }
159 
160  std::vector<LeafNodeType*> nodes;
161  grid.tree().getNodes(nodes);
162 
163  const ValueType offset = ValueType(double(0.5) * grid.transform().voxelSize()[0]);
164 
165  tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
166  OffsetAndMinComp<TreeType>(nodes, filterGrid.tree(), -offset));
167 
168  // Clean up any damanage that was done by the min operation
169  normalizeLevelSet(grid, halfBandWidthInVoxels, interrupt);
170 }
171 
172 } // namespace ttls_internal
173 
174 /// @endcond
175 
176 template<typename GridT, typename InterrupterT>
177 typename GridT::template ValueConverter<float>::Type::Ptr
178 topologyToLevelSet(const GridT& grid, int halfWidth, int closingSteps, int dilation,
179  int smoothingSteps, InterrupterT* interrupt)
180 {
181  using MaskTreeT = typename GridT::TreeType::template ValueConverter<ValueMask>::Type;
182  using FloatTreeT = typename GridT::TreeType::template ValueConverter<float>::Type;
183  using FloatGridT = Grid<FloatTreeT>;
184 
185  // Check inputs
186 
187  halfWidth = std::max(halfWidth, 1);
188  closingSteps = std::max(closingSteps, 0);
189  dilation = std::max(dilation, 0);
190 
191  if (!grid.hasUniformVoxels()) {
192  OPENVDB_THROW(ValueError, "Non-uniform voxels are not supported!");
193  }
194 
195  // Copy the topology into a MaskGrid.
196  MaskTreeT maskTree(grid.tree(), false/*background*/, openvdb::TopologyCopy());
197 
198  // Morphological closing operation.
199  tools::dilateActiveValues(maskTree, closingSteps + dilation, tools::NN_FACE, tools::PRESERVE_TILES);
200  tools::erodeActiveValues(maskTree, /*iterations=*/closingSteps, tools::NN_FACE, tools::PRESERVE_TILES);
201  tools::pruneInactive(maskTree);
202 
203  // Generate a volume with an implicit zero crossing at the boundary
204  // between active and inactive values in the input grid.
205  const float background = float(grid.voxelSize()[0]) * float(halfWidth);
206  typename FloatTreeT::Ptr lsTree(
207  new FloatTreeT(maskTree, /*out=*/background, /*in=*/-background, openvdb::TopologyCopy()));
208 
209  tbb::task_group pool;
210  pool.run([&]() {
211  tools::erodeActiveValues(maskTree, /*iterations=*/halfWidth, tools::NN_FACE, tools::PRESERVE_TILES);
212  tools::pruneInactive(maskTree);
213  });
214  pool.run([&]() {
215  tools::dilateActiveValues(*lsTree, /*iterations=*/halfWidth, tools::NN_FACE, tools::PRESERVE_TILES);
216  });
217  pool.wait();// wait for both tasks to complete
218 
219  lsTree->topologyDifference(maskTree);
220  maskTree.clear();
221 
222  lsTree->voxelizeActiveTiles(); // has to come before ::pruneLevelSet
223  tools::pruneLevelSet(*lsTree, /*threading=*/true);
224 
225  // Create a level set grid from the tree
226  typename FloatGridT::Ptr lsGrid = FloatGridT::create(lsTree);
227  lsGrid->setTransform(grid.transform().copy());
228  lsGrid->setGridClass(openvdb::GRID_LEVEL_SET);
229 
230  // Use a PDE based scheme to propagate distance values from the
231  // implicit zero crossing.
232  ttls_internal::normalizeLevelSet(*lsGrid, 3*halfWidth, interrupt);
233 
234  // Additional filtering
235  if (smoothingSteps > 0) {
236  ttls_internal::smoothLevelSet(*lsGrid, smoothingSteps, halfWidth, interrupt);
237  }
238 
239  return lsGrid;
240 }
241 
242 
243 template<typename GridT>
244 typename GridT::template ValueConverter<float>::Type::Ptr
245 topologyToLevelSet(const GridT& grid, int halfWidth, int closingSteps, int dilation, int smoothingSteps)
246 {
247  util::NullInterrupter interrupt;
248  return topologyToLevelSet(grid, halfWidth, closingSteps, dilation, smoothingSteps, &interrupt);
249 }
250 
251 
252 ////////////////////////////////////////
253 
254 
255 // Explicit Template Instantiation
256 
257 #ifdef OPENVDB_USE_EXPLICIT_INSTANTIATION
258 
259 #ifdef OPENVDB_INSTANTIATE_TOPOLOGYTOLEVELSET
261 #endif
262 
263 #define _FUNCTION(TreeT) \
264  Grid<TreeT>::ValueConverter<float>::Type::Ptr topologyToLevelSet(const Grid<TreeT>&, int, int, int, int, \
265  util::NullInterrupter*)
267 #undef _FUNCTION
268 
269 #endif // OPENVDB_USE_EXPLICIT_INSTANTIATION
270 
271 
272 } // namespace tools
273 } // namespace OPENVDB_VERSION_NAME
274 } // namespace openvdb
275 
276 #endif // OPENVDB_TOOLS_TOPOLOGY_TO_LEVELSET_HAS_BEEN_INCLUDED
277 
Definition: FiniteDifference.h:166
Definition: Exceptions.h:65
void pruneLevelSet(TreeT &tree, bool threaded=true, size_t grainSize=1)
Reduce the memory footprint of a tree by replacing nodes whose values are all inactive with inactive ...
Definition: Prune.h:390
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:74
Base class for interrupters.
Definition: NullInterrupter.h:25
Definition: Morphology.h:81
GridT::template ValueConverter< float >::Type::Ptr topologyToLevelSet(const GridT &grid, int halfWidth=3, int closingSteps=1, int dilation=0, int smoothingSteps=0, InterrupterT *interrupt=nullptr)
Compute the narrow-band signed distance to the interface between active and inactive voxels in the in...
Definition: TopologyToLevelSet.h:178
void erodeActiveValues(TreeOrLeafManagerT &tree, const int iterations=1, const NearestNeighbors nn=NN_FACE, const TilePolicy mode=PRESERVE_TILES, const bool threaded=true)
Topologically erode all active values (i.e. both voxels and tiles) in a tree using one of three neare...
Definition: Morphology.h:1133
Definition: Morphology.h:59
Performs various types of level set deformations with interface tracking. These unrestricted deformat...
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:683
#define OPENVDB_ALL_TREE_INSTANTIATE(Function)
Definition: version.h.in:161
Definition: Types.h:455
Container class that associates a tree with a transform and metadata.
Definition: Grid.h:28
void dilateActiveValues(TreeOrLeafManagerT &tree, const int iterations=1, const NearestNeighbors nn=NN_FACE, const TilePolicy mode=PRESERVE_TILES, const bool threaded=true)
Topologically dilate all active values (i.e. both voxels and tiles) in a tree using one of three near...
Definition: Morphology.h:1056
Implementation of morphological dilation and erosion.
Definition: Exceptions.h:13
GridType
List of types that are currently supported by NanoVDB.
Definition: NanoVDB.h:317
Attribute-owned data structure for points. Point attributes are stored in leaf nodes and ordered by v...
void pruneInactive(TreeT &tree, bool threaded=true, size_t grainSize=1)
Reduce the memory footprint of a tree by replacing with background tiles any nodes whose values are a...
Definition: Prune.h:355
const std::enable_if<!VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:110
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h.in:121
const std::enable_if<!VecTraits< T >::IsVec, T >::type & min(const T &a, const T &b)
Definition: Composite.h:106
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h.in:212