HDK
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 {
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 
void parallel_for(int64_t start, int64_t end, std::function< void(int64_t index)> &&task, parallel_options opt=parallel_options(0, Split_Y, 1))
Definition: parallel.h:127
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.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:1130
Dummy NOOP interrupter class defining interface.
GLdouble GLdouble t
Definition: glew.h:1403
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
GLsizeiptr size
Definition: glcorearb.h:663
SYS_FORCE_INLINE const_iterator end() const
ImageBuf OIIO_API min(Image_or_Const A, Image_or_Const B, ROI roi={}, int nthreads=0)
ImageBuf OIIO_API max(Image_or_Const A, Image_or_Const B, ROI roi={}, int nthreads=0)
GLint GLint GLint GLint GLint GLint GLint GLbitfield GLenum filter
Definition: glcorearb.h:1296
Container class that associates a tree with a transform and metadata.
Definition: Grid.h:28
Propagate the signs of distance values from the active voxels in the narrow band to the inactive valu...
GLdouble n
Definition: glcorearb.h:2007
GLuint GLfloat * val
Definition: glcorearb.h:1607
GridT::template ValueConverter< float >::Type::Ptr topologyToLevelSet(const GridT &grid, int halfWidth=3, int closingSteps=1, int dilation=0, int smoothingSteps=0)
Compute the narrow-band signed distance to the interface between active and inactive voxels in the in...
GA_API const UT_StringHolder N
Implementation of morphological dilation and erosion.
GLenum GLint * range
Definition: glcorearb.h:1924
Performs various types of level set deformations with interface tracking. These unrestricted deformat...
GLintptr offset
Definition: glcorearb.h:664
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
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:114
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
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:74
**Note that the tasks the is the thread number *for the pool
Definition: thread.h:643