ATLAS Offline Software
Loading...
Searching...
No Matches
ISF::SimKernelMT Class Referencefinal

Core Athena algorithm for the Integrated Simulation Framework. More...

#include <SimKernelMT.h>

Inheritance diagram for ISF::SimKernelMT:

Public Member Functions

StatusCode initialize () override
 Check user configuration and configure the algorithm state accordingly.
StatusCode execute () override
 Check the validity of the MC truth event in the current Athena event.
StatusCode finalize () override
 Implements empty finalize (implementation required by AthAlgorithm)
bool isClonable () const override
 AthAlgorithm (const std::string &name, ISvcLocator *pSvcLocator)
 Constructor with parameters:
virtual StatusCode sysInitialize () override
 Override sysInitialize.
virtual const DataObjIDColl & extraOutputDeps () const override
 Return the list of extra output dependencies.
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

Protected Member Functions

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

ISimulatorToolidentifySimulator (const ISF::ISFParticle &particle)
 Returns the simulator to use for the given particle.
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

SG::ReadHandleKey< McEventCollectionm_inputEvgenKey {this, "InputEvgenCollection", "", "Input EVGEN collection."}
 Input Generator Truth collection.
SG::WriteHandleKey< McEventCollectionm_outputTruthKey {this, "OutputTruthCollection", "", "Output Truth collection."}
 Output Simulation Truth collection.
SG::WriteHandleKey< TrackRecordCollectionm_caloEntryLayerKey {this, "CaloEntryLayerKey", "CaloEntryLayer", ""}
 Output TrackRecordCollections.
SG::WriteHandleKey< TrackRecordCollectionm_muonEntryLayerKey {this, "MuonEntryLayerKey", "MuonEntryLayer", ""}
SG::WriteHandleKey< TrackRecordCollectionm_muonExitLayerKey {this, "MuonExitLayerKey", "MuonExitLayer", ""}
BooleanProperty m_useShadowEvent {this, "UseShadowEvent", false, "New approach to selecting particles for simulation" }
Gaudi::Property< bool > m_forceGeoIDSvc {this, "AlwaysUseGeoIDSvc", false, "Force geoID recalculation for each particle" }
 Force geoID recalculation for each particle.
ToolHandle< IGenEventFilterm_truthPreselectionTool {this, "TruthPreselectionTool", "", "Tool for filtering out quasi-stable particle daughters"}
ServiceHandle< IInputConverterm_inputConverter {this, "InputConverter", "", "Input McEventCollection->ISFParticleContainer conversion service."}
 Input converter service (from Generator->ISF particle types)
ServiceHandle< ITruthSvcm_truthRecordSvc {this,"TruthRecordService", "ISF_MC15aPlusTruthService", ""}
 Central truth service.
ServiceHandle< Simulation::IZeroLifetimePatcherm_qspatcher {this, "QuasiStablePatcher", "", "Quasi-Stable Particle Simulation Patcher"}
 Quasi-Stable Particle Simulation Patcher.
ToolHandleArray< ISimulatorToolm_simulationTools {this, "SimulationTools", {}, ""}
 Simulation Tools.
ISimulatorToolm_particleKillerTool {}
 When no appropriate simulator can be found for a given particle, the particle is sent to this "particle killer":
PublicToolHandle< IEntryLayerToolm_entryLayerTool {this, "EntryLayerTool", "ISF_EntryLayerToolMT", ""}
 AthenaTool responsible for writing Calo/Muon Entry/Exit Layer collection.
ServiceHandle< IGeoIDSvcm_geoIDSvc {this, "GeoIDSvc", "", "Service to set particle GeoIDs"}
 Service to set particle GeoIDs.
ToolHandle< IParticleOrderingToolm_orderingTool {this, "ParticleOrderingTool", "", "Tool to set order of particles"}
 Tool to set particle ordering.
std::array< PublicToolHandleArray< ISimulationSelector >, AtlasDetDescr::fNumAtlasRegionsm_simSelectors
 The simulation selectors defining the "routing chain".
std::map< ISF::SimulationFlavor, ISimulatorTool * > m_simToolMap
 Map of the simulation flavours used in this job to the corresponding Simulation Services.
Gaudi::Property< size_t > m_maxParticleVectorSize {this, "MaximumParticleVectorSize", 10240}
 Number of particles simultaneously sent to simulator.
