ATLAS Offline Software
Loading...
Searching...
No Matches
EvtInclusiveDecay Class Reference

#include <EvtInclusiveDecay.h>

Inheritance diagram for EvtInclusiveDecay:
Collaboration diagram for EvtInclusiveDecay:

Classes

struct  ParticleIdCompare

Public Member Functions

 EvtInclusiveDecay (const std::string &name, ISvcLocator *pSvcLocator)
virtual ~EvtInclusiveDecay ()
StatusCode initialize ()
StatusCode execute (const EventContext &ctx)
 Execute method.
StatusCode finalize ()
std::string xmlpath (void)
virtual StatusCode sysInitialize () override
 Override sysInitialize.
virtual bool isClonable () const override
 Specify if the algorithm is clonable.
virtual StatusCode sysExecute (const EventContext &ctx) override
 Execute an algorithm.
virtual const DataObjIDColl & extraOutputDeps () const override
 Return the list of extra output dependencies.
virtual bool filterPassed (const EventContext &ctx) const
 Get filter decision:
virtual void setFilterPassed (bool state, const EventContext &ctx) const
 Set filter decision:
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const
Event collection accessors (const and non-const)
HepMC::GenEvent *event ATLAS_NOT_CONST_THREAD_SAFE ()
 Access the current signal event (first in the McEventCollection).
McEventCollection *events ATLAS_NOT_CONST_THREAD_SAFE ()
 Access the current event's McEventCollection.
const HepMC::GenEventevent_const (const EventContext &ctx) const
 Access the current signal event (const).
const McEventCollectionevents_const (const EventContext &ctx) const
 Access the current event's McEventCollection (const).

Protected Member Functions

virtual bool isReEntrant () const override final
 Legacy algorithms are not thread-safe.
void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.

Private Types

typedef ServiceHandle< StoreGateSvcStoreGateSvc_t

Private Member Functions

Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

DataObjIDColl m_extendedExtraObjects
 Extra output dependency collection, extended by AthAlgorithmDHUpdate to add symlinks.
StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default).
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default).
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared

Features for derived classes to use internally

ServiceHandle< IAthRNGSvcm_rndmSvc {this, "RndmSvc", "AthRNGSvc"}
IntegerProperty m_dsid {this, "Dsid", 999999}
IntegerProperty m_randomSeed {this, "RandomSeed", 1234567, "Random seed for the built-in random engine"}
 Seed for random number engine.
EvtInclusiveAtRndmGenm_evtAtRndmGen {}
EvtGen * m_myEvtGen {}
std::string m_pdtFile
std::string m_decayFile
std::string m_userDecayFile
std::string m_randomStreamName
std::string m_inputKeyName
std::string m_outputKeyName
bool m_readExisting
bool m_prohibitFinalStateDecay
bool m_prohibitReDecay
bool m_prohibitUnDecay
bool m_prohibitRemoveSelfDecay
std::vector< int > m_blackList
std::set< int > m_blackListSet
bool m_allowAllKnownDecays
bool m_allowDefaultBDecays
std::vector< int > m_whiteList
std::set< int > m_whiteListSet
bool m_printHepMCBeforeEvtGen
bool m_printHepMCAfterEvtGen
bool m_printHepMCHighlighted
bool m_printHepMCHighLightTopLevelDecays
bool m_checkDecayTree
bool m_checkDecayChannels
std::map< int, long > m_noDecayChannels
int m_nRepeatedDecays
int m_maxNRepeatedDecays
bool m_applyUserSelection
bool m_userSelRequireOppositeSignedMu
double m_userSelMu1MinPt
double m_userSelMu2MinPt
double m_userSelMu1MaxEta
double m_userSelMu2MaxEta
double m_userSelMinDimuMass
double m_userSelMaxDimuMass
bool m_isfHerwig
bool m_setVMtransversePol
void reseedRandomEngine (const std::string &streamName, const EventContext &ctx)
CLHEP::HepRandomEngine * getRandomEngine (const std::string &streamName, unsigned long int randomSeedOffset, const EventContext &ctx) const
CLHEP::HepRandomEngine * getRandomEngineDuringInitialize (const std::string &streamName, unsigned long int randomSeedOffset, unsigned int conditionsRun=1, unsigned int lbn=1) const
StatusCode traverseDecayTree (HepMC::GenParticlePtr p, bool isToBeRemoved, std::set< HepMC::GenVertexPtr > &visited, std::set< HepMC::GenParticlePtr, ParticleIdCompare > &toBeDecayed)
void removeDecayTree (HepMC::GenEvent *hepMC, HepMC::GenParticlePtr p)
void decayParticle (HepMC::GenEvent *hepMC, HepMC::GenParticlePtr p)
void addEvtGenDecayTree (HepMC::GenEvent *hepMC, HepMC::GenParticlePtr part, EvtParticle *evtPart, EvtVector4R treeStart, double momentumScaleFactor=1.0)
bool isToBeDecayed (HepMC::ConstGenParticlePtr p, bool doCrossChecks)
bool isDefaultB (const int pId) const
bool passesUserSelection (HepMC::GenEvent *hepMC)
double invMass (HepMC::ConstGenParticlePtr p1, HepMC::ConstGenParticlePtr p2)
void printHepMC (HepMC::GenEvent *hepMC, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)
unsigned int printTree (HepMC::GenParticlePtr p, std::set< HepMC::GenVertexPtr > &visited, int level, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)
std::string pdgName (HepMC::ConstGenParticlePtr p, bool statusHighlighting=false, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)

Properties

SG::ReadHandleKey< McEventCollectionm_mcevents_const { this, "McEventKey", "GEN_EVENT", "StoreGate key of the MC event collection" }
 Const handle to the MC event collection.
std::string m_mcEventKey {}
 StoreGate key for the MC event collection (defaults to GEN_EVENT).
BooleanProperty m_mkMcEvent {this, "MakeMcEvent", false, "Create a new MC event collection if it doesn't exist"}
 Flag to determine if a new MC event collection should be made if it doesn't exist.

Detailed Description

Definition at line 54 of file EvtInclusiveDecay.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ EvtInclusiveDecay()

EvtInclusiveDecay::EvtInclusiveDecay ( const std::string & name,
ISvcLocator * pSvcLocator )

Definition at line 54 of file EvtInclusiveDecay.cxx.

