OpenVDB  3.2.0
ParticleAtlas.h
Go to the documentation of this file.
1 //
3 // Copyright (c) 2012-2016 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 //
51 
52 
53 #ifndef OPENVDB_TOOLS_PARTICLE_ATLAS_HAS_BEEN_INCLUDED
54 #define OPENVDB_TOOLS_PARTICLE_ATLAS_HAS_BEEN_INCLUDED
55 
56 #include "PointIndexGrid.h"
57 
58 #include <openvdb/Grid.h>
59 #include <openvdb/Types.h>
60 #include <openvdb/math/Transform.h>
61 #include <openvdb/tree/Tree.h>
62 #include <openvdb/tree/LeafNode.h>
63 
64 #include <boost/scoped_array.hpp>
65 #include <tbb/blocked_range.h>
66 #include <tbb/parallel_for.h>
67 #include <tbb/parallel_reduce.h>
68 #include <deque>
69 
70 
71 namespace openvdb {
73 namespace OPENVDB_VERSION_NAME {
74 namespace tools {
75 
76 
78 
79 
105 template<typename PointIndexGridType = PointIndexGrid>
107 {
108  typedef boost::shared_ptr<ParticleAtlas> Ptr;
109  typedef boost::shared_ptr<const ParticleAtlas> ConstPtr;
110 
111  typedef typename PointIndexGridType::Ptr PointIndexGridPtr;
112  typedef typename PointIndexGridType::ValueType IndexType;
113 
114  struct Iterator;
115 
117 
118  ParticleAtlas() : mIndexGridArray(), mMinRadiusArray(), mMaxRadiusArray() {}
119 
125  template<typename ParticleArrayType>
126  void construct(const ParticleArrayType& particles, double minVoxelSize, size_t maxLevels = 50);
127 
133  template<typename ParticleArrayType>
134  static Ptr create(const ParticleArrayType& particles, double minVoxelSize, size_t maxLevels = 50);
135 
137  size_t levels() const { return mIndexGridArray.size(); }
139  bool empty() const { return mIndexGridArray.empty(); }
140 
142  double minRadius(size_t n) const { return mMinRadiusArray[n]; }
144  double maxRadius(size_t n) const { return mMaxRadiusArray[n]; }
145 
147  PointIndexGridType& pointIndexGrid(size_t n) { return *mIndexGridArray[n]; }
149  const PointIndexGridType& pointIndexGrid(size_t n) const { return *mIndexGridArray[n]; }
150 
151 private:
152  // Disallow copying
154  ParticleAtlas& operator=(const ParticleAtlas&);
155 
156  std::vector<PointIndexGridPtr> mIndexGridArray;
157  std::vector<double> mMinRadiusArray, mMaxRadiusArray;
158 }; // struct ParticleAtlas
159 
160 
162 
163 
165 
166 
172 template<typename PointIndexGridType>
173 struct ParticleAtlas<PointIndexGridType>::Iterator
174 {
175  typedef typename PointIndexGridType::TreeType TreeType;
177  typedef boost::scoped_ptr<ConstAccessor> ConstAccessorPtr;
178 
180 
182  explicit Iterator(const ParticleAtlas& atlas);
183 
189  template<typename ParticleArrayType>
190  void worldSpaceSearchAndUpdate(const Vec3d& center, double radius, const ParticleArrayType& particles);
191 
196  template<typename ParticleArrayType>
197  void worldSpaceSearchAndUpdate(const BBoxd& bbox, const ParticleArrayType& particles);
198 
200  size_t levels() const { return mAtlas->levels(); }
201 
204  void updateFromLevel(size_t level);
205 
207  void reset();
208 
210  const IndexType& operator*() const { return *mRange.first; }
211 
214  bool test() const { return mRange.first < mRange.second || mIter != mRangeList.end(); }
215  operator bool() const { return this->test(); }
217 
219  void increment();
220 
222  void operator++() { this->increment(); }
223 
226  bool next();
227 
229  size_t size() const;
230 
232  bool operator==(const Iterator& p) const { return mRange.first == p.mRange.first; }
233  bool operator!=(const Iterator& p) const { return !this->operator==(p); }
234 
235 private:
236  Iterator(const Iterator& rhs);
237  Iterator& operator=(const Iterator& rhs);
238 
239  void clear();
240 
241  typedef std::pair<const IndexType*, const IndexType*> Range;
242  typedef std::deque<Range> RangeDeque;
243  typedef typename RangeDeque::const_iterator RangeDequeCIter;
244  typedef boost::scoped_array<IndexType> IndexArray;
245 
246  ParticleAtlas const * const mAtlas;
247  boost::scoped_array<ConstAccessorPtr> mAccessorList;
248 
249  // Primary index collection
250  Range mRange;
251  RangeDeque mRangeList;
252  RangeDequeCIter mIter;
253  // Secondary index collection
254  IndexArray mIndexArray;
255  size_t mIndexArraySize, mAccessorListSize;
256 }; // struct ParticleAtlas::Iterator
257 
258 
260 
261 // Internal operators and implementation details
262 
263 
264 namespace particle_atlas_internal {
265 
266 
267 template<typename ParticleArrayT>
269  typedef typename ParticleArrayT::PosType PosType;
270  typedef typename PosType::value_type ScalarType;
271 
272  ComputeExtremas(const ParticleArrayT& particles)
273  : particleArray(&particles)
274  , minRadius(std::numeric_limits<ScalarType>::max())
275  , maxRadius(-std::numeric_limits<ScalarType>::max())
276  {
277  }
278 
280  : particleArray(rhs.particleArray)
281  , minRadius(std::numeric_limits<ScalarType>::max())
282  , maxRadius(-std::numeric_limits<ScalarType>::max())
283  {
284  }
285 
286  void operator()(const tbb::blocked_range<size_t>& range) {
287 
288  ScalarType radius, tmpMin = minRadius, tmpMax = maxRadius;
289 
290  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
291  particleArray->getRadius(n, radius);
292  tmpMin = std::min(radius, tmpMin);
293  tmpMax = std::max(radius, tmpMax);
294  }
295 
296  minRadius = std::min(minRadius, tmpMin);
297  maxRadius = std::max(maxRadius, tmpMax);
298  }
299 
300  void join(const ComputeExtremas& rhs) {
301  minRadius = std::min(minRadius, rhs.minRadius);
302  maxRadius = std::max(maxRadius, rhs.maxRadius);
303  }
304 
305  ParticleArrayT const * const particleArray;
306  ScalarType minRadius, maxRadius;
307 }; // struct ComputeExtremas
308 
309 
310 template<typename ParticleArrayT, typename PointIndex>
312 {
313  typedef boost::shared_ptr<SplittableParticleArray> Ptr;
314  typedef boost::shared_ptr<const SplittableParticleArray> ConstPtr;
315  typedef ParticleArrayT ParticleArray;
316 
317  typedef typename ParticleArray::PosType PosType;
318  typedef typename PosType::value_type ScalarType;
319 
320  SplittableParticleArray(const ParticleArrayT& particles)
321  : mIndexMap(), mParticleArray(&particles), mSize(particles.size())
322  {
323  updateExtremas();
324  }
325 
326  SplittableParticleArray(const ParticleArrayT& particles, double minR, double maxR)
327  : mIndexMap(), mParticleArray(&particles), mSize(particles.size())
328  {
329  mMinRadius = ScalarType(minR);
330  mMaxRadius = ScalarType(maxR);
331  }
332 
333  const ParticleArrayT& particleArray() const { return *mParticleArray; }
334 
335  size_t size() const { return mSize; }
336 
337  void getPos(size_t n, PosType& xyz) const { return mParticleArray->getPos(getGlobalIndex(n), xyz); }
338  void getRadius(size_t n, ScalarType& radius) const { return mParticleArray->getRadius(getGlobalIndex(n), radius); }
339 
340  ScalarType minRadius() const { return mMinRadius; }
341  ScalarType maxRadius() const { return mMaxRadius; }
342 
343  size_t getGlobalIndex(size_t n) const { return mIndexMap ? size_t(mIndexMap[n]) : n; }
344 
347  Ptr split(ScalarType maxRadiusLimit) {
348 
349  if (mMaxRadius < maxRadiusLimit) return Ptr();
350 
351  boost::scoped_array<bool> mask(new bool[mSize]);
352 
353  tbb::parallel_for(tbb::blocked_range<size_t>(0, mSize),
354  MaskParticles(*this, mask, maxRadiusLimit));
355 
356  Ptr output(new SplittableParticleArray(*this, mask));
357  if (output->size() == 0) return Ptr();
358 
359  size_t newSize = 0;
360  for (size_t n = 0, N = mSize; n < N; ++n) {
361  newSize += size_t(!mask[n]);
362  }
363 
364  boost::scoped_array<PointIndex> newIndexMap(new PointIndex[newSize]);
365 
366  setIndexMap(newIndexMap, mask, false);
367 
368  mSize = newSize;
369  mIndexMap.swap(newIndexMap);
370  updateExtremas();
371 
372  return output;
373  }
374 
375 
376 private:
377  // Disallow copying
380 
381  // Masked copy constructor
382  SplittableParticleArray(const SplittableParticleArray& other, const boost::scoped_array<bool>& mask)
383  : mIndexMap(), mParticleArray(&other.particleArray()), mSize(0)
384  {
385  for (size_t n = 0, N = other.size(); n < N; ++n) {
386  mSize += size_t(mask[n]);
387  }
388 
389  if (mSize != 0) {
390  mIndexMap.reset(new PointIndex[mSize]);
391  other.setIndexMap(mIndexMap, mask, true);
392  }
393 
394  updateExtremas();
395  }
396 
397  struct MaskParticles {
398  MaskParticles(const SplittableParticleArray& particles,
399  const boost::scoped_array<bool>& mask, ScalarType radius)
400  : particleArray(&particles)
401  , particleMask(mask.get())
402  , radiusLimit(radius)
403  {
404  }
405 
406  void operator()(const tbb::blocked_range<size_t>& range) const {
407  const ScalarType maxRadius = radiusLimit;
408  ScalarType radius;
409  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
410  particleArray->getRadius(n, radius);
411  particleMask[n] = !(radius < maxRadius);
412  }
413  }
414 
415  SplittableParticleArray const * const particleArray;
416  bool * const particleMask;
417  ScalarType const radiusLimit;
418  }; // struct MaskParticles
419 
420  inline void updateExtremas() {
422  tbb::parallel_reduce(tbb::blocked_range<size_t>(0, mSize), op);
423  mMinRadius = op.minRadius;
424  mMaxRadius = op.maxRadius;
425  }
426 
427  void setIndexMap(boost::scoped_array<PointIndex>& newIndexMap,
428  const boost::scoped_array<bool>& mask, bool maskValue) const
429  {
430  if (mIndexMap.get() != NULL) {
431  const PointIndex* indices = mIndexMap.get();
432  for (size_t idx = 0, n = 0, N = mSize; n < N; ++n) {
433  if (mask[n] == maskValue) newIndexMap[idx++] = indices[n];
434  }
435  } else {
436  for (size_t idx = 0, n = 0, N = mSize; n < N; ++n) {
437  if (mask[n] == maskValue) newIndexMap[idx++] = PointIndex(n);
438  }
439  }
440  }
441 
442 
444 
445  boost::scoped_array<PointIndex> mIndexMap;
446  ParticleArrayT const * const mParticleArray;
447  size_t mSize;
448  ScalarType mMinRadius, mMaxRadius;
449 }; // struct SplittableParticleArray
450 
451 
452 template<typename ParticleArrayType, typename PointIndexLeafNodeType>
453 struct RemapIndices {
454 
455  RemapIndices(const ParticleArrayType& particles, std::vector<PointIndexLeafNodeType*> nodes)
456  : mParticles(&particles)
457  , mNodes(nodes.empty() ? NULL : &nodes.front())
458  {
459  }
460 
461  void operator()(const tbb::blocked_range<size_t>& range) const {
462  typedef typename PointIndexLeafNodeType::ValueType PointIndexType;
463  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
464 
465  PointIndexLeafNodeType& node = *mNodes[n];
466  const size_t numIndices = node.indices().size();
467 
468  if (numIndices > 0) {
469  PointIndexType* begin = &node.indices().front();
470  const PointIndexType* end = begin + numIndices;
471 
472  while (begin < end) {
473  *begin = PointIndexType(mParticles->getGlobalIndex(*begin));
474  ++begin;
475  }
476  }
477  }
478  }
479 
480  ParticleArrayType const * const mParticles;
481  PointIndexLeafNodeType * const * const mNodes;
482 }; // struct RemapIndices
483 
484 
485 template<typename ParticleArrayType, typename IndexT>
487 {
488  typedef typename ParticleArrayType::PosType PosType;
489  typedef typename PosType::value_type ScalarType;
490 
491  typedef std::pair<const IndexT*, const IndexT*> Range;
492  typedef std::deque<Range> RangeDeque;
493  typedef std::deque<IndexT> IndexDeque;
494 
495  RadialRangeFilter(RangeDeque& ranges, IndexDeque& indices, const PosType& xyz,
496  ScalarType radius, const ParticleArrayType& particles, bool hasUniformRadius = false)
497  : mRanges(ranges)
498  , mIndices(indices)
499  , mCenter(xyz)
500  , mRadius(radius)
501  , mParticles(&particles)
502  , mHasUniformRadius(hasUniformRadius)
503  {
504  if (mHasUniformRadius) {
505  ScalarType uniformRadius;
506  mParticles->getRadius(0, uniformRadius);
507  mRadius = mRadius + uniformRadius;
508  mRadius *= mRadius;
509  }
510  }
511 
512  template <typename LeafNodeType>
513  void filterLeafNode(const LeafNodeType& leaf)
514  {
515  const size_t numIndices = leaf.indices().size();
516  if (numIndices > 0) {
517  const IndexT* begin = &leaf.indices().front();
518  filterVoxel(leaf.origin(), begin, begin + numIndices);
519  }
520  }
521 
522  void filterVoxel(const Coord&, const IndexT* begin, const IndexT* end)
523  {
524  PosType pos;
525 
526  if (mHasUniformRadius) {
527 
528  const ScalarType searchRadiusSqr = mRadius;
529 
530  while (begin < end) {
531  mParticles->getPos(size_t(*begin), pos);
532  const ScalarType distSqr = (mCenter - pos).lengthSqr();
533  if (distSqr < searchRadiusSqr) {
534  mIndices.push_back(*begin);
535  }
536  ++begin;
537  }
538  } else {
539  while (begin < end) {
540  const size_t idx = size_t(*begin);
541  mParticles->getPos(idx, pos);
542 
543  ScalarType radius;
544  mParticles->getRadius(idx, radius);
545 
546  ScalarType searchRadiusSqr = mRadius + radius;
547  searchRadiusSqr *= searchRadiusSqr;
548 
549  const ScalarType distSqr = (mCenter - pos).lengthSqr();
550 
551  if (distSqr < searchRadiusSqr) {
552  mIndices.push_back(*begin);
553  }
554 
555  ++begin;
556  }
557  }
558  }
559 
560 private:
562  RadialRangeFilter& operator=(const RadialRangeFilter&);
563 
564  RangeDeque& mRanges;
565  IndexDeque& mIndices;
566  PosType const mCenter;
567  ScalarType mRadius;
568  ParticleArrayType const * const mParticles;
569  bool const mHasUniformRadius;
570 }; // struct RadialRangeFilter
571 
572 
573 template<typename ParticleArrayType, typename IndexT>
575 {
576  typedef typename ParticleArrayType::PosType PosType;
577  typedef typename PosType::value_type ScalarType;
578 
579  typedef std::pair<const IndexT*, const IndexT*> Range;
580  typedef std::deque<Range> RangeDeque;
581  typedef std::deque<IndexT> IndexDeque;
582 
583  BBoxFilter(RangeDeque& ranges, IndexDeque& indices,
584  const BBoxd& bbox, const ParticleArrayType& particles, bool hasUniformRadius = false)
585  : mRanges(ranges)
586  , mIndices(indices)
587  , mBBox(PosType(bbox.min()), PosType(bbox.max()))
588  , mCenter(mBBox.getCenter())
589  , mParticles(&particles)
590  , mHasUniformRadius(hasUniformRadius)
591  , mUniformRadiusSqr(ScalarType(0.0))
592  {
593  if (mHasUniformRadius) {
594  mParticles->getRadius(0, mUniformRadiusSqr);
595  mUniformRadiusSqr *= mUniformRadiusSqr;
596  }
597  }
598 
599  template <typename LeafNodeType>
600  void filterLeafNode(const LeafNodeType& leaf)
601  {
602  const size_t numIndices = leaf.indices().size();
603  if (numIndices > 0) {
604  const IndexT* begin = &leaf.indices().front();
605  filterVoxel(leaf.origin(), begin, begin + numIndices);
606  }
607  }
608 
609  void filterVoxel(const Coord&, const IndexT* begin, const IndexT* end)
610  {
611  PosType pos;
612 
613  if (mHasUniformRadius) {
614  const ScalarType radiusSqr = mUniformRadiusSqr;
615 
616  while (begin < end) {
617 
618  mParticles->getPos(size_t(*begin), pos);
619  if (mBBox.isInside(pos)) {
620  mIndices.push_back(*begin++);
621  continue;
622  }
623 
624  const ScalarType distSqr = pointToBBoxDistSqr(pos);
625  if (!(distSqr > radiusSqr)) {
626  mIndices.push_back(*begin);
627  }
628 
629  ++begin;
630  }
631 
632  } else {
633  while (begin < end) {
634 
635  const size_t idx = size_t(*begin);
636  mParticles->getPos(idx, pos);
637  if (mBBox.isInside(pos)) {
638  mIndices.push_back(*begin++);
639  continue;
640  }
641 
642  ScalarType radius;
643  mParticles->getRadius(idx, radius);
644  const ScalarType distSqr = pointToBBoxDistSqr(pos);
645  if (!(distSqr > (radius * radius))) {
646  mIndices.push_back(*begin);
647  }
648 
649  ++begin;
650  }
651  }
652  }
653 
654 private:
655  BBoxFilter(const BBoxFilter&);
656  BBoxFilter& operator=(const BBoxFilter&);
657 
658  ScalarType pointToBBoxDistSqr(const PosType& pos) const
659  {
660  ScalarType distSqr = ScalarType(0.0);
661 
662  for (int i = 0; i < 3; ++i) {
663 
664  const ScalarType a = pos[i];
665 
666  ScalarType b = mBBox.min()[i];
667  if (a < b) {
668  ScalarType delta = b - a;
669  distSqr += delta * delta;
670  }
671 
672  b = mBBox.max()[i];
673  if (a > b) {
674  ScalarType delta = a - b;
675  distSqr += delta * delta;
676  }
677  }
678 
679  return distSqr;
680  }
681 
682  RangeDeque& mRanges;
683  IndexDeque& mIndices;
684  math::BBox<PosType> const mBBox;
685  PosType const mCenter;
686  ParticleArrayType const * const mParticles;
687  bool const mHasUniformRadius;
688  ScalarType mUniformRadiusSqr;
689 }; // struct BBoxFilter
690 
691 
692 } // namespace particle_atlas_internal
693 
694 
696 
697 
698 template<typename PointIndexGridType>
699 template<typename ParticleArrayType>
700 inline void
702  const ParticleArrayType& particles, double minVoxelSize, size_t maxLevels)
703 {
704  typedef typename particle_atlas_internal::
705  SplittableParticleArray<ParticleArrayType, IndexType> SplittableParticleArray;
706  typedef typename SplittableParticleArray::Ptr SplittableParticleArrayPtr;
707  typedef typename ParticleArrayType::ScalarType ScalarType;
708 
710 
712  tbb::parallel_reduce(tbb::blocked_range<size_t>(0, particles.size()), extremas);
713  const double firstMin = extremas.minRadius;
714  const double firstMax = extremas.maxRadius;
715  const double firstVoxelSize = std::max(minVoxelSize, firstMin);
716 
717  if (!(firstMax < (firstVoxelSize * double(2.0))) && maxLevels > 1) {
718 
719  std::vector<SplittableParticleArrayPtr> levels;
720  levels.push_back(SplittableParticleArrayPtr(
721  new SplittableParticleArray(particles, firstMin, firstMax)));
722 
723  std::vector<double> voxelSizeArray;
724  voxelSizeArray.push_back(firstVoxelSize);
725 
726  for (size_t n = 0; n < maxLevels; ++n) {
727 
728  const double maxParticleRadius = double(levels.back()->maxRadius());
729  const double particleRadiusLimit = voxelSizeArray.back() * double(2.0);
730  if (maxParticleRadius < particleRadiusLimit) break;
731 
732  SplittableParticleArrayPtr newLevel = levels.back()->split(ScalarType(particleRadiusLimit));
733  if (!newLevel) break;
734 
735  levels.push_back(newLevel);
736  voxelSizeArray.push_back(double(newLevel->minRadius()));
737  }
738 
739  size_t numPoints = 0;
740 
741  typedef typename PointIndexGridType::TreeType PointIndexTreeType;
742  typedef typename PointIndexTreeType::LeafNodeType PointIndexLeafNodeType;
743 
744  std::vector<PointIndexLeafNodeType*> nodes;
745 
746  for (size_t n = 0, N = levels.size(); n < N; ++n) {
747 
748  const SplittableParticleArray& particleArray = *levels[n];
749 
750  numPoints += particleArray.size();
751 
752  mMinRadiusArray.push_back(double(particleArray.minRadius()));
753  mMaxRadiusArray.push_back(double(particleArray.maxRadius()));
754 
755  PointIndexGridPtr grid = createPointIndexGrid<PointIndexGridType>(particleArray, voxelSizeArray[n]);
756 
757  nodes.clear();
758  grid->tree().getNodes(nodes);
759 
760  tbb::parallel_for(tbb::blocked_range<size_t>(0, nodes.size()),
762  (particleArray, nodes));
763 
764  mIndexGridArray.push_back(grid);
765  }
766 
767  } else {
768  mMinRadiusArray.push_back(firstMin);
769  mMaxRadiusArray.push_back(firstMax);
770  mIndexGridArray.push_back(
771  createPointIndexGrid<PointIndexGridType>(particles, firstVoxelSize));
772  }
773 }
774 
775 
776 template<typename PointIndexGridType>
777 template<typename ParticleArrayType>
780  const ParticleArrayType& particles, double minVoxelSize, size_t maxLevels)
781 {
782  Ptr ret(new ParticleAtlas());
783  ret->construct(particles, minVoxelSize, maxLevels);
784  return ret;
785 }
786 
787 
789 
790 // ParticleAtlas::Iterator implementation
791 
792 template<typename PointIndexGridType>
793 inline
795  : mAtlas(&atlas)
796  , mAccessorList()
797  , mRange(static_cast<IndexType*>(NULL), static_cast<IndexType*>(NULL))
798  , mRangeList()
799  , mIter(mRangeList.begin())
800  , mIndexArray()
801  , mIndexArraySize(0)
802  , mAccessorListSize(atlas.levels())
803 {
804  if (mAccessorListSize > 0) {
805  mAccessorList.reset(new ConstAccessorPtr[mAccessorListSize]);
806  for (size_t n = 0, N = mAccessorListSize; n < N; ++n) {
807  mAccessorList[n].reset(new ConstAccessor(atlas.pointIndexGrid(n).tree()));
808  }
809  }
810 }
811 
812 
813 template<typename PointIndexGridType>
814 inline void
816 {
817  mIter = mRangeList.begin();
818  if (!mRangeList.empty()) {
819  mRange = mRangeList.front();
820  } else if (mIndexArray) {
821  mRange.first = mIndexArray.get();
822  mRange.second = mRange.first + mIndexArraySize;
823  } else {
824  mRange.first = static_cast<IndexType*>(NULL);
825  mRange.second = static_cast<IndexType*>(NULL);
826  }
827 }
828 
829 
830 template<typename PointIndexGridType>
831 inline void
833 {
834  ++mRange.first;
835  if (mRange.first >= mRange.second && mIter != mRangeList.end()) {
836  ++mIter;
837  if (mIter != mRangeList.end()) {
838  mRange = *mIter;
839  } else if (mIndexArray) {
840  mRange.first = mIndexArray.get();
841  mRange.second = mRange.first + mIndexArraySize;
842  }
843  }
844 }
845 
846 
847 template<typename PointIndexGridType>
848 inline bool
850 {
851  if (!this->test()) return false;
852  this->increment();
853  return this->test();
854 }
855 
856 
857 template<typename PointIndexGridType>
858 inline size_t
860 {
861  size_t count = 0;
862  typename RangeDeque::const_iterator it =
863  mRangeList.begin(), end = mRangeList.end();
864 
865  for ( ; it != end; ++it) {
866  count += it->second - it->first;
867  }
868 
869  return count + mIndexArraySize;
870 }
871 
872 
873 template<typename PointIndexGridType>
874 inline void
876 {
877  mRange.first = static_cast<IndexType*>(NULL);
878  mRange.second = static_cast<IndexType*>(NULL);
879  mRangeList.clear();
880  mIter = mRangeList.end();
881  mIndexArray.reset();
882  mIndexArraySize = 0;
883 }
884 
885 
886 template<typename PointIndexGridType>
887 inline void
889 {
890  typedef typename PointIndexGridType::TreeType TreeType;
891  typedef typename TreeType::LeafNodeType LeafNodeType;
892 
893  this->clear();
894 
895  if (mAccessorListSize > 0) {
896  const size_t levelIdx = std::min(mAccessorListSize - 1, level);
897 
898  const TreeType& tree = mAtlas->pointIndexGrid(levelIdx).tree();
899 
900 
901  std::vector<const LeafNodeType*> nodes;
902  tree.getNodes(nodes);
903 
904  for (size_t n = 0, N = nodes.size(); n < N; ++n) {
905 
906  const LeafNodeType& node = *nodes[n];
907  const size_t numIndices = node.indices().size();
908 
909  if (numIndices > 0) {
910  const IndexType* begin = &node.indices().front();
911  mRangeList.push_back(Range(begin, (begin + numIndices)));
912  }
913  }
914  }
915 
916  this->reset();
917 }
918 
919 
920 template<typename PointIndexGridType>
921 template<typename ParticleArrayType>
922 inline void
924  const Vec3d& center, double radius, const ParticleArrayType& particles)
925 {
926  typedef typename ParticleArrayType::PosType PosType;
927  typedef typename ParticleArrayType::ScalarType ScalarType;
928 
930 
931  this->clear();
932 
933  std::deque<IndexType> filteredIndices;
934  std::vector<CoordBBox> searchRegions;
935 
936  const double iRadius = radius * double(1.0 / std::sqrt(3.0));
937 
938  const Vec3d ibMin(center[0] - iRadius, center[1] - iRadius, center[2] - iRadius);
939  const Vec3d ibMax(center[0] + iRadius, center[1] + iRadius, center[2] + iRadius);
940 
941  const Vec3d bMin(center[0] - radius, center[1] - radius, center[2] - radius);
942  const Vec3d bMax(center[0] + radius, center[1] + radius, center[2] + radius);
943 
944  const PosType pos = PosType(center);
945  const ScalarType dist = ScalarType(radius);
946 
947  for (size_t n = 0, N = mAccessorListSize; n < N; ++n) {
948 
949  const double maxRadius = mAtlas->maxRadius(n);
950  const bool uniformRadius = math::isApproxEqual(mAtlas->minRadius(n), maxRadius);
951 
952  const openvdb::math::Transform& xform = mAtlas->pointIndexGrid(n).transform();
953 
954  ConstAccessor& acc = *mAccessorList[n];
955 
956  openvdb::CoordBBox inscribedRegion(
957  xform.worldToIndexCellCentered(ibMin),
958  xform.worldToIndexCellCentered(ibMax));
959 
960  inscribedRegion.expand(-1); // erode by one voxel
961 
962  // collect indices that don't need to be tested
963  point_index_grid_internal::pointIndexSearch(mRangeList, acc, inscribedRegion);
964 
965  searchRegions.clear();
966 
967  const openvdb::CoordBBox region(
968  xform.worldToIndexCellCentered(bMin - maxRadius),
969  xform.worldToIndexCellCentered(bMax + maxRadius));
970 
971  inscribedRegion.expand(1);
972  point_index_grid_internal::constructExclusiveRegions(searchRegions, region, inscribedRegion);
973 
975  FilterType filter(mRangeList, filteredIndices, pos, dist, particles, uniformRadius);
976 
977  for (size_t i = 0, I = searchRegions.size(); i < I; ++i) {
978  point_index_grid_internal::filteredPointIndexSearch(filter, acc, searchRegions[i]);
979  }
980  }
981 
982  point_index_grid_internal::dequeToArray(filteredIndices, mIndexArray, mIndexArraySize);
983 
984  this->reset();
985 }
986 
987 
988 template<typename PointIndexGridType>
989 template<typename ParticleArrayType>
990 inline void
992  const BBoxd& bbox, const ParticleArrayType& particles)
993 {
994  typedef typename ParticleArrayType::PosType PosType;
995  typedef typename ParticleArrayType::ScalarType ScalarType;
996 
998 
999  this->clear();
1000 
1001  std::deque<IndexType> filteredIndices;
1002  std::vector<CoordBBox> searchRegions;
1003 
1004  for (size_t n = 0, N = mAccessorListSize; n < N; ++n) {
1005 
1006  const double maxRadius = mAtlas->maxRadius(n);
1007  const bool uniformRadius = math::isApproxEqual(mAtlas->minRadius(n), maxRadius);
1008  const openvdb::math::Transform& xform = mAtlas->pointIndexGrid(n).transform();
1009 
1010  ConstAccessor& acc = *mAccessorList[n];
1011 
1012  openvdb::CoordBBox inscribedRegion(
1013  xform.worldToIndexCellCentered(bbox.min()),
1014  xform.worldToIndexCellCentered(bbox.max()));
1015 
1016  inscribedRegion.expand(-1); // erode by one voxel
1017 
1018  // collect indices that don't need to be tested
1019  point_index_grid_internal::pointIndexSearch(mRangeList, acc, inscribedRegion);
1020 
1021  searchRegions.clear();
1022 
1023  const openvdb::CoordBBox region(
1024  xform.worldToIndexCellCentered(bbox.min() - maxRadius),
1025  xform.worldToIndexCellCentered(bbox.max() + maxRadius));
1026 
1027  inscribedRegion.expand(1);
1028  point_index_grid_internal::constructExclusiveRegions(searchRegions, region, inscribedRegion);
1029 
1031  FilterType filter(mRangeList, filteredIndices, bbox, particles, uniformRadius);
1032 
1033  for (size_t i = 0, I = searchRegions.size(); i < I; ++i) {
1034  point_index_grid_internal::filteredPointIndexSearch(filter, acc, searchRegions[i]);
1035  }
1036  }
1037 
1038  point_index_grid_internal::dequeToArray(filteredIndices, mIndexArray, mIndexArraySize);
1039 
1040  this->reset();
1041 }
1042 
1043 
1044 } // namespace tools
1045 } // namespace OPENVDB_VERSION_NAME
1046 } // namespace openvdb
1047 
1048 #endif // OPENVDB_TOOLS_PARTICLE_ATLAS_HAS_BEEN_INCLUDED
1049 
1050 // Copyright (c) 2012-2016 DreamWorks Animation LLC
1051 // All rights reserved. This software is distributed under the
1052 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
ComputeExtremas(ComputeExtremas &rhs, tbb::split)
Definition: ParticleAtlas.h:279
std::pair< const IndexT *, const IndexT * > Range
Definition: ParticleAtlas.h:491
PosType::value_type ScalarType
Definition: ParticleAtlas.h:270
bool test() const
Return true if this iterator is not yet exhausted.
Definition: ParticleAtlas.h:214
void reset()
Reset the iterator to point to the first item.
Definition: ParticleAtlas.h:815
ScalarType minRadius() const
Definition: ParticleAtlas.h:340
boost::shared_ptr< const SplittableParticleArray > ConstPtr
Definition: ParticleAtlas.h:314
ParticleArray::PosType PosType
Definition: ParticleAtlas.h:317
bool operator==(const Iterator &p) const
Return true if both iterators point to the same element.
Definition: ParticleAtlas.h:232
Space-partitioning acceleration structure for points. Partitions the points into voxels to accelerate...
SplittableParticleArray(const ParticleArrayT &particles, double minR, double maxR)
Definition: ParticleAtlas.h:326
ComputeExtremas(const ParticleArrayT &particles)
Definition: ParticleAtlas.h:272
void filterLeafNode(const LeafNodeType &leaf)
Definition: ParticleAtlas.h:600
boost::shared_ptr< ParticleAtlas > Ptr
Definition: ParticleAtlas.h:108
ParticleArrayType const *const mParticles
Definition: ParticleAtlas.h:480
std::deque< IndexT > IndexDeque
Definition: ParticleAtlas.h:493
std::deque< Range > RangeDeque
Definition: ParticleAtlas.h:580
const ParticleArrayT & particleArray() const
Definition: ParticleAtlas.h:333
RadialRangeFilter(RangeDeque &ranges, IndexDeque &indices, const PosType &xyz, ScalarType radius, const ParticleArrayType &particles, bool hasUniformRadius=false)
Definition: ParticleAtlas.h:495
void increment()
Advance iterator to next item.
Definition: ParticleAtlas.h:832
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:47
ParticleArrayType::PosType PosType
Definition: ParticleAtlas.h:576
const boost::disable_if_c< VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:132
void filterVoxel(const Coord &, const IndexT *begin, const IndexT *end)
Definition: ParticleAtlas.h:522
RemapIndices(const ParticleArrayType &particles, std::vector< PointIndexLeafNodeType *> nodes)
Definition: ParticleAtlas.h:455
const IndexType & operator*() const
Return a const reference to the item to which this iterator is pointing.
Definition: ParticleAtlas.h:210
void updateFromLevel(size_t level)
Clear the iterator and update it with all particles that reside at the given resolution level...
Definition: ParticleAtlas.h:888
bool isApproxEqual(const Type &a, const Type &b)
Return true if a is equal to b to within the default floating-point comparison tolerance.
Definition: Math.h:370
PosType::value_type ScalarType
Definition: ParticleAtlas.h:489
Selectively extract and filter point data using a custom filter operator.
ParticleArrayT ParticleArray
Definition: ParticleAtlas.h:315
const PointIndexGridType & pointIndexGrid(size_t n) const
Returns the PointIndexGrid that represents the given level n.
Definition: ParticleAtlas.h:149
boost::shared_ptr< const ParticleAtlas > ConstPtr
Definition: ParticleAtlas.h:109
#define OPENVDB_VERSION_NAME
Definition: version.h:43
PosType::value_type ScalarType
Definition: ParticleAtlas.h:577
Vec3< double > Vec3d
Definition: Vec3.h:651
void expand(ValueType padding)
Pad this bounding box with the specified padding.
Definition: Coord.h:401
bool next()
Advance iterator to next item.
Definition: ParticleAtlas.h:849
size_t size() const
Return the number of point indices in the iterator range.
Definition: ParticleAtlas.h:859
std::deque< Range > RangeDeque
Definition: ParticleAtlas.h:492
Provides accelerated range and nearest-neighbor searches for particles that are partitioned using the...
Definition: ParticleAtlas.h:173
boost::scoped_ptr< ConstAccessor > ConstAccessorPtr
Definition: ParticleAtlas.h:177
ParticleArrayT::PosType PosType
Definition: ParticleAtlas.h:269
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: ParticleAtlas.h:461
Definition: Exceptions.h:39
ScalarType minRadius
Definition: ParticleAtlas.h:306
boost::shared_ptr< SplittableParticleArray > Ptr
Definition: ParticleAtlas.h:313
ParticleArrayT const *const particleArray
Definition: ParticleAtlas.h:305
PointIndexGridType::Ptr PointIndexGridPtr
Definition: ParticleAtlas.h:111
ScalarType maxRadius
Definition: ParticleAtlas.h:306
const Vec3T & min() const
Return a const reference to the minimum point of the BBox.
Definition: BBox.h:81
bool operator!=(const Iterator &p) const
Definition: ParticleAtlas.h:233
void filteredPointIndexSearch(RangeFilterType &filter, ConstAccessor &acc, const CoordBBox &bbox)
Definition: PointIndexGrid.h:816
ParticleArrayType::PosType PosType
Definition: ParticleAtlas.h:488
void join(const ComputeExtremas &rhs)
Definition: ParticleAtlas.h:300
const Vec3T & max() const
Return a const reference to the maximum point of the BBox.
Definition: BBox.h:84
size_t getGlobalIndex(size_t n) const
Definition: ParticleAtlas.h:343
Ptr split(ScalarType maxRadiusLimit)
Definition: ParticleAtlas.h:347
void getPos(size_t n, PosType &xyz) const
Definition: ParticleAtlas.h:337
tree::ValueAccessor< const TreeType > ConstAccessor
Definition: ParticleAtlas.h:176
std::deque< IndexT > IndexDeque
Definition: ParticleAtlas.h:581
bool operator==(const Vec3< T0 > &v0, const Vec3< T1 > &v1)
Equality operator, does exact floating point comparisons.
Definition: Vec3.h:450
void filterLeafNode(const LeafNodeType &leaf)
Definition: ParticleAtlas.h:513
size_t levels() const
Returns the number of resolution levels.
Definition: ParticleAtlas.h:137
ParticleAtlas< PointIndexGrid > ParticleIndexAtlas
Definition: ParticleAtlas.h:161
Axis-aligned bounding box of signed integer coordinates.
Definition: Coord.h:259
double minRadius(size_t n) const
Returns minimum particle radius for level n.
Definition: ParticleAtlas.h:142
size_t levels() const
Returns the total number of resolution levels.
Definition: ParticleAtlas.h:200
Integer wrapper, required to distinguish PointIndexGrid and PointDataGrid from Int32Grid and Int64Gri...
Definition: Types.h:119
void filterVoxel(const Coord &, const IndexT *begin, const IndexT *end)
Definition: ParticleAtlas.h:609
Definition: ParticleAtlas.h:106
void getRadius(size_t n, ScalarType &radius) const
Definition: ParticleAtlas.h:338
PointIndexGridType::ValueType IndexType
Definition: ParticleAtlas.h:112
ParticleAtlas()
Definition: ParticleAtlas.h:118
PointIndexLeafNodeType *const *const mNodes
Definition: ParticleAtlas.h:481
BBoxFilter(RangeDeque &ranges, IndexDeque &indices, const BBoxd &bbox, const ParticleArrayType &particles, bool hasUniformRadius=false)
Definition: ParticleAtlas.h:583
void dequeToArray(const std::deque< T > &d, boost::scoped_array< T > &a, size_t &size)
Definition: PointIndexGrid.h:548
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:71
void operator()(const tbb::blocked_range< size_t > &range)
Definition: ParticleAtlas.h:286
void pointIndexSearch(RangeDeque &rangeList, ConstAccessor &acc, const CoordBBox &bbox)
Definition: PointIndexGrid.h:892
PointIndexGridType & pointIndexGrid(size_t n)
Returns the PointIndexGrid that represents the given level n.
Definition: ParticleAtlas.h:147
void worldSpaceSearchAndUpdate(const Vec3d &center, double radius, const ParticleArrayType &particles)
Clear the iterator and update it with the result of the given world-space radial query.
Definition: ParticleAtlas.h:923
PointIndexGridType::TreeType TreeType
Definition: ParticleAtlas.h:175
void constructExclusiveRegions(std::vector< CoordBBox > &regions, const CoordBBox &bbox, const CoordBBox &ibox)
Definition: PointIndexGrid.h:559
std::pair< const IndexT *, const IndexT * > Range
Definition: ParticleAtlas.h:579
double maxRadius(size_t n) const
Returns maximum particle radius for level n.
Definition: ParticleAtlas.h:144
SplittableParticleArray(const ParticleArrayT &particles)
Definition: ParticleAtlas.h:320
ScalarType maxRadius() const
Definition: ParticleAtlas.h:341
PosType::value_type ScalarType
Definition: ParticleAtlas.h:318
bool empty() const
true if the container size is 0, false otherwise.
Definition: ParticleAtlas.h:139
void operator++()
Advance iterator to next item.
Definition: ParticleAtlas.h:222
const boost::disable_if_c< VecTraits< T >::IsVec, T >::type & min(const T &a, const T &b)
Definition: Composite.h:128