DataObjIDColl m_extendedExtraObjects
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

Friends

class ISFTesting::SimKernelMT_test
 Allow the test class access to all methods.

Detailed Description

Core Athena algorithm for the Integrated Simulation Framework.

This Athena algorithm dispatches all selected particles from the input collection to the attached simulators.

Definition at line 61 of file SimKernelMT.h.

Member Typedef Documentation

◆ StoreGateSvc_t

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

Definition at line 388 of file AthCommonDataStore.h.

Member Function Documentation

◆ AthAlgorithm()

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

Constructor with parameters:

Definition at line 51 of file AthAlgorithm.cxx.

25 :
27{
28 // Set up to run AthAlgorithmDHUpdate in sysInitialize before
29 // merging dependency lists. This extends the output dependency
30 // list with any symlinks implied by inheritance relations.
31 m_updateDataHandles =
32 std::make_unique<AthenaBaseComps::AthAlgorithmDHUpdate>
34 std::move (m_updateDataHandles));
35}
DataObjIDColl m_extendedExtraObjects
AthCommonDataStore(const std::string &name, T... args)

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< 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 }
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< 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< Algorithm > >::detStore ( ) const
inlineinherited

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

Definition at line 95 of file AthCommonDataStore.h.

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< 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 ISF::SimKernelMT::execute ( )
override

Check the validity of the MC truth event in the current Athena event.

Return values
StatusCode::SUCCESSall HepMC::GenEvents of the current Athena event are valid
StatusCode::FAILUREthe MC truth collection could not be read
StatusCode::FAILUREthe MC truth collection does not have a valid tree structure

Definition at line 115 of file SimKernelMT.cxx.

