Skip to content

Commit 067e370

Browse files
MarkGillespieMarkGillespie
committed
Simplification & subdivision
1 parent 648c71e commit 067e370

File tree

8 files changed

+393
-48
lines changed

8 files changed

+393
-48
lines changed

‎include/geometrycentral/surface/mutation_manager.h‎

Lines changed: 15 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -128,20 +128,23 @@ class MutationManager {
128128
// to `tSplit`.
129129
// In general, both the `tSplit` and the `newVertexPosition` are used; `tSplit` is necessary to allow callbacks to
130130
// interpolate data; if called with either the other will be inferred.
131-
void splitEdge(Edge e, double tSplit);
132-
void splitEdge(Edge e, Vector3 newVertexPosition);
133-
void splitEdge(Edge e, double tSplit, Vector3 newVertexPosition);
131+
// Returns the new halfedge which points away from the new vertex (so he.vertex() is the new vertex), and is the same
132+
// direction as e.halfedge() on the original edge. The halfedge direction of the other part of the new split edge is
133+
// also preserved.
134+
Halfedge splitEdge(Edge e, double tSplit);
135+
Halfedge splitEdge(Edge e, Vector3 newVertexPosition);
136+
Halfedge splitEdge(Edge e, double tSplit, Vector3 newVertexPosition);
134137

135138
// Collapse an edge.
136-
// Returns true if the edge could actually be collapsed.
137-
bool collapseEdge(Edge e, double tCollapse);
138-
bool collapseEdge(Edge e, Vector3 newVertexPosition);
139-
bool collapseEdge(Edge e, double tCollapse, Vector3 newVertexPosition);
140-
141-
// Split a face (i.e. insert a vertex into the face)
142-
void splitFace(Face f, const std::vector<double>& bSplit);
143-
void splitFace(Face f, Vector3 newVertexPosition);
144-
void splitFace(Face f, const std::vector<double>& bSplit, Vector3 newVertexPosition);
139+
// Returns the new vertex if the edge could be collapsed, and Vertex() otherwise
140+
Vertex collapseEdge(Edge e, double tCollapse);
141+
Vertex collapseEdge(Edge e, Vector3 newVertexPosition);
142+
Vertex collapseEdge(Edge e, double tCollapse, Vector3 newVertexPosition);
143+
144+
// Split a face (i.e. insert a vertex into the face), and return the new vertex
145+
Vertex splitFace(Face f, const std::vector<double>& bSplit);
146+
Vertex splitFace(Face f, Vector3 newVertexPosition);
147+
Vertex splitFace(Face f, const std::vector<double>& bSplit, Vector3 newVertexPosition);
145148

146149
// ======================================================
147150
// ======== High-level mutations

‎include/geometrycentral/surface/quadric_error_simplification.h‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -4,6 +4,8 @@
44
#include "geometrycentral/surface/mutation_manager.h"
55
#include "geometrycentral/surface/vertex_position_geometry.h"
66

7+
#include "geometrycentral/utilities/elementary_geometry.h"
8+
79
#include <queue>
810

911
namespace geometrycentral {
@@ -29,8 +31,7 @@ class Quadric {
2931
Quadric operator+(const Quadric& Q1, const Quadric& Q2);
3032

3133
void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, double tol = 0.05);
32-
void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, double tol = 0.05,
33-
MutationManager& mm);
34+
void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, double tol, MutationManager& mm);
3435

3536
} // namespace surface
3637
} // namespace geometrycentral

‎include/geometrycentral/surface/subdivide.h‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
#pragma once
2+
3+
#include "geometrycentral/surface/manifold_surface_mesh.h"
4+
#include "geometrycentral/surface/mutation_manager.h"
5+
#include "geometrycentral/surface/vertex_position_geometry.h"
6+
7+
namespace geometrycentral {
8+
namespace surface {
9+
10+
// Subdivide to form quad mesh
11+
// Note: these do not take MutationManager arguments since MutationManager only supports triangular meshes at the moment
12+
void linearSubdivide(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo);
13+
void catmullClarkSubdivide(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo);
14+
15+
void loopSubdivide(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo);
16+
void loopSubdivide(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, MutationManager& mm);
17+
18+
} // namespace surface
19+
} // namespace geometrycentral

