Fork me on GitHub

source: git/external/fastjet/internal/LazyTiling25.hh@ 5e4c40d

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 5e4c40d was 35cdc46, checked in by Pavel Demin <demin@…>, 10 years ago

upgrade FastJet to version 3.1.0-beta.1, upgrade Nsubjettiness to version 2.1.0, add SoftKiller version 1.0.0

  • Property mode set to 100644
File size: 6.3 KB
Line 
1#ifndef __FASTJET_LAZYTILING25_HH__
2#define __FASTJET_LAZYTILING25_HH__
3
4// #define INSTRUMENT2 1
5
6//FJSTARTHEADER
7// $Id: LazyTiling25.hh 3477 2014-07-29 14:34:39Z salam $
8//
9// Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
10//
11//----------------------------------------------------------------------
12// This file is part of FastJet.
13//
14// FastJet is free software; you can redistribute it and/or modify
15// it under the terms of the GNU General Public License as published by
16// the Free Software Foundation; either version 2 of the License, or
17// (at your option) any later version.
18//
19// The algorithms that underlie FastJet have required considerable
20// development. They are described in the original FastJet paper,
21// hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
22// FastJet as part of work towards a scientific publication, please
23// quote the version you use and include a citation to the manual and
24// optionally also to hep-ph/0512210.
25//
26// FastJet is distributed in the hope that it will be useful,
27// but WITHOUT ANY WARRANTY; without even the implied warranty of
28// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
29// GNU General Public License for more details.
30//
31// You should have received a copy of the GNU General Public License
32// along with FastJet. If not, see <http://www.gnu.org/licenses/>.
33//----------------------------------------------------------------------
34//FJENDHEADER
35
36//#include "fastjet/PseudoJet.hh"
37#include "fastjet/internal/MinHeap.hh"
38#include "fastjet/ClusterSequence.hh"
39#include "fastjet/internal/LazyTiling9Alt.hh"
40#include "fastjet/internal/LazyTiling9.hh"
41
42
43
44FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
45
46typedef Tile2Base<25> Tile25;
47
48
49// class Tile25 {
50// public:
51// /// pointers to neighbouring tiles, including self
52// Tile25 * begin_tiles[25];
53// /// neighbouring tiles, excluding self
54// Tile25 ** surrounding_tiles;
55// /// half of neighbouring tiles, no self
56// Tile25 ** RH_tiles;
57// /// just beyond end of tiles
58// Tile25 ** end_tiles;
59// /// start of list of BriefJets contained in this tile
60// TiledJet * head;
61// /// sometimes useful to be able to tag a tile
62// bool tagged;
63// /// for all particles in the tile, this stores the largest of the
64// /// (squared) nearest-neighbour distances.
65// double max_NN_dist;
66// double eta_centre, phi_centre;
67//
68// bool is_near_zero_phi(double tile_size_phi) const {
69// return phi_centre < 2*tile_size_phi || (twopi-phi_centre) < 2*tile_size_phi;
70// }
71// };
72
73
74//----------------------------------------------------------------------
75class LazyTiling25 {
76public:
77 LazyTiling25(ClusterSequence & cs);
78
79 void run();
80
81protected:
82 ClusterSequence & _cs;
83 const std::vector<PseudoJet> & _jets;
84 std::vector<Tile25> _tiles;
85
86#ifdef INSTRUMENT2
87 int _ncall; // GPS tmp
88 int _ncall_dtt; // GPS tmp
89#endif // INSTRUMENT2
90
91 double _Rparam, _R2, _invR2;
92 double _tiles_eta_min, _tiles_eta_max;
93 double _tile_size_eta, _tile_size_phi;
94 double _tile_half_size_eta, _tile_half_size_phi;
95 int _n_tiles_phi,_tiles_ieta_min,_tiles_ieta_max;
96
97 std::vector<TiledJet *> _jets_for_minheap;
98
99 //MinHeap _minheap;
100
101 void _initialise_tiles();
102
103 // reasonably robust return of tile index given ieta and iphi, in particular
104 // it works even if iphi is negative
105 inline int _tile_index (int ieta, int iphi) const {
106 // note that (-1)%n = -1 so that we have to add _n_tiles_phi
107 // before performing modulo operation
108 return (ieta-_tiles_ieta_min)*_n_tiles_phi
109 + (iphi+_n_tiles_phi) % _n_tiles_phi;
110 }
111
112 void _bj_remove_from_tiles(TiledJet * const jet);
113
114 /// returns the tile index given the eta and phi values of a jet
115 int _tile_index(const double eta, const double phi) const;
116
117 // sets up information regarding the tiling of the given jet
118 void _tj_set_jetinfo(TiledJet * const jet, const int _jets_index);
119
120 void _print_tiles(TiledJet * briefjets ) const;
121 void _add_neighbours_to_tile_union(const int tile_index,
122 std::vector<int> & tile_union, int & n_near_tiles) const;
123 void _add_untagged_neighbours_to_tile_union(const int tile_index,
124 std::vector<int> & tile_union, int & n_near_tiles);
125 void _add_untagged_neighbours_to_tile_union_using_max_info(const TiledJet * const jet,
126 std::vector<int> & tile_union, int & n_near_tiles);
127 double _distance_to_tile(const TiledJet * bj, const Tile25 *)
128#ifdef INSTRUMENT2
129 ;
130#else
131 const;
132#endif
133 void _update_jetX_jetI_NN(TiledJet * jetX, TiledJet * jetI, std::vector<TiledJet *> & jets_for_minheap);
134
135 void _set_NN(TiledJet * jetI, std::vector<TiledJet *> & jets_for_minheap);
136
137 // return the diJ (multiplied by _R2) for this jet assuming its NN
138 // info is correct
139 template <class J> double _bj_diJ(const J * const jet) const {
140 double kt2 = jet->kt2;
141 if (jet->NN != NULL) {if (jet->NN->kt2 < kt2) {kt2 = jet->NN->kt2;}}
142 return jet->NN_dist * kt2;
143 }
144
145
146 //----------------------------------------------------------------------
147 template <class J> inline void _bj_set_jetinfo(
148 J * const jetA, const int _jets_index) const {
149 jetA->eta = _jets[_jets_index].rap();
150 jetA->phi = _jets[_jets_index].phi_02pi();
151 jetA->kt2 = _cs.jet_scale_for_algorithm(_jets[_jets_index]);
152 jetA->_jets_index = _jets_index;
153 // initialise NN info as well
154 jetA->NN_dist = _R2;
155 jetA->NN = NULL;
156 }
157
158
159 //----------------------------------------------------------------------
160 template <class J> inline double _bj_dist(
161 const J * const jetA, const J * const jetB)
162#ifdef INSTRUMENT2
163 {
164 _ncall++; // GPS tmp
165#else
166 const {
167#endif
168 double dphi = std::abs(jetA->phi - jetB->phi);
169 double deta = (jetA->eta - jetB->eta);
170 if (dphi > pi) {dphi = twopi - dphi;}
171 return dphi*dphi + deta*deta;
172 }
173
174
175 //----------------------------------------------------------------------
176 template <class J> inline double _bj_dist_not_periodic(
177 const J * const jetA, const J * const jetB)
178#ifdef INSTRUMENT2
179 {
180 _ncall++; // GPS tmp
181#else
182 const {
183#endif
184 //_ncall++; // GPS tmp
185 double dphi = jetA->phi - jetB->phi;
186 double deta = (jetA->eta - jetB->eta);
187 return dphi*dphi + deta*deta;
188 }
189
190};
191
192
193FASTJET_END_NAMESPACE
194
195#endif // __FASTJET_LAZYTILING25_HH__
Note: See TracBrowser for help on using the repository browser.