PeriDyno 1.0.0
Loading...
Searching...
No Matches
VolumeToTriangleSet.cpp
Go to the documentation of this file.
2
4
5namespace dyno
6{
7 template<typename TDataType>
13
14 template<typename TDataType>
19
20 template<typename TDataType>
22 {
23 Real iso = this->varIsoValue()->getValue();
24
25 auto levelset = this->ioVolume()->constDataPtr();
26
27 auto& sdf = levelset->getSDF();
28
29 Coord lowerBound = sdf.lowerBound();
30 Coord upperBound = sdf.upperBound();
31
32 Real h = sdf.getGridSpacing();
33
34 if (h < EPSILON)
35 return false;
36
37 int nx = (upperBound[0] - lowerBound[0]) / h;
38 int ny = (upperBound[1] - lowerBound[1]) / h;
39 int nz = (upperBound[2] - lowerBound[2]) / h;
40
41 DArray3D<Real> distances(nx + 1, ny + 1, nz + 1);
42 DArray<int> voxelVertNum(nx * ny * nz);
43
45 distances,
46 lowerBound,
47 h,
48 sdf);
49
51 voxelVertNum,
52 distances,
53 iso);
54
55 Reduction<int> reduce;
56 int totalVNum = reduce.accumulate(voxelVertNum.begin(), voxelVertNum.size());
57
58 Scan<int> scan;
59 scan.exclusive(voxelVertNum.begin(), voxelVertNum.size());
60
61 DArray<Coord> vertices(totalVNum);
62
63 DArray<TopologyModule::Triangle> triangles(totalVNum / 3);
64
66 vertices,
67 triangles,
68 voxelVertNum,
69 distances,
70 lowerBound,
71 iso,
72 h);
73
74 if (this->outTriangleSet()->isEmpty()) {
75 this->outTriangleSet()->allocate();
76 }
77
78 auto triSet = this->outTriangleSet()->getDataPtr();
79 triSet->setPoints(vertices);
80 triSet->setTriangles(triangles);
81 triSet->update();
82
83 distances.clear();
84 voxelVertNum.clear();
85 vertices.clear();
86 triangles.clear();
87
88 return true;
89 }
90
92}
#define DEFINE_CLASS(name)
Definition Object.h:140
static void reconstructSDF(DArray3D< Real > &distances, Coord origin, Real h, DistanceField3D< TDataType > &sdf)
static void countVerticeNumber(DArray< int > &num, DArray3D< Real > &distances, Real isoValue)
static void constructTriangles(DArray< Coord > &vertices, DArray< TopologyModule::Triangle > &triangles, DArray< int > &vertNum, DArray3D< Real > &distances, Coord origin, Real isoValue, Real h)
T accumulate(const T *val, const uint num)
void exclusive(T *output, const T *input, size_t length, bool bcao=true)
This is an implementation of AdditiveCCD based on peridyno.
Definition Array.h:25
Array< T, DeviceType::GPU > DArray
Definition Array.inl:89
constexpr Real EPSILON
Definition Typedef.inl:42
Array3D< T, DeviceType::GPU > DArray3D
Definition Array3D.inl:90