54 :
55 GenBase( name, pSvcLocator ),
57
58 // Basic EvtGen configuration: decay and particle definition files, random number stream
59 declareProperty("pdtFile", m_pdtFile = "inclusive.pdt");
60 declareProperty("decayFile", m_decayFile = "2014inclusive.dec");
61 declareProperty("userDecayFile", m_userDecayFile = "");
62 declareProperty("randomStreamName", m_randomStreamName = "EVTGEN");
63 declareProperty("inputKeyName", m_inputKeyName = "GEN_EVENT");
64 declareProperty("outputKeyName",m_outputKeyName = "GEN_EVENT_EVTGEN");
65 declareProperty("readExisting",m_readExisting=false);
66
67 // Selection of particles to be decayed
68 declareProperty("prohibitFinalStateDecay", m_prohibitFinalStateDecay=false);
69 declareProperty("prohibitReDecay", m_prohibitReDecay=false);
70 declareProperty("prohibitUnDecay", m_prohibitUnDecay=true);
71 declareProperty("prohibitRemoveSelfDecay", m_prohibitRemoveSelfDecay=false);
72 declareProperty("blackList",m_blackList);
73 declareProperty("allowAllKnownDecays", m_allowAllKnownDecays=true);
74 declareProperty("allowDefaultBDecays", m_allowDefaultBDecays=true);
75 declareProperty("whiteList",m_whiteList);
76
77 // Level of output
78 declareProperty("printHepMCBeforeEvtGen", m_printHepMCBeforeEvtGen=false);
79 declareProperty("printHepMCAfterEvtGen", m_printHepMCAfterEvtGen=false);
80 declareProperty("printHepMCHighlighted", m_printHepMCHighlighted=true);
81 declareProperty("printHepMCHighLightTopLevelDecays", m_printHepMCHighLightTopLevelDecays=true);
82
83 // Optional checks
84 declareProperty("checkDecayTree", m_checkDecayTree=false);
85 declareProperty("checkDecayChannels", m_checkDecayChannels=false);
86
87 // Repeated decays
88 declareProperty("maxNRepeatedDecays", m_maxNRepeatedDecays=1);
89
90 // User selection
91 declareProperty("applyUserSelection", m_applyUserSelection=false);
92 declareProperty("userSelRequireOppositeSignedMu", m_userSelRequireOppositeSignedMu=true);
93 declareProperty("userSelMu1MinPt", m_userSelMu1MinPt=0.);
94 declareProperty("userSelMu2MinPt", m_userSelMu2MinPt=0.);
95 declareProperty("userSelMu1MaxEta", m_userSelMu1MaxEta=102.5);
96 declareProperty("userSelMu2MaxEta", m_userSelMu2MaxEta=102.5);
97 declareProperty("userSelMinDimuMass", m_userSelMinDimuMass=0.);
98 declareProperty("userSelMaxDimuMass", m_userSelMaxDimuMass=-1.); // set to negative to not apply cut
99 declareProperty("isfHerwig", m_isfHerwig=false);
100 declareProperty("setVMtransversePol", m_setVMtransversePol=false);
101
102 // We have decided to blacklist Tau decays because we are not sure whether the polarization
103 // would be properly passed to EvtGen
104 m_blackList.push_back(15);
105}
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
std::string m_randomStreamName
std::string m_userDecayFile
std::vector< int > m_blackList
std::string m_outputKeyName
std::vector< int > m_whiteList
GenBase(const std::string &name, ISvcLocator *pSvcLocator)
Constructor.
Definition GenBase.cxx:11

◆ ~EvtInclusiveDecay()

EvtInclusiveDecay::~EvtInclusiveDecay ( )
virtual

Definition at line 109 of file EvtInclusiveDecay.cxx.

109 {
110 delete m_myEvtGen;
111 delete m_evtAtRndmGen;
112}
EvtInclusiveAtRndmGen * m_evtAtRndmGen

Member Function Documentation

◆ addEvtGenDecayTree()

void EvtInclusiveDecay::addEvtGenDecayTree ( HepMC::GenEvent * hepMC,
HepMC::GenParticlePtr part,
EvtParticle * evtPart,
EvtVector4R treeStart,
double momentumScaleFactor = 1.0 )
private

Definition at line 455 of file EvtInclusiveDecay.cxx.

456 {
457 if(evtPart->getNDaug()!=0) {
458 // Add decay vertex, starting from production vertex of particle
459 double ct=(evtPart->getDaug(0)->get4Pos()).get(0)+treeStart.get(0);
460 double x=(evtPart->getDaug(0)->get4Pos()).get(1)+treeStart.get(1);
461 double y=(evtPart->getDaug(0)->get4Pos()).get(2)+treeStart.get(2);
462 double z=(evtPart->getDaug(0)->get4Pos()).get(3)+treeStart.get(3);
463
465
466 hepMC->add_vertex(end_vtx);
467 end_vtx->add_particle_in(std::move(part));
468
469 // Add decay daughter with their own decay trees
470 for(uint it=0; it<evtPart->getNDaug(); it++) {
471 double e=(evtPart->getDaug(it)->getP4Lab()).get(0) * momentumScaleFactor;
472 double px=(evtPart->getDaug(it)->getP4Lab()).get(1) * momentumScaleFactor;
473 double py=(evtPart->getDaug(it)->getP4Lab()).get(2) * momentumScaleFactor;
474 double pz=(evtPart->getDaug(it)->getP4Lab()).get(3) * momentumScaleFactor;
475 int id=EvtPDL::getStdHep(evtPart->getDaug(it)->getId());
476 int status=1;
477 if(evtPart->getDaug(it)->getNDaug() != 0) status=2;
478 HepMC::GenParticlePtr daughter = HepMC::newGenParticlePtr(HepMC::FourVector(px,py,pz,e),id,status);
479 end_vtx->add_particle_out(daughter);
480 addEvtGenDecayTree(hepMC, std::move(daughter), evtPart->getDaug(it), treeStart, momentumScaleFactor);
481 }
482 }
483}
unsigned int uint
#define y
#define x
#define z
void addEvtGenDecayTree(HepMC::GenEvent *hepMC, HepMC::GenParticlePtr part, EvtParticle *evtPart, EvtVector4R treeStart, double momentumScaleFactor=1.0)
HepMC3::FourVector FourVector
GenParticlePtr newGenParticlePtr(const HepMC3::FourVector &mom=HepMC3::FourVector::ZERO_VECTOR(), int pid=0, int status=0)
Definition GenParticle.h:21
HepMC3::GenParticlePtr GenParticlePtr
Definition GenParticle.h:19
GenVertexPtr newGenVertexPtr(const HepMC3::FourVector &pos=HepMC3::FourVector::ZERO_VECTOR(), const int i=0)
Definition GenVertex.h:25
HepMC3::GenVertexPtr GenVertexPtr
Definition GenVertex.h:23
status
Definition merge.py:16

◆ ATLAS_NOT_CONST_THREAD_SAFE() [1/2]

McEventCollection *events GenBase::ATLAS_NOT_CONST_THREAD_SAFE ( )
inherited

Access the current event's McEventCollection.

Note
This function will make a new McEventCollection if there is not already a valid one and MakeMcEvent=True.

◆ ATLAS_NOT_CONST_THREAD_SAFE() [2/2]

HepMC::GenEvent *event GenBase::ATLAS_NOT_CONST_THREAD_SAFE ( )
inlineinherited

Access the current signal event (first in the McEventCollection).

Note
This function will make a new McEventCollection if there is not already a valid one and MakeMcEvent=True.

Definition at line 72 of file GenBase.h.

72 {
73 if (events()->empty())
74 ATH_MSG_ERROR("McEventCollection is empty during first event access");
75 return *(events()->begin());
76 }
#define ATH_MSG_ERROR(x,...)
static const Attributes_t empty

◆ decayParticle()

void EvtInclusiveDecay::decayParticle ( HepMC::GenEvent * hepMC,
HepMC::GenParticlePtr p )
private

Definition at line 417 of file EvtInclusiveDecay.cxx.

417 {
418 ATH_MSG_DEBUG("Decaying particle " << pdgName(part) << " " << part);
419 if (msgLvl(MSG::VERBOSE)) HepMC::Print::line(std::cout,part);
420
421 // Remove existing decay tree, if any, and flag particle as being decayed by EvtGen
422 removeDecayTree(hepMC,part);
423 part->set_status(0);
424
425 // Create EvtGen version of part and have EvtGen decay it.
426 // Since EvtGen uses GeV, convert particles momentum from MeV to GeV.
427 int id = part->pdg_id();
428 EvtId evtId=EvtPDL::evtIdFromStdHep(id);
429 double en =(part->momentum()).e()/1000.;
430 double px=(part->momentum()).px()/1000.;
431 double py=(part->momentum()).py()/1000.;
432 double pz=(part->momentum()).pz()/1000.;
433 EvtVector4R evtP(en,px,py,pz);
434 EvtParticle* evtPart = EvtParticleFactory::particleFactory(evtId,evtP);
435
436 // set transverse polarization to vector mesons (relevant for coherent production of J/Psi etc in UPC)
437 if(m_setVMtransversePol && (id==113 || id== 443 || id==100443 || id==553 || id==100553 || id==200553) )evtPart->setVectorSpinDensity();
438
439 m_myEvtGen->generateDecay(evtPart);
440 if (msgLvl(MSG::VERBOSE)) evtPart->printTree();
441 double ct_s = part->production_vertex()->position().t();
442 double x_s = part->production_vertex()->position().x();
443 double y_s = part->production_vertex()->position().y();
444 double z_s = part->production_vertex()->position().z();
445
446 EvtVector4R treeStart(ct_s,x_s,y_s,z_s);
447 // Add new decay tree to hepMC, converting back from GeV to MeV.
448 addEvtGenDecayTree(hepMC, part, evtPart, treeStart, 1000.);
449 if(evtPart->getNDaug() !=0) part->set_status(2);
450 evtPart->deleteTree();
451}
#define ATH_MSG_DEBUG(x,...)
bool msgLvl(const MSG::Level lvl) const
void removeDecayTree(HepMC::GenEvent *hepMC, HepMC::GenParticlePtr p)
std::string pdgName(HepMC::ConstGenParticlePtr p, bool statusHighlighting=false, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
160 hndl.value(),
161 hndl.documentation());
162
163 }

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
148 }
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

◆ event_const()

const HepMC::GenEvent * GenBase::event_const ( const EventContext & ctx) const
inlineinherited

Access the current signal event (const).

Definition at line 79 of file GenBase.h.

79 {
80 const McEventCollection* coll = events_const(ctx);
81 if (coll->empty())
82 ATH_MSG_ERROR("Const McEventCollection is empty during first event access");
83 return *(coll->begin());
84 }
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
bool empty() const noexcept
Returns true if the collection is empty.
const McEventCollection * events_const(const EventContext &ctx) const
Access the current event's McEventCollection (const).
Definition GenBase.h:93

◆ events_const()

const McEventCollection * GenBase::events_const ( const EventContext & ctx) const
inlineinherited

Access the current event's McEventCollection (const).

Definition at line 93 of file GenBase.h.

93 {
94 SG::ReadHandle<McEventCollection> ret = SG::makeHandle(m_mcevents_const, ctx);
95 if (!ret.isValid())
96 ATH_MSG_ERROR("No McEventCollection found in StoreGate with key " << m_mcevents_const.key());
97 return ret.cptr();
98 }
SG::ReadHandleKey< McEventCollection > m_mcevents_const
Const handle to the MC event collection.
Definition GenBase.h:117
virtual bool isValid() override final
Can the handle be successfully dereferenced?
const_pointer_type cptr()
Dereference the pointer.
SG::ReadCondHandle< T > makeHandle(const SG::ReadCondHandleKey< T > &key, const EventContext &ctx=Gaudi::Hive::currentContext())

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

◆ execute()

StatusCode EvtInclusiveDecay::execute ( const EventContext & ctx)
virtual

Execute method.

New collection that we made.

Collection on which we're operating.

Reimplemented from GenBase.

Definition at line 215 of file EvtInclusiveDecay.cxx.

215 {
216 ATH_MSG_DEBUG("EvtInclusiveDecay executing");
217
219
220 std::string key = m_inputKeyName;
221 // retrieve event from Transient Store (Storegate)
222
223 // Load HepMC info
224 // FIXME should be using Read/WriteHandles here
225
227 std::unique_ptr<McEventCollection> mcEvtCollPtr;
228
230 McEventCollection* mcEvtColl = nullptr;
231
232 if(m_readExisting) {
233 const McEventCollection* oldmcEvtColl = nullptr;
234 CHECK(evtStore()->retrieve(oldmcEvtColl, key));
235 // Fill the new McEventCollection with a copy of the initial HepMC::GenEvent
236 mcEvtCollPtr = std::make_unique<McEventCollection>(*oldmcEvtColl);
237 mcEvtColl = mcEvtCollPtr.get();
238 }
239 else {CHECK(evtStore()->retrieve(mcEvtColl, key));}
240
242 for( mcItr = mcEvtColl->begin(); mcItr != mcEvtColl->end(); ++mcItr ) {
243 HepMC::GenEvent* hepMC = *mcItr;
244
245 // Search HepMC record for particles to be decayed by EvtGen
246 // NOTE: In order to ensure repeatability, we use customized a std::set to obtain
247 // an ordered list of particles to be decayed by EvtGen.
248 std::set<HepMC::GenVertexPtr> visited;
249 std::set<HepMC::GenParticlePtr,ParticleIdCompare> toBeDecayed;
250 for (auto p: *hepMC) {
251 if ( (!p->production_vertex()) ||
252 (p->production_vertex()->particles_in_size() == 0) ) {
253 StatusCode sc = traverseDecayTree(std::move(p),false,visited,toBeDecayed);
254 if (sc.isFailure())
255 return StatusCode::FAILURE;
256 }
257 }
258 // Print HepMC in tree format if desired (before doing anything)
260 msg(MSG::INFO) << "Printing HepMC record at " << hepMC << " BEFORE running EvtGen:" << endmsg;
262 printHepMC(hepMC,&toBeDecayed);
263 else
264 printHepMC(hepMC);
265 }
266
267
268 // Decay selected particles
269 bool eventPassesCuts(false);
270 int loopCounter(0);
271 while( !eventPassesCuts && loopCounter < m_maxNRepeatedDecays ) {
272
273 for (auto p: toBeDecayed) {
274 if (p == 0) {
275 msg(MSG::ERROR ) << "Overlapping decay tree for particle" << p <<endmsg;
276 return StatusCode::FAILURE;
277 }
278 decayParticle(hepMC,std::move(p));
280 }
281
283 eventPassesCuts = passesUserSelection(hepMC);
284 else
285 eventPassesCuts = true;
286
288 loopCounter++;
289 }
290
291 // Store the number of decay attempts in event weights std::map, only if repeated decays enabled
292
293 if(m_maxNRepeatedDecays > 1) {
294 hepMC->weight("nEvtGenDecayAttempts") = loopCounter;
295 }
296 // Print HepMC in tree format if desired (after finishing all EvtGen decays)
298 msg(MSG::INFO) << "Printing HepMC record at " << hepMC << " AFTER running EvtGen:" << endmsg;
300 printHepMC(hepMC,&toBeDecayed);
301 else
302 printHepMC(hepMC);
303 }
304 }
305
306 if (mcEvtCollPtr) {
307 if(m_outputKeyName!=key) {
308 ATH_CHECK(evtStore()->record(std::move(mcEvtCollPtr),m_outputKeyName));
309 }
310 else {
311 ATH_CHECK(evtStore()->overwrite(std::move(mcEvtCollPtr), m_outputKeyName, true));
312 }
313 }
314
315 return StatusCode::SUCCESS;
316}
#define endmsg
#define ATH_CHECK
Evaluate an expression and check for errors.
#define CHECK(...)
Evaluate an expression and check for errors.
static Double_t sc
ServiceHandle< StoreGateSvc > & evtStore()
DataModel_detail::iterator< DataVector > iterator
Definition DataVector.h:842
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
void printHepMC(HepMC::GenEvent *hepMC, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)
void reseedRandomEngine(const std::string &streamName, const EventContext &ctx)
StatusCode traverseDecayTree(HepMC::GenParticlePtr p, bool isToBeRemoved, std::set< HepMC::GenVertexPtr > &visited, std::set< HepMC::GenParticlePtr, ParticleIdCompare > &toBeDecayed)
bool passesUserSelection(HepMC::GenEvent *hepMC)
void decayParticle(HepMC::GenEvent *hepMC, HepMC::GenParticlePtr p)
::StatusCode StatusCode
StatusCode definition for legacy code.
void fillBarcodesAttribute(GenEvent *e)
Definition GenEvent.h:393
HepMC3::GenEvent GenEvent
Definition GenEvent.h:39

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ extraOutputDeps()

