ATLAS Offline Software
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  declareProperty( "MaxDeltaR", m_dMaxDeltaR = 0.2);
37 }
38 
39 //______________________________________________________________________________
40 DiTauTruthMatchingTool::~DiTauTruthMatchingTool( )
41 {
42 
43 }
44 
45 //______________________________________________________________________________
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 
54  {
55  ATH_MSG_FATAL("Failed initializing BuildTruthTaus");
56  return StatusCode::FAILURE;
57  }
58  return StatusCode::SUCCESS;
59 }
60 
61 //______________________________________________________________________________
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  {
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 
189  decTruthParticleLinks ("truthParticleLinks");
190  decTruthParticleLinks(xDiTau) = vTruthLinks;
191  if (!m_bTruthTauAvailable)
192  {
194  decTruthTaus ("TruthTaus");
195  decTruthTaus(xDiTau) = vTruthLinks;
196  }
197 
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 
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  static const SG::Decorator<float> decTruthLeadPt("TruthVisLeadPt");
237  static const SG::Decorator<float> decTruthLeadEta("TruthVisLeadEta");
238  static const SG::Decorator<float> decTruthLeadPhi("TruthVisLeadPhi");
239  static const SG::Decorator<float> decTruthLeadM("TruthVisLeadM");
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> decTruthDeltaR("TruthVisDeltaR");
245  static const SG::Decorator<float> decTruthMass("TruthVisMass");
246 
247  // the ditau candidate should have at least 2 subjets to be truth matched
248  if ( accNSubjets(xDiTau) < 2) {
249  decIsTruthMatched(xDiTau) = (char)false;
250  decIsTruthHadronic(xDiTau) = (char)false;
251  decTruthLeadPt(xDiTau) = -1234.;
252  decTruthLeadEta(xDiTau) = -1234.;
253  decTruthLeadPhi(xDiTau) = -1234.;
254  decTruthLeadM(xDiTau) = -1234.;
255  decTruthSubleadPt(xDiTau) = -1234.;
256  decTruthSubleadEta(xDiTau) = -1234.;
257  decTruthSubleadPhi(xDiTau) = -1234.;
258  decTruthSubleadM(xDiTau) = -1234.;
259  decTruthDeltaR(xDiTau) = -1234.;
260  decTruthMass(xDiTau) = -1234.;
261  return StatusCode::SUCCESS;
262  }
263 
264  decIsTruthMatched(xDiTau) = (char)bTruthMatched;
265  if (bTruthMatched)
266  {
267  // ditau is hadronic if two leading subjets are truth matched with hadronic decay
268  decIsTruthHadronic(xDiTau) = (char)(vTruthMatchedParticleType[0]==TruthHadronicTau
269  && vTruthMatchedParticleType[1]==TruthHadronicTau );
270  }
271  else
272  decIsTruthHadronic(xDiTau) = (char)false;
273 
274  if (accIsTruthHadronic(xDiTau))
275  {
276  TLorentzVector tlvTruthTau1;
277  TLorentzVector tlvTruthTau2;
278  tlvTruthTau1.SetPtEtaPhiM(m_accPtVis(*(*vTruthLinks.at(0))),
279  m_accEtaVis(*(*vTruthLinks.at(0))),
280  m_accPhiVis(*(*vTruthLinks.at(0))),
281  m_accMVis(*(*vTruthLinks.at(0))));
282  tlvTruthTau2.SetPtEtaPhiM(m_accPtVis(*(*vTruthLinks.at(1))),
283  m_accEtaVis(*(*vTruthLinks.at(1))),
284  m_accPhiVis(*(*vTruthLinks.at(1))),
285  m_accMVis(*(*vTruthLinks.at(1))));
286 
287  decTruthLeadPt(xDiTau) = std::max(tlvTruthTau1.Pt(), tlvTruthTau2.Pt());
288  decTruthLeadEta(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.Eta() : tlvTruthTau2.Eta();
289  decTruthLeadPhi(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.Phi() : tlvTruthTau2.Phi();
290  decTruthLeadM(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau1.M() : tlvTruthTau2.M();
291  decTruthSubleadPt(xDiTau) = std::min(tlvTruthTau1.Pt(), tlvTruthTau2.Pt());
292  decTruthSubleadEta(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.Eta() : tlvTruthTau1.Eta();
293  decTruthSubleadPhi(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.Phi() : tlvTruthTau1.Phi();
294  decTruthSubleadM(xDiTau) = (tlvTruthTau1.Pt() > tlvTruthTau2.Pt()) ? tlvTruthTau2.M() : tlvTruthTau1.M();
295  decTruthDeltaR(xDiTau) = tlvTruthTau1.DeltaR(tlvTruthTau2);
296  decTruthMass(xDiTau) = (tlvTruthTau1 + tlvTruthTau2).M();
297  }
298  else {
299  // set to a default value
300  decTruthLeadPt(xDiTau) = -1234.;
301  decTruthLeadEta(xDiTau) = -1234.;
302  decTruthLeadPhi(xDiTau) = -1234.;
303  decTruthLeadM(xDiTau) = -1234.;
304  decTruthSubleadPt(xDiTau) = -1234.;
305  decTruthSubleadEta(xDiTau) = -1234.;
306  decTruthSubleadPhi(xDiTau) = -1234.;
307  decTruthSubleadM(xDiTau) = -1234.;
308  decTruthDeltaR(xDiTau) = -1234.;
309  decTruthMass(xDiTau) = -1234.;
310  }
311 
312  return StatusCode::SUCCESS;
313 }
314 
315 //______________________________________________________________________________
316 ElementLink<xAOD::TruthParticleContainer> DiTauTruthMatchingTool::checkTruthLepton(const xAOD::IParticle* pLepton) const {
318  static const SG::ConstAccessor<ElementLink<xAOD::TruthParticleContainer>> accTruthParticleLink("truthParticleLink");
319  if(!accTruthParticleLink.isAvailable(*pLepton)){
320  return truthParticleLink;
321  }
322  truthParticleLink = accTruthParticleLink(*pLepton);
323  return truthParticleLink;
324 }
325 
326 //______________________________________________________________________________
327 StatusCode DiTauTruthMatchingTool::truthMatch(const TLorentzVector& vSubjetTLV,
328  const xAOD::TruthParticleContainer& xTruthTauContainer,
329  const xAOD::TruthParticle* &xTruthMatch,
330  TruthMatchedParticleType &eTruthMatchedParticleType) const
331 {
332  for (auto xTruthTauIt : xTruthTauContainer)
333  {
334  TLorentzVector vTruthVisTLV;
335  vTruthVisTLV.SetPtEtaPhiM(m_accPtVis(*xTruthTauIt),
336  m_accEtaVis(*xTruthTauIt),
337  m_accPhiVis(*xTruthTauIt),
338  m_accMVis(*xTruthTauIt));
339  if (vSubjetTLV.DeltaR(vTruthVisTLV) <= m_dMaxDeltaR)
340  {
341  static const SG::ConstAccessor<char> accIsHadronicTau("IsHadronicTau");
342  if ((bool)accIsHadronicTau(*xTruthTauIt))
343  eTruthMatchedParticleType = TruthHadronicTau;
344  else
345  continue; // don't let leptonic taus steal truthmatch just by chance
346 
347  xTruthMatch = xTruthTauIt;
348  break;
349  }
350  }
351 
352  if (!xTruthMatch and m_truthTausEvent.m_xTruthMuonContainerConst)
353  {
354  double dPtMax = 0.;
355  for (auto xTruthMuonIt : *m_truthTausEvent.m_xTruthMuonContainerConst)
356  {
357  if (vSubjetTLV.DeltaR(xTruthMuonIt->p4()) <= m_dMaxDeltaR)
358  {
359  if (xTruthMuonIt->pt()<dPtMax)
360  continue;
361  eTruthMatchedParticleType = TruthMuon;
362 
363  xTruthMatch = xTruthMuonIt;
364  dPtMax = xTruthMuonIt->pt();
365  }
366  }
367  }
368 
369  if (!xTruthMatch and m_truthTausEvent.m_xTruthElectronContainerConst)
370  {
371  double dPtMax = 0.;
372  for (auto xTruthElectronIt : *m_truthTausEvent.m_xTruthElectronContainerConst)
373  {
374  if (vSubjetTLV.DeltaR(xTruthElectronIt->p4()) <= m_dMaxDeltaR)
375  {
376  if (xTruthElectronIt->pt()<dPtMax)
377  continue;
378  eTruthMatchedParticleType = TruthElectron;
379  xTruthMatch = xTruthElectronIt;
380  dPtMax = xTruthElectronIt->pt();
381  }
382  }
383  }
384 
385  return StatusCode::SUCCESS;
386 }
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
AthCommonDataStore< AthCommonMsg< Algorithm > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
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: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: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:228
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.