31 #ifndef OPENVDB_TOOLS_VOLUME_TO_MESH_HAS_BEEN_INCLUDED
32 #define OPENVDB_TOOLS_VOLUME_TO_MESH_HAS_BEEN_INCLUDED
34 #include <openvdb/openvdb.h>
35 #include <openvdb/tree/ValueAccessor.h>
36 #include <openvdb/tools/ValueTransformer.h>
37 #include <openvdb/util/Util.h>
38 #include <openvdb/math/Operators.h>
40 #include <openvdb/tools/Morphology.h>
41 #include <boost/scoped_array.hpp>
42 #include <boost/intrusive/slist.hpp>
44 #include <tbb/mutex.h>
45 #include <tbb/tick_count.h>
46 #include <tbb/blocked_range.h>
50 #include <boost/scoped_ptr.hpp>
78 , mTriangleFlags(NULL)
82 void resetQuads(
size_t size)
86 mQuadFlags.reset(
new char[mNumQuads]);
93 mQuadFlags.reset(NULL);
96 void resetTriangles(
size_t size)
100 mTriangleFlags.reset(
new char[mNumTriangles]);
103 void clearTriangles()
106 mTriangles.reset(NULL);
107 mTriangleFlags.reset(NULL);
110 const size_t&
numQuads()
const {
return mNumQuads; }
122 const char&
quadFlags(
size_t n)
const {
return mQuadFlags[n]; }
128 size_t mNumQuads, mNumTriangles;
129 boost::scoped_array<openvdb::Vec4I> mQuads;
130 boost::scoped_array<openvdb::Vec3I> mTriangles;
131 boost::scoped_array<char> mQuadFlags, mTriangleFlags;
152 VolumeToMesh(
double isovalue = 0,
double adaptivity = 0);
155 const size_t& pointListSize()
const;
159 const size_t& polygonPoolListSize()
const;
163 template<
typename Gr
idT>
164 void operator()(
const GridT&);
195 size_t mPointListSize, mPolygonPoolListSize;
196 double mIsovalue, mPrimAdaptivity, mSecAdaptivity;
267 return double(
hash(x)) * 2.3283064e-10;
275 template<
class AccessorT>
277 typename AccessorT::ValueType isovalue,
int dim)
281 if (accessor.getValue(coord) < isovalue) signs |= 1u;
283 if (accessor.getValue(coord) < isovalue) signs |= 2u;
285 if (accessor.getValue(coord) < isovalue) signs |= 4u;
287 if (accessor.getValue(coord) < isovalue) signs |= 8u;
288 coord[1] += dim; coord[2] = ijk[2];
289 if (accessor.getValue(coord) < isovalue) signs |= 16u;
291 if (accessor.getValue(coord) < isovalue) signs |= 32u;
293 if (accessor.getValue(coord) < isovalue) signs |= 64u;
295 if (accessor.getValue(coord) < isovalue) signs |= 128u;
300 template<
class AccessorT>
303 typename AccessorT::ValueType isovalue,
const int dim)
309 p[0] = accessor.getValue(coord) < isovalue;
311 p[1] = accessor.getValue(coord) < isovalue;
313 p[2] = accessor.getValue(coord) < isovalue;
315 p[3] = accessor.getValue(coord) < isovalue;
316 coord[1] += dim; coord[2] = ijk[2];
317 p[4] = accessor.getValue(coord) < isovalue;
319 p[5] = accessor.getValue(coord) < isovalue;
321 p[6] = accessor.getValue(coord) < isovalue;
323 p[7] = accessor.getValue(coord) < isovalue;
327 if (p[0]) signs |= 1u;
328 if (p[1]) signs |= 2u;
329 if (p[2]) signs |= 4u;
330 if (p[3]) signs |= 8u;
331 if (p[4]) signs |= 16u;
332 if (p[5]) signs |= 32u;
333 if (p[6]) signs |= 64u;
334 if (p[7]) signs |= 128u;
340 int i = ijk[0], ip = ijk[0] + hDim, ipp = ijk[0] + dim;
341 int j = ijk[1], jp = ijk[1] + hDim, jpp = ijk[1] + dim;
342 int k = ijk[2], kp = ijk[2] + hDim, kpp = ijk[2] + dim;
345 coord.
reset(ip, j, k);
346 m = accessor.getValue(coord) < isovalue;
347 if (p[0] != m && p[1] != m)
return true;
350 coord.
reset(ipp, j, kp);
351 m = accessor.getValue(coord) < isovalue;
352 if (p[1] != m && p[2] != m)
return true;
355 coord.
reset(ip, j, kpp);
356 m = accessor.getValue(coord) < isovalue;
357 if (p[2] != m && p[3] != m)
return true;
360 coord.
reset(i, j, kp);
361 m = accessor.getValue(coord) < isovalue;
362 if (p[0] != m && p[3] != m)
return true;
365 coord.
reset(ip, jpp, k);
366 m = accessor.getValue(coord) < isovalue;
367 if (p[4] != m && p[5] != m)
return true;
370 coord.
reset(ipp, jpp, kp);
371 m = accessor.getValue(coord) < isovalue;
372 if (p[5] != m && p[6] != m)
return true;
375 coord.
reset(ip, jpp, kpp);
376 m = accessor.getValue(coord) < isovalue;
377 if (p[6] != m && p[7] != m)
return true;
380 coord.
reset(i, jpp, kp);
381 m = accessor.getValue(coord) < isovalue;
382 if (p[7] != m && p[4] != m)
return true;
385 coord.
reset(i, jp, k);
386 m = accessor.getValue(coord) < isovalue;
387 if (p[0] != m && p[4] != m)
return true;
390 coord.
reset(ipp, jp, k);
391 m = accessor.getValue(coord) < isovalue;
392 if (p[1] != m && p[5] != m)
return true;
395 coord.
reset(ipp, jp, kpp);
396 m = accessor.getValue(coord) < isovalue;
397 if (p[2] != m && p[6] != m)
return true;
401 coord.
reset(i, jp, kpp);
402 m = accessor.getValue(coord) < isovalue;
403 if (p[3] != m && p[7] != m)
return true;
409 coord.
reset(ip, jp, k);
410 m = accessor.getValue(coord) < isovalue;
411 if (p[0] != m && p[1] != m && p[4] != m && p[5] != m)
return true;
414 coord.
reset(ipp, jp, kp);
415 m = accessor.getValue(coord) < isovalue;
416 if (p[1] != m && p[2] != m && p[5] != m && p[6] != m)
return true;
419 coord.
reset(ip, jp, kpp);
420 m = accessor.getValue(coord) < isovalue;
421 if (p[2] != m && p[3] != m && p[6] != m && p[7] != m)
return true;
424 coord.
reset(i, jp, kp);
425 m = accessor.getValue(coord) < isovalue;
426 if (p[0] != m && p[3] != m && p[4] != m && p[7] != m)
return true;
429 coord.
reset(ip, j, kp);
430 m = accessor.getValue(coord) < isovalue;
431 if (p[0] != m && p[1] != m && p[2] != m && p[3] != m)
return true;
434 coord.
reset(ip, jpp, kp);
435 m = accessor.getValue(coord) < isovalue;
436 if (p[4] != m && p[5] != m && p[6] != m && p[7] != m)
return true;
439 coord.
reset(ip, jp, kp);
440 m = accessor.getValue(coord) < isovalue;
441 if (p[0] != m && p[1] != m && p[2] != m && p[3] != m &&
442 p[4] != m && p[5] != m && p[6] != m && p[7] != m)
return true;
451 template <
class LeafType>
455 Coord ijk, end = start;
460 for (ijk[0] = start[0]; ijk[0] < end[0]; ++ijk[0]) {
461 for (ijk[1] = start[1]; ijk[1] < end[1]; ++ijk[1]) {
462 for (ijk[2] = start[2]; ijk[2] < end[2]; ++ijk[2]) {
463 if(leaf.isValueOn(ijk)) leaf.setValue(ijk, regionId);
472 template <
class LeafType>
475 typename LeafType::ValueType::value_type adaptivity)
477 if (adaptivity < 1e-5)
return false;
479 typedef typename LeafType::ValueType VecT;
480 Coord ijk, end = start;
485 std::vector<VecT> norms;
486 for (ijk[0] = start[0]; ijk[0] < end[0]; ++ijk[0]) {
487 for (ijk[1] = start[1]; ijk[1] < end[1]; ++ijk[1]) {
488 for (ijk[2] = start[2]; ijk[2] < end[2]; ++ijk[2]) {
490 if(!leaf.isValueOn(ijk))
continue;
491 norms.push_back(leaf.getValue(ijk));
496 size_t N = norms.size();
497 for (
size_t ni = 0; ni < N; ++ni) {
498 VecT n_i = norms[ni];
499 for (
size_t nj = 0; nj < N; ++nj) {
500 VecT n_j = norms[nj];
501 if ((1.0 - n_i.dot(n_j)) > adaptivity)
return false;
511 template <
class TreeT>
515 typedef std::vector<const typename TreeT::LeafNodeType *>
ListT;
519 mLeafNodes.reserve(tree.leafCount());
520 typename TreeT::LeafCIter iter = tree.cbeginLeaf();
521 for ( ; iter; ++iter) mLeafNodes.push_back(iter.getLeaf());
524 size_t size()
const {
return mLeafNodes.size(); }
526 const typename TreeT::LeafNodeType* operator[](
size_t n)
const
527 {
return mLeafNodes[n]; }
529 tbb::blocked_range<size_t> getRange()
const
530 {
return tbb::blocked_range<size_t>(0, mLeafNodes.size()); }
539 template <
class TreeT>
543 typedef std::vector<typename TreeT::LeafNodeType *>
ListT;
547 mLeafNodes.reserve(tree.leafCount());
548 typename TreeT::LeafIter iter = tree.beginLeaf();
549 for ( ; iter; ++iter) mLeafNodes.push_back(iter.getLeaf());
552 size_t size()
const {
return mLeafNodes.size(); }
554 typename TreeT::LeafNodeType* operator[](
size_t n)
const
555 {
return mLeafNodes[n]; }
557 tbb::blocked_range<size_t> getRange()
const
558 {
return tbb::blocked_range<size_t>(0, mLeafNodes.size()); }
570 template<
typename DistTreeT>
574 typedef typename DistTreeT::template ValueConverter<char>::Type
CharTreeT;
575 typedef typename DistTreeT::template ValueConverter<bool>::Type
BoolTreeT;
580 , mTopologyMaskTree(NULL)
581 , mSeamMaskTree(typename
BoolTreeT::Ptr())
588 return mDistTree && mEdgeTree && mTopologyMaskTree && mSeamMaskTree;
602 template <
class DistTreeT>
612 inline void operator()(
const tbb::blocked_range<size_t>&)
const;
615 std::vector<size_t>& mLeafRegionCount;
619 template <
class DistTreeT>
622 std::vector<size_t>& leafRegionCount)
624 , mLeafRegionCount(leafRegionCount)
629 template <
class DistTreeT>
632 : mLeafNodes(rhs.mLeafNodes)
633 , mLeafRegionCount(rhs.mLeafRegionCount)
638 template <
class DistTreeT>
642 tbb::parallel_for(mLeafNodes.getRange(), *
this);
646 template <
class DistTreeT>
650 (*this)(mLeafNodes.getRange());
654 template <
class DistTreeT>
658 for (
size_t n = range.begin(); n != range.end(); ++n) {
659 mLeafRegionCount[n] = size_t(mLeafNodes[n]->onVoxelCount());
667 template <
class DistTreeT>
672 typedef typename DistTreeT::template ValueConverter<bool>::Type
BoolTreeT;
673 typedef typename DistTreeT::template ValueConverter<int>::Type
IntTreeT;
676 const DistTreeT& distTree,
678 std::vector<size_t>& leafRegionCount,
689 void operator()(
const tbb::blocked_range<size_t>&)
const;
692 const DistTreeT& mDistTree;
694 std::vector<size_t>& mLeafRegionCount;
700 template <
class DistTreeT>
702 const DistTreeT& distTree,
704 std::vector<size_t>& leafRegionCount,
707 : mDistTree(distTree)
708 , mAuxLeafs(auxLeafs)
709 , mLeafRegionCount(leafRegionCount)
711 , mAdaptivity(adaptivity)
717 template <
class DistTreeT>
720 : mDistTree(rhs.mDistTree)
721 , mAuxLeafs(rhs.mAuxLeafs)
722 , mLeafRegionCount(rhs.mLeafRegionCount)
723 , mIsovalue(rhs.mIsovalue)
724 , mAdaptivity(rhs.mAdaptivity)
725 , mRefData(rhs.mRefData)
730 template <
class DistTreeT>
734 tbb::parallel_for(mAuxLeafs.getRange(), *
this);
738 template <
class DistTreeT>
742 (*this)(mAuxLeafs.getRange());
745 template <
class DistTreeT>
752 template <
class DistTreeT>
757 typedef typename BoolTreeT::LeafNodeType BoolLeafT;
758 typedef typename IntTreeT::LeafNodeType IntLeafT;
759 typedef typename BoolLeafT::template ValueConverter<Vec3T>::Type Vec3LeafT;
761 typedef typename IntLeafT::ValueOnIter IntIterT;
762 typedef typename BoolLeafT::ValueOnCIter BoolCIterT;
767 boost::scoped_ptr<BoolTreeAccessorT> seamMaskAcc;
768 boost::scoped_ptr<BoolTreeCAccessorT> topologyMaskAcc;
769 if (mRefData && mRefData->isValid()) {
770 seamMaskAcc.reset(
new BoolTreeAccessorT(*mRefData->mSeamMaskTree.get()));
771 topologyMaskAcc.reset(
new BoolTreeCAccessorT(*mRefData->mTopologyMaskTree));
773 const bool hasRefData = seamMaskAcc && topologyMaskAcc;
775 const int LeafDim = BoolLeafT::DIM;
780 Vec3LeafT gradientBuffer;
781 Coord ijk, coord, end;
783 for (
size_t n = range.begin(); n != range.end(); ++n) {
786 IntLeafT& auxLeaf = *mAuxLeafs[n];
788 const Coord& origin = auxLeaf.getOrigin();
789 end[0] = origin[0] + LeafDim;
790 end[1] = origin[1] + LeafDim;
791 end[2] = origin[2] + LeafDim;
797 const BoolLeafT* seamMask = seamMaskAcc->probeConstLeaf(origin);
798 if (seamMask != NULL) {
799 for (BoolCIterT it = seamMask->cbeginValueOn(); it; ++it) {
801 coord[0] = ijk[0] - (ijk[0] % 2);
802 coord[1] = ijk[1] - (ijk[1] % 2);
803 coord[2] = ijk[2] - (ijk[2] % 2);
804 mask.setActiveState(coord,
true);
807 if (topologyMaskAcc->probeConstLeaf(origin) == NULL) {
808 adaptivity = mRefData->mInternalAdaptivity;
813 for (IntIterT it = auxLeaf.beginValueOn(); it; ++it) {
815 coord[0] = ijk[0] - (ijk[0] % 2);
816 coord[1] = ijk[1] - (ijk[1] % 2);
817 coord[2] = ijk[2] - (ijk[2] % 2);
818 if(mask.isValueOn(coord))
continue;
819 mask.setActiveState(coord,
isAmbiguous(distAcc, ijk, mIsovalue, 1));
824 for (ijk[0] = origin[0]; ijk[0] < end[0]; ijk[0] += dim) {
825 for (ijk[1] = origin[1]; ijk[1] < end[1]; ijk[1] += dim) {
826 for (ijk[2] = origin[2]; ijk[2] < end[2]; ijk[2] += dim) {
828 mask.setActiveState(ijk,
true);
835 gradientBuffer.setValuesOff();
837 for (IntIterT it = auxLeaf.beginValueOn(); it; ++it) {
840 coord[0] = ijk[0] - (ijk[0] % dim);
841 coord[1] = ijk[1] - (ijk[1] % dim);
842 coord[2] = ijk[2] - (ijk[2] % dim);
843 if(mask.isValueOn(coord))
continue;
851 gradientBuffer.setValue(ijk, norm);
854 int regionId = 1, next_dim = dim << 1;
857 for (ijk[0] = 0; ijk[0] < LeafDim; ijk[0] += dim) {
858 coord[0] = ijk[0] - (ijk[0] % next_dim);
859 for (ijk[1] = 0; ijk[1] < LeafDim; ijk[1] += dim) {
860 coord[1] = ijk[1] - (ijk[1] % next_dim);
861 for (ijk[2] = 0; ijk[2] < LeafDim; ijk[2] += dim) {
862 coord[2] = ijk[2] - (ijk[2] % next_dim);
863 if(mask.isValueOn(ijk) || !
isMergable(gradientBuffer, ijk, dim, adaptivity)) {
864 mask.setActiveState(coord,
true);
874 for (dim = 4; dim < LeafDim; dim = dim << 1) {
876 coord[0] = ijk[0] - (ijk[0] % next_dim);
877 for (ijk[0] = origin[0]; ijk[0] < end[0]; ijk[0] += dim) {
878 coord[1] = ijk[1] - (ijk[1] % next_dim);
879 for (ijk[1] = origin[1]; ijk[1] < end[1]; ijk[1] += dim) {
880 coord[2] = ijk[2] - (ijk[2] % next_dim);
881 for (ijk[2] = origin[2]; ijk[2] < end[2]; ijk[2] += dim) {
883 if (mask.isValueOn(ijk) ||
isNonManifold(distAcc, ijk, mIsovalue, dim) ||
884 !
isMergable(gradientBuffer, ijk, dim, adaptivity)) {
885 mask.setActiveState(coord,
true);
895 if (!(mask.isValueOn(origin) ||
isNonManifold(distAcc, origin, mIsovalue, LeafDim))
896 &&
isMergable(gradientBuffer, origin, LeafDim, adaptivity)) {
902 size_t numVoxels = 0;
903 IntLeafT tmpLeaf(auxLeaf);
904 for (IntIterT it = tmpLeaf.beginValueOn(); it; ++it) {
905 if(it.getValue() == 0) {
911 while(tmpLeaf.onVoxelCount() > 0) {
913 IntIterT it = tmpLeaf.beginValueOn();
914 const int regionId = it.getValue();
916 if(it.getValue() == regionId) it.setValueOff();
920 mLeafRegionCount[n] = numVoxels;
930 template <
class DistTreeT>
936 typedef typename DistTreeT::template ValueConverter<int>::Type
IntTreeT;
939 const DistTreeT& distTree,
941 std::vector<size_t>& leafIndices,
942 const openvdb::math::Transform& xform,
953 void operator()(
const tbb::blocked_range<size_t>&)
const;
956 const DistTreeT& mDistTree;
958 const std::vector<size_t>& mLeafIndices;
959 const openvdb::math::Transform& mTransform;
961 const double mIsovalue;
964 double root(
double v0,
double v1)
const {
return (mIsovalue - v0) / (v1 - v0); }
970 template <
class DistTreeT>
972 const DistTreeT& distTree,
974 std::vector<size_t>& leafIndices,
975 const openvdb::math::Transform& xform,
978 : mDistTree(distTree)
979 , mAuxLeafs(auxLeafs)
980 , mLeafIndices(leafIndices)
989 template <
class DistTreeT>
991 : mDistTree(rhs.mDistTree)
992 , mAuxLeafs(rhs.mAuxLeafs)
993 , mLeafIndices(rhs.mLeafIndices)
994 , mTransform(rhs.mTransform)
995 , mPoints(rhs.mPoints)
996 , mIsovalue(rhs.mIsovalue)
997 , mRefData(rhs.mRefData)
1002 template <
class DistTreeT>
1007 mRefData = &refData;
1010 template <
class DistTreeT>
1014 tbb::parallel_for(mAuxLeafs.getRange(), *
this);
1018 template <
class DistTreeT>
1022 (*this)(mAuxLeafs.getRange());
1026 template <
class DistTreeT>
1030 typedef typename DistTreeT::template ValueConverter<bool>::Type BoolTreeT;
1032 typedef typename BoolTreeT::LeafNodeType BoolLeafT;
1033 typedef typename IntTreeT::LeafNodeType IntLeafT;
1036 boost::scoped_ptr<DistTreeAccessorT> refDistAcc;
1037 boost::scoped_ptr<BoolTreeAccessorT> refMaskAcc;
1038 if (mRefData && mRefData->isValid()) {
1040 refMaskAcc.reset(
new BoolTreeAccessorT(*mRefData->mTopologyMaskTree));
1042 const bool hasRefData = refDistAcc && refMaskAcc;
1044 typename IntTreeT::LeafNodeType::ValueOnIter auxIter;
1050 for (
size_t n = range.begin(); n != range.end(); ++n) {
1052 size_t idx = mLeafIndices[n];
1053 IntLeafT& auxLeaf = *mAuxLeafs[n];
1055 const BoolLeafT* maskLeaf = NULL;
1057 maskLeaf = refMaskAcc->probeConstLeaf(auxLeaf.getOrigin());
1060 for (auxIter = auxLeaf.beginValueOn(); auxIter; ++auxIter) {
1062 if(auxIter.getValue() == 0) {
1064 auxIter.setValue(idx);
1065 auxIter.setValueOff();
1066 ijk = auxIter.getCoord();
1068 if (hasRefData && maskLeaf && maskLeaf->isValueOn(ijk)) {
1070 calcAvgPoint(*refDistAcc.get(), ijk, avg);
1072 const int signChangeFlags = calcAvgPoint(*refDistAcc.get(), ijk, avg);
1075 const int flags = ~signChangeFlags & calcEdgeSigns(distAcc, ijk);
1077 const double r =
twist(ijk[0], ijk[1], ijk[2]);
1078 if (flags &
XEDGE) avg[0] = r;
1079 if (flags &
YEDGE) avg[1] = r;
1080 if (flags &
ZEDGE) avg[2] = r;
1085 calcAvgPoint(distAcc, ijk, avg);
1089 avg[0] += double(ijk[0]);
1090 avg[1] += double(ijk[1]);
1091 avg[2] += double(ijk[2]);
1093 avg = mTransform.indexToWorld(avg);
1096 ptn[0] = float(avg[0]);
1097 ptn[1] = float(avg[1]);
1098 ptn[2] = float(avg[2]);
1104 while(auxLeaf.onVoxelCount() > 0) {
1110 auxIter = auxLeaf.beginValueOn();
1111 int regionId = auxIter.getValue(), points = 0;
1113 for (; auxIter; ++auxIter) {
1114 if(auxIter.getValue() == regionId) {
1116 auxIter.setValue(idx);
1117 auxIter.setValueOff();
1118 ijk = auxIter.getCoord();
1120 if (hasRefData && maskLeaf && maskLeaf->isValueOn(ijk)) {
1121 calcAvgPoint(*refDistAcc.get(), ijk, tmp);
1123 calcAvgPoint(distAcc, ijk, tmp);
1127 tmp[0] += double(ijk[0]);
1128 tmp[1] += double(ijk[1]);
1129 tmp[2] += double(ijk[2]);
1131 avg += mTransform.indexToWorld(tmp);
1137 double w = 1.0 / double(points);
1144 ptn[0] = float(avg[0]);
1145 ptn[1] = float(avg[1]);
1146 ptn[2] = float(avg[2]);
1153 template <
class DistTreeT>
1162 signMask[0] = double(acc.getValue(coord)) < mIsovalue;
1165 signMask[1] = double(acc.getValue(coord)) < mIsovalue;
1168 signMask[2] = double(acc.getValue(coord)) < mIsovalue;
1171 signMask[3] = double(acc.getValue(coord)) < mIsovalue;
1173 coord[1] += 1; coord[2] = ijk[2];
1174 signMask[4] = double(acc.getValue(coord)) < mIsovalue;
1177 signMask[5] = double(acc.getValue(coord)) < mIsovalue;
1180 signMask[6] = double(acc.getValue(coord)) < mIsovalue;
1183 signMask[7] = double(acc.getValue(coord)) < mIsovalue;
1187 if (signMask[0] != signMask[1]) {
1191 if (signMask[1] != signMask[2]) {
1195 if (signMask[3] != signMask[2]) {
1199 if (signMask[0] != signMask[3]) {
1203 if (signMask[4] != signMask[5]) {
1207 if (signMask[5] != signMask[6]) {
1211 if (signMask[7] != signMask[6]) {
1215 if (signMask[4] != signMask[7]) {
1219 if (signMask[0] != signMask[4]) {
1223 if (signMask[1] != signMask[5]) {
1227 if (signMask[2] != signMask[6]) {
1231 if (signMask[3] != signMask[7]) {
1239 template <
class DistTreeT>
1241 PointGen<DistTreeT>::calcAvgPoint(DistTreeAccessorT& acc,
1250 values[0] = double(acc.getValue(coord));
1253 values[1] = double(acc.getValue(coord));
1256 values[2] = double(acc.getValue(coord));
1259 values[3] = double(acc.getValue(coord));
1261 coord[1] += 1; coord[2] = ijk[2];
1262 values[4] = double(acc.getValue(coord));
1265 values[5] = double(acc.getValue(coord));
1268 values[6] = double(acc.getValue(coord));
1271 values[7] = double(acc.getValue(coord));
1274 for (
int n = 0; n < 8; ++n) signMask[n] = (values[n] < mIsovalue);
1276 int samples = 0, edgeFlags = 0;
1281 if (signMask[0] != signMask[1]) {
1282 avg[0] += root(values[0], values[1]);
1287 if (signMask[1] != signMask[2]) {
1289 avg[2] += root(values[1], values[2]);
1294 if (signMask[3] != signMask[2]) {
1295 avg[0] += root(values[3], values[2]);
1301 if (signMask[0] != signMask[3]) {
1302 avg[2] += root(values[0], values[3]);
1307 if (signMask[4] != signMask[5]) {
1308 avg[0] += root(values[4], values[5]);
1314 if (signMask[5] != signMask[6]) {
1317 avg[2] += root(values[5], values[6]);
1322 if (signMask[7] != signMask[6]) {
1323 avg[0] += root(values[7], values[6]);
1330 if (signMask[4] != signMask[7]) {
1332 avg[2] += root(values[4], values[7]);
1337 if (signMask[0] != signMask[4]) {
1338 avg[1] += root(values[0], values[4]);
1343 if (signMask[1] != signMask[5]) {
1345 avg[1] += root(values[1], values[5]);
1350 if (signMask[2] != signMask[6]) {
1352 avg[1] += root(values[2], values[6]);
1358 if (signMask[3] != signMask[7]) {
1359 avg[1] += root(values[3], values[7]);
1366 double w = 1.0 / double(samples);
1385 mPolygonPool = &quadPool;
1393 mPolygonPool->
quad(mIdx) = verts;
1415 AdaptiveMeshOp(): mQuadIdx(0), mTriangleIdx(0), mPolygonPool(NULL), mTmpPolygonPool() {}
1419 mPolygonPool = &polygonPool;
1430 if (verts[0] != verts[1] && verts[0] != verts[2] && verts[0] != verts[3]
1431 && verts[1] != verts[2] && verts[1] != verts[3] && verts[2] != verts[3]) {
1432 mTmpPolygonPool.
quadFlags(mQuadIdx) = flags;
1433 addQuad(verts, reverse);
1435 verts[0] == verts[3] &&
1436 verts[1] != verts[2] &&
1437 verts[1] != verts[0] &&
1438 verts[2] != verts[0]) {
1440 addTriangle(verts[0], verts[1], verts[2], reverse);
1442 verts[1] == verts[2] &&
1443 verts[0] != verts[3] &&
1444 verts[0] != verts[1] &&
1445 verts[3] != verts[1]) {
1447 addTriangle(verts[0], verts[1], verts[3], reverse);
1449 verts[0] == verts[1] &&
1450 verts[2] != verts[3] &&
1451 verts[2] != verts[0] &&
1452 verts[3] != verts[0]) {
1454 addTriangle(verts[0], verts[2], verts[3], reverse);
1456 verts[2] == verts[3] &&
1457 verts[0] != verts[1] &&
1458 verts[0] != verts[2] &&
1459 verts[1] != verts[2]) {
1461 addTriangle(verts[0], verts[1], verts[2], reverse);
1469 for (
size_t i = 0; i < mQuadIdx; ++i) {
1470 mPolygonPool->
quad(i) = mTmpPolygonPool.
quad(i);
1476 for (
size_t i = 0; i < mTriangleIdx; ++i) {
1485 void addQuad(
const Vec4I& verts,
bool reverse)
1488 mTmpPolygonPool.
quad(mQuadIdx) = verts;
1490 Vec4I& quad = mTmpPolygonPool.
quad(mQuadIdx);
1499 void addTriangle(
unsigned v0,
unsigned v1,
unsigned v2,
bool reverse)
1515 size_t mQuadIdx, mTriangleIdx;
1516 PolygonPool *mPolygonPool;
1517 PolygonPool mTmpPolygonPool;
1524 template<
class DistTreeT,
class MeshingOp>
1528 typedef typename DistTreeT::template ValueConverter<char>::Type
CharTreeT;
1529 typedef typename DistTreeT::template ValueConverter<int>::Type
IntTreeT;
1539 void operator()(
const tbb::blocked_range<size_t>&)
const;
1550 template<
class DistTreeT,
class MeshingOp>
1553 : mEdgeLeafs(edgeLeafs)
1555 , mPolygonPoolList(polygonPoolList)
1561 template<
class DistTreeT,
class MeshingOp>
1563 : mEdgeLeafs(rhs.mEdgeLeafs)
1564 , mAuxTree(rhs.mAuxTree)
1565 , mPolygonPoolList(rhs.mPolygonPoolList)
1566 , mRefData(rhs.mRefData)
1571 template<
class DistTreeT,
class MeshingOp>
1576 mRefData = &refData;
1580 template<
class DistTreeT,
class MeshingOp>
1584 tbb::parallel_for(mEdgeLeafs.getRange(), *
this);
1588 template<
class DistTreeT,
class MeshingOp>
1592 (*this)(mEdgeLeafs.getRange());
1596 template<
class DistTreeT,
class MeshingOp>
1599 const tbb::blocked_range<size_t>& range)
const
1601 typedef typename DistTreeT::template ValueConverter<bool>::Type BoolTreeT;
1602 typedef typename BoolTreeT::LeafNodeType BoolLeafT;
1603 typedef typename CharTreeT::LeafNodeType CharLeafT;
1605 typedef openvdb::tree::ValueAccessor<const CharTreeT> CharTreeAccessorT;
1606 typedef openvdb::tree::ValueAccessor<const IntTreeT> IntTreeAccessorT;
1607 typedef openvdb::tree::ValueAccessor<const BoolTreeT> BoolTreeAccessorT;
1609 boost::scoped_ptr<CharTreeAccessorT> refEdgeAcc;
1610 boost::scoped_ptr<BoolTreeAccessorT> refMaskAcc;
1611 if (mRefData && mRefData->isValid()) {
1612 refEdgeAcc.reset(
new CharTreeAccessorT(*mRefData->mEdgeTree));
1613 refMaskAcc.reset(
new BoolTreeAccessorT(*mRefData->mSeamMaskTree.get()));
1615 const bool hasRefData = refEdgeAcc && refMaskAcc;
1618 typename CharTreeT::LeafNodeType::ValueOnCIter iter;
1619 IntTreeAccessorT auxAcc(mAuxTree);
1622 char refEdgeFlags, isSemLinePoly;
1629 for (
size_t n = range.begin(); n != range.end(); ++n) {
1631 const Coord origin = mEdgeLeafs[n]->getOrigin();
1635 iter = mEdgeLeafs[n]->cbeginValueOn();
1636 for (; iter; ++iter) {
1637 char edgeFlags = iter.getValue() >> 1;
1638 edgeCount += edgeFlags & 0x1;
1640 edgeFlags = edgeFlags >> 1;
1641 edgeCount += edgeFlags & 0x1;
1643 edgeFlags = edgeFlags >> 1;
1644 edgeCount += edgeFlags & 0x1;
1647 mesher.init(edgeCount, mPolygonPoolList[n]);
1649 const CharLeafT* refEdgeLeaf = NULL;
1650 const BoolLeafT* refMaskLeaf = NULL;
1653 refEdgeLeaf = refEdgeAcc->probeConstLeaf(origin);
1654 refMaskLeaf = refMaskAcc->probeConstLeaf(origin);
1657 iter = mEdgeLeafs[n]->cbeginValueOn();
1658 for (; iter; ++iter) {
1659 ijk = iter.getCoord();
1660 const char& edgeFlags = iter.getValue();
1665 if(refEdgeLeaf) refEdgeFlags = refEdgeLeaf->getValue(ijk);
1666 if (refMaskLeaf && refMaskLeaf->isValueOn(ijk)) {
1671 const bool isInside = edgeFlags &
INSIDE;
1672 int v0 = auxAcc.getValue(ijk);
1674 if (edgeFlags &
XEDGE) {
1677 coord[0] = ijk[0]; coord[1] = ijk[1]-1; coord[2] = ijk[2];
1678 quad[1] = auxAcc.getValue(coord);
1680 quad[2] = auxAcc.getValue(coord);
1682 quad[3] = auxAcc.getValue(coord);
1684 mesher.addPrim(quad, isInside,
1685 (isSemLinePoly | isExteriorPoly[
bool(refEdgeFlags & XEDGE)]));
1689 if (edgeFlags &
YEDGE) {
1692 coord[0] = ijk[0]; coord[1] = ijk[1]; coord[2] = ijk[2]-1;
1693 quad[1] = auxAcc.getValue(coord);
1695 quad[2] = auxAcc.getValue(coord);
1697 quad[3] = auxAcc.getValue(coord);
1699 mesher.addPrim(quad, isInside,
1700 (isSemLinePoly | isExteriorPoly[
bool(refEdgeFlags & YEDGE)]));
1703 if (edgeFlags &
ZEDGE) {
1706 coord[0] = ijk[0]; coord[1] = ijk[1]-1; coord[2] = ijk[2];
1707 quad[1] = auxAcc.getValue(coord);
1709 quad[2] = auxAcc.getValue(coord);
1711 quad[3] = auxAcc.getValue(coord);
1713 mesher.addPrim(quad, !isInside,
1714 (isSemLinePoly | isExteriorPoly[
bool(refEdgeFlags & ZEDGE)]));
1726 template<
class DistTreeT,
class AuxDataT =
int>
1731 typedef typename DistTreeT::ValueType
ValueT;
1733 typedef typename DistTreeT::template ValueConverter<char>::Type
CharTreeT;
1736 typedef typename DistTreeT::template ValueConverter<AuxDataT>::Type
AuxTreeT;
1740 double iso = 0.0,
bool extraCheck =
false);
1746 typename CharTreeT::Ptr
edgeTree()
const {
return mEdgeTree; }
1747 typename AuxTreeT::Ptr
auxTree()
const {
return mAuxTree; }
1749 void operator()(
const tbb::blocked_range<size_t>&);
1753 mEdgeTree->merge(*rhs.mEdgeTree);
1754 mAuxTree->merge(*rhs.mAuxTree);
1759 const DistTreeT& mSourceTree;
1762 typename CharTreeT::Ptr mEdgeTree;
1765 typename AuxTreeT::Ptr mAuxTree;
1768 const double mIsovalue;
1769 const bool mExtraCheck;
1771 int edgeCheck(
const Coord& ijk,
const bool thisInside);
1774 template<
class DistTreeT,
class AuxDataT>
1777 : mLeafNodes(leafNodes)
1779 , mSourceAccessor(mSourceTree)
1781 , mEdgeAccessor(*mEdgeTree)
1782 , mAuxTree(new
AuxTreeT(AuxDataT(0)))
1783 , mAuxAccessor(*mAuxTree)
1785 , mExtraCheck(extraCheck)
1789 template<
class DistTreeT,
class AuxDataT>
1791 : mLeafNodes(rhs.mLeafNodes)
1792 , mSourceTree(rhs.mSourceTree)
1793 , mSourceAccessor(mSourceTree)
1795 , mEdgeAccessor(*mEdgeTree)
1796 , mAuxTree(new
AuxTreeT(AuxDataT(0)))
1797 , mAuxAccessor(*mAuxTree)
1798 , mIsovalue(rhs.mIsovalue)
1799 , mExtraCheck(rhs.mExtraCheck)
1805 template<
class DistTreeT,
typename AuxDataT>
1809 tbb::parallel_reduce(mLeafNodes.getRange(), *
this);
1812 template<
class DistTreeT,
typename AuxDataT>
1816 (*this)(mLeafNodes.getRange());
1819 template<
class DistTreeT,
typename AuxDataT>
1823 typename DistTreeT::LeafNodeType::ValueOnCIter iter;
1830 for (
size_t n = range.begin(); n != range.end(); ++n) {
1831 for (iter = mLeafNodes[n]->cbeginValueOn(); iter; ++iter) {
1832 ijk = iter.getCoord();
1833 thisInside = iter.getValue() < mIsovalue;
1834 edgeFlags = edgeCheck(ijk, thisInside);
1836 if (edgeFlags != 0) {
1837 edgeFlags |= int(thisInside);
1838 mEdgeAccessor.setValue(ijk,
char(edgeFlags));
1843 for (
size_t n = range.begin(); n != range.end(); ++n) {
1844 for (iter = mLeafNodes[n]->cbeginValueOn(); iter; ++iter) {
1846 ijk = iter.getCoord();
1847 thisInside = iter.getValue() < mIsovalue;
1848 edgeFlags = edgeCheck(ijk, thisInside);
1850 if (edgeFlags != 0) {
1851 edgeFlags |= int(thisInside);
1852 mEdgeAccessor.setValue(ijk,
char(edgeFlags));
1856 if (!mSourceAccessor.probeValue(ijk, val)) {
1857 thisInside = val < mIsovalue;
1858 edgeFlags = edgeCheck(ijk, thisInside);
1860 if (edgeFlags != 0) {
1861 edgeFlags |= int(thisInside);
1862 mEdgeAccessor.setValue(ijk,
char(edgeFlags));
1868 if (!mSourceAccessor.probeValue(ijk, val)) {
1869 thisInside = val < mIsovalue;
1870 edgeFlags = edgeCheck(ijk, thisInside);
1872 if (edgeFlags != 0) {
1873 edgeFlags |= int(thisInside);
1874 mEdgeAccessor.setValue(ijk,
char(edgeFlags));
1880 if (!mSourceAccessor.probeValue(ijk, val)) {
1881 thisInside = val < mIsovalue;
1882 edgeFlags = edgeCheck(ijk, thisInside);
1884 if (edgeFlags != 0) {
1885 edgeFlags |= int(thisInside);
1886 mEdgeAccessor.setValue(ijk,
char(edgeFlags));
1894 template<
class DistTreeT,
typename AuxDataT>
1902 n_ijk = ijk; ++n_ijk[0];
1903 bool otherInside = (mSourceAccessor.getValue(n_ijk) < mIsovalue);
1904 if (otherInside != thisInside) {
1908 mAuxAccessor.setActiveState(ijk,
true);
1910 coord[0] = ijk[0]; coord[1] = ijk[1]-1; coord[2] = ijk[2];
1911 mAuxAccessor.setActiveState(coord,
true);
1913 coord[0] = ijk[0]; coord[1] = ijk[1]-1; coord[2] = ijk[2]-1;
1914 mAuxAccessor.setActiveState(coord,
true);
1916 coord[0] = ijk[0]; coord[1] = ijk[1]; coord[2] = ijk[2]-1;
1917 mAuxAccessor.setActiveState(coord,
true);
1921 n_ijk[0] = ijk[0]; ++n_ijk[1];
1922 otherInside = (mSourceAccessor.getValue(n_ijk) < mIsovalue);
1923 if (otherInside != thisInside) {
1927 mAuxAccessor.setActiveState(ijk,
true);
1929 coord[0] = ijk[0]; coord[1] = ijk[1]; coord[2] = ijk[2]-1;
1930 mAuxAccessor.setActiveState(coord,
true);
1932 coord[0] = ijk[0]-1; coord[1] = ijk[1]; coord[2] = ijk[2];
1933 mAuxAccessor.setActiveState(coord,
true);
1935 coord[0] = ijk[0]-1; coord[1] = ijk[1]; coord[2] = ijk[2]-1;
1936 mAuxAccessor.setActiveState(coord,
true);
1940 n_ijk[1] = ijk[1]; ++n_ijk[2];
1941 otherInside = (mSourceAccessor.getValue(n_ijk) < mIsovalue);
1942 if (otherInside != thisInside) {
1946 mAuxAccessor.setActiveState(ijk,
true);
1948 coord[0] = ijk[0]; coord[1] = ijk[1]-1; coord[2] = ijk[2];
1949 mAuxAccessor.setActiveState(coord,
true);
1951 coord[0] = ijk[0]-1; coord[1] = ijk[1]; coord[2] = ijk[2];
1952 mAuxAccessor.setActiveState(coord,
true);
1954 coord[0] = ijk[0]-1; coord[1] = ijk[1]-1; coord[2] = ijk[2];
1955 mAuxAccessor.setActiveState(coord,
true);
1964 template <
class DistTreeT>
1968 typedef typename DistTreeT::template ValueConverter<bool>::Type
BoolTreeT;
1969 typedef typename DistTreeT::template ValueConverter<int>::Type
IntTreeT;
1981 void operator()(
const tbb::blocked_range<size_t>&)
const;
1992 template <
class DistTreeT>
1995 : mSeamMaskLeafs(seamMaskLeafs)
1996 , mTopologyMaskTree(topologyMaskTree)
1997 , mTopologyMaskAcc(mTopologyMaskTree)
2003 template <
class DistTreeT>
2005 : mSeamMaskLeafs(rhs.mSeamMaskLeafs)
2006 , mTopologyMaskTree(rhs.mTopologyMaskTree)
2007 , mTopologyMaskAcc(mTopologyMaskTree)
2008 , mAuxTree(rhs.mAuxTree)
2013 template <
class DistTreeT>
2017 tbb::parallel_for(mSeamMaskLeafs.getRange(), *
this);
2020 template <
class DistTreeT>
2024 (*this)(mSeamMaskLeafs.getRange());
2027 template <
class DistTreeT>
2031 typedef typename BoolTreeT::LeafNodeType::ValueOnIter ValueOnIterT;
2033 for (
size_t n = range.begin(); n != range.end(); ++n) {
2034 ValueOnIterT it = mSeamMaskLeafs[n]->beginValueOn();
2036 ijk = it.getCoord();
2037 if (!mTopologyMaskAcc.isValueOn(ijk)) {
2040 bool turnOff =
true;
2041 for (
size_t n = 0; n < 6; ++n) {
2043 if (!mAuxTree.isValueOn(n_ijk) && mTopologyMaskAcc.isValueOn(n_ijk)) {
2048 if (turnOff) it.setValueOff();
2063 , mPolygonPoolListSize(0)
2064 , mIsovalue(1e-4 + isovalue)
2065 , mPrimAdaptivity(adaptivity)
2066 , mSecAdaptivity(0.0)
2069 , mRefTopologyMaskTree(
TreeBase::Ptr())
2079 inline const size_t&
2082 return mPointListSize;
2100 inline const size_t&
2103 return mPolygonPoolListSize;
2111 mSecAdaptivity = secAdaptivity;
2118 template<
typename Gr
idT>
2122 typedef typename GridT::TreeType DistTreeT;
2123 typedef typename DistTreeT::ValueType DistValueT;
2125 typedef typename DistTreeT::template ValueConverter<char>::Type CharTreeT;
2126 typedef typename DistTreeT::template ValueConverter<int>::Type IntTreeT;
2127 typedef typename DistTreeT::template ValueConverter<bool>::Type BoolTreeT;
2129 const bool noAdaptivity = mPrimAdaptivity < 1e-6 && mSecAdaptivity < 1e-6;
2131 const openvdb::math::Transform& transform = distGrid.transform();
2132 const DistTreeT& distTree = distGrid.tree();
2133 typename CharTreeT::Ptr edgeTree;
2134 typename IntTreeT::Ptr auxTree;
2136 const bool extraSignCheck = std::abs(mIsovalue -
double(distGrid.background()))
2137 < 2.0 * double(transform.voxelSize()[0]);
2151 const GridT* refGrid =
static_cast<const GridT*
>(mRefGrid.get());
2152 if (refGrid && refGrid->activeVoxelCount() > 0) {
2158 if (!mRefEdgeTree && !mRefTopologyMaskTree) {
2161 *refData.
mDistTree, leafs, mIsovalue, extraSignCheck);
2164 mRefTopologyMaskTree = op.
auxTree();
2167 if (mRefEdgeTree && mRefTopologyMaskTree) {
2168 refData.
mEdgeTree =
static_cast<CharTreeT*
>(mRefEdgeTree.get());
2169 refData.
mTopologyMaskTree =
static_cast<BoolTreeT*
>(mRefTopologyMaskTree.get());
2170 refData.
mSeamMaskTree =
typename BoolTreeT::Ptr(
new BoolTreeT(
false));
2192 std::vector<size_t> regions(auxLeafs.
size(), 0);
2200 regions, DistValueT(mIsovalue), DistValueT(mPrimAdaptivity));
2207 for (
size_t n = 0, N = regions.size(); n < N; ++n) {
2209 regions[n] = mPointListSize;
2210 mPointListSize += tmp;
2218 pointGen(distTree, auxLeafs, regions, transform, mPoints, mIsovalue);
2226 mPolygonPoolListSize = edgeLeafs.
size();
2227 mPolygons.reset(
new PolygonPool[mPolygonPoolListSize]);
2231 meshGen(edgeLeafs, *auxTree, mPolygons);
2236 meshGen(edgeLeafs, *auxTree, mPolygons);
2247 #endif // OPENVDB_TOOLS_VOLUME_TO_MESH_HAS_BEEN_INCLUDED