LArSoft  v06_85_00
Liquid Argon Software toolkit - http://larsoft.org/
DCEL.h
Go to the documentation of this file.
1 
10 #ifndef DCEL2D_h
11 #define DCEL2D_h
12 
13 // std includes
14 #include <vector>
15 #include <list>
16 #include <algorithm>
17 
18 // Eigen
19 #include <Eigen/Dense>
20 //------------------------------------------------------------------------------------------------------------------------------------------
21 
22 namespace reco
23 {
24  class ClusterHit3D;
25 }
26 
27 namespace dcel2d
28 {
29 class HalfEdge; // Forward declaration
30 
34 using Point = std::tuple<double,double,const reco::ClusterHit3D*>;
35 using PointList = std::list<Point>;
36 using Coords = Eigen::Vector3f; //std::pair<double,double>;
37 
38 class Vertex
39 {
45 public:
49  Vertex() : fCoords(0.,0.,0.), fHalfEdge(NULL) {}
50 
51  Vertex(const double* coords, HalfEdge* half) : fHalfEdge(half)
52  {
53  setCoords(coords);
54  }
55 
56  Vertex(const Coords& coords, HalfEdge* half) : fCoords(coords), fHalfEdge(half)
57  {}
58 
59  ~Vertex() {}
60 
61  const Coords& getCoords() const {return fCoords;}
62  const HalfEdge* getHalfEdge() const {return fHalfEdge;}
63 
64  void setCoords(const double* coords)
65  {
66  fCoords[0] = coords[0];
67  fCoords[1] = coords[1];
68  fCoords[2] = coords[2];
69  }
70 
71  void setCoords(const Coords& coords) {fCoords = coords;}
72 
73  void setHalfEdge(HalfEdge* half) {fHalfEdge = half;}
74 
75 private:
76  Coords fCoords; // x,y coordinates of this vertex
77  HalfEdge* fHalfEdge; // pointer to one of the half edges
78 };
79 
80 class Face
81 {
88 public:
92 // Face() : fHalfEdge(NULL), fPoint(0.,0.,NULL) {}
93 
94  Face(HalfEdge* half, const Coords& coords, const reco::ClusterHit3D* clusterHit3D) :
95  fHalfEdge(half),
96  fConvexHull(false),
97  fCoords(coords),
98  fFaceArea(0.),
99  fClusterHit3D(clusterHit3D)
100  {}
101 
102  ~Face() {}
103 
104  const HalfEdge* getHalfEdge() const {return fHalfEdge;}
105  const bool onConvexHull() const {return fConvexHull;}
106  const Coords& getCoords() const {return fCoords;}
107  const double getFaceArea() const {return fFaceArea;}
108  const reco::ClusterHit3D* getClusterHit3D() const {return fClusterHit3D;}
109 
110  void setHalfEdge(HalfEdge* half) {fHalfEdge = half;}
111  void setOnConvexHull() {fConvexHull = true;}
112  void setFaceArea(double area) {fFaceArea = area;}
113 
114 private:
115  HalfEdge* fHalfEdge; // pointer to one of the half edges
116  mutable bool fConvexHull; // This face on convex hull
117  Coords fCoords; // projected coordinates of associated point
118  double fFaceArea; // The area of the face once constructed
119  const reco::ClusterHit3D* fClusterHit3D; // The physical 3D hit this corresponds to
120 };
121 
122 class HalfEdge
123 {
131 public:
136  m_targetVertex(NULL),
137  m_face(NULL),
138  m_twinHalfEdge(NULL),
139  m_nextHalfEdge(NULL),
140  m_lastHalfEdge(NULL)
141  {}
142 
144  Face* face,
145  HalfEdge* twin,
146  HalfEdge* next,
147  HalfEdge* last) :
148  m_targetVertex(vertex),
149  m_face(face),
150  m_twinHalfEdge(twin),
151  m_nextHalfEdge(next),
152  m_lastHalfEdge(last)
153  {}
154 
156 
157  Vertex* getTargetVertex() const {return m_targetVertex;}
158  Face* getFace() const {return m_face;}
159  HalfEdge* getTwinHalfEdge() const {return m_twinHalfEdge;}
160  HalfEdge* getNextHalfEdge() const {return m_nextHalfEdge;}
161  HalfEdge* getLastHalfEdge() const {return m_lastHalfEdge;}
162 
163  void setTargetVertex(Vertex* vertex) {m_targetVertex = vertex;}
164  void setFace(Face* face) {m_face = face;}
165  void setTwinHalfEdge(HalfEdge* twin) {m_twinHalfEdge = twin;}
166  void setNextHalfEdge(HalfEdge* next) {m_nextHalfEdge = next;}
167  void setLastHalfEdge(HalfEdge* last) {m_lastHalfEdge = last;}
168 
169 private:
170  Vertex* m_targetVertex; // Pointer to the vertex we point to
171  Face* m_face; // Pointer to the face we are associated with
172  HalfEdge* m_twinHalfEdge; // Pointer to the twin half edge (pointing opposite)
173  HalfEdge* m_nextHalfEdge; // Pointer ot the next half edge
174  HalfEdge* m_lastHalfEdge; // Pointer to the previous half edge
175 };
176 
177 // Define containers to hold the above objects
178 using VertexList = std::list<Vertex>;
179 using FaceList = std::list<Face>;
180 using HalfEdgeList = std::list<HalfEdge>;
181 
182 } // namespace lar_cluster3d
183 #endif
void setFaceArea(double area)
Definition: DCEL.h:112
const reco::ClusterHit3D * getClusterHit3D() const
Definition: DCEL.h:108
Vertex()
Vertex class definition for use in a doubly connected edge list a Vertex will contain the coordinates...
Definition: DCEL.h:49
HalfEdge * getLastHalfEdge() const
Definition: DCEL.h:161
void setCoords(const Coords &coords)
Definition: DCEL.h:71
Definition: DCEL.h:27
HalfEdge * getNextHalfEdge() const
Definition: DCEL.h:160
std::list< HalfEdge > HalfEdgeList
Definition: DCEL.h:180
void setHalfEdge(HalfEdge *half)
Definition: DCEL.h:73
const bool onConvexHull() const
Definition: DCEL.h:105
Face(HalfEdge *half, const Coords &coords, const reco::ClusterHit3D *clusterHit3D)
Face class definition for use in a doubly connected edge list A Face represents the area enclosed by ...
Definition: DCEL.h:94
HalfEdge * getTwinHalfEdge() const
Definition: DCEL.h:159
Face * getFace() const
Definition: DCEL.h:158
HalfEdge * fHalfEdge
Definition: DCEL.h:115
HalfEdge * m_lastHalfEdge
Definition: DCEL.h:174
HalfEdge * m_nextHalfEdge
Definition: DCEL.h:173
Vertex(const double *coords, HalfEdge *half)
Definition: DCEL.h:51
double fFaceArea
Definition: DCEL.h:118
Vertex(const Coords &coords, HalfEdge *half)
Definition: DCEL.h:56
Vertex * m_targetVertex
Definition: DCEL.h:170
Coords fCoords
Definition: DCEL.h:76
void setLastHalfEdge(HalfEdge *last)
Definition: DCEL.h:167
void setTwinHalfEdge(HalfEdge *twin)
Definition: DCEL.h:165
HalfEdge()
HalfEdge class definition for use in a doubly connected edge list The half edge class represents one ...
Definition: DCEL.h:135
std::list< Face > FaceList
Definition: DCEL.h:179
std::tuple< double, double, const reco::ClusterHit3D * > Point
Definitions used by the VoronoiDiagram algorithm.
Definition: DCEL.h:34
HalfEdge * fHalfEdge
Definition: DCEL.h:77
Coords fCoords
Definition: DCEL.h:117
const Coords & getCoords() const
Definition: DCEL.h:106
void setFace(Face *face)
Definition: DCEL.h:164
bool fConvexHull
Definition: DCEL.h:116
const HalfEdge * getHalfEdge() const
Definition: DCEL.h:104
void setNextHalfEdge(HalfEdge *next)
Definition: DCEL.h:166
void setCoords(const double *coords)
Definition: DCEL.h:64
HalfEdge(Vertex *vertex, Face *face, HalfEdge *twin, HalfEdge *next, HalfEdge *last)
Definition: DCEL.h:143
void setHalfEdge(HalfEdge *half)
Definition: DCEL.h:110
Definition of utility objects for use in the 3D clustering for LArSoft.
Definition: Cluster3D.cxx:17
Face * m_face
Definition: DCEL.h:171
HalfEdge * m_twinHalfEdge
Definition: DCEL.h:172
std::list< Point > PointList
Definition: DCEL.h:35
Vertex * getTargetVertex() const
Definition: DCEL.h:157
Eigen::Vector3f Coords
Definition: DCEL.h:36
void setOnConvexHull()
Definition: DCEL.h:111
const HalfEdge * getHalfEdge() const
Definition: DCEL.h:62
const double getFaceArea() const
Definition: DCEL.h:107
const Coords & getCoords() const
Definition: DCEL.h:61
const reco::ClusterHit3D * fClusterHit3D
Definition: DCEL.h:119
std::list< Vertex > VertexList
Definition: DCEL.h:178
void setTargetVertex(Vertex *vertex)
Definition: DCEL.h:163
vertex reconstruction