31 #ifndef OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
32 #define OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
34 #include <openvdb/Types.h>
35 #include <openvdb/math/FiniteDifference.h>
36 #include <openvdb/math/Operators.h>
37 #include <openvdb/math/Proximity.h>
38 #include <openvdb/tools/Morphology.h>
39 #include <openvdb/util/NullInterrupter.h>
40 #include <openvdb/util/Util.h>
42 #include <tbb/blocked_range.h>
43 #include <tbb/parallel_for.h>
44 #include <tbb/parallel_reduce.h>
75 template<
typename Gr
idType>
76 inline typename GridType::Ptr
78 const openvdb::math::Transform& xform,
79 const std::vector<Vec3s>& points,
80 const std::vector<Vec3I>& triangles,
99 template<
typename Gr
idType>
100 inline typename GridType::Ptr
102 const openvdb::math::Transform& xform,
103 const std::vector<Vec3s>& points,
104 const std::vector<Vec4I>& quads,
124 template<
typename Gr
idType>
125 inline typename GridType::Ptr
127 const openvdb::math::Transform& xform,
128 const std::vector<Vec3s>& points,
129 const std::vector<Vec3I>& triangles,
130 const std::vector<Vec4I>& quads,
152 template<
typename Gr
idType>
153 inline typename GridType::Ptr
155 const openvdb::math::Transform& xform,
156 const std::vector<Vec3s>& points,
157 const std::vector<Vec3I>& triangles,
158 const std::vector<Vec4I>& quads,
177 template<
typename Gr
idType>
178 inline typename GridType::Ptr
180 const openvdb::math::Transform& xform,
181 const std::vector<Vec3s>& points,
182 const std::vector<Vec3I>& triangles,
183 const std::vector<Vec4I>& quads,
195 template<
typename FloatGr
idT,
typename InterruptT = util::NullInterrupter>
203 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
206 MeshToVolume(openvdb::math::Transform::Ptr&,
int conversionFlags = 0,
207 InterruptT *interrupter = NULL,
int signSweeps = 1);
220 void convertToLevelSet(
221 const std::vector<Vec3s>& pointList,
222 const std::vector<Vec4I>& polygonList,
234 void convertToUnsignedDistanceField(
const std::vector<Vec3s>& pointList,
235 const std::vector<Vec4I>& polygonList,
FloatValueT exBandWidth);
240 typename FloatGridT::Ptr
distGridPtr()
const {
return mDistGrid; }
250 void doConvert(
const std::vector<Vec3s>&,
const std::vector<Vec4I>&,
251 FloatValueT exBandWidth, FloatValueT inBandWidth,
bool unsignedDistField =
false);
254 return mInterrupter && mInterrupter->wasInterrupted(percent);
257 openvdb::math::Transform::Ptr mTransform;
258 int mConversionFlags, mSignSweeps;
260 typename FloatGridT::Ptr mDistGrid;
261 typename IntGridT::Ptr mIndexGrid;
262 typename BoolGridT::Ptr mIntersectingVoxelsGrid;
264 InterruptT *mInterrupter;
281 : mXDist(dist), mYDist(dist), mZDist(dist)
324 void convert(
const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList);
329 void getEdgeData(
Accessor& acc,
const Coord& ijk,
330 std::vector<Vec3d>& points, std::vector<Index32>& primitives);
355 PointTransform(
const std::vector<Vec3s>& pointsIn, std::vector<Vec3s>& pointsOut,
357 : mPointsIn(pointsIn)
358 , mPointsOut(&pointsOut)
363 void run(
bool threaded =
true)
366 tbb::parallel_for(tbb::blocked_range<size_t>(0, mPointsOut->size()), *
this);
368 (*this)(tbb::blocked_range<size_t>(0, mPointsOut->size()));
372 inline void operator()(
const tbb::blocked_range<size_t>& range)
const
374 for (
size_t n = range.begin(); n < range.end(); ++n) {
375 (*mPointsOut)[n] = mXform.worldToIndex(mPointsIn[n]);
380 const std::vector<Vec3s>& mPointsIn;
381 std::vector<Vec3s> *
const mPointsOut;
386 template<
typename ValueType>
389 static ValueType
epsilon() {
return ValueType(1e-7); }
394 template<
typename FloatTreeT,
typename IntTreeT>
396 combine(FloatTreeT& lhsDist, IntTreeT& lhsIndex, FloatTreeT& rhsDist, IntTreeT& rhsIndex)
398 typedef typename FloatTreeT::ValueType FloatValueT;
402 typename FloatTreeT::LeafCIter iter = rhsDist.cbeginLeaf();
404 FloatValueT rhsValue;
407 for ( ; iter; ++iter) {
408 typename FloatTreeT::LeafNodeType::ValueOnCIter it = iter->cbeginValueOn();
413 rhsValue = it.getValue();
414 FloatValueT& lhsValue =
const_cast<FloatValueT&
>(lhsDistAccessor.
getValue(ijk));
416 if (-rhsValue < std::abs(lhsValue)) {
435 template<
typename FloatTreeT,
typename InterruptT = util::NullInterrupter>
445 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
451 const std::vector<Vec4I>& polygonList, InterruptT *interrupter = NULL);
455 void run(
bool threaded =
true);
458 void operator() (
const tbb::blocked_range<size_t> &range);
468 bool wasInterrupted()
const {
return mInterrupter && mInterrupter->wasInterrupted(); }
470 bool evalVoxel(
const Coord& ijk,
const Int32 polyIdx);
472 const std::vector<Vec3s>& mPointList;
473 const std::vector<Vec4I>& mPolygonList;
475 FloatTreeT mSqrDistTree;
476 FloatAccessorT mSqrDistAccessor;
478 IntTreeT mPrimIndexTree;
479 IntAccessorT mPrimIndexAccessor;
481 BoolTreeT mIntersectionTree;
482 BoolAccessorT mIntersectionAccessor;
485 IntTreeT mLastPrimTree;
486 IntAccessorT mLastPrimAccessor;
488 InterruptT *mInterrupter;
491 struct Primitive {
Vec3d a, b, c, d;
Int32 index; };
493 template<
bool IsQuad>
494 bool evalPrimitive(
const Coord&,
const Primitive&);
496 template<
bool IsQuad>
497 void voxelize(
const Primitive&);
501 template<
typename FloatTreeT,
typename InterruptT>
506 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, mPolygonList.size()), *
this);
508 (*this)(tbb::blocked_range<size_t>(0, mPolygonList.size()));
512 template<
typename FloatTreeT,
typename InterruptT>
514 const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList,
515 InterruptT *interrupter)
516 : mPointList(pointList)
517 , mPolygonList(polygonList)
519 , mSqrDistAccessor(mSqrDistTree)
521 , mPrimIndexAccessor(mPrimIndexTree)
522 , mIntersectionTree(false)
523 , mIntersectionAccessor(mIntersectionTree)
525 , mLastPrimAccessor(mLastPrimTree)
526 , mInterrupter(interrupter)
530 template<
typename FloatTreeT,
typename InterruptT>
533 : mPointList(rhs.mPointList)
534 , mPolygonList(rhs.mPolygonList)
536 , mSqrDistAccessor(mSqrDistTree)
538 , mPrimIndexAccessor(mPrimIndexTree)
539 , mIntersectionTree(false)
540 , mIntersectionAccessor(mIntersectionTree)
542 , mLastPrimAccessor(mLastPrimTree)
543 , mInterrupter(rhs.mInterrupter)
548 template<
typename FloatTreeT,
typename InterruptT>
554 for (
size_t n = range.begin(); n < range.end(); ++n) {
556 if (mInterrupter && mInterrupter->wasInterrupted()) {
557 tbb::task::self().cancel_group_execution();
561 const Vec4I& verts = mPolygonList[n];
563 prim.index =
Int32(n);
564 prim.a =
Vec3d(mPointList[verts[0]]);
565 prim.b =
Vec3d(mPointList[verts[1]]);
566 prim.c =
Vec3d(mPointList[verts[2]]);
569 prim.d =
Vec3d(mPointList[verts[3]]);
570 voxelize<true>(prim);
572 voxelize<false>(prim);
578 template<
typename FloatTreeT,
typename InterruptT>
579 template<
bool IsQuad>
583 std::deque<Coord> coordList;
587 coordList.push_back(ijk);
589 evalPrimitive<IsQuad>(ijk, prim);
591 while (!coordList.empty()) {
594 ijk = coordList.back();
595 coordList.pop_back();
597 mIntersectionAccessor.setActiveState(ijk,
true);
599 for (
Int32 i = 0; i < 26; ++i) {
601 if (prim.index != mLastPrimAccessor.getValue(nijk)) {
602 mLastPrimAccessor.setValue(nijk, prim.index);
603 if(evalPrimitive<IsQuad>(nijk, prim)) coordList.push_back(nijk);
610 template<
typename FloatTreeT,
typename InterruptT>
611 template<
bool IsQuad>
613 MeshVoxelizer<FloatTreeT, InterruptT>::evalPrimitive(
const Coord& ijk,
const Primitive& prim)
615 Vec3d uvw, voxelCenter(ijk[0], ijk[1], ijk[2]);
618 FloatValueT dist = FloatValueT((voxelCenter -
623 FloatValueT secondDist = FloatValueT((voxelCenter -
626 if (secondDist < dist) dist = secondDist;
629 FloatValueT oldDist = std::abs(mSqrDistAccessor.getValue(ijk));
631 if (dist < oldDist) {
632 mSqrDistAccessor.setValue(ijk, -dist);
633 mPrimIndexAccessor.setValue(ijk, prim.index);
637 mPrimIndexAccessor.setValue(ijk,
std::min(prim.index, mPrimIndexAccessor.getValue(ijk)));
640 return (dist < 0.86602540378443861);
644 template<
typename FloatTreeT,
typename InterruptT>
648 typedef typename FloatTreeT::RootNodeType FloatRootNodeT;
649 typedef typename FloatRootNodeT::NodeChainType FloatNodeChainT;
650 BOOST_STATIC_ASSERT(boost::mpl::size<FloatNodeChainT>::value > 1);
651 typedef typename boost::mpl::at<FloatNodeChainT, boost::mpl::int_<1> >::type FloatInternalNodeT;
653 typedef typename IntTreeT::RootNodeType IntRootNodeT;
654 typedef typename IntRootNodeT::NodeChainType IntNodeChainT;
655 BOOST_STATIC_ASSERT(boost::mpl::size<IntNodeChainT>::value > 1);
656 typedef typename boost::mpl::at<IntNodeChainT, boost::mpl::int_<1> >::type IntInternalNodeT;
663 rhs.mSqrDistTree.clearAllAccessors();
664 rhs.mPrimIndexTree.clearAllAccessors();
666 typename FloatTreeT::LeafIter leafIt = rhs.mSqrDistTree.beginLeaf();
667 for ( ; leafIt; ++leafIt) {
669 ijk = leafIt->origin();
670 FloatLeafT* lhsDistLeafPt = mSqrDistAccessor.probeLeaf(ijk);
672 if (!lhsDistLeafPt) {
678 mSqrDistAccessor.addLeaf(rhs.mSqrDistAccessor.probeLeaf(ijk));
679 FloatInternalNodeT* floatNode =
680 rhs.mSqrDistAccessor.template getNode<FloatInternalNodeT>();
681 floatNode->template stealNode<FloatLeafT>(ijk, background,
false);
683 mPrimIndexAccessor.addLeaf(rhs.mPrimIndexAccessor.probeLeaf(ijk));
684 IntInternalNodeT* intNode =
685 rhs.mPrimIndexAccessor.template getNode<IntInternalNodeT>();
690 IntLeafT* lhsIdxLeafPt = mPrimIndexAccessor.probeLeaf(ijk);
691 IntLeafT* rhsIdxLeafPt = rhs.mPrimIndexAccessor.probeLeaf(ijk);
694 typename FloatLeafT::ValueOnCIter it = leafIt->cbeginValueOn();
699 lhsValue = std::abs(lhsDistLeafPt->getValue(offset));
700 rhsValue = std::abs(it.getValue());
702 if (rhsValue < lhsValue) {
703 lhsDistLeafPt->setValueOn(offset, it.getValue());
704 lhsIdxLeafPt->setValueOn(offset, rhsIdxLeafPt->getValue(offset));
706 lhsIdxLeafPt->setValueOn(offset,
707 std::min(lhsIdxLeafPt->getValue(offset), rhsIdxLeafPt->getValue(offset)));
713 mIntersectionTree.merge(rhs.mIntersectionTree);
715 rhs.mSqrDistTree.clear();
716 rhs.mPrimIndexTree.clear();
717 rhs.mIntersectionTree.clear();
727 template<
typename FloatTreeT,
typename InterruptT = util::NullInterrupter>
733 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
739 void run(
bool threaded =
true);
742 void operator()(
const tbb::blocked_range<int> &range)
const;
747 int sparseScan(
int slice)
const;
749 FloatTreeT& mDistTree;
758 std::vector<Index> mStepSize;
760 InterruptT *mInterrupter;
764 template<
typename FloatTreeT,
typename InterruptT>
769 tbb::parallel_for(tbb::blocked_range<int>(mBBox.min()[0], mBBox.max()[0]+1), *
this);
771 (*this)(tbb::blocked_range<int>(mBBox.min()[0], mBBox.max()[0]+1));
776 template<
typename FloatTreeT,
typename InterruptT>
778 FloatTreeT& distTree,
const BoolTreeT& intersectionTree, InterruptT *interrupter)
779 : mDistTree(distTree)
780 , mDistAccessor(mDistTree)
781 , mIntersectionTree(intersectionTree)
782 , mIntersectionAccessor(mIntersectionTree)
785 , mInterrupter(interrupter)
788 std::vector<Index> dims;
789 mDistTree.getNodeLog2Dims(dims);
791 mStepSize.resize(dims.size()+1, 1);
793 for (
int idx = static_cast<int>(dims.size()) - 1; idx > -1; --idx) {
794 exponent += dims[idx];
795 mStepSize[idx] = 1 << exponent;
798 mDistTree.evalLeafBoundingBox(mBBox);
801 const int tileDim = mStepSize[0];
803 for (
size_t i = 0; i < 3; ++i) {
806 double diff = std::abs(
double(mBBox.min()[i])) / double(tileDim);
808 if (mBBox.min()[i] <= tileDim) {
809 n = int(std::ceil(diff));
810 mBBox.min()[i] = - n * tileDim;
812 n = int(std::floor(diff));
813 mBBox.min()[i] = n * tileDim;
816 n = int(std::ceil(std::abs(
double(mBBox.max()[i] - mBBox.min()[i])) / double(tileDim)));
817 mBBox.max()[i] = mBBox.min()[i] + n * tileDim;
822 template<
typename FloatTreeT,
typename InterruptT>
825 : mDistTree(rhs.mDistTree)
826 , mDistAccessor(mDistTree)
827 , mIntersectionTree(rhs.mIntersectionTree)
828 , mIntersectionAccessor(mIntersectionTree)
830 , mStepSize(rhs.mStepSize)
831 , mInterrupter(rhs.mInterrupter)
836 template<
typename FloatTreeT,
typename InterruptT>
842 for (
int n = range.begin(); n < range.end(); n += iStep) {
844 if (mInterrupter && mInterrupter->wasInterrupted()) {
845 tbb::task::self().cancel_group_execution();
849 iStep = sparseScan(n);
854 template<
typename FloatTreeT,
typename InterruptT>
858 bool lastVoxelWasOut =
true;
859 int last_k = mBBox.min()[2];
861 Coord ijk(slice, mBBox.min()[1], mBBox.min()[2]);
862 Coord step(mStepSize[mDistAccessor.getValueDepth(ijk) + 1]);
865 for (ijk[1] = mBBox.min()[1]; ijk[1] <= mBBox.max()[1]; ijk[1] += step[1]) {
867 if (mInterrupter && mInterrupter->wasInterrupted()) {
871 step[1] = mStepSize[mDistAccessor.getValueDepth(ijk) + 1];
872 step[0] =
std::min(step[0], step[1]);
874 for (ijk[2] = mBBox.min()[2]; ijk[2] <= mBBox.max()[2]; ijk[2] += step[2]) {
876 step[2] = mStepSize[mDistAccessor.getValueDepth(ijk) + 1];
877 step[1] =
std::min(step[1], step[2]);
878 step[0] =
std::min(step[0], step[2]);
881 if (mDistAccessor.isValueOn(ijk)) {
884 if (mIntersectionAccessor.isValueOn(ijk)) {
886 lastVoxelWasOut =
false;
889 }
else if (lastVoxelWasOut) {
891 FloatValueT& val =
const_cast<FloatValueT&
>(mDistAccessor.getValue(ijk));
897 for (
Int32 n = 3; n < 6; n += 2) {
898 n_ijk = ijk + util::COORD_OFFSETS[n];
900 if (mDistAccessor.probeValue(n_ijk, val) && val > 0) {
901 lastVoxelWasOut =
true;
906 if (lastVoxelWasOut) {
908 FloatValueT& v =
const_cast<FloatValueT&
>(mDistAccessor.getValue(ijk));
911 const int tmp_k = ijk[2];
914 for (--ijk[2]; ijk[2] >= last_k; --ijk[2]) {
915 if (mIntersectionAccessor.isValueOn(ijk))
break;
917 const_cast<FloatValueT&
>(mDistAccessor.getValue(ijk));
918 if (vb < FloatValueT(0.0)) vb = -vb;
941 template<
typename FloatTreeT,
typename InterruptT = util::NullInterrupter>
949 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
956 InterruptT *interrupter = NULL);
960 void run(
bool threaded =
true);
963 void operator() (
const tbb::blocked_range<size_t> &range);
971 bool wasInterrupted()
const {
return mInterrupter && mInterrupter->wasInterrupted(); }
974 const FloatTreeT& mDistTree;
979 InterruptT *mInterrupter;
983 template<
typename FloatTreeT,
typename InterruptT>
986 const BoolTreeT& intersectionTree, InterruptT *interrupter)
987 : mDistLeafs(distLeafs)
988 , mDistTree(distTree)
989 , mIntersectionTree(intersectionTree)
990 , mSignMaskTree(false)
991 , mInterrupter(interrupter)
996 template<
typename FloatTreeT,
typename InterruptT>
999 : mDistLeafs(rhs.mDistLeafs)
1000 , mDistTree(rhs.mDistTree)
1001 , mIntersectionTree(rhs.mIntersectionTree)
1002 , mSignMaskTree(false)
1003 , mInterrupter(rhs.mInterrupter)
1008 template<
typename FloatTreeT,
typename InterruptT>
1012 if (threaded) tbb::parallel_reduce(mDistLeafs.getRange(), *
this);
1013 else (*
this)(mDistLeafs.getRange());
1017 template<
typename FloatTreeT,
typename InterruptT>
1027 Coord& maxCoord = bbox.max();
1028 Coord& minCoord = bbox.min();
1030 const int extent = BoolLeafT::DIM - 1;
1032 for (
size_t n = range.begin(); n < range.end(); ++n) {
1034 const FloatLeafT& distLeaf = mDistLeafs.leaf(n);
1036 minCoord = distLeaf.origin();
1037 maxCoord[0] = minCoord[0] + extent;
1038 maxCoord[1] = minCoord[1] + extent;
1039 maxCoord[2] = minCoord[2] + extent;
1045 bool addLeaf =
false;
1049 typename FloatLeafT::ValueOnCIter it = distLeaf.cbeginValueOn();
1051 if (intersectionLeaf && intersectionLeaf->isValueOn(it.pos()))
continue;
1053 ijk = it.getCoord();
1054 if (bbox.isInside(ijk)) {
1055 for (
size_t i = 0; i < 6; ++i) {
1056 if (distLeaf.probeValue(ijk+util::COORD_OFFSETS[i], value) && value>0.0) {
1057 maskLeaf.setValueOn(ijk);
1063 for (
size_t i = 0; i < 6; ++i) {
1064 if (distAcc.
probeValue(ijk+util::COORD_OFFSETS[i], value) && value>0.0) {
1065 maskLeaf.setValueOn(ijk);
1074 if (addLeaf) maskAcc.
addLeaf(maskLeafPt);
1075 else delete maskLeafPt;
1080 template<
typename FloatTreeT,
typename InterruptT>
1084 mSignMaskTree.merge(rhs.mSignMaskTree);
1092 template<
typename FloatTreeT,
typename InterruptT = util::NullInterrupter>
1099 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1109 void run(
bool threaded =
true);
1112 void operator() (
const tbb::blocked_range<size_t> &range);
1120 bool wasInterrupted()
const {
return mInterrupter && mInterrupter->wasInterrupted(); }
1123 FloatTreeT& mDistTree;
1127 InterruptT *mInterrupter;
1131 template<
typename FloatTreeT,
typename InterruptT>
1133 FloatTreeT& distTree,
const BoolTreeT& intersectionTree, InterruptT *interrupter)
1134 : mOldSignMaskLeafs(signMaskLeafs)
1135 , mDistTree(distTree)
1136 , mIntersectionTree(intersectionTree)
1137 , mSignMaskTree(false)
1138 , mInterrupter(interrupter)
1143 template<
typename FloatTreeT,
typename InterruptT>
1146 : mOldSignMaskLeafs(rhs.mOldSignMaskLeafs)
1147 , mDistTree(rhs.mDistTree)
1148 , mIntersectionTree(rhs.mIntersectionTree)
1149 , mSignMaskTree(false)
1150 , mInterrupter(rhs.mInterrupter)
1155 template<
typename FloatTreeT,
typename InterruptT>
1159 if (threaded) tbb::parallel_reduce(mOldSignMaskLeafs.getRange(), *
this);
1160 else (*
this)(mOldSignMaskLeafs.getRange());
1164 template<
typename FloatTreeT,
typename InterruptT>
1172 std::deque<Coord> coordList;
1176 Coord& maxCoord = bbox.max();
1177 Coord& minCoord = bbox.min();
1179 const int extent = BoolLeafT::DIM - 1;
1181 for (
size_t n = range.begin(); n < range.end(); ++n) {
1182 BoolLeafT& oldMaskLeaf = mOldSignMaskLeafs.leaf(n);
1184 minCoord = oldMaskLeaf.origin();
1185 maxCoord[0] = minCoord[0] + extent;
1186 maxCoord[1] = minCoord[1] + extent;
1187 maxCoord[2] = minCoord[2] + extent;
1192 typename BoolLeafT::ValueOnCIter it = oldMaskLeaf.cbeginValueOn();
1194 coordList.push_back(it.getCoord());
1196 while (!coordList.empty()) {
1198 ijk = coordList.back();
1199 coordList.pop_back();
1205 for (
size_t i = 0; i < 6; ++i) {
1206 nijk = ijk + util::COORD_OFFSETS[i];
1207 if (bbox.isInside(nijk)) {
1208 if (intersectionLeaf && intersectionLeaf->isValueOn(nijk))
continue;
1210 if (distLeaf.probeValue(nijk, value) && value < 0.0) {
1211 coordList.push_back(nijk);
1216 distAcc.
probeValue(nijk, value) && value < 0.0) {
1228 template<
typename FloatTreeT,
typename InterruptT>
1232 mSignMaskTree.merge(rhs.mSignMaskTree);
1243 template<
typename FloatTreeT>
1251 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1256 const std::vector<Vec3s>& pointList,
1257 const std::vector<Vec4I>& polygonList,
1258 FloatTreeT& distTree,
1265 void run(
bool threaded =
true);
1268 void operator()(
const tbb::blocked_range<size_t>&)
const;
1273 Vec3d getClosestPoint(
const Coord& ijk,
const Vec4I& prim)
const;
1275 std::vector<Vec3s>
const *
const mPointList;
1276 std::vector<Vec4I>
const *
const mPolygonList;
1278 FloatTreeT& mDistTree;
1286 template<
typename FloatTreeT>
1290 if (threaded) tbb::parallel_for(mLeafs.getRange(), *
this);
1291 else (*
this)(mLeafs.getRange());
1295 template<
typename FloatTreeT>
1297 const std::vector<Vec3s>& pointList,
1298 const std::vector<Vec4I>& polygonList,
1299 FloatTreeT& distTree,
1303 : mPointList(&pointList)
1304 , mPolygonList(&polygonList)
1305 , mDistTree(distTree)
1306 , mIndexTree(indexTree)
1307 , mIntersectionTree(intersectionTree)
1313 template<
typename FloatTreeT>
1316 : mPointList(rhs.mPointList)
1317 , mPolygonList(rhs.mPolygonList)
1318 , mDistTree(rhs.mDistTree)
1319 , mIndexTree(rhs.mIndexTree)
1320 , mIntersectionTree(rhs.mIntersectionTree)
1321 , mLeafs(rhs.mLeafs)
1326 template<
typename FloatTreeT>
1329 const tbb::blocked_range<size_t>& range)
const
1340 typename BoolTreeT::LeafNodeType::ValueOnCIter iter;
1341 for (
size_t n = range.begin(); n < range.end(); ++n) {
1342 iter = mLeafs.leaf(n).cbeginValueOn();
1343 for (; iter; ++iter) {
1345 ijk = iter.getCoord();
1351 center =
Vec3d(ijk[0], ijk[1], ijk[2]);
1353 for (
Int32 i = 0; i < 26; ++i) {
1354 nijk = ijk + util::COORD_OFFSETS[i];
1361 cpt = getClosestPoint(nijk, prim);
1363 dir1 = center -
cpt;
1366 dir2 =
Vec3d(nijk[0], nijk[1], nijk[2]) -
cpt;
1369 if (dir2.dot(dir1) > 0.0) {
1380 template<
typename FloatTreeT>
1384 Vec3d voxelCenter(ijk[0], ijk[1], ijk[2]);
1387 const Vec3d a((*mPointList)[prim[0]]);
1388 const Vec3d b((*mPointList)[prim[1]]);
1389 const Vec3d c((*mPointList)[prim[2]]);
1397 Vec3d diff1 = voxelCenter - cpt1;
1399 const Vec3d d((*mPointList)[prim[3]]);
1402 Vec3d diff2 = voxelCenter - cpt2;
1404 if (diff2.lengthSqr() < diff1.lengthSqr()) {
1419 template<
typename FloatTreeT>
1429 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1439 void run(
bool threaded =
true);
1442 void operator()(
const tbb::blocked_range<size_t>&)
const;
1447 FloatTreeT& mDistTree;
1454 template<
typename FloatTreeT>
1458 if (threaded) tbb::parallel_for(mLeafs.getRange(), *
this);
1459 else (*
this)(mLeafs.getRange());
1461 mIntersectionTree.pruneInactive();
1465 template<
typename FloatTreeT>
1467 FloatTreeT& distTree,
1471 : mDistTree(distTree)
1472 , mIndexTree(indexTree)
1473 , mIntersectionTree(intersectionTree)
1479 template<
typename FloatTreeT>
1482 : mDistTree(rhs.mDistTree)
1483 , mIndexTree(rhs.mIndexTree)
1484 , mIntersectionTree(rhs.mIntersectionTree)
1485 , mLeafs(rhs.mLeafs)
1490 template<
typename FloatTreeT>
1493 const tbb::blocked_range<size_t>& range)
const
1500 typename BoolLeafT::ValueOnCIter iter;
1506 for (
size_t n = range.begin(); n < range.end(); ++n) {
1510 ijk = maskLeaf.origin();
1514 iter = maskLeaf.cbeginValueOn();
1515 for (; iter; ++iter) {
1517 offset = iter.pos();
1519 if(distLeaf->getValue(offset) > 0.0)
continue;
1521 ijk = iter.getCoord();
1523 for (
Int32 m = 0; m < 26; ++m) {
1524 m_ijk = ijk + util::COORD_OFFSETS[m];
1534 maskLeaf.setValueOff(offset);
1535 distLeaf->setValueOn(offset, -0.86602540378443861);
1549 template<
typename FloatTreeT>
1560 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1569 void run(
bool threaded =
true);
1572 void operator()(
const tbb::blocked_range<size_t>&)
const;
1577 FloatTreeT& mDistTree;
1584 template<
typename FloatTreeT>
1588 if (threaded) tbb::parallel_for(mLeafs.getRange(), *
this);
1589 else (*
this)(mLeafs.getRange());
1591 mDistTree.pruneInactive();
1592 mIndexTree.pruneInactive();
1596 template<
typename FloatTreeT>
1598 FloatTreeT& distTree,
1602 : mDistTree(distTree)
1604 , mIndexTree(indexTree)
1605 , mIntersectionTree(intersectionTree)
1610 template<
typename FloatTreeT>
1613 : mDistTree(rhs.mDistTree)
1614 , mLeafs(rhs.mLeafs)
1615 , mIndexTree(rhs.mIndexTree)
1616 , mIntersectionTree(rhs.mIntersectionTree)
1621 template<
typename FloatTreeT>
1624 const tbb::blocked_range<size_t>& range)
const
1631 typename DistLeafT::ValueOnCIter iter;
1632 const FloatValueT distBG = mDistTree.background();
1633 const Int32 indexBG = mIntersectionTree.background();
1640 for (
size_t n = range.begin(); n < range.end(); ++n) {
1644 ijk = distLeaf.origin();
1649 iter = distLeaf.cbeginValueOn();
1650 for (; iter; ++iter) {
1652 value = iter.getValue();
1653 if(value > 0.0)
continue;
1655 offset = iter.pos();
1656 if (maskLeaf && maskLeaf->isValueOn(offset))
continue;
1658 ijk = iter.getCoord();
1660 for (
Int32 m = 0; m < 26; ++m) {
1661 m_ijk = ijk + util::COORD_OFFSETS[m];
1669 distLeaf.setValueOff(offset, distBG);
1670 indexLeaf.setValueOff(offset, indexBG);
1680 template<
typename TreeType>
1687 template <
typename LeafNodeType>
1690 LeafNodeType* rhsLeaf =
const_cast<LeafNodeType*
>(mAcc.
probeLeaf(leaf.origin()));
1691 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
1692 for (; iter; ++iter) {
1693 rhsLeaf->setValueOnly(iter.pos(), iter.getValue());
1705 template<
typename FloatTreeT>
1715 typedef typename FloatTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1723 const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList);
1725 void run(
bool threaded =
true);
1727 void operator()(
const tbb::blocked_range<size_t>&);
1741 double closestPrimDist(
const Coord&, std::vector<Int32>&,
Int32&)
const;
1745 FloatTreeT& mDistTree;
1749 const FloatValueT mExteriorBandWidth, mInteriorBandWidth, mVoxelSize;
1750 const std::vector<Vec3s>& mPointList;
1751 const std::vector<Vec4I>& mPolygonList;
1753 FloatTreeT mNewDistTree;
1759 template<
typename FloatTreeT>
1762 FloatTreeT& distTree,
1768 const std::vector<Vec3s>& pointList,
1769 const std::vector<Vec4I>& polygonList)
1771 , mDistTree(distTree)
1772 , mIndexTree(indexTree)
1773 , mMaskTree(maskTree)
1774 , mExteriorBandWidth(exteriorBandWidth)
1775 , mInteriorBandWidth(interiorBandWidth)
1776 , mVoxelSize(voxelSize)
1777 , mPointList(pointList)
1778 , mPolygonList(polygonList)
1781 , mNewMaskTree(false)
1786 template<
typename FloatTreeT>
1788 : mMaskLeafs(rhs.mMaskLeafs)
1789 , mDistTree(rhs.mDistTree)
1790 , mIndexTree(rhs.mIndexTree)
1791 , mMaskTree(rhs.mMaskTree)
1792 , mExteriorBandWidth(rhs.mExteriorBandWidth)
1793 , mInteriorBandWidth(rhs.mInteriorBandWidth)
1794 , mVoxelSize(rhs.mVoxelSize)
1795 , mPointList(rhs.mPointList)
1796 , mPolygonList(rhs.mPolygonList)
1799 , mNewMaskTree(false)
1804 template<
typename FloatTreeT>
1808 if (threaded) tbb::parallel_reduce(mMaskLeafs.getRange(), *
this);
1809 else (*
this)(mMaskLeafs.getRange());
1813 mDistTree.topologyUnion(mNewDistTree);
1817 mIndexTree.merge(mNewIndexTree);
1820 mMaskTree.merge(mNewMaskTree);
1825 template<
typename FloatTreeT>
1830 Int32 closestPrim = 0;
1844 std::vector<Int32> primitives(18);
1846 for (
size_t n = range.begin(); n < range.end(); ++n) {
1848 BoolLeafT& maskLeaf = mMaskLeafs.leaf(n);
1850 if (maskLeaf.isEmpty())
continue;
1852 ijk = maskLeaf.origin();
1858 newDistAcc.
addLeaf(distLeafPt);
1862 if (!indexLeafPt) indexLeafPt = newIndexAcc.
touchLeaf(ijk);
1864 bbox = maskLeaf.getNodeBoundingBox();
1867 typename BoolLeafT::ValueOnIter iter = maskLeaf.beginValueOn();
1868 for (; iter; ++iter) {
1870 ijk = iter.getCoord();
1872 if (bbox.isInside(ijk)) {
1873 distance = evalVoxelDist(ijk, *distLeafPt, *indexLeafPt, maskLeaf,
1874 primitives, closestPrim);
1876 distance = evalVoxelDist(ijk, distAcc, indexAcc, maskAcc,
1877 primitives, closestPrim);
1882 inside = distLeafPt->getValue(pos) <
FloatValueT(0.0);
1884 if (!inside && distance < mExteriorBandWidth) {
1885 distLeafPt->setValueOn(pos, distance);
1886 indexLeafPt->setValueOn(pos, closestPrim);
1887 }
else if (inside && distance < mInteriorBandWidth) {
1888 distLeafPt->setValueOn(pos, -distance);
1889 indexLeafPt->setValueOn(pos, closestPrim);
1894 for (
Int32 i = 0; i < 6; ++i) {
1895 newMaskAcc.
setValueOn(ijk + util::COORD_OFFSETS[i]);
1902 template<
typename FloatTreeT>
1906 FloatAccessorT& distAcc,
1907 IntAccessorT& indexAcc,
1908 BoolAccessorT& maskAcc,
1909 std::vector<Int32>& prims,
1910 Int32& closestPrim)
const
1917 for (
Int32 n = 0; n < 18; ++n) {
1918 n_ijk = ijk + util::COORD_OFFSETS[n];
1919 if (!maskAcc.isValueOn(n_ijk) && distAcc.probeValue(n_ijk, tmpDist)) {
1920 prims.push_back(indexAcc.getValue(n_ijk));
1921 tmpDist = std::abs(tmpDist);
1922 if (tmpDist < minDist) minDist = tmpDist;
1927 tmpDist = FloatValueT(closestPrimDist(ijk, prims, closestPrim));
1931 return tmpDist > minDist ? tmpDist : minDist + mVoxelSize;
1936 template<
typename FloatTreeT>
1938 ExpandNB<FloatTreeT>::evalVoxelDist(
1940 FloatLeafT& distLeaf,
1941 IntLeafT& indexLeaf,
1942 BoolLeafT& maskLeaf,
1943 std::vector<Int32>& prims,
1944 Int32& closestPrim)
const
1950 for (
Int32 n = 0; n < 18; ++n) {
1951 pos = FloatLeafT::coordToOffset(ijk + util::COORD_OFFSETS[n]);
1952 if (!maskLeaf.isValueOn(pos) && distLeaf.probeValue(pos, tmpDist)) {
1953 prims.push_back(indexLeaf.getValue(pos));
1954 tmpDist = std::abs(tmpDist);
1955 if (tmpDist < minDist) minDist = tmpDist;
1959 tmpDist = FloatValueT(closestPrimDist(ijk, prims, closestPrim));
1960 return tmpDist > minDist ? tmpDist : minDist + mVoxelSize;
1964 template<
typename FloatTreeT>
1966 ExpandNB<FloatTreeT>::closestPrimDist(
const Coord& ijk,
1967 std::vector<Int32>& prims,
Int32& closestPrim)
const
1969 std::sort(prims.begin(), prims.end());
1971 Int32 lastPrim = -1;
1972 Vec3d uvw, voxelCenter(ijk[0], ijk[1], ijk[2]);
1975 for (
size_t n = 0, N = prims.size(); n < N; ++n) {
1976 if (prims[n] == lastPrim)
continue;
1978 lastPrim = prims[n];
1980 const Vec4I& verts = mPolygonList[lastPrim];
1983 const Vec3d a(mPointList[verts[0]]);
1984 const Vec3d b(mPointList[verts[1]]);
1985 const Vec3d c(mPointList[verts[2]]);
1987 primDist = (voxelCenter -
1992 const Vec3d d(mPointList[verts[3]]);
1994 tmpDist = (voxelCenter -
1997 if (tmpDist < primDist) primDist = tmpDist;
2000 if (primDist < dist) {
2002 closestPrim = lastPrim;
2006 return std::sqrt(dist) * double(mVoxelSize);
2010 template<
typename FloatTreeT>
2014 mNewDistTree.merge(rhs.mNewDistTree);
2015 mNewIndexTree.merge(rhs.mNewIndexTree);
2016 mNewMaskTree.merge(rhs.mNewMaskTree);
2023 template<
typename ValueType>
2027 : mVoxelSize(voxelSize)
2028 , mUnsigned(unsignedDist)
2032 template <
typename LeafNodeType>
2039 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2040 for (; iter; ++iter) {
2041 ValueType& val =
const_cast<ValueType&
>(iter.getValue());
2042 val = w[!mUnsigned && int(val < ValueType(0.0))] * std::sqrt(std::abs(val));
2047 ValueType mVoxelSize;
2048 const bool mUnsigned;
2052 template<
typename ValueType>
2056 : mExBandWidth(exBandWidth)
2057 , mInBandWidth(inBandWidth)
2061 template <
typename LeafNodeType>
2064 ValueType bgValues[2];
2065 bgValues[0] = mExBandWidth;
2066 bgValues[1] = -mInBandWidth;
2068 typename LeafNodeType::ValueOffIter iter = leaf.beginValueOff();
2070 for (; iter; ++iter) {
2071 ValueType& val =
const_cast<ValueType&
>(iter.getValue());
2072 val = bgValues[int(val < ValueType(0.0))];
2077 ValueType mExBandWidth, mInBandWidth;
2081 template<
typename ValueType>
2084 TrimOp(ValueType exBandWidth, ValueType inBandWidth)
2085 : mExBandWidth(exBandWidth)
2086 , mInBandWidth(inBandWidth)
2090 template <
typename LeafNodeType>
2093 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2095 for (; iter; ++iter) {
2096 ValueType& val =
const_cast<ValueType&
>(iter.getValue());
2097 const bool inside = val < ValueType(0.0);
2099 if (inside && !(val > -mInBandWidth)) {
2100 val = -mInBandWidth;
2102 }
else if (!inside && !(val < mExBandWidth)) {
2110 ValueType mExBandWidth, mInBandWidth;
2114 template<
typename ValueType>
2121 template <
typename LeafNodeType>
2124 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2125 for (; iter; ++iter) {
2126 ValueType& val =
const_cast<ValueType&
>(iter.getValue());
2136 template<
typename Gr
idType,
typename ValueType>
2140 typedef typename Scheme::template ISStencil<GridType>::StencilType
Stencil;
2147 , mVoxelSize(voxelSize)
2154 template <
typename LeafNodeType>
2157 const ValueType dt = mCFL * mVoxelSize, one(1.0), invDx = one / mVoxelSize;
2161 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2162 for (; iter; ++iter) {
2163 stencil.moveTo(iter);
2165 const ValueType normSqGradPhi =
2168 const ValueType phi0 = iter.getValue();
2169 const ValueType diff =
math::Sqrt(normSqGradPhi) * invDx - one;
2172 buffer.setValue(iter.pos(), phi0 - dt * S * diff);
2179 ValueType mVoxelSize, mCFL;
2183 template<
typename TreeType,
typename ValueType>
2191 template <
typename LeafNodeType>
2195 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2197 for (; iter; ++iter) {
2198 ValueType& val =
const_cast<ValueType&
>(iter.getValue());
2199 val =
std::min(val, buffer.getValue(iter.pos()));
2208 template<
typename TreeType,
typename ValueType>
2216 , mBufferIndex(bufferIndex)
2220 template <
typename LeafNodeType>
2224 typename LeafNodeType::ValueOnIter iter = leaf.beginValueOn();
2227 for (; iter; ++iter) {
2228 offset = iter.pos();
2229 leaf.setValueOnly(offset, buffer.getValue(offset));
2235 const size_t mBufferIndex;
2239 template<
typename TreeType>
2247 template <
typename LeafNodeType>
2251 if (rhsLeaf) leaf.topologyDifference(*rhsLeaf,
false);
2259 template<
typename ValueType>
2268 template <
typename ChildType>
2271 child.pruneOp(*
this);
2272 if (!child.isInactive())
return false;
2291 template<
typename FloatGr
idT,
typename InterruptT>
2293 openvdb::math::Transform::Ptr& transform,
int conversionFlags,
2294 InterruptT *interrupter,
int signSweeps)
2295 : mTransform(transform)
2296 , mConversionFlags(conversionFlags)
2297 , mSignSweeps(signSweeps)
2298 , mInterrupter(interrupter)
2301 mSignSweeps =
std::min(mSignSweeps, 1);
2305 template<
typename FloatGr
idT,
typename InterruptT>
2311 mIntersectingVoxelsGrid = BoolGridT::create(
false);
2315 template<
typename FloatGr
idT,
typename InterruptT>
2318 const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList,
2324 const FloatValueT vs = mTransform->voxelSize()[0];
2325 doConvert(pointList, polygonList, vs * exBandWidth, vs * inBandWidth);
2330 template<
typename FloatGr
idT,
typename InterruptT>
2333 const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList,
2338 const FloatValueT vs = mTransform->voxelSize()[0];
2339 doConvert(pointList, polygonList, vs * exBandWidth, 0.0,
true);
2344 template<
typename FloatGr
idT,
typename InterruptT>
2347 const std::vector<Vec3s>& pointList,
const std::vector<Vec4I>& polygonList,
2348 FloatValueT exBandWidth, FloatValueT inBandWidth,
bool unsignedDistField)
2350 mDistGrid->setTransform(mTransform);
2351 mIndexGrid->setTransform(mTransform);
2363 voxelizer(pointList, polygonList, mInterrupter);
2369 mDistGrid->tree().merge(voxelizer.sqrDistTree());
2370 mIndexGrid->tree().merge(voxelizer.primIndexTree());
2371 mIntersectingVoxelsGrid->tree().merge(voxelizer.intersectionTree());
2374 if (!unsignedDistField) {
2378 internal::ContourTracer<FloatTreeT, InterruptT> trace(
2379 mDistGrid->tree(), mIntersectingVoxelsGrid->tree(), mInterrupter);
2380 for (
int i = 0; i < mSignSweeps; ++i) {
2389 BoolTreeT signMaskTree(
false);
2392 internal::SignMask<FloatTreeT, InterruptT> signMaskOp(leafs,
2393 mDistGrid->tree(), mIntersectingVoxelsGrid->tree(), mInterrupter);
2395 signMaskTree.merge(signMaskOp.signMaskTree());
2402 if(leafs.leafCount() == 0)
break;
2404 internal::PropagateSign<FloatTreeT, InterruptT> sign(leafs,
2405 mDistGrid->tree(), mIntersectingVoxelsGrid->tree(), mInterrupter);
2409 signMaskTree.clear();
2410 signMaskTree.merge(sign.signMaskTree());
2418 tree::LeafManager<BoolTreeT> leafs(mIntersectingVoxelsGrid->tree());
2421 internal::IntersectingVoxelSign<FloatTreeT> sign(pointList, polygonList,
2422 mDistGrid->tree(), mIndexGrid->tree(), mIntersectingVoxelsGrid->tree(), leafs);
2430 internal::IntersectingVoxelCleaner<FloatTreeT> cleaner(mDistGrid->tree(),
2431 mIndexGrid->tree(), mIntersectingVoxelsGrid->tree(), leafs);
2438 tree::LeafManager<FloatTreeT> leafs(mDistGrid->tree());
2440 internal::ShellVoxelCleaner<FloatTreeT> cleaner(mDistGrid->tree(),
2441 leafs, mIndexGrid->tree(), mIntersectingVoxelsGrid->tree());
2449 inBandWidth = FloatValueT(0.0);
2452 if (mDistGrid->activeVoxelCount() == 0)
return;
2454 mIntersectingVoxelsGrid->clear();
2455 const FloatValueT voxelSize(mTransform->voxelSize()[0]);
2458 tree::LeafManager<FloatTreeT> leafs(mDistGrid->tree());
2459 leafs.foreach(internal::SqrtAndScaleOp<FloatValueT>(voxelSize, unsignedDistField));
2464 if (!unsignedDistField) {
2465 mDistGrid->tree().root().setBackground(exBandWidth,
false);
2466 mDistGrid->tree().signedFloodFill(exBandWidth, -inBandWidth);
2472 const FloatValueT minWidth = voxelSize * 2.0;
2473 if (inBandWidth > minWidth || exBandWidth > minWidth) {
2476 BoolTreeT maskTree(
false);
2477 maskTree.topologyUnion(mDistGrid->tree());
2481 internal::LeafTopologyDiffOp<FloatTreeT> diffOp(mDistGrid->tree());
2485 float progress = 48, step = 0.0;
2488 2.0 * std::ceil((
std::max(inBandWidth, exBandWidth) - minWidth) / voxelSize);
2489 if (estimated <
double(maxIterations)) {
2490 maxIterations = unsigned(estimated);
2491 step = 42.0 / float(maxIterations);
2499 tree::LeafManager<BoolTreeT> leafs(maskTree);
2501 if (leafs.leafCount() == 0)
break;
2503 leafs.foreach(diffOp);
2505 internal::ExpandNB<FloatTreeT> expand(
2506 leafs, mDistGrid->tree(), mIndexGrid->tree(), maskTree,
2507 exBandWidth, inBandWidth, voxelSize, pointList, polygonList);
2511 if ((++count) >= maxIterations)
break;
2522 if (!unsignedDistField && !rawData) {
2524 mDistGrid->tree().pruneLevelSet();
2525 tree::LeafManager<FloatTreeT> leafs(mDistGrid->tree(), 1);
2527 const FloatValueT offset = 0.8 * voxelSize;
2530 internal::OffsetOp<FloatValueT> offsetOp(-offset);
2532 leafs.foreach(offsetOp);
2536 leafs.foreach(internal::RenormOp<FloatGridT, FloatValueT>(*mDistGrid, leafs, voxelSize));
2538 leafs.foreach(internal::MinOp<FloatTreeT, FloatValueT>(leafs));
2542 offsetOp.resetOffset(offset - internal::Tolerance<FloatValueT>::epsilon());
2543 leafs.foreach(offsetOp);
2548 const FloatValueT minTrimWidth = voxelSize * 4.0;
2549 if (inBandWidth < minTrimWidth || exBandWidth < minTrimWidth) {
2555 tree::LeafManager<FloatTreeT> leafs(mDistGrid->tree());
2556 leafs.foreach(internal::TrimOp<FloatValueT>(
2557 exBandWidth, unsignedDistField ? exBandWidth : inBandWidth));
2559 internal::SDFPrune<FloatValueT> sdfPrune(exBandWidth, -inBandWidth);
2560 mDistGrid->tree().pruneOp(sdfPrune);
2569 template<
typename Gr
idType>
2570 inline typename boost::enable_if<boost::is_floating_point<typename GridType::ValueType>,
2571 typename GridType::Ptr>::type
2573 const openvdb::math::Transform& xform,
2574 const std::vector<Vec3s>& points,
2575 const std::vector<Vec3I>& triangles,
2576 const std::vector<Vec4I>& quads,
2579 bool unsignedDistanceField =
false)
2581 std::vector<Vec3s> indexSpacePoints(points.size());
2589 std::vector<Vec4I> primitives(triangles.size() + quads.size());
2591 for (
size_t n = 0, N = triangles.size(); n < N; ++n) {
2592 Vec4I& prim = primitives[n];
2593 const Vec3I& triangle = triangles[n];
2594 prim[0] = triangle[0];
2595 prim[1] = triangle[1];
2596 prim[2] = triangle[2];
2600 for (
size_t n = 0, N = quads.size(); n < N; ++n) {
2601 primitives[n + triangles.size()] = quads[n];
2604 typename GridType::ValueType exWidth(exBandWidth);
2605 typename GridType::ValueType inWidth(inBandWidth);
2611 if (!unsignedDistanceField) {
2623 template<
typename Gr
idType>
2624 inline typename boost::disable_if<boost::is_floating_point<typename GridType::ValueType>,
2625 typename GridType::Ptr>::type
2628 const std::vector<Vec3s>& ,
2629 const std::vector<Vec3I>& ,
2630 const std::vector<Vec4I>& ,
2633 bool unsignedDistanceField =
false)
2636 "mesh to volume conversion is supported only for scalar, floating-point grids");
2643 template<
typename Gr
idType>
2644 inline typename GridType::Ptr
2646 const openvdb::math::Transform& xform,
2647 const std::vector<Vec3s>& points,
2648 const std::vector<Vec3I>& triangles,
2651 std::vector<Vec4I> quads(0);
2652 return doMeshConversion<GridType>(xform, points, triangles, quads,
2653 halfWidth, halfWidth);
2657 template<
typename Gr
idType>
2658 inline typename GridType::Ptr
2660 const openvdb::math::Transform& xform,
2661 const std::vector<Vec3s>& points,
2662 const std::vector<Vec4I>& quads,
2665 std::vector<Vec3I> triangles(0);
2666 return doMeshConversion<GridType>(xform, points, triangles, quads,
2667 halfWidth, halfWidth);
2671 template<
typename Gr
idType>
2672 inline typename GridType::Ptr
2674 const openvdb::math::Transform& xform,
2675 const std::vector<Vec3s>& points,
2676 const std::vector<Vec3I>& triangles,
2677 const std::vector<Vec4I>& quads,
2680 return doMeshConversion<GridType>(xform, points, triangles, quads,
2681 halfWidth, halfWidth);
2685 template<
typename Gr
idType>
2686 inline typename GridType::Ptr
2688 const openvdb::math::Transform& xform,
2689 const std::vector<Vec3s>& points,
2690 const std::vector<Vec3I>& triangles,
2691 const std::vector<Vec4I>& quads,
2695 return doMeshConversion<GridType>(xform, points, triangles,
2696 quads, exBandWidth, inBandWidth);
2700 template<
typename Gr
idType>
2701 inline typename GridType::Ptr
2703 const openvdb::math::Transform& xform,
2704 const std::vector<Vec3s>& points,
2705 const std::vector<Vec3I>& triangles,
2706 const std::vector<Vec4I>& quads,
2709 return doMeshConversion<GridType>(xform, points, triangles, quads,
2710 bandWidth, bandWidth,
true);
2718 inline std::ostream&
2721 ostr <<
"{[ " << rhs.
mXPrim <<
", " << rhs.
mXDist <<
"]";
2722 ostr <<
" [ " << rhs.
mYPrim <<
", " << rhs.
mYDist <<
"]";
2723 ostr <<
" [ " << rhs.
mZPrim <<
", " << rhs.
mZDist <<
"]}";
2728 inline MeshToVoxelEdgeData::EdgeData
2743 const std::vector<Vec3s>& pointList,
2744 const std::vector<Vec4I>& polygonList);
2746 void run(
bool threaded =
true);
2749 inline void operator() (
const tbb::blocked_range<size_t> &range);
2757 struct Primitive {
Vec3d a, b, c, d;
Int32 index; };
2759 template<
bool IsQuad>
2760 inline void voxelize(
const Primitive&);
2762 template<
bool IsQuad>
2763 inline bool evalPrimitive(
const Coord&,
const Primitive&);
2765 inline bool rayTriangleIntersection(
const Vec3d& origin,
const Vec3d& dir,
2772 const std::vector<Vec3s>& mPointList;
2773 const std::vector<Vec4I>& mPolygonList;
2777 IntTreeT mLastPrimTree;
2784 const std::vector<Vec3s>& pointList,
2785 const std::vector<Vec4I>& polygonList)
2788 , mPointList(pointList)
2789 , mPolygonList(polygonList)
2791 , mLastPrimAccessor(mLastPrimTree)
2800 , mPointList(rhs.mPointList)
2801 , mPolygonList(rhs.mPolygonList)
2803 , mLastPrimAccessor(mLastPrimTree)
2812 tbb::parallel_reduce(tbb::blocked_range<size_t>(0, mPolygonList.size()), *
this);
2814 (*this)(tbb::blocked_range<size_t>(0, mPolygonList.size()));
2823 typedef RootNodeType::NodeChainType NodeChainType;
2824 BOOST_STATIC_ASSERT(boost::mpl::size<NodeChainType>::value > 1);
2825 typedef boost::mpl::at<NodeChainType, boost::mpl::int_<1> >::type InternalNodeType;
2833 for ( ; leafIt; ++leafIt) {
2834 ijk = leafIt->origin();
2840 mAccessor.addLeaf(rhs.mAccessor.
probeLeaf(ijk));
2841 InternalNodeType* node = rhs.mAccessor.
getNode<InternalNodeType>();
2843 rhs.mAccessor.
clear();
2853 if (!lhsLeafPt->isValueOn(offset)) {
2854 lhsLeafPt->setValueOn(offset, rhsValue);
2886 for (
size_t n = range.begin(); n < range.end(); ++n) {
2888 const Vec4I& verts = mPolygonList[n];
2890 prim.index =
Int32(n);
2891 prim.a =
Vec3d(mPointList[verts[0]]);
2892 prim.b =
Vec3d(mPointList[verts[1]]);
2893 prim.c =
Vec3d(mPointList[verts[2]]);
2896 prim.d =
Vec3d(mPointList[verts[3]]);
2897 voxelize<true>(prim);
2899 voxelize<false>(prim);
2905 template<
bool IsQuad>
2907 MeshToVoxelEdgeData::GenEdgeData::voxelize(
const Primitive& prim)
2909 std::deque<Coord> coordList;
2913 coordList.push_back(ijk);
2915 evalPrimitive<IsQuad>(ijk, prim);
2917 while (!coordList.empty()) {
2919 ijk = coordList.back();
2920 coordList.pop_back();
2922 for (
Int32 i = 0; i < 26; ++i) {
2923 nijk = ijk + util::COORD_OFFSETS[i];
2925 if (prim.index != mLastPrimAccessor.getValue(nijk)) {
2926 mLastPrimAccessor.setValue(nijk, prim.index);
2927 if(evalPrimitive<IsQuad>(nijk, prim)) coordList.push_back(nijk);
2934 template<
bool IsQuad>
2936 MeshToVoxelEdgeData::GenEdgeData::evalPrimitive(
const Coord& ijk,
const Primitive& prim)
2938 Vec3d uvw, org(ijk[0], ijk[1], ijk[2]);
2939 bool intersecting =
false;
2943 mAccessor.probeValue(ijk, edgeData);
2946 double dist = (org -
2949 if (rayTriangleIntersection(org,
Vec3d(1.0, 0.0, 0.0), prim.a, prim.c, prim.b, t)) {
2950 if (t < edgeData.mXDist) {
2951 edgeData.mXDist = t;
2952 edgeData.mXPrim = prim.index;
2953 intersecting =
true;
2957 if (rayTriangleIntersection(org,
Vec3d(0.0, 1.0, 0.0), prim.a, prim.c, prim.b, t)) {
2958 if (t < edgeData.mYDist) {
2959 edgeData.mYDist = t;
2960 edgeData.mYPrim = prim.index;
2961 intersecting =
true;
2965 if (rayTriangleIntersection(org,
Vec3d(0.0, 0.0, 1.0), prim.a, prim.c, prim.b, t)) {
2966 if (t < edgeData.mZDist) {
2967 edgeData.mZDist = t;
2968 edgeData.mZPrim = prim.index;
2969 intersecting =
true;
2975 double secondDist = (org -
2978 if (secondDist < dist) dist = secondDist;
2980 if (rayTriangleIntersection(org,
Vec3d(1.0, 0.0, 0.0), prim.a, prim.d, prim.c, t)) {
2981 if (t < edgeData.mXDist) {
2982 edgeData.mXDist = t;
2983 edgeData.mXPrim = prim.index;
2984 intersecting =
true;
2988 if (rayTriangleIntersection(org,
Vec3d(0.0, 1.0, 0.0), prim.a, prim.d, prim.c, t)) {
2989 if (t < edgeData.mYDist) {
2990 edgeData.mYDist = t;
2991 edgeData.mYPrim = prim.index;
2992 intersecting =
true;
2996 if (rayTriangleIntersection(org,
Vec3d(0.0, 0.0, 1.0), prim.a, prim.d, prim.c, t)) {
2997 if (t < edgeData.mZDist) {
2998 edgeData.mZDist = t;
2999 edgeData.mZPrim = prim.index;
3000 intersecting =
true;
3005 if (intersecting) mAccessor.setValue(ijk, edgeData);
3007 return (dist < 0.86602540378443861);
3012 MeshToVoxelEdgeData::GenEdgeData::rayTriangleIntersection(
3023 double divisor = s1.
dot(e1);
3024 if (!(std::abs(divisor) > 0.0))
return false;
3028 double inv_divisor = 1.0 / divisor;
3029 Vec3d d = origin - a;
3030 double b1 = d.
dot(s1) * inv_divisor;
3032 if (b1 < 0.0 || b1 > 1.0)
return false;
3035 double b2 = dir.
dot(s2) * inv_divisor;
3037 if (b2 < 0.0 || (b1 + b2) > 1.0)
return false;
3041 t = e2.dot(s2) * inv_divisor;
3042 return (t < 0.0) ?
false :
true;
3058 const std::vector<Vec3s>& pointList,
3059 const std::vector<Vec4I>& polygonList)
3073 std::vector<Vec3d>& points,
3074 std::vector<Index32>& primitives)
3084 point[0] = double(coord[0]) + data.
mXDist;
3085 point[1] = double(coord[1]);
3086 point[2] = double(coord[2]);
3088 points.push_back(point);
3089 primitives.push_back(data.
mXPrim);
3093 point[0] = double(coord[0]);
3094 point[1] = double(coord[1]) + data.
mYDist;
3095 point[2] = double(coord[2]);
3097 points.push_back(point);
3098 primitives.push_back(data.
mYPrim);
3102 point[0] = double(coord[0]);
3103 point[1] = double(coord[1]);
3104 point[2] = double(coord[2]) + data.
mZDist;
3106 points.push_back(point);
3107 primitives.push_back(data.
mZPrim);
3117 point[0] = double(coord[0]);
3118 point[1] = double(coord[1]) + data.
mYDist;
3119 point[2] = double(coord[2]);
3121 points.push_back(point);
3122 primitives.push_back(data.
mYPrim);
3126 point[0] = double(coord[0]);
3127 point[1] = double(coord[1]);
3128 point[2] = double(coord[2]) + data.
mZDist;
3130 points.push_back(point);
3131 primitives.push_back(data.
mZPrim);
3139 point[0] = double(coord[0]);
3140 point[1] = double(coord[1]) + data.
mYDist;
3141 point[2] = double(coord[2]);
3143 points.push_back(point);
3144 primitives.push_back(data.
mYPrim);
3153 point[0] = double(coord[0]) + data.
mXDist;
3154 point[1] = double(coord[1]);
3155 point[2] = double(coord[2]);
3157 points.push_back(point);
3158 primitives.push_back(data.
mXPrim);
3162 point[0] = double(coord[0]);
3163 point[1] = double(coord[1]) + data.
mYDist;
3164 point[2] = double(coord[2]);
3166 points.push_back(point);
3167 primitives.push_back(data.
mYPrim);
3177 point[0] = double(coord[0]) + data.
mXDist;
3178 point[1] = double(coord[1]);
3179 point[2] = double(coord[2]);
3181 points.push_back(point);
3182 primitives.push_back(data.
mXPrim);
3191 point[0] = double(coord[0]) + data.
mXDist;
3192 point[1] = double(coord[1]);
3193 point[2] = double(coord[2]);
3195 points.push_back(point);
3196 primitives.push_back(data.
mXPrim);
3200 point[0] = double(coord[0]);
3201 point[1] = double(coord[1]);
3202 point[2] = double(coord[2]) + data.
mZDist;
3204 points.push_back(point);
3205 primitives.push_back(data.
mZPrim);
3214 point[0] = double(coord[0]);
3215 point[1] = double(coord[1]);
3216 point[2] = double(coord[2]) + data.
mZDist;
3218 points.push_back(point);
3219 primitives.push_back(data.
mZPrim);
3229 #endif // OPENVDB_TOOLS_MESH_TO_VOLUME_HAS_BEEN_INCLUDED
OPENVDB_API Hermite min(const Hermite &, const Hermite &)
min and max operations done directly on the compressed data.
void addLeaf(LeafNodeT *leaf)
Add the specified leaf to this tree, possibly creating a child branch in the process. If the leaf node already exists, replace it.
Definition: ValueAccessor.h:328
void clear()
Remove all tiles from this tree and all nodes other than the root node.
Definition: Tree.h:533
NodeType * getNode()
Return the cached node of type NodeType. [Mainly for internal use].
Definition: ValueAccessor.h:303
void setValue(const Coord &xyz, const ValueType &value)
Set the value of the voxel at the given coordinates and mark the voxel as active. ...
Definition: ValueAccessor.h:241
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:653
void clearAllAccessors()
Clear all registered accessors.
Definition: Tree.h:1370
Index32 Index
Definition: Types.h:57
virtual void clear()
Remove all nodes from this cache, then reinsert the root node.
Definition: ValueAccessor.h:392
LeafNodeT * probeLeaf(const Coord &xyz)
Return a pointer to the leaf node that contains voxel (x, y, z), or NULL if no such node exists...
Definition: ValueAccessor.h:378
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:97
OPENVDB_API const Index32 INVALID_IDX
LeafIter beginLeaf()
Return an iterator over all leaf nodes in this tree.
Definition: Tree.h:1064
bool wasInterrupted(T *i, int percent=-1)
Definition: NullInterrupter.h:76
Type Pow2(Type x)
Return .
Definition: Math.h:458
static const Real LEVEL_SET_HALF_WIDTH
Definition: Types.h:176
const ValueT & getValue() const
Return the tile or voxel value to which this iterator is currently pointing.
Definition: TreeIterator.h:733
OPENVDB_API Vec3d closestPointOnTriangleToPoint(const Vec3d &a, const Vec3d &b, const Vec3d &c, const Vec3d &p, Vec3d &uvw)
Closest Point on Triangle to Point. Given a triangle abc and a point p, return the point on abc close...
LeafNodeType * probeLeaf(const Coord &xyz)
Return a pointer to the leaf node that contains voxel (x, y, z). If no such node exists, return NULL.
Definition: Tree.h:1602
boost::shared_ptr< Grid > Ptr
Definition: Grid.h:461
const ValueType & getValue(const Coord &xyz) const
Return the value of the voxel at the given coordinates.
Definition: ValueAccessor.h:210
Container class that associates a tree with a transform and metadata.
Definition: Grid.h:54
Vec3< T > cross(const Vec3< T > &v) const
Return the cross product of "this" vector and v;.
Definition: Vec3.h:228
Definition: Exceptions.h:87
#define OPENVDB_VERSION_NAME
Definition: version.h:45
void merge(Tree &other, MergePolicy=MERGE_ACTIVE_STATES)
Efficiently merge another tree into this tree using one of several schemes.
Definition: Tree.h:1693
T dot(const Vec3< T > &v) const
Dot product.
Definition: Vec3.h:199
bool probeValue(const Coord &xyz, ValueType &value) const
Return the active state of the voxel as well as its value.
Definition: ValueAccessor.h:220
static Accessor::ValueType result(const Accessor &grid, const Coord &ijk)
Definition: Operators.h:260
OPENVDB_API Hermite max(const Hermite &, const Hermite &)
min and max operations done directly on the compressed data.
bool isExactlyEqual(const T0 &a, const T1 &b)
Return true if a is exactly equal to b.
Definition: Math.h:351
Coord nearestCoord(const Vec3d &voxelCoord)
Return voxelCoord rounded to the closest integer coordinates.
Definition: util/Util.h:55
const LeafNodeT * probeConstLeaf(const Coord &xyz) const
Return a pointer to the leaf node that contains voxel (x, y, z), or NULL if no such node exists...
Definition: ValueAccessor.h:383
RootNodeType::LeafNodeType LeafNodeType
Definition: Tree.h:186
bool operator<(const Tuple< SIZE, T0 > &t0, const Tuple< SIZE, T1 > &t1)
Definition: Tuple.h:158
Definition: Operators.h:152
uint32_t Index32
Definition: Types.h:55
int32_t Int32
Definition: Types.h:59
math::Vec4< Index32 > Vec4I
Definition: Types.h:89
Base class for tree-traversal iterators over tile and voxel values.
Definition: TreeIterator.h:658
Base class for tree-traversal iterators over all leaf nodes (but not leaf voxels) ...
Definition: TreeIterator.h:1211
void foreach(const LeafOp &op, bool threaded=true, size_t grainSize=1)
Threaded method that applies a user-supplied functor to each leaf node in the LeafManager.
Definition: LeafManager.h:463
Vec3< double > Vec3d
Definition: Vec3.h:625
_RootNodeType RootNodeType
Definition: Tree.h:184
CopyConstness< TreeType, NonConstBufferType >::Type BufferType
Definition: LeafManager.h:117
OPENVDB_API const Coord COORD_OFFSETS[26]
coordinate offset table for neighboring voxels
bool operator>(const Tuple< SIZE, T0 > &t0, const Tuple< SIZE, T1 > &t1)
Definition: Tuple.h:170
bool isValueOn(const Coord &xyz) const
Return the active state of the voxel at the given coordinates.
Definition: ValueAccessor.h:217
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:67
bool isNegative(const Type &x)
Return true if x is less than zero.
Definition: Math.h:304
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:109
Tree< typename RootNodeType::template ValueConverter< Int32 >::Type > Type
Definition: Tree.h:198
void setValueOn(const Coord &xyz, const ValueType &value)
Set the value of the voxel at the given coordinates and mark the voxel as active. ...
Definition: ValueAccessor.h:246
BufferType & getBuffer(size_t leafIdx, size_t bufferIdx) const
Return the leaf or auxiliary buffer for the leaf node at index leafIdx. If bufferIdx is zero...
Definition: LeafManager.h:317
LeafNodeT * touchLeaf(const Coord &xyz)
Return a pointer to the leaf node that contains voxel (x, y, z). If no such node exists, create one, but preserve the values and active states of all voxels.
Definition: ValueAccessor.h:347