OpenVDB  11.0.0
Filter.h
Go to the documentation of this file.
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3 //
4 /// @author Ken Museth
5 ///
6 /// @file tools/Filter.h
7 ///
8 /// @brief Filtering of VDB volumes. All operations can optionally be masked
9 /// with another grid that acts as an alpha-mask. By default, filtering
10 /// operations do not modify the topology of the input tree and thus do
11 /// not process active tiles. However Filter::setProcessTiles can be
12 /// used to process active tiles, densifying them on demand when necessary.
13 
14 #ifndef OPENVDB_TOOLS_FILTER_HAS_BEEN_INCLUDED
15 #define OPENVDB_TOOLS_FILTER_HAS_BEEN_INCLUDED
16 
17 #include <openvdb/Types.h>
18 #include <openvdb/Grid.h>
19 #include <openvdb/math/Math.h>
20 #include <openvdb/math/Stencils.h>
21 #include <openvdb/math/Transform.h>
25 #include <openvdb/util/Util.h>
26 #include <openvdb/thread/Threading.h>
27 #include "Interpolation.h"
28 
29 #include <tbb/parallel_for.h>
30 #include <tbb/concurrent_vector.h>
31 
32 #include <algorithm> // for std::max()
33 #include <functional>
34 #include <type_traits>
35 
36 namespace openvdb {
38 namespace OPENVDB_VERSION_NAME {
39 namespace tools {
40 
41 /// @brief Volume filtering (e.g., diffusion) with optional alpha masking
42 template<typename GridT,
43  typename MaskT = typename GridT::template ValueConverter<float>::Type,
44  typename InterruptT = util::NullInterrupter>
45 class Filter
46 {
47 public:
48  using GridType = GridT;
49  using MaskType = MaskT;
50  using TreeType = typename GridType::TreeType;
51  using LeafType = typename TreeType::LeafNodeType;
52  using ValueType = typename GridType::ValueType;
53  using AlphaType = typename MaskType::ValueType;
55  using RangeType = typename LeafManagerType::LeafRange;
56  using BufferType = typename LeafManagerType::BufferType;
57  static_assert(std::is_floating_point<AlphaType>::value,
58  "openvdb::tools::Filter requires a mask grid with floating-point values");
59 
60  /// Constructor
61  /// @param grid Grid to be filtered.
62  /// @param interrupt Optional interrupter.
63  Filter(GridT& grid, InterruptT* interrupt = nullptr)
64  : mGrid(&grid)
65  , mTask(nullptr)
66  , mInterrupter(interrupt)
67  , mMask(nullptr)
68  , mGrainSize(1)
69  , mMinMask(0)
70  , mMaxMask(1)
71  , mInvertMask(false)
72  , mTiles(false) {}
73 
74  /// @brief Shallow copy constructor called by tbb::parallel_for()
75  /// threads during filtering.
76  /// @param other The other Filter from which to copy.
77  Filter(const Filter& other)
78  : mGrid(other.mGrid)
79  , mTask(other.mTask)
80  , mInterrupter(other.mInterrupter)
81  , mMask(other.mMask)
82  , mGrainSize(other.mGrainSize)
83  , mMinMask(other.mMinMask)
84  , mMaxMask(other.mMaxMask)
85  , mInvertMask(other.mInvertMask)
86  , mTiles(other.mTiles) {}
87 
88  /// @return the grain-size used for multi-threading
89  int getGrainSize() const { return mGrainSize; }
90  /// @brief Set the grain-size used for multi-threading.
91  /// @note A grain size of 0 or less disables multi-threading!
92  void setGrainSize(int grainsize) { mGrainSize = grainsize; }
93 
94  /// @return whether active tiles are being processed
95  bool getProcessTiles() const { return mTiles; }
96  /// @brief Set whether active tiles should also be processed.
97  /// @note If true, some tiles may become voxelized
98  /// @warning If using with a mask, ensure that the mask topology matches the
99  /// tile topology of the filter grid as tiles will not respect overlapping
100  /// mask values at tree levels finer than themselves e.g. a leaf level tile
101  /// will only use the corresponding tile ijk value in the mask grid
102  void setProcessTiles(bool flag) { mTiles = flag; }
103 
104  /// @brief Return the minimum value of the mask to be used for the
105  /// derivation of a smooth alpha value.
106  AlphaType minMask() const { return mMinMask; }
107  /// @brief Return the maximum value of the mask to be used for the
108  /// derivation of a smooth alpha value.
109  AlphaType maxMask() const { return mMaxMask; }
110  /// @brief Define the range for the (optional) scalar mask.
111  /// @param min Minimum value of the range.
112  /// @param max Maximum value of the range.
113  /// @details Mask values outside the range are clamped to zero or one, and
114  /// values inside the range map smoothly to 0->1 (unless the mask is inverted).
115  /// @throw ValueError if @a min is not smaller than @a max.
117  {
118  if (!(min < max)) OPENVDB_THROW(ValueError, "Invalid mask range (expects min < max)");
119  mMinMask = min;
120  mMaxMask = max;
121  }
122 
123  /// @brief Return true if the mask is inverted, i.e. min->max in the
124  /// original mask maps to 1->0 in the inverted alpha mask.
125  bool isMaskInverted() const { return mInvertMask; }
126  /// @brief Invert the optional mask, i.e. min->max in the original
127  /// mask maps to 1->0 in the inverted alpha mask.
128  void invertMask(bool invert=true) { mInvertMask = invert; }
129 
130  /// @brief One iteration of a fast separable mean-value (i.e. box) filter.
131  /// @param width The width of the mean-value filter is 2*width+1 voxels.
132  /// @param iterations Number of times the mean-value filter is applied.
133  /// @param mask Optional alpha mask.
134  void mean(int width = 1, int iterations = 1, const MaskType* mask = nullptr);
135 
136  /// @brief One iteration of a fast separable Gaussian filter.
137  ///
138  /// @note This is approximated as 4 iterations of a separable mean filter
139  /// which typically leads an approximation that's better than 95%!
140  /// @param width The width of the mean-value filter is 2*width+1 voxels.
141  /// @param iterations Number of times the mean-value filter is applied.
142  /// @param mask Optional alpha mask.
143  void gaussian(int width = 1, int iterations = 1, const MaskType* mask = nullptr);
144 
145  /// @brief One iteration of a median-value filter
146  ///
147  /// @note This filter is not separable and is hence relatively slow!
148  /// @param width The width of the mean-value filter is 2*width+1 voxels.
149  /// @param iterations Number of times the mean-value filter is applied.
150  /// @param mask Optional alpha mask.
151  void median(int width = 1, int iterations = 1, const MaskType* mask = nullptr);
152 
153  /// Offsets (i.e. adds) a constant value to all active voxels.
154  /// @param offset Offset in the same units as the grid.
155  /// @param mask Optional alpha mask.
156  void offset(ValueType offset, const MaskType* mask = nullptr);
157 
158  /// @brief Used internally by tbb::parallel_for()
159  /// @param range Range of LeafNodes over which to multi-thread.
160  ///
161  /// @warning Never call this method directly!
162  void operator()(const RangeType& range) const
163  {
164  if (mTask) mTask(const_cast<Filter*>(this), range);
165  else OPENVDB_THROW(ValueError, "task is undefined - call median(), mean(), etc.");
166  }
167 
168 private:
169  using LeafT = typename TreeType::LeafNodeType;
170  using VoxelIterT = typename LeafT::ValueOnIter;
171  using VoxelCIterT = typename LeafT::ValueOnCIter;
172  using BufferT = typename tree::LeafManager<TreeType>::BufferType;
173  using LeafIterT = typename RangeType::Iterator;
175 
176  void cook(LeafManagerType& leafs);
177 
178  template<size_t Axis>
179  struct Avg {
180  Avg(const GridT* grid, Int32 w): acc(grid->tree()), width(w), frac(1.f/float(2*w+1)) {}
181  inline ValueType operator()(Coord xyz);
182  typename GridT::ConstAccessor acc;
183  const Int32 width;
184  const float frac;
185  };
186 
187  // Private filter methods called by tbb::parallel_for threads
188  template <typename AvgT>
189  void doBox(const RangeType& r, Int32 w);
190  void doBoxX(const RangeType& r, Int32 w) { this->doBox<Avg<0> >(r,w); }
191  void doBoxY(const RangeType& r, Int32 w) { this->doBox<Avg<1> >(r,w); }
192  void doBoxZ(const RangeType& r, Int32 w) { this->doBox<Avg<2> >(r,w); }
193  void doMedian(const RangeType&, int);
194  void doOffset(const RangeType&, ValueType);
195  /// @return true if the process was interrupted
196  bool wasInterrupted();
197 
198  GridType* mGrid;
199  typename std::function<void (Filter*, const RangeType&)> mTask;
200  InterruptT* mInterrupter;
201  const MaskType* mMask;
202  int mGrainSize;
203  AlphaType mMinMask, mMaxMask;
204  bool mInvertMask;
205  bool mTiles;
206 }; // end of Filter class
207 
208 
209 ////////////////////////////////////////
210 
211 /// @cond OPENVDB_DOCS_INTERNAL
212 
213 namespace filter_internal {
214 
215 template<typename TreeT>
216 struct Voxelizer
217 {
218  // NodeManager for processing internal/root node values
219  // @note Should not cache leaf nodes
220  using NodeManagerT = tree::NodeManager<TreeT, TreeT::RootNodeType::LEVEL-1>;
221  using MaskT = typename TreeT::template ValueConverter<ValueMask>::Type;
222 
223  Voxelizer(TreeT& tree, const bool allNeighbors, const size_t grainSize)
224  : mVoxelTopology()
225  , mManager(nullptr)
226  , mGrainSize(grainSize)
227  , mOp(tree, mVoxelTopology, allNeighbors ? 26 : 6) {}
228 
229  /// @brief Convert tiles to leaf nodes that exist at a particular
230  /// voxel distance away
231  /// @param width distance in voxels to seach for tiles from each leaf
232  /// @return Returns how many search iterations were performed, which
233  /// also represents how many leaf node neighbors may have been
234  /// created. Returns 0 if the tree is already entirely voxelized
235  int run(const int width)
236  {
237  if (!mOp.tree().hasActiveTiles()) return 0;
238  this->init();
239  int count = 0;
240  for (int i = 0; i < width; i += int(TreeT::LeafNodeType::DIM), ++count) {
241  if (i > 0) mManager->rebuild();
242  mManager->foreachBottomUp(mOp, mGrainSize > 0, mGrainSize);
243  mOp.tree().topologyUnion(mVoxelTopology);
244  }
245  return count;
246  }
247 
248 private:
249  void init()
250  {
251  if (mManager) {
252  mManager->rebuild();
253  }
254  else {
255  // @note We don't actually need the leaf topology here, just the
256  // internal node structure so that we can generate leaf nodes in parallel
257  mVoxelTopology.topologyUnion(mOp.tree());
258  mManager.reset(new NodeManagerT(mOp.tree()));
259  }
260  }
261 
262  struct CreateVoxelMask
263  {
264  using LeafT = typename TreeT::LeafNodeType;
265  using RootT = typename TreeT::RootNodeType;
266 
267  CreateVoxelMask(TreeT& tree, MaskT& mask, const size_t NN)
268  : mTree(tree), mVoxelTopology(mask), mNeighbors(NN) {}
269 
270  TreeT& tree() { return mTree; }
271 
272  // do nothing for leaf nodes. They shouldn't even be cached as
273  // part of the NodeManager used with this method.
274  void operator()(const LeafT&) const { assert(false); }
275 
276  void operator()(const RootT& node) const
277  {
278  using ChildT = typename RootT::ChildNodeType;
279  static constexpr Int32 CHILDDIM = Int32(ChildT::DIM);
280  static constexpr Int32 LEAFDIM = Int32(LeafT::DIM);
281  const Tester op(mTree, mNeighbors);
282 
283  auto step =
284  [&](const Coord& ijk,
285  const size_t axis1,
286  const size_t axis2,
287  const auto& val)
288  {
289  Coord offset(0);
290  Int32& a = offset[axis1];
291  Int32& b = offset[axis2];
292  for (a = 0; a < CHILDDIM; a+=LEAFDIM) {
293  for (b = 0; b < CHILDDIM; b+=LEAFDIM) {
294  const Coord childijk = ijk + offset;
295  if (op.test(childijk, val)) {
296  mVoxelTopology.touchLeaf(childijk);
297  }
298  }
299  }
300 
301  offset.reset(CHILDDIM-1);
302  for (a = 0; a < CHILDDIM; a+=LEAFDIM) {
303  for (b = 0; b < CHILDDIM; b+=LEAFDIM) {
304  const Coord childijk = ijk + offset;
305  if (op.test(childijk, val)) {
306  mVoxelTopology.touchLeaf(childijk);
307  }
308  }
309  }
310  };
311 
312  for (auto iter = node.cbeginValueOn(); iter; ++iter) {
313  const Coord& ijk = iter.getCoord();
314  // @todo step only needs to search if a given direction
315  // depending on the face
316  step(ijk, 0, 1, *iter);
317  step(ijk, 0, 2, *iter);
318  step(ijk, 1, 2, *iter);
319  }
320  }
321 
322  template<typename NodeT>
323  void operator()(const NodeT& node) const
324  {
325  using ChildT = typename NodeT::ChildNodeType;
326  static constexpr Int32 CHILDDIM = Int32(ChildT::DIM);
327  static constexpr Int32 LEAFDIM = Int32(LeafT::DIM);
328 
329  static auto step =
330  [](const Tester& op,
331  const Coord& ijk,
332  const size_t axis1,
333  const size_t axis2,
334  const auto& val,
335  std::vector<Coord>& coords)
336  {
337  Coord offset(0);
338  Int32& a = offset[axis1];
339  Int32& b = offset[axis2];
340  for (a = 0; a < CHILDDIM; a+=LEAFDIM) {
341  for (b = 0; b < CHILDDIM; b+=LEAFDIM) {
342  const Coord childijk = ijk + offset;
343  if (op.test(childijk, val)) {
344  coords.emplace_back(childijk);
345  }
346  }
347  }
348 
349  offset.reset(CHILDDIM-1);
350  for (a = 0; a < CHILDDIM; a+=LEAFDIM) {
351  for (b = 0; b < CHILDDIM; b+=LEAFDIM) {
352  const Coord childijk = ijk + offset;
353  if (op.test(childijk, val)) {
354  coords.emplace_back(childijk);
355  }
356  }
357  }
358  };
359 
360  /// Two types of algorithms here
361  /// 1) For the case where this node is the direct parent of leaf nodes
362  /// 2) For all other node types
363  ///
364  /// In general, given a tile's ijk, search its faces/edges/corners for
365  /// values which differ from its own or leaf level topology. When a
366  /// difference is detected, mask topology is generated which can be used
367  /// with topologyUnion to ensure valid voxel values exist in the source
368  /// grid.
369  ///
370  /// This operator handles all internal node types. For example, for the
371  /// lowest level internal node (which contains leaf nodes as children)
372  /// each tile is at the leaf level (a single tile represents an 8x8x8
373  /// node). CHILDDIM is this case will match the valid of LEAFDIM, as we
374  /// only need to check each tiles immediate neighbors. For higher level
375  /// internal nodes (and the root node) each child tile will have a
376  /// significantly larger CHILDDIM than the grid's LEAFDIM. We
377  /// consistently probe values along the LEAFDIM stride to ensure no
378  /// changes are missed.
379 
380  if (CHILDDIM == LEAFDIM) {
381  // If the current node is the parent of leaf nodes, search each
382  // neighbor directly and use a flag buffer to test offsets in
383  // this node which need converting to leaf level topology.
384  // This is faster than the more general method which steps across
385  // faces (unecessary due to CHILDDIM == LEAFDIM) and provides
386  // a simpler way of tracking new topology
387 
388  std::vector<char> flags(NodeT::NUM_VALUES, char(0));
389  tbb::parallel_for(tbb::blocked_range<size_t>(0, NodeT::NUM_VALUES),
390  [&](const tbb::blocked_range<size_t>& range) {
391  const Tester op(mTree, mNeighbors);
392  for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
393  if (node.isValueMaskOn(Index(n))) {
394  // if index is a tile, search its neighbors
395  const Coord ijk = node.offsetToGlobalCoord(Index(n));
396  flags[n] = op.test(ijk, node.getValue(ijk));
397  }
398  }
399  });
400 
401  // create leaf level topology in this internal node
402  Index idx = 0;
403  for (auto iter = flags.begin(); iter != flags.end(); ++iter, ++idx) {
404  if (*iter) mVoxelTopology.touchLeaf(node.offsetToGlobalCoord(idx));
405  }
406  }
407  else {
408  // If this is a higher level internal node, we only need to search its
409  // face/edge/vertex neighbors for values which differ or leaf level
410  // topology. When a difference is detected, store the coordinate which
411  // needs to be voxelized.
412  // @todo investigate better threaded impl
413 
414  tbb::concurrent_vector<Coord> nodes;
415  tbb::parallel_for(tbb::blocked_range<size_t>(0, NodeT::NUM_VALUES),
416  [&](const tbb::blocked_range<size_t>& range)
417  {
418  const Tester op(mTree, mNeighbors);
419  std::vector<Coord> coords;
420 
421  for (size_t n = range.begin(), N = range.end(); n < N; ++n) {
422  if (!node.isValueMaskOn(Index(n))) continue;
423 
424  const Coord ijk = node.offsetToGlobalCoord(Index(n));
425  const auto& val = node.getValue(ijk);
426  // @todo step only needs to search if a given direction
427  // depending on the face
428  step(op, ijk, 0, 1, val, coords);
429  step(op, ijk, 0, 2, val, coords);
430  step(op, ijk, 1, 2, val, coords);
431  }
432 
433  if (!coords.empty()) {
434  std::copy(coords.begin(), coords.end(),
435  nodes.grow_by(coords.size()));
436  }
437  });
438 
439  // create leaf level topology in this internal node
440  // @note nodes may contain duplicate coords
441  for (const auto& coord : nodes) {
442  mVoxelTopology.touchLeaf(coord);
443  }
444  }
445  }
446 
447  private:
448  struct Tester
449  {
450  Tester(const TreeT& tree, const size_t NN)
451  : mAcc(tree), mNeighbors(NN) {}
452 
453  inline bool test(const Coord& ijk,
454  const typename TreeT::ValueType& val) const
455  {
456  static constexpr Int32 LEAFDIM = Int32(LeafT::DIM);
457  const Coord* NN = util::COORD_OFFSETS;
458  for (size_t i = 0; i < mNeighbors; ++i, ++NN) {
459  Coord neighbor(*NN);
460  neighbor.x() *= LEAFDIM;
461  neighbor.y() *= LEAFDIM;
462  neighbor.z() *= LEAFDIM;
463  neighbor += ijk;
464  // if a leaf exists, assume its buffer is not constant
465  if (mAcc.getValue(neighbor) != val ||
466  mAcc.probeConstLeaf(neighbor)) {
467  return true;
468  }
469  }
470  return false;
471  }
472  private:
474  const size_t mNeighbors;
475  };
476 
477  private:
478  TreeT& mTree;
479  MaskT& mVoxelTopology;
480  const size_t mNeighbors;
481  };// CreateVoxelMask
482 
483 private:
484  MaskT mVoxelTopology;
485  std::unique_ptr<NodeManagerT> mManager;
486  const size_t mGrainSize;
487  CreateVoxelMask mOp;
488 };
489 
490 // Helper function for Filter::Avg::operator()
491 template<typename T> static inline void accum(T& sum, T addend) { sum += addend; }
492 // Overload for bool ValueType
493 inline void accum(bool& sum, bool addend) { sum = sum || addend; }
494 
495 } // namespace filter_internal
496 
497 /// @endcond
498 
499 ////////////////////////////////////////
500 
501 
502 template<typename GridT, typename MaskT, typename InterruptT>
503 template<size_t Axis>
504 inline typename GridT::ValueType
506 {
507  ValueType sum = zeroVal<ValueType>();
508  Int32 &i = xyz[Axis], j = i + width;
509  for (i -= width; i <= j; ++i) filter_internal::accum(sum, acc.getValue(xyz));
510  if constexpr(std::is_same<ValueType, bool>::value) {
511  return sum && frac > 0.0f;
512  } else {
514  ValueType value = static_cast<ValueType>(sum * frac);
516  return value;
517  }
518 }
519 
520 
521 ////////////////////////////////////////
522 
523 
524 template<typename GridT, typename MaskT, typename InterruptT>
525 void
526 Filter<GridT, MaskT, InterruptT>::mean(int width, int iterations, const MaskType* mask)
527 {
528  if (iterations <= 0) return;
529  mMask = mask;
530  const int w = std::max(1, width);
531  const bool serial = mGrainSize == 0;
532 
533  if (mInterrupter) mInterrupter->start("Applying mean filter");
534 
535  std::unique_ptr<filter_internal::Voxelizer<TreeType>> voxelizer;
536  if (this->getProcessTiles()) {
537  // if performing multiple iterations, also search edge/vertex
538  // neighbors for difference topology.
539  const bool allNeighbors = iterations > 1;
540  // If processing tiles, create a voxelizer and run a single
541  // width based search for tiles that need to be voxelized
542  voxelizer.reset(new filter_internal::Voxelizer<TreeType>
543  (mGrid->tree(), allNeighbors, mGrainSize));
544  if (!voxelizer->run(w)) voxelizer.reset();
545  }
546 
547  LeafManagerType leafs(mGrid->tree(), 1, serial);
548 
549  int iter = 1; // num of leaf level neighbor based searches performed
550  int dist = w; // kernel distance of the current iteration
551  for (int i=0; i<iterations && !this->wasInterrupted(); ++i, dist+=w) {
552  if (i > 0 && voxelizer) {
553  // the total influence distance in voxels of this iteration
554  // minus how far we've already accounted for
555  const int remain = dist - iter * int(TreeType::LeafNodeType::DIM);
556  if (remain > 0) {
557  const int searches = voxelizer->run(remain);
558  if (searches == 0) voxelizer.reset();
559  else leafs.rebuild(serial);
560  iter += searches;
561  }
562  }
563 
564  mTask = std::bind(&Filter::doBoxX, std::placeholders::_1, std::placeholders::_2, w);
565  this->cook(leafs);
566  // note that the order of the YZ passes are flipped to maintain backwards-compatibility
567  // with an indexing typo in the original logic
568  mTask = std::bind(&Filter::doBoxZ, std::placeholders::_1, std::placeholders::_2, w);
569  this->cook(leafs);
570  mTask = std::bind(&Filter::doBoxY, std::placeholders::_1, std::placeholders::_2, w);
571  this->cook(leafs);
572  }
573 
574  if (mInterrupter) mInterrupter->end();
575 }
576 
577 
578 template<typename GridT, typename MaskT, typename InterruptT>
579 void
580 Filter<GridT, MaskT, InterruptT>::gaussian(int width, int iterations, const MaskType* mask)
581 {
582  if (iterations <= 0) return;
583  mMask = mask;
584  const int w = std::max(1, width);
585  const bool serial = mGrainSize == 0;
586 
587  if (mInterrupter) mInterrupter->start("Applying Gaussian filter");
588 
589  std::unique_ptr<filter_internal::Voxelizer<TreeType>> voxelizer;
590  if (this->getProcessTiles()) {
591  // if performing multiple iterations, also search edge/vertex
592  // neighbors for difference topology.
593  const bool allNeighbors = iterations > 1;
594  // If processing tiles, create a voxelizer and run a single
595  // width based search for tiles that need to be voxelized
596  // @note account for sub iteration due to gaussian filter
597  voxelizer.reset(new filter_internal::Voxelizer<TreeType>
598  (mGrid->tree(), allNeighbors, mGrainSize));
599  if (!voxelizer->run(w*4)) voxelizer.reset();
600  }
601 
602  LeafManagerType leafs(mGrid->tree(), 1, serial);
603 
604  int iter = 1; // num of leaf level neighbor based searches performed
605  int dist = w*4; // kernel distance of the current iteration
606  for (int i=0; i<iterations; ++i, dist+=(w*4)) {
607  if (i > 0 && voxelizer) {
608  // the total influence distance in voxels of this iteration
609  // minus how far we've already accounted for
610  const int remain = dist - iter * int(TreeType::LeafNodeType::DIM);
611  if (remain > 0) {
612  const int searches = voxelizer->run(remain);
613  if (searches == 0) voxelizer.reset();
614  else leafs.rebuild(serial);
615  iter += searches;
616  }
617  }
618 
619  for (int n=0; n<4 && !this->wasInterrupted(); ++n) {
620  mTask = std::bind(&Filter::doBoxX, std::placeholders::_1, std::placeholders::_2, w);
621  this->cook(leafs);
622  // note that the order of the YZ passes are flipped to maintain backwards-compatibility
623  // with an indexing typo in the original logic
624  mTask = std::bind(&Filter::doBoxZ, std::placeholders::_1, std::placeholders::_2, w);
625  this->cook(leafs);
626  mTask = std::bind(&Filter::doBoxY, std::placeholders::_1, std::placeholders::_2, w);
627  this->cook(leafs);
628  }
629  }
630 
631  if (mInterrupter) mInterrupter->end();
632 }
633 
634 
635 template<typename GridT, typename MaskT, typename InterruptT>
636 void
637 Filter<GridT, MaskT, InterruptT>::median(int width, int iterations, const MaskType* mask)
638 {
639  if (iterations <= 0) return;
640  mMask = mask;
641  const int w = std::max(1, width);
642  const bool serial = mGrainSize == 0;
643 
644  if (mInterrupter) mInterrupter->start("Applying median filter");
645 
646  std::unique_ptr<filter_internal::Voxelizer<TreeType>> voxelizer;
647  if (this->getProcessTiles()) {
648  // If processing tiles, create a voxelizer and run a single
649  // width based search for tiles that need to be voxelized
650  voxelizer.reset(new filter_internal::Voxelizer<TreeType>
651  (mGrid->tree(), /*allNeighbors*/true, mGrainSize));
652  if (!voxelizer->run(w)) voxelizer.reset();
653  }
654 
655  LeafManagerType leafs(mGrid->tree(), 1, serial);
656 
657  mTask = std::bind(&Filter::doMedian, std::placeholders::_1, std::placeholders::_2, w);
658 
659  int iter = 1; // num of leaf level neighbor based searches performed
660  int dist = w; // kernel distance of the current iteration
661  for (int i=0; i<iterations && !this->wasInterrupted(); ++i, dist+=w) {
662  if (i > 0 && voxelizer) {
663  // the total influence distance in voxels of this iteration
664  // minus how far we've already accounted for
665  const int remain = dist - iter * int(TreeType::LeafNodeType::DIM);
666  if (remain > 0) {
667  const int searches = voxelizer->run(remain);
668  if (searches == 0) voxelizer.reset();
669  else leafs.rebuild(serial);
670  iter += searches;
671  }
672  }
673 
674  this->cook(leafs);
675  }
676 
677  if (mInterrupter) mInterrupter->end();
678 }
679 
680 
681 template<typename GridT, typename MaskT, typename InterruptT>
682 void
684 {
685  mMask = mask;
686 
687  if (mInterrupter) mInterrupter->start("Applying offset");
688 
689  if (this->getProcessTiles()) {
690  // Don't process leaf nodes with the node manager - we'll do them
691  // separately to allow for cleaner branching
692  using NodeManagerT = tree::NodeManager<TreeType, TreeType::RootNodeType::LEVEL-1>;
693  NodeManagerT manager(mGrid->tree());
694 
695  if (mask) {
696  manager.foreachBottomUp([&](auto& node) {
697  this->wasInterrupted();
698  AlphaMaskT alpha(*mGrid, *mMask, mMinMask, mMaxMask, mInvertMask);
699  typename AlphaMaskT::FloatType a, b;
700  for (auto iter = node.beginValueOn(); iter; ++iter) {
701  if (!alpha(iter.getCoord(), a, b)) continue;
703  iter.modifyValue([&](ValueType& v) { v += a*value; });
705  }
706  });
707  }
708  else {
709  manager.foreachBottomUp([&](auto& node) {
710  this->wasInterrupted();
711  for (auto iter = node.beginValueOn(); iter; ++iter) {
712  iter.modifyValue([&](ValueType& v) { v += value; });
713  }
714  });
715  }
716  }
717 
718  LeafManagerType leafs(mGrid->tree(), 0, mGrainSize==0);
719  mTask = std::bind(&Filter::doOffset, std::placeholders::_1, std::placeholders::_2, value);
720  this->cook(leafs);
721 
722  if (mInterrupter) mInterrupter->end();
723 }
724 
725 
726 ////////////////////////////////////////
727 
728 
729 /// Private method to perform the task (serial or threaded) and
730 /// subsequently swap the leaf buffers.
731 template<typename GridT, typename MaskT, typename InterruptT>
732 void
734 {
735  if (mGrainSize>0) {
736  tbb::parallel_for(leafs.leafRange(mGrainSize), *this);
737  } else {
738  (*this)(leafs.leafRange());
739  }
740  leafs.swapLeafBuffer(1, mGrainSize==0);
741 }
742 
743 
744 /// One dimensional convolution of a separable box filter
745 template<typename GridT, typename MaskT, typename InterruptT>
746 template <typename AvgT>
747 void
749 {
750  this->wasInterrupted();
751  AvgT avg(mGrid, w);
752  if (mMask) {
753  typename AlphaMaskT::FloatType a, b;
754  AlphaMaskT alpha(*mGrid, *mMask, mMinMask, mMaxMask, mInvertMask);
755  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
756  BufferT& buffer = leafIter.buffer(1);
757  for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
758  const Coord xyz = iter.getCoord();
759  if (alpha(xyz, a, b)) {
761  const ValueType value(b*(*iter) + a*avg(xyz));
763  buffer.setValue(iter.pos(), value);
764  }
765  }
766  }
767  } else {
768  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
769  BufferT& buffer = leafIter.buffer(1);
770  for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
771  buffer.setValue(iter.pos(), avg(iter.getCoord()));
772  }
773  }
774  }
775 }
776 
777 
778 /// Performs simple but slow median-value diffusion
779 template<typename GridT, typename MaskT, typename InterruptT>
780 void
782 {
783  this->wasInterrupted();
784  typename math::DenseStencil<GridType> stencil(*mGrid, width);//creates local cache!
785  if (mMask) {
786  typename AlphaMaskT::FloatType a, b;
787  AlphaMaskT alpha(*mGrid, *mMask, mMinMask, mMaxMask, mInvertMask);
788  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
789  BufferT& buffer = leafIter.buffer(1);
790  for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
791  if (alpha(iter.getCoord(), a, b)) {
792  stencil.moveTo(iter);
794  ValueType value(b*(*iter) + a*stencil.median());
796  buffer.setValue(iter.pos(), value);
797  }
798  }
799  }
800  } else {
801  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
802  BufferT& buffer = leafIter.buffer(1);
803  for (VoxelCIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
804  stencil.moveTo(iter);
805  buffer.setValue(iter.pos(), stencil.median());
806  }
807  }
808  }
809 }
810 
811 
812 /// Offsets the values by a constant
813 template<typename GridT, typename MaskT, typename InterruptT>
814 void
816 {
817  this->wasInterrupted();
818  if (mMask) {
819  typename AlphaMaskT::FloatType a, b;
820  AlphaMaskT alpha(*mGrid, *mMask, mMinMask, mMaxMask, mInvertMask);
821  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
822  for (VoxelIterT iter = leafIter->beginValueOn(); iter; ++iter) {
823  if (alpha(iter.getCoord(), a, b)) {
825  ValueType value(*iter + a*offset);
827  iter.setValue(value);
828  }
829  }
830  }
831  } else {
832  for (LeafIterT leafIter=range.begin(); leafIter; ++leafIter) {
833  for (VoxelIterT iter = leafIter->beginValueOn(); iter; ++iter) {
834  iter.setValue(*iter + offset);
835  }
836  }
837  }
838 }
839 
840 
841 template<typename GridT, typename MaskT, typename InterruptT>
842 inline bool
844 {
845  if (util::wasInterrupted(mInterrupter)) {
846  thread::cancelGroupExecution();
847  return true;
848  }
849  return false;
850 }
851 
852 
853 ////////////////////////////////////////
854 
855 
856 // Explicit Template Instantiation
857 
858 #ifdef OPENVDB_USE_EXPLICIT_INSTANTIATION
859 
860 #ifdef OPENVDB_INSTANTIATE_FILTER
862 #endif
863 
866 
867 #endif // OPENVDB_USE_EXPLICIT_INSTANTIATION
868 
869 
870 } // namespace tools
871 } // namespace OPENVDB_VERSION_NAME
872 } // namespace openvdb
873 
874 #endif // OPENVDB_TOOLS_FILTER_HAS_BEEN_INCLUDED
#define OPENVDB_NO_TYPE_CONVERSION_WARNING_BEGIN
Bracket code with OPENVDB_NO_TYPE_CONVERSION_WARNING_BEGIN/_END, to inhibit warnings about type conve...
Definition: Platform.h:227
To facilitate threading over the nodes of a tree, cache node pointers in linear arrays, one for each level of the tree.
Definition: NodeManager.h:30
Definition: Exceptions.h:65
void invertMask(bool invert=true)
Invert the optional mask, i.e. min->max in the original mask maps to 1->0 in the inverted alpha mask...
Definition: Filter.h:128
typename CopyConstness< TreeType, NonConstBufferType >::Type BufferType
Definition: LeafManager.h:95
typename LeafManagerType::BufferType BufferType
Definition: Filter.h:56
bool wasInterrupted(T *i, int percent=-1)
Definition: NullInterrupter.h:49
The Value Accessor Implementation and API methods. The majoirty of the API matches the API of a compa...
Definition: ValueAccessor.h:68
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:74
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
AlphaType minMask() const
Return the minimum value of the mask to be used for the derivation of a smooth alpha value...
Definition: Filter.h:106
void moveTo(const Coord &ijk)
Initialize the stencil buffer with the values of voxel (x, y, z) and its neighbors.
Definition: Stencils.h:1784
void setProcessTiles(bool flag)
Set whether active tiles should also be processed.
Definition: Filter.h:102
Filter(GridT &grid, InterruptT *interrupt=nullptr)
Definition: Filter.h:63
typename MaskType::ValueType AlphaType
Definition: Filter.h:53
Volume filtering (e.g., diffusion) with optional alpha masking.
Definition: Filter.h:45
void foreachBottomUp(const NodeOp &op, bool threaded=true, size_t grainSize=1)
Threaded method that applies a user-supplied functor to all the nodes in the tree.
Definition: NodeManager.h:624
typename TreeType::LeafNodeType LeafType
Definition: Filter.h:51
bool getProcessTiles() const
Definition: Filter.h:95
GridT GridType
Definition: Filter.h:48
#define OPENVDB_NO_TYPE_CONVERSION_WARNING_END
Definition: Platform.h:228
FloatT FloatType
Definition: Interpolation.h:551
Definition: Interpolation.h:543
constexpr Coord COORD_OFFSETS[26]
coordinate offset table for neighboring voxels
Definition: Util.h:22
bool isMaskInverted() const
Return true if the mask is inverted, i.e. min->max in the original mask maps to 1->0 in the inverted ...
Definition: Filter.h:125
typename GridType::TreeType TreeType
Definition: Filter.h:50
void setMaskRange(AlphaType min, AlphaType max)
Define the range for the (optional) scalar mask.
Definition: Filter.h:116
MaskT MaskType
Definition: Filter.h:49
typename tree::LeafManager< TreeType > LeafManagerType
Definition: Filter.h:54
int getGrainSize() const
Definition: Filter.h:89
Axis
Definition: Math.h:901
Definition: Exceptions.h:13
ValueType median() const
Return the median value of the current stencil.
Definition: Stencils.h:121
typename LeafManagerType::LeafRange RangeType
Definition: Filter.h:55
Filter(const Filter &other)
Shallow copy constructor called by tbb::parallel_for() threads during filtering.
Definition: Filter.h:77
Index32 Index
Definition: Types.h:54
Dense stencil of a given width.
Definition: Stencils.h:1764
#define OPENVDB_INSTANTIATE_CLASS
Definition: version.h.in:153
AlphaType maxMask() const
Return the maximum value of the mask to be used for the derivation of a smooth alpha value...
Definition: Filter.h:109
typename GridType::ValueType ValueType
Definition: Filter.h:52
OPENVDB_AX_API void run(const char *ax, openvdb::GridBase &grid, const AttributeBindings &bindings={})
Run a full AX pipeline (parse, compile and execute) on a single OpenVDB Grid.
void operator()(const RangeType &range) const
Used internally by tbb::parallel_for()
Definition: Filter.h:162
const std::enable_if<!VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:110
A LeafManager manages a linear array of pointers to a given tree&#39;s leaf nodes, as well as optional au...
NodeManager produces linear arrays of all tree nodes allowing for efficient threading and bottom-up p...
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h.in:121
void setGrainSize(int grainsize)
Set the grain-size used for multi-threading.
Definition: Filter.h:92
int32_t Int32
Definition: Types.h:56
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