34 , m_metaStore(
"MetaDataStore",
n )
36 declareInterface<DerivationFramework::IAugmentationTool>(
this);
51 ATH_MSG_INFO(
"Using " << m_particlesKey.key() <<
" as the input truth container key");
54 ATH_CHECK(m_outputParticlesKey.initialize());
55 ATH_MSG_INFO(
"New truth particles container key: " << m_outputParticlesKey.key() );
57 if (m_partString.empty()) {
59 return StatusCode::FAILURE;
60 }
else {
ATH_MSG_INFO(
"Truth particle selection string: " << m_partString );}
63 if (!m_partString.empty()) {
64 ATH_CHECK( initializeParser(m_partString) );
68 ATH_CHECK(m_originReadDecorKey.initialize());
69 ATH_CHECK(m_typeReadDecorKey.initialize());
70 ATH_CHECK(m_outcomeReadDecorKey.initialize());
71 ATH_CHECK(m_classificationReadDecorKey.initialize());
74 ATH_CHECK(m_linkDecoratorKey.initialize());
75 ATH_CHECK(m_originDecoratorKey.initialize());
76 ATH_CHECK(m_typeDecoratorKey.initialize());
77 ATH_CHECK(m_outcomeDecoratorKey.initialize());
78 ATH_CHECK(m_classificationDecoratorKey.initialize());
79 ATH_CHECK(m_motherIDDecoratorKey.initialize());
80 ATH_CHECK(m_daughterIDDecoratorKey.initialize());
81 ATH_CHECK(m_hadronOriginDecoratorKey.initialize());
83 return StatusCode::SUCCESS;
90 ATH_MSG_INFO(
"Processed "<< m_ntotpart <<
" truth particles, "<< m_npasspart <<
" were retained ");
92 return StatusCode::SUCCESS;
99 const EventContext& ctx = Gaudi::Hive::currentContext();
103 static const bool is_sherpa = [
this]() {
104 bool is_sherpa =
false;
112 if (m_metaStore->retrieve(truthMetaData).isSuccess() && !truthMetaData->
empty()){
114 const std::string gens = truthMetaData->
at(0)->generators();
115 is_sherpa = (gens.find(
"sherpa")==std::string::npos &&
116 gens.find(
"Sherpa")==std::string::npos &&
117 gens.find(
"SHERPA")==std::string::npos) ?
false :
true;
120 ATH_MSG_WARNING(
"Found xAODTruthMetaDataContainer empty! Configuring to be NOT Sherpa.");
122 ATH_MSG_INFO(
"From metadata configured: Sherpa? " << is_sherpa);
124 ATH_MSG_WARNING(
"Could not find metadata container in storegate; assuming NOT Sherpa");
131 if (!truthParticles.
isValid()) {
132 ATH_MSG_ERROR(
"Couldn't retrieve TruthParticle collection with name " << m_particlesKey);
133 return StatusCode::FAILURE;
138 ATH_CHECK(newParticlesWriteHandle.
record(std::make_unique<xAOD::TruthParticleContainer>(),
139 std::make_unique<xAOD::TruthParticleAuxContainer>()));
140 ATH_MSG_DEBUG(
"Recorded new TruthParticleContainer with key: " << (m_outputParticlesKey.key()));
143 unsigned int nParticles = truthParticles->
size();
144 m_ntotpart += nParticles;
163 if (!m_partString.empty()) {
164 std::vector<int>
entries = m_parser->evaluateAsVector();
168 ATH_MSG_ERROR(
"Sizes incompatible! Are you sure your selection string used TruthParticles?");
169 return StatusCode::FAILURE;
178 bool SherpaW =
false;
179 bool SherpaZ =
false;
180 if (m_partString.value().find(
"24") != std::string::npos && m_do_sherpa && is_sherpa) {
183 if (m_partString.value().find(
"23") != std::string::npos && m_do_sherpa && is_sherpa) {
186 if ((SherpaW or SherpaZ) && is_sherpa){
193 if ((SherpaW || SherpaZ) && is_sherpa){
195 std::vector<const xAOD::TruthParticle*> status20, status3;
196 for (
unsigned int i=0;
i<nParticles; ++
i) {
198 if (!truthParticles->
at(
i))
continue;
202 if (truthParticles->
at(
i)->status()==20) status20.push_back( truthParticles->
at(
i) );
203 if (truthParticles->
at(
i)->status()== 3) status3.push_back( truthParticles->
at(
i) );
207 if (!status20.empty()){
208 status3.swap(status20);
211 if ((status3.size()==2 || status3.size()==4 || status3.size()==6) && (SherpaZ || SherpaW)){
213 int gens[3] = {0,0,0};
214 for (
size_t i=0;
i<status3.size();++
i){
215 if (status3[
i]->absPdgId()<13) gens[0]++;
216 else if (status3[
i]->absPdgId()<15) gens[1]++;
220 if (gens[0]>2 || gens[0]==1 || gens[1]>2 || gens[1]==1 || gens[2]>2 || gens[2]==1){
222 ATH_MSG_VERBOSE(
"Too many leptons of one generation. Cannot make bosons. Give up");
223 return StatusCode::SUCCESS;
225 std::vector<const xAOD::TruthParticle*> boson;
226 for (
size_t i=0;
i<status3.size();++
i){
227 if (status3[
i]->absPdgId()<13) boson.push_back(status3[
i]);
228 else if (gens[0]==0 && status3[
i]->absPdgId()<15) boson.push_back(status3[
i]);
229 else if (gens[0]==0 && gens[1]==0) boson.push_back(status3[
i]);
230 if (boson.size()==2){
234 if (boson[0]->
pdgId()==-boson[1]->
pdgId()) pdg_id=23;
235 else if (std::abs(boson[0]->
pdgId()+boson[1]->
pdgId())!=1){
240 pdg_id=24*(boson[0]->pdgId()+boson[1]->pdgId());
242 if ( (SherpaW && std::abs(pdg_id)==24) ||
243 (SherpaZ && pdg_id==23) ){
246 newParticlesWriteHandle->
push_back( xTruthParticle );
247 if (m_keep_navigation_info){
248 if (boson[0]->hasProdVtx()) {
249 if ((boson[0]->prodVtx()->nIncomingParticles() > 0) && (boson[0]->prodVtx()->incomingParticle(0)!=
nullptr)) {
250 motherIDDecorator(*xTruthParticle) = boson[0]->prodVtx()->incomingParticle(0)->pdgId();
251 }
else {motherIDDecorator(*xTruthParticle) = 0;}
252 }
else {motherIDDecorator(*xTruthParticle) = 0;}
253 if (boson[0]->hasDecayVtx()) {
254 if ((boson[0]->decayVtx()->nOutgoingParticles() > 0) && (boson[0]->decayVtx()->outgoingParticle(0)!=
nullptr)) {
255 daughterIDDecorator(*xTruthParticle) = boson[0]->decayVtx()->outgoingParticle(0)->pdgId();
256 }
else {daughterIDDecorator(*xTruthParticle) = 0;}
257 }
else {daughterIDDecorator(*xTruthParticle) = 0;}
266 xTruthParticle->
setM(new_mom.M());
267 xTruthParticle->
setPx(new_mom.Px());
268 xTruthParticle->
setPy(new_mom.Py());
269 xTruthParticle->
setPz(new_mom.Pz());
270 xTruthParticle->
setE(new_mom.E());
275 if (gens[0]==0 && gens[1]==0) gens[2]=0;
276 else if (gens[0]==0) gens[1]=0;
281 if (status3.size()==1 || status3.size()==3 || status3.size()==5 || status3.size()>6){
282 ATH_MSG_WARNING(status3.size() <<
" leptons found in the Sherpa event record. Not sure how to deal with this.");
284 return StatusCode::SUCCESS;
286 for (
unsigned int i=0;
i<nParticles; ++
i) {
293 bool same_as_mother =
false;
294 bool same_as_daughter =
false;
298 same_as_mother =
true;
305 same_as_daughter =
true;
309 if (same_as_mother && same_as_daughter){
315 newParticlesWriteHandle->
push_back( xTruthParticle );
316 if (m_keep_navigation_info){
320 }
else {motherIDDecorator(*xTruthParticle) = 0;}
321 }
else {motherIDDecorator(*xTruthParticle) = 0;}
325 }
else {daughterIDDecorator(*xTruthParticle) = 0;}
326 }
else {daughterIDDecorator(*xTruthParticle) = 0;}
329 *xTruthParticle=*theParticle;
331 typeDecorator(*xTruthParticle) =
334 originDecorator(*xTruthParticle) =
337 outcomeDecorator(*xTruthParticle) =
340 classificationDecorator(*xTruthParticle) =
341 classificationReadDecor.
withDefault(*theParticle, 0);
343 if (m_outputParticlesKey.key()==
"TruthHFHadrons"){
345 hadronOriginDecorator(*xTruthParticle) =
346 TopHadronOriginFlagAcc.
withDefault (*theParticle, 0);
349 if(m_keep_navigation_info) linkDecorator(*xTruthParticle) = eltp;
354 for (
unsigned int i=0;
i<nParticles; ++
i)
if (
entries[
i]) ++m_npasspart;
357 return StatusCode::SUCCESS;