ATLAS Offline Software
Loading...
Searching...
No Matches
BuildTruthTaus.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2026 CERN for the benefit of the ATLAS collaboration
3*/
4
5// Local include(s)
11
12// Core include(s):
13#include "AthLinks/ElementLink.h"
15
16// EDM include(s):
19
21
22using namespace ADD_CMODE(TauAnalysisTools);
23
24//=================================PUBLIC-PART==================================
25//______________________________________________________________________________
26BuildTruthTaus::BuildTruthTaus( const std::string& name )
27 : AsgMetadataTool(name)
28{
29}
30
31//______________________________________________________________________________
32StatusCode BuildTruthTaus::initialize()
33{
34 if (!m_truthMatchingMode) {
35 ATH_MSG_INFO( "Initializing BuildTruthTaus, will generate " << m_truthTauOutputContainer.key() << " from " << m_truthParticleContainer.key() << " container" );
36 }
37 else {
38 ATH_MSG_INFO( "Initializing BuildTruthTaus in truth matching mode, using input container " << m_truthTauInputContainer.key() );
39 }
40
41 // input containers
42 ATH_CHECK( m_truthTauInputContainer.initialize(m_truthMatchingMode) );
43 ATH_CHECK( m_truthParticleContainer.initialize(SG::AllowEmpty) );
44 ATH_CHECK( m_truthElectronContainer.initialize(SG::AllowEmpty) );
45 ATH_CHECK( m_truthMuonContainer.initialize(SG::AllowEmpty) );
46 ATH_CHECK( m_truthJetContainer.initialize(SG::AllowEmpty) );
47
48 // in truth matching mode, truth electron/muon/jet containers are optional but expected
49 if (m_truthMatchingMode) {
50 if (m_truthElectronContainer.empty()) {
51 ATH_MSG_WARNING("Truth electron container is not available, won't perform matching to truth electrons");
52 }
53 if (m_truthMuonContainer.empty()) {
54 ATH_MSG_WARNING("Truth muon container is not available, won't perform matching to truth muons");
55 }
56 if (m_truthJetContainer.empty()) {
57 ATH_MSG_WARNING("Truth jet container is not available, won't perform matching to truth jets");
58 }
59 }
60
61 // output container
62 ATH_CHECK( m_truthTauOutputContainer.initialize(!m_truthMatchingMode) );
63
64 // input decorations
65 ATH_CHECK( m_originReadDecorKey.initialize(!m_truthMatchingMode) );
66 ATH_CHECK( m_typeReadDecorKey.initialize(!m_truthMatchingMode) );
67 ATH_CHECK( m_outcomeReadDecorKey.initialize(!m_truthMatchingMode) );
68 ATH_CHECK( m_classificationReadDecorKey.initialize(!m_truthMatchingMode) );
69 // output decorations
70 ATH_CHECK( m_linkDecoratorKey.initialize(!m_truthMatchingMode) );
71 ATH_CHECK( m_originDecoratorKey.initialize(!m_truthMatchingMode) );
72 ATH_CHECK( m_typeDecoratorKey.initialize(!m_truthMatchingMode) );
73 ATH_CHECK( m_outcomeDecoratorKey.initialize(!m_truthMatchingMode) );
74 ATH_CHECK( m_classificationDecoratorKey.initialize(!m_truthMatchingMode) );
75
76 // drop at earliest occasion
77 m_bTruthTauAvailable = !m_truthTauInputContainer.empty();
78
79 return StatusCode::SUCCESS;
80}
81
82
83StatusCode BuildTruthTaus::retrieveTruthTaus()
84{
85 const EventContext& ctx = Gaudi::Hive::currentContext();
86 return retrieveTruthTausImpl (m_truthTausEvent, ctx);
87}
88
89
90StatusCode BuildTruthTaus::retrieveTruthTaus(ITruthTausEvent& truthTausEvent, const EventContext& ctx) const
91{
92 return retrieveTruthTausImpl (dynamic_cast<TruthTausEvent&> (truthTausEvent), ctx);
93}
94
95
96StatusCode BuildTruthTaus::retrieveTruthTausImpl(TruthTausEvent& truthTausEvent, const EventContext& ctx) const
97{
98 // truth matching mode
99 if (m_truthMatchingMode) {
100 if (!m_truthElectronContainer.empty()) {
101 SG::ReadHandle<xAOD::TruthParticleContainer> truthElectronsHandle(m_truthElectronContainer, ctx);
102 if (!truthElectronsHandle.isValid()) {
103 ATH_MSG_ERROR ("Could not retrieve " << truthElectronsHandle.key());
104 return StatusCode::FAILURE;
105 }
106 truthTausEvent.m_xTruthElectronContainerConst = truthElectronsHandle.cptr();
107 }
108
109 if (!m_truthMuonContainer.empty()) {
110 SG::ReadHandle<xAOD::TruthParticleContainer> truthMuonsHandle(m_truthMuonContainer, ctx);
111 if (!truthMuonsHandle.isValid()) {
112 ATH_MSG_ERROR ("Could not retrieve " << truthMuonsHandle.key());
113 return StatusCode::FAILURE;
114 }
115 truthTausEvent.m_xTruthMuonContainerConst = truthMuonsHandle.cptr();
116 }
117
118 if (!m_truthJetContainer.empty()) {
119 SG::ReadHandle<xAOD::JetContainer> truthJetsHandle(m_truthJetContainer, ctx);
120 if (!truthJetsHandle.isValid()) {
121 ATH_MSG_ERROR ("Could not retrieve " << truthJetsHandle.key());
122 return StatusCode::FAILURE;
123 }
124 truthTausEvent.m_xTruthJetContainerConst = truthJetsHandle.cptr();
125 }
126
127 SG::ReadHandle<xAOD::TruthParticleContainer> truthTausHandle(m_truthTauInputContainer, ctx);
128 if (!truthTausHandle.isValid()) {
129 ATH_MSG_ERROR ("Could not retrieve " << truthTausHandle.key());
130 return StatusCode::FAILURE;
131 }
132 truthTausEvent.m_xTruthTauContainerConst = truthTausHandle.cptr();
133 }
134 // truth tau building mode
135 else {
136 SG::ReadHandle<xAOD::TruthParticleContainer> truthParticlesHandle(m_truthParticleContainer, ctx);
137 if (!truthParticlesHandle.isValid()) {
138 ATH_MSG_ERROR ("Could not retrieve " << truthParticlesHandle.key());
139 return StatusCode::FAILURE;
140 }
141 truthTausEvent.m_xTruthParticleContainer = truthParticlesHandle.cptr();
142
143 auto truthTausOutput = std::make_unique<xAOD::TruthParticleContainer>();
144 auto truthTausOutputAux = std::make_unique<xAOD::TruthParticleAuxContainer>();
145 truthTausOutput->setStore(truthTausOutputAux.get());
146 truthTausEvent.m_xTruthTauContainer = truthTausOutput.get();
147
148 auto writeHandle = SG::makeHandle(m_truthTauOutputContainer, ctx);
149 ATH_CHECK(writeHandle.record(std::move(truthTausOutput), std::move(truthTausOutputAux)));
150
151 ATH_CHECK( buildTruthTausFromTruthParticles(truthTausEvent, ctx) );
152 }
153
154 return StatusCode::SUCCESS;
155}
156
157//=================================PRIVATE-PART=================================
158//______________________________________________________________________________
159//______________________________________________________________________________
160StatusCode
161BuildTruthTaus::buildTruthTausFromTruthParticles(TruthTausEvent& truthTausEvent, const EventContext& ctx) const
162{
163 SG::ReadDecorHandle<xAOD::TruthParticleContainer, unsigned int > originReadDecor(m_originReadDecorKey, ctx);
164 SG::ReadDecorHandle<xAOD::TruthParticleContainer, unsigned int > typeReadDecor(m_typeReadDecorKey, ctx);
165 SG::ReadDecorHandle<xAOD::TruthParticleContainer, unsigned int > outcomeReadDecor(m_outcomeReadDecorKey, ctx);
166 SG::ReadDecorHandle<xAOD::TruthParticleContainer, unsigned int > classificationReadDecor(m_classificationReadDecorKey, ctx);
167
169 SG::WriteDecorHandle<xAOD::TruthParticleContainer, unsigned int > originDecorator(m_originDecoratorKey, ctx);
170 SG::WriteDecorHandle<xAOD::TruthParticleContainer, unsigned int > typeDecorator(m_typeDecoratorKey, ctx);
171 SG::WriteDecorHandle<xAOD::TruthParticleContainer, unsigned int > outcomeDecorator(m_outcomeDecoratorKey, ctx);
172 SG::WriteDecorHandle<xAOD::TruthParticleContainer, unsigned int > classificationDecorator(m_classificationDecoratorKey, ctx);
173
174 for (auto xTruthParticle : *truthTausEvent.m_xTruthParticleContainer)
175 {
176 if ( xTruthParticle->isTau() )
177 {
178 auto xTruthTau = std::make_unique<xAOD::TruthParticle>();
179 xTruthTau->makePrivateStore( *xTruthParticle );
180
181 if ( examineTruthTau(*xTruthTau).isFailure() )
182 {
183 continue;
184 }
185
186 // adding decorations onto the truth tau requires that it is in the final Aux store
187 truthTausEvent.m_xTruthTauContainer->push_back(std::move(xTruthTau));
188 xAOD::TruthParticle* truthTau = truthTausEvent.m_xTruthTauContainer->back();
189
190 // propagate MCTruthClassifier decorations
191 typeDecorator(*truthTau) = typeReadDecor.isAvailable() ? typeReadDecor(*xTruthParticle) : -1234;
192 originDecorator(*truthTau) = originReadDecor.isAvailable() ? originReadDecor(*xTruthParticle) : -1234;
193 outcomeDecorator(*truthTau) = outcomeReadDecor.isAvailable() ? outcomeReadDecor(*xTruthParticle) : -1234;
194 classificationDecorator(*truthTau) = classificationReadDecor.isAvailable() ? classificationReadDecor(*xTruthParticle) : -1234;
195
196 // create link to the original TruthParticle
197 ElementLink < xAOD::TruthParticleContainer > lTruthParticleLink(xTruthParticle, *truthTausEvent.m_xTruthParticleContainer, ctx);
198 linkDecorator(*truthTau) = lTruthParticleLink;
199 }
200 }
201 return StatusCode::SUCCESS;
202}
203
204//______________________________________________________________________________
205StatusCode BuildTruthTaus::examineTruthTau(const xAOD::TruthParticle& xTruthParticle) const
206{
207 // skip this tau if it has no decay vertex, should not happen
208 if ( !xTruthParticle.hasDecayVtx() )
209 return StatusCode::FAILURE;
210
211 ATH_MSG_VERBOSE("looking for charged daughters of a truth tau");
212
213 TauTruthInfo truthInfo;
214
215 const xAOD::TruthVertex* xDecayVertex = xTruthParticle.decayVtx();
216 if (xDecayVertex == nullptr)
217 return StatusCode::FAILURE;
218 for ( size_t iOutgoingParticle = 0; iOutgoingParticle < xDecayVertex->nOutgoingParticles(); ++iOutgoingParticle )
219 {
220 const xAOD::TruthParticle* xTruthDaughter = xDecayVertex->outgoingParticle(iOutgoingParticle);
221 if (xTruthDaughter == nullptr)
222 {
223 ATH_MSG_ERROR("examineTruthTau: Truth daughter of tau decay was not found in "<< m_truthParticleContainer.key() <<" container. Please ensure that this container has the full tau decay information or produce the TruthTaus container in AtlasDerivation.");
224 return StatusCode::FAILURE;
225 }
226
227 // if tau decays into tau this is not a proper tau decay
228 if ( xTruthDaughter->isTau() )
229 {
230 ATH_MSG_VERBOSE("Tau decays into a tau itself -> skip this decay");
231 return StatusCode::FAILURE;
232 }
233 }
234
235 examineTruthTauDecay(xTruthParticle, truthInfo).ignore();
236
237 if (truthInfo.m_bIsHadronicTau)
238 ATH_MSG_VERBOSE(truthInfo.m_iNChargedDaughters << " prong hadronic truth tau was found with uniqueID "<<HepMC::uniqueID(xTruthParticle));
239 else
240 ATH_MSG_VERBOSE(truthInfo.m_iNChargedDaughters << " prong leptonic truth tau was found with uniqueID "<<HepMC::uniqueID(xTruthParticle));
241 if ( truthInfo.m_iNChargedDaughters%2 == 0 )
242 {
243 ATH_MSG_WARNING("found tau with even multiplicity: " << truthInfo.m_iNChargedDaughters);
244 printDecay(xTruthParticle);
245 }
246
247 static const SG::Decorator<double> decPtVis("pt_vis");
248 static const SG::Decorator<double> decEtaVis("eta_vis");
249 static const SG::Decorator<double> decPhiVis("phi_vis");
250 static const SG::Decorator<double> decMVis("m_vis");
251
252 static const SG::Decorator<size_t> decNumCharged("numCharged");
253 static const SG::Decorator<size_t> decNumChargedPion("numChargedPion");
254 static const SG::Decorator<size_t> decNumNeutral("numNeutral");
255 static const SG::Decorator<size_t> decNumNeutralPion("numNeutralPion");
256
257 decPtVis(xTruthParticle) = truthInfo.m_vTruthVisTLV.Pt();
258 decEtaVis(xTruthParticle) = truthInfo.m_vTruthVisTLV.Eta();
259 decPhiVis(xTruthParticle) = truthInfo.m_vTruthVisTLV.Phi();
260 decMVis(xTruthParticle) = truthInfo.m_vTruthVisTLV.M();
261
262 decNumCharged(xTruthParticle) = truthInfo.m_iNChargedDaughters;
263 decNumChargedPion(xTruthParticle) = truthInfo.m_iNChargedPions;
264 decNumNeutral(xTruthParticle) = truthInfo.m_iNNeutralPions+truthInfo.m_iNNeutralOthers;
265 decNumNeutralPion(xTruthParticle) = truthInfo.m_iNNeutralPions;
266
267 static const SG::Decorator<char> decIsHadronicTau("IsHadronicTau");
268 decIsHadronicTau(xTruthParticle) = (char)truthInfo.m_bIsHadronicTau;
269
270 if ( m_bWriteInvisibleFourMomentum )
271 {
272 TLorentzVector vTruthInvisTLV = xTruthParticle.p4() - truthInfo.m_vTruthVisTLV;
273 static const SG::Decorator<double> decPtInvis("pt_invis");
274 static const SG::Decorator<double> decEtaInvis("eta_invis");
275 static const SG::Decorator<double> decPhiInvis("phi_invis");
276 static const SG::Decorator<double> decMInvis("m_invis");
277 decPtInvis(xTruthParticle) = vTruthInvisTLV.Pt();
278 decEtaInvis(xTruthParticle) = vTruthInvisTLV.Eta();
279 decPhiInvis(xTruthParticle) = vTruthInvisTLV.Phi();
280 decMInvis(xTruthParticle) = vTruthInvisTLV.M();
281 }
282
283 if ( m_bWriteVisibleChargedFourMomentum )
284 {
285 static const SG::Decorator<double> decPtVisCharged("pt_vis_charged");
286 static const SG::Decorator<double> decEtaVisCharged("eta_vis_charged");
287 static const SG::Decorator<double> decPhiVisCharged("phi_vis_charged");
288 static const SG::Decorator<double> decMVisCharged("m_vis_charged");
289 decPtVisCharged(xTruthParticle) = truthInfo.m_vTruthVisTLVCharged.Pt();
290 decEtaVisCharged(xTruthParticle) = truthInfo.m_vTruthVisTLVCharged.Eta();
291 decPhiVisCharged(xTruthParticle) = truthInfo.m_vTruthVisTLVCharged.Phi();
292 decMVisCharged(xTruthParticle) = truthInfo.m_vTruthVisTLVCharged.M();
293 }
294
295 if ( m_bWriteVisibleChargedFourMomentumComponent )
296 {
297 static const SG::Decorator<std::vector<double>> decCharged_pt("pt_charged"); std::vector<double> vec_ptcharged;
298 static const SG::Decorator<std::vector<double>> decCharged_eta("eta_charged"); std::vector<double> vec_etacharged;
299 static const SG::Decorator<std::vector<double>> decCharged_phi("phi_charged"); std::vector<double> vec_phicharged;
300 static const SG::Decorator<std::vector<double>> decCharged_m("m_charged"); std::vector<double> vec_mcharged;
301 for(size_t iProducts=0; iProducts < truthInfo.m_TLVCharged.size(); ++iProducts){
302 TLorentzVector iTLV_c = truthInfo.m_TLVCharged.at(iProducts);
303 vec_ptcharged.push_back(iTLV_c.Pt());
304 vec_etacharged.push_back(iTLV_c.Eta());
305 vec_phicharged.push_back(iTLV_c.Phi());
306 vec_mcharged.push_back(iTLV_c.M());
307 }
308 decCharged_pt(xTruthParticle) = std::move(vec_ptcharged);
309 decCharged_eta(xTruthParticle) = std::move(vec_etacharged);
310 decCharged_phi(xTruthParticle) = std::move(vec_phicharged);
311 decCharged_m(xTruthParticle) = std::move(vec_mcharged);
312 }
313
314 if ( m_bWriteVisibleNeutralFourMomentum )
315 {
316 static const SG::Decorator<double> decPtVisNeutral("pt_vis_neutral");
317 static const SG::Decorator<double> decEtaVisNeutral("eta_vis_neutral");
318 static const SG::Decorator<double> decPhiVisNeutral("phi_vis_neutral");
319 static const SG::Decorator<double> decMVisNeutral("m_vis_neutral");
320 decPtVisNeutral(xTruthParticle) = truthInfo.m_vTruthVisTLVNeutral.Pt();
321 decEtaVisNeutral(xTruthParticle) = truthInfo.m_vTruthVisTLVNeutral.Eta();
322 decPhiVisNeutral(xTruthParticle) = truthInfo.m_vTruthVisTLVNeutral.Phi();
323 decMVisNeutral(xTruthParticle) = truthInfo.m_vTruthVisTLVNeutral.M();
324 }
325
326 if ( m_bWriteVisibleNeutralFourMomentumComponent )
327 {
328 static const SG::Decorator<std::vector<double>> decNeutral_pt("pt_neutral"); std::vector<double> vec_ptneutral;
329 static const SG::Decorator<std::vector<double>> decNeutral_eta("eta_neutral"); std::vector<double> vec_etaneutral;
330 static const SG::Decorator<std::vector<double>> decNeutral_phi("phi_neutral"); std::vector<double> vec_phineutral;
331 static const SG::Decorator<std::vector<double>> decNeutral_m("m_neutral"); std::vector<double> vec_mneutral;
332 for(size_t iProducts=0; iProducts < truthInfo.m_TLVNeutral.size(); ++iProducts){
333 TLorentzVector iTLV_n = truthInfo.m_TLVNeutral.at(iProducts);
334 vec_ptneutral.push_back(iTLV_n.Pt());
335 vec_etaneutral.push_back(iTLV_n.Eta());
336 vec_phineutral.push_back(iTLV_n.Phi());
337 vec_mneutral.push_back(iTLV_n.M());
338 }
339 decNeutral_pt(xTruthParticle) = std::move(vec_ptneutral);
340 decNeutral_eta(xTruthParticle) = std::move(vec_etaneutral);
341 decNeutral_phi(xTruthParticle) = std::move(vec_phineutral);
342 decNeutral_m(xTruthParticle) = std::move(vec_mneutral);
343 }
344
345 if ( m_bWriteDecayModeVector )
346 {
347 static const SG::Decorator<std::vector<int> > decDecayModeVector("DecayModeVector");
348 decDecayModeVector(xTruthParticle) = truthInfo.m_vDecayMode;
349 }
350
351 if ( m_bWriteVertices )
352 {
353 // tau decay vertex
354 static const SG::Decorator<float> decDecayVertexX("decay_vertex_x");
355 static const SG::Decorator<float> decDecayVertexY("decay_vertex_y");
356 static const SG::Decorator<float> decDecayVertexZ("decay_vertex_z");
357
358 decDecayVertexX(xTruthParticle) = truthInfo.m_vDecayVertex.X();
359 decDecayVertexY(xTruthParticle) = truthInfo.m_vDecayVertex.Y();
360 decDecayVertexZ(xTruthParticle) = truthInfo.m_vDecayVertex.Z();
361
362 // tau production vertex
363 static const SG::Decorator<float> decProdVertexX("prod_vertex_x");
364 static const SG::Decorator<float> decProdVertexY("prod_vertex_y");
365 static const SG::Decorator<float> decProdVertexZ("prod_vertex_z");
366
367 decProdVertexX(xTruthParticle) = truthInfo.m_vProdVertex.X();
368 decProdVertexY(xTruthParticle) = truthInfo.m_vProdVertex.Y();
369 decProdVertexZ(xTruthParticle) = truthInfo.m_vProdVertex.Z();
370 }
371
372 return StatusCode::SUCCESS;
373}
374
375//______________________________________________________________________________
376StatusCode BuildTruthTaus::examineTruthTauDecay (const xAOD::TruthParticle& xTruthParticle,
377 TauTruthInfo& truthInfo) const
378{
379 // get vertex and check if it exists
380 const xAOD::TruthVertex* xDecayVertex = xTruthParticle.decayVtx();
381 if (!xDecayVertex)
382 return StatusCode::SUCCESS;
383
384 truthInfo.m_vDecayVertex.SetXYZ(xDecayVertex->x(),xDecayVertex->y(),xDecayVertex->z());
385
386 if (xTruthParticle.hasProdVtx() ) {
387 const xAOD::TruthVertex* xProdVertex = xTruthParticle.prodVtx();
388 truthInfo.m_vProdVertex.SetXYZ(xProdVertex->x(),xProdVertex->y(),xProdVertex->z());
389 } else {
390 truthInfo.m_vProdVertex.SetXYZ(-1234,-1234,-1234);
391 }
392
393 for ( size_t iOutgoingParticle = 0; iOutgoingParticle < xDecayVertex->nOutgoingParticles(); ++iOutgoingParticle )
394 {
395 const xAOD::TruthParticle* xTruthDaughter = xDecayVertex->outgoingParticle(iOutgoingParticle);
396 if (xTruthDaughter == nullptr)
397 {
398 ATH_MSG_ERROR("examineTruthTauDecay: Truth daughter of tau decay was not found in "<< m_truthParticleContainer.key() <<" container. Please ensure that this container has the full tau decay information or produce the TruthTaus container in AtlasDerivation.");
399 return StatusCode::FAILURE;
400 }
401
402 int iAbsPdgId = xTruthDaughter->absPdgId();
403 int iPdgId = xTruthDaughter->pdgId();
404
405 // look at decay of unstable particles
406 if (MC::isDecayed(xTruthDaughter) || !MC::isPhysical(xTruthDaughter))
407 {
408 if ( iAbsPdgId != 111 && iAbsPdgId != 311 && iAbsPdgId != 310 && iAbsPdgId != 130 )
409 {
410 examineTruthTauDecay(*xTruthDaughter, truthInfo).ignore();
411 continue;
412 }
413 }
414
415 // only process stable particles
416 if (!MC::isStable(xTruthDaughter) && !MC::isDecayed(xTruthDaughter))
417 continue;
418
419 // add pdgID to vector for decay mode classification
420 truthInfo.m_vDecayMode.push_back(iPdgId);
421
422 // if tau decays leptonically, indicated by an electron/muon neutrino then
423 // it is not a hadronic decay
424 if ( xTruthDaughter->isHadron() )
425 truthInfo.m_bIsHadronicTau = true;
426
427 // ignore neutrinos for further progress
428 if ( xTruthDaughter->isNeutrino() )
429 {
430 ATH_MSG_VERBOSE("found neutrino decay particle with PdgId "<<iPdgId);
431 continue;
432 }
433
434 // add momentum of non-neutrino particle to visible momentum
435 truthInfo.m_vTruthVisTLV += xTruthDaughter->p4();
436 if ( m_bWriteVisibleChargedFourMomentum )
437 if ( xTruthDaughter->isCharged() )
438 truthInfo.m_vTruthVisTLVCharged += xTruthDaughter->p4();
439 if ( m_bWriteVisibleChargedFourMomentumComponent)
440 if ( xTruthDaughter->isCharged() )
441 truthInfo.m_TLVCharged.push_back(xTruthDaughter->p4());
442 if ( m_bWriteVisibleNeutralFourMomentum )
443 if ( xTruthDaughter->isNeutral() )
444 truthInfo.m_vTruthVisTLVNeutral += xTruthDaughter->p4();
445 if ( m_bWriteVisibleNeutralFourMomentumComponent)
446 if ( xTruthDaughter->isNeutral() )
447 truthInfo.m_TLVNeutral.push_back(xTruthDaughter->p4());
448
449 // only count charged decay particles
450 if ( xTruthDaughter->isCharged() )
451 {
452 ATH_MSG_VERBOSE("found charged decay particle with PdgId "<<iPdgId);
453 truthInfo.m_iNChargedDaughters++;
454 // count charged pions
455 if (iAbsPdgId==211) truthInfo.m_iNChargedPions++;
456 else truthInfo.m_iNChargedOthers++;
457 }
458 else
459 {
460 ATH_MSG_VERBOSE("found neutral decay particle with PdgId "<<iPdgId);
461 // count neutral pions
462 if (iAbsPdgId==111) truthInfo.m_iNNeutralPions++;
463 else truthInfo.m_iNNeutralOthers++;
464 }
465 }
466 return StatusCode::SUCCESS;
467}
468
469void BuildTruthTaus::printDecay(const xAOD::TruthParticle& xTruthParticle, int depth) const
470{
471 // loop over all decay particles, print their kinematic and other properties
472
473 const xAOD::TruthVertex* xDecayVertex = xTruthParticle.decayVtx();
474 if (xDecayVertex == nullptr)
475 return;
476
477 for ( size_t iOutgoingParticle = 0; iOutgoingParticle < xDecayVertex->nOutgoingParticles(); ++iOutgoingParticle )
478 {
479 const xAOD::TruthParticle* xTruthDaughter = xDecayVertex->outgoingParticle(iOutgoingParticle);
480 if (xTruthDaughter == nullptr)
481 {
482 ATH_MSG_WARNING("Truth daughter of tau decay was not found in "<< m_truthParticleContainer.key() <<" container. Please ensure that this container has the full tau decay information or produce the TruthTaus container in AtlasDerivation.");
483 return;
484 }
485 ATH_MSG_WARNING("depth "<<depth
486 <<" e "<<xTruthDaughter->e()
487 <<" eta "<<xTruthDaughter->p4().Eta()
488 <<" phi "<<xTruthDaughter->p4().Phi()
489 <<" pdgid "<<xTruthDaughter->pdgId()
490 <<" status "<<xTruthDaughter->status()
491 <<" uniqueID "<<HepMC::uniqueID(xTruthDaughter));
492 printDecay(*xTruthDaughter, depth+1);
493 }
494}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_ERROR(x,...)
#define ATH_MSG_WARNING(x,...)
#define ATH_MSG_VERBOSE(x,...)
#define ATH_MSG_INFO(x,...)
Handle class for reading a decoration on an object.
Handle class for reading from StoreGate.
Handle class for adding a decoration to an object.
Handle class for recording to StoreGate.
#define ADD_CMODE(name)
Helper class to provide type-safe access to aux data.
Handle class for reading a decoration on an object.
Handle class for adding a decoration to an object.
bool isNeutrino() const
Whether the particle is a neutrino (or antineutrino).
int status() const
Status code.
bool isHadron() const
Whether the particle is a hadron.
int pdgId() const
PDG ID code.
const TruthVertex_v1 * decayVtx() const
The decay vertex of this particle.
bool isCharged() const
Whether the particle is electrically charged.
int absPdgId() const
Absolute PDG ID code (often useful).
bool hasProdVtx() const
Check for a production vertex on this particle.
bool hasDecayVtx() const
Check for a decay vertex on this particle.
virtual double e() const override final
The total energy of the particle.
bool isNeutral() const
Whether the particle is electrically neutral.
const TruthVertex_v1 * prodVtx() const
The production vertex of this particle.
virtual FourMom_t p4() const override final
The full 4-momentum of the particle.
bool isTau() const
Whether the particle is a tau (or antitau).
const TruthParticle_v1 * outgoingParticle(size_t index) const
Get one of the outgoing particles.
float z() const
Vertex longitudinal distance along the beam line form the origin.
float y() const
Vertex y displacement.
size_t nOutgoingParticles() const
Get the number of outgoing particles.
float x() const
Vertex x displacement.
std::string depth
tag string for intendation
Definition fastadd.cxx:46
int uniqueID(const T &p)
bool isStable(const T &p)
Identify if the particle is stable, i.e. has not decayed.
bool isDecayed(const T &p)
Identify if the particle decayed.
bool isPhysical(const T &p)
Identify if the particle is physical, i.e. is stable or decayed.
SG::Decorator< T, ALLOC > Decorator
Helper class to provide type-safe access to aux data, specialized for JaggedVecElt.
Definition AuxElement.h:576
SG::ReadCondHandle< T > makeHandle(const SG::ReadCondHandleKey< T > &key, const EventContext &ctx=Gaudi::Hive::currentContext())
TruthVertex_v1 TruthVertex
Typedef to implementation.
Definition TruthVertex.h:15
TruthParticle_v1 TruthParticle
Typedef to implementation.