Fork me on GitHub

source: git/external/fastjet/internal/Dnn3piCylinder.hh@ b4a3c55

Timing
Last change on this file since b4a3c55 was b7b836a, checked in by Pavel Demin <pavel-demin@…>, 6 years ago

update FastJet library to 3.3.1 and FastJet Contrib library to 1.036

  • Property mode set to 100644
File size: 10.0 KB
Line 
1//FJSTARTHEADER
2// $Id: Dnn3piCylinder.hh 4354 2018-04-22 07:12:37Z salam $
3//
4// Copyright (c) 2005-2018, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5//
6//----------------------------------------------------------------------
7// This file is part of FastJet.
8//
9// FastJet is free software; you can redistribute it and/or modify
10// it under the terms of the GNU General Public License as published by
11// the Free Software Foundation; either version 2 of the License, or
12// (at your option) any later version.
13//
14// The algorithms that underlie FastJet have required considerable
15// development. They are described in the original FastJet paper,
16// hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17// FastJet as part of work towards a scientific publication, please
18// quote the version you use and include a citation to the manual and
19// optionally also to hep-ph/0512210.
20//
21// FastJet is distributed in the hope that it will be useful,
22// but WITHOUT ANY WARRANTY; without even the implied warranty of
23// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24// GNU General Public License for more details.
25//
26// You should have received a copy of the GNU General Public License
27// along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28//----------------------------------------------------------------------
29//FJENDHEADER
30
31
32#ifndef DROP_CGAL // in case we do not have the code for CGAL
33#ifndef __FASTJET_DNN3PICYLINDER_HH__
34#define __FASTJET_DNN3PICYLINDER_HH__
35
36#include "fastjet/internal/DynamicNearestNeighbours.hh"
37#include "fastjet/internal/DnnPlane.hh"
38#include "fastjet/internal/numconsts.hh"
39
40FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
41
42/// \if internal_doc
43/// @ingroup internal
44/// \class Dnn3piCylinder
45/// class derived from DynamicNearestNeighbours that provides an
46/// implementation for the surface of cylinder (using one
47/// DnnPlane object spanning 0--3pi).
48/// \endif
49class Dnn3piCylinder : public DynamicNearestNeighbours {
50 public:
51 /// empty initaliser
52 Dnn3piCylinder() {}
53
54 /// Initialiser from a set of points on an Eta-Phi plane, where
55 /// eta can have an arbitrary ranges and phi must be in range
56 /// 0 <= phi < 2pi;
57 ///
58 /// NB: this class is more efficient than the plain Dnn4piCylinder
59 /// class, but can give wrong answers when the nearest neighbour is
60 /// further away than 2pi (in this case a point's nearest neighbour
61 /// becomes itself, because it is considered to be a distance 2pi
62 /// away). For the kt-algorithm (e.g.) this is actually not a
63 /// problem (the distance need only be accurate when it is less than
64 /// R), so we can tell the routine to ignore this problem --
65 /// alternatively the routine will crash if it detects it occurring
66 /// (only when finding the nearest neighbour index, not its
67 /// distance).
68 Dnn3piCylinder(const std::vector<EtaPhi> &,
69 const bool & ignore_nearest_is_mirror = false,
70 const bool & verbose = false );
71
72 /// Returns the index of the nearest neighbour of point labelled
73 /// by ii (assumes ii is valid)
74 int NearestNeighbourIndex(const int ii) const ;
75
76 /// Returns the distance to the nearest neighbour of point labelled
77 /// by index ii (assumes ii is valid)
78 double NearestNeighbourDistance(const int ii) const ;
79
80 /// Returns true iff the given index corresponds to a point that
81 /// exists in the DNN structure (meaning that it has been added, and
82 /// not removed in the meantime)
83 bool Valid(const int index) const;
84
85 void RemoveAndAddPoints(const std::vector<int> & indices_to_remove,
86 const std::vector<EtaPhi> & points_to_add,
87 std::vector<int> & indices_added,
88 std::vector<int> & indices_of_updated_neighbours);
89
90 ~Dnn3piCylinder();
91
92 private:
93
94 // our extras to help us navigate, find distance, etc.
95 const static int INEXISTENT_VERTEX=-3;
96
97 bool _verbose;
98
99 bool _ignore_nearest_is_mirror;
100
101 /// Picture of how things will work... Copy 0--pi part of the 0--2pi
102 /// cylinder into a region 2pi--3pi of a Euclidean plane. Below we
103 /// show points labelled by + that have a mirror image in this
104 /// manner, while points labelled by * do not have a mirror image.
105 ///
106 /// | . |
107 /// | . |
108 /// | . |
109 /// | . |
110 /// | 2 . |
111 /// | * . |
112 /// | + . + |
113 /// | 0 . 1 |
114 /// | . |
115 /// 0 2pi 3pi
116 ///
117 /// Each "true" point has its true "cylinder" index (the index that
118 /// is known externally to this class) as well as euclidean plane
119 /// indices (main_index and mirror index in the MirrorVertexInfo
120 /// structure), which are private concepts of this class.
121 ///
122 /// In above picture our structures would hold the following info
123 /// (the picture shows the euclidean-plane numbering)
124 ///
125 /// _mirror_info[cylinder_index = 0] = (0, 1)
126 /// _mirror_info[cylinder_index = 1] = (2, INEXISTENT_VERTEX)
127 ///
128 /// We also need to be able to go from the euclidean plane indices
129 /// back to the "true" cylinder index, and for this purpose we use
130 /// the vector _cylinder_index_of_plane_vertex[...], which in the above example has
131 /// the following contents
132 ///
133 /// _cylinder_index_of_plane_vertex[0] = 0
134 /// _cylinder_index_of_plane_vertex[1] = 0
135 /// _cylinder_index_of_plane_vertex[2] = 1
136 ///
137
138 ///
139 struct MirrorVertexInfo {
140 /// index of the given point (appearing in the range 0--2pi) in the
141 /// 0--3pi euclidean plane structure (position will coincide with
142 /// that on the 0--2pi cylinder, but index labelling it will be
143 /// different)
144 int main_index;
145 /// index of the mirror point (appearing in the range 2pi--3pi) in the
146 /// 0--3pi euclidean plane structure
147 int mirror_index;
148 };
149
150 // for each "true" vertex we have reference to indices in the euclidean
151 // plane structure
152 std::vector<MirrorVertexInfo> _mirror_info;
153 // for each index in the euclidean 0--3pi plane structure we want to
154 // be able to get back to the "true" vertex index on the overall
155 // 0--2pi cylinder structure
156 std::vector<int> _cylinder_index_of_plane_vertex;
157
158 // NB: we define POINTERS here because the initialisation gave
159 // us problems (things crashed!), perhaps because in practice
160 // we were making a copy without being careful and defining
161 // a proper copy constructor.
162 DnnPlane * _DNN;
163
164 /// given a phi value in the 0--2pi range return one
165 /// in the pi--3pi range.
166 inline EtaPhi _remap_phi(const EtaPhi & point) {
167 double phi = point.second;
168 if (phi < pi) { phi += twopi ;}
169 return EtaPhi(point.first, phi);}
170
171
172 //----------------------------------------------------------------------
173 /// What on earth does this do?
174 ///
175 /// Example: last true "cylinder" index was 15
176 /// last plane index was 23
177 ///
178 /// Then: _cylinder_index_of_plane_vertex.size() = 24 and
179 /// _mirror_info.size() = 16
180 ///
181 /// IF cylinder_point's phi < pi then
182 /// create: _mirror_info[16] = (main_index = 24, mirror_index=25)
183 /// _cylinder_index_of_plane_vertex[24] = 16
184 /// _cylinder_index_of_plane_vertex[25] = 16
185 /// ELSE
186 /// create: _mirror_info[16] = (main_index = 24, mirror_index=INEXISTENT..)
187 /// _cylinder_index_of_plane_vertex[24] = 16
188 ///
189 /// ADDITIONALLY push the cylinder_point (and if it exists the mirror
190 /// copy) onto the vector plane_points.
191 void _RegisterCylinderPoint (const EtaPhi & cylinder_point,
192 std::vector<EtaPhi> & plane_points);
193};
194
195
196// here follow some inline implementations of the simpler of the
197// functions defined above
198
199//----------------------------------------------------------------------
200/// Note: one of the difficulties of the 0--3pi mapping is that
201/// a point may have its mirror copy as its own nearest neighbour
202/// (if no other point is within a distance of 2pi). This does
203/// not matter for the kt_algorithm with
204/// reasonable values of radius, but might matter for a general use
205/// of this algorithm -- depending on whether or not the user has
206/// initialised the class with instructions to ignore this problem the
207/// program will detect and ignore it, or crash.
208inline int Dnn3piCylinder::NearestNeighbourIndex(const int current) const {
209 int main_index = _mirror_info[current].main_index;
210 int mirror_index = _mirror_info[current].mirror_index;
211 int plane_index;
212 if (mirror_index == INEXISTENT_VERTEX ) {
213 plane_index = _DNN->NearestNeighbourIndex(main_index);
214 } else {
215 plane_index = (
216 _DNN->NearestNeighbourDistance(main_index) <
217 _DNN->NearestNeighbourDistance(mirror_index)) ?
218 _DNN->NearestNeighbourIndex(main_index) :
219 _DNN->NearestNeighbourIndex(mirror_index) ;
220 }
221 int this_cylinder_index = _cylinder_index_of_plane_vertex[plane_index];
222 // either the user has acknowledged the fact that they may get the
223 // mirror copy as the closest point, or crash if it should occur
224 // that mirror copy is the closest point.
225 assert(_ignore_nearest_is_mirror || this_cylinder_index != current);
226 //if (this_cylinder_index == current) {
227 // std::cerr << "WARNING point "<<current<<
228 // " has its mirror copy as its own nearest neighbour"<<endl;
229 //}
230 return this_cylinder_index;
231}
232
233inline double Dnn3piCylinder::NearestNeighbourDistance(const int current) const {
234 int main_index = _mirror_info[current].main_index;
235 int mirror_index = _mirror_info[current].mirror_index;
236 if (mirror_index == INEXISTENT_VERTEX ) {
237 return _DNN->NearestNeighbourDistance(main_index);
238 } else {
239 return (
240 _DNN->NearestNeighbourDistance(main_index) <
241 _DNN->NearestNeighbourDistance(mirror_index)) ?
242 _DNN->NearestNeighbourDistance(main_index) :
243 _DNN->NearestNeighbourDistance(mirror_index) ;
244 }
245
246}
247
248inline bool Dnn3piCylinder::Valid(const int index) const {
249 return (_DNN->Valid(_mirror_info[index].main_index));
250}
251
252
253inline Dnn3piCylinder::~Dnn3piCylinder() {
254 delete _DNN;
255}
256
257
258FASTJET_END_NAMESPACE
259
260#endif // __FASTJET_DNN3PICYLINDER_HH__
261#endif // DROP_CGAL
Note: See TracBrowser for help on using the repository browser.