OpenVDB  0.104.0
LevelSetUtil.h
Go to the documentation of this file.
1 
2 //
3 // Copyright (c) 2012 DreamWorks Animation LLC
4 //
5 // All rights reserved. This software is distributed under the
6 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
7 //
8 // Redistributions of source code must retain the above copyright
9 // and license notice and the following restrictions and disclaimer.
10 //
11 // * Neither the name of DreamWorks Animation nor the names of
12 // its contributors may be used to endorse or promote products derived
13 // from this software without specific prior written permission.
14 //
15 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
16 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
17 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
18 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
19 // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY INDIRECT, INCIDENTAL,
20 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
21 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
22 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
23 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 // IN NO EVENT SHALL THE COPYRIGHT HOLDERS' AND CONTRIBUTORS' AGGREGATE
27 // LIABILITY FOR ALL CLAIMS REGARDLESS OF THEIR BASIS EXCEED US$250.00.
28 //
30 
31 #ifndef OPENVDB_TOOLS_LEVELSETUTIL_HAS_BEEN_INCLUDED
32 #define OPENVDB_TOOLS_LEVELSETUTIL_HAS_BEEN_INCLUDED
33 
34 #include <openvdb/Grid.h>
35 #include <openvdb/tools/Composite.h>
36 #include <openvdb/tools/GridTransformer.h>
37 #include <openvdb/tree/LeafManager.h>
38 #include <openvdb/util/NullInterrupter.h>
39 #include <openvdb/util/Util.h>
40 
41 #include <tbb/blocked_range.h>
42 #include <tbb/parallel_for.h>
43 #include <tbb/parallel_reduce.h>
44 
45 #include <boost/random.hpp>
46 #include <boost/generator_iterator.hpp>
47 #include <boost/math/constants/constants.hpp>
48 
49 #include <limits>
50 #include <list>
51 #include <deque> // used by segment()
52 
53 
54 namespace openvdb {
56 namespace OPENVDB_VERSION_NAME {
57 namespace tools {
58 
59 // MS Visual C++ requires this extra level of indirection in order to compile
60 // THIS MUST EXIST IN AN UNNAMED NAMESPACE IN ORDER TO COMPILE ON WINDOWS
61 namespace {
62 
63 template<typename GridType>
64 inline typename GridType::ValueType lsutilGridMax()
65 {
67 }
68 
69 template<typename GridType>
70 inline typename GridType::ValueType lsutilGridZero()
71 {
72  return zeroVal<typename GridType::ValueType>();
73 }
74 
75 } // unnamed namespace
76 
77 
89 template<class GridType>
90 inline void
92  GridType& grid,
93  typename GridType::ValueType cutOffDistance = lsutilGridMax<GridType>());
94 
95 
97 
98 
108 template<class GridType>
109 inline typename Grid<typename GridType::TreeType::template ValueConverter<bool>::Type>::Ptr
111  const GridType& grid,
112  typename GridType::ValueType iso = lsutilGridZero<GridType>());
113 
114 
116 
117 
119 template<class TreeType>
121 {
122 public:
124  typedef typename TreeType::ValueType ValueType;
125 
128 
129  void runParallel();
130  void runSerial();
131 
132  const ValueType& minVoxel() const { return mMin; }
133  const ValueType& maxVoxel() const { return mMax; }
134 
135 
136  inline MinMaxVoxel(const MinMaxVoxel<TreeType>&, tbb::split);
137  inline void operator()(const tbb::blocked_range<size_t>&);
138  inline void join(const MinMaxVoxel<TreeType>&);
139 
140 private:
141  LeafArray& mLeafArray;
142  ValueType mMin, mMax;
143 };
144 
145 
147 
148 
150 template<class TreeType, class LeafOp>
152 {
153 public:
155  typedef typename TreeType::ValueType ValueType;
156 
158  LeafTransformer(LeafArray&, LeafOp&);
159 
160  void runParallel();
161  void runSerial();
162 
163 
164  inline void operator()(const tbb::blocked_range<size_t>&) const;
166 
167 private:
168  LeafArray& mLeafArray;
169  LeafOp& mLeafOp;
170 };
171 
172 
174 
175 
177 template<class GridType, class InterruptType = util::NullInterrupter>
179 {
180 public:
181  typedef std::vector<Vec3s> PointVec;
182  typedef std::list<typename GridType::Ptr> GridPtrList;
183  typedef typename GridPtrList::iterator GridPtrListIter;
184 
185 
189  LevelSetFracture(InterruptType* interrupter = NULL);
190 
191 
207  void fracture(GridPtrList& grids, const GridType& cutter, bool segment = false,
208  const PointVec* points = NULL, bool randomizeRotation = true, bool cutterOverlap = true);
209 
210 
212  GridPtrList& fragments() { return mFragments; }
213 
214 
216  void clear() { mFragments.clear(); }
217 
218 
219 private:
220  // disallow copy by assignment
221  void operator=(const LevelSetFracture<GridType, InterruptType>&) {}
222 
223  bool wasInterrupted(int percent = -1) const {
224  return mInterrupter && mInterrupter->wasInterrupted(percent);
225  }
226 
227  bool isValidFragment(GridType&) const;
228  void segmentFragments(GridPtrList&) const;
229  void process(GridPtrList&, const GridType& cutter);
230 
231  InterruptType* mInterrupter;
232  GridPtrList mFragments;
233 };
234 
235 
238 
239 // Internal utility objects and implementation details
240 
241 
242 namespace internal {
243 
244 
245 template<typename ValueType>
247 {
248  FogVolumeOp(ValueType cutOffDistance)
249  : mWeight(ValueType(1.0) / cutOffDistance)
250  {
251  }
252 
253  // cutOff has to be < 0.0
254  template <typename LeafNodeType>
255  void operator()(LeafNodeType &leaf, size_t/*leafIndex*/) const
256  {
257  const ValueType zero = zeroVal<ValueType>();
258 
259  for (typename LeafNodeType::ValueAllIter iter = leaf.beginValueAll(); iter; ++iter) {
260  ValueType& value = const_cast<ValueType&>(iter.getValue());
261  if (value > zero) {
262  value = zero;
263  iter.setValueOff();
264  } else {
265  value = std::min(ValueType(1.0), value * mWeight);
266  iter.setValueOn(value > zero);
267  }
268  }
269  }
270 
271 private:
272  ValueType mWeight;
273 };
274 
275 
276 template<typename TreeType>
278 {
279  InteriorMaskOp(const TreeType& tree, typename TreeType::ValueType iso)
280  : mTree(tree)
281  , mIso(iso)
282  {
283  }
284 
285  template <typename LeafNodeType>
286  void operator()(LeafNodeType &leaf, size_t/*leafIndex*/) const
287  {
288  const Coord origin = leaf.getOrigin();
289  const typename TreeType::LeafNodeType* refLeafPt = mTree.probeConstLeaf(origin);
290 
291  if (refLeafPt != NULL) {
292 
293  const typename TreeType::LeafNodeType& refLeaf = *refLeafPt;
294  typename LeafNodeType::ValueAllIter iter = leaf.beginValueAll();
295 
296  for (; iter; ++iter) {
297  if (refLeaf.getValue(iter.pos()) < mIso) {
298  iter.setValueOn();
299  } else {
300  iter.setValueOff();
301  }
302  }
303  }
304  }
305 
306 private:
307  const TreeType& mTree;
308  typename TreeType::ValueType mIso;
309 };
310 
311 
314 template<typename GridType, typename InterruptType>
315 inline
316 std::vector<typename GridType::Ptr> segment(GridType& grid, InterruptType* interrupter = NULL)
317 {
318  typedef typename GridType::Ptr GridPtr;
319  typedef typename GridType::TreeType TreeType;
320  typedef typename TreeType::Ptr TreePtr;
321  typedef typename TreeType::ValueType ValueType;
322 
323  std::vector<GridPtr> segments;
324 
325  while (grid.activeVoxelCount() > 0) {
326 
327  if (interrupter && interrupter->wasInterrupted()) break;
328 
329  // Deep copy the grid's metadata (tree and transform are shared)
330  GridPtr segment(new GridType(grid, ShallowCopy()));
331  // Make the transform unique and insert an empty tree
332  segment->setTransform(grid.transform().copy());
333  TreePtr tree(new TreeType(grid.background()));
334  segment->setTree(tree);
335 
336  std::deque<Coord> coordList;
337  coordList.push_back(grid.tree().beginLeaf()->beginValueOn().getCoord());
338 
339  Coord ijk, n_ijk;
340  ValueType value;
341 
342  typename tree::ValueAccessor<TreeType> sourceAcc(grid.tree());
343  typename tree::ValueAccessor<TreeType> targetAcc(segment->tree());
344 
345  while (!coordList.empty()) {
346 
347  if (interrupter && interrupter->wasInterrupted()) break;
348 
349  ijk = coordList.back();
350  coordList.pop_back();
351 
352  if (!sourceAcc.probeValue(ijk, value)) continue;
353  if (targetAcc.isValueOn(ijk)) continue;
354 
355  targetAcc.setValue(ijk, value);
356  sourceAcc.setValueOff(ijk);
357 
358  for (int n = 0; n < 6; n++) {
359  n_ijk = ijk + util::COORD_OFFSETS[n];
360  if (!targetAcc.isValueOn(n_ijk) && sourceAcc.isValueOn(n_ijk)) {
361  coordList.push_back(n_ijk);
362  }
363  }
364  }
365 
366  grid.tree().pruneInactive();
367  segment->tree().signedFloodFill();
368  segments.push_back(segment);
369  }
370 
371 
372 
373  return segments;
374 }
375 
376 } // namespace internal
377 
378 
380 
381 
382 template <class TreeType>
384  : mLeafArray(leafs)
385  , mMin(std::numeric_limits<ValueType>::max())
386  , mMax(-mMin)
387 {
388 }
389 
390 template <class TreeType>
391 inline
393  const MinMaxVoxel<TreeType>& rhs, tbb::split)
394  : mLeafArray(rhs.mLeafArray)
395  , mMin(rhs.mMin)
396  , mMax(rhs.mMax)
397 {
398 }
399 
400 template <class TreeType>
401 void
403 {
404  tbb::parallel_reduce(mLeafArray.getRange(), *this);
405 }
406 
407 
408 template <class TreeType>
409 void
411 {
412  (*this)(mLeafArray.getRange());
413 }
414 
415 
416 template <class TreeType>
417 inline void
418 MinMaxVoxel<TreeType>::operator()(const tbb::blocked_range<size_t>& range)
419 {
420  typename TreeType::LeafNodeType::ValueOnCIter iter;
421 
422  for (size_t n = range.begin(); n < range.end(); ++n) {
423  iter = mLeafArray.leaf(n).cbeginValueOn();
424  for (; iter; ++iter) {
425  const ValueType value = iter.getValue();
426  mMin = std::min(mMin, value);
427  mMax = std::max(mMax, value);
428  }
429  }
430 }
431 
432 template <class TreeType>
433 inline void
435 {
436  mMin = std::min(mMin, rhs.mMin);
437  mMax = std::max(mMax, rhs.mMax);
438 }
439 
440 
442 
443 
444 template <class TreeType, class LeafOp>
446  LeafTransformer(LeafArray& leafs, LeafOp& leafOp)
447  : mLeafArray(leafs)
448  , mLeafOp(leafOp)
449 {
450 }
451 
452 template <class TreeType, class LeafOp>
453 inline
456  : mLeafArray(rhs.mLeafArray)
457  , mLeafOp(rhs.mLeafOp)
458 {
459 }
460 
461 template <class TreeType, class LeafOp>
462 void
464 {
465  tbb::parallel_for(mLeafArray.getRange(), *this);
466 }
467 
468 template <class TreeType, class LeafOp>
469 void
471 {
472  (*this)(mLeafArray.getRange());
473 }
474 
475 template <class TreeType, class LeafOp>
476 inline void
477 LeafTransformer<TreeType, LeafOp>::operator()(const tbb::blocked_range<size_t>& range) const
478 {
479  for (size_t n = range.begin(); n < range.end(); ++n) {
480  mLeafOp(mLeafArray.leaf(n), n);
481  }
482 }
483 
484 
486 
487 
488 template <class GridType>
489 inline void
490 sdfToFogVolume(GridType& grid, typename GridType::ValueType cutOffDistance)
491 {
492  typedef typename GridType::TreeType TreeType;
493  typedef typename GridType::ValueType ValueType;
494 
495  cutOffDistance = -std::abs(cutOffDistance);
496 
497  TreeType& tree = const_cast<TreeType&>(grid.tree());
498 
499  { // Transform all voxels (parallel, over leaf nodes)
500  tree::LeafManager<TreeType> leafs(tree);
501 
502  MinMaxVoxel<TreeType> minmax(leafs);
503  minmax.runParallel();
504 
505  // Clamp to the interior band width.
506  if (minmax.minVoxel() > cutOffDistance) {
507  cutOffDistance = minmax.minVoxel();
508  }
509 
510  internal::FogVolumeOp<ValueType> op(cutOffDistance);
512  transform.runParallel();
513  }
514 
515  // Transform all tile values (serial, but the iteration
516  // is constrained from descending into leaf nodes)
517  const ValueType zero = zeroVal<ValueType>();
518  typename TreeType::ValueAllIter iter(tree);
519  iter.setMaxDepth(TreeType::ValueAllIter::LEAF_DEPTH - 1);
520 
521  for ( ; iter; ++iter) {
522  ValueType& value = const_cast<ValueType&>(iter.getValue());
523 
524  if (value > zero) {
525  value = zero;
526  iter.setValueOff();
527  } else {
528  value = ValueType(1.0);
529  iter.setActiveState(true);
530  }
531  }
532 
533  // Update the tree background value.
534 
535  typename TreeType::Ptr newTree(new TreeType(/*background=*/zero));
536  newTree->merge(tree);
537  // This is faster than calling Tree::setBackground, since we only need
538  // to update the value that is returned for coordinates that don't fall
539  // inside an allocated node. All inactive tiles and voxels have already
540  // been updated in the previous step so the Tree::setBackground method
541  // will in this case do a redundant traversal of the tree to update the
542  // inactive values once more.
543 
544  //newTree->pruneInactive();
545  grid.setTree(newTree);
546 
547  grid.setGridClass(GRID_FOG_VOLUME);
548 }
549 
550 
552 
553 
554 template <class GridType>
556 sdfInteriorMask(const GridType& grid, typename GridType::ValueType iso)
557 {
558  typedef typename GridType::TreeType::template ValueConverter<bool>::Type BoolTreeType;
559  typedef Grid<BoolTreeType> BoolGridType;
560 
561  typename BoolGridType::Ptr maskGrid(BoolGridType::create(false));
562  maskGrid->setTransform(grid.transform().copy());
563  BoolTreeType& maskTree = maskGrid->tree();
564 
565  maskTree.topologyUnion(grid.tree());
566 
567  { // Evaluate voxels (parallel, over leaf nodes)
568 
569  tree::LeafManager<BoolTreeType> leafs(maskTree);
570 
572  InteriorMaskOp op(grid.tree(), iso);
573 
575  transform.runParallel();
576  }
577 
578  // Evaluate tile values (serial, but the iteration
579  // is constrained from descending into leaf nodes)
580 
582  typename BoolTreeType::ValueAllIter iter(maskTree);
583  iter.setMaxDepth(BoolTreeType::ValueAllIter::LEAF_DEPTH - 1);
584 
585  for ( ; iter; ++iter) {
586  iter.setActiveState(acc.getValue(iter.getCoord()) < iso);
587  }
588 
589  maskTree.pruneInactive();
590 
591  return maskGrid;
592 }
593 
594 
596 
597 
598 template<class GridType, class InterruptType>
600  : mInterrupter(interrupter)
601  , mFragments()
602 {
603 }
604 
605 
606 template<class GridType, class InterruptType>
607 void
609  bool segment, const PointVec* points, bool randomizeRotation, bool cutterOverlap)
610 {
611  if (points && points->size() != 0) {
612 
613 
614 
615  math::Transform::Ptr originalTransform = cutter.transform().copy();
616  const Vec3s cutterCenter = originalTransform->indexToWorld(
617  cutter.evalActiveVoxelBoundingBox().getCenter());
618 
619  boost::mt19937 rng(1);
620  boost::uniform_int<int> range(0, 359);
621  boost::variate_generator< boost::mt19937, boost::uniform_int<int> > randNr(rng, range);
622  const float deg2rad = boost::math::constants::pi<float>() / 180.0;
623 
624  GridType cutterGrid(cutter, ShallowCopy());
625 
626  // for each instance point..
627  for (size_t p = 0, P = points->size(); p < P; ++p) {
628 
629  int percent = int((float(p) / float(P)) * 100.0);
630  if (wasInterrupted(percent)) break;
631 
632  GridType instCutterGrid;
633  instCutterGrid.setTransform(originalTransform->copy());
634 
635  math::Transform::Ptr xform = originalTransform->copy();
636  xform->postTranslate(-cutterCenter);
637 
638  if (randomizeRotation) {
639 
640  const float xRot = float(randNr()) * deg2rad;
641  const float yRot = float(randNr()) * deg2rad;
642  const float zRot = float(randNr()) * deg2rad;
643 
644  xform->preRotate(xRot, openvdb::math::X_AXIS);
645  xform->preRotate(yRot, openvdb::math::Y_AXIS);
646  xform->preRotate(zRot, openvdb::math::Z_AXIS);
647 
648  math::Transform::Ptr T = originalTransform->copy();
649 
650  T->preRotate(-zRot, openvdb::math::Z_AXIS);
651  T->preRotate(-yRot, openvdb::math::Y_AXIS);
652  T->preRotate(-xRot, openvdb::math::X_AXIS);
653 
654  xform->postTranslate(T->indexToWorld((*points)[p]));
655 
656  } else {
657  xform->postTranslate(originalTransform->indexToWorld((*points)[p]));
658  }
659 
660  cutterGrid.setTransform(xform);
661 
662  if (wasInterrupted()) break;
663 
664  resampleToMatch<openvdb::tools::BoxSampler>(cutterGrid, instCutterGrid);
665 
666  if (cutterOverlap) process(mFragments, instCutterGrid);
667  process(grids, instCutterGrid);
668  }
669 
670  } else {
671  // use cutter in place
672  process(grids, cutter);
673  }
674 
675  if (segment) {
676  segmentFragments(mFragments);
677  segmentFragments(grids);
678  }
679 }
680 
681 
682 template<class GridType, class InterruptType>
683 bool
685 {
686  typedef typename GridType::TreeType TreeType;
687  if (grid.activeVoxelCount() < 27) return false;
688 
689  // Check if valid level-set
690  {
691  tree::LeafManager<TreeType> leafs(grid.tree());
692  MinMaxVoxel<TreeType> minmax(leafs);
693  minmax.runParallel();
694 
695  if ((minmax.minVoxel() < 0) == (minmax.maxVoxel() < 0)) return false;
696  }
697 
698  return true;
699 }
700 
701 
702 template<class GridType, class InterruptType>
703 void
704 LevelSetFracture<GridType, InterruptType>::segmentFragments(GridPtrList& grids) const
705 {
706  GridPtrList newFragments;
707  for (GridPtrListIter it = grids.begin(); it != grids.end(); ++it) {
708 
709  if (wasInterrupted()) break;
710 
711  std::vector<typename GridType::Ptr> segments = internal::segment(*(*it), mInterrupter);
712 
713  for (size_t n = 0, N = segments.size(); n < N; ++n) {
714 
715  if (wasInterrupted()) break;
716 
717  if (isValidFragment(*segments[n])) {
718  newFragments.push_back(segments[n]);
719  }
720  }
721  }
722  grids.swap(newFragments);
723 }
724 
725 
726 template<class GridType, class InterruptType>
727 void
728 LevelSetFracture<GridType, InterruptType>::process(
729  GridPtrList& grids, const GridType& cutter)
730 {
731  typedef typename GridType::Ptr GridPtr;
732 
733  GridPtrList newFragments;
734 
735  for (GridPtrListIter it = grids.begin(); it != grids.end(); ++it) {
736 
737  if (wasInterrupted()) break;
738 
739  GridPtr grid = *it;
740 
741  // gen new fragment
742  GridPtr fragment = grid->deepCopy();
743  openvdb::tools::csgIntersection(*fragment, *cutter.deepCopy());
744 
745  if (wasInterrupted()) break;
746 
747  if (!isValidFragment(*fragment)) continue;
748 
749  // update residual
750  GridPtr residual = grid->deepCopy();
751  openvdb::tools::csgDifference(*residual, *cutter.deepCopy());
752 
753  if (wasInterrupted()) break;
754 
755  if (!isValidFragment(*residual)) continue;
756 
757  newFragments.push_back(fragment);
758 
759  grid->tree().clear();
760  grid->tree().merge(residual->tree());
761  }
762 
763  if(!newFragments.empty()) {
764  mFragments.splice(mFragments.end(), newFragments);
765  }
766 }
767 
768 } // namespace tools
769 } // namespace OPENVDB_VERSION_NAME
770 } // namespace openvdb
771 
772 #endif // OPENVDB_TOOLS_LEVELSETUTIL_HAS_BEEN_INCLUDED
773 
774 // Copyright (c) 2012 DreamWorks Animation LLC
775 // All rights reserved. This software is distributed under the
776 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )