ISF::TruthSvc Class Reference

#include <TruthSvc.h>

Inheritance diagram for ISF::TruthSvc:

Collaboration diagram for ISF::TruthSvc:

Public Member Functions
	TruthSvc (const std::string &name, ISvcLocator *pSvcLocator)
	Constructor. More...
virtual	~TruthSvc ()=default
	Destructor. More...
StatusCode	initialize () override
	Athena algorithm's interface method initialize() More...
StatusCode	finalize () override
	Athena algorithm's interface method finalize() More...
void	registerTruthIncident (ITruthIncident &truthincident, bool saveAllChildren=false) const override
	Register a truth incident. More...
StatusCode	initializeTruthCollection (int largestGeneratedParticleBC=0, int largestGeneratedVertexBC=0) override
	Initialize the Truth Svc at the beginning of each event. More...
StatusCode	releaseEvent () override
	Finalize the Truth Svc at the end of each event. More...

Private Member Functions
void	recordIncidentToMCTruth (ITruthIncident &truthincident, bool passWholeVertex) const
	Record the given truth incident to the MC Truth. More...
HepMC::GenVertexPtr	createGenVertexFromTruthIncident (ITruthIncident &truthincident) const
	Record and end vertex to the MC Truth for the parent particle. More...

Private Attributes
ServiceHandle< Barcode::IBarcodeSvc >	m_barcodeSvc {this, "BarcodeSvc", "BarcodeSvc", ""}
	The Barcode service. More...
ToolHandleArray< ITruthStrategy >	m_truthStrategies {this, "TruthStrategies", {}, ""}
	the truth strategies applied (as AthenaToolHandle Array) More...
ITruthStrategy **	m_geoStrategies [AtlasDetDescr::fNumAtlasRegions]
	for faster access: using an internal pointer to the actual ITruthStrategy instances More...
unsigned short	m_numStrategies [AtlasDetDescr::fNumAtlasRegions]
Gaudi::Property< bool >	m_skipIfNoChildren {this, "SkipIfNoChildren", true, ""}
	MCTruth steering. More...
Gaudi::Property< bool >	m_skipIfNoParentId {this, "SkipIfNoParentId", true, ""}
	do not record if parentId == HepMC::UNDEFINED_ID More...
Gaudi::Property< bool >	m_ignoreUndefinedBarcodes {this, "IgnoreUndefinedBarcodes", false, ""}
	do/don't abort if retrieve an undefined barcode More...
Gaudi::Property< bool >	m_passWholeVertex {this, "PassWholeVertices", true, ""}
Gaudi::Property< std::vector< unsigned int > >	m_forceEndVtxRegionsVec {this, "ForceEndVtxInRegions", {}, ""}
	property containing AtlasRegions for which More...
std::array< bool, AtlasDetDescr::fNumAtlasRegions >	m_forceEndVtx
	attach end vertex to More...
Gaudi::Property< bool >	m_quasiStableParticlesIncluded {this, "QuasiStableParticlesIncluded", false, ""}
	does this job simulate quasi-stable particles. More...
Gaudi::Property< bool >	m_quasiStableParticleOverwrite {this, "QuasiStableParticleOverwrite", true, ""}
	should the job be overwriting QS particles More...

Detailed Description

HepMC based version of the ISF::ITruthSvc, currently it takes an ITruthIncident base class

Author

Andreas.Salzburger -at- cern.ch , Elmar.Ritsch -at- cern.ch

Definition at line 50 of file TruthSvc.h.

Constructor & Destructor Documentation

TruthSvc()

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

Constructor.

Definition at line 49 of file TruthSvc.cxx.

                                                            :
  base_class(name,svc),
  m_geoStrategies(),
  m_numStrategies()
{
}

~TruthSvc()

virtual ISF::TruthSvc::~TruthSvc ( )

virtualdefault

Destructor.

Member Function Documentation

createGenVertexFromTruthIncident()

HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident ( ISF::ITruthIncident &  ti ) const

private

Record and end vertex to the MC Truth for the parent particle.

Record the given truth incident to the MC Truth.

Definition at line 300 of file TruthSvc.cxx.

                                                                                              {
 
  int processCode = ti.physicsProcessCode();
  int       parentBC = ti.parentBarcode();
 
  std::vector<double> weights(1);
  int primaryBC = parentBC % HepMC::SIM_REGENERATION_INCREMENT;
  weights[0] = static_cast<double>( primaryBC ); // FIXME vertex weights should not be used to encode other info.
 
  // Check for a previous end vertex on this particle.  If one existed, then we should put down next to this
  //  a new copy of the particle.  This is the agreed upon version of the quasi-stable particle truth, where
  //  the vertex at which we start Q-S simulation no longer conserves energy, but we keep both copies of the
  //  truth particles
  HepMC::GenParticlePtr  parent = ti.parentParticle();
  if (!parent) {
    ATH_MSG_ERROR("Unable to write particle interaction to MC truth due to missing parent HepMC::GenParticle instance");
    abort();
  }
  HepMC::GenEvent *mcEvent = parent->parent_event();
  if (!mcEvent) {
    ATH_MSG_ERROR("Unable to write particle interaction to MC truth due to missing parent HepMC::GenEvent instance");
    abort();
  }
 
  // generate vertex
  int vtxbcode = m_barcodeSvc->newSimulationVertex(); // TODO replace barcodeSvc
  if ( vtxbcode == HepMC::UNDEFINED_ID) {
    if (m_ignoreUndefinedBarcodes) {
      ATH_MSG_WARNING("Unable to generate new Truth Vertex Barcode. Continuing due to 'IgnoreUndefinedBarcodes'==True");
    } else {
      ATH_MSG_ERROR("Unable to generate new Truth Vertex Barcode. Aborting");
      abort();
    }
  }
  const int vtxStatus = 1000 + static_cast<int>(processCode) + HepMC::SIM_STATUS_THRESHOLD;
#ifdef HEPMC3
  auto newVtx = HepMC::newGenVertexPtr( ti.position(),vtxStatus);
#else
  // NB In HepMC2 there is no GenVertex status, so we set the GenVertex ID.
  std::unique_ptr<HepMC::GenVertex> newVtx = std::make_unique<HepMC::GenVertex>( ti.position(), vtxStatus, weights );
  HepMC::suggest_barcode( newVtx.get(), vtxbcode );
#endif
 
  if (parent->end_vertex()){
      ATH_MSG_ERROR ("createGVfromTI: Parent particle found with an end vertex attached.  This should not happen!");
      ATH_MSG_ERROR ("createGVfromTI: Parent 1: " << parent << ", barcode: " << HepMC::barcode(parent));
      ATH_MSG_ERROR ( "createGVfromTI: parent->end_vertex(): " << parent->end_vertex() << ", barcode: " << HepMC::barcode(parent->end_vertex()) );
     abort();
  } else { // Normal simulation
#ifdef DEBUG_TRUTHSVC
    ATH_MSG_VERBOSE ("createGVfromTI Parent 1: " << parent << ", barcode: " << HepMC::barcode(parent));
#endif
    // add parent particle to newVtx
    newVtx->add_particle_in( parent );
#ifdef DEBUG_TRUTHSVC
    ATH_MSG_VERBOSE ( "createGVfromTI End Vertex representing process: " << processCode << ", for parent with barcode "<<parentBC<<". Creating." );
    ATH_MSG_VERBOSE ( "createGVfromTI Parent 2: " << parent << ", barcode: " << HepMC::barcode(parent));
#endif
#ifdef HEPMC3
    mcEvent->add_vertex(newVtx);
    HepMC::suggest_barcode( newVtx, vtxbcode );
    newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights));
#else
    mcEvent->add_vertex( newVtx.release() );
#endif
  }
 
  return parent->end_vertex();
}

finalize()

StatusCode ISF::TruthSvc::finalize ( )

override

Athena algorithm's interface method finalize()

framework methods

Definition at line 105 of file TruthSvc.cxx.

{
  ATH_MSG_VERBOSE("Finalizing ...");
  return StatusCode::SUCCESS;
}

initialize()

StatusCode ISF::TruthSvc::initialize ( )

override

Athena algorithm's interface method initialize()

framework methods

Definition at line 57 of file TruthSvc.cxx.

{
  ATH_MSG_VERBOSE( "initialize()" );
 
  // Screen output
  ATH_MSG_DEBUG("--------------------------------------------------------");
 
  // retrieve BarcodeSvc
  if ( m_barcodeSvc.retrieve().isFailure() ) {
    ATH_MSG_FATAL("Could not retrieve BarcodeService. Abort.");
    return StatusCode::FAILURE;
  }
 
  // copy a pointer to the strategy instance to the local
  // array of pointers (for faster access)
  ATH_CHECK(m_truthStrategies.retrieve());
  // Would be nicer to make m_geoStrategies a vector of vectors
  for ( unsigned short geoID=AtlasDetDescr::fFirstAtlasRegion; geoID<AtlasDetDescr::fNumAtlasRegions; ++geoID) {
    m_numStrategies[geoID] = 0;
    for ( const auto &truthStrategy : m_truthStrategies) {
      if(truthStrategy->appliesToRegion(geoID)) {
        ++m_numStrategies[geoID];
      }
    }
  }
  for ( unsigned short geoID=AtlasDetDescr::fFirstAtlasRegion; geoID<AtlasDetDescr::fNumAtlasRegions; ++geoID) {
    m_geoStrategies[geoID] = new ISF::ITruthStrategy*[m_numStrategies[geoID]];
    unsigned short curNumStrategies = m_truthStrategies.size();
    unsigned short nStrat(0);
    for ( unsigned short i = 0; i < curNumStrategies; ++i) {
      if(m_truthStrategies[i]->appliesToRegion(geoID)) {
        m_geoStrategies[geoID][nStrat++] = &(*m_truthStrategies[i]);
      }
    }
 
    // setup whether we want to write end-vertices in this region whenever a truth particle dies
    // create an end-vertex for all truth particles ending in the current AtlasRegion?
    bool forceEndVtx = std::find( m_forceEndVtxRegionsVec.begin(),
                                  m_forceEndVtxRegionsVec.end(),
                                  geoID ) != m_forceEndVtxRegionsVec.end();
    m_forceEndVtx[geoID] = forceEndVtx;
  }
  ATH_MSG_VERBOSE("initialize() successful");
  return StatusCode::SUCCESS;
}

initializeTruthCollection()

StatusCode ISF::TruthSvc::initializeTruthCollection ( int  largestGeneratedParticleBC = 0, int  largestGeneratedVertexBC = 0  )

override

Initialize the Truth Svc at the beginning of each event.

Initialize the TruthSvc and the truthSvc.

Definition at line 113 of file TruthSvc.cxx.

{
  ATH_CHECK( m_barcodeSvc->initializeBarcodes(largestGeneratedParticleBC, largestGeneratedVertexBC) );
  return StatusCode::SUCCESS;
}

recordIncidentToMCTruth()

void ISF::TruthSvc::recordIncidentToMCTruth ( ISF::ITruthIncident &  ti, bool  passWholeVertex  ) const

private

Record the given truth incident to the MC Truth.

Definition at line 199 of file TruthSvc.cxx.

                                                                                            {
#ifdef  DEBUG_TRUTHSVC
  ATH_MSG_INFO("Starting recordIncidentToMCTruth(...)");
#endif
  int       parentBC = ti.parentBarcode();
 
  if (ti.parentParticle()->end_vertex()) {
    ATH_MSG_WARNING ("Attempting to record a TruthIncident for a particle which has already decayed!");
    ATH_MSG_WARNING ("No action will be taken - please fix client code.");
    return;
  }
 
  // record the GenVertex
  HepMC::GenVertexPtr  vtxFromTI = createGenVertexFromTruthIncident( ti );
  const ISF::InteractionClass_t classification = ti.interactionClassification();
#ifdef DEBUG_TRUTHSVC
  const std::string survival = (ti.parentSurvivesIncident()) ? "parent survives" : "parent destroyed";
  const std::string vtxType = (ti.interactionClassification()==ISF::STD_VTX) ? "Normal" : "Quasi-stable";
  ATH_MSG_INFO("TruthSvc: " << vtxType << " vertex + " << survival
               << ", TI Class: " << ti.interactionClassification()
               << ", ProcessType: " << ti.physicsProcessCategory()
               << ", ProcessSubType: " << ti.physicsProcessCode());
#endif
 
  ATH_MSG_VERBOSE ( "Outgoing particles:" );
  // update parent barcode and add it to the vertex as outgoing particle
  int newPrimaryBC = HepMC::UNDEFINED_ID;
  if (classification == ISF::QS_SURV_VTX) {
    // Special case when a particle with a pre-defined decay interacts
    // and survives.
    // Set the barcode to the next available value below the simulation
    // barcode offset.
    newPrimaryBC = m_barcodeSvc->newGeneratedParticle(parentBC);
  }
  else {
    newPrimaryBC = parentBC + HepMC::SIM_REGENERATION_INCREMENT;
  }
 
  HepMC::GenParticlePtr  parentAfterIncident = ti.parentParticleAfterIncident( newPrimaryBC ); // This call changes ti.parentParticle() output
  if(parentAfterIncident) {
    if (classification==ISF::QS_SURV_VTX) {
      // Special case when a particle with a pre-defined decay
      // interacts and survives.
      // As the parentParticleAfterIncident has a pre-defined decay
      // its status should end in 2, but we should flag that it has
      // survived an interaction.
      if (!MC::isDecayed(parentAfterIncident)) {
        ATH_MSG_WARNING ( "recordIncidentToMCTruth - check parentAfterIncident: " << parentAfterIncident );
      }
    }
    vtxFromTI->add_particle_out( parentAfterIncident );
#ifdef HEPMC3
    HepMC::suggest_barcode( parentAfterIncident, newPrimaryBC ); // TODO check this works correctly
#endif
    // NB For ISFTruthIncident the m_parent ISFParticle still needs
    // its id and particleLink properties to be properly updated at
    // this point.
    ATH_MSG_VERBOSE ( "Parent After Incident: " << parentAfterIncident << ", barcode: " << HepMC::barcode(parentAfterIncident));
  }
 
  const bool isQuasiStableVertex = (classification == ISF::QS_PREDEF_VTX); // QS_DEST_VTX and QS_SURV_VTX should be treated as normal from now on.
  // add child particles to the vertex
  const unsigned short numSec = ti.numberOfChildren();
  for ( unsigned short i=0; i<numSec; ++i) {
    bool writeOutChild = isQuasiStableVertex || passWholeVertex || ti.childPassedFilters(i);
    if (writeOutChild) {
      HepMC::GenParticlePtr  p = nullptr;
      // generate a new barcode for the child particle
      int secondaryParticleBC = (isQuasiStableVertex) ?
        m_barcodeSvc->newGeneratedParticle(parentBC) : m_barcodeSvc->newSecondaryParticle(parentBC); // TODO replace m_barcodeSvc
      if ( secondaryParticleBC == HepMC::UNDEFINED_ID) {
        if (m_ignoreUndefinedBarcodes)
          ATH_MSG_WARNING("Unable to generate new Secondary Particle Barcode. Continuing due to 'IgnoreUndefinedBarcodes'==True");
        else {
          ATH_MSG_ERROR("Unable to generate new Secondary Particle Barcode. Aborting");
          abort();
        }
      }
      p = ti.childParticle( i, secondaryParticleBC ); // potentially overrides secondaryParticleBC
      if (p) {
        // add particle to vertex
        vtxFromTI->add_particle_out( p);
#ifdef HEPMC3
        int secondaryParticleBCFromTI = ti.childBarcode(i);
        HepMC::suggest_barcode( p, secondaryParticleBCFromTI ? secondaryParticleBCFromTI : secondaryParticleBC );
        // NB For ISFTruthIncident the current child ISFParticle still needs
        // its id and particleLink properties to be properly updated at
        // this point.
#endif
      }
      ATH_MSG_VERBOSE ( "Writing out " << i << "th child particle: " << p << ", barcode: " << HepMC::barcode(p));
    } // <-- if write out child particle
    else {
      ATH_MSG_VERBOSE ( "Not writing out " << i << "th child particle." );
    }
 
  } // <-- loop over all child particles
  ATH_MSG_VERBOSE("--------------------------------------------------------");
}

registerTruthIncident()

void ISF::TruthSvc::registerTruthIncident ( ISF::ITruthIncident &  ti, bool  saveAllChildren = false  ) const

override

Register a truth incident.

Definition at line 125 of file TruthSvc.cxx.

                                                                                          {
 
  const bool passWholeVertex = m_passWholeVertex || saveAllChildren;
  // pass whole vertex or individual child particles
  ti.setPassWholeVertices(passWholeVertex);
 
  // the GeoID
  AtlasDetDescr::AtlasRegion geoID = ti.geoID();
 
  // check geoID assigned to the TruthIncident
  if ( !validAtlasRegion(geoID) ) {
    const auto& position = ti.position();
    ATH_MSG_ERROR("Unable to register truth incident with unknown SimGeoID="<< geoID
                  << " at position z=" << position.z() << " r=" << position.perp());
    return;
  }
 
  ATH_MSG_VERBOSE( "Registering TruthIncident for SimGeoID="
                   << AtlasDetDescr::AtlasRegionHelper::getName(geoID) );
 
  // number of child particles
  const unsigned short numSec = ti.numberOfChildren();
  if ( m_skipIfNoChildren && (numSec==0) ) {
    ATH_MSG_VERBOSE( "No child particles present in the TruthIncident,"
                     << " will not record this TruthIncident.");
    return;
  }
 
  // the parent particle -> get its id
  const int parentID = ti.parentUniqueID();
  if ( m_skipIfNoParentId && (parentID == HepMC::UNDEFINED_ID) ) {
    ATH_MSG_VERBOSE( "Parent particle in TruthIncident does not have an id,"
                     << " will not record this TruthIncident.");
    return;
  }
 
  // loop over registered truth strategies for given geoID
  bool pass = false;
  for ( unsigned short stratID=0; (!pass) && (stratID<m_numStrategies[geoID]); stratID++) {
    // (*) test if given TruthIncident passes current strategy
    pass = m_geoStrategies[geoID][stratID]->pass(ti);
  }
 
  if (pass) {
    ATH_MSG_VERBOSE("At least one TruthStrategy passed.");
    // at least one truth strategy returned true
    //  -> record incident
    recordIncidentToMCTruth(ti, passWholeVertex);
 
  } else {
    // none of the truth strategies returned true
    //  -> child particles will NOT be added to the TruthEvent collection
    // attach parent particle end vertex if it gets killed by this interaction
    if ( m_forceEndVtx[geoID] && !ti.parentSurvivesIncident() ) {
      ATH_MSG_VERBOSE("No TruthStrategies passed and parent destroyed - create end vertex.");
      this->createGenVertexFromTruthIncident( ti );
 
#ifdef DEBUG_TRUTHSVC
      const std::string survival = (ti.parentSurvivesIncident()) ? "parent survives" : "parent destroyed";
      const std::string vtxType = (ti.interactionClassification()==ISF::STD_VTX) ? "Normal" : "Quasi-stable";
      ATH_MSG_INFO("TruthSvc: " << vtxType << " vertex + " << survival
                   << ", TI Class: " << ti.interactionClassification()
                   << ", ProcessType: " << ti.physicsProcessCategory()
                   << ", ProcessSubType: " << ti.physicsProcessCode());
#endif
 
    }
 
  }
 
  return;
}

releaseEvent()

StatusCode ISF::TruthSvc::releaseEvent ( )

override

Finalize the Truth Svc at the end of each event.

Definition at line 119 of file TruthSvc.cxx.

                                     {
  return StatusCode::SUCCESS;
}

Member Data Documentation

m_barcodeSvc

ServiceHandle<Barcode::IBarcodeSvc> ISF::TruthSvc::m_barcodeSvc {this, "BarcodeSvc", "BarcodeSvc", ""}

private

The Barcode service.

Definition at line 83 of file TruthSvc.h.

m_forceEndVtx

std::array<bool, AtlasDetDescr::fNumAtlasRegions> ISF::TruthSvc::m_forceEndVtx

private

attach end vertex to

Definition at line 100 of file TruthSvc.h.

m_forceEndVtxRegionsVec

Gaudi::Property<std::vector<unsigned int> > ISF::TruthSvc::m_forceEndVtxRegionsVec {this, "ForceEndVtxInRegions", {}, ""}

private

property containing AtlasRegions for which

Definition at line 98 of file TruthSvc.h.

m_geoStrategies

ITruthStrategy** ISF::TruthSvc::m_geoStrategies[AtlasDetDescr::fNumAtlasRegions]

private

for faster access: using an internal pointer to the actual ITruthStrategy instances

Definition at line 88 of file TruthSvc.h.

m_ignoreUndefinedBarcodes

Gaudi::Property<bool> ISF::TruthSvc::m_ignoreUndefinedBarcodes {this, "IgnoreUndefinedBarcodes", false, ""}

private

do/don't abort if retrieve an undefined barcode

Definition at line 94 of file TruthSvc.h.

m_numStrategies

unsigned short ISF::TruthSvc::m_numStrategies[AtlasDetDescr::fNumAtlasRegions]

private

Definition at line 89 of file TruthSvc.h.

m_passWholeVertex

Gaudi::Property<bool> ISF::TruthSvc::m_passWholeVertex {this, "PassWholeVertices", true, ""}

private

Definition at line 96 of file TruthSvc.h.

m_quasiStableParticleOverwrite

Gaudi::Property<bool> ISF::TruthSvc::m_quasiStableParticleOverwrite {this, "QuasiStableParticleOverwrite", true, ""}

private

should the job be overwriting QS particles

Definition at line 104 of file TruthSvc.h.

m_quasiStableParticlesIncluded

Gaudi::Property<bool> ISF::TruthSvc::m_quasiStableParticlesIncluded {this, "QuasiStableParticlesIncluded", false, ""}

private

does this job simulate quasi-stable particles.

Definition at line 103 of file TruthSvc.h.

m_skipIfNoChildren

Gaudi::Property<bool> ISF::TruthSvc::m_skipIfNoChildren {this, "SkipIfNoChildren", true, ""}

private

MCTruth steering.

do not record incident if numChildren==0

Definition at line 92 of file TruthSvc.h.

m_skipIfNoParentId

Gaudi::Property<bool> ISF::TruthSvc::m_skipIfNoParentId {this, "SkipIfNoParentId", true, ""}

private

do not record if parentId == HepMC::UNDEFINED_ID

Definition at line 93 of file TruthSvc.h.

m_truthStrategies

ToolHandleArray<ITruthStrategy> ISF::TruthSvc::m_truthStrategies {this, "TruthStrategies", {}, ""}

private

the truth strategies applied (as AthenaToolHandle Array)

Definition at line 86 of file TruthSvc.h.

---

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

• TruthSvc.h
• TruthSvc.cxx