FastJet  3.0.6
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Voronoi.hh
1 #ifndef __FASTJET__VORONOI_H__
2 #define __FASTJET__VORONOI_H__
3 
4 //STARTHEADER
5 // $Id: Voronoi.hh 2686 2011-11-14 09:28:22Z soyez $
6 //
7 // Copyright (c) 1994 by AT&T Bell Laboratories (see below)
8 //
9 //
10 //----------------------------------------------------------------------
11 // This file is included as part of FastJet but was mostly written by
12 // S. Fortune in C, put into C++ with memory management by S
13 // O'Sullivan, and with further interface and memeory management
14 // modifications by Gregory Soyez.
15 //
16 // Permission to use, copy, modify, and distribute this software for
17 // any purpose without fee is hereby granted, provided that this
18 // entire notice is included in all copies of any software which is or
19 // includes a copy or modification of this software and in all copies
20 // of the supporting documentation for such software. THIS SOFTWARE IS
21 // BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED WARRANTY.
22 // IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY REPRESENTATION
23 // OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY OF THIS
24 // SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
25 //
26 //----------------------------------------------------------------------
27 //ENDHEADER
28 
29 
30 /*
31 * The author of this software is Steven Fortune.
32 * Copyright (c) 1994 by AT&T Bell Laboratories.
33 * Permission to use, copy, modify, and distribute this software for any
34 * purpose without fee is hereby granted, provided that this entire notice
35 * is included in all copies of any software which is or includes a copy
36 * or modification of this software and in all copies of the supporting
37 * documentation for such software.
38 * THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
39 * WARRANTY. IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY
40 * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY
41 * OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
42 */
43 
44 /*
45 * This code was originally written by Stephan Fortune in C code. I,
46 * Shane O'Sullivan, have since modified it, encapsulating it in a C++
47 * class and, fixing memory leaks and adding accessors to the Voronoi
48 * Edges. Permission to use, copy, modify, and distribute this
49 * software for any purpose without fee is hereby granted, provided
50 * that this entire notice is included in all copies of any software
51 * which is or includes a copy or modification of this software and in
52 * all copies of the supporting documentation for such software. THIS
53 * SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
54 * WARRANTY. IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY
55 * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE
56 * MERCHANTABILITY OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR
57 * PURPOSE.
58 */
59 
60 /*
61  * This code, included in the FastJet distribution, was originally
62  * written by Stephan Fortune in C and adapted to C++ by Shane
63  * O'Sullivan under the terms repported above.
64  *
65  * Below are the list of changes implemented by the FastJet authors:
66  *
67  * 2011-11-14 Gregory Soyez <soyez@fastjet.fr>
68  *
69  * * removed 'plot' and 'triangulate' (were always 0)
70  * * removed unused plot functions (openpl, circle, range,
71  * out_bisector, out_ep, out_vertex, out_site, out_triple)
72  * * removed unused 'VPoint p' in 'intersect'
73  *
74  *
75  * 2011-07-22 Gregory Soyez <soyez@fastjet.fr>
76  *
77  * * replaced Point by VPoint (to avoid any potential conflict
78  * with an already existing class Point in FastJet
79  *
80  *
81  * 2008-04-01 Gregory Soyez <soyez@fastjet.fr>
82  *
83  * * declared ystar volatile in HalfEdge (apparently fixes a bug
84  * related to VD computations with points on a grid)
85  *
86  *
87  * 2007-05-07 Gregory Soyez <soyez@fastjet.fr>
88  *
89  * * put the code in the fastjet namespace
90  *
91  * * replaced float by double
92  *
93  * * generateVoronoi() takes a vector of Point instead of 2
94  * pointers
95  *
96  * * added info about the parent sites to GraphEdge
97  *
98  * * removed condition on minimal distance between sites
99  *
100  */
101 
102 #include "fastjet/LimitedWarning.hh"
103 #include <vector>
104 #include <math.h>
105 #include <stdlib.h>
106 #include <string.h>
107 
108 #define DELETED -2
109 #define le 0
110 #define re 1
111 
112 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
113 
114 /**
115  * \if internal_doc
116  * @ingroup internal
117  * \class VPoint
118  * class to handle a 2d point
119  * \endif
120  */
121 class VPoint{
122 public:
123  /// defailt ctor
124  VPoint() : x(0.0), y(0.0) {}
125 
126  /// ctor with initialisation
127  VPoint(double _x, double _y) : x(_x), y(_y) {}
128 
129  /// addition
130  inline VPoint operator + (const VPoint &p) const{
131  return VPoint(x+p.x, y+p.y);
132  }
133 
134  /// subtraction
135  inline VPoint operator - (const VPoint &p) const{
136  return VPoint(x-p.x, y-p.y);
137  }
138 
139  /// scalar multiplication
140  inline VPoint operator * (const double t) const{
141  return VPoint(x*t, y*t);
142  }
143 
144  /// vector coordinates
145  double x,y;
146 };
147 
148 
149 /// norm of a vector
150 inline double norm(const VPoint p){
151  return p.x*p.x+p.y*p.y;
152 }
153 
154 
155 /// 2D vector product
156 inline double vector_product(const VPoint &p1, const VPoint &p2){
157  return p1.x*p2.y-p1.y*p2.x;
158 }
159 
160 
161 /// scalar product
162 inline double scalar_product(const VPoint &p1, const VPoint &p2){
163  return p1.x*p2.x+p1.y*p2.y;
164 }
165 
166 
167 /**
168  * \if internal_doc
169  * @ingroup internal
170  * \class GraphEdge
171  * handle an edge of the Voronoi Diagram.
172  * \endif
173  */
174 class GraphEdge{
175 public:
176  /// coordinates of the extreme points
177  double x1,y1,x2,y2;
178 
179  /// indices of the parent sites that define the edge
180  int point1, point2;
181 
182  /// pointer to the next edge
183  GraphEdge* next;
184 };
185 
186 
187 /**
188  * \if internal_doc
189  * @ingroup internal
190  * \class Site
191  * structure used both for particle sites and for vertices.
192  * \endif
193  */
194 class Site{
195  public:
196  VPoint coord;
197  int sitenbr;
198  int refcnt;
199 };
200 
201 
202 
203 class Freenode{
204 public:
205  Freenode *nextfree;
206 };
207 
208 
209 class FreeNodeArrayList{
210 public:
211  Freenode* memory;
212  FreeNodeArrayList* next;
213 };
214 
215 
216 class Freelist{
217 public:
218  Freenode *head;
219  int nodesize;
220 };
221 
222 class Edge{
223 public:
224  double a,b,c;
225  Site *ep[2];
226  Site *reg[2];
227  int edgenbr;
228 };
229 
230 
231 class Halfedge{
232 public:
233  Halfedge *ELleft, *ELright;
234  Edge *ELedge;
235  int ELrefcnt;
236  char ELpm;
237  Site *vertex;
238  volatile double ystar;
239  Halfedge *PQnext;
240 };
241 
242 /**
243  * \if internal_doc
244  * @ingroup internal
245  * \class VoronoiDiagramGenerator
246  * Shane O'Sullivan C++ version of Stephan Fortune Voronoi diagram
247  * generator
248  * \endif
249  */
250 class VoronoiDiagramGenerator{
251 public:
252  VoronoiDiagramGenerator();
253  ~VoronoiDiagramGenerator();
254 
255  bool generateVoronoi(std::vector<VPoint> *_parent_sites,
256  double minX, double maxX, double minY, double maxY,
257  double minDist=0);
258 
259  inline void resetIterator(){
260  iteratorEdges = allEdges;
261  }
262 
263  bool getNext(GraphEdge **e){
264  if(iteratorEdges == 0)
265  return false;
266 
267  *e = iteratorEdges;
268  iteratorEdges = iteratorEdges->next;
269  return true;
270  }
271 
272  std::vector<VPoint> *parent_sites;
273  int n_parent_sites;
274 
275 private:
276  void cleanup();
277  void cleanupEdges();
278  char *getfree(Freelist *fl);
279  Halfedge *PQfind();
280  int PQempty();
281 
282  Halfedge **ELhash;
283  Halfedge *HEcreate(), *ELleft(), *ELright(), *ELleftbnd();
284  Halfedge *HEcreate(Edge *e,int pm);
285 
286  VPoint PQ_min();
287  Halfedge *PQextractmin();
288  void freeinit(Freelist *fl,int size);
289  void makefree(Freenode *curr,Freelist *fl);
290  void geominit();
291  void plotinit();
292 
293  // GS: removed the unused (always ==0) argument
294  bool voronoi(/*int triangulate*/);
295  void ref(Site *v);
296  void deref(Site *v);
297  void endpoint(Edge *e,int lr,Site * s);
298 
299  void ELdelete(Halfedge *he);
300  Halfedge *ELleftbnd(VPoint *p);
301  Halfedge *ELright(Halfedge *he);
302  void makevertex(Site *v);
303 
304  void PQinsert(Halfedge *he,Site * v, double offset);
305  void PQdelete(Halfedge *he);
306  bool ELinitialize();
307  void ELinsert(Halfedge *lb, Halfedge *newHe);
308  Halfedge * ELgethash(int b);
309  Halfedge *ELleft(Halfedge *he);
310  Site *leftreg(Halfedge *he);
311  bool PQinitialize();
312  int PQbucket(Halfedge *he);
313  void clip_line(Edge *e);
314  char *myalloc(unsigned n);
315  int right_of(Halfedge *el,VPoint *p);
316 
317  Site *rightreg(Halfedge *he);
318  Edge *bisect(Site *s1, Site *s2);
319  double dist(Site *s,Site *t);
320 
321  // GS: 'p' is unused and always ==0 (see also comment by
322  // S. O'Sullivan in the source file), so we remove it
323  Site *intersect(Halfedge *el1, Halfedge *el2 /*, VPoint *p=0*/);
324 
325  Site *nextone();
326 
327  void pushGraphEdge(double x1, double y1, double x2, double y2,
328  Site *s1, Site *s2);
329 
330  // Gregory Soyez: unused plotting methods
331  // void openpl();
332  // void circle(double x, double y, double radius);
333  // void range(double minX, double minY, double maxX, double maxY);
334  //
335  // void out_bisector(Edge *e);
336  // void out_ep(Edge *e);
337  // void out_vertex(Site *v);
338  // void out_site(Site *s);
339  //
340  // void out_triple(Site *s1, Site *s2,Site * s3);
341 
342  Freelist hfl;
343  Halfedge *ELleftend, *ELrightend;
344  int ELhashsize;
345 
346  int sorted, debug;
347  double xmin, xmax, ymin, ymax, deltax, deltay;
348 
349  Site *sites;
350  int nsites;
351  int siteidx;
352  int sqrt_nsites;
353  int nvertices;
354  Freelist sfl;
355  Site *bottomsite;
356 
357  int nedges;
358  Freelist efl;
359  int PQhashsize;
360  Halfedge *PQhash;
361  int PQcount;
362  int PQmin;
363 
364  int ntry, totalsearch;
365  double pxmin, pxmax, pymin, pymax, cradius;
366  int total_alloc;
367 
368  double borderMinX, borderMaxX, borderMinY, borderMaxY;
369 
370  FreeNodeArrayList* allMemoryList;
371  FreeNodeArrayList* currentMemoryBlock;
372 
373  GraphEdge* allEdges;
374  GraphEdge* iteratorEdges;
375 
376  double minDistanceBetweenSites;
377 
378  static LimitedWarning _warning_degeneracy;
379 };
380 
381 int scomp(const void *p1,const void *p2);
382 
383 
384 FASTJET_END_NAMESPACE
385 
386 #endif
double norm(const VPoint p)
norm of a vector
Definition: Voronoi.hh:150
Selector operator*(const Selector &s1, const Selector &s2)
successive application of 2 selectors
Definition: Selector.cc:507
double scalar_product(const VPoint &p1, const VPoint &p2)
scalar product
Definition: Voronoi.hh:162
double vector_product(const VPoint &p1, const VPoint &p2)
2D vector product
Definition: Voronoi.hh:156