115 {
116
117 const EventContext& ctx = Gaudi::Hive::currentContext();
118 auto hitCollections = std::make_shared<HitCollectionMap>();
119 // Call setupEvent for all simulators (TODO: make the tools do this)
120 for (auto& curSimTool: m_simulationTools) {
121 if ( curSimTool ) {
122 if (auto* g4sim =
123 dynamic_cast<ISF::BaseSimulatorG4Tool*>(curSimTool.get())) {
124 ATH_CHECK(g4sim->setupEvent(ctx, *hitCollections));
125 } else {
126 ATH_CHECK(curSimTool->setupEvent(ctx));
127 }
128 ATH_MSG_DEBUG( "Event setup done for " << curSimTool->name() );
129 }
130 }
131
132 SG::ReadHandle<McEventCollection> inputEvgen(m_inputEvgenKey);
133 if (!inputEvgen.isValid()) {
134 ATH_MSG_FATAL("Unable to read input GenEvent collection '" << inputEvgen.key() << "'");
135 return StatusCode::FAILURE;
136 }
137
138 // create output Truth collection
139 SG::WriteHandle<McEventCollection> outputTruth(m_outputTruthKey, ctx);
140 std::unique_ptr<McEventCollection> shadowTruth{};
141 if (m_useShadowEvent) {
142 outputTruth = std::make_unique<McEventCollection>();
143 // copy input Evgen collection to shadow Truth collection
144 shadowTruth = std::make_unique<McEventCollection>(*inputEvgen);
145 for (HepMC::GenEvent* currentGenEvent : *shadowTruth ) {
146 // Apply QS patch if required
147 if ( not m_qspatcher.empty() ) {
148 ATH_CHECK(m_qspatcher->applyWorkaround(*currentGenEvent));
149 }
150 // Copy GenEvent and remove daughters of quasi-stable particles to be simulated
151 std::unique_ptr<HepMC::GenEvent> outputEvent = m_truthPreselectionTool->filterGenEvent(*currentGenEvent);
152 outputTruth->push_back(outputEvent.release());
153 }
154 }
155 else {
156 // copy input Evgen collection to output Truth collection
157 outputTruth = std::make_unique<McEventCollection>(*inputEvgen);
158 // Apply QS patch if required
159 if ( not m_qspatcher.empty() ) {
160 for (HepMC::GenEvent* currentGenEvent : *outputTruth ) {
161 ATH_CHECK(m_qspatcher->applyWorkaround(*currentGenEvent));
162 }
163 }
164 }
165
166 const int largestGeneratedParticleBC = (outputTruth->empty()) ? HepMC::UNDEFINED_ID
167 : HepMC::maxGeneratedParticleBarcode(outputTruth->at(0)); // TODO make this more robust
168 const int largestGeneratedVertexBC = (outputTruth->empty()) ? HepMC::UNDEFINED_ID
169 : HepMC::maxGeneratedVertexBarcode(outputTruth->at(0)); // TODO make this more robust
170 // tell TruthService we're starting a new event
171 ATH_CHECK( m_truthRecordSvc->initializeTruthCollection(largestGeneratedParticleBC, largestGeneratedVertexBC) );
172
173 // Create TrackRecordCollections and pass them to the entryLayerTool
174 SG::WriteHandle<TrackRecordCollection> caloEntryLayer(m_caloEntryLayerKey, ctx);
175 caloEntryLayer = std::make_unique<TrackRecordCollection>(caloEntryLayer.name());
176 ATH_CHECK(m_entryLayerTool->registerTrackRecordCollection(caloEntryLayer.ptr(), fAtlasCaloEntry));
177 SG::WriteHandle<TrackRecordCollection> muonEntryLayer(m_muonEntryLayerKey, ctx);
178 muonEntryLayer = std::make_unique<TrackRecordCollection>(muonEntryLayer.name());
179 ATH_CHECK(m_entryLayerTool->registerTrackRecordCollection(muonEntryLayer.ptr(), fAtlasMuonEntry));
180 SG::WriteHandle<TrackRecordCollection> muonExitLayer(m_muonExitLayerKey, ctx);
181 muonExitLayer = std::make_unique<TrackRecordCollection>(muonExitLayer.name());
182 ATH_CHECK(m_entryLayerTool->registerTrackRecordCollection(muonExitLayer.ptr(), fAtlasMuonExit));
183
184 // read and convert input
185 ISFParticleContainer simParticles;
186 ATH_CHECK( m_inputConverter->convert(*outputTruth, simParticles) );
187
188 // create an ordered queue of particles
189 ISFParticleOrderedQueue particleQueue;
190 for ( auto& particle : simParticles ) {
191
192 if ( m_orderingTool.empty() ) {
193 // Without a defined ordering, preserve old FIFO behaviour
194 particle->setOrder( -particleQueue.size() );
195 } else {
196 m_orderingTool->setOrder( *particle );
197 }
198
199 particleQueue.push( particle );
200 }
201
202 unsigned int loopCounter{0};
203 // loop until there are no more particles to simulate
205 ISimulatorTool* lastSimulator{};
206 ISFParticleContainer newSecondaries{};
207 while ( particleQueue.size() ) {
208 ++loopCounter;
209 ATH_MSG_VERBOSE("Main Loop pass no. " << loopCounter);
210 ATH_MSG_VERBOSE("Queue starts with " << particleQueue.size() << " particles.");
211 // Create a vector of particles with the same simulator
212 ISFParticleOrderedQueue tempQueue;
213 while ( particleQueue.size() ) {
214 auto particlePtr = particleQueue.top();
215 ISFParticle& curParticle( *particlePtr );
216 particleQueue.pop();
217
218 // Get the geo ID for the particle
219 if ( m_forceGeoIDSvc || !validAtlasRegion( curParticle.nextGeoID() ) ) {
220 m_geoIDSvc->identifyAndRegNextGeoID( curParticle );
221 }
222
223 // Get the simulator using the GeoID
224 auto& simTool = identifySimulator(curParticle);
225
226 // Fill the vector
227 if ( particles.empty() ) {
228 // First particle in the vector defines the simulator
229 particles.push_back(particlePtr);
230 lastSimulator=&simTool;
231 }
232 else if (&simTool!=lastSimulator || particles.size() >= m_maxParticleVectorSize.value() ) {
233 // Change of simulator, end the current vector
234 tempQueue.push(particlePtr);
235 }
236 else {
237 particles.push_back(particlePtr);
238 }
239 }
240 particleQueue = std::move(tempQueue);
241 #ifdef ISFDEBUG
242 if (loopCounter>100 && particles.size()<3) {
243 ATH_MSG_INFO("Main Loop pass no. " << loopCounter);
244 ATH_MSG_INFO("Selected " << particles.size() << " particles to be processed by " << lastSimulator->name());
245 for ( const ISFParticle *particle : particles ) {
246 ATH_MSG_INFO(*particle);
247 }
248 }
249 #endif // ISFDEBUG
250
251 ATH_MSG_VERBOSE("Selected " << particles.size() << " particles to be processed by " << lastSimulator->name());
252 // Run the simulation
253 if (auto* g4sim = dynamic_cast<ISF::BaseSimulatorG4Tool*>(lastSimulator)) {
254 ATH_CHECK(g4sim->simulateVector(ctx, particles, newSecondaries,
255 outputTruth.ptr(), hitCollections,
256 shadowTruth.get()));
257 } else {
258 ATH_CHECK(lastSimulator->simulateVector(ctx, particles, newSecondaries,
259 outputTruth.ptr(),
260 shadowTruth.get()));
261 }
262
263 ATH_MSG_VERBOSE(lastSimulator->name() << " returned " << newSecondaries.size() << " new particles to be added to the queue." );
264 // Register returned particles with the entry layer tool, set their order and enqueue them
265 for ( auto* secondary : newSecondaries ) {
266
267 // Set up particle in ISF
268 m_entryLayerTool->registerParticle( *secondary );
269 if ( m_forceGeoIDSvc || !validAtlasRegion( secondary->nextGeoID() ) ) {
270 m_geoIDSvc->identifyAndRegNextGeoID( *secondary );
271 }
272
273 if ( m_orderingTool.empty() ) {
274 // Without a defined ordering, preserve old FIFO behaviour
275 secondary->setOrder( -particleQueue.size() );
276 } else {
277 m_orderingTool->setOrder( *secondary );
278 }
279
280 particleQueue.push( secondary );
281 }
282 newSecondaries.clear();
283
284 // Delete simulated particles
285 for ( auto usedParticle : particles ) {
286 delete usedParticle;
287 }
288 particles.clear();
289 }
290 ATH_MSG_VERBOSE("Final status: queue contains " << particleQueue.size() << " particles.");
291
292 // Release the event from all simulators (TODO: make the tools do this)
293 for (auto& curSimTool: m_simulationTools) {
294 if ( curSimTool ) {
295 if (auto* g4sim =
296 dynamic_cast<ISF::BaseSimulatorG4Tool*>(curSimTool.get())) {
297 ATH_CHECK(g4sim->releaseEvent(ctx, *hitCollections));
298 } else {
299 ATH_CHECK(curSimTool->releaseEvent(ctx));
300 }
301 ATH_MSG_DEBUG( "releaseEvent() completed for " << curSimTool->name() );
302 }
303 }
304
305 // Remove QS patch if required
306 if(!m_qspatcher.empty()) {
307 for (HepMC::GenEvent* currentGenEvent : *outputTruth ) {
308 ATH_CHECK(m_qspatcher->removeWorkaround(*currentGenEvent));
309 }
310 }
311
312 return StatusCode::SUCCESS;
313}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_FATAL(x)
#define ATH_MSG_INFO(x)
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_DEBUG(x)
bool validAtlasRegion(AtlasDetDescr::AtlasRegion region)
Check a given AtlasRegion for its validity.
Definition AtlasRegion.h:40
SG::WriteHandleKey< McEventCollection > m_outputTruthKey
Output Simulation Truth collection.
Definition SimKernelMT.h:91
ToolHandleArray< ISimulatorTool > m_simulationTools
Simulation Tools.
SG::ReadHandleKey< McEventCollection > m_inputEvgenKey
Input Generator Truth collection.
Definition SimKernelMT.h:89
BooleanProperty m_useShadowEvent
Definition SimKernelMT.h:97
SG::WriteHandleKey< TrackRecordCollection > m_muonEntryLayerKey
Definition SimKernelMT.h:94
ToolHandle< IGenEventFilter > m_truthPreselectionTool
ToolHandle< IParticleOrderingTool > m_orderingTool
Tool to set particle ordering.
Gaudi::Property< size_t > m_maxParticleVectorSize
Number of particles simultaneously sent to simulator.
ISimulatorTool & identifySimulator(const ISF::ISFParticle &particle)
Returns the simulator to use for the given particle.
ServiceHandle< IGeoIDSvc > m_geoIDSvc
Service to set particle GeoIDs.
SG::WriteHandleKey< TrackRecordCollection > m_caloEntryLayerKey
Output TrackRecordCollections.
Definition SimKernelMT.h:93
ServiceHandle< Simulation::IZeroLifetimePatcher > m_qspatcher
Quasi-Stable Particle Simulation Patcher.
ServiceHandle< IInputConverter > m_inputConverter
Input converter service (from Generator->ISF particle types)
SG::WriteHandleKey< TrackRecordCollection > m_muonExitLayerKey
Definition SimKernelMT.h:95
Gaudi::Property< bool > m_forceGeoIDSvc
Force geoID recalculation for each particle.
Definition SimKernelMT.h:99
PublicToolHandle< IEntryLayerTool > m_entryLayerTool
AthenaTool responsible for writing Calo/Muon Entry/Exit Layer collection.
ServiceHandle< ITruthSvc > m_truthRecordSvc
Central truth service.
int maxGeneratedVertexBarcode(const HepMC::GenEvent *genEvent)
Get the maximal absolute value of barcode of vertex present in the event. Returns a negative number.
constexpr int UNDEFINED_ID
int maxGeneratedParticleBarcode(const HepMC::GenEvent *genEvent)
Get the maximal value of barcode of particle present in the event.
std::priority_queue< ISF::ISFParticle *, ISF::ISFParticleVector, ISF::ISFParticleOrdering > ISFParticleOrderedQueue
the actual particle priority_queue
@ fAtlasMuonExit
Definition EntryLayer.h:39
@ fAtlasCaloEntry
Definition EntryLayer.h:37
@ fAtlasMuonEntry
Definition EntryLayer.h:38
std::list< ISF::ISFParticle * > ISFParticleContainer
generic ISFParticle container (not necessarily a std::list!)
std::vector< ISF::ISFParticle * > ISFParticleVector
ISFParticle vector.
constexpr ParticleHypothesis particle[PARTICLEHYPOTHESES]
the array of masses

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< 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 & AthAlgorithm::extraOutputDeps ( ) const
overridevirtualinherited