const DataObjIDColl & AthCommonAlgorithm< Gaudi::Algorithm >::extraOutputDeps ( ) const
overridevirtualinherited

Return the list of extra output dependencies.

This list is extended to include symlinks implied by inheritance relations.

Definition at line 89 of file AthCommonAlgorithm.cxx.

54{
55 // If we didn't find any symlinks to add, just return the collection
56 // from the base class. Otherwise, return the extended collection.
57 if (!m_extendedExtraObjects.empty()) {
59 }
61}
Common base class for algorithms.

◆ filterPassed()

virtual bool AthCommonAlgorithm< Gaudi::Algorithm >::filterPassed ( const EventContext & ctx) const
inlinevirtualinherited

Get filter decision:

Definition at line 93 of file AthCommonAlgorithm.h.

93 {
94 return execState( ctx ).filterPassed();
95 }
virtual bool filterPassed(const EventContext &ctx) const
Get filter decision:

◆ finalize()

StatusCode EvtInclusiveDecay::finalize ( )

Definition at line 319 of file EvtInclusiveDecay.cxx.

319 {
320
322 ATH_MSG_INFO("The following particles were checked and didn't have any decay channels:");
323 if (msgLvl(MSG::INFO)) {
324 std::cout << std::endl;
325 std::cout << " Particle code # Occurences" << std::endl;
326 std::cout << "----------------------------------" << std::endl;
327 for (std::map<int,long>::iterator p = m_noDecayChannels.begin(); p!=m_noDecayChannels.end(); ++p) {
328 int id = p->first;
329 int count = p->second;
330 std::cout << std::setw(14) << id
331 << std::setw(20) << count
332 << std::endl;
333 }
334 std::cout << std::endl;
335 }
336 }
337 ATH_MSG_INFO("Total number of repeated decays: " << m_nRepeatedDecays);
338 ATH_MSG_INFO("EvtInclusiveDecay finalized");
339 return StatusCode::SUCCESS;
340}
#define ATH_MSG_INFO(x,...)
std::map< int, long > m_noDecayChannels
int count(std::string s, const std::string &regx)
count how many occurances of a regx are in a string
Definition hcg.cxx:148

◆ getRandomEngine()

CLHEP::HepRandomEngine * EvtInclusiveDecay::getRandomEngine ( const std::string & streamName,
unsigned long int randomSeedOffset,
const EventContext & ctx ) const
private

Definition at line 189 of file EvtInclusiveDecay.cxx.

191{
192 ATHRNG::RNGWrapper* rngWrapper = m_rndmSvc->getEngine(this, streamName);
193 rngWrapper->setSeed( streamName, ctx.slot(), randomSeedOffset, ctx.eventID().run_number() );
194 return rngWrapper->getEngine(ctx);
195}
void setSeed(const std::string &algName, const EventContext &ctx)
Set the random seed using a string (e.g.
Definition RNGWrapper.h:154
CLHEP::HepRandomEngine * getEngine(const EventContext &ctx) const
Retrieve the random engine corresponding to the provided EventContext.
Definition RNGWrapper.h:108
ServiceHandle< IAthRNGSvc > m_rndmSvc

◆ getRandomEngineDuringInitialize()

CLHEP::HepRandomEngine * EvtInclusiveDecay::getRandomEngineDuringInitialize ( const std::string & streamName,
unsigned long int randomSeedOffset,
unsigned int conditionsRun = 1,
unsigned int lbn = 1 ) const
private

Definition at line 198 of file EvtInclusiveDecay.cxx.

199{
200 const size_t slot=0;
201 EventContext ctx;
202 ctx.setSlot( slot );
203 ctx.setEventID (EventIDBase (conditionsRun,
204 EventIDBase::UNDEFEVT, // event
205 EventIDBase::UNDEFNUM, // timestamp
206 EventIDBase::UNDEFNUM, // timestamp ns
207 lbn));
209 Atlas::ExtendedEventContext( evtStore()->hiveProxyDict(),
210 conditionsRun) );
211 return getRandomEngine(streamName, randomSeedOffset, ctx);
212}
CLHEP::HepRandomEngine * getRandomEngine(const std::string &streamName, unsigned long int randomSeedOffset, const EventContext &ctx) const
void setExtendedEventContext(EventContext &ctx, ExtendedEventContext &&ectx)
Move an extended context into a context object.

◆ initialize()

StatusCode EvtInclusiveDecay::initialize ( )
virtual

Reimplemented from GenBase.

Definition at line 116 of file EvtInclusiveDecay.cxx.

116 {
117
119 // Get the random number service
120 CHECK(m_rndmSvc.retrieve());
121
122 msg(MSG::INFO) << "EvtInclusiveDecay initialize" << endmsg;
123 msg(MSG::INFO) << "Particle properties definition file = " << m_pdtFile << endmsg;
124 msg(MSG::INFO) << "Main decay file = " << m_decayFile << endmsg;
125 msg(MSG::INFO) << "User decay file = " << m_userDecayFile << endmsg;
126 msg(MSG::INFO) << "Max number of repeated decays = " << m_maxNRepeatedDecays << endmsg;
127 msg(MSG::INFO) << "EvtInclusiveDecay selection parameters:" << endmsg;
128 msg(MSG::INFO) << "* prohibitFinalStateDecay = " << m_prohibitFinalStateDecay << endmsg;
129 msg(MSG::INFO) << "* prohibitReDecay = " << m_prohibitReDecay << endmsg;
130 msg(MSG::INFO) << "* prohibitUnDecay = " << m_prohibitUnDecay << endmsg;
131 msg(MSG::INFO) << "* prohibitRemoveSelfDecay = " << m_prohibitRemoveSelfDecay << endmsg;
132 msg(MSG::INFO) << "* allowAllKnownDecays = " << m_allowAllKnownDecays << endmsg;
133 msg(MSG::INFO) << "* allowDefaultBDecays = " << m_allowDefaultBDecays << endmsg;
134 msg(MSG::INFO) << "User selection parameters:" << endmsg;
135 msg(MSG::INFO) << "* applyUserSelection = " << m_applyUserSelection << endmsg;
136 msg(MSG::INFO) << "* userSelRequireOppositeSignedMu = " << m_userSelRequireOppositeSignedMu << endmsg;
137 msg(MSG::INFO) << "* userSelMu1MinPt = " << m_userSelMu1MinPt << endmsg;
138 msg(MSG::INFO) << "* userSelMu2MinPt = " << m_userSelMu2MinPt << endmsg;
139 msg(MSG::INFO) << "* userSelMu1MaxEta = " << m_userSelMu1MaxEta << endmsg;
140 msg(MSG::INFO) << "* userSelMu2MaxEta = " << m_userSelMu2MaxEta << endmsg;
141 msg(MSG::INFO) << "* userSelMinDimuMass = " << m_userSelMinDimuMass << endmsg;
142 msg(MSG::INFO) << "* userSelMaxDimuMass = " << m_userSelMaxDimuMass << endmsg;
143
144 // Initialize and print blackList
145 m_blackListSet.insert(m_blackList.begin(),m_blackList.end());
146 msg(MSG::INFO) << "* blackList; = ";
147 for (std::set<int>::iterator i = m_blackListSet.begin(); i!=m_blackListSet.end(); ++i)
148 msg(MSG::INFO) << (*i) << " ";
149 msg(MSG::INFO)<< endmsg;
150
151 // Initialize and print whiteList
152 m_whiteListSet.insert(m_whiteList.begin(),m_whiteList.end());
153 msg(MSG::INFO) << "* whiteList = ";
154 for (std::set<int>::iterator i = m_whiteListSet.begin(); i!=m_whiteListSet.end(); ++i)
155 msg(MSG::INFO) << (*i) << " ";
156 msg(MSG::INFO) << endmsg;
157
158 CLHEP::HepRandomEngine* rndmEngine = getRandomEngineDuringInitialize(m_randomStreamName, m_randomSeed, m_dsid);
159 // Obtain random number generator for EvtGen
160 m_evtAtRndmGen = new EvtInclusiveAtRndmGen(rndmEngine);
161
162 // Create an instance of EvtGen and read particle properties and decay files
163 EvtExternalGenList genList(true,xmlpath(),"gamma");
164 EvtAbsRadCorr* radCorrEngine = genList.getPhotosModel();
165 std::list<EvtDecayBase*> extraModels = genList.getListOfModels();
166
167 // Create the EvtGen generator object
168 // EvtGen myGenerator("decayFile.dec", "evt.pdl", randomEnginePointer,
169 // radCorrEngine, &extraModels);
170
171
172 m_myEvtGen = new EvtGen( m_decayFile.c_str(), m_pdtFile.c_str(), m_evtAtRndmGen, radCorrEngine, &extraModels);
173 if(!m_userDecayFile.empty())
174 m_myEvtGen->readUDecay(m_userDecayFile.c_str());
175
176 return StatusCode::SUCCESS;
177}
std::set< int > m_blackListSet
CLHEP::HepRandomEngine * getRandomEngineDuringInitialize(const std::string &streamName, unsigned long int randomSeedOffset, unsigned int conditionsRun=1, unsigned int lbn=1) const
IntegerProperty m_dsid
std::string xmlpath(void)
std::set< int > m_whiteListSet
IntegerProperty m_randomSeed
Seed for random number engine.
virtual StatusCode initialize() override
Definition GenBase.cxx:17

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ invMass()

