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 
11 
13 
14 #include <vector>
15 
16 namespace l1t {
17 
18  int deltaGctPhi(const CaloRegion& region, const CaloRegion& neighbor) {
19  int phi1 = region.hwPhi();
20  int phi2 = neighbor.hwPhi();
21  int diff = phi1 - phi2;
22  if (std::abs(phi1 - phi2) == L1CaloRegionDetId::N_PHI - 1) { //18 regions in phi
23  diff = -diff / std::abs(diff);
24  }
25  return diff;
26  }
27 
28  // turn each central region into a jet
29  void passThroughJets(const std::vector<l1t::CaloRegion>* regions, std::vector<l1t::Jet>* uncalibjets) {
30  for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
31  int jetQual = 0;
32  if (region->hwEta() < 4 || region->hwEta() > 17)
33  jetQual = 2;
34  int jetET = region->hwPt();
35  int jetEta = region->hwEta();
36  int jetPhi = region->hwPhi();
37 
38  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
39  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
40  uncalibjets->push_back(theJet);
41  }
42  }
43 
45  const std::vector<l1t::CaloRegion>* regions,
46  std::vector<l1t::Jet>* uncalibjets) {
47  for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
48  int regionET = region->hwPt(); //regionPhysicalEt(*region);
49  if (regionET <= jetSeedThreshold)
50  continue;
51  int neighborN_et = 0;
52  int neighborS_et = 0;
53  int neighborE_et = 0;
54  int neighborW_et = 0;
55  int neighborNE_et = 0;
56  int neighborSW_et = 0;
57  int neighborNW_et = 0;
58  int neighborSE_et = 0;
59  unsigned int nNeighbors = 0;
60  for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
61  neighbor++) {
62  int neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
63  if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
64  neighborN_et = neighborET;
65  nNeighbors++;
66  continue;
67  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
68  neighborS_et = neighborET;
69  nNeighbors++;
70  continue;
71  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
72  neighborE_et = neighborET;
73  nNeighbors++;
74  continue;
75  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
76  neighborW_et = neighborET;
77  nNeighbors++;
78  continue;
79  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
80  neighborNE_et = neighborET;
81  nNeighbors++;
82  continue;
83  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
84  neighborSW_et = neighborET;
85  nNeighbors++;
86  continue;
87  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
88  neighborNW_et = neighborET;
89  nNeighbors++;
90  continue;
91  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
92  neighborSE_et = neighborET;
93  nNeighbors++;
94  continue;
95  }
96  }
97  if (regionET > neighborN_et && regionET > neighborNW_et && regionET > neighborW_et && regionET > neighborSW_et &&
98  regionET >= neighborNE_et && regionET >= neighborE_et && regionET >= neighborSE_et &&
99  regionET >= neighborS_et) {
100  unsigned int jetET = regionET + neighborN_et + neighborS_et + neighborE_et + neighborW_et + neighborNE_et +
101  neighborSW_et + neighborSE_et + neighborNW_et;
102 
103  int jetPhi = region->hwPhi();
104  int jetEta = region->hwEta();
105 
106  bool neighborCheck = (nNeighbors == 8);
107  // On the eta edge we only expect 5 neighbors
108  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
109  neighborCheck = true;
110 
111  if (!neighborCheck) {
112  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
113  assert(false);
114  }
115 
116  //first iteration, eta cut defines forward
117  //const bool forward = (jetEta <= 4 || jetEta >= 17);
118  const bool forward = (jetEta < 4 || jetEta > 17);
119  int jetQual = 0;
120  if (forward)
121  jetQual |= 0x2;
122 
123  // check for input overflow regions
124  if (forward && regionET == 255) {
125  jetET = 1023; // 10 bit max
126  } else if (!forward && regionET == 1023) {
127  jetET = 1023; // 10 bit max
128  } else if (region->hwEta() == 17) {
129  if (neighborNE_et == 255 || neighborE_et == 255 || neighborSE_et == 255)
130  jetET = 1023; // 10 bit max
131  } else if (region->hwEta() == 4) {
132  if (neighborNW_et == 255 || neighborW_et == 255 || neighborSW_et == 255)
133  jetET = 1023; // 10 bit max
134  }
135 
136  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
137  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
138  //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
139 
140  uncalibjets->push_back(theJet);
141  }
142  }
143  }
144 
146  const std::vector<l1t::CaloRegion>* regions,
147  std::vector<l1t::Jet>* uncalibjets) {
148  for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
149  int regionET = region->hwPt(); //regionPhysicalEt(*region);
150  if (regionET < jetSeedThreshold)
151  continue;
152  int neighborN_et = 0;
153  int neighborS_et = 0;
154  int neighborE_et = 0;
155  int neighborW_et = 0;
156  int neighborNE_et = 0;
157  int neighborSW_et = 0;
158  int neighborNW_et = 0;
159  int neighborSE_et = 0;
160  unsigned int nNeighbors = 0;
161  for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
162  neighbor++) {
163  int neighborET = neighbor->hwPt(); //regionPhysicalEt(*neighbor);
164  if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
165  neighborN_et = neighborET;
166  nNeighbors++;
167  continue;
168  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
169  neighborS_et = neighborET;
170  nNeighbors++;
171  continue;
172  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
173  neighborE_et = neighborET;
174  nNeighbors++;
175  continue;
176  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
177  neighborW_et = neighborET;
178  nNeighbors++;
179  continue;
180  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
181  neighborNE_et = neighborET;
182  nNeighbors++;
183  continue;
184  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
185  neighborSW_et = neighborET;
186  nNeighbors++;
187  continue;
188  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
189  neighborNW_et = neighborET;
190  nNeighbors++;
191  continue;
192  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
193  neighborSE_et = neighborET;
194  nNeighbors++;
195  continue;
196  }
197  }
198  unsigned int jetET = regionET + neighborN_et + neighborS_et + neighborE_et + neighborW_et + neighborNE_et +
199  neighborSW_et + neighborSE_et + neighborNW_et;
200 
201  int jetPhi = region->hwPhi();
202  int jetEta = region->hwEta();
203 
204  bool neighborCheck = (nNeighbors == 8);
205  // On the eta edge we only expect 5 neighbors
206  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
207  neighborCheck = true;
208 
209  if (!neighborCheck) {
210  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
211  assert(false);
212  }
213 
214  //first iteration, eta cut defines forward
215  //const bool forward = (jetEta <= 4 || jetEta >= 17);
216  const bool forward = (jetEta < 4 || jetEta > 17);
217  int jetQual = 0;
218  if (forward)
219  jetQual |= 0x2;
220 
221  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
222  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
223  //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
224 
225  uncalibjets->push_back(theJet);
226  }
227  }
228 
230  const int etaMask,
231  const std::vector<l1t::CaloRegion>* regions,
232  std::vector<l1t::Jet>* uncalibjets) {
233  for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
234  int regionET = region->hwPt();
235  if (regionET <= jetSeedThreshold)
236  continue;
237  int subEta = region->hwEta();
238  if ((etaMask & (1 << subEta)) >> subEta)
239  regionET = 0;
240  int neighborN_et = 0;
241  int neighborS_et = 0;
242  int neighborE_et = 0;
243  int neighborW_et = 0;
244  int neighborNE_et = 0;
245  int neighborSW_et = 0;
246  int neighborNW_et = 0;
247  int neighborSE_et = 0;
248  unsigned int nNeighbors = 0;
249  for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
250  neighbor++) {
251  int neighborET = neighbor->hwPt();
252  int subEta2 = neighbor->hwEta();
253  if ((etaMask & (1 << subEta2)) >> subEta2)
254  neighborET = 0;
255 
256  if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
257  neighborN_et = neighborET;
258  nNeighbors++;
259  continue;
260  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
261  neighborS_et = neighborET;
262  nNeighbors++;
263  continue;
264  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
265  neighborE_et = neighborET;
266  nNeighbors++;
267  continue;
268  } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
269  neighborW_et = neighborET;
270  nNeighbors++;
271  continue;
272  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
273  neighborNE_et = neighborET;
274  nNeighbors++;
275  continue;
276  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
277  neighborSW_et = neighborET;
278  nNeighbors++;
279  continue;
280  } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
281  neighborNW_et = neighborET;
282  nNeighbors++;
283  continue;
284  } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
285  neighborSE_et = neighborET;
286  nNeighbors++;
287  continue;
288  }
289  }
290  if (regionET > neighborN_et && regionET > neighborNW_et && regionET > neighborW_et && regionET > neighborSW_et &&
291  regionET >= neighborNE_et && regionET >= neighborE_et && regionET >= neighborSE_et &&
292  regionET >= neighborS_et) {
293  // use the highest-pT 2x2 jet inside this 3x3
294  unsigned int jetET_NW;
295  unsigned int jetET_NE;
296  unsigned int jetET_SW;
297  unsigned int jetET_SE;
298 
299  jetET_NW = regionET + neighborW_et + neighborNW_et + neighborN_et;
300  jetET_NE = regionET + neighborE_et + neighborNE_et + neighborN_et;
301  jetET_SW = regionET + neighborS_et + neighborSW_et + neighborW_et;
302  jetET_SE = regionET + neighborS_et + neighborSE_et + neighborE_et;
303 
304  unsigned int jetET = std::max(jetET_NW, jetET_NE);
305  jetET = std::max(jetET, jetET_SW);
306  jetET = std::max(jetET, jetET_SE);
307 
308  int jetPhi = region->hwPhi();
309  int jetEta = region->hwEta();
310 
311  bool neighborCheck = (nNeighbors == 8);
312  // On the eta edge we only expect 5 neighbor
313  if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
314  neighborCheck = true;
315 
316  if (!neighborCheck) {
317  std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
318  assert(false);
319  }
320 
321  //first iteration, eta cut defines forward
322  const bool forward = (jetEta < 4 || jetEta > 17);
323  int jetQual = 0;
324  if (forward)
325  jetQual |= 0x2;
326 
327  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
328  l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
329  uncalibjets->push_back(theJet);
330  }
331  }
332  }
333 } // namespace l1t
int hwPhi() const
Definition: L1Candidate.h:37
delete x;
Definition: CaloConfig.h:22
void TwelveByTwelveFinder(const int, const std::vector< l1t::CaloRegion > *regions, std::vector< l1t::Jet > *uncalibjets)
assert(be >=bs)
Definition: Jet.h:20
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)
static const unsigned N_PHI