18 #ifndef OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED
19 #define OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED
24 #include <boost/integer.hpp>
25 #include <boost/scoped_array.hpp>
27 #include <tbb/blocked_range.h>
28 #include <tbb/parallel_for.h>
29 #include <tbb/task_scheduler_init.h>
78 template<
typename Po
intIndexType = u
int32_t, Index BucketLog2Dim = 3>
82 enum { LOG2DIM = BucketLog2Dim };
106 template<
typename Po
intArray>
108 bool voxelOrder =
false,
bool recordVoxelOffsets =
false,
109 bool cellCenteredTransform =
true);
121 template<
typename Po
intArray>
123 bool voxelOrder =
false,
bool recordVoxelOffsets =
false,
124 bool cellCenteredTransform =
true);
128 size_t size()
const {
return mPageCount; }
131 bool empty()
const {
return mPageCount == 0; }
140 IndexIterator indices(
size_t n)
const;
144 return CoordBBox::createCube(mPageCoordinates[n], (1u << BucketLog2Dim));
148 const Coord&
origin(
size_t n)
const {
return mPageCoordinates[n]; }
164 boost::scoped_array<IndexType> mPointIndices;
165 VoxelOffsetArray mVoxelOffsets;
167 boost::scoped_array<IndexType> mPageOffsets;
168 boost::scoped_array<Coord> mPageCoordinates;
169 IndexType mPageCount;
170 bool mUsingCellCenteredTransform;
177 template<
typename Po
intIndexType, Index BucketLog2Dim>
184 : mBegin(begin), mEnd(end), mItem(begin) {}
190 size_t size()
const {
return mEnd - mBegin; }
197 operator bool()
const {
return mItem < mEnd; }
198 bool test()
const {
return mItem < mEnd; }
204 bool next() { this->operator++();
return this->test(); }
223 namespace point_partitioner_internal {
226 template<
typename Po
intIndexType>
230 const PointIndexType* bucketCounters,
const PointIndexType* bucketOffsets)
231 : mPointOrder(pointOrder)
232 , mBucketCounters(bucketCounters)
233 , mBucketOffsets(bucketOffsets)
237 void operator()(
const tbb::blocked_range<size_t>& range)
const {
238 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
239 mPointOrder[n] += mBucketCounters[mBucketOffsets[n]];
249 template<
typename Po
intIndexType>
253 const PointIndexType* pointOrder,
const PointIndexType* indices)
254 : mOrderedIndexArray(orderedIndexArray)
255 , mPointOrder(pointOrder)
260 void operator()(
const tbb::blocked_range<size_t>& range)
const {
261 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
262 mOrderedIndexArray[mPointOrder[n]] = mIndices[n];
272 template<
typename Po
intIndexType, Index BucketLog2Dim>
280 : mIndices(indices.get())
281 , mPages(pages.get())
282 , mVoxelOffsets(offsets.get())
286 void operator()(
const tbb::blocked_range<size_t>& range)
const {
289 for (
size_t n(range.begin()), N(range.end()); n != N; ++n) {
293 const PointIndexType voxelCount = 1 << (3 * BucketLog2Dim);
297 boost::scoped_array<PointIndexType> sortedIndices(
new PointIndexType[
pointCount]);
298 boost::scoped_array<PointIndexType>
histogram(
new PointIndexType[voxelCount]);
300 for (
size_t n(range.begin()), N(range.end()); n != N; ++n) {
302 PointIndexType *
const indices = mIndices + mPages[n];
306 for (PointIndexType i = 0; i <
pointCount; ++i) {
307 offsets[i] = mVoxelOffsets[ indices[i] ];
311 memset(&
histogram[0], 0, voxelCount *
sizeof(PointIndexType));
314 for (PointIndexType i = 0; i <
pointCount; ++i) {
318 PointIndexType count = 0, startOffset;
319 for (
int i = 0; i < int(voxelCount); ++i) {
328 for (PointIndexType i = 0; i <
pointCount; ++i) {
329 sortedIndices[
histogram[ offsets[i] ]++ ] = indices[i];
332 memcpy(&indices[0], &sortedIndices[0],
sizeof(PointIndexType) *
pointCount);
342 template<
typename Po
intArray,
typename Po
intIndexType>
351 : mCoordinates(coordinates.get())
352 , mIndices(indices.get())
353 , mPages(pages.get())
357 , mCellCenteredTransform(cellCenteredTransform)
361 void operator()(
const tbb::blocked_range<size_t>& range)
const {
363 using PosType =
typename PointArray::PosType;
365 const bool cellCentered = mCellCenteredTransform;
366 const int mask = ~((1 << mLog2Dim) - 1);
370 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
372 mPoints->getPos(mIndices[mPages[n]], pos);
374 if (std::isfinite(pos[0]) && std::isfinite(pos[1]) && std::isfinite(pos[2])) {
375 ijk = cellCentered ? mXForm.worldToIndexCellCentered(pos) :
376 mXForm.worldToIndexNodeCentered(pos);
382 mCoordinates[n] = ijk;
405 Array(
size_t size) : mSize(size), mData(new T[size]) { }
407 size_t size()
const {
return mSize; }
409 T*
data() {
return mData.get(); }
410 const T*
data()
const {
return mData.get(); }
412 void clear() { mSize = 0; mData.reset(); }
416 boost::scoped_array<T> mData;
420 template<
typename Po
intIndexType>
427 : mIndexLists(&indexLists[0]), mSegments(segments)
431 void operator()(
const tbb::blocked_range<size_t>& range)
const {
432 for (
size_t n(range.begin()), N(range.end()); n != N; ++n) {
433 PointIndexType* indices = mIndexLists[n];
436 tbb::parallel_for(tbb::blocked_range<size_t>(0, segment->size()),
437 CopyData(indices, segment->data()));
447 CopyData(PointIndexType* lhs,
const PointIndexType* rhs) : mLhs(lhs), mRhs(rhs) { }
449 void operator()(
const tbb::blocked_range<size_t>& range)
const {
450 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
455 PointIndexType *
const mLhs;
456 PointIndexType
const *
const mRhs;
459 PointIndexType *
const *
const mIndexLists;
460 SegmentPtr *
const mSegments;
464 template<
typename Po
intIndexType>
470 using IndexPair = std::pair<PointIndexType, PointIndexType>;
482 , mIndexSegments(indexSegments)
483 , mOffsetSegments(offsetSegments)
485 , mNumSegments(numSegments)
489 void operator()(
const tbb::blocked_range<size_t>& range)
const {
491 std::vector<IndexPairListPtr*> data;
492 std::vector<PointIndexType> arrayOffsets;
494 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
496 const Coord& ijk = mCoords[n];
497 size_t numIndices = 0;
501 for (
size_t i = 0, I = mNumSegments; i < I; ++i) {
504 typename IndexPairListMap::iterator iter = idxMap.find(ijk);
506 if (iter != idxMap.end() && iter->second) {
509 data.push_back(&idxListPtr);
510 numIndices += idxListPtr->size();
514 if (data.empty() || numIndices == 0)
continue;
517 SegmentPtr& offsetSegment = mOffsetSegments[n];
522 arrayOffsets.clear();
523 arrayOffsets.reserve(data.size());
525 for (
size_t i = 0, count = 0, I = data.size(); i < I; ++i) {
526 arrayOffsets.push_back(PointIndexType(count));
527 count += (*data[i])->size();
530 tbb::parallel_for(tbb::blocked_range<size_t>(0, data.size()),
531 CopyData(&data[0], &arrayOffsets[0], indexSegment->data(), offsetSegment->data()));
539 CopyData(IndexPairListPtr** indexLists,
540 const PointIndexType* arrayOffsets,
541 PointIndexType* indices,
542 PointIndexType* offsets)
543 : mIndexLists(indexLists)
544 , mArrayOffsets(arrayOffsets)
550 void operator()(
const tbb::blocked_range<size_t>& range)
const {
552 using CIter =
typename IndexPairList::const_iterator;
554 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
556 const PointIndexType arrayOffset = mArrayOffsets[n];
557 PointIndexType* indexPtr = &mIndices[arrayOffset];
558 PointIndexType* offsetPtr = &mOffsets[arrayOffset];
560 IndexPairListPtr& list = *mIndexLists[n];
562 for (CIter it = list->begin(), end = list->end(); it != end; ++it) {
564 *indexPtr++ = data.first;
565 *offsetPtr++ = data.second;
572 IndexPairListPtr *
const *
const mIndexLists;
573 PointIndexType
const *
const mArrayOffsets;
574 PointIndexType *
const mIndices;
575 PointIndexType *
const mOffsets;
578 IndexPairListMapPtr *
const mBins;
579 SegmentPtr *
const mIndexSegments;
580 SegmentPtr *
const mOffsetSegments;
581 Coord
const *
const mCoords;
582 size_t const mNumSegments;
586 template<
typename Po
intArray,
typename Po
intIndexType,
typename VoxelOffsetType>
590 using IndexPair = std::pair<PointIndexType, PointIndexType>;
598 VoxelOffsetType* voxelOffsets,
603 bool cellCenteredTransform)
606 , mVoxelOffsets(voxelOffsets)
608 , mBinLog2Dim(binLog2Dim)
609 , mBucketLog2Dim(bucketLog2Dim)
610 , mNumSegments(numSegments)
611 , mCellCenteredTransform(cellCenteredTransform)
615 void operator()(
const tbb::blocked_range<size_t>& range)
const {
617 const Index log2dim = mBucketLog2Dim;
618 const Index log2dim2 = 2 * log2dim;
619 const Index bucketMask = (1u << log2dim) - 1u;
621 const Index binLog2dim = mBinLog2Dim;
622 const Index binLog2dim2 = 2 * binLog2dim;
624 const Index binMask = (1u << (log2dim + binLog2dim)) - 1u;
625 const Index invBinMask = ~binMask;
628 Coord ijk(0, 0, 0), loc(0, 0, 0), binCoord(0, 0, 0), lastBinCoord(1, 2, 3);
631 PointIndexType bucketOffset = 0;
632 VoxelOffsetType voxelOffset = 0;
634 const bool cellCentered = mCellCenteredTransform;
636 const size_t numPoints = mPoints->size();
637 const size_t segmentSize = numPoints / mNumSegments;
639 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
645 const bool isLastSegment = (n + 1) >= mNumSegments;
647 const size_t start = n * segmentSize;
648 const size_t end = isLastSegment ? numPoints : (start + segmentSize);
650 for (
size_t i = start; i != end; ++i) {
652 mPoints->getPos(i, pos);
654 if (std::isfinite(pos[0]) && std::isfinite(pos[1]) && std::isfinite(pos[2])) {
655 ijk = cellCentered ? mXForm.worldToIndexCellCentered(pos) :
656 mXForm.worldToIndexNodeCentered(pos);
659 loc[0] = ijk[0] & bucketMask;
660 loc[1] = ijk[1] & bucketMask;
661 loc[2] = ijk[2] & bucketMask;
662 voxelOffset = VoxelOffsetType(
663 (loc[0] << log2dim2) + (loc[1] << log2dim) + loc[2]);
666 binCoord[0] = ijk[0] & invBinMask;
667 binCoord[1] = ijk[1] & invBinMask;
668 binCoord[2] = ijk[2] & invBinMask;
678 bucketOffset = PointIndexType(
679 (ijk[0] << binLog2dim2) + (ijk[1] << binLog2dim) + ijk[2]);
681 if (lastBinCoord != binCoord) {
682 lastBinCoord = binCoord;
685 idxList = idxListPtr.get();
688 idxList->push_back(
IndexPair(PointIndexType(i), bucketOffset));
689 if (mVoxelOffsets) mVoxelOffsets[i] = voxelOffset;
706 template<
typename Po
intIndexType>
714 : mIndexSegments(indexSegments)
715 , mOffsetSegments(offestSegments)
716 , mPageOffsetArrays(pageOffsetArrays)
717 , mBinVolume(binVolume)
721 void operator()(
const tbb::blocked_range<size_t>& range)
const {
723 const size_t bucketCountersSize = size_t(mBinVolume);
724 IndexArray bucketCounters(
new PointIndexType[bucketCountersSize]);
726 size_t maxSegmentSize = 0;
727 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
728 maxSegmentSize =
std::max(maxSegmentSize, mIndexSegments[n]->size());
731 IndexArray bucketIndices(
new PointIndexType[maxSegmentSize]);
734 for (
size_t n = range.begin(), N = range.end(); n != N; ++n) {
736 memset(bucketCounters.get(), 0,
sizeof(PointIndexType) * bucketCountersSize);
738 const size_t segmentSize = mOffsetSegments[n]->size();
739 PointIndexType* offsets = mOffsetSegments[n]->data();
743 for (
size_t i = 0; i < segmentSize; ++i) {
744 bucketIndices[i] = bucketCounters[offsets[i]]++;
747 PointIndexType nonemptyBucketCount = 0;
748 for (
size_t i = 0; i < bucketCountersSize; ++i) {
749 nonemptyBucketCount +=
static_cast<PointIndexType
>(bucketCounters[i] != 0);
753 IndexArray& pageOffsets = mPageOffsetArrays[n];
754 pageOffsets.reset(
new PointIndexType[nonemptyBucketCount + 1]);
755 pageOffsets[0] = nonemptyBucketCount + 1;
758 PointIndexType count = 0, idx = 1;
759 for (
size_t i = 0; i < bucketCountersSize; ++i) {
760 if (bucketCounters[i] != 0) {
761 pageOffsets[idx] = bucketCounters[i];
762 bucketCounters[i] = count;
763 count += pageOffsets[idx];
768 PointIndexType* indices = mIndexSegments[n]->data();
769 const tbb::blocked_range<size_t> segmentRange(0, segmentSize);
774 bucketIndices.get(), bucketCounters.get(), offsets));
777 offsets, bucketIndices.get(), indices));
779 mIndexSegments[n]->clear();
794 template<
typename Po
intIndexType,
typename VoxelOffsetType,
typename Po
intArray>
800 size_t& segmentCount,
801 const Index binLog2Dim,
802 const Index bucketLog2Dim,
803 VoxelOffsetType* voxelOffsets =
nullptr,
804 bool cellCenteredTransform =
true)
806 using IndexPair = std::pair<PointIndexType, PointIndexType>;
807 using IndexPairList = std::deque<IndexPair>;
809 using IndexPairListMap = std::map<Coord, IndexPairListPtr>;
812 size_t numTasks = 1, numThreads = size_t(tbb::task_scheduler_init::default_num_threads());
813 if (points.size() > (numThreads * 2)) numTasks = numThreads * 2;
814 else if (points.size() > numThreads) numTasks = numThreads;
816 boost::scoped_array<IndexPairListMapPtr> bins(
new IndexPairListMapPtr[numTasks]);
820 tbb::parallel_for(tbb::blocked_range<size_t>(0, numTasks),
821 BinOp(bins.get(), points, voxelOffsets, xform, binLog2Dim, bucketLog2Dim,
822 numTasks, cellCenteredTransform));
824 std::set<Coord> uniqueCoords;
826 for (
size_t i = 0; i < numTasks; ++i) {
827 IndexPairListMap& idxMap = *bins[i];
828 for (
typename IndexPairListMap::iterator it = idxMap.begin(); it != idxMap.end(); ++it) {
829 uniqueCoords.insert(it->first);
833 std::vector<Coord> coords(uniqueCoords.begin(), uniqueCoords.end());
834 uniqueCoords.clear();
836 segmentCount = coords.size();
840 indexSegments.reset(
new SegmentPtr[segmentCount]);
841 offsetSegments.reset(
new SegmentPtr[segmentCount]);
845 tbb::parallel_for(tbb::blocked_range<size_t>(0, segmentCount),
846 MergeOp(bins.get(), indexSegments.get(), offsetSegments.get(), &coords[0], numTasks));
850 template<
typename Po
intIndexType,
typename VoxelOffsetType,
typename Po
intArray>
854 const Index bucketLog2Dim,
855 boost::scoped_array<PointIndexType>& pointIndices,
856 boost::scoped_array<PointIndexType>& pageOffsets,
857 PointIndexType& pageCount,
858 boost::scoped_array<VoxelOffsetType>& voxelOffsets,
859 bool recordVoxelOffsets,
860 bool cellCenteredTransform)
862 if (recordVoxelOffsets) voxelOffsets.reset(
new VoxelOffsetType[points.size()]);
863 else voxelOffsets.reset();
865 const Index binLog2Dim = 5u;
871 size_t numSegments = 0;
873 boost::scoped_array<typename Array<PointIndexType>::Ptr> indexSegments;
874 boost::scoped_array<typename Array<PointIndexType>::Ptr> offestSegments;
876 binAndSegment<PointIndexType, VoxelOffsetType, PointArray>(points, xform,
877 indexSegments, offestSegments, numSegments, binLog2Dim, bucketLog2Dim,
878 voxelOffsets.get(), cellCenteredTransform);
880 const tbb::blocked_range<size_t> segmentRange(0, numSegments);
882 using IndexArray = boost::scoped_array<PointIndexType>;
883 boost::scoped_array<IndexArray> pageOffsetArrays(
new IndexArray[numSegments]);
885 const Index binVolume = 1u << (3u * binLog2Dim);
888 (indexSegments.get(), offestSegments.get(), pageOffsetArrays.get(), binVolume));
890 indexSegments.reset();
893 for (
size_t n = 0; n < numSegments; ++n) {
894 pageCount += pageOffsetArrays[n][0] - 1;
897 pageOffsets.reset(
new PointIndexType[pageCount + 1]);
899 PointIndexType count = 0;
900 for (
size_t n = 0, idx = 0; n < numSegments; ++n) {
902 PointIndexType* offsets = pageOffsetArrays[n].get();
903 size_t size = size_t(offsets[0]);
905 for (
size_t i = 1; i < size; ++i) {
906 pageOffsets[idx++] = count;
911 pageOffsets[pageCount] = count;
913 pointIndices.reset(
new PointIndexType[points.size()]);
915 std::vector<PointIndexType*> indexArray;
916 indexArray.reserve(numSegments);
918 PointIndexType* index = pointIndices.get();
919 for (
size_t n = 0; n < numSegments; ++n) {
920 indexArray.push_back(index);
921 index += offestSegments[n]->size();
924 tbb::parallel_for(segmentRange,
935 template<
typename Po
intIndexType, Index BucketLog2Dim>
937 : mPointIndices(nullptr)
938 , mVoxelOffsets(nullptr)
939 , mPageOffsets(nullptr)
940 , mPageCoordinates(nullptr)
942 , mUsingCellCenteredTransform(true)
947 template<
typename Po
intIndexType, Index BucketLog2Dim>
952 mUsingCellCenteredTransform =
true;
953 mPointIndices.reset();
954 mVoxelOffsets.reset();
955 mPageOffsets.reset();
956 mPageCoordinates.reset();
960 template<
typename Po
intIndexType, Index BucketLog2Dim>
964 const IndexType tmpLhsPageCount = mPageCount;
965 mPageCount = rhs.mPageCount;
966 rhs.mPageCount = tmpLhsPageCount;
968 mPointIndices.swap(rhs.mPointIndices);
969 mVoxelOffsets.swap(rhs.mVoxelOffsets);
970 mPageOffsets.swap(rhs.mPageOffsets);
971 mPageCoordinates.swap(rhs.mPageCoordinates);
973 bool lhsCellCenteredTransform = mUsingCellCenteredTransform;
974 mUsingCellCenteredTransform = rhs.mUsingCellCenteredTransform;
975 rhs.mUsingCellCenteredTransform = lhsCellCenteredTransform;
979 template<
typename Po
intIndexType, Index BucketLog2Dim>
983 assert(
bool(mPointIndices) &&
bool(mPageCount));
985 mPointIndices.get() + mPageOffsets[n],
986 mPointIndices.get() + mPageOffsets[n + 1]);
990 template<
typename Po
intIndexType, Index BucketLog2Dim>
991 template<
typename Po
intArray>
997 bool recordVoxelOffsets,
998 bool cellCenteredTransform)
1000 mUsingCellCenteredTransform = cellCenteredTransform;
1003 mPointIndices, mPageOffsets, mPageCount, mVoxelOffsets,
1004 (voxelOrder || recordVoxelOffsets), cellCenteredTransform);
1006 const tbb::blocked_range<size_t> pageRange(0, mPageCount);
1007 mPageCoordinates.reset(
new Coord[mPageCount]);
1009 tbb::parallel_for(pageRange,
1011 (mPageCoordinates, mPointIndices, mPageOffsets, points, xform,
1012 BucketLog2Dim, cellCenteredTransform));
1014 if (mVoxelOffsets && voxelOrder) {
1016 IndexType, BucketLog2Dim>(mPointIndices, mPageOffsets, mVoxelOffsets));
1019 if (mVoxelOffsets && !recordVoxelOffsets) {
1020 mVoxelOffsets.reset();
1025 template<
typename Po
intIndexType, Index BucketLog2Dim>
1026 template<
typename Po
intArray>
1032 bool recordVoxelOffsets,
1033 bool cellCenteredTransform)
1036 ret->construct(points, xform, voxelOrder, recordVoxelOffsets, cellCenteredTransform);
1049 #endif // OPENVDB_TOOLS_POINT_PARTITIONER_HAS_BEEN_INCLUDED