Return the list of extra output dependencies.

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

Definition at line 50 of file AthAlgorithm.cxx.

51{
52 // If we didn't find any symlinks to add, just return the collection
53 // from the base class. Otherwise, return the extended collection.
54 if (!m_extendedExtraObjects.empty()) {
56 }
57 return Algorithm::extraOutputDeps();
58}

◆ finalize()

StatusCode ISF::SimKernelMT::finalize ( )
override

Implements empty finalize (implementation required by AthAlgorithm)

Definition at line 316 of file SimKernelMT.cxx.

316 {
317 return StatusCode::SUCCESS;
318}

◆ identifySimulator()

ISF::ISimulatorTool & ISF::SimKernelMT::identifySimulator ( const ISF::ISFParticle & particle)
private

Returns the simulator to use for the given particle.

Definition at line 322 of file SimKernelMT.cxx.

322 {
323 AtlasDetDescr::AtlasRegion geoID = particle.nextGeoID();
324
325 auto& localSelectors = m_simSelectors[geoID];
326 for (auto& selector: localSelectors) {
327 bool selected = selector->selfSelect(particle);
328 if (selected) {
329 return *m_simToolMap.at(selector->simFlavor());
330 }
331 }
332
333 ATH_MSG_WARNING("No simulator found for particle (" << particle << ")."
334 << " Will send it to " << m_particleKillerTool->name());
335 return *m_particleKillerTool;
336}
#define ATH_MSG_WARNING(x)
std::map< ISF::SimulationFlavor, ISimulatorTool * > m_simToolMap
Map of the simulation flavours used in this job to the corresponding Simulation Services.
ISimulatorTool * m_particleKillerTool
When no appropriate simulator can be found for a given particle, the particle is sent to this "partic...
std::array< PublicToolHandleArray< ISimulationSelector >, AtlasDetDescr::fNumAtlasRegions > m_simSelectors
The simulation selectors defining the "routing chain".
AtlasRegion
A simple enum of ATLAS regions and sub-detectors.
Definition AtlasRegion.h:21

