DGtal  1.5.beta
SurfaceMesh.h
1 
17 #pragma once
18 
31 #if defined(SurfaceMesh_RECURSES)
32 #error Recursive header files inclusion detected in SurfaceMesh.h
33 #else // defined(SurfaceMesh_RECURSES)
35 #define SurfaceMesh_RECURSES
36 
37 #if !defined SurfaceMesh_h
39 #define SurfaceMesh_h
40 
42 // Inclusions
43 #include <iostream>
44 #include <sstream>
45 #include <string>
46 #include "DGtal/base/Common.h"
47 #include "DGtal/base/IntegerSequenceIterator.h"
48 #include "DGtal/helpers/StdDefs.h"
49 
50 namespace DGtal
51 {
53  // template class SurfaceMesh
90  template < typename TRealPoint, typename TRealVector >
91  struct SurfaceMesh
92  {
93  typedef TRealPoint RealPoint;
94  typedef TRealVector RealVector;
96 
97  static const Dimension dimension = RealPoint::dimension;
99 
100  typedef typename RealVector::Component Scalar;
101  typedef std::vector<Scalar> Scalars;
103  typedef std::size_t Size;
105  typedef std::size_t Index;
106  typedef Index Face;
107  typedef Index Edge;
108  typedef Index Vertex;
109  typedef std::pair< Edge, Scalar > WeightedEdge;
110  typedef std::pair< Face, Scalar > WeightedFace;
112  typedef std::vector< Vertex > Vertices;
114  typedef std::vector< Edge > Edges;
115  typedef std::vector< WeightedEdge > WeightedEdges;
116  typedef std::vector< Face > Faces;
117  typedef std::vector< WeightedFace > WeightedFaces;
118  typedef std::pair< Vertex, Vertex > VertexPair;
119 
120  // Required by CUndirectedSimpleLocalGraph
121  typedef std::set<Vertex> VertexSet;
122  template <typename Value> struct VertexMap {
123  typedef typename std::map<Vertex, Value> Type;
124  };
125 
126  // Required by CUndirectedSimpleGraph
127 
130 
131  //---------------------------------------------------------------------------
132  public:
135 
137  ~SurfaceMesh() = default;
147  SurfaceMesh() = default;
150  SurfaceMesh( const Self& other ) = default;
153  SurfaceMesh( Self&& other ) = default;
157  Self& operator=( const Self& other ) = default;
158 
178  template <typename RealPointIterator, typename VerticesIterator>
179  SurfaceMesh( RealPointIterator itPos, RealPointIterator itPosEnd,
180  VerticesIterator itVertices, VerticesIterator itVerticesEnd );
181 
203  template <typename RealPointIterator, typename VerticesIterator>
204  bool init( RealPointIterator itPos, RealPointIterator itPosEnd,
205  VerticesIterator itVertices, VerticesIterator itVerticesEnd );
206 
208  void clear();
209 
211 
212  //---------------------------------------------------------------------------
213  public:
216 
219  template <typename RealVectorIterator>
220  bool setVertexNormals( RealVectorIterator itN, RealVectorIterator itNEnd );
221 
224  template <typename RealVectorIterator>
225  bool setFaceNormals( RealVectorIterator itN, RealVectorIterator itNEnd );
226 
232 
239 
244 
249 
255 
259  template <typename AnyRing>
261  ( const std::vector<AnyRing>& vvalues ) const;
262 
266  template <typename AnyRing>
268  ( const std::vector<AnyRing>& fvalues ) const;
269 
273  ( const std::vector<RealVector>& vuvectors ) const;
274 
278  ( const std::vector<RealVector>& fuvectors ) const;
279 
281 
282  //---------------------------------------------------------------------------
283  public:
286 
288  Size nbVertices() const
289  { return myIncidentFaces.size(); }
290 
292  Size nbEdges() const
293  { return myEdgeVertices.size(); }
294 
296  Size nbFaces() const
297  { return myIncidentVertices.size(); }
298 
302  long Euler() const
303  { return nbVertices() - nbEdges() + nbFaces(); }
304 
310  Edge makeEdge( Vertex i, Vertex j ) const;
311 
315  const Vertices& incidentVertices( Face f ) const
316  { return myIncidentVertices[ f ]; }
317 
321  const Faces& incidentFaces( Vertex v ) const
322  { return myIncidentFaces[ v ]; }
323 
326  const Faces& neighborFaces( Face f ) const
327  { return myNeighborFaces[ f ]; }
328 
331  const Vertices& neighborVertices( Vertex v ) const
332  { return myNeighborVertices[ v ]; }
333 
339  { return myEdgeVertices[ e ]; }
340 
344  const Faces& edgeFaces( Edge e ) const
345  { return myEdgeFaces[ e ]; }
346 
355  const Faces& edgeRightFaces( Edge e ) const
356  { return myEdgeRightFaces[ e ]; }
357 
366  const Faces& edgeLeftFaces( Edge e ) const
367  { return myEdgeLeftFaces[ e ]; }
368 
371  const std::vector< Vertices >& allIncidentVertices() const
372  { return myIncidentVertices; }
373 
376  const std::vector< Faces >& allIncidentFaces() const
377  { return myIncidentFaces; }
378 
380  const std::vector< Faces >& allNeighborFaces() const
381  { return myNeighborFaces; }
382 
384  const std::vector< Vertices >& allNeighborVertices() const
385  { return myNeighborVertices; }
386 
391  const std::vector< VertexPair >& allEdgeVertices() const
392  { return myEdgeVertices; }
393 
396  const std::vector< Faces >& allEdgeFaces() const
397  { return myEdgeFaces; }
398 
406  const std::vector< Faces >& allEdgeRightFaces() const
407  { return myEdgeRightFaces; }
408 
416  const std::vector< Faces >& allEdgeLeftFaces() const
417  { return myEdgeLeftFaces; }
418 
420 
421  //---------------------------------------------------------------------------
422  public:
425 
448 
456  bool isBoundariesManifold( bool checkClosed = true ) const
457  {
458  // computes unordered list of boundary vertices
459  std::map<Vertex,Vertices> adjacent;
460  auto MBE = this->computeManifoldBoundaryEdges();
461  if ( MBE.size() == 0 ) return false;
462 
463  for (auto e : MBE)
464  {
465  auto ij = this->edgeVertices(e);
466  adjacent[ ij.first ].push_back( ij.second );
467  if ( adjacent[ ij.first ] .size() > 2 ) return false;
468  adjacent[ ij.second ].push_back( ij.first );
469  if ( adjacent[ ij.second ].size() > 2 ) return false;
470  }
471  //we may check if all curves are closed.
472  if ( checkClosed )
473  for ( const auto &adj : adjacent )
474  if ( adj.second.size() != 2 ) return false;
475  return true;
476  }
477 
483  std::vector<Vertices> computeManifoldBoundaryChains() const
484  {
485  std::vector<Vertices> boundaries;
486  Vertices boundary;
487  auto MBE = this->computeManifoldBoundaryEdges();
488  std::map<Vertex,bool> visited;
489  std::map<Vertex,Vertices> adjacent;
490 
491  ASSERT_MSG(MBE.size()>0,"The surface mesh must have boundary edges");
492  ASSERT_MSG(this->isBoundariesManifold(), "The surface mesh mush have manifold boundaries");
493 
494  //Coompute adjecency relationships
495  for (auto e : MBE)
496  {
497  auto ij = this->edgeVertices(e);
498  visited[ij.first] = false;
499  visited[ij.second] = false;
500  adjacent[ij.first].push_back(ij.second);
501  adjacent[ij.second].push_back(ij.first);
502  }
503 
504  auto boundary_it = visited.begin();
505 
506  while(boundary_it != visited.end() )
507  {
508  Vertex first = (*boundary_it).first;
509  visited[first] = true;
510  boundary.clear();
511  boundary.push_back(first);
512 
513  Vertex current = first;
514  size_t nb_iter = 0;
515  bool pushed=false;
516 
517  while ((!pushed) && (nb_iter < MBE.size()*2))
518  {
519  bool ok = false;
520  for (auto other : adjacent[current])
521  if (!visited[other])
522  {
523  boundary.push_back(other);
524  current = other;
525  visited[other] = true;
526  ok = true;
527  break;
528  }
529  if (!ok)
530  {
531  //all neighboors are visited
532  for (auto other : adjacent[current])
533  if (other == first)
534  {
535  boundaries.push_back(boundary);
536  pushed = true;
537  break;
538  }
539  //if first vertex isn't found then this chain is not
540  //homeomorphic to a circle, hence isn't added to boundaries
541  }
542  nb_iter++;
543  }
544  boundary_it = std::find_if(visited.begin(), visited.end(),
545  []
546  (std::pair<Vertex,bool> x){return !x.second;});
547  }
548  return boundaries;
549  }
551 
552  // ----------------------- Undirected simple graph services ----------------------
553  public:
556 
560  Size size() const
561  { return nbVertices(); }
562 
570  { return 8; }
571 
577  Size degree( const Vertex & v ) const
578  { return myNeighborVertices[ v ].size(); }
579 
590  template <typename OutputIterator>
591  void
592  writeNeighbors( OutputIterator &it ,
593  const Vertex & v ) const
594  {
595  for ( auto&& nv : myNeighborVertices[ v ] )
596  *it++ = nv;
597  }
598 
615  template <typename OutputIterator, typename VertexPredicate>
616  void
617  writeNeighbors( OutputIterator &it ,
618  const Vertex & v,
619  const VertexPredicate & pred) const
620  {
621  for ( auto&& nv : myNeighborVertices[ v ] )
622  if ( pred( nv ) ) *it++ = nv;
623  }
624 
627  { return ConstIterator( 0 ); }
628 
631  { return ConstIterator( nbVertices() ); }
632 
634 
635  //---------------------------------------------------------------------------
636  public:
639 
641  const std::vector< RealPoint >& positions() const
642  { return myPositions; }
643 
648  { return myPositions[ v ]; }
649 
653  const RealPoint& position( Vertex v ) const
654  { return myPositions[ v ]; }
655 
657  const std::vector< RealVector >& vertexNormals() const
658  { return myVertexNormals; }
659 
661  std::vector< RealVector >& vertexNormals()
662  { return myVertexNormals; }
663 
668  { return myVertexNormals[ v ]; }
669 
673  const RealVector& vertexNormal( Vertex v ) const
674  { return myVertexNormals[ v ]; }
675 
677  const std::vector< RealVector >& faceNormals() const
678  { return myFaceNormals; }
679 
681  std::vector< RealVector >& faceNormals()
682  { return myFaceNormals; }
683 
688  { return myFaceNormals[ f ]; }
689 
693  const RealVector& faceNormal( Face f ) const
694  { return myFaceNormals[ f ]; }
695 
698 
702  Scalar distance(const Vertex i, const Vertex j) const
703  {
704  //unrolling for compiler optimization
705  const auto p=this->myPositions[ i ];
706  const auto q=this->myPositions[ j ];
707  return std::sqrt( (p[0]-q[0])*(p[0]-q[0]) + (p[1]-q[1])*(p[1]-q[1])+ (p[2]-q[2])*(p[2]-q[2]));
708  }
709 
713  Scalar localWindow( Face f ) const;
714 
719 
724 
728 
732 
735  Scalar faceArea( Index f ) const;
736 
741 
763 
786 
802  std::tuple< Vertices, WeightedEdges, WeightedFaces >
804 
821  std::tuple< Vertices, WeightedEdges, WeightedFaces >
823 
838 
853 
864 
866 
867  //---------------------------------------------------------------------------
868  public:
871 
882  bool isFlippable( const Edge e ) const;
883 
887  VertexPair otherDiagonal( const Edge e ) const;
888 
933  void flip( const Edge e, bool recompute_face_normals = false );
934 
936 
937  //---------------------------------------------------------------------------
938  public:
941 
945  void computeEdges();
946 
948 
949  // ----------------------- Interface --------------------------------------
950  public:
951 
956  void selfDisplay ( std::ostream & out ) const;
957 
962  bool isValid() const;
963 
964  // ------------------------- Protected Datas ------------------------------
965  protected:
967  std::vector< Vertices > myIncidentVertices;
969  std::vector< Faces > myIncidentFaces;
971  std::vector< RealPoint > myPositions;
973  std::vector< RealVector > myVertexNormals;
975  std::vector< RealVector > myFaceNormals;
977  std::vector< Faces > myNeighborFaces;
979  std::vector< Vertices > myNeighborVertices;
981  std::vector< VertexPair > myEdgeVertices;
983  std::map< VertexPair,Edge > myVertexPairEdge;
985  std::vector< Faces > myEdgeFaces;
991  std::vector< Faces > myEdgeRightFaces;
997  std::vector< Faces > myEdgeLeftFaces;
998 
999  // ------------------------- Private Datas --------------------------------
1000  private:
1001 
1002 
1003  // ------------------------- Internals ------------------------------------
1004  protected:
1005 
1009  void removeIndex( std::vector< Index >& v, Index i )
1010  {
1011  const std::size_t n = v.size();
1012  for ( std::size_t j = 0; j < n; j++ )
1013  if ( v[ j ] == i )
1014  {
1015  std::swap( v[ j ], v.back() );
1016  v.resize( n - 1 );
1017  return;
1018  }
1019  trace.error() << "[SurfaceMesh::removeIndex] Index " << i
1020  << " is not in vector:";
1021  for ( auto e : v ) std::cerr << " " << e;
1022  std::cerr << std::endl;
1023  }
1024 
1029  void replaceIndex( std::vector< Index >& v, Index i, Index ri )
1030  {
1031  const std::size_t n = v.size();
1032  for ( std::size_t j = 0; j < n; j++ )
1033  if ( v[ j ] == i )
1034  {
1035  v[ j ] = ri;
1036  return;
1037  }
1038  trace.error() << "[SurfaceMesh::replaceIndex] Index " << i
1039  << " (subs=" << ri << ") is not in vector:";
1040  for ( auto e : v ) std::cerr << " " << e;
1041  std::cerr << std::endl;
1042  }
1043 
1047  void addIndex( std::vector< Index >& v, Index i )
1048  {
1049  v.push_back( i );
1050  }
1051 
1052 
1054  static Scalar rand01()
1055  { return (Scalar) rand() / (Scalar) RAND_MAX; }
1056 
1057  }; // end of class SurfaceMesh
1058 
1065  template < typename TRealPoint, typename TRealVector >
1066  std::ostream&
1067  operator<< ( std::ostream & out,
1068  const SurfaceMesh<TRealPoint, TRealVector> & object );
1069 
1070 } // namespace DGtal
1071 
1073 // Includes inline functions.
1074 #include "SurfaceMesh.ih"
1075 // //
1077 
1078 #endif // !defined SurfaceMesh_h
1079 
1080 #undef SurfaceMesh_RECURSES
1081 #endif // else defined(SurfaceMesh_RECURSES)
Aim: It is a simple class that mimics a (non mutable) iterator over integers. You can increment it,...
TEuclideanRing Component
Type for Vector elements.
Definition: PointVector.h:614
std::ostream & error()
DGtal is the top-level namespace which contains all DGtal functions and types.
std::ostream & operator<<(std::ostream &out, const ATu0v1< TKSpace, TLinearAlgebra > &object)
DGtal::uint32_t Dimension
Definition: Common.h:136
Trace trace
Definition: Common.h:153
std::map< Vertex, Value > Type
Definition: SurfaceMesh.h:123
Aim: Represents an embedded mesh as faces and a list of vertices. Vertices may be shared among faces ...
Definition: SurfaceMesh.h:92
std::vector< WeightedEdge > WeightedEdges
Definition: SurfaceMesh.h:115
void writeNeighbors(OutputIterator &it, const Vertex &v, const VertexPredicate &pred) const
Definition: SurfaceMesh.h:617
std::vector< RealVector > & faceNormals()
Definition: SurfaceMesh.h:681
WeightedFaces computeFacesInclusionsInBall(Scalar r, Index f, RealPoint p) const
Size size() const
Definition: SurfaceMesh.h:560
const Faces & neighborFaces(Face f) const
Definition: SurfaceMesh.h:326
std::pair< Vertex, Vertex > VertexPair
Definition: SurfaceMesh.h:118
Scalar edgeInclusionRatio(RealPoint p, Scalar r, Index e) const
Edges computeManifoldInnerUnconsistentEdges() const
RealPoint & position(Vertex v)
Definition: SurfaceMesh.h:647
SurfaceMesh(RealPointIterator itPos, RealPointIterator itPosEnd, VerticesIterator itVertices, VerticesIterator itVerticesEnd)
SurfaceMesh(const Self &other)=default
std::vector< RealVector > computeVertexUnitVectorsFromFaceUnitVectors(const std::vector< RealVector > &fuvectors) const
void perturbateWithAdaptiveUniformRandomNoise(Scalar p)
bool setVertexNormals(RealVectorIterator itN, RealVectorIterator itNEnd)
std::vector< Vertex > Vertices
The type that defines a list/range of vertices (e.g. to define faces)
Definition: SurfaceMesh.h:112
Edges computeManifoldInnerEdges() const
std::vector< Faces > myNeighborFaces
For each face, its range of neighbor faces (no particular order)
Definition: SurfaceMesh.h:977
RealPoint faceCentroid(Index f) const
std::vector< Scalar > Scalars
Definition: SurfaceMesh.h:101
static Scalar rand01()
Definition: SurfaceMesh.h:1054
const std::vector< Faces > & allEdgeRightFaces() const
Definition: SurfaceMesh.h:406
long Euler() const
Definition: SurfaceMesh.h:302
void computeEdges()
Computes edge information.
SurfaceMesh< RealPoint, RealVector > Self
Definition: SurfaceMesh.h:95
const std::vector< RealVector > & faceNormals() const
Definition: SurfaceMesh.h:677
static const Dimension dimension
Definition: SurfaceMesh.h:97
Scalars getMaxWeights(Index v) const
std::vector< WeightedFace > WeightedFaces
Definition: SurfaceMesh.h:117
const Faces & edgeFaces(Edge e) const
Definition: SurfaceMesh.h:344
Scalar localWindow(Face f) const
void perturbateWithUniformRandomNoise(Scalar p)
TRealPoint RealPoint
Definition: SurfaceMesh.h:93
void replaceIndex(std::vector< Index > &v, Index i, Index ri)
Definition: SurfaceMesh.h:1029
const std::vector< RealVector > & vertexNormals() const
Definition: SurfaceMesh.h:657
Scalar faceArea(Index f) const
void addIndex(std::vector< Index > &v, Index i)
Definition: SurfaceMesh.h:1047
std::vector< Faces > myIncidentFaces
For each vertex, its range of incident faces.
Definition: SurfaceMesh.h:969
void computeFaceNormalsFromPositions()
std::vector< Faces > myEdgeRightFaces
Definition: SurfaceMesh.h:991
std::vector< AnyRing > computeFaceValuesFromVertexValues(const std::vector< AnyRing > &vvalues) const
const std::vector< RealPoint > & positions() const
Definition: SurfaceMesh.h:641
bool isBoundariesManifold(bool checkClosed=true) const
Definition: SurfaceMesh.h:456
RealVector & vertexNormal(Vertex v)
Definition: SurfaceMesh.h:667
SurfaceMesh()=default
ConstIterator begin() const
Definition: SurfaceMesh.h:626
std::vector< AnyRing > computeVertexValuesFromFaceValues(const std::vector< AnyRing > &fvalues) const
void computeVertexNormalsFromFaceNormals()
RealVector & faceNormal(Face f)
Definition: SurfaceMesh.h:687
std::size_t Index
The type used for numbering vertices and faces.
Definition: SurfaceMesh.h:105
std::vector< Vertices > computeManifoldBoundaryChains() const
Definition: SurfaceMesh.h:483
const Vertices & incidentVertices(Face f) const
Definition: SurfaceMesh.h:315
~SurfaceMesh()=default
Default destructor.
const Faces & edgeRightFaces(Edge e) const
Definition: SurfaceMesh.h:355
IntegerSequenceIterator< Vertex > ConstIterator
Non mutable iterator for visiting vertices.
Definition: SurfaceMesh.h:129
VertexPair otherDiagonal(const Edge e) const
std::vector< Vertices > myNeighborVertices
For each vertex, its range of neighbor vertices (no particular order)
Definition: SurfaceMesh.h:979
BOOST_STATIC_ASSERT((dimension==3))
const std::vector< VertexPair > & allEdgeVertices() const
Definition: SurfaceMesh.h:391
std::vector< Edge > Edges
The type that defines a list/range of edges.
Definition: SurfaceMesh.h:114
std::pair< Face, Scalar > WeightedFace
Definition: SurfaceMesh.h:110
const std::vector< Faces > & allEdgeFaces() const
Definition: SurfaceMesh.h:396
std::map< VertexPair, Edge > myVertexPairEdge
For each vertex pair, its edge index.
Definition: SurfaceMesh.h:983
const std::vector< Vertices > & allNeighborVertices() const
Definition: SurfaceMesh.h:384
const RealVector & faceNormal(Face f) const
Definition: SurfaceMesh.h:693
std::vector< Faces > myEdgeFaces
For each edge, its faces (one, two, or more if non manifold)
Definition: SurfaceMesh.h:985
Scalar distance(const Vertex i, const Vertex j) const
Definition: SurfaceMesh.h:702
void computeVertexNormalsFromFaceNormalsWithMaxWeights()
std::vector< RealVector > & vertexNormals()
Definition: SurfaceMesh.h:661
std::set< Vertex > VertexSet
Definition: SurfaceMesh.h:121
std::vector< RealVector > computeFaceUnitVectorsFromVertexUnitVectors(const std::vector< RealVector > &vuvectors) const
bool setFaceNormals(RealVectorIterator itN, RealVectorIterator itNEnd)
const RealVector & vertexNormal(Vertex v) const
Definition: SurfaceMesh.h:673
SurfaceMesh(Self &&other)=default
Scalar faceInclusionRatio(RealPoint p, Scalar r, Index f) const
std::size_t Size
The type for counting elements.
Definition: SurfaceMesh.h:103
Self & operator=(const Self &other)=default
Size nbFaces() const
Definition: SurfaceMesh.h:296
std::tuple< Vertices, WeightedEdges, WeightedFaces > computeCellsInclusionsInBall(Scalar r, Index f) const
void computeFaceNormalFromPositions(const Face f)
Size nbEdges() const
Definition: SurfaceMesh.h:292
Edges computeManifoldBoundaryEdges() const
Edges computeManifoldInnerConsistentEdges() const
const std::vector< Faces > & allEdgeLeftFaces() const
Definition: SurfaceMesh.h:416
VertexPair edgeVertices(Edge e) const
Definition: SurfaceMesh.h:338
void flip(const Edge e, bool recompute_face_normals=false)
Size bestCapacity() const
Definition: SurfaceMesh.h:569
WeightedFaces computeFacesInclusionsInBall(Scalar r, Index f) const
const std::vector< Vertices > & allIncidentVertices() const
Definition: SurfaceMesh.h:371
const Faces & incidentFaces(Vertex v) const
Definition: SurfaceMesh.h:321
bool init(RealPointIterator itPos, RealPointIterator itPosEnd, VerticesIterator itVertices, VerticesIterator itVerticesEnd)
bool isValid() const
const Faces & edgeLeftFaces(Edge e) const
Definition: SurfaceMesh.h:366
RealVector::Component Scalar
Definition: SurfaceMesh.h:100
Edge makeEdge(Vertex i, Vertex j) const
Size degree(const Vertex &v) const
Definition: SurfaceMesh.h:577
Size nbVertices() const
Definition: SurfaceMesh.h:288
RealPoint edgeCentroid(Index e) const
void computeFaceNormalsFromVertexNormals()
std::vector< VertexPair > myEdgeVertices
For each edge, its two vertices.
Definition: SurfaceMesh.h:981
void clear()
Clears everything. The object is empty.
Edges computeNonManifoldEdges() const
const Vertices & neighborVertices(Vertex v) const
Definition: SurfaceMesh.h:331
std::vector< RealVector > myFaceNormals
For each face, its normal vector.
Definition: SurfaceMesh.h:975
std::tuple< Vertices, WeightedEdges, WeightedFaces > computeCellsInclusionsInBall(Scalar r, Index f, RealPoint p) const
ConstIterator end() const
Definition: SurfaceMesh.h:630
std::vector< Face > Faces
Definition: SurfaceMesh.h:116
std::pair< Edge, Scalar > WeightedEdge
Definition: SurfaceMesh.h:109
void computeNeighbors()
Computes neighboring information.
std::vector< Vertices > myIncidentVertices
For each face, its range of incident vertices.
Definition: SurfaceMesh.h:967
std::vector< Faces > myEdgeLeftFaces
Definition: SurfaceMesh.h:997
std::vector< RealVector > myVertexNormals
For each vertex, its normal vector.
Definition: SurfaceMesh.h:973
const std::vector< Faces > & allNeighborFaces() const
Definition: SurfaceMesh.h:380
TRealVector RealVector
Definition: SurfaceMesh.h:94
void removeIndex(std::vector< Index > &v, Index i)
Definition: SurfaceMesh.h:1009
std::vector< RealPoint > myPositions
For each vertex, its position.
Definition: SurfaceMesh.h:971
void writeNeighbors(OutputIterator &it, const Vertex &v) const
Definition: SurfaceMesh.h:592
const RealPoint & position(Vertex v) const
Definition: SurfaceMesh.h:653
Scalar averageEdgeLength() const
Scalar vertexInclusionRatio(RealPoint p, Scalar r, Index v) const
bool isFlippable(const Edge e) const
const std::vector< Faces > & allIncidentFaces() const
Definition: SurfaceMesh.h:376
void selfDisplay(std::ostream &out) const