double EvtInclusiveDecay::invMass ( HepMC::ConstGenParticlePtr p1,
HepMC::ConstGenParticlePtr p2 )
private

Definition at line 614 of file EvtInclusiveDecay.cxx.

614 {
615 double p1Px = p1->momentum().px();
616 double p1Py = p1->momentum().py();
617 double p1Pz = p1->momentum().pz();
618 double p1E = p1->momentum().e();
619 double p2Px = p2->momentum().px();
620 double p2Py = p2->momentum().py();
621 double p2Pz = p2->momentum().pz();
622 double p2E = p2->momentum().e();
623 double dimuE = p2E + p1E;
624 double dimuPx = p2Px + p1Px;
625 double dimuPy = p2Py + p1Py;
626 double dimuPz = p2Pz + p1Pz;
627 double invMass = std::sqrt(dimuE*dimuE - dimuPx*dimuPx - dimuPy*dimuPy - dimuPz*dimuPz);
628
629 return invMass;
630}
double invMass(HepMC::ConstGenParticlePtr p1, HepMC::ConstGenParticlePtr p2)

◆ isClonable()

virtual bool AthCommonAlgorithm< Gaudi::Algorithm >::isClonable ( ) const
inlineoverridevirtualinherited

Specify if the algorithm is clonable.

Only relevant for non-reentrant algorithms. Actual number of clones needs to be set via the "Cardinality" property.

Reimplemented in AFP_DigiTop, AlgB, AlgT, BCM_Digitization, CscDigitBuilder, CscDigitToCscRDO, G4AtlasAlg, G4RunAlg, HGTD_Digitization, HiveAlgBase, InDet::GNNSeedingTrackMaker, InDet::SCT_Clusterization, InDet::SiSPGNNTrackMaker, InDet::SiSPSeededTrackFinder, InDet::SiTrackerSpacePointFinder, ISF::SimKernelMT, ITk::StripDigitization, ITkPixelCablingAlg, ITkStripCablingAlg, LArHitEMapMaker, LArTTL1Maker, LUCID_DigiTop, LVL1::L1TopoSimulation, MergeCalibHits, MergeGenericMuonSimHitColl, MergeHijingPars, MergeMcEventCollection, MergeTrackRecordCollection, MergeTruthJets, MergeTruthParticles, MuonDigitizer, PileUpMTAlg, PixelDigitization, RoIBResultToxAOD, SCT_ByteStreamErrorsTestAlg, SCT_CablingCondAlgFromCoraCool, SCT_CablingCondAlgFromText, SCT_ConditionsParameterTestAlg, SCT_ConditionsSummaryTestAlg, SCT_ConfigurationConditionsTestAlg, SCT_Digitization, SCT_FlaggedConditionTestAlg, SCT_LinkMaskingTestAlg, SCT_MajorityConditionsTestAlg, SCT_ModuleVetoTestAlg, SCT_MonitorConditionsTestAlg, SCT_PrepDataToxAOD, SCT_RawDataToxAOD, SCT_ReadCalibChipDataTestAlg, SCT_ReadCalibDataTestAlg, SCT_RODVetoTestAlg, SCT_SensorsTestAlg, SCT_SiliconConditionsTestAlg, SCT_StripVetoTestAlg, SCT_TdaqEnabledTestAlg, SCT_TestCablingAlg, SCTEventFlagWriter, SCTRawDataProvider, SCTSiLorentzAngleTestAlg, SCTSiPropertiesTestAlg, SGInputLoader, Simulation::BeamEffectsAlg, TileHitVecToCnt, TileMuonFitter, TilePulseForTileMuonReceiver, TileRawChannelMaker, TRTDigitization, and ZDC_DigiTop.

Definition at line 68 of file AthCommonAlgorithm.h.

68 {
69 return true;
70 }

◆ isDefaultB()

bool EvtInclusiveDecay::isDefaultB ( const int pId) const
private

Definition at line 562 of file EvtInclusiveDecay.cxx.

562 {
563 int id = std::abs(pId);
564 if ( id == 511 ||
565 id == 521 ||
566 id == 531 ||
567 id == 541 ||
568 id == 5122 ||
569 id == 5132 ||
570 id == 5232 ||
571 id == 5112 ||
572 id == 5212 ||
573 id == 5222 )
574 return true;
575 else
576 return false;
577}

◆ isReEntrant()

virtual bool AthAlgorithm::isReEntrant ( ) const
inlinefinaloverrideprotectedvirtualinherited

Legacy algorithms are not thread-safe.

Definition at line 47 of file AthAlgorithm.h.

47{ return false; }

◆ isToBeDecayed()

bool EvtInclusiveDecay::isToBeDecayed ( HepMC::ConstGenParticlePtr p,
bool doCrossChecks )
private

Definition at line 493 of file EvtInclusiveDecay.cxx.