◆ initialize()

StatusCode ISF::SimKernelMT::initialize ( )
override

Check user configuration and configure the algorithm state accordingly.

Core Athena algorithm for the Integrated Simulation Framework.

Author
Elmar Ritsch Elmar.nosp@m..Rit.nosp@m.sch@c.nosp@m.ern..nosp@m.ch
Date
October 2016

Definition at line 27 of file SimKernelMT.cxx.

27 {
28
29 ATH_CHECK( m_simulationTools.retrieve() );
30 for ( auto& curSimTool : m_simulationTools )
31 {
32 if ( curSimTool )
33 {
34 const auto flavor = curSimTool->simFlavor();
35 auto itr = m_simToolMap.find(flavor);
36 if (itr != m_simToolMap.end() )
37 {
38 ATH_MSG_FATAL("Two ISimulatorTool instances (" << itr->second->name() << "," << curSimTool->name() << ") with the same flavor in this job!\n Check your configuration!");
39 return StatusCode::FAILURE;
40 }
41 // New flavour add it to the map.
42 m_simToolMap[flavor] = &*curSimTool;
43 if (flavor == ISF::ParticleKiller) {
44 m_particleKillerTool = &*curSimTool;
45 }
46 }
47 }
48 // Check that a particle killer simulator tool was provided
50 ATH_MSG_FATAL("No fallback ParticleKiller Simulator Tool provided in SimulationTools, the job will bail out now.");
51 return StatusCode::FAILURE;
52 }
53
54 ATH_MSG_INFO("The following Simulators will be used in this job: \t" << m_simulationTools);
55 // retrieve simulation selectors (i.e. the "routing chain")
56 for ( auto& selectorsToolHandleArray: m_simSelectors )
57 {
58 ATH_CHECK( selectorsToolHandleArray.retrieve() );
59 }
60
61 // info screen output
62 ATH_MSG_INFO( "The following routing chains are defined:" );
64 {
65 auto& localSelectors = m_simSelectors[geoID];
66 ATH_MSG_INFO( AtlasDetDescr::AtlasRegionHelper::getName(geoID) << " (GeoID=" << geoID << "): \t" << localSelectors );
67
68 for ( auto& selector : localSelectors )
69 {
70 const auto flavor = selector->simFlavor();
71 auto itr = m_simToolMap.find(flavor);
72 if (itr == m_simToolMap.end() )
73 {
74 ATH_MSG_WARNING( "Map from SimulationFlavor to SimulatorTool:" );
75 for ( auto& entry : m_simToolMap )
76 {
77 ATH_MSG_WARNING( "SimulationFlavor: " << entry.first << ", SimulatorTool Name: " << entry.second->name() );
78 }
79 ATH_MSG_FATAL( "No SimulationTool with flavor " << flavor << " expected by " << selector->name() << " found in this job!\n Check your configuration!" );
80 return StatusCode::FAILURE;
81 }
82 }
83 }
84
85 ATH_CHECK( m_entryLayerTool.retrieve() );
86
87 if ( not m_orderingTool.empty() ) {
88 ATH_CHECK( m_orderingTool.retrieve() );
89 }
90
91 ATH_CHECK( m_inputEvgenKey.initialize() );
92 ATH_CHECK( m_outputTruthKey.initialize() );
93
94 ATH_CHECK( m_caloEntryLayerKey.initialize() );
95 ATH_CHECK( m_muonEntryLayerKey.initialize() );
96 ATH_CHECK( m_muonExitLayerKey.initialize() );
97
98 ATH_CHECK( m_inputConverter.retrieve() );
99 if ( not m_truthPreselectionTool.empty() ) {
101 }
102
103 ATH_CHECK ( m_truthRecordSvc.retrieve() );
104
105 if ( not m_qspatcher.empty() ) {
106 ATH_CHECK( m_qspatcher.retrieve() );
107 }
108
109 ATH_CHECK( m_geoIDSvc.retrieve() );
110
111 return StatusCode::SUCCESS;
112}
static const char * getName(int region)
@ ParticleKiller

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< 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.

