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