OpenVDB  8.2.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 #include "SignedFloodFill.h"
18 
19 #include <openvdb/Grid.h>
20 #include <openvdb/Types.h>
21 #include <openvdb/math/FiniteDifference.h> // for math::BiasedGradientScheme
23 #include <tbb/task_group.h>
24 #include <algorithm> // for std::min(), std::max()
25 #include <vector>
26 
27 
28 namespace openvdb {
30 namespace OPENVDB_VERSION_NAME {
31 namespace tools {
32 
33 
34 /// @brief Compute the narrow-band signed distance to the interface between
35 /// active and inactive voxels in the input grid.
36 ///
37 /// @return A shared pointer to a new sdf / level set grid of type @c float
38 ///
39 /// @param grid Input grid of arbitrary type whose active voxels are used
40 /// in constructing the level set.
41 /// @param halfWidth Half the width of the narrow band in voxel units.
42 /// @param closingSteps Number of morphological closing steps used to fill gaps
43 /// in the active voxel region.
44 /// @param dilation Number of voxels to expand the active voxel region.
45 /// @param smoothingSteps Number of smoothing interations.
46 template<typename GridT>
47 inline typename GridT::template ValueConverter<float>::Type::Ptr
48 topologyToLevelSet(const GridT& grid, int halfWidth = 3, int closingSteps = 1, int dilation = 0,
49  int smoothingSteps = 0);
50 
51 
52 /// @brief Compute the narrow-band signed distance to the interface between
53 /// active and inactive voxels in the input grid.
54 ///
55 /// @return A shared pointer to a new sdf / level set grid of type @c float
56 ///
57 /// @param grid Input grid of arbitrary type whose active voxels are used
58 /// in constructing the level set.
59 /// @param halfWidth Half the width of the narrow band in voxel units.
60 /// @param closingSteps Number of morphological closing steps used to fill gaps
61 /// in the active voxel region.
62 /// @param dilation Number of voxels to expand the active voxel region.
63 /// @param smoothingSteps Number of smoothing interations.
64 /// @param interrupt Optional object adhering to the util::NullInterrupter interface.
65 template<typename GridT, typename InterrupterT>
66 inline typename GridT::template ValueConverter<float>::Type::Ptr
67 topologyToLevelSet(const GridT& grid, int halfWidth = 3, int closingSteps = 1, int dilation = 0,
68  int smoothingSteps = 0, InterrupterT* interrupt = nullptr);
69 
70 
71 ////////////////////////////////////////
72 
73 /// @cond OPENVDB_DOCS_INTERNAL
74 
75 namespace ttls_internal {
76 
77 
78 template<typename TreeT>
79 struct DilateOp
80 {
81  DilateOp(TreeT& t, int n) : tree(&t), size(n) {}
82  void operator()() const {
84  }
85  TreeT* tree;
86  const int size;
87 };
88 
89 
90 template<typename TreeT>
91 struct ErodeOp
92 {
93  ErodeOp(TreeT& t, int n) : tree(&t), size(n) {}
94  void operator()() const {
95  tools::erodeActiveValues(*tree, /*iterations=*/size, tools::NN_FACE, tools::IGNORE_TILES);
96  tools::pruneInactive(*tree);
97  }
98  TreeT* tree;
99  const int size;
100 };
101 
102 
103 template<typename TreeType>
104 struct OffsetAndMinComp
105 {
106  using LeafNodeType = typename TreeType::LeafNodeType;
107  using ValueType = typename TreeType::ValueType;
108 
109  OffsetAndMinComp(std::vector<LeafNodeType*>& lhsNodes,
110  const TreeType& rhsTree, ValueType offset)
111  : mLhsNodes(lhsNodes.empty() ? nullptr : &lhsNodes[0]), mRhsTree(&rhsTree), mOffset(offset)
112  {
113  }
114 
115  void operator()(const tbb::blocked_range<size_t>& range) const
116  {
117  using Iterator = typename LeafNodeType::ValueOnIter;
118 
119  tree::ValueAccessor<const TreeType> rhsAcc(*mRhsTree);
120  const ValueType offset = mOffset;
121 
122  for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
123 
124  LeafNodeType& lhsNode = *mLhsNodes[n];
125  const LeafNodeType * rhsNodePt = rhsAcc.probeConstLeaf(lhsNode.origin());
126  if (!rhsNodePt) continue;
127 
128  for (Iterator it = lhsNode.beginValueOn(); it; ++it) {
129  ValueType& val = const_cast<ValueType&>(it.getValue());
130  val = std::min(val, offset + rhsNodePt->getValue(it.pos()));
131  }
132  }
133  }
134 
135 private:
136  LeafNodeType * * const mLhsNodes;
137  TreeType const * const mRhsTree;
138  ValueType const mOffset;
139 }; // struct OffsetAndMinComp
140 
141 
142 template<typename GridType, typename InterrupterType>
143 inline void
144 normalizeLevelSet(GridType& grid, const int halfWidthInVoxels, InterrupterType* interrupt = nullptr)
145 {
146  LevelSetFilter<GridType, GridType, InterrupterType> filter(grid, interrupt);
147  filter.setSpatialScheme(math::FIRST_BIAS);
148  filter.setNormCount(halfWidthInVoxels);
149  filter.normalize();
150  filter.prune();
151 }
152 
153 
154 template<typename GridType, typename InterrupterType>
155 inline void
156 smoothLevelSet(GridType& grid, int iterations, int halfBandWidthInVoxels,
157  InterrupterType* interrupt = nullptr)
158 {
159  using ValueType = typename GridType::ValueType;
160  using TreeType = typename GridType::TreeType;
161  using LeafNodeType = typename TreeType::LeafNodeType;
162 
163  GridType filterGrid(grid);
164 
165  LevelSetFilter<GridType, GridType, InterrupterType> filter(filterGrid, interrupt);
166  filter.setSpatialScheme(math::FIRST_BIAS);
167 
168  for (int n = 0; n < iterations; ++n) {
169  if (interrupt && interrupt->wasInterrupted()) break;
170  filter.mean(1);
171  }
172 
173  std::vector<LeafNodeType*> nodes;
174  grid.tree().getNodes(nodes);
175 
176  const ValueType offset = ValueType(double(0.5) * grid.transform().voxelSize()[0]);
177 
178  tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
179  OffsetAndMinComp<TreeType>(nodes, filterGrid.tree(), -offset));
180 
181  // Clean up any damanage that was done by the min operation
182  normalizeLevelSet(grid, halfBandWidthInVoxels, interrupt);
183 }
184 
185 
186 } // namespace ttls_internal
187 
188 /// @endcond
189 
190 template<typename GridT, typename InterrupterT>
191 inline typename GridT::template ValueConverter<float>::Type::Ptr
192 topologyToLevelSet(const GridT& grid, int halfWidth, int closingSteps, int dilation,
193  int smoothingSteps, InterrupterT* interrupt)
194 {
195  using MaskTreeT = typename GridT::TreeType::template ValueConverter<ValueMask>::Type;
196  using FloatTreeT = typename GridT::TreeType::template ValueConverter<float>::Type;
197  using FloatGridT = Grid<FloatTreeT>;
198 
199  // Check inputs
200 
201  halfWidth = std::max(halfWidth, 1);
202  closingSteps = std::max(closingSteps, 0);
203  dilation = std::max(dilation, 0);
204 
205  if (!grid.hasUniformVoxels()) {
206  OPENVDB_THROW(ValueError, "Non-uniform voxels are not supported!");
207  }
208 
209  // Copy the topology into a MaskGrid.
210  MaskTreeT maskTree( grid.tree(), false/*background*/, openvdb::TopologyCopy() );
211 
212  // Morphological closing operation.
213  tools::dilateActiveValues(maskTree, closingSteps + dilation, tools::NN_FACE, tools::IGNORE_TILES);
214  tools::erodeActiveValues(maskTree, /*iterations=*/closingSteps, tools::NN_FACE, tools::IGNORE_TILES);
215  tools::pruneInactive(maskTree);
216 
217  // Generate a volume with an implicit zero crossing at the boundary
218  // between active and inactive values in the input grid.
219  const float background = float(grid.voxelSize()[0]) * float(halfWidth);
220  typename FloatTreeT::Ptr lsTree(
221  new FloatTreeT( maskTree, /*out=*/background, /*in=*/-background, openvdb::TopologyCopy() ) );
222 
223  tbb::task_group pool;
224  pool.run( ttls_internal::ErodeOp< MaskTreeT >( maskTree, halfWidth ) );
225  pool.run( ttls_internal::DilateOp<FloatTreeT>( *lsTree , halfWidth ) );
226  pool.wait();// wait for both tasks to complete
227 
228  lsTree->topologyDifference( maskTree );
229  tools::pruneLevelSet( *lsTree, /*threading=*/true);
230 
231  // Create a level set grid from the tree
232  typename FloatGridT::Ptr lsGrid = FloatGridT::create( lsTree );
233  lsGrid->setTransform( grid.transform().copy() );
234  lsGrid->setGridClass( openvdb::GRID_LEVEL_SET );
235 
236  // Use a PDE based scheme to propagate distance values from the
237  // implicit zero crossing.
238  ttls_internal::normalizeLevelSet(*lsGrid, 3*halfWidth, interrupt);
239 
240  // Additional filtering
241  if (smoothingSteps > 0) {
242  ttls_internal::smoothLevelSet(*lsGrid, smoothingSteps, halfWidth, interrupt);
243  }
244 
245  return lsGrid;
246 }
247 
248 
249 template<typename GridT>
250 inline typename GridT::template ValueConverter<float>::Type::Ptr
251 topologyToLevelSet(const GridT& grid, int halfWidth, int closingSteps, int dilation, int smoothingSteps)
252 {
253  util::NullInterrupter interrupt;
254  return topologyToLevelSet(grid, halfWidth, closingSteps, dilation, smoothingSteps, &interrupt);
255 }
256 
257 
258 } // namespace tools
259 } // namespace OPENVDB_VERSION_NAME
260 } // namespace openvdb
261 
262 #endif // OPENVDB_TOOLS_TOPOLOGY_TO_LEVELSET_HAS_BEEN_INCLUDED
263 
Performs various types of level set deformations with interface tracking. These unrestricted deformat...
Implementation of morphological dilation and erosion.
Propagate the signs of distance values from the active voxels in the narrow band to the inactive valu...
Container class that associates a tree with a transform and metadata.
Definition: Grid.h:577
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:564
Definition: Exceptions.h:65
@ FIRST_BIAS
Definition: FiniteDifference.h:167
const std::enable_if<!VecTraits< T >::IsVec, T >::type & min(const T &a, const T &b)
Definition: Composite.h:102
const std::enable_if<!VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:106
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:192
@ NN_FACE
Definition: Morphology.h:56
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:389
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:1130
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:1053
@ IGNORE_TILES
Definition: Morphology.h:78
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:354
@ GRID_LEVEL_SET
Definition: Types.h:337
Definition: Exceptions.h:13
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:74
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:26
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h.in:116
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h.in:180