◆ isClonable()

bool ISF::SimKernelMT::isClonable ( ) const
inlineoverride

Definition at line 82 of file SimKernelMT.h.

82{ return true; }

◆ msg()

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

Definition at line 24 of file AthCommonMsg.h.

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

◆ msgLvl()

bool AthCommonMsg< 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< 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.

◆ 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< 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< Algorithm > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ sysInitialize()

StatusCode AthAlgorithm::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< Algorithm > >.

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

Definition at line 66 of file AthAlgorithm.cxx.

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

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< 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.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< 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 }
std::vector< SG::VarHandleKeyArray * > m_vhka

◆ ISFTesting::SimKernelMT_test

friend class ISFTesting::SimKernelMT_test
friend

Allow the test class access to all methods.

Definition at line 64 of file SimKernelMT.h.

Member Data Documentation

◆ m_caloEntryLayerKey

SG::WriteHandleKey<TrackRecordCollection> ISF::SimKernelMT::m_caloEntryLayerKey {this, "CaloEntryLayerKey", "CaloEntryLayer", ""}
private

Output TrackRecordCollections.

Definition at line 93 of file SimKernelMT.h.

93{this, "CaloEntryLayerKey", "CaloEntryLayer", ""};

◆ m_detStore

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

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_entryLayerTool