493 {
494 int id = p->pdg_id();
495 int nDaughters = 0;
496 auto v = p->end_vertex();
497 if (v) nDaughters = v->particles_out_size();
498
499 // Ignore documentation lines
500 if (p->status() == 3) return false;
501 // And any particles that aren't stable or decayed
502 if(!m_isfHerwig && !MC::isPhysical(p)) return false;
503
504 // Particularly for Herwig, try to ignore particles that really should
505 // be flagged as documentation lines
506 double m2 = p->momentum().m2();
507 if (m2 < -1.0E-3) {
508 ATH_MSG_DEBUG("Ignoring particle " << pdgName(std::move(p)) << " with m^2 = " << m2);
509 return false;
510 }
511
512 // Check whether EvtGen has any decay channels defined for this particle
513 EvtId evtId = EvtPDL::evtIdFromStdHep(id);
514 // std::cout << "EVTID: " << evtId.getId() << " alias " << evtId.getAlias() << std::endl;
515 int nModes = 0;
516 if (evtId.getId()>=0)
517 // nModes = EvtDecayTable::getNMode(evtId.getAlias());
518 nModes = EvtDecayTable::getInstance()->getNMode(evtId.getAlias());
519 if (doCrossChecks) {
520 ATH_MSG_VERBOSE("Checking particle " << pdgName(p)
521 << " (status = " << p->status()
522 <<") -- " << nModes << " decay modes found");
523 if (m_checkDecayChannels && nModes==0) {
524 std::map<int,long>::iterator pos = m_noDecayChannels.find(id);
525 if (pos != m_noDecayChannels.end())
526 (pos->second)++;
527 else
529 }
530 }
531
532 // Check prohibit* settings
533 if (m_prohibitFinalStateDecay && MC::isStable(p)) return false;
534 if (m_prohibitReDecay && nDaughters>0) return false;
535 if (m_prohibitUnDecay && nModes==0) return false;
536 if (m_prohibitRemoveSelfDecay && nDaughters>0) {
537 // For now, check only children - this should be sufficient and checking all
538 // descendants would be very expensive.
539 for (auto itd: *v) {
540 if (std::abs(itd->pdg_id()) == std::abs(id)) return false;
541 }
542 }
543
544 // Check blackList
545 if (m_blackListSet.count(std::abs(id))>0) return false;
546
547 // Check allow* settings
548 if (m_allowAllKnownDecays && nModes>0) return true;
549 if (m_allowDefaultBDecays && isDefaultB(id)) return true;
550
551 // Check whiteList
552 if (m_whiteListSet.count(std::abs(id))>0) return true;
553
554 return false; // Default is NOT to decay through EvtGen
555}
#define ATH_MSG_VERBOSE(x,...)
bool isDefaultB(const int pId) const
bool isStable(const T &p)
Identify if the particle is stable, i.e. has not decayed.
bool isPhysical(const T &p)
Identify if the particle is physical, i.e. is stable or decayed.

◆ msg()

MsgStream & AthCommonMsg< Gaudi::Algorithm >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msgLvl()

bool AthCommonMsg< Gaudi::Algorithm >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ passesUserSelection()

bool EvtInclusiveDecay::passesUserSelection ( HepMC::GenEvent * hepMC)
private

Definition at line 584 of file EvtInclusiveDecay.cxx.

584 {
585 bool passed(false);
586 std::vector<HepMC::GenParticlePtr> *muons = new std::vector<HepMC::GenParticlePtr>;
587
588 for ( const auto& p: *hepMC) {
589 if( std::abs(p->pdg_id()) == 13 )
590 muons->push_back(p);
591 }
592
593 for (auto muItr1 = muons->begin(); muItr1 != muons->end(); ++muItr1) {
594 for (auto muItr2 = muItr1+1; muItr2 != muons->end(); ++muItr2) {
595 if( m_userSelRequireOppositeSignedMu && (*muItr1)->pdg_id() * (*muItr2)->pdg_id() > 0)
596 continue;
597 if( !( (*muItr1)->momentum().perp() > m_userSelMu1MinPt && std::abs((*muItr1)->momentum().pseudoRapidity()) < m_userSelMu1MaxEta &&
598 (*muItr2)->momentum().perp() > m_userSelMu2MinPt && std::abs((*muItr2)->momentum().pseudoRapidity()) < m_userSelMu2MaxEta ) &&
599 !( (*muItr2)->momentum().perp() > m_userSelMu1MinPt && std::abs((*muItr2)->momentum().pseudoRapidity()) < m_userSelMu1MaxEta &&
600 (*muItr1)->momentum().perp() > m_userSelMu2MinPt && std::abs((*muItr1)->momentum().pseudoRapidity()) < m_userSelMu2MaxEta ) )
601 continue;
602 double dimuMass = invMass((*muItr1),(*muItr2));
603 if( !( dimuMass > m_userSelMinDimuMass && (dimuMass < m_userSelMaxDimuMass || m_userSelMaxDimuMass < 0.) ) )
604 continue;
605 passed = true;
606 }
607 }
608
609 delete muons;
610
611 return passed;
612}
bool passed(DecisionID id, const DecisionIDContainer &)
checks if required decision ID is in the set of IDs in the container

◆ pdgName()

std::string EvtInclusiveDecay::pdgName ( HepMC::ConstGenParticlePtr p,
bool statusHighlighting = false,
std::set< HepMC::GenParticlePtr, ParticleIdCompare > * particleSet = nullptr )
private

Definition at line 685 of file EvtInclusiveDecay.cxx.

685 {
686 std::ostringstream buf;
687 bool inlist = false;
688 if (particleSet) for (const auto& pinl: *particleSet) if (pinl&&p) if (pinl.get() == p.get()) inlist=true;
689 if (statusHighlighting) {
690 if ( ((particleSet!=0) && (inlist)) ||
691 ((particleSet==0) && isToBeDecayed(p,false)) )
692 buf << "\033[7m"; // reverse
693 if (p and !MC::isStable(p)) {
694 if (MC::isDecayed(p))
695 buf << "\033[33m"; // yellow
696 else
697 buf << "\033[31m"; // red
698 }
699 }
700 if (p){
701 buf << p->pdg_id();
702 if (statusHighlighting) {
703 buf << "\033[0m"; // revert color attributes
704 }
705 }
706 return buf.str();
707}
bool isToBeDecayed(HepMC::ConstGenParticlePtr p, bool doCrossChecks)
bool isDecayed(const T &p)
Identify if the particle decayed.

◆ printHepMC()

void EvtInclusiveDecay::printHepMC ( HepMC::GenEvent * hepMC,
std::set< HepMC::GenParticlePtr, ParticleIdCompare > * particleSet = nullptr )
private

Definition at line 638 of file EvtInclusiveDecay.cxx.

638 {
639 std::set<HepMC::GenVertexPtr> visited;
640 unsigned int nParticlesFound = 0;
641 unsigned int nTreesFound = 0;
642 for (auto p: *hepMC) {
643 if ( (!p->production_vertex()) ||
644 (p->production_vertex()->particles_in_size() == 0) ) {
645 nTreesFound++;
646 std::cout << "\n Found new partial decay tree:\n" << std::endl;
647 unsigned int nParticlesVisited = printTree(std::move(p),visited,1,particleSet);
648 std::cout << "\n " << nParticlesVisited << " particles in this subtree" << std::endl;
649 nParticlesFound += nParticlesVisited;
650 }
651 }
652 std::cout << "\n Total of " << nParticlesFound << " particles found in "
653 << nTreesFound << " decay subtrees in HepMC event record\n" << std::endl;
654}
unsigned int printTree(HepMC::GenParticlePtr p, std::set< HepMC::GenVertexPtr > &visited, int level, std::set< HepMC::GenParticlePtr, ParticleIdCompare > *particleSet=nullptr)

◆ printTree()

unsigned int EvtInclusiveDecay::printTree ( HepMC::GenParticlePtr p,
std::set< HepMC::GenVertexPtr > & visited,
int level,
std::set< HepMC::GenParticlePtr, ParticleIdCompare > * particleSet = nullptr )
private

Definition at line 656 of file EvtInclusiveDecay.cxx.

656 {
657
658 unsigned int nParticlesVisited = 1;
659 for (int i=0; i<level; i++) std::cout << " ";
660 std::cout << pdgName(p,m_printHepMCHighlighted,particleSet);
661 auto v = p->end_vertex();
662 if (v) {
663 if (v->particles_in_size() > 1)
664 std::cout << " [interaction: " << v->particles_in_size() << " particles, vertex " << v << "] --> ";
665 else
666 std::cout << " --> ";
667 if (visited.insert(v).second) {
668 for (auto itp: *v) {
669 std::cout << pdgName(itp,m_printHepMCHighlighted,particleSet) << " ";
670 }
671 std::cout << std::endl;
672 for (auto itp: *v) {
673 if (itp->end_vertex())
674 nParticlesVisited += printTree(std::move(itp), visited, level+1, particleSet);
675 else
676 nParticlesVisited++;
677 }
678 } else
679 std::cout << "see above" << std::endl;
680 } else
681 std::cout << " no decay vertex\n" << std::endl;
682 return nParticlesVisited;
683}

