ATLAS Offline Software
TTBarElectronJetOverlap.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
3  */
4 
5 #include "TopJetSubstructure/TTBarElectronJetOverlap.h"
6 
7 #include <iostream>
8 #include <cmath>
9 
10 #include "xAODJet/JetContainer.h"
11 #include "xAODEgamma/ElectronContainer.h"
12 
13 
14 using namespace std;
15 
16 void TTBarElectronJetOverlap::Load(xAOD::JetContainer* jets, xAOD::ElectronContainer* electrons,
17  const std::string& leptonDef) {
18  m_jets = jets;
19  m_electrons = electrons;
20  fDebug = false;
21 
22  size_t se = m_electrons->size();
23  fElClTLVs.resize(se);
24  fElTLVs.resize(se);
25  fElGood.resize(se);
26 
27  size_t sj = m_jets->size();
28  fJetTLVs.resize(sj);
29  fOrigJetTLVs.resize(sj);
30  fJetJVFs.resize(sj);
31  fJetD3PDTrkPtSums.resize(sj);
32  fJetD3PDTrkPtPVSums.resize(sj);
33 
34  for (size_t i = 0; i < sj; i++) {
35  fJetTLVs[i].SetPtEtaPhiE(m_jets->at(i)->pt(), m_jets->at(i)->eta(), m_jets->at(i)->phi(), m_jets->at(i)->e());
36  fOrigJetTLVs[i].SetPtEtaPhiE(m_jets->at(i)->pt(), m_jets->at(i)->eta(), m_jets->at(i)->phi(), m_jets->at(i)->e());
37 
38  //std::vector<float> JVF = m_jets->at(i)->getAttribute<std::vector<float> >("JVF");
39  //fJetJVFs[i] = JVF.size() > 0 ? JVF[0] : -1; // (TODO) NOTA BENE: Assumes PV is the first vertex
40  //std::vector<float> trkPtSumVec = m_jets->at(i)->getAttribute<std::vector<float> >("SumPtTrkPt500");
41  //float trkPtSum = trkPtSumVec.size() > 0 ? trkPtSumVec[0] : 0; // (TODO) NOTA BENE: Assumes PV is the first vertex
42  //fJetD3PDTrkPtSums[i] = fJetJVFs[i] ? trkPtSum/JVF[0] : -1;
43  //fJetD3PDTrkPtPVSums[i] = trkPtSum;
44  if (fDebug) {
45  cout << " - (pre-OR) Jet # " << i <<
46  " Pt Eta Phi: " <<
47  fJetTLVs[i].Pt() << " " <<
48  fJetTLVs[i].Eta() << " " <<
49  fJetTLVs[i].Phi() << endl;
50  }
51  }
52 
53  for (size_t i = 0; i < se; i++) {
54  fElGood[i] = m_electrons->at(i)->auxdataConst<char>(leptonDef.c_str());
55 
56  if (m_electrons->at(i)->caloCluster()) {
57  fElClTLVs[i].SetPtEtaPhiM(m_electrons->at(i)->caloCluster()->e() /
58  std::cosh(m_electrons->at(i)->caloCluster()->eta()),
59  m_electrons->at(i)->caloCluster()->eta(), m_electrons->at(
60  i)->caloCluster()->phi(), 0.511);
61  } else {
62  fElClTLVs[i] = m_electrons->at(i)->p4();
63  }
64 
65  if (m_electrons->at(i)->trackParticle() && m_electrons->at(i)->caloCluster()) {
66  fElTLVs[i].SetPtEtaPhiM(m_electrons->at(i)->caloCluster()->e() /
67  std::cosh(m_electrons->at(i)->trackParticle()->eta()),
68  m_electrons->at(i)->trackParticle()->eta(), m_electrons->at(
69  i)->trackParticle()->phi(), 0.511);
70  } else {
71  fElTLVs[i] = m_electrons->at(i)->p4();
72  }
73  if (fDebug) {
74  cout << " - (pre-OR) El # " << i <<
75  " Pt Eta Phi: " <<
76  fElTLVs[i].Pt() << " " <<
77  fElTLVs[i].Eta() << " " <<
78  fElTLVs[i].Phi() << endl;
79  }
80  }
81 }
82 
83 void TTBarElectronJetOverlap::FindAssocEls() {
84  size_t nJets = m_jets->size();
85 
86  fJetAssocElCls = vector<set<int> >(nJets);
87  size_t nEls = m_electrons->size();
88  fElClAssocJet = vector<int>(nEls, -1);
89 
90  // find the associated electron clusters in each jet
91  double drmin, dr;
92  TLorentzVector elcl, jet;
93  for (size_t z = 0; z < nEls; ++z) {
94  if (!fElGood[z]) continue;
95  elcl = fElClTLVs[z];
96  drmin = 0.4;
97  int drmin_idx = -1;
98  for (size_t y = 0; y < nJets; ++y) {
99  jet = fJetTLVs[y];
100  dr = elcl.DeltaR(jet);
101  // attempt to match the jet to the electron.
102  if (dr < drmin) {
103  drmin = dr;
104  drmin_idx = y;
105  }
106  }
107  if (drmin_idx == -1) continue;
108  fJetAssocElCls[drmin_idx].insert(z);
109  fElClAssocJet[z] = drmin_idx;
110  }
111 }
112 
113 void TTBarElectronJetOverlap::SubtractEls() {
114  size_t nJets = m_jets->size();
115  int ElIdx;
116 
117  fSubJets.resize(nJets);
118  for (unsigned int k = 0; k < nJets; ++k) {
119  fSubJets[k] = 0;
120  }
121 
122  for (size_t iJet = 0; iJet < nJets; iJet++) {
123  TLorentzVector originalJet = fOrigJetTLVs[iJet];
124  for (set<int>::iterator iEl = fJetAssocElCls[iJet].begin(); iEl != fJetAssocElCls[iJet].end(); ++iEl) {
125  ElIdx = *iEl;
126  TLorentzVector elcorr = fElTLVs[ElIdx];
127  fJetTLVs[iJet] -= elcorr;
128  fSubJets[iJet] = 1;
129  }
130  }
131 }
132 
133 void TTBarElectronJetOverlap::FindGoodObjects() {
134  size_t nEls = m_electrons->size();
135 
136  fGoodEls = vector<bool>(nEls, true);
137 
138  size_t nJets = fJetTLVs.size();
139  vector<TLorentzVector> TmpJetTLVs = fJetTLVs;
140  vector<TLorentzVector> TmpOrigJetTLVs = fOrigJetTLVs;
141  TLorentzVector el, jet, jetOrig;
142 
143  for (size_t iEl = 0; iEl < nEls; iEl++) {
144  if (!fElGood[iEl]) continue;
145 
146  el = fElTLVs[iEl];
147 
148  if (fDebug) {
149  cout << " - (good) El # " << iEl <<
150  " Pt Eta Phi: " <<
151  el.Pt() << " " <<
152  el.Eta() << " " <<
153  el.Phi() << endl;
154  }
155 
156  for (size_t iJet = 0; iJet < nJets; iJet++) {
157  jet = fJetTLVs[iJet];
158  jetOrig = fOrigJetTLVs[iJet];
159  if (fDebug) {
160  cout << " - (good) Jet # " << iJet <<
161  " Pt Eta Phi M: " <<
162  jet.Pt() << " " <<
163  jet.Eta() << " " <<
164  jet.Phi() << " " <<
165  jet.M() << endl;
166  }
167 
168  //if (jet.Pt() < 25e3)
169  //continue;
170 
171  // New overlap removal procedure:
172  // subtracted ets failing the pT cut are assumed to be a result of electron energy deposits. Continue to the next
173  // jet.
174 
175  if ((jetOrig.Pt() >= 30e3 && jetOrig.Pt() < 50e3) && jet.Pt() < 20e3) {
176  continue;
177  } else if ((jetOrig.Pt() >= 50e3 && jetOrig.Pt() < 100e3) && jet.Pt() < 30e3) {
178  continue;
179  } else if ((jetOrig.Pt() >= 100e3 && jetOrig.Pt() < 400e3) && jet.Pt() < 90e3) {
180  continue;
181  } else if ((jetOrig.Pt() >= 400e3 && jetOrig.Pt() < 800e3) && jet.Pt() < 130e3) {
182  continue;
183  } else if ((jetOrig.Pt() >= 800e3 && jetOrig.Pt() < 1200e3) && jet.Pt() < 230e3) {
184  continue;
185  } else if ((jetOrig.Pt() >= 1200e3 && jetOrig.Pt() < 1600e3) && jet.Pt() < 700e3) {
186  continue;
187  } else if ((jetOrig.Pt() >= 1600e3 && jetOrig.Pt() < 2000e3) && jet.Pt() < 1250e3) {
188  continue;
189  } else if ((jetOrig.Pt() >= 2000e3 && jetOrig.Pt() < 2200e3) && jet.Pt() < 1300e3) {
190  continue;
191  } else if (jetOrig.Pt() > 2200e3 && jet.Pt() < 1350e3) {
192  continue;
193  }
194 
195  // if the electron is too close to a subtracted jet...
196  if (jet.DeltaR(el) < 0.2) {
197  if (fDebug) cout << " El too close to jet. Removing." << endl;
198 
199  // remove from good electrons list.
200  fGoodEls[iEl] = false;
201 
202  if (fElClAssocJet[iEl] >= 0) {
203  // add electron 4-vector back to jet.
204  TmpJetTLVs[fElClAssocJet[iEl]] += el;
205  //remove from subtracted jets list
206  fSubJets[iJet] = 0;
207  // remove this electron from the association.
208  fJetAssocElCls[iJet].erase(iEl);
209  }
210  }
211  }
212  }
213 
214  fJetTLVs = TmpJetTLVs;
215  fOrigJetTLVs = TmpOrigJetTLVs;
216  fGoodJets = vector<bool>(nJets, true);
217 
218  for (size_t iJet = 0; iJet < nJets; iJet++) {
219  jetOrig = fOrigJetTLVs[iJet];
220 
221  if (jetOrig.Pt() >= 30e3 && jetOrig.Pt() < 50e3) {
222  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 20e3;
223  } else if (jetOrig.Pt() >= 50e3 && jetOrig.Pt() < 100e3) {
224  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 30e3;
225  } else if (jetOrig.Pt() >= 100e3 && jetOrig.Pt() < 400e3) {
226  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 90e3;
227  } else if (jetOrig.Pt() >= 400e3 && jetOrig.Pt() < 800e3) {
228  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 130e3;
229  } else if (jetOrig.Pt() >= 800e3 && jetOrig.Pt() < 1200e3) {
230  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 230e3;
231  } else if (jetOrig.Pt() >= 1200e3 && jetOrig.Pt() < 1600e3) {
232  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 700e3;
233  } else if (jetOrig.Pt() >= 1600e3 && jetOrig.Pt() < 2000e3) {
234  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 1250e3;
235  } else if (jetOrig.Pt() >= 2000e3 && jetOrig.Pt() < 2200e3) {
236  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 1300e3;
237  } else if (jetOrig.Pt() > 2200e3) {
238  fGoodJets[iJet] = fJetTLVs[iJet].Pt() > 1350e3;
239  }
240  }
241 }
242 
243 void TTBarElectronJetOverlap::RecalcJVF() {
244  size_t nJets = m_jets->size();
245  int ElIdx;
246  //int TrkIdx;
247  TLorentzVector trk;
248 
249  for (size_t iJet = 0; iJet < nJets; iJet++) {
250  // if jvf == 0 or -1, no need to recalulate...
251  if (!fJetJVFs[iJet] ||
252  fJetJVFs[iJet] < 0) continue;
253 
254  for (set<int>::iterator iEl = fJetAssocElCls[iJet].begin();
255  iEl != fJetAssocElCls[iJet].end(); ++iEl) {
256  ElIdx = *iEl;
257 
258  trk.SetPtEtaPhiE(m_electrons->at(ElIdx)->trackParticle()->pt(),
259  m_electrons->at(ElIdx)->trackParticle()->eta(), m_electrons->at(
260  ElIdx)->trackParticle()->phi(), m_electrons->at(ElIdx)->trackParticle()->e());
261 
262  // not in the associated tracks.
263  bool foundMatch = false;
264  std::vector<const xAOD::TrackParticle*> jetTracks;
265  //bool haveJetTracks = m_jets->at(iJet)->getAssociatedObjects(xAOD::JetAttribute::GhostTrack, jetTracks);
266  for (size_t t = 0; t < jetTracks.size(); ++t) {
267  TLorentzVector jet_trk;
268  jet_trk.SetPtEtaPhiE(jetTracks[t]->pt(), jetTracks[t]->eta(), jetTracks[t]->phi(), jetTracks[t]->e());
269  if (jet_trk.DeltaR(trk) < 0.01) {
270  foundMatch = true;
271  }
272  }
273 
274  if (!foundMatch) continue;
275 
276  // recompute JVF.
277  fJetD3PDTrkPtSums[iJet] -= trk.Pt();
278 
279  if (std::fabs(std::sin(m_electrons->at(ElIdx)->trackParticle()->theta()) *
280  m_electrons->at(ElIdx)->trackParticle()->z0()) < 1.0 &&
281  std::fabs(m_electrons->at(ElIdx)->trackParticle()->d0()) < 1.0) fJetD3PDTrkPtPVSums[iJet] -= trk.Pt();
282  }
283 
284  // we subtracted too much for some reason?
285  if (fJetD3PDTrkPtSums[iJet] < 0 ||
286  fJetD3PDTrkPtPVSums[iJet] < 0) fJetJVFs[iJet] = 0;
287  // no tracks associated with the jet.
288  else if (fJetD3PDTrkPtSums[iJet] == 0) fJetJVFs[iJet] = -1;
289  // all good. recalculate JVF.
290  else fJetJVFs[iJet] = fJetD3PDTrkPtPVSums[iJet] / fJetD3PDTrkPtSums[iJet];
291  }
292 }
293 
294 void TTBarElectronJetOverlap::AnalyzeEvent(const std::string& leptonDef) {
295  FindAssocEls();
296  SubtractEls();
297  FindGoodObjects();
298  //RecalcJVF();
299 
300  // put the results back in place
301  size_t se = m_electrons->size();
302  size_t sj = m_jets->size();
303  for (size_t i = 0; i < sj; i++) {
304  m_jets->at(i)->setJetP4(xAOD::JetFourMom_t(fJetTLVs[i].Perp(), fJetTLVs[i].Eta(), fJetTLVs[i].Phi(),
305  fJetTLVs[i].M()));
306 
307  if (fGoodJets[i]) m_jets->at(i)->auxdata<char>("passesFancyOR") = 1; // Subtracted
308  // jets
309  // that
310  // pass
311  // the
312  // OR
313  // test
314  // and
315  // are
316  // kept.
317  else m_jets->at(i)->auxdata<char>("passesFancyOR") = 0;
318 
319  if (fSubJets[i] == 1) m_jets->at(i)->auxdata<char>("subtractedJet") = 1; // Jets
320  // close
321  // to
322  // electrons
323  // that
324  // undergo
325  // the
326  // OR
327  // test.
328  else m_jets->at(i)->auxdata<char>("subtractedJet") = 0;
329 
330 
331  //std::vector<float> newJVF;
332  //newJVF.push_back(fJetJVFs[i]);
333  //m_jets->at(i)->setAttribute<std::vector<float> >("JVF", newJVF);
334  if (fDebug) {
335  cout << " - (final) Jet # " << i <<
336  " Pt Eta Phi M: " <<
337  m_jets->at(i)->pt() << " " <<
338  m_jets->at(i)->eta() << " " <<
339  m_jets->at(i)->phi() << " " <<
340  m_jets->at(i)->m() << " good " << (int) m_jets->at(i)->auxdata<char>(leptonDef.c_str()) << endl;
341  }
342  }
343 
344 
345  for (size_t j = 0; j < se; j++) {
346  if (!fGoodEls[j] || !fElGood[j]) {
347  m_electrons->at(j)->auxdata<char>(leptonDef.c_str()) = 0;
348  } else {
349  m_electrons->at(j)->auxdata<char>(leptonDef.c_str()) = 1;
350  }
351 
352  if (fDebug) {
353  cout << " - (final) El # " << j <<
354  " Pt Eta Phi M: " <<
355  m_electrons->at(j)->pt() << " " <<
356  m_electrons->at(j)->eta() << " " <<
357  m_electrons->at(j)->phi() << " " <<
358  m_electrons->at(j)->m() << " good " << (int) m_electrons->at(j)->auxdata<char>(leptonDef.c_str()) << endl;
359  }
360  }
361 }
TTBarElectronJetOverlap::FindGoodObjects
void FindGoodObjects()
Definition: TTBarElectronJetOverlap.cxx:133
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
PlotCalibFromCool.begin
begin
Definition: PlotCalibFromCool.py:94
test_pyathena.pt
pt
Definition: test_pyathena.py:11
Phi
@ Phi
Definition: RPCdef.h:8
python.TurnDataReader.dr
dr
Definition: TurnDataReader.py:112
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
keylayer_zslicemap.se
se
Definition: keylayer_zslicemap.py:194
TTBarElectronJetOverlap::Load
void Load(xAOD::JetContainer *jets, xAOD::ElectronContainer *electrons, const std::string &leptonDef)
Definition: TTBarElectronJetOverlap.cxx:16
TTBarElectronJetOverlap::FindAssocEls
void FindAssocEls()
Definition: TTBarElectronJetOverlap.cxx:83
jet
Definition: JetCalibTools_PlotJESFactors.cxx:23
CheckAppliedSFs.e3
e3
Definition: CheckAppliedSFs.py:264
ElectronContainer.h
TTBarElectronJetOverlap.h
lumiFormat.i
int i
Definition: lumiFormat.py:92
z
#define z
TTBarElectronJetOverlap::AnalyzeEvent
void AnalyzeEvent(const std::string &leptonDef)
Definition: TTBarElectronJetOverlap.cxx:294
plotIsoValidation.el
el
Definition: plotIsoValidation.py:197
DataVector
Derived DataVector<T>.
Definition: DataVector.h:581
xAOD::JetFourMom_t
ROOT::Math::LorentzVector< ROOT::Math::PtEtaPhiM4D< double > > JetFourMom_t
Base 4 Momentum type for Jet.
Definition: JetTypes.h:17
TTBarElectronJetOverlap::SubtractEls
void SubtractEls()
Definition: TTBarElectronJetOverlap.cxx:113
DiTauMassTools::MaxHistStrategyV2::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
y
#define y
JetContainer.h
defineDB.jets
list jets
Definition: JetTagCalibration/share/defineDB.py:24
TTBarElectronJetOverlap::RecalcJVF
void RecalcJVF()
Definition: TTBarElectronJetOverlap.cxx:243
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
InDetDD::electrons
@ electrons
Definition: InDetDD_Defs.h:17
Eta
@ Eta
Definition: RPCdef.h:8
fitman.k
k
Definition: fitman.py:528
Generated on Fri Jul 5 2024 21:32:13 for ATLAS Offline Software by   doxygen 1.8.18