PublicToolHandle<IEntryLayerTool> ISF::SimKernelMT::m_entryLayerTool {this, "EntryLayerTool", "ISF_EntryLayerToolMT", ""}
private

AthenaTool responsible for writing Calo/Muon Entry/Exit Layer collection.

Definition at line 119 of file SimKernelMT.h.

119{this, "EntryLayerTool", "ISF_EntryLayerToolMT", ""};

◆ m_evtStore

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

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_extendedExtraObjects

DataObjIDColl AthAlgorithm::m_extendedExtraObjects
privateinherited

Definition at line 79 of file AthAlgorithm.h.

◆ m_forceGeoIDSvc

Gaudi::Property<bool> ISF::SimKernelMT::m_forceGeoIDSvc {this, "AlwaysUseGeoIDSvc", false, "Force geoID recalculation for each particle" }
private

Force geoID recalculation for each particle.

Definition at line 99 of file SimKernelMT.h.

99{this, "AlwaysUseGeoIDSvc", false, "Force geoID recalculation for each particle" };

◆ m_geoIDSvc

ServiceHandle<IGeoIDSvc> ISF::SimKernelMT::m_geoIDSvc {this, "GeoIDSvc", "", "Service to set particle GeoIDs"}
private

Service to set particle GeoIDs.

Definition at line 122 of file SimKernelMT.h.

122{this, "GeoIDSvc", "", "Service to set particle GeoIDs"};

◆ m_inputConverter

ServiceHandle<IInputConverter> ISF::SimKernelMT::m_inputConverter {this, "InputConverter", "", "Input McEventCollection->ISFParticleContainer conversion service."}
private

Input converter service (from Generator->ISF particle types)

Definition at line 104 of file SimKernelMT.h.

104{this, "InputConverter", "", "Input McEventCollection->ISFParticleContainer conversion service."};

◆ m_inputEvgenKey

SG::ReadHandleKey<McEventCollection> ISF::SimKernelMT::m_inputEvgenKey {this, "InputEvgenCollection", "", "Input EVGEN collection."}
private

Input Generator Truth collection.

Definition at line 89 of file SimKernelMT.h.

89{this, "InputEvgenCollection", "", "Input EVGEN collection."};

◆ m_maxParticleVectorSize

Gaudi::Property<size_t> ISF::SimKernelMT::m_maxParticleVectorSize {this, "MaximumParticleVectorSize", 10240}
private

Number of particles simultaneously sent to simulator.

Definition at line 141 of file SimKernelMT.h.

141{this, "MaximumParticleVectorSize", 10240};

◆ m_muonEntryLayerKey

SG::WriteHandleKey<TrackRecordCollection> ISF::SimKernelMT::m_muonEntryLayerKey {this, "MuonEntryLayerKey", "MuonEntryLayer", ""}
private

Definition at line 94 of file SimKernelMT.h.

94{this, "MuonEntryLayerKey", "MuonEntryLayer", ""};

◆ m_muonExitLayerKey

SG::WriteHandleKey<TrackRecordCollection> ISF::SimKernelMT::m_muonExitLayerKey {this, "MuonExitLayerKey", "MuonExitLayer", ""}
private

Definition at line 95 of file SimKernelMT.h.

95{this, "MuonExitLayerKey", "MuonExitLayer", ""};

◆ m_orderingTool