◆ removeDecayTree()

void EvtInclusiveDecay::removeDecayTree ( HepMC::GenEvent * hepMC,
HepMC::GenParticlePtr p )
private

Definition at line 390 of file EvtInclusiveDecay.cxx.

390 {
391 auto v = p->end_vertex();
392 if (v) {
393 //This is recursive in HepMC3. But explicit deletion is allowed as well.
394 hepMC->remove_vertex(std::move(v));
395 p->set_status(1); // For now, flag particle as undecayed (stable)
396 ATH_MSG_DEBUG("Removed existing " << pdgName(p) << " " << p );
397 }
398}

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
384 }
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ reseedRandomEngine()

void EvtInclusiveDecay::reseedRandomEngine ( const std::string & streamName,
const EventContext & ctx )
private

Definition at line 180 of file EvtInclusiveDecay.cxx.

182{
183 long seeds[7];
184 ATHRNG::calculateSeedsMC21(seeds, streamName, ctx.eventID().event_number(), m_dsid, m_randomSeed);
185 m_evtAtRndmGen->getEngine()->setSeeds(seeds, 0); // NOT THREAD-SAFE
186}
void calculateSeedsMC21(long *seeds, const std::string &algName, uint64_t ev, uint64_t run, uint64_t offset=0)
Set the random seed using a string (e.g.

◆ setFilterPassed()

virtual void AthCommonAlgorithm< Gaudi::Algorithm >::setFilterPassed ( bool state,
const EventContext & ctx ) const
inlinevirtualinherited

Set filter decision:

Reimplemented in AthFilterAlgorithm.

Definition at line 99 of file AthCommonAlgorithm.h.

99 {
101 }
virtual void setFilterPassed(bool state, const EventContext &ctx) const
Set filter decision:

◆ sysExecute()

StatusCode AthCommonAlgorithm< Gaudi::Algorithm >::sysExecute ( const EventContext & ctx)
overridevirtualinherited

Execute an algorithm.

We override this in order to work around an issue with the Algorithm base class storing the event context in a member variable that can cause crashes in MT jobs.

Reimplemented in AthAnalysisAlgorithm.

Definition at line 80 of file AthCommonAlgorithm.cxx.

41{
42 return BaseAlg::sysExecute (ctx);
43}

◆ sysInitialize()

StatusCode AthCommonAlgorithm< Gaudi::Algorithm >::sysInitialize ( )
overridevirtualinherited

Override sysInitialize.

Override sysInitialize from the base class.

Loop through all output handles, and if they're WriteCondHandles, automatically register them and this Algorithm with the CondSvc

Scan through all outputHandles, and if they're WriteCondHandles, register them with the CondSvc

Reimplemented from AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >.

Reimplemented in AthAnalysisAlgorithm, AthFilterAlgorithm, AthHistogramAlgorithm, HypoBase, InputMakerBase, and PyAthena::Alg.

Definition at line 60 of file AthCommonAlgorithm.cxx.

71 {
73
74 if (sc.isFailure()) {
75 return sc;
76 }
77
78 ServiceHandle<ICondSvc> cs("CondSvc",name());
79 for (auto h : outputHandles()) {
80 if (h->isCondition() && h->mode() == Gaudi::DataHandle::Writer) {
81 // do this inside the loop so we don't create the CondSvc until needed
82 if ( cs.retrieve().isFailure() ) {
83 ATH_MSG_WARNING("no CondSvc found: won't autoreg WriteCondHandles");
85 }
86 if (cs->regHandle(this,*h).isFailure()) {
88 ATH_MSG_ERROR("unable to register WriteCondHandle " << h->fullKey()
89 << " with CondSvc");
90 }
91 }
92 }
93 return sc;
94}
#define ATH_MSG_WARNING(x,...)
virtual StatusCode sysInitialize() override
virtual std::vector< Gaudi::DataHandle * > outputHandles() const override

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ traverseDecayTree()

StatusCode EvtInclusiveDecay::traverseDecayTree ( HepMC::GenParticlePtr p,
bool isToBeRemoved,
std::set< HepMC::GenVertexPtr > & visited,
std::set< HepMC::GenParticlePtr, ParticleIdCompare > & toBeDecayed )
private

Definition at line 351 of file EvtInclusiveDecay.cxx.

354 {
355 ATH_MSG_VERBOSE("Inspecting: " << pdgName(p) << " " << p);
356 if (!isToBeRemoved) {
357 if (isToBeDecayed(p,true)) {
358 toBeDecayed.insert(p);
359 isToBeRemoved = true;
360 ATH_MSG_VERBOSE("Selected particle for decay: " << pdgName(p) << " " << p );
361
362 // In principle we could stop the recursion here. However, to prevent
363 // pathological cases in certain decay trees (in particular from Herwig),
364 // we continue in order to mark all descendants of this particle
365 // as visited. Thus none of these descendants can be flagged for further
366 // decay, even if it has several mothers.
367 }
368 }
369 auto v = p->end_vertex();
370 if (v) {
371 if (visited.insert(v).second) {
372 if ( isToBeRemoved && (v->particles_in_size()>1) && m_checkDecayTree ) {
373 ATH_MSG_WARNING("Found particle to be decayed with vertex with >1 incoming mother particles in decay tree");
374 ATH_MSG_WARNING( ([&p, &v](){ std::stringstream ss; HepMC::Print::line(ss,p); HepMC::Print::line(ss,v); return ss.str();})());
375 }
376 for (auto itp: *v) {
377 ATH_CHECK(traverseDecayTree(std::move(itp),isToBeRemoved,visited,toBeDecayed) );
378 }
379
380 }
381 }
382 return StatusCode::SUCCESS;
383}
static Double_t ss

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }

◆ xmlpath()

std::string EvtInclusiveDecay::xmlpath ( void )

Definition at line 723 of file EvtInclusiveDecay.cxx.

723 {
724 return PathResolverFindCalibDirectory( "Pythia8/xmldoc" );
725}
std::string PathResolverFindCalibDirectory(const std::string &logical_file_name)

Member Data Documentation

◆ m_allowAllKnownDecays

bool EvtInclusiveDecay::m_allowAllKnownDecays
private

Definition at line 123 of file EvtInclusiveDecay.h.

◆ m_allowDefaultBDecays

bool EvtInclusiveDecay::m_allowDefaultBDecays
private

Definition at line 124 of file EvtInclusiveDecay.h.

◆ m_applyUserSelection

bool EvtInclusiveDecay::m_applyUserSelection
private

Definition at line 141 of file EvtInclusiveDecay.h.

◆ m_blackList

std::vector<int> EvtInclusiveDecay::m_blackList
private

Definition at line 120 of file EvtInclusiveDecay.h.

◆ m_blackListSet

std::set<int> EvtInclusiveDecay::m_blackListSet
private

Definition at line 121 of file EvtInclusiveDecay.h.

◆ m_checkDecayChannels

bool EvtInclusiveDecay::m_checkDecayChannels
private

Definition at line 134 of file EvtInclusiveDecay.h.

◆ m_checkDecayTree

bool EvtInclusiveDecay::m_checkDecayTree
private

Definition at line 133 of file EvtInclusiveDecay.h.

◆ m_decayFile

std::string EvtInclusiveDecay::m_decayFile
private

