Loading [MathJax]/extensions/tex2jax.js
ATLAS Offline Software
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
DiTauTruthMatchingTool.cxx
Go to the documentation of this file.
1 
10 // Local include(s)
12 
13 // Core include(s):
14 #include "AthLinks/ElementLink.h"
18 
19 // EDM include(s):
21 #include "xAODMuon/MuonContainer.h"
22 
24 
25 using namespace TauAnalysisTools;
26 
27 //=================================PUBLIC-PART==================================
28 //______________________________________________________________________________
29 DiTauTruthMatchingTool::DiTauTruthMatchingTool( const std::string& name )
31  , m_accPtVis("pt_vis")
32  , m_accEtaVis("eta_vis")
33  , m_accPhiVis("phi_vis")
34  , m_accMVis("m_vis")
35 {
36 }
37 
38 //______________________________________________________________________________
39 DiTauTruthMatchingTool::~DiTauTruthMatchingTool( )
40 {
41 
42 }
43 
44 //______________________________________________________________________________
46 {
47  ATH_MSG_INFO( "Initializing DiTauTruthMatchingTool" );
48 
49  // configure BuildTruthTaus in truth matching mode, not truth tau building mode
50  DiTauTruthMatchingTool::BuildTruthTaus::setTruthMatchingMode();
51 
53  {
54  ATH_MSG_FATAL("Failed initializing BuildTruthTaus");
55  return StatusCode::FAILURE;
56  }
57  return StatusCode::SUCCESS;
58 }
59 
60 //______________________________________________________________________________
62 {
63  if (retrieveTruthTaus().isFailure())
64  return;
65 
66  if (findTruthTau(xDiTau).isFailure())
67  ATH_MSG_WARNING("There was a failure in finding the matched truth tau");
68 
69  return;
70 }
71 
72 //______________________________________________________________________________
73 void DiTauTruthMatchingTool::getTruth(const std::vector<const xAOD::DiTauJet*>& vDiTaus)
74 {
75  for (auto xDiTau : vDiTaus)
76  getTruth(*xDiTau);
77  return;
78 }
79 
80 
82 // Private Part //
84 
85 //______________________________________________________________________________
86 StatusCode DiTauTruthMatchingTool::findTruthTau(const xAOD::DiTauJet& xDiTau)
87 {
88  // check if decorations were already added to the first passed tau
89  if (!m_bIsTruthMatchedAvailable.isValid()) {
90  static const SG::ConstAccessor<char> accIsTruthMatched("IsTruthMatched");
91  m_bIsTruthMatchedAvailable.set (accIsTruthMatched.isAvailable(xDiTau));
92  }
93  if (*m_bIsTruthMatchedAvailable.ptr())
94  return StatusCode::SUCCESS;
95 
96  if (m_bTruthTauAvailable)
97  return checkTruthMatch(xDiTau, *m_truthTausEvent.m_xTruthTauContainerConst);
98  else
99  return checkTruthMatch(xDiTau, *m_truthTausEvent.m_xTruthTauContainer);
100 }
101 
102 //______________________________________________________________________________
103 StatusCode DiTauTruthMatchingTool::checkTruthMatch (const xAOD::DiTauJet& xDiTau, const xAOD::TruthParticleContainer& xTruthTauContainer) const
104 {
105  std::vector<const xAOD::TruthParticle*> vTruthMatch;
106  std::vector<TruthMatchedParticleType> vTruthMatchedParticleType;
107 
108  xAOD::TruthParticleContainer xRemainingTruthTaus = xTruthTauContainer;
109 
110  static const SG::Decorator<char> decIsTruthMatched("IsTruthMatched");
111  static const SG::Decorator<char> decIsTruthHadronic("IsTruthHadronic");
112  static const SG::Decorator<char> decIsTruthHadMu("IsTruthHadMu");
113  static const SG::Decorator<char> decIsTruthHadEl("IsTruthHadEl");
114  static const SG::ConstAccessor<int> accNSubjets("n_subjets");
115  static const SG::ConstAccessor<char> accIsTruthHadronic("IsTruthHadronic");
116 
117  // set default values for each subjet
118  for (int i = 0; i < accNSubjets(xDiTau); ++i)
119  {
120  const xAOD::TruthParticle* xTruthMatch = nullptr;
121  TruthMatchedParticleType eTruthMatchedParticleType = Unknown;
122 
123  vTruthMatch.push_back(xTruthMatch);
124  vTruthMatchedParticleType.push_back(eTruthMatchedParticleType);
125  }
126 
127  // truthmatching for subjets:
128  for (int i = 0; i < accNSubjets(xDiTau); ++i)
129  {
130  TLorentzVector vSubjetTLV;
131  vSubjetTLV.SetPtEtaPhiE(xDiTau.subjetPt(i),
132  xDiTau.subjetEta(i),
133  xDiTau.subjetPhi(i),
134  xDiTau.subjetE(i));
135  if ( truthMatch(vSubjetTLV,
136  xRemainingTruthTaus,
137  vTruthMatch.at(i),
138  vTruthMatchedParticleType.at(i)).isFailure() )
139  {
140  ATH_MSG_WARNING("There was a failure in matching truth taus with subjet " << i);
141  return StatusCode::FAILURE;
142  }
143  if (vTruthMatch.at(i) && vTruthMatchedParticleType.at(i) == TruthHadronicTau)
144  {
145  xRemainingTruthTaus.erase( std::find(xRemainingTruthTaus.begin(),
146  xRemainingTruthTaus.end(),
147  vTruthMatch.at(i)) );
148  }
149  }
150 
151  bool bTruthMatched = true;
152 
153  // create link to the original TruthParticle
154  std::vector< ElementLink < xAOD::TruthParticleContainer > > vTruthLinks;
155  for (int i = 0; i < accNSubjets(xDiTau); ++i)
156  {
157  const xAOD::TruthParticle* xTruthMatch = vTruthMatch.at(i);
158  TruthMatchedParticleType eTruthMatchedParticleType = vTruthMatchedParticleType.at(i);
159  if (xTruthMatch)
160  {
161  if (eTruthMatchedParticleType == TruthHadronicTau)
162  {
163  ElementLink < xAOD::TruthParticleContainer > lTruthParticleLink(xTruthMatch, xTruthTauContainer);
164  vTruthLinks.push_back(lTruthParticleLink);
165  }
166  else if (eTruthMatchedParticleType == TruthMuon)
167  {
168  ElementLink <xAOD::TruthParticleContainer> lTruthParticleLink(xTruthMatch, *m_truthTausEvent.m_xTruthMuonContainerConst);
169  vTruthLinks.push_back(lTruthParticleLink);
170  }
171  else if (eTruthMatchedParticleType == TruthElectron)
172  {
173  ElementLink <xAOD::TruthParticleContainer> lTruthParticleLink(xTruthMatch, *m_truthTausEvent.m_xTruthElectronContainerConst);
174  vTruthLinks.push_back(lTruthParticleLink);
175  }
176  }
177  else
178  {
180  vTruthLinks.push_back(lTruthParticleLink);
181 
182  // ditau is not truth matched if one of the two leading subjets is not truth matched
183  if (i == 0 || i == 1) bTruthMatched = false;
184  }
185  }
186 
188  decTruthParticleLinks ("truthParticleLinks");
189  decTruthParticleLinks(xDiTau) = vTruthLinks;
190  if (!m_bTruthTauAvailable)
191  {
193  decTruthTaus ("TruthTaus");
194  decTruthTaus(xDiTau) = vTruthLinks;
195  }
196 
198  static const SG::Decorator<unsigned int> decClassifierParticleType("classifierParticleTypeTruthLepton");
199  static const SG::Decorator<unsigned int> decClassifierParticleOrigin("classifierParticleOriginTruthLepton");
200  static const SG::Decorator<ElementLink<xAOD::TruthParticleContainer>> decTruthLeptonLink("truthLeptonLink");
201 
204  static const SG::ConstAccessor<int> accTruthType("truthType");
205  static const SG::ConstAccessor<int> accTruthOrigin("truthOrigin");
206  static const SG::ConstAccessor<ElementLink<xAOD::ElectronContainer>> accElLink("elLink");
207  static const SG::ConstAccessor<ElementLink<xAOD::MuonContainer>> accMuLink("muonLink");
208  if(accElLink.isAvailable(xDiTau) && accMuLink.isAvailable(xDiTau))
209  ATH_MSG_ERROR("Links to reco electron and reco muon available for one ditau candidate.");
210  if(accElLink.isAvailable(xDiTau)){
211  const xAOD::Electron* pElectron = *accElLink(xDiTau);
212  if ((accTruthType.isAvailable(*pElectron) && accTruthOrigin.isAvailable(*pElectron)))
213  {
214  mcTruthType = accTruthType(*pElectron);
215  mcTruthOrigin = accTruthOrigin(*pElectron);
216  }
217  lTruthLeptonLink = checkTruthLepton(pElectron);
218  }
219  if(accMuLink.isAvailable(xDiTau)){
220  const xAOD::Muon* pMuon = *accMuLink(xDiTau);
221  if (accTruthType.isAvailable(*pMuon) && accTruthOrigin.isAvailable(*pMuon))
222  {
223  mcTruthType = accTruthType(*pMuon);
224  mcTruthOrigin = accTruthOrigin(*pMuon);
225  }
226  lTruthLeptonLink = checkTruthLepton(pMuon);
227  }
228 
229  decIsTruthHadEl(xDiTau) = (char)(mcTruthType == MCTruthPartClassifier::ParticleType::IsoElectron && accNSubjets(xDiTau) != 0 && vTruthMatchedParticleType[0] == TruthHadronicTau);
230  decIsTruthHadMu(xDiTau) = (char)(mcTruthType == MCTruthPartClassifier::ParticleType::IsoMuon && accNSubjets(xDiTau) != 0 && vTruthMatchedParticleType[0] == TruthHadronicTau);
231  decClassifierParticleType(xDiTau) = mcTruthType;
232  decClassifierParticleOrigin(xDiTau) = mcTruthOrigin;
233  decTruthLeptonLink(xDiTau) = lTruthLeptonLink;
234 
235  static const SG::Decorator<float> decTruthLeadPt("TruthVisLeadPt");
236  static const SG::Decorator<float> decTruthLeadEta("TruthVisLeadEta");
237  static const SG::Decorator<float> decTruthLeadPhi("TruthVisLeadPhi");
238  static const SG::Decorator<float> decTruthLeadM("TruthVisLeadM");
239  static const SG::Decorator<float> decTruthLeadPdgID("TruthLeadPdgID");
240  static const SG::Decorator<float> decTruthSubleadPt("TruthVisSubleadPt");
241  static const SG::Decorator<float> decTruthSubleadEta("TruthVisSubleadEta");
242  static const SG::Decorator<float> decTruthSubleadPhi("TruthVisSubleadPhi");
243  static const SG::Decorator<float> decTruthSubleadM("TruthVisSubleadM");
244  static const SG::Decorator<float> decTruthSubleadPdgID("TruthSubleadPdgID");
245  static const SG::Decorator<float> decTruthDeltaR("TruthVisDeltaR");
246  static const SG::Decorator<float> decTruthMass("TruthVisMass");
247 
248  // the ditau candidate should have at least 2 subjets to be truth matched
249  if ( accNSubjets(xDiTau) < 2) {
250  decIsTruthMatched(xDiTau) = (char)false;
251  decIsTruthHadronic(xDiTau) = (char)false;
252  decTruthLeadPt(xDiTau) = -1234.;
253  decTruthLeadEta(xDiTau) = -1234.;
254  decTruthLeadPhi(xDiTau) = -1234.;
255  decTruthLeadM(xDiTau) = -1234.;
256  decTruthLeadPdgID(xDiTau) = -1234.;
257  decTruthSubleadPt(xDiTau) = -1234.;
258  decTruthSubleadEta(xDiTau) = -1234.;
259  decTruthSubleadPhi(xDiTau) = -1234.;
260  decTruthSubleadM(xDiTau) = -1234.;
261  decTruthSubleadPdgID(xDiTau) = -1234.;
262  decTruthDeltaR(xDiTau) = -1234.;
263  decTruthMass(xDiTau) = -1234.;
264  return StatusCode::SUCCESS;
265  }
266 
267  decIsTruthMatched(xDiTau) = (char)bTruthMatched;
268  if (bTruthMatched)
269  {
270  // ditau is hadronic if two leading subjets are truth matched with hadronic decay
271  decIsTruthHadronic(xDiTau) = (char)(vTruthMatchedParticleType[0]==TruthHadronicTau
272  && vTruthMatchedParticleType[1]==TruthHadronicTau );
273  }
274  else
275  decIsTruthHadronic(xDiTau) = (char)false;
276 
277  if (accIsTruthHadronic(xDiTau))
278  {
279  TLorentzVector tlvTruthTau1;
280  TLorentzVector tlvTruthTau2;
281  tlvTruthTau1.SetPtEtaPhiM(m_accPtVis(*(*vTruthLinks.at(0))),
282  m_accEtaVis(*(*vTruthLinks.at(0))),
283  m_accPhiVis(*(*vTruthLinks.at(0))),
284  m_accMVis(*(*vTruthLinks.at(0))));
285  tlvTruthTau2.SetPtEtaPhiM(m_accPtVis(*(*vTruthLinks.at(1))),
286  m_accEtaVis(*(*vTruthLinks.at(1))),
287  m_accPhiVis(*(*vTruthLinks.at(1))),
288  m_accMVis(*(*vTruthLinks.at(1))));
289 
290  decTruthLeadPt(xDiTau) = std::max(tlvTruthTau1.Pt(), tlvTruthTau2.Pt());
291  decTruthLeadEta(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.Eta() : tlvTruthTau2.Eta();
292  decTruthLeadPhi(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.Phi() : tlvTruthTau2.Phi();
293  decTruthLeadM(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.M() : tlvTruthTau2.M();
294  decTruthLeadPdgID(xDiTau) = (*vTruthLinks.at(0))->pdgId();
295  decTruthSubleadPt(xDiTau) = std::min(tlvTruthTau1.Pt(), tlvTruthTau2.Pt());
296  decTruthSubleadEta(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.Eta() : tlvTruthTau1.Eta();
297  decTruthSubleadPhi(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.Phi() : tlvTruthTau1.Phi();
298  decTruthSubleadM(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.M() : tlvTruthTau1.M();
299  decTruthSubleadPdgID(xDiTau) = (*vTruthLinks.at(1))->pdgId();
300  decTruthDeltaR(xDiTau) = tlvTruthTau1.DeltaR(tlvTruthTau2);
301  decTruthMass(xDiTau) = (tlvTruthTau1 + tlvTruthTau2).M();
302  }
303  else {
304  // set to a default value
305  decTruthLeadPt(xDiTau) = -1234.;
306  decTruthLeadEta(xDiTau) = -1234.;
307  decTruthLeadPhi(xDiTau) = -1234.;
308  decTruthLeadM(xDiTau) = -1234.;
309  decTruthLeadPdgID(xDiTau) = -1234.;
310  decTruthSubleadPt(xDiTau) = -1234.;
311  decTruthSubleadEta(xDiTau) = -1234.;
312  decTruthSubleadPhi(xDiTau) = -1234.;
313  decTruthSubleadM(xDiTau) = -1234.;
314  decTruthSubleadPdgID(xDiTau) = -1234.;
315  decTruthDeltaR(xDiTau) = -1234.;
316  decTruthMass(xDiTau) = -1234.;
317  }
318 
319  return StatusCode::SUCCESS;
320 }
321 
322 //______________________________________________________________________________
323 ElementLink<xAOD::TruthParticleContainer> DiTauTruthMatchingTool::checkTruthLepton(const xAOD::IParticle* pLepton) const {
325  static const SG::ConstAccessor<ElementLink<xAOD::TruthParticleContainer>> accTruthParticleLink("truthParticleLink");
326  if(!accTruthParticleLink.isAvailable(*pLepton)){
327  return truthParticleLink;
328  }
329  truthParticleLink = accTruthParticleLink(*pLepton);
330  return truthParticleLink;
331 }
332 
333 //______________________________________________________________________________
334 StatusCode DiTauTruthMatchingTool::truthMatch(const TLorentzVector& vSubjetTLV,
335  const xAOD::TruthParticleContainer& xTruthTauContainer,
336  const xAOD::TruthParticle* &xTruthMatch,
337  TruthMatchedParticleType &eTruthMatchedParticleType) const
338 {
339  for (auto xTruthTauIt : xTruthTauContainer)
340  {
341  TLorentzVector vTruthVisTLV;
342  vTruthVisTLV.SetPtEtaPhiM(m_accPtVis(*xTruthTauIt),
343  m_accEtaVis(*xTruthTauIt),
344  m_accPhiVis(*xTruthTauIt),
345  m_accMVis(*xTruthTauIt));
346  if (vSubjetTLV.DeltaR(vTruthVisTLV) <= m_dMaxDeltaR)
347  {
348  static const SG::ConstAccessor<char> accIsHadronicTau("IsHadronicTau");
349  if ((bool)accIsHadronicTau(*xTruthTauIt))
350  eTruthMatchedParticleType = TruthHadronicTau;
351  else
352  continue; // don't let leptonic taus steal truthmatch just by chance
353 
354  xTruthMatch = xTruthTauIt;
355  break;
356  }
357  }
358 
359  if (!xTruthMatch and m_truthTausEvent.m_xTruthMuonContainerConst)
360  {
361  double dPtMax = 0.;
362  for (auto xTruthMuonIt : *m_truthTausEvent.m_xTruthMuonContainerConst)
363  {
364  if (vSubjetTLV.DeltaR(xTruthMuonIt->p4()) <= m_dMaxDeltaR)
365  {
366  if (xTruthMuonIt->pt()<dPtMax)
367  continue;
368  eTruthMatchedParticleType = TruthMuon;
369 
370  xTruthMatch = xTruthMuonIt;
371  dPtMax = xTruthMuonIt->pt();
372  }
373  }
374  }
375 
376  if (!xTruthMatch and m_truthTausEvent.m_xTruthElectronContainerConst)
377  {
378  double dPtMax = 0.;
379  for (auto xTruthElectronIt : *m_truthTausEvent.m_xTruthElectronContainerConst)
380  {
381  if (vSubjetTLV.DeltaR(xTruthElectronIt->p4()) <= m_dMaxDeltaR)
382  {
383  if (xTruthElectronIt->pt()<dPtMax)
384  continue;
385  eTruthMatchedParticleType = TruthElectron;
386  xTruthMatch = xTruthElectronIt;
387  dPtMax = xTruthElectronIt->pt();
388  }
389  }
390  }
391 
392  return StatusCode::SUCCESS;
393 }
TauAnalysisTools
Definition: TruthCollectionMakerTau.h:16
GetLCDefs::Unknown
@ Unknown
Definition: GetLCDefs.h:21
TauAnalysisTools::TruthElectron
@ TruthElectron
Definition: PhysicsAnalysis/TauID/TauAnalysisTools/TauAnalysisTools/Enums.h:100
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
TauAnalysisTools::getTruth
const xAOD::TruthParticle * getTruth(const xAOD::TauJet &xTau)
Definition: PhysicsAnalysis/TauID/TauAnalysisTools/Root/HelperFunctions.cxx:194
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
find
std::string find(const std::string &s)
return a remapped string
Definition: hcg.cxx:135
TauAnalysisTools::TruthHadronicTau
@ TruthHadronicTau
Definition: PhysicsAnalysis/TauID/TauAnalysisTools/TauAnalysisTools/Enums.h:97
max
constexpr double max()
Definition: ap_fixedTest.cxx:33
min
constexpr double min()
Definition: ap_fixedTest.cxx:26
initialize
void initialize()
Definition: run_EoverP.cxx:894
xAOD::char
char
Definition: TrigDecision_v1.cxx:38
xAOD::DiTauJet_v1::subjetPhi
float subjetPhi(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:112
IsoElectron
@ IsoElectron
Definition: TruthClasses.h:11
SG::ConstAccessor< char >
MCTruthClassifier.h
xAOD::IParticle
Class providing the definition of the 4-vector interface.
Definition: Event/xAOD/xAODBase/xAODBase/IParticle.h:41
xAOD::Muon_v1
Class describing a Muon.
Definition: Muon_v1.h:38
xAOD::DiTauJet_v1::subjetE
float subjetE(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:122
NonDefined
@ NonDefined
Definition: TruthClasses.h:52
TauAnalysisTools::BuildTruthTaus
Definition: BuildTruthTaus.h:35
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
SG::Decorator< char >
ElectronContainer.h
lumiFormat.i
int i
Definition: lumiFormat.py:85
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
xAOD::TruthParticle_v1
Class describing a truth particle in the MC record.
Definition: TruthParticle_v1.h:37
IsoMuon
@ IsoMuon
Definition: TruthClasses.h:15
DataVector
Derived DataVector<T>.
Definition: DataVector.h:794
TauAnalysisTools::TruthMuon
@ TruthMuon
Definition: PhysicsAnalysis/TauID/TauAnalysisTools/TauAnalysisTools/Enums.h:99
xAOD::DiTauJet_v1::subjetEta
float subjetEta(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:102
TauAnalysisTools::TruthMatchedParticleType
TruthMatchedParticleType
Definition: PhysicsAnalysis/TauID/TauAnalysisTools/TauAnalysisTools/Enums.h:95
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
xAOD::Electron_v1
Definition: Electron_v1.h:34
MuonContainer.h
DataVector::end
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
DiTauTruthMatchingTool.h
Tau, lepton and jet truth matching for ditau jets.
xAOD::DiTauJet_v1
Definition: DiTauJet_v1.h:31
xAOD::TruthParticle_v1::pt
virtual double pt() const override final
The transverse momentum ( ) of the particle.
Definition: TruthParticle_v1.cxx:166
SG::ConstAccessor::isAvailable
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.
xAOD::DiTauJet_v1::subjetPt
float subjetPt(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:92
ConstAccessor.h
Helper class to provide constant type-safe access to aux data.
DataVector::erase
iterator erase(iterator position)
Remove element at a given position.
Decorator.h
Helper class to provide type-safe access to aux data.
HepMCHelpers.h
DataVector::begin
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.