ToolHandle<IParticleOrderingTool> ISF::SimKernelMT::m_orderingTool {this, "ParticleOrderingTool", "", "Tool to set order of particles"}
private

Tool to set particle ordering.

Definition at line 125 of file SimKernelMT.h.

125{this, "ParticleOrderingTool", "", "Tool to set order of particles"};

◆ m_outputTruthKey

SG::WriteHandleKey<McEventCollection> ISF::SimKernelMT::m_outputTruthKey {this, "OutputTruthCollection", "", "Output Truth collection."}
private

Output Simulation Truth collection.

Definition at line 91 of file SimKernelMT.h.

91{this, "OutputTruthCollection", "", "Output Truth collection."};

◆ m_particleKillerTool

ISimulatorTool* ISF::SimKernelMT::m_particleKillerTool {}
private

When no appropriate simulator can be found for a given particle, the particle is sent to this "particle killer":

Definition at line 116 of file SimKernelMT.h.

116{};

◆ m_qspatcher

ServiceHandle<Simulation::IZeroLifetimePatcher> ISF::SimKernelMT::m_qspatcher {this, "QuasiStablePatcher", "", "Quasi-Stable Particle Simulation Patcher"}
private

Quasi-Stable Particle Simulation Patcher.

Definition at line 110 of file SimKernelMT.h.

110{this, "QuasiStablePatcher", "", "Quasi-Stable Particle Simulation Patcher"};

◆ m_simSelectors

std::array<PublicToolHandleArray<ISimulationSelector>, AtlasDetDescr::fNumAtlasRegions> ISF::SimKernelMT::m_simSelectors
private
Initial value:
{{
{},
{this, "IDSimulationSelectors", {} },
{this, "BeamPipeSimulationSelectors", {} },
{this, "CaloSimulationSelectors", {} },
{this, "MSSimulationSelectors", {} },
{this, "CavernSimulationSelectors", {} }
}}

The simulation selectors defining the "routing chain".

Definition at line 128 of file SimKernelMT.h.

128 {{
129 {}, // fUndefinedAtlasRegion
130 {this, "IDSimulationSelectors", {} }, // fAtlasID
131 {this, "BeamPipeSimulationSelectors", {} }, // fAtlasFoward
132 {this, "CaloSimulationSelectors", {} }, // fAtlasCalo
133 {this, "MSSimulationSelectors", {} }, // fAtlasMS
134 {this, "CavernSimulationSelectors", {} } // fAtlasCavern
135 }};

◆ m_simToolMap

std::map<ISF::SimulationFlavor, ISimulatorTool*> ISF::SimKernelMT::m_simToolMap
private

Map of the simulation flavours used in this job to the corresponding Simulation Services.

Definition at line 138 of file SimKernelMT.h.

◆ m_simulationTools

ToolHandleArray<ISimulatorTool> ISF::SimKernelMT::m_simulationTools {this, "SimulationTools", {}, ""}
private

Simulation Tools.

Definition at line 113 of file SimKernelMT.h.

113{this, "SimulationTools", {}, ""};

◆ m_truthPreselectionTool

ToolHandle<IGenEventFilter> ISF::SimKernelMT::m_truthPreselectionTool {this, "TruthPreselectionTool", "", "Tool for filtering out quasi-stable particle daughters"}
private

Definition at line 101 of file SimKernelMT.h.

101{this, "TruthPreselectionTool", "", "Tool for filtering out quasi-stable particle daughters"};

◆ m_truthRecordSvc

ServiceHandle<ITruthSvc> ISF::SimKernelMT::m_truthRecordSvc {this,"TruthRecordService", "ISF_MC15aPlusTruthService", ""}
private

Central truth service.

Definition at line 107 of file SimKernelMT.h.

107{this,"TruthRecordService", "ISF_MC15aPlusTruthService", ""};

◆ m_useShadowEvent

BooleanProperty ISF::SimKernelMT::m_useShadowEvent {this, "UseShadowEvent", false, "New approach to selecting particles for simulation" }
private

Definition at line 97 of file SimKernelMT.h.

97{this, "UseShadowEvent", false, "New approach to selecting particles for simulation" };

◆ m_varHandleArraysDeclared

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

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

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

Definition at line 398 of file AthCommonDataStore.h.


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