39 #ifndef OPENVDB_TOOLS_LEVEL_SET_TRACKER_HAS_BEEN_INCLUDED
40 #define OPENVDB_TOOLS_LEVEL_SET_TRACKER_HAS_BEEN_INCLUDED
42 #include <tbb/parallel_for.h>
43 #include <boost/bind.hpp>
44 #include <boost/function.hpp>
45 #include <boost/type_traits/is_floating_point.hpp>
66 template<
typename Gr
idT,
typename InterruptT = util::NullInterrupter>
72 typedef typename TreeType::LeafNodeType
LeafType;
77 typedef typename TreeType::template ValueConverter<bool>::Type
BoolMaskType;
78 BOOST_STATIC_ASSERT(boost::is_floating_point<ValueType>::value);
85 : spatialScheme(s), temporalScheme(t), normCount(n), grainSize(g) {}
100 template <
typename MaskType>
104 void normalize() { this->normalize<BoolMaskType>(NULL); }
128 void dilate(
int iterations = 1);
132 void erode(
int iterations = 1);
176 void startInterrupter(
const char* msg);
178 void endInterrupter();
181 bool checkInterrupter();
183 const GridType&
grid()
const {
return *mGrid; }
185 LeafManagerType&
leafs() {
return *mLeafs; }
187 const LeafManagerType&
leafs()
const {
return *mLeafs; }
201 void operator()(
const LeafRange& r)
const;
202 LevelSetTracker& mTracker;
212 typedef typename SchemeT::template ISStencil<GridType>::StencilType StencilT;
213 typedef typename MaskT::LeafNodeType MaskLeafT;
214 typedef typename MaskLeafT::ValueOnCIter MaskIterT;
215 typedef typename LeafType::ValueOnCIter VoxelIterT;
216 Normalizer(LevelSetTracker& tracker,
const MaskT* mask);
218 void operator()(
const LeafRange& r)
const {mTask(const_cast<Normalizer*>(
this), r);}
219 void cook(
int swapBuffer=0);
220 template <
int Nominator,
int Denominator>
221 void euler(
const LeafRange& range,
Index phiBuffer,
Index resultBuffer);
222 inline void euler01(
const LeafRange& r) {this->euler<0,1>(r, 0, 1);}
223 inline void euler12(
const LeafRange& r) {this->euler<1,2>(r, 1, 1);}
224 inline void euler34(
const LeafRange& r) {this->euler<3,4>(r, 1, 2);}
225 inline void euler13(
const LeafRange& r) {this->euler<1,3>(r, 1, 2);}
226 template <
int Nominator,
int Denominator>
227 void eval(StencilT& stencil,
const ValueType* phi, ValueType* result,
Index n)
const;
228 LevelSetTracker& mTracker;
230 const ValueType mDt, mInvDx;
231 typename boost::function<void (Normalizer*, const LeafRange&)> mTask;
234 template<math::BiasedGradientScheme SpatialScheme,
typename MaskT>
235 void normalize1(
const MaskT* mask);
239 void normalize2(
const MaskT* mask);
246 LeafManagerType* mLeafs;
247 InterruptT* mInterrupter;
252 template<
typename Gr
idT,
typename InterruptT>
253 LevelSetTracker<GridT, InterruptT>::
254 LevelSetTracker(GridT& grid, InterruptT* interrupt):
256 mLeafs(new LeafManagerType(grid.tree())),
257 mInterrupter(interrupt),
258 mDx(static_cast<ValueType>(grid.voxelSize()[0])),
261 if ( !grid.hasUniformVoxels() ) {
263 "The transform must have uniform scale for the LevelSetTracker to function");
267 "LevelSetTracker expected a level set, got a grid of class \""
268 + grid.gridClassToString(grid.getGridClass())
269 +
"\" [hint: Grid::setGridClass(openvdb::GRID_LEVEL_SET)]");
273 template<
typename Gr
idT,
typename InterruptT>
278 this->startInterrupter(
"Pruning Level Set");
288 mLeafs->rebuildLeafArray();
289 this->endInterrupter();
292 template<
typename Gr
idT,
typename InterruptT>
307 template<
typename Gr
idT,
typename InterruptT>
312 if (this->getNormCount() == 0) {
313 for (
int i=0; i < iterations; ++i) {
315 mLeafs->rebuildLeafArray();
319 for (
int i=0; i < iterations; ++i) {
320 BoolMaskType mask0(mGrid->tree(),
false,
TopologyCopy());
322 mLeafs->rebuildLeafArray();
324 BoolMaskType mask(mGrid->tree(),
false,
TopologyCopy());
325 mask.topologyDifference(mask0);
331 template<
typename Gr
idT,
typename InterruptT>
337 mLeafs->rebuildLeafArray();
338 const ValueType background = mGrid->background() - iterations*mDx;
342 template<
typename Gr
idT,
typename InterruptT>
347 const int wOld =
static_cast<int>(
math::RoundDown(this->getHalfWidth()));
348 const int wNew =
static_cast<int>(halfWidth);
350 this->dilate(wNew - wOld);
351 }
else if (wOld > wNew) {
352 this->erode(wOld - wNew);
357 template<
typename Gr
idT,
typename InterruptT>
362 if (mInterrupter) mInterrupter->start(msg);
365 template<
typename Gr
idT,
typename InterruptT>
370 if (mInterrupter) mInterrupter->end();
373 template<
typename Gr
idT,
typename InterruptT>
379 tbb::task::self().cancel_group_execution();
385 template<
typename Gr
idT,
typename InterruptT>
386 template<
typename MaskT>
391 switch (this->getSpatialScheme()) {
393 this->normalize1<math::FIRST_BIAS , MaskT>(mask);
break;
395 this->normalize1<math::SECOND_BIAS, MaskT>(mask);
break;
397 this->normalize1<math::THIRD_BIAS, MaskT>(mask);
break;
399 this->normalize1<math::WENO5_BIAS, MaskT>(mask);
break;
401 this->normalize1<math::HJWENO5_BIAS, MaskT>(mask);
break;
407 template<
typename Gr
idT,
typename InterruptT>
408 template<math::BiasedGradientScheme SpatialScheme,
typename MaskT>
413 switch (this->getTemporalScheme()) {
415 this->normalize2<SpatialScheme, math::TVD_RK1, MaskT>(mask);
break;
417 this->normalize2<SpatialScheme, math::TVD_RK2, MaskT>(mask);
break;
419 this->normalize2<SpatialScheme, math::TVD_RK3, MaskT>(mask);
break;
425 template<
typename Gr
idT,
typename InterruptT>
430 LevelSetTracker<GridT, InterruptT>::
431 normalize2(
const MaskT* mask)
433 Normalizer<SpatialScheme, TemporalScheme, MaskT> tmp(*
this, mask);
439 template<
typename Gr
idT,
typename InterruptT>
441 LevelSetTracker<GridT, InterruptT>::
444 const int grainSize = mTracker.getGrainSize();
445 const LeafRange range = mTracker.leafs().leafRange(grainSize);
448 tbb::parallel_for(range, *
this);
455 template<
typename Gr
idT,
typename InterruptT>
457 LevelSetTracker<GridT, InterruptT>::
458 Trim::operator()(
const LeafRange& range)
const
460 typedef typename LeafType::ValueOnIter VoxelIterT;
461 mTracker.checkInterrupter();
462 const ValueType gamma = mTracker.mGrid->background();
464 for (
typename LeafRange::Iterator leafIter = range.begin(); leafIter; ++leafIter) {
465 LeafType &leaf = *leafIter;
466 for (VoxelIterT iter = leaf.beginValueOn(); iter; ++iter) {
467 const ValueType val = *iter;
469 leaf.setValueOff(iter.pos(), -gamma);
470 else if (val >= gamma)
471 leaf.setValueOff(iter.pos(), gamma);
478 template<
typename Gr
idT,
typename InterruptT>
483 LevelSetTracker<GridT, InterruptT>::
484 Normalizer<SpatialScheme, TemporalScheme, MaskT>::
485 Normalizer(LevelSetTracker& tracker,
const MaskT* mask)
488 , mDt(tracker.voxelSize()*(TemporalScheme == math::
TVD_RK1 ? 0.3f :
489 TemporalScheme == math::
TVD_RK2 ? 0.9f : 1.0f))
490 , mInvDx(1.0f/tracker.voxelSize())
495 template<
typename Gr
idT,
typename InterruptT>
505 mTracker.mLeafs->rebuildAuxBuffers(TemporalScheme ==
math::TVD_RK3 ? 2 : 1);
507 for (
int n=0, e=mTracker.getNormCount(); n < e; ++n) {
510 switch(TemporalScheme) {
514 mTask = boost::bind(&Normalizer::euler01, _1, _2);
522 mTask = boost::bind(&Normalizer::euler01, _1, _2);
529 mTask = boost::bind(&Normalizer::euler12, _1, _2);
537 mTask = boost::bind(&Normalizer::euler01, _1, _2);
544 mTask = boost::bind(&Normalizer::euler34, _1, _2);
551 mTask = boost::bind(&Normalizer::euler13, _1, _2);
557 OPENVDB_THROW(ValueError,
"Temporal integration scheme not supported!");
561 mTracker.mLeafs->removeAuxBuffers();
566 template<
typename Gr
idT,
typename InterruptT>
571 LevelSetTracker<GridT, InterruptT>::
572 Normalizer<SpatialScheme, TemporalScheme, MaskT>::
575 mTracker.startInterrupter(
"Normalizing Level Set");
577 const int grainSize = mTracker.getGrainSize();
578 const LeafRange range = mTracker.leafs().leafRange(grainSize);
581 tbb::parallel_for(range, *
this);
586 mTracker.leafs().swapLeafBuffer(swapBuffer, grainSize==0);
588 mTracker.endInterrupter();
591 template<
typename Gr
idT,
typename InterruptT>
595 template <
int Nominator,
int Denominator>
597 LevelSetTracker<GridT, InterruptT>::
598 Normalizer<SpatialScheme, TemporalScheme, MaskT>::
599 eval(StencilT& stencil,
const ValueType* phi, ValueType* result,
Index n)
const
601 typedef typename math::ISGradientNormSqrd<SpatialScheme> GradientT;
602 static const ValueType alpha = ValueType(Nominator)/ValueType(Denominator);
603 static const ValueType beta = ValueType(1) - alpha;
605 const ValueType normSqGradPhi = GradientT::result(stencil);
606 const ValueType phi0 = stencil.getValue();
609 v = phi0 - mDt * v * (
math::Sqrt(normSqGradPhi) * mInvDx - 1.0f);
610 result[n] = Nominator ? alpha * phi[n] + beta * v : v;
613 template<
typename Gr
idT,
typename InterruptT>
617 template <
int Nominator,
int Denominator>
619 LevelSetTracker<GridT,InterruptT>::
620 Normalizer<SpatialScheme, TemporalScheme, MaskT>::
621 euler(
const LeafRange& range,
Index phiBuffer,
Index resultBuffer)
623 typedef typename LeafType::ValueOnCIter VoxelIterT;
625 mTracker.checkInterrupter();
627 StencilT stencil(mTracker.grid());
629 for (
typename LeafRange::Iterator leafIter = range.begin(); leafIter; ++leafIter) {
630 const ValueType* phi = leafIter.buffer(phiBuffer).data();
631 ValueType* result = leafIter.buffer(resultBuffer).data();
633 for (VoxelIterT iter = leafIter->cbeginValueOn(); iter; ++iter) {
634 stencil.moveTo(iter);
635 this->eval<Nominator, Denominator>(stencil, phi, result, iter.pos());
637 }
else if (
const MaskLeafT* mask = mMask->probeLeaf(leafIter->origin())) {
638 const ValueType* phi0 = leafIter->buffer().data();
639 for (MaskIterT iter = mask->cbeginValueOn(); iter; ++iter) {
640 const Index i = iter.pos();
641 stencil.moveTo(iter.getCoord(), phi0[i]);
642 this->eval<Nominator, Denominator>(stencil, phi, result, i);
652 #endif // OPENVDB_TOOLS_LEVEL_SET_TRACKER_HAS_BEEN_INCLUDED
static T value()
Definition: Math.h:125
Index32 Index
Definition: Types.h:58
bool wasInterrupted(T *i, int percent=-1)
Definition: NullInterrupter.h:76
Definition: FiniteDifference.h:198
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
Definition: FiniteDifference.h:263
Definition: Exceptions.h:88
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:97
Efficient multi-threaded replacement of the background values in tree.
Definition: FiniteDifference.h:197
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:114
Definition: Operators.h:152
float RoundDown(float x)
Return x rounded down to the nearest integer.
Definition: Math.h:751
Definition: FiniteDifference.h:195
Defined various multi-threaded utility functions for trees.
#define OPENVDB_VERSION_NAME
Definition: version.h:43
Definition: Exceptions.h:86
Implementation of morphological dilation and erosion.
Definition: FiniteDifference.h:196
Definition: Exceptions.h:39
CopyConstness< TreeType, NonConstBufferType >::Type BufferType
Definition: LeafManager.h:125
Type Pow2(Type x)
Return .
Definition: Math.h:514
Definition: FiniteDifference.h:194
BiasedGradientScheme
Biased Gradients are limited to non-centered differences.
Definition: FiniteDifference.h:192
Definition: FiniteDifference.h:265
Defines various finite difference stencils by means of the "curiously recurring template pattern" on ...
Definition: FiniteDifference.h:264
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:71
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:709
A LeafManager manages a linear array of pointers to a given tree's leaf nodes, as well as optional au...
static const Real LEVEL_SET_HALF_WIDTH
Definition: Types.h:212
TemporalIntegrationScheme
Temporal integration schemes.
Definition: FiniteDifference.h:261
Definition: LeafManager.h:131