OpenVDB  6.2.0
PointPartitioner.h
Go to the documentation of this file.
1 //
3 // Copyright (c) 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 
44 
45 #ifndef OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED
46 #define OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED
47 
48 #include <openvdb/Types.h>
49 #include <openvdb/math/Transform.h>
50 
51 #include <boost/integer.hpp> // boost::int_t<N>::least
52 #include <boost/scoped_array.hpp>
53 
54 #include <tbb/blocked_range.h>
55 #include <tbb/parallel_for.h>
56 #include <tbb/task_scheduler_init.h>
57 
58 #include <algorithm>
59 #include <cmath> // for std::isfinite()
60 #include <deque>
61 #include <map>
62 #include <set>
63 #include <utility> // std::pair
64 #include <vector>
65 
66 
67 namespace openvdb {
69 namespace OPENVDB_VERSION_NAME {
70 namespace tools {
71 
72 
74 
75 
105 template<typename PointIndexType = uint32_t, Index BucketLog2Dim = 3>
107 {
108 public:
109  enum { LOG2DIM = BucketLog2Dim };
110 
113 
114  using IndexType = PointIndexType;
115  using VoxelOffsetType = typename boost::int_t<1 + (3 * BucketLog2Dim)>::least;
116  using VoxelOffsetArray = boost::scoped_array<VoxelOffsetType>;
117 
118  class IndexIterator;
119 
121 
123 
133  template<typename PointArray>
134  void construct(const PointArray& points, const math::Transform& xform,
135  bool voxelOrder = false, bool recordVoxelOffsets = false,
136  bool cellCenteredTransform = true);
137 
138 
148  template<typename PointArray>
149  static Ptr create(const PointArray& points, const math::Transform& xform,
150  bool voxelOrder = false, bool recordVoxelOffsets = false,
151  bool cellCenteredTransform = true);
152 
153 
155  size_t size() const { return mPageCount; }
156 
158  bool empty() const { return mPageCount == 0; }
159 
161  void clear();
162 
164  void swap(PointPartitioner&);
165 
167  IndexIterator indices(size_t n) const;
168 
170  CoordBBox getBBox(size_t n) const {
171  return CoordBBox::createCube(mPageCoordinates[n], (1u << BucketLog2Dim));
172  }
173 
175  const Coord& origin(size_t n) const { return mPageCoordinates[n]; }
176 
179  const VoxelOffsetArray& voxelOffsets() const { return mVoxelOffsets; }
180 
184  bool usingCellCenteredTransform() const { return mUsingCellCenteredTransform; }
185 
186 private:
187  // Disallow copying
189  PointPartitioner& operator=(const PointPartitioner&);
190 
191  boost::scoped_array<IndexType> mPointIndices;
192  VoxelOffsetArray mVoxelOffsets;
193 
194  boost::scoped_array<IndexType> mPageOffsets;
195  boost::scoped_array<Coord> mPageCoordinates;
196  IndexType mPageCount;
197  bool mUsingCellCenteredTransform;
198 }; // class PointPartitioner
199 
200 
202 
203 
204 template<typename PointIndexType, Index BucketLog2Dim>
205 class PointPartitioner<PointIndexType, BucketLog2Dim>::IndexIterator
206 {
207 public:
208  using IndexType = PointIndexType;
209 
210  IndexIterator(IndexType* begin = nullptr, IndexType* end = nullptr)
211  : mBegin(begin), mEnd(end), mItem(begin) {}
212 
214  void reset() { mItem = mBegin; }
215 
217  size_t size() const { return mEnd - mBegin; }
218 
220  IndexType& operator*() { assert(mItem != nullptr); return *mItem; }
221  const IndexType& operator*() const { assert(mItem != nullptr); return *mItem; }
222 
224  operator bool() const { return mItem < mEnd; }
225  bool test() const { return mItem < mEnd; }
226 
228  IndexIterator& operator++() { assert(this->test()); ++mItem; return *this; }
229 
231  bool next() { this->operator++(); return this->test(); }
232  bool increment() { this->next(); return this->test(); }
233 
235  bool operator==(const IndexIterator& other) const { return mItem == other.mItem; }
236  bool operator!=(const IndexIterator& other) const { return !this->operator==(other); }
237 
238 private:
239  IndexType * const mBegin, * const mEnd;
240  IndexType * mItem;
241 }; // class PointPartitioner::IndexIterator
242 
243 
246 
247 // Implementation details
248 
249 
250 namespace point_partitioner_internal {
251 
252 
253 template<typename PointIndexType>
255 {
256  ComputePointOrderOp(PointIndexType* pointOrder,
257  const PointIndexType* bucketCounters, const PointIndexType* bucketOffsets)
258  : mPointOrder(pointOrder)
259  , mBucketCounters(bucketCounters)
260  , mBucketOffsets(bucketOffsets)
261  {
262  }
263 
264  void operator()(const tbb::blocked_range<size_t>& range) const {
265  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
266  mPointOrder[n] += mBucketCounters[mBucketOffsets[n]];
267  }
268  }
269 
270  PointIndexType * const mPointOrder;
271  PointIndexType const * const mBucketCounters;
272  PointIndexType const * const mBucketOffsets;
273 }; // struct ComputePointOrderOp
274 
275 
276 template<typename PointIndexType>
278 {
279  CreateOrderedPointIndexArrayOp(PointIndexType* orderedIndexArray,
280  const PointIndexType* pointOrder, const PointIndexType* indices)
281  : mOrderedIndexArray(orderedIndexArray)
282  , mPointOrder(pointOrder)
283  , mIndices(indices)
284  {
285  }
286 
287  void operator()(const tbb::blocked_range<size_t>& range) const {
288  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
289  mOrderedIndexArray[mPointOrder[n]] = mIndices[n];
290  }
291  }
292 
293  PointIndexType * const mOrderedIndexArray;
294  PointIndexType const * const mPointOrder;
295  PointIndexType const * const mIndices;
296 }; // struct CreateOrderedPointIndexArrayOp
297 
298 
299 template<typename PointIndexType, Index BucketLog2Dim>
301 {
302  using VoxelOffsetType = typename boost::int_t<1 + (3 * BucketLog2Dim)>::least;
303  using VoxelOffsetArray = boost::scoped_array<VoxelOffsetType>;
304  using IndexArray = boost::scoped_array<PointIndexType>;
305 
306  VoxelOrderOp(IndexArray& indices, const IndexArray& pages,const VoxelOffsetArray& offsets)
307  : mIndices(indices.get())
308  , mPages(pages.get())
309  , mVoxelOffsets(offsets.get())
310  {
311  }
312 
313  void operator()(const tbb::blocked_range<size_t>& range) const {
314 
315  PointIndexType pointCount = 0;
316  for (size_t n(range.begin()), N(range.end()); n != N; ++n) {
317  pointCount = std::max(pointCount, (mPages[n + 1] - mPages[n]));
318  }
319 
320  const PointIndexType voxelCount = 1 << (3 * BucketLog2Dim);
321 
322  // allocate histogram buffers
323  boost::scoped_array<VoxelOffsetType> offsets(new VoxelOffsetType[pointCount]);
324  boost::scoped_array<PointIndexType> sortedIndices(new PointIndexType[pointCount]);
325  boost::scoped_array<PointIndexType> histogram(new PointIndexType[voxelCount]);
326 
327  for (size_t n(range.begin()), N(range.end()); n != N; ++n) {
328 
329  PointIndexType * const indices = mIndices + mPages[n];
330  pointCount = mPages[n + 1] - mPages[n];
331 
332  // local copy of voxel offsets.
333  for (PointIndexType i = 0; i < pointCount; ++i) {
334  offsets[i] = mVoxelOffsets[ indices[i] ];
335  }
336 
337  // reset histogram
338  memset(&histogram[0], 0, voxelCount * sizeof(PointIndexType));
339 
340  // compute histogram
341  for (PointIndexType i = 0; i < pointCount; ++i) {
342  ++histogram[ offsets[i] ];
343  }
344 
345  PointIndexType count = 0, startOffset;
346  for (int i = 0; i < int(voxelCount); ++i) {
347  if (histogram[i] > 0) {
348  startOffset = count;
349  count += histogram[i];
350  histogram[i] = startOffset;
351  }
352  }
353 
354  // sort indices based on voxel offset
355  for (PointIndexType i = 0; i < pointCount; ++i) {
356  sortedIndices[ histogram[ offsets[i] ]++ ] = indices[i];
357  }
358 
359  memcpy(&indices[0], &sortedIndices[0], sizeof(PointIndexType) * pointCount);
360  }
361  }
362 
363  PointIndexType * const mIndices;
364  PointIndexType const * const mPages;
366 }; // struct VoxelOrderOp
367 
368 
369 template<typename PointArray, typename PointIndexType>
371 {
372  using IndexArray = boost::scoped_array<PointIndexType>;
373  using CoordArray = boost::scoped_array<Coord>;
374 
376  const IndexArray& indices, const IndexArray& pages,
377  const PointArray& points, const math::Transform& m, int log2dim, bool cellCenteredTransform)
378  : mCoordinates(coordinates.get())
379  , mIndices(indices.get())
380  , mPages(pages.get())
381  , mPoints(&points)
382  , mXForm(m)
383  , mLog2Dim(log2dim)
384  , mCellCenteredTransform(cellCenteredTransform)
385  {
386  }
387 
388  void operator()(const tbb::blocked_range<size_t>& range) const {
389 
390  using PosType = typename PointArray::PosType;
391 
392  const bool cellCentered = mCellCenteredTransform;
393  const int mask = ~((1 << mLog2Dim) - 1);
394  Coord ijk;
395  PosType pos;
396 
397  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
398 
399  mPoints->getPos(mIndices[mPages[n]], pos);
400 
401  if (std::isfinite(pos[0]) && std::isfinite(pos[1]) && std::isfinite(pos[2])) {
402  ijk = cellCentered ? mXForm.worldToIndexCellCentered(pos) :
403  mXForm.worldToIndexNodeCentered(pos);
404 
405  ijk[0] &= mask;
406  ijk[1] &= mask;
407  ijk[2] &= mask;
408 
409  mCoordinates[n] = ijk;
410  }
411  }
412  }
413 
414  Coord * const mCoordinates;
415  PointIndexType const * const mIndices;
416  PointIndexType const * const mPages;
417  PointArray const * const mPoints;
419  int const mLog2Dim;
421 }; // struct LeafNodeOriginOp
422 
423 
425 
426 
427 template<typename T>
428 struct Array
429 {
431 
432  Array(size_t size) : mSize(size), mData(new T[size]) { }
433 
434  size_t size() const { return mSize; }
435 
436  T* data() { return mData.get(); }
437  const T* data() const { return mData.get(); }
438 
439  void clear() { mSize = 0; mData.reset(); }
440 
441 private:
442  size_t mSize;
443  boost::scoped_array<T> mData;
444 }; // struct Array
445 
446 
447 template<typename PointIndexType>
449 {
451  using SegmentPtr = typename Segment::Ptr;
452 
453  MoveSegmentDataOp(std::vector<PointIndexType*>& indexLists, SegmentPtr* segments)
454  : mIndexLists(&indexLists[0]), mSegments(segments)
455  {
456  }
457 
458  void operator()(const tbb::blocked_range<size_t>& range) const {
459  for (size_t n(range.begin()), N(range.end()); n != N; ++n) {
460  PointIndexType* indices = mIndexLists[n];
461  SegmentPtr& segment = mSegments[n];
462 
463  tbb::parallel_for(tbb::blocked_range<size_t>(0, segment->size()),
464  CopyData(indices, segment->data()));
465 
466  segment.reset(); // clear data
467  }
468  }
469 
470 private:
471 
472  struct CopyData
473  {
474  CopyData(PointIndexType* lhs, const PointIndexType* rhs) : mLhs(lhs), mRhs(rhs) { }
475 
476  void operator()(const tbb::blocked_range<size_t>& range) const {
477  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
478  mLhs[n] = mRhs[n];
479  }
480  }
481 
482  PointIndexType * const mLhs;
483  PointIndexType const * const mRhs;
484  };
485 
486  PointIndexType * const * const mIndexLists;
487  SegmentPtr * const mSegments;
488 }; // struct MoveSegmentDataOp
489 
490 
491 template<typename PointIndexType>
493 {
495  using SegmentPtr = typename Segment::Ptr;
496 
497  using IndexPair = std::pair<PointIndexType, PointIndexType>;
498  using IndexPairList = std::deque<IndexPair>;
500  using IndexPairListMap = std::map<Coord, IndexPairListPtr>;
502 
504  SegmentPtr* indexSegments,
505  SegmentPtr* offsetSegments,
506  Coord* coords,
507  size_t numSegments)
508  : mBins(bins)
509  , mIndexSegments(indexSegments)
510  , mOffsetSegments(offsetSegments)
511  , mCoords(coords)
512  , mNumSegments(numSegments)
513  {
514  }
515 
516  void operator()(const tbb::blocked_range<size_t>& range) const {
517 
518  std::vector<IndexPairListPtr*> data;
519  std::vector<PointIndexType> arrayOffsets;
520 
521  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
522 
523  const Coord& ijk = mCoords[n];
524  size_t numIndices = 0;
525 
526  data.clear();
527 
528  for (size_t i = 0, I = mNumSegments; i < I; ++i) {
529 
530  IndexPairListMap& idxMap = *mBins[i];
531  typename IndexPairListMap::iterator iter = idxMap.find(ijk);
532 
533  if (iter != idxMap.end() && iter->second) {
534  IndexPairListPtr& idxListPtr = iter->second;
535 
536  data.push_back(&idxListPtr);
537  numIndices += idxListPtr->size();
538  }
539  }
540 
541  if (data.empty() || numIndices == 0) continue;
542 
543  SegmentPtr& indexSegment = mIndexSegments[n];
544  SegmentPtr& offsetSegment = mOffsetSegments[n];
545 
546  indexSegment.reset(new Segment(numIndices));
547  offsetSegment.reset(new Segment(numIndices));
548 
549  arrayOffsets.clear();
550  arrayOffsets.reserve(data.size());
551 
552  for (size_t i = 0, count = 0, I = data.size(); i < I; ++i) {
553  arrayOffsets.push_back(PointIndexType(count));
554  count += (*data[i])->size();
555  }
556 
557  tbb::parallel_for(tbb::blocked_range<size_t>(0, data.size()),
558  CopyData(&data[0], &arrayOffsets[0], indexSegment->data(), offsetSegment->data()));
559  }
560  }
561 
562 private:
563 
564  struct CopyData
565  {
566  CopyData(IndexPairListPtr** indexLists,
567  const PointIndexType* arrayOffsets,
568  PointIndexType* indices,
569  PointIndexType* offsets)
570  : mIndexLists(indexLists)
571  , mArrayOffsets(arrayOffsets)
572  , mIndices(indices)
573  , mOffsets(offsets)
574  {
575  }
576 
577  void operator()(const tbb::blocked_range<size_t>& range) const {
578 
579  using CIter = typename IndexPairList::const_iterator;
580 
581  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
582 
583  const PointIndexType arrayOffset = mArrayOffsets[n];
584  PointIndexType* indexPtr = &mIndices[arrayOffset];
585  PointIndexType* offsetPtr = &mOffsets[arrayOffset];
586 
587  IndexPairListPtr& list = *mIndexLists[n];
588 
589  for (CIter it = list->begin(), end = list->end(); it != end; ++it) {
590  const IndexPair& data = *it;
591  *indexPtr++ = data.first;
592  *offsetPtr++ = data.second;
593  }
594 
595  list.reset(); // clear data
596  }
597  }
598 
599  IndexPairListPtr * const * const mIndexLists;
600  PointIndexType const * const mArrayOffsets;
601  PointIndexType * const mIndices;
602  PointIndexType * const mOffsets;
603  }; // struct CopyData
604 
605  IndexPairListMapPtr * const mBins;
606  SegmentPtr * const mIndexSegments;
607  SegmentPtr * const mOffsetSegments;
608  Coord const * const mCoords;
609  size_t const mNumSegments;
610 }; // struct MergeBinsOp
611 
612 
613 template<typename PointArray, typename PointIndexType, typename VoxelOffsetType>
615 {
616  using PosType = typename PointArray::PosType;
617  using IndexPair = std::pair<PointIndexType, PointIndexType>;
618  using IndexPairList = std::deque<IndexPair>;
620  using IndexPairListMap = std::map<Coord, IndexPairListPtr>;
622 
624  const PointArray& points,
625  VoxelOffsetType* voxelOffsets,
626  const math::Transform& m,
627  Index binLog2Dim,
628  Index bucketLog2Dim,
629  size_t numSegments,
630  bool cellCenteredTransform)
631  : mData(data)
632  , mPoints(&points)
633  , mVoxelOffsets(voxelOffsets)
634  , mXForm(m)
635  , mBinLog2Dim(binLog2Dim)
636  , mBucketLog2Dim(bucketLog2Dim)
637  , mNumSegments(numSegments)
638  , mCellCenteredTransform(cellCenteredTransform)
639  {
640  }
641 
642  void operator()(const tbb::blocked_range<size_t>& range) const {
643 
644  const Index log2dim = mBucketLog2Dim;
645  const Index log2dim2 = 2 * log2dim;
646  const Index bucketMask = (1u << log2dim) - 1u;
647 
648  const Index binLog2dim = mBinLog2Dim;
649  const Index binLog2dim2 = 2 * binLog2dim;
650 
651  const Index binMask = (1u << (log2dim + binLog2dim)) - 1u;
652  const Index invBinMask = ~binMask;
653 
654  IndexPairList * idxList = nullptr;
655  Coord ijk(0, 0, 0), loc(0, 0, 0), binCoord(0, 0, 0), lastBinCoord(1, 2, 3);
656  PosType pos;
657 
658  PointIndexType bucketOffset = 0;
659  VoxelOffsetType voxelOffset = 0;
660 
661  const bool cellCentered = mCellCenteredTransform;
662 
663  const size_t numPoints = mPoints->size();
664  const size_t segmentSize = numPoints / mNumSegments;
665 
666  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
667 
668  IndexPairListMapPtr& dataPtr = mData[n];
669  if (!dataPtr) dataPtr.reset(new IndexPairListMap());
670  IndexPairListMap& idxMap = *dataPtr;
671 
672  const bool isLastSegment = (n + 1) >= mNumSegments;
673 
674  const size_t start = n * segmentSize;
675  const size_t end = isLastSegment ? numPoints : (start + segmentSize);
676 
677  for (size_t i = start; i != end; ++i) {
678 
679  mPoints->getPos(i, pos);
680 
681  if (std::isfinite(pos[0]) && std::isfinite(pos[1]) && std::isfinite(pos[2])) {
682  ijk = cellCentered ? mXForm.worldToIndexCellCentered(pos) :
683  mXForm.worldToIndexNodeCentered(pos);
684 
685  if (mVoxelOffsets) {
686  loc[0] = ijk[0] & bucketMask;
687  loc[1] = ijk[1] & bucketMask;
688  loc[2] = ijk[2] & bucketMask;
689  voxelOffset = VoxelOffsetType(
690  (loc[0] << log2dim2) + (loc[1] << log2dim) + loc[2]);
691  }
692 
693  binCoord[0] = ijk[0] & invBinMask;
694  binCoord[1] = ijk[1] & invBinMask;
695  binCoord[2] = ijk[2] & invBinMask;
696 
697  ijk[0] &= binMask;
698  ijk[1] &= binMask;
699  ijk[2] &= binMask;
700 
701  ijk[0] >>= log2dim;
702  ijk[1] >>= log2dim;
703  ijk[2] >>= log2dim;
704 
705  bucketOffset = PointIndexType(
706  (ijk[0] << binLog2dim2) + (ijk[1] << binLog2dim) + ijk[2]);
707 
708  if (lastBinCoord != binCoord) {
709  lastBinCoord = binCoord;
710  IndexPairListPtr& idxListPtr = idxMap[lastBinCoord];
711  if (!idxListPtr) idxListPtr.reset(new IndexPairList());
712  idxList = idxListPtr.get();
713  }
714 
715  idxList->push_back(IndexPair(PointIndexType(i), bucketOffset));
716  if (mVoxelOffsets) mVoxelOffsets[i] = voxelOffset;
717  }
718  }
719  }
720  }
721 
723  PointArray const * const mPoints;
724  VoxelOffsetType * const mVoxelOffsets;
728  size_t const mNumSegments;
730 }; // struct BinPointIndicesOp
731 
732 
733 template<typename PointIndexType>
735 {
736  using IndexArray = boost::scoped_array<PointIndexType>;
738 
739  OrderSegmentsOp(SegmentPtr* indexSegments, SegmentPtr* offestSegments,
740  IndexArray* pageOffsetArrays, Index binVolume)
741  : mIndexSegments(indexSegments)
742  , mOffsetSegments(offestSegments)
743  , mPageOffsetArrays(pageOffsetArrays)
744  , mBinVolume(binVolume)
745  {
746  }
747 
748  void operator()(const tbb::blocked_range<size_t>& range) const {
749 
750  const size_t bucketCountersSize = size_t(mBinVolume);
751  IndexArray bucketCounters(new PointIndexType[bucketCountersSize]);
752 
753  size_t maxSegmentSize = 0;
754  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
755  maxSegmentSize = std::max(maxSegmentSize, mIndexSegments[n]->size());
756  }
757 
758  IndexArray bucketIndices(new PointIndexType[maxSegmentSize]);
759 
760 
761  for (size_t n = range.begin(), N = range.end(); n != N; ++n) {
762 
763  memset(bucketCounters.get(), 0, sizeof(PointIndexType) * bucketCountersSize);
764 
765  const size_t segmentSize = mOffsetSegments[n]->size();
766  PointIndexType* offsets = mOffsetSegments[n]->data();
767 
768  // Count the number of points per bucket and assign a local bucket index
769  // to each point.
770  for (size_t i = 0; i < segmentSize; ++i) {
771  bucketIndices[i] = bucketCounters[offsets[i]]++;
772  }
773 
774  PointIndexType nonemptyBucketCount = 0;
775  for (size_t i = 0; i < bucketCountersSize; ++i) {
776  nonemptyBucketCount += static_cast<PointIndexType>(bucketCounters[i] != 0);
777  }
778 
779 
780  IndexArray& pageOffsets = mPageOffsetArrays[n];
781  pageOffsets.reset(new PointIndexType[nonemptyBucketCount + 1]);
782  pageOffsets[0] = nonemptyBucketCount + 1; // stores array size in first element
783 
784  // Compute bucket counter prefix sum
785  PointIndexType count = 0, idx = 1;
786  for (size_t i = 0; i < bucketCountersSize; ++i) {
787  if (bucketCounters[i] != 0) {
788  pageOffsets[idx] = bucketCounters[i];
789  bucketCounters[i] = count;
790  count += pageOffsets[idx];
791  ++idx;
792  }
793  }
794 
795  PointIndexType* indices = mIndexSegments[n]->data();
796  const tbb::blocked_range<size_t> segmentRange(0, segmentSize);
797 
798  // Compute final point order by incrementing the local bucket point index
799  // with the prefix sum offset.
800  tbb::parallel_for(segmentRange, ComputePointOrderOp<PointIndexType>(
801  bucketIndices.get(), bucketCounters.get(), offsets));
802 
803  tbb::parallel_for(segmentRange, CreateOrderedPointIndexArrayOp<PointIndexType>(
804  offsets, bucketIndices.get(), indices));
805 
806  mIndexSegments[n]->clear(); // clear data
807  }
808  }
809 
814 }; // struct OrderSegmentsOp
815 
816 
818 
819 
821 template<typename PointIndexType, typename VoxelOffsetType, typename PointArray>
822 inline void binAndSegment(
823  const PointArray& points,
824  const math::Transform& xform,
825  boost::scoped_array<typename Array<PointIndexType>::Ptr>& indexSegments,
826  boost::scoped_array<typename Array<PointIndexType>::Ptr>& offsetSegments,
827  size_t& segmentCount,
828  const Index binLog2Dim,
829  const Index bucketLog2Dim,
830  VoxelOffsetType* voxelOffsets = nullptr,
831  bool cellCenteredTransform = true)
832 {
833  using IndexPair = std::pair<PointIndexType, PointIndexType>;
834  using IndexPairList = std::deque<IndexPair>;
835  using IndexPairListPtr = SharedPtr<IndexPairList>;
836  using IndexPairListMap = std::map<Coord, IndexPairListPtr>;
837  using IndexPairListMapPtr = SharedPtr<IndexPairListMap>;
838 
839  size_t numTasks = 1, numThreads = size_t(tbb::task_scheduler_init::default_num_threads());
840  if (points.size() > (numThreads * 2)) numTasks = numThreads * 2;
841  else if (points.size() > numThreads) numTasks = numThreads;
842 
843  boost::scoped_array<IndexPairListMapPtr> bins(new IndexPairListMapPtr[numTasks]);
844 
846 
847  tbb::parallel_for(tbb::blocked_range<size_t>(0, numTasks),
848  BinOp(bins.get(), points, voxelOffsets, xform, binLog2Dim, bucketLog2Dim,
849  numTasks, cellCenteredTransform));
850 
851  std::set<Coord> uniqueCoords;
852 
853  for (size_t i = 0; i < numTasks; ++i) {
854  IndexPairListMap& idxMap = *bins[i];
855  for (typename IndexPairListMap::iterator it = idxMap.begin(); it != idxMap.end(); ++it) {
856  uniqueCoords.insert(it->first);
857  }
858  }
859 
860  std::vector<Coord> coords(uniqueCoords.begin(), uniqueCoords.end());
861  uniqueCoords.clear();
862 
863  segmentCount = coords.size();
864 
865  using SegmentPtr = typename Array<PointIndexType>::Ptr;
866 
867  indexSegments.reset(new SegmentPtr[segmentCount]);
868  offsetSegments.reset(new SegmentPtr[segmentCount]);
869 
870  using MergeOp = MergeBinsOp<PointIndexType>;
871 
872  tbb::parallel_for(tbb::blocked_range<size_t>(0, segmentCount),
873  MergeOp(bins.get(), indexSegments.get(), offsetSegments.get(), &coords[0], numTasks));
874 }
875 
876 
877 template<typename PointIndexType, typename VoxelOffsetType, typename PointArray>
878 inline void partition(
879  const PointArray& points,
880  const math::Transform& xform,
881  const Index bucketLog2Dim,
882  boost::scoped_array<PointIndexType>& pointIndices,
883  boost::scoped_array<PointIndexType>& pageOffsets,
884  PointIndexType& pageCount,
885  boost::scoped_array<VoxelOffsetType>& voxelOffsets,
886  bool recordVoxelOffsets,
887  bool cellCenteredTransform)
888 {
889  if (recordVoxelOffsets) voxelOffsets.reset(new VoxelOffsetType[points.size()]);
890  else voxelOffsets.reset();
891 
892  const Index binLog2Dim = 5u;
893  // note: Bins span a (2^(binLog2Dim + bucketLog2Dim))^3 voxel region,
894  // i.e. bucketLog2Dim = 3 and binLog2Dim = 5 corresponds to a
895  // (2^8)^3 = 256^3 voxel region.
896 
897 
898  size_t numSegments = 0;
899 
900  boost::scoped_array<typename Array<PointIndexType>::Ptr> indexSegments;
901  boost::scoped_array<typename Array<PointIndexType>::Ptr> offestSegments;
902 
903  binAndSegment<PointIndexType, VoxelOffsetType, PointArray>(points, xform,
904  indexSegments, offestSegments, numSegments, binLog2Dim, bucketLog2Dim,
905  voxelOffsets.get(), cellCenteredTransform);
906 
907  const tbb::blocked_range<size_t> segmentRange(0, numSegments);
908 
909  using IndexArray = boost::scoped_array<PointIndexType>;
910  boost::scoped_array<IndexArray> pageOffsetArrays(new IndexArray[numSegments]);
911 
912  const Index binVolume = 1u << (3u * binLog2Dim);
913 
914  tbb::parallel_for(segmentRange, OrderSegmentsOp<PointIndexType>
915  (indexSegments.get(), offestSegments.get(), pageOffsetArrays.get(), binVolume));
916 
917  indexSegments.reset();
918 
919  pageCount = 0;
920  for (size_t n = 0; n < numSegments; ++n) {
921  pageCount += pageOffsetArrays[n][0] - 1;
922  }
923 
924  pageOffsets.reset(new PointIndexType[pageCount + 1]);
925 
926  PointIndexType count = 0;
927  for (size_t n = 0, idx = 0; n < numSegments; ++n) {
928 
929  PointIndexType* offsets = pageOffsetArrays[n].get();
930  size_t size = size_t(offsets[0]);
931 
932  for (size_t i = 1; i < size; ++i) {
933  pageOffsets[idx++] = count;
934  count += offsets[i];
935  }
936  }
937 
938  pageOffsets[pageCount] = count;
939 
940  pointIndices.reset(new PointIndexType[points.size()]);
941 
942  std::vector<PointIndexType*> indexArray;
943  indexArray.reserve(numSegments);
944 
945  PointIndexType* index = pointIndices.get();
946  for (size_t n = 0; n < numSegments; ++n) {
947  indexArray.push_back(index);
948  index += offestSegments[n]->size();
949  }
950 
951  tbb::parallel_for(segmentRange,
952  MoveSegmentDataOp<PointIndexType>(indexArray, offestSegments.get()));
953 }
954 
955 
956 } // namespace point_partitioner_internal
957 
958 
960 
961 
962 template<typename PointIndexType, Index BucketLog2Dim>
964  : mPointIndices(nullptr)
965  , mVoxelOffsets(nullptr)
966  , mPageOffsets(nullptr)
967  , mPageCoordinates(nullptr)
968  , mPageCount(0)
969  , mUsingCellCenteredTransform(true)
970 {
971 }
972 
973 
974 template<typename PointIndexType, Index BucketLog2Dim>
975 inline void
977 {
978  mPageCount = 0;
979  mUsingCellCenteredTransform = true;
980  mPointIndices.reset();
981  mVoxelOffsets.reset();
982  mPageOffsets.reset();
983  mPageCoordinates.reset();
984 }
985 
986 
987 template<typename PointIndexType, Index BucketLog2Dim>
988 inline void
990 {
991  const IndexType tmpLhsPageCount = mPageCount;
992  mPageCount = rhs.mPageCount;
993  rhs.mPageCount = tmpLhsPageCount;
994 
995  mPointIndices.swap(rhs.mPointIndices);
996  mVoxelOffsets.swap(rhs.mVoxelOffsets);
997  mPageOffsets.swap(rhs.mPageOffsets);
998  mPageCoordinates.swap(rhs.mPageCoordinates);
999 
1000  bool lhsCellCenteredTransform = mUsingCellCenteredTransform;
1001  mUsingCellCenteredTransform = rhs.mUsingCellCenteredTransform;
1002  rhs.mUsingCellCenteredTransform = lhsCellCenteredTransform;
1003 }
1004 
1005 
1006 template<typename PointIndexType, Index BucketLog2Dim>
1009 {
1010  assert(bool(mPointIndices) && bool(mPageCount));
1011  return IndexIterator(
1012  mPointIndices.get() + mPageOffsets[n],
1013  mPointIndices.get() + mPageOffsets[n + 1]);
1014 }
1015 
1016 
1017 template<typename PointIndexType, Index BucketLog2Dim>
1018 template<typename PointArray>
1019 inline void
1021  const PointArray& points,
1022  const math::Transform& xform,
1023  bool voxelOrder,
1024  bool recordVoxelOffsets,
1025  bool cellCenteredTransform)
1026 {
1027  mUsingCellCenteredTransform = cellCenteredTransform;
1028 
1029  point_partitioner_internal::partition(points, xform, BucketLog2Dim,
1030  mPointIndices, mPageOffsets, mPageCount, mVoxelOffsets,
1031  (voxelOrder || recordVoxelOffsets), cellCenteredTransform);
1032 
1033  const tbb::blocked_range<size_t> pageRange(0, mPageCount);
1034  mPageCoordinates.reset(new Coord[mPageCount]);
1035 
1036  tbb::parallel_for(pageRange,
1038  (mPageCoordinates, mPointIndices, mPageOffsets, points, xform,
1039  BucketLog2Dim, cellCenteredTransform));
1040 
1041  if (mVoxelOffsets && voxelOrder) {
1042  tbb::parallel_for(pageRange, point_partitioner_internal::VoxelOrderOp<
1043  IndexType, BucketLog2Dim>(mPointIndices, mPageOffsets, mVoxelOffsets));
1044  }
1045 
1046  if (mVoxelOffsets && !recordVoxelOffsets) {
1047  mVoxelOffsets.reset();
1048  }
1049 }
1050 
1051 
1052 template<typename PointIndexType, Index BucketLog2Dim>
1053 template<typename PointArray>
1056  const PointArray& points,
1057  const math::Transform& xform,
1058  bool voxelOrder,
1059  bool recordVoxelOffsets,
1060  bool cellCenteredTransform)
1061 {
1062  Ptr ret(new PointPartitioner());
1063  ret->construct(points, xform, voxelOrder, recordVoxelOffsets, cellCenteredTransform);
1064  return ret;
1065 }
1066 
1067 
1069 
1070 
1071 } // namespace tools
1072 } // namespace OPENVDB_VERSION_NAME
1073 } // namespace openvdb
1074 
1075 
1076 #endif // OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED
1077 
1078 // Copyright (c) DreamWorks Animation LLC
1079 // All rights reserved. This software is distributed under the
1080 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
SharedPtr< Array > Ptr
Definition: PointPartitioner.h:430
OrderSegmentsOp(SegmentPtr *indexSegments, SegmentPtr *offestSegments, IndexArray *pageOffsetArrays, Index binVolume)
Definition: PointPartitioner.h:739
PointArray const *const mPoints
Definition: PointPartitioner.h:723
Definition: PointPartitioner.h:106
typename PointArray::PosType PosType
Definition: PointPartitioner.h:616
size_t size() const
Definition: PointPartitioner.h:434
PointIndexType *const mOrderedIndexArray
Definition: PointPartitioner.h:293
bool empty() const
true if the container size is 0, false otherwise.
Definition: PointPartitioner.h:158
static Ptr create(const PointArray &points, const math::Transform &xform, bool voxelOrder=false, bool recordVoxelOffsets=false, bool cellCenteredTransform=true)
Partitions point indices into BucketLog2Dim aligned buckets.
Definition: PointPartitioner.h:1055
SegmentPtr *const mIndexSegments
Definition: PointPartitioner.h:810
Coord *const mCoordinates
Definition: PointPartitioner.h:414
math::Transform const mXForm
Definition: PointPartitioner.h:418
size_t size() const
Number of point indices in the iterator range.
Definition: PointPartitioner.h:217
std::pair< PointIndexType, PointIndexType > IndexPair
Definition: PointPartitioner.h:617
bool next()
Advance to the next item.
Definition: PointPartitioner.h:231
bool increment()
Definition: PointPartitioner.h:232
PointIndexType const *const mPointOrder
Definition: PointPartitioner.h:294
PointIndexType *const mIndices
Definition: PointPartitioner.h:363
void swap(PointPartitioner &)
Exchanges the content of the container by another.
Definition: PointPartitioner.h:989
void clear()
Definition: PointPartitioner.h:439
typename boost::int_t< 1+(3 *BucketLog2Dim)>::least VoxelOffsetType
Definition: PointPartitioner.h:115
PointIndexType IndexType
Definition: PointPartitioner.h:114
PointIndexType const *const mIndices
Definition: PointPartitioner.h:415
const std::enable_if<!VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:133
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:458
T * data()
Definition: PointPartitioner.h:436
IndexPairListMapPtr *const mData
Definition: PointPartitioner.h:722
CoordBBox getBBox(size_t n) const
Returns the coordinate-aligned bounding box for bucket n.
Definition: PointPartitioner.h:170
bool operator==(const IndexIterator &other) const
Equality operators.
Definition: PointPartitioner.h:235
IndexIterator & operator++()
Advance to the next item.
Definition: PointPartitioner.h:228
std::pair< PointIndexType, PointIndexType > IndexPair
Definition: PointPartitioner.h:497
int const mLog2Dim
Definition: PointPartitioner.h:419
bool operator!=(const IndexIterator &other) const
Definition: PointPartitioner.h:236
std::map< Coord, IndexPairListPtr > IndexPairListMap
Definition: PointPartitioner.h:500
boost::scoped_array< PointIndexType > IndexArray
Definition: PointPartitioner.h:304
Index32 Index
Definition: Types.h:61
math::Histogram histogram(const IterT &iter, double minVal, double maxVal, size_t numBins=10, bool threaded=true)
Iterate over a scalar grid and compute a histogram of the values of the voxels that are visited...
Definition: Statistics.h:368
std::shared_ptr< T > SharedPtr
Definition: Types.h:139
IndexIterator indices(size_t n) const
Returns the point indices for bucket n.
Definition: PointPartitioner.h:1008
IndexIterator(IndexType *begin=nullptr, IndexType *end=nullptr)
Definition: PointPartitioner.h:210
boost::scoped_array< VoxelOffsetType > VoxelOffsetArray
Definition: PointPartitioner.h:303
LeafNodeOriginOp(CoordArray &coordinates, const IndexArray &indices, const IndexArray &pages, const PointArray &points, const math::Transform &m, int log2dim, bool cellCenteredTransform)
Definition: PointPartitioner.h:375
VoxelOrderOp(IndexArray &indices, const IndexArray &pages, const VoxelOffsetArray &offsets)
Definition: PointPartitioner.h:306
PointIndexType const *const mPages
Definition: PointPartitioner.h:364
PointIndexType *const mPointOrder
Definition: PointPartitioner.h:270
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:287
const Coord & origin(size_t n) const
Returns the origin coordinate for bucket n.
Definition: PointPartitioner.h:175
SharedPtr< IndexPairListMap > IndexPairListMapPtr
Definition: PointPartitioner.h:501
VoxelOffsetType const *const mVoxelOffsets
Definition: PointPartitioner.h:365
bool const mCellCenteredTransform
Definition: PointPartitioner.h:729
boost::scoped_array< PointIndexType > IndexArray
Definition: PointPartitioner.h:736
math::Transform const mXForm
Definition: PointPartitioner.h:725
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:516
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:128
ComputePointOrderOp(PointIndexType *pointOrder, const PointIndexType *bucketCounters, const PointIndexType *bucketOffsets)
Definition: PointPartitioner.h:256
bool test() const
Definition: PointPartitioner.h:225
bool usingCellCenteredTransform() const
Returns true if this point partitioning was constructed using a cell-centered transform.
Definition: PointPartitioner.h:184
Index const mBucketLog2Dim
Definition: PointPartitioner.h:727
PointIndexType const *const mBucketCounters
Definition: PointPartitioner.h:271
PointIndexType IndexType
Definition: PointPartitioner.h:208
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:264
Definition: Exceptions.h:40
Index const mBinLog2Dim
Definition: PointPartitioner.h:726
typename Segment::Ptr SegmentPtr
Definition: PointPartitioner.h:451
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:748
const T * data() const
Definition: PointPartitioner.h:437
typename Segment::Ptr SegmentPtr
Definition: PointPartitioner.h:495
BinPointIndicesOp(IndexPairListMapPtr *data, const PointArray &points, VoxelOffsetType *voxelOffsets, const math::Transform &m, Index binLog2Dim, Index bucketLog2Dim, size_t numSegments, bool cellCenteredTransform)
Definition: PointPartitioner.h:623
IndexArray *const mPageOffsetArrays
Definition: PointPartitioner.h:812
PointIndexType const *const mIndices
Definition: PointPartitioner.h:295
void partition(const PointArray &points, const math::Transform &xform, const Index bucketLog2Dim, boost::scoped_array< PointIndexType > &pointIndices, boost::scoped_array< PointIndexType > &pageOffsets, PointIndexType &pageCount, boost::scoped_array< VoxelOffsetType > &voxelOffsets, bool recordVoxelOffsets, bool cellCenteredTransform)
Definition: PointPartitioner.h:878
SharedPtr< IndexPairListMap > IndexPairListMapPtr
Definition: PointPartitioner.h:621
SharedPtr< PointPartitioner > Ptr
Definition: PointPartitioner.h:111
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:642
std::pair< Index, Index > IndexPair
Definition: PointMove.h:195
IndexType & operator*()
Returns the item to which this iterator is currently pointing.
Definition: PointPartitioner.h:220
Index const mBinVolume
Definition: PointPartitioner.h:813
boost::scoped_array< PointIndexType > IndexArray
Definition: PointPartitioner.h:372
void construct(const PointArray &points, const math::Transform &xform, bool voxelOrder=false, bool recordVoxelOffsets=false, bool cellCenteredTransform=true)
Partitions point indices into BucketLog2Dim aligned buckets.
Definition: PointPartitioner.h:1020
SegmentPtr *const mOffsetSegments
Definition: PointPartitioner.h:811
Index64 pointCount(const PointDataTreeT &tree, const FilterT &filter=NullFilter(), const bool inCoreOnly=false, const bool threaded=true)
Count the total number of points in a PointDataTree.
Definition: PointCount.h:115
std::deque< IndexPair > IndexPairList
Definition: PointPartitioner.h:498
size_t const mNumSegments
Definition: PointPartitioner.h:728
Array(size_t size)
Definition: PointPartitioner.h:432
PointPartitioner()
Definition: PointPartitioner.h:963
boost::scoped_array< Coord > CoordArray
Definition: PointPartitioner.h:373
const IndexType & operator*() const
Definition: PointPartitioner.h:221
bool operator==(const Vec3< T0 > &v0, const Vec3< T1 > &v1)
Equality operator, does exact floating point comparisons.
Definition: Vec3.h:498
bool const mCellCenteredTransform
Definition: PointPartitioner.h:420
Definition: Transform.h:66
boost::scoped_array< VoxelOffsetType > VoxelOffsetArray
Definition: PointPartitioner.h:116
SharedPtr< IndexPairList > IndexPairListPtr
Definition: PointPartitioner.h:619
typename boost::int_t< 1+(3 *BucketLog2Dim)>::least VoxelOffsetType
Definition: PointPartitioner.h:302
void clear()
Removes all data and frees up memory.
Definition: PointPartitioner.h:976
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:313
MergeBinsOp(IndexPairListMapPtr *bins, SegmentPtr *indexSegments, SegmentPtr *offsetSegments, Coord *coords, size_t numSegments)
Definition: PointPartitioner.h:503
CreateOrderedPointIndexArrayOp(PointIndexType *orderedIndexArray, const PointIndexType *pointOrder, const PointIndexType *indices)
Definition: PointPartitioner.h:279
MoveSegmentDataOp(std::vector< PointIndexType * > &indexLists, SegmentPtr *segments)
Definition: PointPartitioner.h:453
std::deque< IndexPair > IndexPairList
Definition: PointPartitioner.h:618
size_t size() const
Returns the number of buckets.
Definition: PointPartitioner.h:155
const VoxelOffsetArray & voxelOffsets() const
Returns a list of LeafNode voxel offsets for the points.
Definition: PointPartitioner.h:179
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:180
std::map< Coord, IndexPairListPtr > IndexPairListMap
Definition: PointPartitioner.h:620
PointArray const *const mPoints
Definition: PointPartitioner.h:417
void binAndSegment(const PointArray &points, const math::Transform &xform, boost::scoped_array< typename Array< PointIndexType >::Ptr > &indexSegments, boost::scoped_array< typename Array< PointIndexType >::Ptr > &offsetSegments, size_t &segmentCount, const Index binLog2Dim, const Index bucketLog2Dim, VoxelOffsetType *voxelOffsets=nullptr, bool cellCenteredTransform=true)
Segment points using one level of least significant digit radix bins.
Definition: PointPartitioner.h:822
PointIndexType const *const mPages
Definition: PointPartitioner.h:416
VoxelOffsetType *const mVoxelOffsets
Definition: PointPartitioner.h:724
SharedPtr< IndexPairList > IndexPairListPtr
Definition: PointPartitioner.h:499
void operator()(const tbb::blocked_range< size_t > &range) const
Definition: PointPartitioner.h:388
std::vector< Index > IndexArray
Definition: PointMove.h:188
typename Array< PointIndexType >::Ptr SegmentPtr
Definition: PointPartitioner.h:737
void reset()
Rewind to first item.
Definition: PointPartitioner.h:214
PointIndexType const *const mBucketOffsets
Definition: PointPartitioner.h:272
SharedPtr< const PointPartitioner > ConstPtr
Definition: PointPartitioner.h:112