CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
JetFinderMethods.cc
Go to the documentation of this file.
1 // JetFinderMethods.cc
2 // Author: Alex Barbieri
3 //
4 // This file should contain the different algorithms used to find jets.
5 // Currently the standard is the sliding window method, used by both
6 // HI and PP.
7 
11 
13 
14 #include <vector>
15 
16 namespace l1t {
17 
18  int deltaGctPhi(const CaloRegion & region, const CaloRegion & neighbor)
19  {
20  int phi1 = region.hwPhi();
21  int phi2 = neighbor.hwPhi();
22  int diff = phi1 - phi2;
23  if (std::abs(phi1 - phi2) == L1CaloRegionDetId::N_PHI-1) { //18 regions in phi
24  diff = -diff/std::abs(diff);
25  }
26  return diff;
27  }
28 
29  // turn each central region into a jet
30  void passThroughJets(const std::vector<l1t::CaloRegion> * regions,
31  std::vector<l1t::Jet> * uncalibjets)
32  {
33  for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
34  int jetQual = 0;
35  if( region->hwEta() < 4 || region->hwEta() > 17)
36  jetQual = 2;
37  int jetET = region->hwPt();
38  int jetEta = region->hwEta();
39  int jetPhi = region->hwPhi();
40 
41  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0,0,0,0);
42  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
43  uncalibjets->push_back(theJet);
44  }
45  }
46 
47  void slidingWindowJetFinder(const int jetSeedThreshold,
48  const std::vector<l1t::CaloRegion> * regions,
49  std::vector<l1t::Jet> * uncalibjets)
50  {
51  // std::cout << "Jet Seed: " << jetSeedThreshold << std::endl;
52  for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
53  int regionET = region->hwPt(); //regionPhysicalEt(*region);
54  if (regionET <= jetSeedThreshold) continue;
55  int neighborN_et = 0;
56  int neighborS_et = 0;
57  int neighborE_et = 0;
58  int neighborW_et = 0;
59  int neighborNE_et = 0;
60  int neighborSW_et = 0;
61  int neighborNW_et = 0;
62  int neighborSE_et = 0;
63  unsigned int nNeighbors = 0;
64  for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
65  int neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
66  if(deltaGctPhi(*region, *neighbor) == 1 &&
67  (region->hwEta() ) == neighbor->hwEta()) {
68  neighborN_et = neighborET;
69  nNeighbors++;
70  continue;
71  }
72  else if(deltaGctPhi(*region, *neighbor) == -1 &&
73  (region->hwEta() ) == neighbor->hwEta()) {
74  neighborS_et = neighborET;
75  nNeighbors++;
76  continue;
77  }
78  else if(deltaGctPhi(*region, *neighbor) == 0 &&
79  (region->hwEta() + 1) == neighbor->hwEta()) {
80  neighborE_et = neighborET;
81  nNeighbors++;
82  continue;
83  }
84  else if(deltaGctPhi(*region, *neighbor) == 0 &&
85  (region->hwEta() - 1) == neighbor->hwEta()) {
86  neighborW_et = neighborET;
87  nNeighbors++;
88  continue;
89  }
90  else if(deltaGctPhi(*region, *neighbor) == 1 &&
91  (region->hwEta() + 1) == neighbor->hwEta()) {
92  neighborNE_et = neighborET;
93  nNeighbors++;
94  continue;
95  }
96  else if(deltaGctPhi(*region, *neighbor) == -1 &&
97  (region->hwEta() - 1) == neighbor->hwEta()) {
98  neighborSW_et = neighborET;
99  nNeighbors++;
100  continue;
101  }
102  else if(deltaGctPhi(*region, *neighbor) == 1 &&
103  (region->hwEta() - 1) == neighbor->hwEta()) {
104  neighborNW_et = neighborET;
105  nNeighbors++;
106  continue;
107  }
108  else if(deltaGctPhi(*region, *neighbor) == -1 &&
109  (region->hwEta() + 1) == neighbor->hwEta()) {
110  neighborSE_et = neighborET;
111  nNeighbors++;
112  continue;
113  }
114  }
115  if(regionET > neighborN_et &&
116  regionET > neighborNW_et &&
117  regionET > neighborW_et &&
118  regionET > neighborSW_et &&
119  regionET >= neighborNE_et &&
120  regionET >= neighborE_et &&
121  regionET >= neighborSE_et &&
122  regionET >= neighborS_et) {
123  unsigned int jetET = regionET +
124  neighborN_et + neighborS_et + neighborE_et + neighborW_et +
125  neighborNE_et + neighborSW_et + neighborSE_et + neighborNW_et;
126 
127  int jetPhi = region->hwPhi();
128  int jetEta = region->hwEta();
129 
130  bool neighborCheck = (nNeighbors == 8);
131  // On the eta edge we only expect 5 neighbors
132  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
133  neighborCheck = true;
134 
135  if (!neighborCheck) {
136  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
137  assert(false);
138  }
139 
140  //first iteration, eta cut defines forward
141  //const bool forward = (jetEta <= 4 || jetEta >= 17);
142  const bool forward = (jetEta < 4 || jetEta > 17);
143  int jetQual = 0;
144  if(forward)
145  jetQual |= 0x2;
146 
147  // check for input overflow regions
148  if(forward && regionET == 255) {
149  jetET = 1023; // 10 bit max
150  } else if(!forward && regionET == 1023) {
151  jetET = 1023; // 10 bit max
152  } else if(region->hwEta() == 17) {
153  if(neighborNE_et == 255 || neighborE_et == 255 || neighborSE_et == 255)
154  jetET = 1023; // 10 bit max
155  } else if(region->hwEta() == 4) {
156  if(neighborNW_et == 255 || neighborW_et == 255 || neighborSW_et == 255)
157  jetET = 1023; // 10 bit max
158  }
159 
160  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0,0,0,0);
161  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
162  //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
163 
164  uncalibjets->push_back(theJet);
165  }
166  }
167  }
168 
169  void TwelveByTwelveFinder(const int jetSeedThreshold,
170  const std::vector<l1t::CaloRegion> * regions,
171  std::vector<l1t::Jet> * uncalibjets)
172  {
173  // std::cout << "Jet Seed: " << jetSeedThreshold << std::endl;
174  for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
175  int regionET = region->hwPt(); //regionPhysicalEt(*region);
176  if (regionET < jetSeedThreshold) continue;
177  int neighborN_et = 0;
178  int neighborS_et = 0;
179  int neighborE_et = 0;
180  int neighborW_et = 0;
181  int neighborNE_et = 0;
182  int neighborSW_et = 0;
183  int neighborNW_et = 0;
184  int neighborSE_et = 0;
185  unsigned int nNeighbors = 0;
186  for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
187  int neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
188  if(deltaGctPhi(*region, *neighbor) == 1 &&
189  (region->hwEta() ) == neighbor->hwEta()) {
190  neighborN_et = neighborET;
191  nNeighbors++;
192  continue;
193  }
194  else if(deltaGctPhi(*region, *neighbor) == -1 &&
195  (region->hwEta() ) == neighbor->hwEta()) {
196  neighborS_et = neighborET;
197  nNeighbors++;
198  continue;
199  }
200  else if(deltaGctPhi(*region, *neighbor) == 0 &&
201  (region->hwEta() + 1) == neighbor->hwEta()) {
202  neighborE_et = neighborET;
203  nNeighbors++;
204  continue;
205  }
206  else if(deltaGctPhi(*region, *neighbor) == 0 &&
207  (region->hwEta() - 1) == neighbor->hwEta()) {
208  neighborW_et = neighborET;
209  nNeighbors++;
210  continue;
211  }
212  else if(deltaGctPhi(*region, *neighbor) == 1 &&
213  (region->hwEta() + 1) == neighbor->hwEta()) {
214  neighborNE_et = neighborET;
215  nNeighbors++;
216  continue;
217  }
218  else if(deltaGctPhi(*region, *neighbor) == -1 &&
219  (region->hwEta() - 1) == neighbor->hwEta()) {
220  neighborSW_et = neighborET;
221  nNeighbors++;
222  continue;
223  }
224  else if(deltaGctPhi(*region, *neighbor) == 1 &&
225  (region->hwEta() - 1) == neighbor->hwEta()) {
226  neighborNW_et = neighborET;
227  nNeighbors++;
228  continue;
229  }
230  else if(deltaGctPhi(*region, *neighbor) == -1 &&
231  (region->hwEta() + 1) == neighbor->hwEta()) {
232  neighborSE_et = neighborET;
233  nNeighbors++;
234  continue;
235  }
236  }
237  unsigned int jetET = regionET +
238  neighborN_et + neighborS_et + neighborE_et + neighborW_et +
239  neighborNE_et + neighborSW_et + neighborSE_et + neighborNW_et;
240 
241  int jetPhi = region->hwPhi();
242  int jetEta = region->hwEta();
243 
244  bool neighborCheck = (nNeighbors == 8);
245  // On the eta edge we only expect 5 neighbors
246  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
247  neighborCheck = true;
248 
249  if (!neighborCheck) {
250  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
251  assert(false);
252  }
253 
254  //first iteration, eta cut defines forward
255  //const bool forward = (jetEta <= 4 || jetEta >= 17);
256  const bool forward = (jetEta < 4 || jetEta > 17);
257  int jetQual = 0;
258  if(forward)
259  jetQual |= 0x2;
260 
261  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0,0,0,0);
262  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
263  //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
264 
265  uncalibjets->push_back(theJet);
266  }
267  }
268 
269  void TwoByTwoFinder(const int jetSeedThreshold,
270  const int etaMask,
271  const std::vector<l1t::CaloRegion> * regions,
272  std::vector<l1t::Jet> * uncalibjets)
273  {
274  for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
275  int regionET = region->hwPt();
276  if (regionET <= jetSeedThreshold) continue;
277  int subEta = region->hwEta();
278  if((etaMask & (1<<subEta))>>subEta) regionET = 0;
279  int neighborN_et = 0;
280  int neighborS_et = 0;
281  int neighborE_et = 0;
282  int neighborW_et = 0;
283  int neighborNE_et = 0;
284  int neighborSW_et = 0;
285  int neighborNW_et = 0;
286  int neighborSE_et = 0;
287  unsigned int nNeighbors = 0;
288  for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
289  int neighborET = neighbor->hwPt();
290  int subEta2 = neighbor->hwEta();
291  if((etaMask & (1<<subEta2))>>subEta2) neighborET = 0;
292 
293  if(deltaGctPhi(*region, *neighbor) == 1 &&
294  (region->hwEta() ) == neighbor->hwEta()) {
295  neighborN_et = neighborET;
296  nNeighbors++;
297  continue;
298  }
299  else if(deltaGctPhi(*region, *neighbor) == -1 &&
300  (region->hwEta() ) == neighbor->hwEta()) {
301  neighborS_et = neighborET;
302  nNeighbors++;
303  continue;
304  }
305  else if(deltaGctPhi(*region, *neighbor) == 0 &&
306  (region->hwEta() + 1) == neighbor->hwEta()) {
307  neighborE_et = neighborET;
308  nNeighbors++;
309  continue;
310  }
311  else if(deltaGctPhi(*region, *neighbor) == 0 &&
312  (region->hwEta() - 1) == neighbor->hwEta()) {
313  neighborW_et = neighborET;
314  nNeighbors++;
315  continue;
316  }
317  else if(deltaGctPhi(*region, *neighbor) == 1 &&
318  (region->hwEta() + 1) == neighbor->hwEta()) {
319  neighborNE_et = neighborET;
320  nNeighbors++;
321  continue;
322  }
323  else if(deltaGctPhi(*region, *neighbor) == -1 &&
324  (region->hwEta() - 1) == neighbor->hwEta()) {
325  neighborSW_et = neighborET;
326  nNeighbors++;
327  continue;
328  }
329  else if(deltaGctPhi(*region, *neighbor) == 1 &&
330  (region->hwEta() - 1) == neighbor->hwEta()) {
331  neighborNW_et = neighborET;
332  nNeighbors++;
333  continue;
334  }
335  else if(deltaGctPhi(*region, *neighbor) == -1 &&
336  (region->hwEta() + 1) == neighbor->hwEta()) {
337  neighborSE_et = neighborET;
338  nNeighbors++;
339  continue;
340  }
341 
342  }
343  if(regionET > neighborN_et &&
344  regionET > neighborNW_et &&
345  regionET > neighborW_et &&
346  regionET > neighborSW_et &&
347  regionET >= neighborNE_et &&
348  regionET >= neighborE_et &&
349  regionET >= neighborSE_et &&
350  regionET >= neighborS_et) {
351 
352  // use the highest-pT 2x2 jet inside this 3x3
353  unsigned int jetET_NW;
354  unsigned int jetET_NE;
355  unsigned int jetET_SW;
356  unsigned int jetET_SE;
357 
358  jetET_NW = regionET + neighborW_et + neighborNW_et + neighborN_et;
359  jetET_NE = regionET + neighborE_et + neighborNE_et + neighborN_et;
360  jetET_SW = regionET + neighborS_et + neighborSW_et + neighborW_et;
361  jetET_SE = regionET + neighborS_et + neighborSE_et + neighborE_et;
362 
363  unsigned int jetET = std::max(jetET_NW, jetET_NE);
364  jetET = std::max(jetET, jetET_SW);
365  jetET = std::max(jetET, jetET_SE);
366 
367  int jetPhi = region->hwPhi();
368  int jetEta = region->hwEta();
369 
370  bool neighborCheck = (nNeighbors == 8);
371  // On the eta edge we only expect 5 neighbor
372  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
373  neighborCheck = true;
374 
375  if (!neighborCheck) {
376  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
377  assert(false);
378  }
379 
380  //first iteration, eta cut defines forward
381  const bool forward = (jetEta < 4 || jetEta > 17);
382  int jetQual = 0;
383  if(forward)
384  jetQual |= 0x2;
385 
386  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0,0,0,0);
387  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
388  uncalibjets->push_back(theJet);
389  }
390  }
391  }
392 }
assert(m_qm.get())
int hwPhi() const
Definition: L1Candidate.cc:79
void TwelveByTwelveFinder(const int, const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
Definition: Jet.h:13
void TwoByTwoFinder(const int, const int, const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void passThroughJets(const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
int deltaGctPhi(const CaloRegion &region, const CaloRegion &neighbor)
void slidingWindowJetFinder(const int, const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
tuple cout
Definition: gather_cfg.py:121
static const unsigned N_PHI