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