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 // The sorting of the jets in pT order is handled here.
8 
9 #include "DataFormats/L1TCalorimeter/interface/CaloRegion.h"
10 #include "DataFormats/L1Trigger/interface/Jet.h"
11 #include "L1Trigger/L1TCalorimeter/interface/JetFinderMethods.h"
12 
13 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegionDetId.h"
14 
15 #include <vector>
16 
17 namespace l1t {
18 
19  int deltaGctPhi(const CaloRegion & region, const CaloRegion & neighbor)
20  {
21  int phi1 = region.hwPhi();
22  int phi2 = neighbor.hwPhi();
23  int diff = phi1 - phi2;
24  if (std::abs(phi1 - phi2) == L1CaloRegionDetId::N_PHI-1) { //18 regions in phi
25  diff = -diff/std::abs(diff);
26  }
27  return diff;
28  }
29 
30  bool compareJets (l1t::Jet i, l1t::Jet j){
31  return (i.hwPt() < j.hwPt() );
32  }
33 
34  void slidingWindowJetFinder(const int jetSeedThreshold, const std::vector<l1t::CaloRegion> * regions,
35  std::vector<l1t::Jet> * uncalibjets)
36  {
37  // std::cout << "Jet Seed: " << jetSeedThreshold << std::endl;
38  for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
39  double regionET = region->hwPt(); //regionPhysicalEt(*region);
40  if (regionET <= jetSeedThreshold) continue;
41  double neighborN_et = 0;
42  double neighborS_et = 0;
43  double neighborE_et = 0;
44  double neighborW_et = 0;
45  double neighborNE_et = 0;
46  double neighborSW_et = 0;
47  double neighborNW_et = 0;
48  double neighborSE_et = 0;
49  unsigned int nNeighbors = 0;
50  for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
51  double neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
52  if(deltaGctPhi(*region, *neighbor) == 1 &&
53  (region->hwEta() ) == neighbor->hwEta()) {
54  neighborN_et = neighborET;
55  nNeighbors++;
56  continue;
57  }
58  else if(deltaGctPhi(*region, *neighbor) == -1 &&
59  (region->hwEta() ) == neighbor->hwEta()) {
60  neighborS_et = neighborET;
61  nNeighbors++;
62  continue;
63  }
64  else if(deltaGctPhi(*region, *neighbor) == 0 &&
65  (region->hwEta() + 1) == neighbor->hwEta()) {
66  neighborE_et = neighborET;
67  nNeighbors++;
68  continue;
69  }
70  else if(deltaGctPhi(*region, *neighbor) == 0 &&
71  (region->hwEta() - 1) == neighbor->hwEta()) {
72  neighborW_et = neighborET;
73  nNeighbors++;
74  continue;
75  }
76  else if(deltaGctPhi(*region, *neighbor) == 1 &&
77  (region->hwEta() + 1) == neighbor->hwEta()) {
78  neighborNE_et = neighborET;
79  nNeighbors++;
80  continue;
81  }
82  else if(deltaGctPhi(*region, *neighbor) == -1 &&
83  (region->hwEta() - 1) == neighbor->hwEta()) {
84  neighborSW_et = neighborET;
85  nNeighbors++;
86  continue;
87  }
88  else if(deltaGctPhi(*region, *neighbor) == 1 &&
89  (region->hwEta() - 1) == neighbor->hwEta()) {
90  neighborNW_et = neighborET;
91  nNeighbors++;
92  continue;
93  }
94  else if(deltaGctPhi(*region, *neighbor) == -1 &&
95  (region->hwEta() + 1) == neighbor->hwEta()) {
96  neighborSE_et = neighborET;
97  nNeighbors++;
98  continue;
99  }
100  }
101  if(regionET > neighborN_et &&
102  regionET > neighborNW_et &&
103  regionET > neighborW_et &&
104  regionET > neighborSW_et &&
105  regionET >= neighborNE_et &&
106  regionET >= neighborE_et &&
107  regionET >= neighborSE_et &&
108  regionET >= neighborS_et) {
109  unsigned int jetET = regionET +
110  neighborN_et + neighborS_et + neighborE_et + neighborW_et +
111  neighborNE_et + neighborSW_et + neighborSE_et + neighborNW_et;
112  /*
113  int jetPhi = region->hwPhi() * 4 +
114  ( - 2 * (neighborS_et + neighborSE_et + neighborSW_et)
115  + 2 * (neighborN_et + neighborNE_et + neighborNW_et) ) / jetET;
116  if(jetPhi < 0) {
117 
118  }
119  else if(jetPhi >= ((int) N_JET_PHI)) {
120  jetPhi -= N_JET_PHI;
121  }
122  int jetEta = region->hwEta() * 4 +
123  ( - 2 * (neighborW_et + neighborNW_et + neighborSW_et)
124  + 2 * (neighborE_et + neighborNE_et + neighborSE_et) ) / jetET;
125  if(jetEta < 0) jetEta = 0;
126  if(jetEta >= ((int) N_JET_ETA)) jetEta = N_JET_ETA - 1;
127  */
128  // Temporarily use the region granularity -- we will try to improve as above when code is debugged
129  int jetPhi = region->hwPhi();
130  int jetEta = region->hwEta();
131 
132  bool neighborCheck = (nNeighbors == 8);
133  // On the eta edge we only expect 5 neighbors
134  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
135  neighborCheck = true;
136 
137  if (!neighborCheck) {
138  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
139  assert(false);
140  }
141 
142  //first iteration, eta cut defines forward
143  //const bool forward = (jetEta <= 4 || jetEta >= 17);
144  const bool forward = (jetEta < 4 || jetEta > 17);
145  int jetQual = 0;
146  if(forward)
147  jetQual |= 0x2;
148 
149  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > *jetLorentz =
150  new ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> >();
151  l1t::Jet theJet(*jetLorentz, jetET, jetEta, jetPhi, jetQual);
152  //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
153 
154  uncalibjets->push_back(theJet);
155  }
156  }
157 
158  // // separate loops for the central jets at the edges of HF, composed of 6 regions only.
159  // for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
160  // double regionET = region->hwPt();
161  // // look at only the left eta wall
162  // if (region->hwEta() != 4) continue;
163  // if (regionET < jetSeedThreshold) continue;
164  // double neighborN_et = 0;
165  // double neighborS_et = 0;
166  // double neighborE_et = 0;
167  // //double neighborW_et = 0;
168  // double neighborNE_et = 0;
169  // //double neighborSW_et = 0;
170  // //double neighborNW_et = 0;
171  // double neighborSE_et = 0;
172  // unsigned int nNeighbors = 0;
173  // for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
174 
175  // double neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
176  // if(deltaGctPhi(*region, *neighbor) == 1 &&
177  // (region->hwEta() ) == neighbor->hwEta()) {
178  // neighborN_et = neighborET;
179  // nNeighbors++;
180  // continue;
181  // }
182  // else if(deltaGctPhi(*region, *neighbor) == -1 &&
183  // (region->hwEta() ) == neighbor->hwEta()) {
184  // neighborS_et = neighborET;
185  // nNeighbors++;
186  // continue;
187  // }
188  // else if(deltaGctPhi(*region, *neighbor) == 0 &&
189  // (region->hwEta() + 1) == neighbor->hwEta()) {
190  // neighborE_et = neighborET;
191  // nNeighbors++;
192  // continue;
193  // }
194  // // else if(deltaGctPhi(*region, *neighbor) == 0 &&
195  // // (region->hwEta() - 1) == neighbor->hwEta()) {
196  // // neighborW_et = neighborET;
197  // // nNeighbors++;
198  // // continue;
199  // // }
200  // else if(deltaGctPhi(*region, *neighbor) == 1 &&
201  // (region->hwEta() + 1) == neighbor->hwEta()) {
202  // neighborNE_et = neighborET;
203  // nNeighbors++;
204  // continue;
205  // }
206  // // else if(deltaGctPhi(*region, *neighbor) == -1 &&
207  // // (region->hwEta() - 1) == neighbor->hwEta()) {
208  // // neighborSW_et = neighborET;
209  // // nNeighbors++;
210  // // continue;
211  // // }
212  // // else if(deltaGctPhi(*region, *neighbor) == 1 &&
213  // // (region->hwEta() - 1) == neighbor->hwEta()) {
214  // // neighborNW_et = neighborET;
215  // // nNeighbors++;
216  // // continue;
217  // // }
218  // else if(deltaGctPhi(*region, *neighbor) == -1 &&
219  // (region->hwEta() + 1) == neighbor->hwEta()) {
220  // neighborSE_et = neighborET;
221  // nNeighbors++;
222  // continue;
223  // }
224  // }
225  // if(regionET > neighborN_et &&
226  // //regionET > neighborNW_et &&
227  // //regionET > neighborW_et &&
228  // //regionET > neighborSW_et &&
229  // regionET >= neighborNE_et &&
230  // regionET >= neighborE_et &&
231  // regionET >= neighborSE_et &&
232  // regionET >= neighborS_et) {
233  // unsigned int jetET = regionET +
234  // neighborN_et + neighborS_et + neighborE_et + /*neighborW_et +*/
235  // neighborNE_et + /*neighborSW_et +*/ neighborSE_et;// + neighborNW_et;
236  // /*
237  // int jetPhi = region->hwPhi() * 4 +
238  // ( - 2 * (neighborS_et + neighborSE_et + neighborSW_et)
239  // + 2 * (neighborN_et + neighborNE_et + neighborNW_et) ) / jetET;
240  // if(jetPhi < 0) {
241 
242  // }
243  // else if(jetPhi >= ((int) N_JET_PHI)) {
244  // jetPhi -= N_JET_PHI;
245  // }
246  // int jetEta = region->hwEta() * 4 +
247  // ( - 2 * (neighborW_et + neighborNW_et + neighborSW_et)
248  // + 2 * (neighborE_et + neighborNE_et + neighborSE_et) ) / jetET;
249  // if(jetEta < 0) jetEta = 0;
250  // if(jetEta >= ((int) N_JET_ETA)) jetEta = N_JET_ETA - 1;
251  // */
252  // // Temporarily use the region granularity -- we will try to improve as above when code is debugged
253  // int jetPhi = region->hwPhi();
254  // int jetEta = region->hwEta();
255 
256  // bool neighborCheck = (nNeighbors == 5);
257  // // On the eta edge we only expect 5 neighbors
258  // // if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
259  // // neighborCheck = true;
260 
261  // if (!neighborCheck) {
262  // std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
263  // assert(false);
264  // }
265 
266  // const bool forward = false; //by definition
267  // int jetQual = 0;
268  // if(forward)
269  // jetQual |= 0x2;
270 
271  // ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > *jetLorentz =
272  // new ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> >();
273  // l1t::Jet theJet(*jetLorentz, jetET, jetEta, jetPhi, jetQual);
274  // //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
275 
276  // jets->push_back(theJet);
277  // }
278  // }
279 
280  // // separate loops for the central jets at the edges of HF, composed of 6 regions only.
281  // for(std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
282  // double regionET = region->hwPt();
283  // // look at only the right eta wall
284  // if (region->hwEta() != 17) continue;
285  // if (regionET < jetSeedThreshold) continue;
286  // double neighborN_et = 0;
287  // double neighborS_et = 0;
288  // //double neighborE_et = 0;
289  // double neighborW_et = 0;
290  // //double neighborNE_et = 0;
291  // double neighborSW_et = 0;
292  // double neighborNW_et = 0;
293  // //double neighborSE_et = 0;
294  // unsigned int nNeighbors = 0;
295  // for(std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end(); neighbor++) {
296 
297  // double neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
298  // if(deltaGctPhi(*region, *neighbor) == 1 &&
299  // (region->hwEta() ) == neighbor->hwEta()) {
300  // neighborN_et = neighborET;
301  // nNeighbors++;
302  // continue;
303  // }
304  // else if(deltaGctPhi(*region, *neighbor) == -1 &&
305  // (region->hwEta() ) == neighbor->hwEta()) {
306  // neighborS_et = neighborET;
307  // nNeighbors++;
308  // continue;
309  // }
310  // // else if(deltaGctPhi(*region, *neighbor) == 0 &&
311  // // (region->hwEta() + 1) == neighbor->hwEta()) {
312  // // neighborE_et = neighborET;
313  // // nNeighbors++;
314  // // continue;
315  // // }
316  // else if(deltaGctPhi(*region, *neighbor) == 0 &&
317  // (region->hwEta() - 1) == neighbor->hwEta()) {
318  // neighborW_et = neighborET;
319  // nNeighbors++;
320  // continue;
321  // }
322  // // else if(deltaGctPhi(*region, *neighbor) == 1 &&
323  // // (region->hwEta() + 1) == neighbor->hwEta()) {
324  // // neighborNE_et = neighborET;
325  // // nNeighbors++;
326  // // continue;
327  // // }
328  // else if(deltaGctPhi(*region, *neighbor) == -1 &&
329  // (region->hwEta() - 1) == neighbor->hwEta()) {
330  // neighborSW_et = neighborET;
331  // nNeighbors++;
332  // continue;
333  // }
334  // else if(deltaGctPhi(*region, *neighbor) == 1 &&
335  // (region->hwEta() - 1) == neighbor->hwEta()) {
336  // neighborNW_et = neighborET;
337  // nNeighbors++;
338  // continue;
339  // }
340  // // else if(deltaGctPhi(*region, *neighbor) == -1 &&
341  // // (region->hwEta() + 1) == neighbor->hwEta()) {
342  // // neighborSE_et = neighborET;
343  // // nNeighbors++;
344  // // continue;
345  // // }
346  // }
347  // if(//regionET > neighborN_et &&
348  // regionET > neighborNW_et &&
349  // regionET > neighborW_et &&
350  // regionET > neighborSW_et &&
351  // //regionET >= neighborNE_et &&
352  // //regionET >= neighborE_et &&
353  // //regionET >= neighborSE_et &&
354  // regionET >= neighborS_et) {
355  // unsigned int jetET = regionET +
356  // neighborN_et + neighborS_et + /*neighborE_et + */neighborW_et +
357  // /*neighborNE_et +*/ neighborSW_et + /*neighborSE_et + */neighborNW_et;
358  // /*
359  // int jetPhi = region->hwPhi() * 4 +
360  // ( - 2 * (neighborS_et + neighborSE_et + neighborSW_et)
361  // + 2 * (neighborN_et + neighborNE_et + neighborNW_et) ) / jetET;
362  // if(jetPhi < 0) {
363 
364  // }
365  // else if(jetPhi >= ((int) N_JET_PHI)) {
366  // jetPhi -= N_JET_PHI;
367  // }
368  // int jetEta = region->hwEta() * 4 +
369  // ( - 2 * (neighborW_et + neighborNW_et + neighborSW_et)
370  // + 2 * (neighborE_et + neighborNE_et + neighborSE_et) ) / jetET;
371  // if(jetEta < 0) jetEta = 0;
372  // if(jetEta >= ((int) N_JET_ETA)) jetEta = N_JET_ETA - 1;
373  // */
374  // // Temporarily use the region granularity -- we will try to improve as above when code is debugged
375  // int jetPhi = region->hwPhi();
376  // int jetEta = region->hwEta();
377 
378  // bool neighborCheck = (nNeighbors == 5);
379  // // On the eta edge we only expect 5 neighbors
380  // // if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
381  // // neighborCheck = true;
382 
383  // if (!neighborCheck) {
384  // std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
385  // assert(false);
386  // }
387 
388  // const bool forward = false; //by definition
389  // int jetQual = 0;
390  // if(forward)
391  // jetQual |= 0x2;
392 
393  // ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > *jetLorentz =
394  // new ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> >();
395  // l1t::Jet theJet(*jetLorentz, jetET, jetEta, jetPhi, jetQual);
396  // //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
397 
398  // jets->push_back(theJet);
399  // }
400  // }
401 
402  //the jets should be sorted, highest pT first.
403  // do not truncate the jet list, GT converter handles that
404  std::sort(uncalibjets->begin(), uncalibjets->end(), compareJets);
405  std::reverse(uncalibjets->begin(), uncalibjets->end());
406  }
407 }
i
int i
Definition: DBlmapReader.cc:9
l1t::L1Candidate::hwPhi
int hwPhi() const
Definition: L1Candidate.cc:79
l1t::Jet
Definition: Jet.h:13
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
j
int j
Definition: DBlmapReader.cc:9
l1t::deltaGctPhi
int deltaGctPhi(const CaloRegion &region, const CaloRegion &neighbor)
Definition: JetFinderMethods.cc:19
l1t::L1Candidate::hwPt
int hwPt() const
Definition: L1Candidate.cc:69
l1t::slidingWindowJetFinder
void slidingWindowJetFinder(const int, const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
Definition: JetFinderMethods.cc:34
l1t::compareJets
bool compareJets(l1t::Jet i, l1t::Jet j)
Definition: JetFinderMethods.cc:30
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
L1CaloRegionDetId::N_PHI
static const unsigned N_PHI
Definition: L1CaloRegionDetId.h:22
Generated for CMSSW Reference Manual by doxygen 1.8.5