BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSEventAction.hh
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#ifndef BDSEVENTACTION_H
20#define BDSEVENTACTION_H
21
22#include "BDSHitEnergyDeposition.hh"
23#include "BDSHitSampler.hh"
24#include "BDSTrajectoryFilter.hh"
25
26#include "globals.hh" // geant4 types / globals
27#include "G4UserEventAction.hh"
28
29#include <bitset>
30#include <ctime>
31#include <map>
32#include <string>
33#include <vector>
34
35class BDSEventInfo;
36class BDSOutput;
38class BDSTrajectory;
40class G4Event;
41class G4PrimaryVertex;
42
47class BDSEventAction: public G4UserEventAction
48{
49public:
50 explicit BDSEventAction(BDSOutput* outputIn);
51 virtual ~BDSEventAction();
52
53 virtual void BeginOfEventAction(const G4Event*);
54 virtual void EndOfEventAction(const G4Event*);
55
56 void StoreSeedState(G4String seedState) {seedStateAtStart = seedState;}
57 G4int CurrentEventIndex() const {return currentEventIndex;}
58
60 void SetPrimaryAbsorbedInCollimator(G4bool stoppedIn) {primaryAbsorbedInCollimator = stoppedIn;}
61
63 void SetSamplerIDsForTrajectories(const std::vector<G4int>& samplerIDsIn) {trajectorySamplerID = samplerIDsIn;}
64
67
69 void RegisterPrimaryTrajectory(const BDSTrajectoryPrimary* trajectoryIn);
70
71protected:
74 G4bool verboseThisEvent,
76 BDSHitsCollectionEnergyDeposition* eCounterFullHits,
77 const std::vector<BDSHitsCollectionSampler*>& allSamplerHits,
78 G4int nChar = 50) const;
79
82 void ConnectTrajectory(std::map<BDSTrajectory*, bool>& interestingTraj,
83 BDSTrajectory* trajectoryToConnect,
84 std::map<BDSTrajectory*, std::bitset<BDS::NTrajectoryFilters> >& trajectoryFilters) const;
85
86private:
88 G4bool verboseEventBDSIM;
89 G4int verboseEventStart;
90 G4int verboseEventStop;
93 G4bool storeTrajectorySecondary;
94 G4int printModulo;
95
99 G4int eCounterID;
109 std::map<G4String, G4int> scorerCollectionIDs;
110 std::map<G4String, G4int> extraSamplerCollectionIDs;
111 std::map<G4String, G4int> extraSamplerCylinderCollectionIDs;
112 std::map<G4String, G4int> extraSamplerSphereCollectionIDs;
113
114 time_t startTime;
115 time_t stopTime;
116
117 G4double starts;
118 G4double stops;
119
120 // Note that individual calls to std::clock have no meaning, only
121 // the differences between them, and therefore this should not be
122 // written to the output.
123 std::clock_t cpuStartTime;
124
126
135 std::vector<int> trajParticleIDIntToStore;
137 std::vector<int> trajectorySamplerID;
138 std::vector<std::pair<double,double>> trajSRangeToStore;
139 std::bitset<BDS::NTrajectoryFilters> trajFiltersSet;
141
142 std::string seedStateAtStart;
143 G4int currentEventIndex;
144
148
149 long long int nTracks;
150
155 std::map<G4int, const BDSTrajectoryPrimary*> primaryTrajectoriesCache;// Cache of primary trajectories as constructed
156};
157
158#endif
159
Process information at the event level.
G4String trajParticleIDToStore
Cache of variable from global constants.
void IncrementNTracks()
Interface for tracking action to increment the number of tracks in each event.
void ConnectTrajectory(std::map< BDSTrajectory *, bool > &interestingTraj, BDSTrajectory *trajectoryToConnect, std::map< BDSTrajectory *, std::bitset< BDS::NTrajectoryFilters > > &trajectoryFilters) const
G4int eCounterVacuumID
Collection ID for the vacuum energy deposition hits.
G4bool primaryAbsorbedInCollimator
Whether primary stopped in a collimator.
G4double trajectoryCutZ
Cache of variable from global constants.
G4int collimatorCollID
Collection ID for the collimator hits.
G4int eCounterFullID
Collection ID for general energy deposition full hits.
G4int worldExitCollID
Collection ID for the world exit hits.
G4int thinThingCollID
Collection ID for the thin thing hits.
std::vector< std::pair< double, double > > trajSRangeToStore
Cache of variable from global constants.
void SetPrimaryAbsorbedInCollimator(G4bool stoppedIn)
Flag that the primary was absorbed in a collimator - can be done externally to this class.
G4bool storeTrajectory
Cache of whether to store trajectories or not.
G4bool trajectoryFilterLogicAND
Cache of variable from global constants.
G4int eCounterWorldID
Collection ID for the world energy deposition hits.
std::bitset< BDS::NTrajectoryFilters > trajFiltersSet
Cache of variable from global constants.
std::vector< int > trajParticleIDIntToStore
Cache of variable from global constants.
G4int samplerCollID_cylin
Collection ID for cylindrical sampler hits.
std::map< G4String, G4int > extraSamplerCylinderCollectionIDs
Collection IDs for extra samplers.
G4int eCounterWorldContentsID
Collection ID for the world energy deposition hits.
BDSOutput * output
Cache of output instance. Not owned by this class.
BDSTrajectoriesToStore * IdentifyTrajectoriesForStorage(const G4Event *evt, G4bool verboseThisEvent, BDSHitsCollectionEnergyDeposition *eCounterHits, BDSHitsCollectionEnergyDeposition *eCounterFullHits, const std::vector< BDSHitsCollectionSampler * > &allSamplerHits, G4int nChar=50) const
Sift through all trajectories (if any) and mark for storage.
G4int samplerCollID_plane
Collection ID for plane sampler hits.
long long int nTracks
Accumulated number of tracks for the event.
std::map< G4String, G4int > extraSamplerSphereCollectionIDs
Collection IDs for extra samplers.
time_t startTime
Time at the start of the event.
std::map< G4String, G4int > scorerCollectionIDs
Collection IDs for all scorers.
G4bool trajConnect
Cache of variable from global constants.
BDSEventInfo * eventInfo
G4int apertureCollID
Collection ID for the aperture hits.
std::vector< int > trajectorySamplerID
Cache of variable from global constants.
std::map< G4String, G4int > extraSamplerCollectionIDs
Collection IDs for extra samplers.
G4double trajectoryCutR
Cache of variable from global constants.
G4String trajParticleNameToStore
Cache of variable from global constants.
G4double stops
Precise stop time in seconds.
std::clock_t cpuStartTime
CPU time at the start of the event.
G4int trajDepth
Cache of variable from global constants.
std::map< G4int, const BDSTrajectoryPrimary * > primaryTrajectoriesCache
G4int eCounterID
Collection ID for general energy deposition hits.
void SetSamplerIDsForTrajectories(const std::vector< G4int > &samplerIDsIn)
Update the vector of sampler IDs to match for trajectories.
void RegisterPrimaryTrajectory(const BDSTrajectoryPrimary *trajectoryIn)
Append this trajectory to vector of primaries we keep to avoid sifting at the end of event.
G4int samplerCollID_sphere
Collection ID for spherical sampler hits.
G4bool storeTrajectoryAll
Store all trajectories irrespective of filters.
time_t stopTime
Time at the end of the event.
std::string seedStateAtStart
Seed state at start of the event.
G4double starts
Precise start time in seconds.
G4int eCounterTunnelID
Collection ID for the tunnel energy deposition hits.
G4double trajectoryEnergyThreshold
Cache of variable from global constants.
Interface to store event information use G4 hooks.
Definition: BDSEventInfo.hh:42
Output base class that defines interface for all output types.
Definition: BDSOutput.hh:73
Double map of trajectories to bitset of which filters matched whether to store them.
Trajectory information for only the primary.
Trajectory information from track including last scatter etc.