ATLAS Offline Software
GenericTruthThinning.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 // GenericTruthThinning.cxx, (c) ATLAS Detector software
8 // Author: James Catmore (James.Catmore@cern.ch)
9 // Removes all truth particles/vertices which do not pass a user-defined cut
10 
16 #include "GaudiKernel/ThreadLocalContext.h"
17 
19 
20 #include <vector>
21 #include <string>
22 
23 // Constructor
25  const std::string& n,
26  const IInterface* p ) :
27 base_class(t,n,p),
28 m_ntotvtx(0),
29 m_ntotpart(0),
30 m_npassvtx(0),
31 m_npasspart(0),
32 m_partString(""),
33 //m_vtxString(""),
34 m_preserveDescendants(false),
35 m_preserveGeneratorDescendants(false),
36 m_preserveAncestors(false),
37 m_tauHandling(true)
38 {
39  declareProperty("ParticleSelectionString", m_partString);
40  //declareProperty("VertexSelectionString", m_vtxString);
41  declareProperty("PreserveDescendants", m_preserveDescendants);
42  declareProperty("PreserveGeneratorDescendants", m_preserveGeneratorDescendants);
43  declareProperty("PreserveAncestors", m_preserveAncestors);
44  declareProperty("TauHandling", m_tauHandling);
45 }
46 
47 // Destructor
49 }
50 
51 // Athena initialize and finalize
53 {
54  ATH_MSG_VERBOSE("initialize() ...");
55  ATH_CHECK( m_particlesKey.initialize (m_streamName) );
56  ATH_CHECK( m_verticesKey.initialize (m_streamName) );
57  ATH_MSG_INFO("Using " << m_particlesKey.key() << " and "<< m_verticesKey.key() << " as the source collections for truth thinning");
58 
59  if (m_partString.empty()/* && m_vtxString==""*/) {
60  ATH_MSG_FATAL("No selection string provided either for vertices or particles!");
61  return StatusCode::FAILURE;
62  } else {ATH_MSG_INFO("Truth thinning selection strings: " << m_partString /*<< " " << m_vtxString*/);}
63 
64  if (m_preserveDescendants && m_preserveGeneratorDescendants) {
65  ATH_MSG_FATAL("You are asking to keep both all descendants, and only those from the event generator. Please check your job options.");
66  return StatusCode::FAILURE;
67  }
68  // Set up the text-parsing machinery for thinning the truth directly according to user cuts
69  if (!m_partString.empty()) {
70  ATH_CHECK( initializeParser(m_partString) );
71  }
72  return StatusCode::SUCCESS;
73 }
74 
76 {
77  ATH_MSG_VERBOSE("finalize() ...");
78  ATH_MSG_INFO("Processed "<< m_ntotvtx <<" truth vertices, "<< m_npassvtx << " were retained ");
79  ATH_MSG_INFO("Processed "<< m_ntotpart <<" truth particles, "<< m_npasspart << " were retained ");
80  ATH_CHECK(finalizeParser());
81  return StatusCode::SUCCESS;
82 }
83 
84 // The thinning itself
86 {
87  const EventContext& ctx = Gaudi::Hive::currentContext();
88 
89  // Retrieve truth collections
91  (m_particlesKey, ctx);
93  (m_verticesKey, ctx);
94 
95  // Set up a mask with the same entries as the full collections
96  unsigned int nParticles = importedTruthParticles->size();
97  unsigned int nVertices = importedTruthVertices->size();
98  std::vector<bool> partMask, vertMask;
99  partMask.assign(nParticles,false); // default: don't keep any truth items
100  vertMask.assign(nVertices,false);
101  m_ntotvtx += nVertices; m_ntotpart += nParticles;
102 
103  // Execute the text parsers and update the mask
104  if (!m_partString.empty()) {
105  std::vector<int> entries = m_parser->evaluateAsVector();
106  unsigned int nEntries = entries.size();
107  // check the sizes are compatible
108  if (nParticles != nEntries ) {
109  ATH_MSG_ERROR("Sizes incompatible! Are you sure your selection string used TruthParticles?");
110  return StatusCode::FAILURE;
111  } else {
112  // set mask
113  for (unsigned int i=0; i<nParticles; ++i) if (entries[i]==1) partMask[i]=true;
114  }
115  }
116 
117  // Special treatment of taus such that only the last one in the chain is kept
118  // Needs another run over the particle collection
119  if (m_tauHandling) {
121  for (unsigned int i=0; i<nParticles; ++i) {
122  const xAOD::TruthParticle* particle = (*importedTruthParticles)[i];
123  if ( MC::isTau(particle) ) { // This is a tau
124  bool last(true);
125  std::vector<int> tauDecayProducts; // all decay products of the tau
126  std::unordered_set<int> tauDecayEncounteredUniqueIDs; // loop checking
127  tauDecayHelper.descendants(particle,tauDecayProducts,tauDecayEncounteredUniqueIDs); // recursive
128  for (unsigned int tauDecIt=0; tauDecIt<tauDecayProducts.size(); ++tauDecIt) {
129  if (std::abs(tauDecayProducts[tauDecIt])==15) { // any taus in the decay products?
130  last = false;
131  break;
132  }
133  }
134  if (!last) partMask[i]=false;
135  } // end of code for tau
136  } // end of loop over particles for tau checking
137  } // end of tau handling option
138 
139  // If user requested preservation of descendants/ancestors:
140  // - loop over the masks and work out which particles need to be descended/ascended from
141  // - do the recursive loop
142  // - update the masks including the descendants/ancestors
143  // To ensure graph completeness, this over-rides anything set by the special treatment
144  // of taus in the section above
146  std::unordered_set<int> encounteredUniqueIDs; // to enable loop handling
147  if (m_preserveDescendants || m_preserveGeneratorDescendants || m_preserveAncestors) {
148  for (unsigned int i=0; i<nParticles; ++i) {
149  bool toKeep = partMask[i];
150  if (!toKeep) continue;
151  const xAOD::TruthParticle* particle = (*importedTruthParticles)[i];
152  encounteredUniqueIDs.clear();
153  if (m_preserveDescendants) decayHelper.descendants(particle,partMask,vertMask,encounteredUniqueIDs,true);
154  encounteredUniqueIDs.clear();
155  if (m_preserveGeneratorDescendants) decayHelper.descendants(particle,partMask,vertMask,encounteredUniqueIDs,false);
156  encounteredUniqueIDs.clear();
157  if (m_preserveAncestors) decayHelper.ancestors(particle,partMask,vertMask,encounteredUniqueIDs);
158  encounteredUniqueIDs.clear();
159  }
160  }
161  //for (unsigned int i=0; i<nVertices; ++i) {
162  // bool toKeep = vertMask[i];
163  // if (!toKeep) continue;
164  // const xAOD::TruthVertex* vertex = (*importedTruthVertices)[i];
165  // decayHelper.descend(vertex,partMask,vertMask);
166  //}
167 
168  // Count the masks
169  m_npasspart += std::count (partMask.begin(), partMask.end(), true);
170  m_npassvtx += std::count (vertMask.begin(), vertMask.end(), true);
171 
172  // Execute the thinning service based on the mask. Finish.
173  importedTruthParticles.keep (partMask);
174  importedTruthVertices.keep (vertMask);
175 
176  return StatusCode::SUCCESS;
177 }
178 
DerivationFramework::GenericTruthThinning::m_tauHandling
bool m_tauHandling
Definition: GenericTruthThinning.h:49
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
Trk::ParticleSwitcher::particle
constexpr ParticleHypothesis particle[PARTICLEHYPOTHESES]
the array of masses
Definition: ParticleHypothesis.h:79
NSWL1::nVertices
int nVertices(const Polygon &p)
Definition: GeoUtils.cxx:35
DerivationFramework::GenericTruthThinning::m_preserveAncestors
bool m_preserveAncestors
Definition: GenericTruthThinning.h:48
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
ThinningHandle.h
Handle for requesting thinning for a data object.
TruthVertexContainer.h
TruthParticleContainer.h
DerivationFramework::GenericTruthThinning::m_partString
std::string m_partString
Definition: GenericTruthThinning.h:44
DerivationFramework::GenericTruthThinning::finalize
virtual StatusCode finalize() override
Definition: GenericTruthThinning.cxx:75
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
XMLtoHeader.count
count
Definition: XMLtoHeader.py:84
SG::ThinningHandle
Handle for requesting thinning for a data object.
Definition: ThinningHandle.h:84
DerivationFramework::DecayGraphHelper::descendants
void descendants(const xAOD::TruthParticle *pHead, std::vector< int > &particleList, std::unordered_set< int > &encounteredUniqueIDs)
Definition: DecayGraphHelper.h:116
SG::ThinningHandleBase::keep
void keep(size_t ndx)
Mark that index ndx in the container should be kept (not thinned away).
Definition: ThinningHandleBase.cxx:75
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:209
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
lumiFormat.i
int i
Definition: lumiFormat.py:85
beamspotman.n
n
Definition: beamspotman.py:729
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
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
GenericTruthThinning.h
DerivationFramework::GenericTruthThinning::m_preserveGeneratorDescendants
bool m_preserveGeneratorDescendants
Definition: GenericTruthThinning.h:47
isTau
bool isTau(const T &p)
Definition: AtlasPID.h:205
DerivationFramework::GenericTruthThinning::initialize
virtual StatusCode initialize() override
Definition: GenericTruthThinning.cxx:52
DerivationFramework::DecayGraphHelper::ancestors
void ancestors(const xAOD::TruthParticle *pHead, std::vector< bool > &particleMask, std::vector< bool > &vertexMask, std::unordered_set< int > &encounteredUniqueIDs)
Definition: DecayGraphHelper.h:186
DerivationFramework::DecayGraphHelper
Definition: DecayGraphHelper.h:26
DerivationFramework::GenericTruthThinning::GenericTruthThinning
GenericTruthThinning(const std::string &t, const std::string &n, const IInterface *p)
Definition: GenericTruthThinning.cxx:24
DerivationFramework::GenericTruthThinning::doThinning
virtual StatusCode doThinning() const override
Definition: GenericTruthThinning.cxx:85
DerivationFramework::GenericTruthThinning::m_preserveDescendants
bool m_preserveDescendants
Definition: GenericTruthThinning.h:46
DerivationFramework::GenericTruthThinning::~GenericTruthThinning
virtual ~GenericTruthThinning()
Definition: GenericTruthThinning.cxx:48
entries
double entries
Definition: listroot.cxx:49
TruthEventContainer.h
dqBeamSpot.nEntries
int nEntries
Definition: dqBeamSpot.py:72
HepMCHelpers.h