Definition at line 109 of file EvtInclusiveDecay.h.

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default).

Definition at line 393 of file AthCommonDataStore.h.

◆ m_dsid

IntegerProperty EvtInclusiveDecay::m_dsid {this, "Dsid", 999999}
private

Definition at line 98 of file EvtInclusiveDecay.h.

98{this, "Dsid", 999999};

◆ m_evtAtRndmGen

EvtInclusiveAtRndmGen* EvtInclusiveDecay::m_evtAtRndmGen {}
private

Definition at line 104 of file EvtInclusiveDecay.h.

104{};

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default).

Definition at line 390 of file AthCommonDataStore.h.

◆ m_extendedExtraObjects

DataObjIDColl AthCommonAlgorithm< Gaudi::Algorithm >::m_extendedExtraObjects
privateinherited

Extra output dependency collection, extended by AthAlgorithmDHUpdate to add symlinks.

Empty if no symlinks were found.

Definition at line 108 of file AthCommonAlgorithm.h.

◆ m_inputKeyName

std::string EvtInclusiveDecay::m_inputKeyName
private

Definition at line 112 of file EvtInclusiveDecay.h.

◆ m_isfHerwig

bool EvtInclusiveDecay::m_isfHerwig
private

Definition at line 149 of file EvtInclusiveDecay.h.

◆ m_maxNRepeatedDecays

int EvtInclusiveDecay::m_maxNRepeatedDecays
private

Definition at line 139 of file EvtInclusiveDecay.h.

◆ m_mcEventKey

std::string GenBase::m_mcEventKey {}
protectedinherited

StoreGate key for the MC event collection (defaults to GEN_EVENT).

Definition at line 108 of file GenBase.h.

108{};

◆ m_mcevents_const

SG::ReadHandleKey<McEventCollection> GenBase::m_mcevents_const { this, "McEventKey", "GEN_EVENT", "StoreGate key of the MC event collection" }
privateinherited

Const handle to the MC event collection.

Definition at line 117 of file GenBase.h.

117{ this, "McEventKey", "GEN_EVENT", "StoreGate key of the MC event collection" };

◆ m_mkMcEvent

BooleanProperty GenBase::m_mkMcEvent {this, "MakeMcEvent", false, "Create a new MC event collection if it doesn't exist"}
protectedinherited

Flag to determine if a new MC event collection should be made if it doesn't exist.

Definition at line 110 of file GenBase.h.

110{this, "MakeMcEvent", false, "Create a new MC event collection if it doesn't exist"};

◆ m_myEvtGen

EvtGen* EvtInclusiveDecay::m_myEvtGen {}
private

Definition at line 105 of file EvtInclusiveDecay.h.

105{};

◆ m_noDecayChannels

std::map<int,long> EvtInclusiveDecay::m_noDecayChannels
private

Definition at line 135 of file EvtInclusiveDecay.h.

◆ m_nRepeatedDecays

int EvtInclusiveDecay::m_nRepeatedDecays
private

Definition at line 137 of file EvtInclusiveDecay.h.

◆ m_outputKeyName

std::string EvtInclusiveDecay::m_outputKeyName
private

Definition at line 113 of file EvtInclusiveDecay.h.

◆ m_pdtFile

std::string EvtInclusiveDecay::m_pdtFile
private

Definition at line 108 of file EvtInclusiveDecay.h.

◆ m_printHepMCAfterEvtGen

bool EvtInclusiveDecay::m_printHepMCAfterEvtGen
private

Definition at line 129 of file EvtInclusiveDecay.h.

◆ m_printHepMCBeforeEvtGen

bool EvtInclusiveDecay::m_printHepMCBeforeEvtGen
private

Definition at line 128 of file EvtInclusiveDecay.h.

◆ m_printHepMCHighlighted

bool EvtInclusiveDecay::m_printHepMCHighlighted
private

Definition at line 130 of file EvtInclusiveDecay.h.

◆ m_printHepMCHighLightTopLevelDecays

bool EvtInclusiveDecay::m_printHepMCHighLightTopLevelDecays
private

Definition at line 131 of file EvtInclusiveDecay.h.

◆ m_prohibitFinalStateDecay

bool EvtInclusiveDecay::m_prohibitFinalStateDecay
private

Definition at line 116 of file EvtInclusiveDecay.h.

◆ m_prohibitReDecay

bool EvtInclusiveDecay::m_prohibitReDecay
private

Definition at line 117 of file EvtInclusiveDecay.h.

◆ m_prohibitRemoveSelfDecay

bool EvtInclusiveDecay::m_prohibitRemoveSelfDecay
private

Definition at line 119 of file EvtInclusiveDecay.h.

◆ m_prohibitUnDecay

bool EvtInclusiveDecay::m_prohibitUnDecay
private

Definition at line 118 of file EvtInclusiveDecay.h.

◆ m_randomSeed

IntegerProperty EvtInclusiveDecay::m_randomSeed {this, "RandomSeed", 1234567, "Random seed for the built-in random engine"}
private

Seed for random number engine.

Definition at line 101 of file EvtInclusiveDecay.h.

101{this, "RandomSeed", 1234567, "Random seed for the built-in random engine"}; // FIXME make this into an unsigned long int?

◆ m_randomStreamName

std::string EvtInclusiveDecay::m_randomStreamName
private

Definition at line 111 of file EvtInclusiveDecay.h.

◆ m_readExisting

bool EvtInclusiveDecay::m_readExisting
private

Definition at line 115 of file EvtInclusiveDecay.h.

◆ m_rndmSvc

ServiceHandle<IAthRNGSvc> EvtInclusiveDecay::m_rndmSvc {this, "RndmSvc", "AthRNGSvc"}
private

Definition at line 95 of file EvtInclusiveDecay.h.

95{this, "RndmSvc", "AthRNGSvc"};

◆ m_setVMtransversePol

bool EvtInclusiveDecay::m_setVMtransversePol
private

Definition at line 150 of file EvtInclusiveDecay.h.

◆ m_userDecayFile

std::string EvtInclusiveDecay::m_userDecayFile
private

Definition at line 110 of file EvtInclusiveDecay.h.

◆ m_userSelMaxDimuMass

double EvtInclusiveDecay::m_userSelMaxDimuMass
private

Definition at line 148 of file EvtInclusiveDecay.h.

◆ m_userSelMinDimuMass

double EvtInclusiveDecay::m_userSelMinDimuMass
private

Definition at line 147 of file EvtInclusiveDecay.h.

◆ m_userSelMu1MaxEta

double EvtInclusiveDecay::m_userSelMu1MaxEta
private

Definition at line 145 of file EvtInclusiveDecay.h.

◆ m_userSelMu1MinPt

double EvtInclusiveDecay::m_userSelMu1MinPt
private

Definition at line 143 of file EvtInclusiveDecay.h.

◆ m_userSelMu2MaxEta

double EvtInclusiveDecay::m_userSelMu2MaxEta
private

Definition at line 146 of file EvtInclusiveDecay.h.

◆ m_userSelMu2MinPt

double EvtInclusiveDecay::m_userSelMu2MinPt
private

Definition at line 144 of file EvtInclusiveDecay.h.

◆ m_userSelRequireOppositeSignedMu

bool EvtInclusiveDecay::m_userSelRequireOppositeSignedMu
private

Definition at line 142 of file EvtInclusiveDecay.h.

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ m_whiteList

std::vector<int> EvtInclusiveDecay::m_whiteList
private

Definition at line 125 of file EvtInclusiveDecay.h.

◆ m_whiteListSet

std::set<int> EvtInclusiveDecay::m_whiteListSet
private

Definition at line 126 of file EvtInclusiveDecay.h.


The documentation for this class was generated from the following files: