Euclid
Geometry Processing and Shape Analysis in C++
TriMeshGeometry

Discrete (differential) geometry on a triangle mesh. More...

Enumerations

enum  Euclid::VertexNormal { Euclid::VertexNormal::uniform, Euclid::VertexNormal::face_area, Euclid::VertexNormal::incident_angle }
 Strategies to compute vertex normal. More...
 
enum  Euclid::VertexArea { Euclid::VertexArea::barycentric, Euclid::VertexArea::voronoi, Euclid::VertexArea::mixed_voronoi }
 Strategies to compute vertex area. More...
 

Functions

template<typename Mesh >
Vector_3_t< Mesh > Euclid::vertex_normal (vertex_t< Mesh > v, const Mesh &mesh, const VertexNormal &weight=VertexNormal::incident_angle)
 Normal vector of a vertex on the mesh. More...
 
template<typename Mesh >
Vector_3_t< Mesh > Euclid::vertex_normal (vertex_t< Mesh > v, const Mesh &mesh, const std::vector< Vector_3_t< Mesh >> &face_normals, const VertexNormal &weight=VertexNormal::incident_angle)
 Normal vector of a vertex on the mesh. More...
 
template<typename Mesh >
std::vector< Vector_3_t< Mesh > > Euclid::vertex_normals (const Mesh &mesh, const std::vector< Vector_3_t< Mesh >> &face_normals, const VertexNormal &weight=VertexNormal::incident_angle)
 Normal vectors of all vertices on the mesh. More...
 
template<typename Mesh >
FT_t< Mesh > Euclid::vertex_area (vertex_t< Mesh > v, const Mesh &mesh, const VertexArea &method=VertexArea::mixed_voronoi)
 Area of a vertex on the mesh. More...
 
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::vertex_areas (const Mesh &mesh, const VertexArea &method=VertexArea::mixed_voronoi)
 Areas of all vertices on the mesh. More...
 
template<typename Mesh >
FT_t< Mesh > Euclid::edge_length (halfedge_t< Mesh > he, const Mesh &mesh)
 Edge length.
 
template<typename Mesh >
FT_t< Mesh > Euclid::edge_length (edge_t< Mesh > e, const Mesh &mesh)
 Edge length.
 
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::edge_lengths (const Mesh &mesh)
 Edge lengths.
 
template<typename Mesh >
FT_t< Mesh > Euclid::squared_edge_length (halfedge_t< Mesh > he, const Mesh &mesh)
 Squared edge length.
 
template<typename Mesh >
FT_t< Mesh > Euclid::squared_edge_length (edge_t< Mesh > e, const Mesh &mesh)
 Squared edge length.
 
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::squared_edge_lengths (const Mesh &mesh)
 Squared edge lengths.
 
template<typename Mesh >
FT_t< Mesh > Euclid::dihedral_angle (edge_t< Mesh > e, const Mesh &mesh)
 Dihedral angle between adjacent faces.
 
template<typename Mesh >
FT_t< Mesh > Euclid::corner_angle (halfedge_t< Mesh > h, const Mesh &mesh)
 Corner angle.
 
template<typename Mesh >
Vector_3_t< Mesh > Euclid::face_normal (face_t< Mesh > f, const Mesh &mesh)
 Normal of a face on the mesh.
 
template<typename Mesh >
std::vector< Vector_3_t< Mesh > > Euclid::face_normals (const Mesh &mesh)
 Normals of all faces on the mesh.
 
template<typename Mesh >
FT_t< Mesh > Euclid::face_area (face_t< Mesh > f, const Mesh &mesh)
 Area of a face on the mesh.
 
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::face_areas (const Mesh &mesh)
 Areas of all faces on the mesh.
 
template<typename Mesh >
Point_3_t< Mesh > Euclid::barycenter (face_t< Mesh > f, const Mesh &mesh)
 Barycenter/centroid of a face on the mesh.
 
template<typename Mesh >
std::vector< Point_3_t< Mesh > > Euclid::barycenters (const Mesh &mesh)
 Barycenters/centroids of all faces on the mesh.
 
template<typename Mesh >
FT_t< Mesh > Euclid::gaussian_curvature (vertex_t< Mesh > v, const Mesh &mesh)
 Gaussian curvature of a vertex on the mesh. More...
 
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::gaussian_curvatures (const Mesh &mesh)
 Gaussian curvatures of all vertices on the mesh. More...
 
template<typename Mesh >
std::tuple< Eigen::SparseMatrix< FT_t< Mesh > >, Eigen::SparseMatrix< FT_t< Mesh > > > Euclid::adjacency_matrix (const Mesh &mesh)
 Adjacency matrix of the mesh. More...
 
template<typename Mesh >
Eigen::SparseMatrix< FT_t< Mesh > > Euclid::graph_laplacian_matrix (const Mesh &mesh)
 Graph laplacian matrix of the mesh.
 
template<typename Mesh >
FT_t< Mesh > Euclid::cotangent_weight (halfedge_t< Mesh > he, const Mesh &mesh)
 Cotangent weight associated with an edge.
 
template<typename Mesh >
FT_t< Mesh > Euclid::cotangent_weight (edge_t< Mesh > e, const Mesh &mesh)
 Cotangent weight associated with an edge.
 
template<typename Mesh >
Eigen::SparseMatrix< FT_t< Mesh > > Euclid::cotangent_matrix (const Mesh &mesh)
 Cotangent matrix of the mesh. More...
 
template<typename Mesh >
Eigen::SparseMatrix< FT_t< Mesh > > Euclid::mass_matrix (const Mesh &mesh, const VertexArea &method=VertexArea::mixed_voronoi)
 Mass matrix of the mesh. More...
 

Detailed Description

This package contains functions to compute geometric properties and differential operators on a triangle mesh.

Note that, edge_length() and squared_edge_length() are partial specialized for CGAL::Dual<Mesh> which approximate the length with barycenters.

References

[1] Botsch, M., Kobbelt L., Pauly M., et al. Polygon Mesh Processing. CRC press, 2010.

[2] Crane, K. Discrete Differential Geometry: an Applied Introduction.

[3] Meyer, M., Schroder, P., Barrr, A. H. Discrete Differential-Geometry Operators for Triangulated 2-Manifolds. Visualization and Mathematics III, 2003.

[4] Zhang, H., Van Kaick, O., Dyer, R. Spectral Mesh Processing. Computer Graphics Forum, 2010.

Enumeration Type Documentation

enum Euclid::VertexArea
strong
See also
vertex_area()
http://www.alecjacobson.com/weblog/?p=1146
Enumerator
barycentric 

Barycentric cell.

Uses the barycenters of incident triangles as the vertices of the local cell. In other words the area is 1/3 of the area of the incident triangles. This area is fast and positive but sensitive to mesh triangulations.

voronoi 

Voronoi cell.

Uses the circumcenters of incident triangles as the vertices of the local cell, which gives tighter error bounds. However the circumcenters might appear outside the triangles, and the signed area might be negative.

mixed_voronoi 

A Mixed voronoi cell.

Ensures the vertices of the local voronoi cell are within the triangles by replacing circumcenters which are outside of the triangles by the midpoints of the opposite edges of the obtuse vertices. The signed area is guaranteed to be positive, but not as smooth as pure voronoi cells.

enum Euclid::VertexNormal
strong
See also
vertex_normal()
Enumerator
uniform 

Assign equal weights to neighboring face normals.

This is the fastest yet most inaccurate strategy.

face_area 

Uses face area as weight for face normals.

This gives good balance of speed and quality.

incident_angle 

Uses the incident angle of a face to the vertex as weight for the face normal.

It often gives the best result but is computationally expensive.

Function Documentation

template<typename Mesh >
std::tuple< Eigen::SparseMatrix< FT_t< Mesh > >, Eigen::SparseMatrix< FT_t< Mesh > > > Euclid::adjacency_matrix ( const Mesh &  mesh)

Return the unweighted adjacency matrix as well as the degree matrix of a mesh.

template<typename Mesh >
Eigen::SparseMatrix< FT_t< Mesh > > Euclid::cotangent_matrix ( const Mesh &  mesh)

The cotangent matrix is a discretization of the Laplacian-Beltrami operator. It is a symemtric, positive semi-definitive matrix. Assume that the cotangent matrix is \(L\), and the mass matrix is \(A\), here're some variants of discrete LB operator.

  • \(L\), alone is not invariant to meshing
  • \(A^{-1}L\), accurate discretization but not symmetric
  • \(A^{-1/2}LA^{-1/2}\), symmetric and mesh independent

Note

It's kinda arbitary in the literature when people refer to the Laplacian matrix, especially the sign. In this particular implementation, the diagonal elements are positive and the others are negative, thus forming a positive smei-definitive matrix.

See also
mass_matrix
template<typename Mesh >
FT_t< Mesh > Euclid::gaussian_curvature ( vertex_t< Mesh >  v,
const Mesh &  mesh 
)

Discrete Gaussian curvature using the angle deficit method.

template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::gaussian_curvatures ( const Mesh &  mesh)

Discrete Gaussian curvature using the angle deficit method.

template<typename Mesh >
Eigen::SparseMatrix< FT_t< Mesh > > Euclid::mass_matrix ( const Mesh &  mesh,
const VertexArea method = VertexArea::mixed_voronoi 
)

The mass matrix is simply the vertex areas of all the vertices of a mesh written in a diagonal matrix, which serves as the local averaging region commonly used in discrete differential geometry.

See also
VertexArea, vertex_area, cotangent_matrix
template<typename Mesh >
FT_t< Mesh > Euclid::vertex_area ( vertex_t< Mesh >  v,
const Mesh &  mesh,
const VertexArea method = VertexArea::mixed_voronoi 
)

Compute a small area around the vertex as the local average of the integral of the differential area dA. In the discrete settings, this involves constructing a small cell around the vertex and use its area as the local averaging region. Choose one type of region as specified in VertexArea. Only 1-ring neighborhood is considered in this function.

See also
VertexArea
template<typename Mesh >
std::vector< FT_t< Mesh > > Euclid::vertex_areas ( const Mesh &  mesh,
const VertexArea method = VertexArea::mixed_voronoi 
)

Compute a small area around the vertex as the local average of the integral of the differential area dA. In the discrete settings, this involves constructing a small cell around the vertex and use its area as the local averaging region. Choose one type of region as specified in VertexArea. Only 1-ring neighborhood is considered in this function.

See also
VertexArea
template<typename Mesh >
Vector_3_t< Mesh > Euclid::vertex_normal ( vertex_t< Mesh >  v,
const Mesh &  mesh,
const VertexNormal weight = VertexNormal::incident_angle 
)

Compute vertex normal from the incident face normals. Choose a weighting strategy of face normals from VertexNormal.

Note

This function will compute face normals for all the incident triangles of vertex v. To avoid wasting computation, provide face normals using the other overload method.

See also
VertexNormal
template<typename Mesh >
Vector_3_t< Mesh > Euclid::vertex_normal ( vertex_t< Mesh >  v,
const Mesh &  mesh,
const std::vector< Vector_3_t< Mesh >> &  face_normals,
const VertexNormal weight = VertexNormal::incident_angle 
)

Compute vertex normal from the incident face normals. Choose a weighting strategy of face normals from VertexNormal.

See also
VertexNormal
template<typename Mesh >
std::vector< Vector_3_t< Mesh > > Euclid::vertex_normals ( const Mesh &  mesh,
const std::vector< Vector_3_t< Mesh >> &  face_normals,
const VertexNormal weight = VertexNormal::incident_angle 
)

Compute vertex normal from the incident face normals. Choose a weighting strategy of face normals from VertexNormal.

See also
VertexNormal