/
basicmesh.h
193 lines (162 loc) · 4.97 KB
/
basicmesh.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
#ifndef BASICMESH_H
#define BASICMESH_H
#ifndef MKL_BLAS
#define MKL_BLAS MKL_DOMAIN_BLAS
#endif
#define EIGEN_USE_MKL_ALL
#include <eigen3/Eigen/Dense>
#include "common.h"
#include <random>
#include <OpenMesh/Core/IO/MeshIO.hh>
#include <OpenMesh/Core/Mesh/TriMesh_ArrayKernelT.hh>
typedef OpenMesh::TriMesh_ArrayKernelT<> HalfEdgeMesh;
using namespace std;
using namespace Eigen;
class BasicMesh
{
public:
BasicMesh() {}
BasicMesh(const string& filename);
void set_vertex(int i, const Vector3d& v) {
verts.row(i) = v;
}
Vector3d vertex(int i) const {
return verts.row(i);
}
Vector3i face(int i) const {
return faces.row(i);
}
Vector3i face_texture(int i) const {
return face_tex_index.row(i);
}
Vector3d normal(int i) const {
return norms.row(i);
}
Vector3d vertex_normal(int i) const {
return vertex_norms.row(i);
}
Vector2d texture_coords(int i) const {
return texcoords.row(i);
}
const MatrixX3d& vertices() const {
return verts;
}
MatrixX3d& vertices() {
return verts;
}
int NumVertices() const { return static_cast<int>(verts.rows()); }
int NumFaces() const { return static_cast<int>(faces.rows()); }
bool LoadOBJMesh(const string &filename);
void ComputeNormals();
void UpdateVertices(const VectorXd& vertices);
void Subdivide();
vector<int> GetNeighbors() const;
void Write(const string& filename) const;
template <typename Pred>
vector<int> filterFaces(Pred p) {
vector<int> v;
for(int i=0;i<NumFaces();++i) {
if( p(faces.row(i)) ) {
v.push_back(i);
}
}
return v;
}
template <typename Pred>
vector<int> filterVertices(Pred p) {
vector<int> v;
for(int i=0;i<NumVertices();++i) {
if( p(verts.row(i)) ) {
v.push_back(i);
}
}
return v;
}
MatrixX3d samplePoints(int points_per_face, double zcutoff) const {
int npoints = 0;
vector<int> validfaces;
for (int i = 0; i < NumFaces(); ++i) {
// sample 8 points per face
int v1 = faces(i, 0), v2 = faces(i, 1), v3 = faces(i, 2);
double z1 = verts(v1, 2), z2 = verts(v2, 2), z3 = verts(v3, 2);
double zc = (z1 + z2 + z3) / 3.0;
if (zc > zcutoff) {
npoints += points_per_face;
validfaces.push_back(i);
}
}
return samplePoints(points_per_face, set<int>(validfaces.begin(), validfaces.end()));
}
MatrixX3d samplePoints(int points_per_face, const set<int>& faces_to_sample) const {
int npoints = faces_to_sample.size() * points_per_face;
cout << "npoints = " << npoints << endl;
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_real_distribution<> dis(0, 1);
vector<int> validfaces(faces_to_sample.begin(), faces_to_sample.end());
MatrixXd P(npoints, 3);
for (size_t i = 0, offset=0; i < validfaces.size(); ++i) {
int fidx = validfaces[i];
#if 0
int v1 = faces(fidx, 0), v2 = faces(fidx, 1), v3 = faces(fidx, 2);
double x1 = verts(v1, 0), x2 = verts(v2, 0), x3 = verts(v3, 0);
double y1 = verts(v1, 1), y2 = verts(v2, 1), y3 = verts(v3, 1);
double z1 = verts(v1, 2), z2 = verts(v2, 2), z3 = verts(v3, 2);
for(int j=0;j<points_per_face;++j) {
// sample a point
double alpha = rand()/(double)RAND_MAX,
beta = rand()/(double)RAND_MAX * (1-alpha),
gamma = 1.0 - alpha - beta;
double xij = x1*alpha + x2*beta + x3*gamma;
double yij = y1*alpha + y2*beta + y3*gamma;
double zij = z1*alpha + z2*beta + z3*gamma;
P.row(offset) = Vector3d(xij, yij, zij); ++offset;
}
#else
Vector3d v0 = vertex(faces(fidx, 0));
Vector3d v1 = vertex(faces(fidx, 1));
Vector3d v2 = vertex(faces(fidx, 2));
Vector3d v0v1 = v1 - v0;
Vector3d v0v2 = v2 - v0;
for(int j=0;j<points_per_face;++j) {
Vector3d v;
while(true) {
double alpha = dis(gen);
double beta = dis(gen);
v = alpha * v0v1 + beta * v0v2;
if (v.dot(v0v1) <= v0v1.norm() && v.dot(v0v2) <= v0v2.norm()) {
v = v + v0;
break;
}
}
P.row(offset) = v; ++offset;
}
#endif
}
cout << "points sampled." << endl;
return P;
}
void BuildHalfEdgeMesh();
private:
unordered_map<int, int> vert_face_map;
MatrixX3d verts;
MatrixX3i faces, face_tex_index;
MatrixX2d texcoords;
// Per-face normal vector
MatrixX3d norms;
MatrixX3d vertex_norms;
template <typename Handle>
struct HandleHasher {
std::size_t operator()(const Handle& h) const {
return h.idx();
}
};
// This mesh stored in half edge data structure
// Makes mesh manipulation easier
vector<HalfEdgeMesh::VertexHandle> vhandles;
unordered_map<HalfEdgeMesh::VertexHandle, int, HandleHasher<HalfEdgeMesh::VertexHandle>> vhandles_map;
vector<HalfEdgeMesh::FaceHandle> fhandles;
unordered_map<HalfEdgeMesh::FaceHandle, int, HandleHasher<HalfEdgeMesh::FaceHandle>> fhandles_map;
HalfEdgeMesh hemesh;
};
#endif // BASICMESH_H