97 #ifndef OPENVDB_TOOLS_PARTICLES_TO_LEVELSET_HAS_BEEN_INCLUDED
98 #define OPENVDB_TOOLS_PARTICLES_TO_LEVELSET_HAS_BEEN_INCLUDED
100 #include <tbb/parallel_reduce.h>
101 #include <tbb/blocked_range.h>
102 #include <boost/bind.hpp>
103 #include <boost/function.hpp>
104 #include <boost/type_traits/is_floating_point.hpp>
105 #include <boost/utility/enable_if.hpp>
106 #include <boost/mpl/if.hpp>
107 #include <openvdb/Types.h>
108 #include <openvdb/Grid.h>
109 #include <openvdb/math/Math.h>
110 #include <openvdb/math/Transform.h>
111 #include <openvdb/util/NullInterrupter.h>
122 namespace p2ls_internal {
126 template<
typename VisibleT,
typename BlindT>
class BlindData;
130 template<
typename SdfGridT,
131 typename AttributeT = void,
136 typedef typename boost::is_void<AttributeT>::type
DisableT;
142 typedef typename boost::mpl::if_<DisableT, size_t, AttributeT>::type
AttType;
143 typedef typename SdfGridT::template ValueConverter<AttType>::Type
AttGridType;
145 BOOST_STATIC_ASSERT(boost::is_floating_point<SdfType>::value);
181 void finalize(
bool prune =
false);
223 template <
typename ParticleListT>
224 void rasterizeSpheres(
const ParticleListT& pa);
231 template <
typename ParticleListT>
232 void rasterizeSpheres(
const ParticleListT& pa,
Real radius);
250 template <
typename ParticleListT>
251 void rasterizeTrails(
const ParticleListT& pa,
Real delta=1.0);
255 typedef typename SdfGridT::template ValueConverter<BlindType>::Type BlindGridType;
258 template<
typename ParticleListT,
typename Gr
idT>
struct Raster;
260 SdfGridType* mSdfGrid;
261 typename AttGridType::Ptr mAttGrid;
262 BlindGridType* mBlindGrid;
263 InterrupterT* mInterrupter;
264 Real mDx, mHalfWidth;
266 size_t mMinCount, mMaxCount;
271 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
272 inline ParticlesToLevelSet<SdfGridT, AttributeT, InterrupterT>::
273 ParticlesToLevelSet(SdfGridT& grid, InterrupterT* interrupter) :
276 mInterrupter(interrupter),
277 mDx(grid.voxelSize()[0]),
278 mHalfWidth(grid.background()/mDx),
285 if (!mSdfGrid->hasUniformVoxels() ) {
287 "ParticlesToLevelSet only supports uniform voxels!");
291 "ParticlesToLevelSet only supports level sets!"
292 "\nUse Grid::setGridClass(openvdb::GRID_LEVEL_SET)");
295 if (!DisableT::value) {
296 mBlindGrid =
new BlindGridType(
BlindType(grid.background()));
297 mBlindGrid->setTransform(mSdfGrid->transform().copy());
301 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
302 template <
typename ParticleListT>
306 if (DisableT::value) {
307 Raster<ParticleListT, SdfGridT> r(*
this, mSdfGrid, pa);
308 r.rasterizeSpheres();
310 Raster<ParticleListT, BlindGridType> r(*
this, mBlindGrid, pa);
311 r.rasterizeSpheres();
315 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
316 template <
typename ParticleListT>
320 if (DisableT::value) {
321 Raster<ParticleListT, SdfGridT> r(*
this, mSdfGrid, pa);
322 r.rasterizeSpheres(radius/mDx);
324 Raster<ParticleListT, BlindGridType> r(*
this, mBlindGrid, pa);
325 r.rasterizeSpheres(radius/mDx);
329 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
330 template <
typename ParticleListT>
334 if (DisableT::value) {
335 Raster<ParticleListT, SdfGridT> r(*
this, mSdfGrid, pa);
336 r.rasterizeTrails(delta);
338 Raster<ParticleListT, BlindGridType> r(*
this, mBlindGrid, pa);
339 r.rasterizeTrails(delta);
343 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
347 if (mBlindGrid==NULL) {
354 typedef typename SdfGridType::TreeType SdfTreeT;
355 typedef typename AttGridType::TreeType AttTreeT;
356 typedef typename BlindGridType::TreeType BlindTreeT;
358 const BlindTreeT& tree = mBlindGrid->tree();
361 typename SdfTreeT::Ptr sdfTree(
new SdfTreeT(
365 typename AttTreeT::Ptr attTree(
new AttTreeT(
367 mAttGrid =
typename AttGridType::Ptr(
new AttGridType(attTree));
368 mAttGrid->setTransform(mBlindGrid->transform().copy());
373 typedef typename BlindTreeT::LeafCIter LeafIterT;
374 typedef typename BlindTreeT::LeafNodeType LeafT;
375 typedef typename SdfTreeT::LeafNodeType SdfLeafT;
376 typedef typename AttTreeT::LeafNodeType AttLeafT;
377 for (LeafIterT n = tree.cbeginLeaf(); n; ++n) {
378 const LeafT& leaf = *n;
379 const openvdb::Coord xyz = leaf.origin();
381 SdfLeafT* sdfLeaf = sdfTree->probeLeaf(xyz);
382 AttLeafT* attLeaf = attTree->probeLeaf(xyz);
384 typename LeafT::ValueOnCIter m=leaf.cbeginValueOn();
388 sdfLeaf->setValueOnly(k, v.
visible());
389 attLeaf->setValueOnly(k, v.
blind());
395 sdfLeaf->setValueOnly(k, v.
visible());
396 attLeaf->setValueOnly(k, v.
blind());
403 if (mSdfGrid->empty()) {
404 mSdfGrid->setTree(sdfTree);
412 template<
typename SdfGr
idT,
typename AttributeT,
typename InterrupterT>
413 template<
typename ParticleListT,
typename Gr
idT>
416 typedef typename boost::is_void<AttributeT>::type
DisableT;
418 typedef typename ParticlesToLevelSetT::SdfType SdfT;
419 typedef typename ParticlesToLevelSetT::AttType AttT;
420 typedef typename GridT::ValueType ValueT;
421 typedef typename GridT::Accessor AccessorT;
422 typedef typename GridT::TreeType TreeT;
423 typedef typename TreeT::LeafNodeType LeafNodeT;
428 Raster(ParticlesToLevelSetT& parent, GridT* grid,
const ParticleListT& particles)
430 , mParticles(particles)
432 , mMap(*(mGrid->transform().baseMap()))
437 mPointPartitioner =
new PointPartitionerT();
438 mPointPartitioner->construct(particles, mGrid->transform());
442 Raster(Raster& other, tbb::split)
443 : mParent(other.mParent)
444 , mParticles(other.mParticles)
451 , mPointPartitioner(other.mPointPartitioner)
463 delete mPointPartitioner;
471 mMinCount = mMaxCount = 0;
472 if (mParent.mInterrupter) {
473 mParent.mInterrupter->start(
"Rasterizing particles to level set using spheres");
475 mTask = boost::bind(&Raster::rasterSpheres, _1, _2);
477 if (mParent.mInterrupter) mParent.mInterrupter->end();
484 mMinCount = radius < mParent.mRmin ? mParticles.size() : 0;
485 mMaxCount = radius > mParent.mRmax ? mParticles.size() : 0;
486 if (mMinCount>0 || mMaxCount>0) {
487 mParent.mMinCount = mMinCount;
488 mParent.mMaxCount = mMaxCount;
490 if (mParent.mInterrupter) {
491 mParent.mInterrupter->start(
492 "Rasterizing particles to level set using const spheres");
494 mTask = boost::bind(&Raster::rasterFixedSpheres, _1, _2, SdfT(radius));
496 if (mParent.mInterrupter) mParent.mInterrupter->end();
515 mMinCount = mMaxCount = 0;
516 if (mParent.mInterrupter) {
517 mParent.mInterrupter->start(
"Rasterizing particles to level set using trails");
519 mTask = boost::bind(&Raster::rasterTrails, _1, _2, SdfT(delta));
521 if (mParent.mInterrupter) mParent.mInterrupter->end();
525 void operator()(
const tbb::blocked_range<size_t>& r)
529 mParent.mMinCount = mMinCount;
530 mParent.mMaxCount = mMaxCount;
534 void join(Raster& other)
537 mMinCount += other.mMinCount;
538 mMaxCount += other.mMaxCount;
542 Raster& operator=(
const Raster& other) {
return *
this; }
545 bool ignoreParticle(SdfT R)
547 if (R < mParent.mRmin) {
551 if (R > mParent.mRmax) {
561 void rasterSpheres(
const tbb::blocked_range<size_t>& r)
563 AccessorT acc = mGrid->getAccessor();
565 const SdfT invDx = SdfT(1/mParent.mDx);
571 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
573 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
574 for ( ; run && iter; ++iter) {
576 mParticles.getPosRad(
id, pos, rad);
577 const SdfT R = SdfT(invDx * rad);
578 if (this->ignoreParticle(R))
continue;
579 const Vec3R P = mMap.applyInverseMap(pos);
580 this->getAtt<DisableT>(id, att);
581 run = this->makeSphere(P, R, att, acc);
590 void rasterFixedSpheres(
const tbb::blocked_range<size_t>& r, SdfT R)
593 dx =
static_cast<SdfT
>(mParent.mDx),
594 w = static_cast<SdfT>(mParent.mHalfWidth);
595 AccessorT acc = mGrid->getAccessor();
596 const ValueT inside = -mGrid->background();
597 const SdfT
max = R + w;
606 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
608 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
609 for ( ; iter; ++iter) {
611 this->getAtt<DisableT>(id, att);
612 mParticles.getPos(
id, pos);
613 const Vec3R P = mMap.applyInverseMap(pos);
616 for (Coord c = a; c.x() <= b.x(); ++c.x()) {
619 tbb::task::self().cancel_group_execution();
622 SdfT x2 =
static_cast<SdfT
>(
math::Pow2(c.x() - P[0]));
623 for (c.y() = a.y(); c.y() <= b.y(); ++c.y()) {
624 SdfT x2y2 =
static_cast<SdfT
>(x2 +
math::Pow2(c.y() - P[1]));
625 for (c.z() = a.z(); c.z() <= b.z(); ++c.z()) {
626 SdfT x2y2z2 =
static_cast<SdfT
>(
628 if (x2y2z2 >= max2 || (!acc.probeValue(c,v) && v<ValueT(0)))
630 if (x2y2z2 <= min2) {
631 acc.setValueOff(c, inside);
635 const ValueT d=Merge(dx*(
math::Sqrt(x2y2z2) - R), att);
636 if (d < v) acc.setValue(c, d);
648 void rasterTrails(
const tbb::blocked_range<size_t>& r, SdfT delta)
650 AccessorT acc = mGrid->getAccessor();
655 const Vec3R origin = mMap.applyInverseMap(
Vec3R(0,0,0));
656 const SdfT Rmin = SdfT(mParent.mRmin), invDx = SdfT(1/mParent.mDx);
659 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
661 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
662 for ( ; run && iter; ++iter) {
664 mParticles.getPosRadVel(
id, pos, rad, vel);
665 const SdfT R0 = SdfT(invDx*rad);
666 if (this->ignoreParticle(R0))
continue;
667 this->getAtt<DisableT>(id, att);
668 const Vec3R P0 = mMap.applyInverseMap(pos);
669 const Vec3R V = mMap.applyInverseMap(vel) - origin;
670 const SdfT speed = SdfT(V.length()), inv_speed = SdfT(1.0/speed);
671 const Vec3R Nrml = -V*inv_speed;
674 for (
size_t m=0; run && d <= speed ; ++m) {
675 run = this->makeSphere(P, R, att, acc);
676 P += 0.5*delta*R*Nrml;
677 d = SdfT((P-P0).length());
678 R = R0-(R0-Rmin)*d*inv_speed;
689 if (mParent.mGrainSize>0) {
690 tbb::parallel_reduce(
691 tbb::blocked_range<size_t>(0, bucketCount, mParent.mGrainSize), *
this);
693 (*this)(tbb::blocked_range<size_t>(0, bucketCount));
712 bool makeSphere(
const Vec3R &P, SdfT R,
const AttT& att, AccessorT& acc)
714 const ValueT inside = -mGrid->background();
715 const SdfT dx = SdfT(mParent.mDx), w = SdfT(mParent.mHalfWidth);
716 const SdfT max = R + w;
723 for ( Coord c = a; c.x() <= b.x(); ++c.x() ) {
726 tbb::task::self().cancel_group_execution();
730 for (c.y() = a.y(); c.y() <= b.y(); ++c.y()) {
731 SdfT x2y2 = SdfT(x2 +
math::Pow2(c.y() - P[1]));
732 for (c.z() = a.z(); c.z() <= b.z(); ++c.z()) {
733 SdfT x2y2z2 = SdfT(x2y2 +
math::Pow2(c.z()-P[2]));
734 if (x2y2z2 >= max2 || (!acc.probeValue(c,v) && v<ValueT(0)))
736 if (x2y2z2 <= min2) {
737 acc.setValueOff(c, inside);
742 const ValueT d=Merge(dx*(
math::Sqrt(x2y2z2) - R), att);
743 if (d < v) acc.setValue(c, d);
749 typedef typename boost::function<void (Raster*, const tbb::blocked_range<size_t>&)> FuncType;
751 template <
typename DisableType>
752 typename boost::enable_if<DisableType>::type
753 getAtt(
size_t, AttT&)
const {;}
755 template <
typename DisableType>
756 typename boost::disable_if<DisableType>::type
757 getAtt(
size_t n, AttT& a)
const { mParticles.getAtt(n, a); }
759 template <
typename T>
760 typename boost::enable_if<boost::is_same<T,ValueT>, ValueT>::type
761 Merge(T s,
const AttT&)
const {
return s; }
763 template <
typename T>
764 typename boost::disable_if<boost::is_same<T,ValueT>, ValueT>::type
765 Merge(T s,
const AttT& a)
const {
return ValueT(s,a); }
767 ParticlesToLevelSetT& mParent;
768 const ParticleListT& mParticles;
770 const math::MapBase& mMap;
771 size_t mMinCount, mMaxCount;
774 PointPartitionerT* mPointPartitioner;
780 namespace p2ls_internal {
785 template<
typename VisibleT,
typename BlindT>
818 template<
typename VisibleT,
typename BlindT>
819 inline std::ostream& operator<<(std::ostream& ostr, const BlindData<VisibleT, BlindT>& rhs)
821 ostr << rhs.visible();
826 template<
typename VisibleT,
typename BlindT>
840 #endif // OPENVDB_TOOLS_PARTICLES_TO_LEVELSET_HAS_BEEN_INCLUDED
int32_t Abs(int32_t i)
Return the absolute value of the given quantity.
Definition: Math.h:284
~ParticlesToLevelSet()
Destructor.
Definition: ParticlesToLevelSet.h:171
Real getRmin() const
Return the smallest radius allowed in voxel units.
Definition: ParticlesToLevelSet.h:197
Functions to efficiently perform various compositing operations on grids.
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:97
bool ignoredParticles() const
Return true if any particles were ignored due to their size.
Definition: ParticlesToLevelSet.h:202
Type Pow2(Type x)
Return .
Definition: Math.h:498
Real getVoxelSize() const
Return the size of a voxel in world units.
Definition: ParticlesToLevelSet.h:191
double Real
Definition: Types.h:65
int getGrainSize() const
Rreturn the grain-size used for multi-threading.
Definition: ParticlesToLevelSet.h:214
Multi-threaded space-partitioning scheme for points.
Defined various multi-threaded utility functions for trees.
size_t getMaxCount() const
Return number of large particles that were ignore due to Rmax.
Definition: ParticlesToLevelSet.h:206
void rasterizeTrails(const ParticleListT &pa, Real delta=1.0)
Rasterize a trail per particle derived from their position, radius and velocity. Each trail is genera...
Definition: ParticlesToLevelSet.h:332
bool wasInterrupted(T *i, int percent=-1)
Definition: NullInterrupter.h:76
SdfGridT::ValueType SdfType
Definition: ParticlesToLevelSet.h:140
void setRmax(Real Rmax)
set the largest radius allowed in voxel units
Definition: ParticlesToLevelSet.h:211
Definition: Exceptions.h:86
void rasterizeSpheres(const ParticleListT &pa)
Rasterize a sphere per particle derived from their position and radius. All spheres are CSG unioned...
Definition: ParticlesToLevelSet.h:304
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:693
#define OPENVDB_VERSION_NAME
Definition: version.h:43
Index32 Index
Definition: Types.h:59
Propagates the sign of distance values from the active voxels in the narrow band to the inactive valu...
math::Vec3< Real > Vec3R
Definition: Types.h:77
AttGridType::Ptr attributeGrid()
Return a shared pointer to the grid containing the (optional) attribute.
Definition: ParticlesToLevelSet.h:188
Definition: Exceptions.h:39
void setGrainSize(int grainSize)
Set the grain-size used for multi-threading.
Definition: ParticlesToLevelSet.h:217
SdfGridT SdfGridType
Definition: ParticlesToLevelSet.h:139
const Type & Max(const Type &a, const Type &b)
Return the maximum of two values.
Definition: Math.h:545
InterrupterT InterrupterType
Definition: ParticlesToLevelSet.h:137
Real getRmax() const
Return the largest radius allowed in voxel units.
Definition: ParticlesToLevelSet.h:199
void setRmin(Real Rmin)
set the smallest radius allowed in voxel units
Definition: ParticlesToLevelSet.h:209
int Floor(float x)
Return the floor of x.
Definition: Math.h:780
int Ceil(float x)
Return the ceiling of x.
Definition: Math.h:788
Real getHalfWidth() const
Return the half-width of the narrow band in voxel units.
Definition: ParticlesToLevelSet.h:194
OPENVDB_API Hermite max(const Hermite &, const Hermite &)
min and max operations done directly on the compressed data.
boost::mpl::if_< DisableT, size_t, AttributeT >::type AttType
Definition: ParticlesToLevelSet.h:142
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
uint32_t Index32
Definition: Types.h:57
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:71
boost::is_void< AttributeT >::type DisableT
Definition: ParticlesToLevelSet.h:136
T zeroVal()
Return the value of type T that corresponds to zero.
Definition: Math.h:85
size_t getMinCount() const
Return number of small particles that were ignore due to Rmin.
Definition: ParticlesToLevelSet.h:204
SdfGridT::template ValueConverter< AttType >::Type AttGridType
Definition: ParticlesToLevelSet.h:143
Definition: ParticlesToLevelSet.h:133
void finalize(bool prune=false)
This methods syncs up the level set and attribute grids and therefore needs to be called before any o...
Definition: ParticlesToLevelSet.h:345