‎src/CMakeLists.txt‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,8 @@ SET(SRCS
4646
surface/tufted_laplacian.cpp
4747
surface/flip_geodesics.cpp
4848
surface/transfer_functions.cpp
49+
surface/quadric_error_simplification.cpp
50+
surface/subdivide.cpp
4951
#surface/detect_symmetry.cpp
5052
#surface/mesh_ray_tracer.cpp
5153

@@ -105,13 +107,15 @@ SET(HEADERS
105107
${INCLUDE_ROOT}/surface/mesh_graph_algorithms.h
106108
${INCLUDE_ROOT}/surface/mesh_ray_tracer.h
107109
${INCLUDE_ROOT}/surface/parameterize.h
110+
${INCLUDE_ROOT}/surface/quadric_error_simplification.h
108111
${INCLUDE_ROOT}/surface/rich_surface_mesh_data.h
109112
${INCLUDE_ROOT}/surface/rich_surface_mesh_data.ipp
110113
${INCLUDE_ROOT}/surface/polygon_soup_mesh.h
111114
${INCLUDE_ROOT}/surface/signpost_intrinsic_triangulation.h
112115
${INCLUDE_ROOT}/surface/signpost_intrinsic_triangulation.ipp
113116
${INCLUDE_ROOT}/surface/simple_idt.h
114117
${INCLUDE_ROOT}/surface/simple_polygon_mesh.h
118+
${INCLUDE_ROOT}/surface/subdivide.h
115119
${INCLUDE_ROOT}/surface/surface_centers.h
116120
${INCLUDE_ROOT}/surface/surface_mesh.h
117121
${INCLUDE_ROOT}/surface/surface_mesh.ipp

‎src/surface/manifold_surface_mesh.cpp‎

Lines changed: 58 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1263,6 +1263,64 @@ Vertex ManifoldSurfaceMesh::collapseEdgeTriangular(Edge e) {
12631263
modificationTick++;
12641264
}
12651265

1266+
Face ManifoldSurfaceMesh::removeEdge(Edge e) {
1267+
if (e.isBoundary()) {
1268+
throw std::runtime_error("not implemented");
1269+
}
1270+
1271+
// Halfedges/edges/faces that will be removed
1272+
// (except first face)
1273+
std::vector<Halfedge> toRemove{e.halfedge(), e.halfedge().twin()};
1274+
std::vector<Halfedge> ringHalfedges;
1275+
for (Halfedge heStart : toRemove) {
1276+
Halfedge he = heStart.next();
1277+
while (he != heStart) {
1278+
// The one-ring must not contain any other copies of e, or we cannot remove the edge
1279+
if (he.edge() == e) {
1280+
return Face();
1281+
}
1282+
ringHalfedges.push_back(he);
1283+
he = he.next();
1284+
}
1285+
}
1286+
1287+
// If both faces are the same, we cannot remove the edge. This should have been caught above
1288+
if (toRemove[0].face() == toRemove[1].face()) {
1289+
return Face();
1290+
}
1291+
Face keepFace = toRemove[0].face();
1292+
1293+
1294+
// Record these before we break pointers
1295+
Vertex src = e.halfedge().vertex();
1296+
Vertex dst = e.halfedge().twin().vertex();
1297+
Halfedge altSrcHedge = e.halfedge().twin().next();
1298+
Halfedge altDstHedge = e.halfedge().next();
1299+
1300+
// Hook up next and face refs for the halfedges along the ring
1301+
size_t N = ringHalfedges.size();
1302+
for (size_t i = 0; i < N; i++) {
1303+
heNextArr[ringHalfedges[i].getIndex()] = ringHalfedges[(i + 1) % N].getIndex();
1304+
heFaceArr[ringHalfedges[i].getIndex()] = keepFace.getIndex();
1305+
}
1306+
1307+
// only update vHalfedgeArr if needed to avoid disturbing boundary halfedges
1308+
if (src.halfedge().edge() == e) {
1309+
vHalfedgeArr[src.getIndex()] = altSrcHedge.getIndex();
1310+
}
1311+
if (dst.halfedge().edge() == e) {
1312+
vHalfedgeArr[dst.getIndex()] = altDstHedge.getIndex();
1313+
}
1314+
1315+
fHalfedgeArr[keepFace.getIndex()] = ringHalfedges[0].getIndex();
1316+
1317+
// Actually delete all of the elements
1318+
deleteElement(toRemove[1].face());
1319+
deleteEdgeBundle(e);
1320+
1321+
modificationTick++;
1322+
return keepFace;
1323+
}
12661324

12671325

12681326
bool ManifoldSurfaceMesh::removeFaceAlongBoundary(Face f) {

‎src/surface/mutation_manager.cpp‎

Lines changed: 21 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -59,16 +59,16 @@ bool MutationManager::flipEdge(Edge e) {
5959
return true;
6060
}
6161

62-
void MutationManager::splitEdge(Edge e, double tSplit) {
62+
Halfedge MutationManager::splitEdge(Edge e, double tSplit) {
6363
Vector3 newPos{0., 0., 0.};
6464
if (geometry) {
6565
VertexData<Vector3>& pos = geometry->vertexPositions;
6666
newPos = (1. - tSplit) * pos[e.halfedge().tailVertex()] + tSplit * pos[e.halfedge().tipVertex()];
6767
}
68-
splitEdge(e, tSplit, newPos);
68+
return splitEdge(e, tSplit, newPos);
6969
}
7070

71-
void MutationManager::splitEdge(Edge e, Vector3 newVertexPosition) {
71+
Halfedge MutationManager::splitEdge(Edge e, Vector3 newVertexPosition) {
7272

7373
double tSplit = -1;
7474
GC_SAFETY_ASSERT(geometry, "must have geometry to split by position");
@@ -80,10 +80,10 @@ void MutationManager::splitEdge(Edge e, Vector3 newVertexPosition) {
8080
tSplit = pointLineSegmentNeaestLocation(newVertexPosition, posTail, posTip);
8181
}
8282

83-
splitEdge(e, tSplit, newVertexPosition);
83+
return splitEdge(e, tSplit, newVertexPosition);
8484
}
8585

86-
void MutationManager::splitEdge(Edge e, double tSplit, Vector3 newVertexPosition) {
86+
Halfedge MutationManager::splitEdge(Edge e, double tSplit, Vector3 newVertexPosition) {
8787

8888
// Invoke before callbacks
8989
for (EdgeSplitPolicy* policy : edgeSplitPolicies) {
@@ -102,11 +102,13 @@ void MutationManager::splitEdge(Edge e, double tSplit, Vector3 newVertexPosition
102102
for (EdgeSplitPolicy* policy : edgeSplitPolicies) {
103103
policy->afterEdgeSplit(newHeFront, newHeBack, tSplit);
104104
}
105+
106+
return newHeFront;
105107
}
106108

107109
// Collapse an edge.
108-
// Returns true if the edge could actually be collapsed.
109-
bool MutationManager::collapseEdge(Edge e, double tCollapse) {
110+
// Returns the new vertex if the edge could be collapsed, and Vertex() otherwise
111+
Vertex MutationManager::collapseEdge(Edge e, double tCollapse) {
110112
Vector3 newPos{0., 0., 0.};
111113
if (geometry) {
112114
// Find the nearest tCoord
@@ -116,7 +118,7 @@ bool MutationManager::collapseEdge(Edge e, double tCollapse) {
116118
return collapseEdge(e, tCollapse, newPos);
117119
}
118120

119-
bool MutationManager::collapseEdge(Edge e, Vector3 newVertexPosition) {
121+
Vertex MutationManager::collapseEdge(Edge e, Vector3 newVertexPosition) {
120122

121123
double tCollapse = -1;
122124
GC_SAFETY_ASSERT(geometry, "must have geometry to split by position");
@@ -130,8 +132,7 @@ bool MutationManager::collapseEdge(Edge e, Vector3 newVertexPosition) {
130132

131133
return collapseEdge(e, tCollapse, newVertexPosition);
132134
}
133-
134-
bool MutationManager::collapseEdge(Edge e, double tCollapse, Vector3 newVertexPosition) {
135+
Vertex MutationManager::collapseEdge(Edge e, double tCollapse, Vector3 newVertexPosition) {
135136

136137
// Invoke before callbacks
137138
// TODO need to handle possiblity that collapse fails -- check before calling
@@ -140,7 +141,7 @@ bool MutationManager::collapseEdge(Edge e, double tCollapse, Vector3 newVertexPo
140141
}
141142

142143
Vertex newV = mesh.collapseEdgeTriangular(e);
143-
if (newV == Vertex()) return false;
144+
if (newV == Vertex()) return Vertex();
144145

145146
if (geometry) {
146147
VertexData<Vector3>& pos = geometry->vertexPositions;
@@ -152,11 +153,11 @@ bool MutationManager::collapseEdge(Edge e, double tCollapse, Vector3 newVertexPo
152153
policy->afterEdgeCollapse(newV, tCollapse);
153154
}
154155

155-
return true;
156+
return newV;
156157
}
157158

158159
// Split a face (i.e. insert a vertex into the face)
159-
void MutationManager::splitFace(Face f, const std::vector<double>& bSplit) {
160+
Vertex MutationManager::splitFace(Face f, const std::vector<double>& bSplit) {
160161
Vector3 newPos = Vector3::zero();
161162
if (geometry) {
162163
size_t iV = 0;
@@ -167,16 +168,17 @@ void MutationManager::splitFace(Face f, const std::vector<double>& bSplit) {
167168
}
168169
}
169170

170-
splitFace(f, bSplit, newPos);
171+
return splitFace(f, bSplit, newPos);
171172
}
172173

173-
void MutationManager::splitFace(Face f, Vector3 newVertexPosition) {
174+
Vertex MutationManager::splitFace(Face f, Vector3 newVertexPosition) {
174175
// TODO
175176
throw std::runtime_error("Face split based on vertex position not implemented yet");
177+
return Vertex();
176178
}
177179

178180

179-
void MutationManager::splitFace(Face f, const std::vector<double>& bSplit, Vector3 newVertexPosition) {
181+
Vertex MutationManager::splitFace(Face f, const std::vector<double>& bSplit, Vector3 newVertexPosition) {
180182
// Invoke before callbacks
181183
for (FaceSplitPolicy* policy : faceSplitPolicies) {
182184
policy->beforeFaceSplit(f, bSplit);
@@ -192,6 +194,8 @@ void MutationManager::splitFace(Face f, const std::vector<double>& bSplit, Vecto
192194
for (FaceSplitPolicy* policy : faceSplitPolicies) {
193195
policy->afterFaceSplit(newV, bSplit);
194196
}
197+
198+
return newV;
195199
}
196200

197201

@@ -290,7 +294,7 @@ MutationPolicyHandle MutationManager::registerPolicy(MutationPolicy* policyObjec
290294
}
291295

292296
void MutationManager::removePolicy(const MutationPolicyHandle& toRemove) {
293-
// Remove from all lists
297+
// Remove from all lists
294298
removeFromVector(vertexRepositionPolicies, toRemove.policy);
295299
removeFromVector(edgeFlipPolicies, toRemove.policy);
296300
removeFromVector(edgeSplitPolicies, toRemove.policy);

‎src/surface/quadric_error_simplification.cpp‎

Lines changed: 18 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
#include "quadric_error_simplification.h"
1+
#include "geometrycentral/surface/quadric_error_simplification.h"
22

33
namespace geometrycentral {
44
namespace surface {
@@ -32,12 +32,19 @@ Quadric Quadric::operator+=(const Quadric& Q) {
3232
Quadric operator+(const Quadric& Q1, const Quadric& Q2) { return Quadric(Q1.A + Q2.A, Q1.b + Q2.b, Q1.c + Q2.c); }
3333

3434
void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, double tol) {
35-
MutationManager mm(mesh);
35+
MutationManager mm(mesh, geo);
3636
quadricErrorSimplify(mesh, geo, tol, mm);
3737
}
3838

3939
void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo, double tol, MutationManager& mm) {
40-
VertexData<QES::Quadric> Q(mesh, QES::Quadric());
40+
auto toEigen = [](Vector3 v) -> Eigen::Vector3d {
41+
Eigen::Vector3d ret;
42+
ret << v.x, v.y, v.z;
43+
return ret;
44+
};
45+
auto fromEigen = [](Eigen::Vector3d v) -> Vector3 { return Vector3{v(0), v(1), v(2)}; };
46+
47+
VertexData<Quadric> Q(mesh, Quadric());
4148

4249
geo.requireFaceNormals();
4350

@@ -48,7 +55,7 @@ void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo
4855
Eigen::Vector3d q = toEigen(geo.inputVertexPositions[v]);
4956
double d = -n.dot(q);
5057

51-
Q[v] += QES::Quadric(M, d * n, d * d);
58+
Q[v] += Quadric(M, d * n, d * d);
5259
}
5360
}
5461

@@ -59,7 +66,7 @@ void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo
5966
std::priority_queue<PotentialEdge, std::vector<PotentialEdge>, decltype(cmp)> edgesToCheck(cmp);
6067

6168
for (Edge e : mesh.edges()) {
62-
QES::Quadric Qe = Q[e.src()] + Q[e.dst()];
69+
Quadric Qe = Q[e.halfedge().tailVertex()] + Q[e.halfedge().tipVertex()];
6370
Eigen::Vector3d q = Qe.optimalPoint();
6471
double cost = Qe.cost(q);
6572
edgesToCheck.push(std::make_tuple(cost, e));
@@ -76,29 +83,23 @@ void quadricErrorSimplify(ManifoldSurfaceMesh& mesh, VertexPositionGeometry& geo
7683
Edge e = std::get<1>(best);
7784
if (!e.isDead()) {
7885

86+
Vertex v1 = e.halfedge().tailVertex();
87+
Vertex v2 = e.halfedge().tipVertex();
88+
7989
// Get edge quadric
80-
QES::Quadric Qe(Q[e.src()], Q[e.dst()]);
90+
Quadric Qe(Q[v1], Q[v2]);
8191
Eigen::Vector3d q = Qe.optimalPoint();
8292

8393
// If either vertex has been collapsed since the edge was pushed
8494
// onto the queue, the old cost was wrong. In that case, give up
8595
if (abs(cost - Qe.cost(q)) > 1e-8) continue;
8696

87-
Vertex v;
88-
try {
89-
v = mm.collapseEdge(e);
90-
} catch (std::runtime_error& e) {
91-
}
92-
93-
// v == Vertex() if collapseEdge threw an exception, or if
94-
// collapseEdge failed and returned Vertex(). In either case, we
95-
// give up
97+
Vertex v = mm.collapseEdge(e, fromEigen(q));
9698
if (v == Vertex()) continue;
97-
geo.vertexPositions[v] = fromEigen(q);
9899
Q[v] = Qe;
99100

100101
for (Edge f : v.adjacentEdges()) {
101-
QES::Quadric Qf(Q[f.src()], Q[f.dst()]);
102+
Quadric Qf(Q[f.halfedge().tailVertex()], Q[f.halfedge().tipVertex()]);
102103
Eigen::Vector3d q = Qf.optimalPoint();
103104
double cost = Qf.cost(q);
104105
edgesToCheck.push(std::make_tuple(cost, f));

0 commit comments

Comments
 (0)