98 #ifndef OPENVDB_TOOLS_PARTICLES_TO_LEVELSET_HAS_BEEN_INCLUDED 99 #define OPENVDB_TOOLS_PARTICLES_TO_LEVELSET_HAS_BEEN_INCLUDED 101 #include <tbb/parallel_reduce.h> 102 #include <tbb/blocked_range.h> 113 #include <type_traits> 121 namespace p2ls_internal {
125 template<
typename VisibleT,
typename BlindT>
class BlindData;
129 template<
typename SdfGridT,
130 typename AttributeT = void,
135 using DisableT =
typename std::is_void<AttributeT>::type;
141 using AttType =
typename std::conditional<DisableT::value, size_t, AttributeT>::type;
142 using AttGridType =
typename SdfGridT::template ValueConverter<AttType>::Type;
144 static_assert(std::is_floating_point<SdfType>::value,
145 "ParticlesToLevelSet requires an SDF grid with floating-point values");
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 using BlindGridType =
typename SdfGridT::template ValueConverter<BlindType>::Type;
258 template<
typename ParticleListT,
typename Gr
idT>
struct Raster;
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>
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 ==
nullptr) {
354 using SdfTreeT =
typename SdfGridType::TreeType;
355 using AttTreeT =
typename AttGridType::TreeType;
356 using BlindTreeT =
typename BlindGridType::TreeType;
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 using LeafIterT =
typename BlindTreeT::LeafCIter;
374 using LeafT =
typename BlindTreeT::LeafNodeType;
375 using SdfLeafT =
typename SdfTreeT::LeafNodeType;
376 using AttLeafT =
typename AttTreeT::LeafNodeType;
377 for (LeafIterT n = tree.cbeginLeaf(); n; ++n) {
378 const LeafT& leaf = *n;
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 using DisableT =
typename std::is_void<AttributeT>::type;
418 using SdfT =
typename ParticlesToLevelSetT::SdfType;
419 using AttT =
typename ParticlesToLevelSetT::AttType;
420 using ValueT =
typename GridT::ValueType;
421 using AccessorT =
typename GridT::Accessor;
422 using TreeT =
typename GridT::TreeType;
423 using LeafNodeT =
typename TreeT::LeafNodeType;
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 = std::bind(&Raster::rasterSpheres, std::placeholders::_1, std::placeholders::_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 = std::bind(&Raster::rasterFixedSpheres,
495 std::placeholders::_1, std::placeholders::_2, SdfT(radius));
497 if (mParent.mInterrupter) mParent.mInterrupter->end();
516 mMinCount = mMaxCount = 0;
517 if (mParent.mInterrupter) {
518 mParent.mInterrupter->start(
"Rasterizing particles to level set using trails");
520 mTask = std::bind(&Raster::rasterTrails,
521 std::placeholders::_1, std::placeholders::_2, SdfT(delta));
523 if (mParent.mInterrupter) mParent.mInterrupter->end();
527 void operator()(
const tbb::blocked_range<size_t>& r)
531 mParent.mMinCount = mMinCount;
532 mParent.mMaxCount = mMaxCount;
536 void join(Raster& other)
539 mMinCount += other.mMinCount;
540 mMaxCount += other.mMaxCount;
544 Raster& operator=(
const Raster&) {
return *
this; }
547 bool ignoreParticle(SdfT R)
549 if (R < mParent.mRmin) {
553 if (R > mParent.mRmax) {
563 void rasterSpheres(
const tbb::blocked_range<size_t>& r)
565 AccessorT acc = mGrid->getAccessor();
567 const SdfT invDx = SdfT(1/mParent.mDx);
573 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
575 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
576 for ( ; run && iter; ++iter) {
578 mParticles.getPosRad(
id, pos, rad);
579 const SdfT R = SdfT(invDx * rad);
580 if (this->ignoreParticle(R))
continue;
581 const Vec3R P = mMap.applyInverseMap(pos);
582 this->getAtt<DisableT>(id, att);
583 run = this->makeSphere(P, R, att, acc);
592 void rasterFixedSpheres(
const tbb::blocked_range<size_t>& r, SdfT R)
595 dx =
static_cast<SdfT
>(mParent.mDx),
596 w = static_cast<SdfT>(mParent.mHalfWidth);
597 AccessorT acc = mGrid->getAccessor();
598 const ValueT inside = -mGrid->background();
599 const SdfT
max = R + w;
608 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
610 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
611 for ( ; iter; ++iter) {
613 this->getAtt<DisableT>(id, att);
614 mParticles.getPos(
id, pos);
615 const Vec3R P = mMap.applyInverseMap(pos);
618 for (
Coord c = a; c.
x() <= b.
x(); ++c.x()) {
621 tbb::task::self().cancel_group_execution();
624 SdfT x2 =
static_cast<SdfT
>(
math::Pow2(c.x() - P[0]));
625 for (c.y() = a.
y(); c.y() <= b.
y(); ++c.y()) {
626 SdfT x2y2 =
static_cast<SdfT
>(x2 +
math::Pow2(c.y() - P[1]));
627 for (c.z() = a.
z(); c.z() <= b.
z(); ++c.z()) {
628 SdfT x2y2z2 =
static_cast<SdfT
>(
630 if (x2y2z2 >= max2 || (!acc.probeValue(c,v) && v<ValueT(0)))
632 if (x2y2z2 <= min2) {
633 acc.setValueOff(c, inside);
637 const ValueT d=Merge(dx*(
math::Sqrt(x2y2z2) - R), att);
638 if (d < v) acc.setValue(c, d);
650 void rasterTrails(
const tbb::blocked_range<size_t>& r, SdfT delta)
652 AccessorT acc = mGrid->getAccessor();
657 const Vec3R origin = mMap.applyInverseMap(
Vec3R(0,0,0));
658 const SdfT Rmin = SdfT(mParent.mRmin), invDx = SdfT(1/mParent.mDx);
661 for (
size_t n = r.begin(), N = r.end(); n < N; ++n) {
663 typename PointPartitionerT::IndexIterator iter = mPointPartitioner->indices(n);
664 for ( ; run && iter; ++iter) {
666 mParticles.getPosRadVel(
id, pos, rad, vel);
667 const SdfT R0 = SdfT(invDx*rad);
668 if (this->ignoreParticle(R0))
continue;
669 this->getAtt<DisableT>(id, att);
670 const Vec3R P0 = mMap.applyInverseMap(pos);
671 const Vec3R V = mMap.applyInverseMap(vel) - origin;
672 const SdfT speed = SdfT(V.
length()), inv_speed = SdfT(1.0/speed);
673 const Vec3R Nrml = -V*inv_speed;
676 for (
size_t m=0; run && d <= speed ; ++m) {
677 run = this->makeSphere(P, R, att, acc);
678 P += 0.5*delta*R*Nrml;
679 d = SdfT((P-P0).length());
680 R = R0-(R0-Rmin)*d*inv_speed;
691 if (mParent.mGrainSize>0) {
692 tbb::parallel_reduce(
693 tbb::blocked_range<size_t>(0, bucketCount, mParent.mGrainSize), *
this);
695 (*this)(tbb::blocked_range<size_t>(0, bucketCount));
714 bool makeSphere(
const Vec3R &P, SdfT R,
const AttT& att, AccessorT& acc)
716 const ValueT inside = -mGrid->background();
717 const SdfT dx = SdfT(mParent.mDx), w = SdfT(mParent.mHalfWidth);
718 const SdfT
max = R + w;
725 for (
Coord c = a; c.
x() <= b.
x(); ++c.x() ) {
728 tbb::task::self().cancel_group_execution();
732 for (c.y() = a.
y(); c.y() <= b.
y(); ++c.y()) {
733 SdfT x2y2 = SdfT(x2 +
math::Pow2(c.y() - P[1]));
734 for (c.z() = a.
z(); c.z() <= b.
z(); ++c.z()) {
735 SdfT x2y2z2 = SdfT(x2y2 +
math::Pow2(c.z()-P[2]));
736 if (x2y2z2 >= max2 || (!acc.probeValue(c,v) && v<ValueT(0)))
738 if (x2y2z2 <= min2) {
739 acc.setValueOff(c, inside);
744 const ValueT d=Merge(dx*(
math::Sqrt(x2y2z2) - R), att);
745 if (d < v) acc.setValue(c, d);
751 using FuncType =
typename std::function<void (Raster*, const tbb::blocked_range<size_t>&)>;
753 template<
typename DisableType>
754 typename std::enable_if<DisableType::value>::type
755 getAtt(
size_t, AttT&)
const {}
757 template<
typename DisableType>
758 typename std::enable_if<!DisableType::value>::type
759 getAtt(
size_t n, AttT& a)
const { mParticles.getAtt(n, a); }
762 typename std::enable_if<std::is_same<T, ValueT>::value, ValueT>::type
763 Merge(T s,
const AttT&)
const {
return s; }
766 typename std::enable_if<!std::is_same<T, ValueT>::value, ValueT>::type
767 Merge(T s,
const AttT& a)
const {
return ValueT(s,a); }
769 ParticlesToLevelSetT& mParent;
770 const ParticleListT& mParticles;
773 size_t mMinCount, mMaxCount;
776 PointPartitionerT* mPointPartitioner;
782 namespace p2ls_internal {
787 template<
typename VisibleT,
typename BlindT>
797 BlindData(VisibleT v, BlindT b) : mVisible(v), mBlind(b) {}
800 const VisibleT&
visible()
const {
return mVisible; }
801 const BlindT&
blind()
const {
return mBlind; }
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 void setRmin(Real Rmin)
set the smallest radius allowed in voxel units
Definition: ParticlesToLevelSet.h:209
~ParticlesToLevelSet()
Destructor.
Definition: ParticlesToLevelSet.h:171
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
SdfGridT SdfGridType
Definition: ParticlesToLevelSet.h:138
const Type & Max(const Type &a, const Type &b)
Return the maximum of two values.
Definition: Math.h:545
bool operator<(const Tuple< SIZE, T0 > &t0, const Tuple< SIZE, T1 > &t1)
Definition: Tuple.h:206
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
int Ceil(float x)
Return the ceiling of x.
Definition: Math.h:806
void setGrainSize(int grainSize)
Set the grain-size used for multi-threading.
Definition: ParticlesToLevelSet.h:217
Functions to efficiently perform various compositing operations on grids.
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:108
T length() const
Length of the vector.
Definition: Vec3.h:225
Type Pow2(Type x)
Return x2.
Definition: Math.h:498
math::Vec3< Real > Vec3R
Definition: Types.h:78
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:51
Coord Abs(const Coord &xyz)
Definition: Coord.h:513
typename std::conditional< DisableT::value, size_t, AttributeT >::type AttType
Definition: ParticlesToLevelSet.h:141
typename SdfGridT::template ValueConverter< AttType >::Type AttGridType
Definition: ParticlesToLevelSet.h:142
bool operator>(const Tuple< SIZE, T0 > &t0, const Tuple< SIZE, T1 > &t1)
Definition: Tuple.h:218
bool ignoredParticles() const
Return true if any particles were ignored due to their size.
Definition: ParticlesToLevelSet.h:202
Spatially partitions points using a parallel radix-based sorting algorithm.
Defined various multi-threaded utility functions for trees.
T zeroVal()
Return the value of type T that corresponds to zero.
Definition: Math.h:86
Definition: Exceptions.h:89
Real getRmax() const
Return the largest radius allowed in voxel units.
Definition: ParticlesToLevelSet.h:199
Real getRmin() const
Return the smallest radius allowed in voxel units.
Definition: ParticlesToLevelSet.h:197
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
uint32_t Index32
Definition: Types.h:58
Index32 Index
Definition: Types.h:60
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:136
typename std::is_void< AttributeT >::type DisableT
Definition: ParticlesToLevelSet.h:135
Propagate the signs of distance values from the active voxels in the narrow band to the inactive valu...
void setRmax(Real Rmax)
set the largest radius allowed in voxel units
Definition: ParticlesToLevelSet.h:211
Int32 z() const
Definition: Coord.h:159
Tag dispatch class that distinguishes shallow copy constructors from deep copy constructors.
Definition: Types.h:505
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:711
Definition: Exceptions.h:39
InterrupterT InterrupterType
Definition: ParticlesToLevelSet.h:136
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
double Real
Definition: Types.h:66
Real getHalfWidth() const
Return the half-width of the narrow band in voxel units.
Definition: ParticlesToLevelSet.h:194
Int32 y() const
Definition: Coord.h:158
#define OPENVDB_NO_FP_EQUALITY_WARNING_END
Definition: Math.h:75
AttGridType::Ptr attributeGrid()
Return a shared pointer to the grid containing the (optional) attribute.
Definition: ParticlesToLevelSet.h:188
Definition: ParticlesToLevelSet.h:132
int getGrainSize() const
Returns the grain-size used for multi-threading.
Definition: ParticlesToLevelSet.h:214
size_t getMinCount() const
Return number of small particles that were ignore due to Rmin.
Definition: ParticlesToLevelSet.h:204
bool wasInterrupted(T *i, int percent=-1)
Definition: NullInterrupter.h:76
#define OPENVDB_NO_FP_EQUALITY_WARNING_BEGIN
Definition: Math.h:74
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:508
Abstract base class for maps.
Definition: Maps.h:161
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
int Floor(float x)
Return the floor of x.
Definition: Math.h:798
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:188
typename SdfGridT::ValueType SdfType
Definition: ParticlesToLevelSet.h:139
Real getVoxelSize() const
Return the size of a voxel in world units.
Definition: ParticlesToLevelSet.h:191
Int32 x() const
Definition: Coord.h:157
size_t getMaxCount() const
Return number of large particles that were ignore due to Rmax.
Definition: ParticlesToLevelSet.h:206