FastJet 3.0.2
Voronoi.hh
00001 #ifndef __FASTJET__VORONOI_H__
00002 #define __FASTJET__VORONOI_H__
00003 
00004 //STARTHEADER
00005 // $Id: Voronoi.hh 2686 2011-11-14 09:28:22Z soyez $
00006 //
00007 // Copyright (c) 1994 by AT&T Bell Laboratories (see below)
00008 //
00009 //
00010 //----------------------------------------------------------------------
00011 // This file is included as part of FastJet but was mostly written by
00012 // S. Fortune in C, put into C++ with memory management by S
00013 // O'Sullivan, and with further interface and memeory management
00014 // modifications by Gregory Soyez.
00015 //
00016 // Permission to use, copy, modify, and distribute this software for
00017 // any purpose without fee is hereby granted, provided that this
00018 // entire notice is included in all copies of any software which is or
00019 // includes a copy or modification of this software and in all copies
00020 // of the supporting documentation for such software. THIS SOFTWARE IS
00021 // BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED WARRANTY.
00022 // IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY REPRESENTATION
00023 // OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY OF THIS
00024 // SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
00025 //
00026 //----------------------------------------------------------------------
00027 //ENDHEADER
00028 
00029 
00030 /*
00031 * The author of this software is Steven Fortune.  
00032 * Copyright (c) 1994 by AT&T Bell Laboratories.
00033 * Permission to use, copy, modify, and distribute this software for any
00034 * purpose without fee is hereby granted, provided that this entire notice
00035 * is included in all copies of any software which is or includes a copy
00036 * or modification of this software and in all copies of the supporting
00037 * documentation for such software.
00038 * THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
00039 * WARRANTY.  IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY
00040 * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY
00041 * OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
00042 */
00043 
00044 /* 
00045 * This code was originally written by Stephan Fortune in C code.  I,
00046 * Shane O'Sullivan, have since modified it, encapsulating it in a C++
00047 * class and, fixing memory leaks and adding accessors to the Voronoi
00048 * Edges.  Permission to use, copy, modify, and distribute this
00049 * software for any purpose without fee is hereby granted, provided
00050 * that this entire notice is included in all copies of any software
00051 * which is or includes a copy or modification of this software and in
00052 * all copies of the supporting documentation for such software.  THIS
00053 * SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
00054 * WARRANTY.  IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY
00055 * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE
00056 * MERCHANTABILITY OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR
00057 * PURPOSE.
00058 */
00059 
00060 /*
00061  * This code, included in the FastJet distribution, was originally
00062  * written by Stephan Fortune in C and adapted to C++ by Shane
00063  * O'Sullivan under the terms repported above.
00064  *
00065  * Below are the list of changes implemented by the FastJet authors:
00066  *
00067  * 2011-11-14  Gregory Soyez  <soyez@fastjet.fr>
00068  * 
00069  *      * removed 'plot' and 'triangulate' (were always 0)
00070  *      * removed unused plot functions (openpl, circle, range, 
00071  *        out_bisector, out_ep, out_vertex, out_site, out_triple)
00072  *      * removed unused 'VPoint p' in 'intersect'
00073  * 
00074  * 
00075  * 2011-07-22  Gregory Soyez  <soyez@fastjet.fr>
00076  * 
00077  *      * replaced Point by VPoint (to avoid any potential conflict
00078  *        with an already existing class Point in FastJet
00079  * 
00080  * 
00081  * 2008-04-01  Gregory Soyez  <soyez@fastjet.fr>
00082  * 
00083  *      * declared ystar volatile in HalfEdge (apparently fixes a bug
00084  *        related to VD computations with points on a grid)
00085  * 
00086  * 
00087  * 2007-05-07  Gregory Soyez  <soyez@fastjet.fr>
00088  * 
00089  *      * put the code in the fastjet namespace
00090  * 
00091  *      * replaced float by double
00092  * 
00093  *      * generateVoronoi() takes a vector of Point instead of 2
00094  *        pointers
00095  * 
00096  *      * added info about the parent sites to GraphEdge
00097  * 
00098  *      * removed condition on minimal distance between sites
00099  * 
00100  */
00101 
00102 #include "fastjet/LimitedWarning.hh"
00103 #include <vector>
00104 #include <math.h>
00105 #include <stdlib.h>
00106 #include <string.h>
00107 
00108 #define DELETED -2
00109 #define le 0
00110 #define re 1
00111 
00112 FASTJET_BEGIN_NAMESPACE      // defined in fastjet/internal/base.hh
00113 
00114 /**
00115  * \if internal_doc
00116  * @ingroup internal
00117  * \class VPoint
00118  * class to handle a 2d point
00119  * \endif
00120  */
00121 class VPoint{
00122 public:
00123   /// defailt ctor
00124   VPoint() : x(0.0), y(0.0) {}
00125 
00126   /// ctor with initialisation
00127   VPoint(double _x, double _y) : x(_x), y(_y) {}
00128 
00129   /// addition
00130   inline VPoint operator + (const VPoint &p) const{
00131     return VPoint(x+p.x, y+p.y);
00132   }
00133 
00134   /// subtraction
00135   inline VPoint operator - (const VPoint &p) const{
00136     return VPoint(x-p.x, y-p.y);
00137   }
00138 
00139   /// scalar multiplication
00140   inline VPoint operator * (const double t) const{
00141     return VPoint(x*t, y*t);
00142   }
00143 
00144   /// vector coordinates
00145   double x,y;
00146 };
00147 
00148 
00149 /// norm of a vector
00150 inline double norm(const VPoint p){
00151   return p.x*p.x+p.y*p.y;
00152 }
00153 
00154 
00155 /// 2D vector product
00156 inline double vector_product(const VPoint &p1, const VPoint &p2){
00157   return p1.x*p2.y-p1.y*p2.x;
00158 }
00159 
00160 
00161 /// scalar product
00162 inline double scalar_product(const VPoint &p1, const VPoint &p2){
00163   return p1.x*p2.x+p1.y*p2.y;
00164 }
00165 
00166 
00167 /**
00168  * \if internal_doc
00169  * @ingroup internal
00170  * \class GraphEdge
00171  * handle an edge of the Voronoi Diagram.
00172  * \endif
00173  */
00174 class GraphEdge{
00175 public:
00176   /// coordinates of the extreme points
00177   double x1,y1,x2,y2;
00178 
00179   /// indices of the parent sites that define the edge
00180   int point1, point2;
00181 
00182   /// pointer to the next edge
00183   GraphEdge* next;
00184 };
00185 
00186 
00187 /**
00188  * \if internal_doc
00189  * @ingroup internal
00190  * \class Site
00191  * structure used both for particle sites and for vertices.
00192  * \endif
00193  */
00194 class Site{
00195  public:
00196   VPoint        coord;
00197   int sitenbr;
00198   int refcnt;
00199 };
00200 
00201 
00202 
00203 class Freenode{
00204 public:
00205   Freenode *nextfree;
00206 };
00207 
00208 
00209 class FreeNodeArrayList{
00210 public:
00211   Freenode* memory;
00212   FreeNodeArrayList* next;
00213 };
00214 
00215 
00216 class Freelist{
00217 public:
00218   Freenode *head;
00219   int nodesize;
00220 };
00221 
00222 class Edge{
00223 public:
00224   double a,b,c;
00225   Site *ep[2];
00226   Site *reg[2];
00227   int edgenbr;
00228 };
00229 
00230 
00231 class Halfedge{
00232 public:
00233   Halfedge *ELleft, *ELright;
00234   Edge *ELedge;
00235   int ELrefcnt;
00236   char ELpm;
00237   Site *vertex;
00238   volatile double ystar;
00239   Halfedge *PQnext;
00240 };
00241 
00242 /**
00243  * \if internal_doc
00244  * @ingroup internal
00245  * \class VoronoiDiagramGenerator
00246  * Shane O'Sullivan C++ version of Stephan Fortune Voronoi diagram
00247  * generator
00248  * \endif
00249  */
00250 class VoronoiDiagramGenerator{
00251 public:
00252   VoronoiDiagramGenerator();
00253   ~VoronoiDiagramGenerator();
00254 
00255   bool generateVoronoi(std::vector<VPoint> *_parent_sites,
00256                        double minX, double maxX, double minY, double maxY, 
00257                        double minDist=0);
00258 
00259   inline void resetIterator(){
00260     iteratorEdges = allEdges;
00261   }
00262 
00263   bool getNext(GraphEdge **e){
00264     if(iteratorEdges == 0)
00265       return false;
00266     
00267     *e = iteratorEdges;
00268     iteratorEdges = iteratorEdges->next;
00269     return true;
00270   }
00271   
00272   std::vector<VPoint> *parent_sites;
00273   int n_parent_sites;
00274 
00275 private:
00276   void cleanup();
00277   void cleanupEdges();
00278   char *getfree(Freelist *fl);  
00279   Halfedge *PQfind();
00280   int PQempty();
00281         
00282   Halfedge **ELhash;
00283   Halfedge *HEcreate(), *ELleft(), *ELright(), *ELleftbnd();
00284   Halfedge *HEcreate(Edge *e,int pm);
00285   
00286   VPoint PQ_min();
00287   Halfedge *PQextractmin();     
00288   void freeinit(Freelist *fl,int size);
00289   void makefree(Freenode *curr,Freelist *fl);
00290   void geominit();
00291   void plotinit();
00292 
00293   // GS: removed the unused (always ==0) argument
00294   bool voronoi(/*int triangulate*/);
00295   void ref(Site *v);
00296   void deref(Site *v);
00297   void endpoint(Edge *e,int lr,Site * s);
00298 
00299   void ELdelete(Halfedge *he);
00300   Halfedge *ELleftbnd(VPoint *p);
00301   Halfedge *ELright(Halfedge *he);
00302   void makevertex(Site *v);
00303   
00304   void PQinsert(Halfedge *he,Site * v, double offset);
00305   void PQdelete(Halfedge *he);
00306   bool ELinitialize();
00307   void ELinsert(Halfedge *lb, Halfedge *newHe);
00308   Halfedge * ELgethash(int b);
00309   Halfedge *ELleft(Halfedge *he);
00310   Site *leftreg(Halfedge *he);
00311   bool PQinitialize();
00312   int PQbucket(Halfedge *he);
00313   void clip_line(Edge *e);
00314   char *myalloc(unsigned n);
00315   int right_of(Halfedge *el,VPoint *p);
00316 
00317   Site *rightreg(Halfedge *he);
00318   Edge *bisect(Site *s1, Site *s2);
00319   double dist(Site *s,Site *t);
00320 
00321   // GS: 'p' is unused and always ==0 (see also comment by
00322   //     S. O'Sullivan in the source file), so we remove it
00323   Site *intersect(Halfedge *el1, Halfedge *el2 /*, VPoint *p=0*/);
00324 
00325   Site *nextone();
00326 
00327   void pushGraphEdge(double x1, double y1, double x2, double y2, 
00328                      Site *s1, Site *s2);
00329 
00330   // Gregory Soyez: unused plotting methods
00331   // void openpl();
00332   // void circle(double x, double y, double radius);
00333   // void range(double minX, double minY, double maxX, double maxY);
00334   // 
00335   // void out_bisector(Edge *e);
00336   // void out_ep(Edge *e);
00337   // void out_vertex(Site *v);
00338   // void out_site(Site *s);
00339   // 
00340   // void out_triple(Site *s1, Site *s2,Site * s3);
00341 
00342   Freelist hfl;
00343   Halfedge *ELleftend, *ELrightend;
00344   int ELhashsize;
00345   
00346   int sorted, debug;
00347   double xmin, xmax, ymin, ymax, deltax, deltay;
00348   
00349   Site *sites;
00350   int nsites;
00351   int siteidx;
00352   int sqrt_nsites;
00353   int nvertices;
00354   Freelist sfl;
00355   Site *bottomsite;
00356   
00357   int nedges;
00358   Freelist efl;
00359   int PQhashsize;
00360   Halfedge *PQhash;
00361   int PQcount;
00362   int PQmin;
00363   
00364   int ntry, totalsearch;
00365   double pxmin, pxmax, pymin, pymax, cradius;
00366   int total_alloc;
00367   
00368   double borderMinX, borderMaxX, borderMinY, borderMaxY;
00369   
00370   FreeNodeArrayList* allMemoryList;
00371   FreeNodeArrayList* currentMemoryBlock;
00372   
00373   GraphEdge* allEdges;
00374   GraphEdge* iteratorEdges;
00375   
00376   double minDistanceBetweenSites;
00377 
00378   static LimitedWarning _warning_degeneracy;
00379 };
00380 
00381 int scomp(const void *p1,const void *p2);
00382 
00383 
00384 FASTJET_END_NAMESPACE
00385 
00386 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends