Fork me on GitHub

source: svn/trunk/external/fastjet/internal/Voronoi.hh@ 1204

Last change on this file since 1204 was 859, checked in by Pavel Demin, 12 years ago

update fastjet to version 3.0.3

File size: 9.4 KB
Line 
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
112FASTJET_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 */
121class VPoint{
122public:
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
150inline double norm(const VPoint p){
151 return p.x*p.x+p.y*p.y;
152}
153
154
155/// 2D vector product
156inline 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
162inline 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 */
174class GraphEdge{
175public:
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 */
194class Site{
195 public:
196 VPoint coord;
197 int sitenbr;
198 int refcnt;
199};
200
201
202
203class Freenode{
204public:
205 Freenode *nextfree;
206};
207
208
209class FreeNodeArrayList{
210public:
211 Freenode* memory;
212 FreeNodeArrayList* next;
213};
214
215
216class Freelist{
217public:
218 Freenode *head;
219 int nodesize;
220};
221
222class Edge{
223public:
224 double a,b,c;
225 Site *ep[2];
226 Site *reg[2];
227 int edgenbr;
228};
229
230
231class Halfedge{
232public:
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 */
250class VoronoiDiagramGenerator{
251public:
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
275private:
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
381int scomp(const void *p1,const void *p2);
382
383
384FASTJET_END_NAMESPACE
385
386#endif
Note: See TracBrowser for help on using the repository browser.