Belle II Software development
Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
SteppingAction Class Reference

The Class for the stepping action. More...

#include <SteppingAction.h>

Inheritance diagram for SteppingAction:

Public Member Functions

 SteppingAction ()
 Constructor.
 
virtual ~SteppingAction ()
 Destructor.
 
void setMaxNumberSteps (int maxSteps)
 Sets the maximum number of steps before a track is stopped and killed.
 
void setStoreTrajectories (bool store)
 Sets the trajectory option to enable storing of the simulated particle trajectories.
 
void setAbsorbersR (const std::vector< float > &vec)
 Sets the radii of absorbers for killing tracks across them.
 
virtual void UserSteppingAction (const G4Step *step)
 The method will be called at each step during simulation.
 

Protected Attributes

int m_maxNumberSteps
 The maximum number of steps before the track transportation is stopped and the track is killed.
 
bool m_storeTrajectories = false
 if true, check if the track has attached trajectory info and append step information if necessary
 
std::vector< float > m_absorbers
 The absorbers defined at given radii where tracks across them will be destroyed.
 

Private Member Functions

void writeVREventStep (const G4Step *, const G4Track *)
 Method to write (almost) each G4Step to the VR event file.
 

Private Attributes

bool m_writeSimSteps {false}
 Flag for writing out the simulation steps.
 

Detailed Description

The Class for the stepping action.

This class is called for each step during the Geant4 transportation. It implements protection mechanisms to remove unreasonable tracks.

Definition at line 30 of file SteppingAction.h.

Constructor & Destructor Documentation

◆ SteppingAction()

SteppingAction ( )

Constructor.

Definition at line 27 of file SteppingAction.cc.

28{
29 //Default value for the maximum number of steps
30 m_maxNumberSteps = 100000;
31 m_writeSimSteps = Environment::Instance().getWriteSimSteps();
32}
Belle2::Environment::getWriteSimSteps
bool getWriteSimSteps() const
Get the flag for writing the simulation steps into an output csv file.
Definition Environment.h:220
Belle2::Environment::Instance
static Environment & Instance()
Static method to get a reference to the Environment instance.
Definition Environment.cc:28

◆ ~SteppingAction()

~SteppingAction ( )
virtual

Destructor.

Definition at line 35 of file SteppingAction.cc.

36{
37
38}

Member Function Documentation

◆ setAbsorbersR()

void setAbsorbersR ( const std::vector< float > & vec)
inline

Sets the radii of absorbers for killing tracks across them.

This set is used in the PXD only simulation for PXD gain calibration.

Parameters
vecThe c++ vector with the radii of absorbers in cm

Definition at line 56 of file SteppingAction.h.

56{ m_absorbers = vec; };

◆ setMaxNumberSteps()

void setMaxNumberSteps ( int maxSteps)
inline

Sets the maximum number of steps before a track is stopped and killed.

Parameters
maxStepsThe maximum number of steps.

Definition at line 47 of file SteppingAction.h.

47{ m_maxNumberSteps = maxSteps; };

◆ setStoreTrajectories()

void setStoreTrajectories ( bool store)
inline

Sets the trajectory option to enable storing of the simulated particle trajectories.

Definition at line 50 of file SteppingAction.h.

50{ m_storeTrajectories = store; }

◆ UserSteppingAction()

void UserSteppingAction ( const G4Step * step)
virtual

The method will be called at each step during simulation.

Parameters
stepThe pointer of current step.

Definition at line 41 of file SteppingAction.cc.

42{
43 G4Track* track = step->GetTrack();
44
45 //------------------------------
46 // Check for NULL world volume
47 //------------------------------
48 if (track->GetVolume() == NULL) {
49 B2WARNING("SteppingAction: Track in NULL volume, terminating!\n"
50 << "step_no=" << track->GetCurrentStepNumber() << " type=" << track->GetDefinition()->GetParticleName()
51 << "\n position=" << G4BestUnit(track->GetPosition(), "Length") << " momentum=" << G4BestUnit(track->GetMomentum(), "Energy"));
52 track->SetTrackStatus(fStopAndKill);
53 return;
54 }
55
56 //------------------------------
57 // Check for absorbers
58 //------------------------------
59 for (const auto& rAbsorber : m_absorbers) {
60 const G4ThreeVector stepPrePos = step->GetPreStepPoint()->GetPosition() / CLHEP::mm * Unit::mm;
61 const G4ThreeVector stepPostPos = step->GetPostStepPoint()->GetPosition() / CLHEP::mm * Unit::mm;
62 if (stepPrePos.perp() < (rAbsorber * Unit::cm) && stepPostPos.perp() > (rAbsorber * Unit::cm)) {
63 //B2WARNING("SteppingAction: Track across absorbers, terminating!\n"
64 //<< "step_no=" << track->GetCurrentStepNumber() << " type=" << track->GetDefinition()->GetParticleName()
65 //<< "\n position=" << G4BestUnit(track->GetPosition(), "Length") << " momentum=" << G4BestUnit(track->GetMomentum(), "Energy") << "\n PrePos.perp=" << stepPrePos.perp() << ", PostPos.perp=" << stepPostPos.perp() << " cm" );
66 track->SetTrackStatus(fStopAndKill);
67 return;
68 }
69 }
70
71 // If we are running this job for producing virtual reality events, let's run the relevant method
72 if (m_writeSimSteps) {
73 writeVREventStep(step, track);
74 }
75
76 //---------------------------------------
77 // Check for very high number of steps.
78 //---------------------------------------
79 if (track->GetCurrentStepNumber() > m_maxNumberSteps) {
80 B2WARNING("SteppingAction: Too many steps for this track, terminating!\n"
81 << "step_no=" << track->GetCurrentStepNumber() << "type=" << track->GetDefinition()->GetParticleName()
82 << "\n position=" << G4BestUnit(track->GetPosition(), "Length") << " momentum=" << G4BestUnit(track->GetMomentum(), "Energy"));
83 track->SetTrackStatus(fStopAndKill);
84 return;
85 }
86
87 //-----------------------------------------------------------
88 // Check if there is an attached trajectory. If so, fill it.
89 //-----------------------------------------------------------
90 if (m_storeTrajectories) {
91 TrackInfo* info = dynamic_cast<TrackInfo*>(track->GetUserInformation());
92 if (info && info->getTrajectory()) {
93 MCParticleTrajectory& trajectory = *(info->getTrajectory());
94 if (trajectory.empty()) {
95 const G4ThreeVector stepPos = step->GetPreStepPoint()->GetPosition() / CLHEP::mm * Unit::mm;
96 const G4ThreeVector stepMom = step->GetPreStepPoint()->GetMomentum() / CLHEP::MeV * Unit::MeV;
97 trajectory.addPoint(
98 stepPos.x(), stepPos.y(), stepPos.z(),
99 stepMom.x(), stepMom.y(), stepMom.z()
100 );
101 }
102 const G4ThreeVector stepPos = step->GetPostStepPoint()->GetPosition() / CLHEP::mm * Unit::mm;
103 const G4ThreeVector stepMom = step->GetPostStepPoint()->GetMomentum() / CLHEP::MeV * Unit::MeV;
104 trajectory.addPoint(
105 stepPos.x(), stepPos.y(), stepPos.z(),
106 stepMom.x(), stepMom.y(), stepMom.z()
107 );
108 }
109 }
110}
Belle2::MCParticleTrajectory::addPoint
void addPoint(float x, float y, float z, float px, float py, float pz)
Add a point to the trajectory.
Definition MCParticleTrajectory.h:56
Belle2::MCParticleTrajectory::empty
bool empty() const
return true if size()==0
Definition MCParticleTrajectory.h:41
Belle2::Unit::mm
static const double mm
[millimeters]
Definition Unit.h:70
Belle2::Unit::MeV
static const double MeV
[megaelectronvolt]
Definition Unit.h:114
Belle2::Unit::cm
static const double cm
Standard units with the value = 1.
Definition Unit.h:47

◆ writeVREventStep()

void writeVREventStep ( const G4Step * step,
const G4Track * track )
private

Method to write (almost) each G4Step to the VR event file.

Definition at line 113 of file SteppingAction.cc.

114{
115 // There must be an open output file (opened in simulation's EventAction)
116 RunManager& runManager = RunManager::Instance();
117 const EventAction* eventAction = static_cast<const Belle2::Simulation::EventAction*>(runManager.GetUserEventAction());
118 if (eventAction == nullptr)
119 return;
120 std::ofstream* output = eventAction->getVREventStream();
121 if (output == nullptr)
122 return;
123 if (!output->is_open())
124 return;
125 // Ignore (hyper)nuclei heavier than alphas
126 if (abs(track->GetDefinition()->GetPDGEncoding()) > 1000020040)
127 return;
128 // Limit the VR event history to 100 ns
129 G4StepPoint* postStepPoint = step->GetPostStepPoint();
130 if (postStepPoint->GetGlobalTime() > 100.0)
131 return;
132 // Discard soft particles (KE < 500 keV) unless they are optical photons (for which KE=0)
133 G4StepPoint* preStepPoint = step->GetPreStepPoint();
134 double KE = preStepPoint->GetTotalEnergy() - track->GetDefinition()->GetPDGMass();
135 if ((KE < 0.0005) && (track->GetDefinition()->GetParticleName() != "opticalphoton"))
136 return;
137
138 // We will definitely write one record to the VR event file
139 G4String pVolName = track->GetVolume()->GetName();
140 G4String sensitiveDetectorName = "";
141 if (pVolName.compare(0, 4, "PXD.") == 0) {
142 if (pVolName.find(".Active") != std::string::npos) { sensitiveDetectorName = "PXD"; }
143 } else if (pVolName.compare(0, 4, "SVD.") == 0) {
144 if (pVolName.find(".Active") != std::string::npos) { sensitiveDetectorName = "SVD"; }
145 } else if (pVolName.compare(0, 20, "physicalSD_CDCLayer_") == 0) {
146 sensitiveDetectorName = "CDC";
147 } else if (pVolName.compare(0, 19, "TOP.moduleSensitive") == 0) {
148 sensitiveDetectorName = "TOP";
149 } else if (pVolName.compare(0, 23, "av_1_impr_1_cuttest_pv_") == 0) {
150 sensitiveDetectorName = "TOP";
151 } else if (pVolName.compare(0, 12, "moduleWindow") == 0) {
152 sensitiveDetectorName = "ARICH";
153 } else if (pVolName.compare(0, 25, "ARICH.AerogelSupportPlate") == 0) {
154 sensitiveDetectorName = "ARICH";
155 } else if (pVolName.compare(0, 25, "eclBarrelCrystalPhysical_") == 0) {
156 sensitiveDetectorName = "ECL";
157 } else if (pVolName.compare(0, 22, "eclFwdCrystalPhysical_") == 0) {
158 sensitiveDetectorName = "ECL";
159 } else if (pVolName.compare(0, 22, "eclBwdCrystalPhysical_") == 0) {
160 sensitiveDetectorName = "ECL";
161 } else if (pVolName.compare(0, 20, "BKLM.ScintActiveType") == 0) {
162 sensitiveDetectorName = "BKLM";
163 } else if (pVolName.compare(0, 10, "BKLM.Layer") == 0) {
164 if (pVolName.find("GasPhysical") != std::string::npos) {
165 sensitiveDetectorName = "BKLM";
166 }
167 } else if (pVolName.compare(0, 15, "StripSensitive_") == 0) {
168 sensitiveDetectorName = "EKLM";
169 }
170// Content of each record:
171// TrackID,ParentID,ParticleName,Mass,Charge,StepNumber,Status,VolumeName,
172// MaterialName,IsFirstStepInVolume,IsLastStepInVolume,EnergyDeposit,
173// ProcessType,ProcessName,PrePointX,PrePointY,PrePointZ,PrePointT,
174// PrePointPX,PrePointPY,PrePointPZ,PrePointE,PostPointX,PostPointY,
175// PostPointZ,PostPointT,PostPointPX,PostPointPY,PostPointPZ,PostPointE
176 (*output) << std::fixed << std::setprecision(4)
177 << track->GetTrackID() << ","
178 << track->GetParentID() << ","
179 << track->GetDefinition()->GetParticleName() << ","
180 << track->GetDefinition()->GetPDGMass() << ","
181 << int(track->GetDefinition()->GetPDGCharge()) << ","
182 << track->GetCurrentStepNumber() << ","
183 << track->GetTrackStatus() << ","
184 << pVolName << ","
185 << sensitiveDetectorName << ","
186 << track->GetMaterial()->GetName() << ","
187 << step->IsFirstStepInVolume() << ","
188 << step->IsLastStepInVolume() << ","
189 << step->GetTotalEnergyDeposit() << ","
190 << postStepPoint->GetProcessDefinedStep()->GetProcessType() << ","
191 << postStepPoint->GetProcessDefinedStep()->GetProcessName() << ","
192 << preStepPoint->GetPosition().x() << ","
193 << preStepPoint->GetPosition().y() << ","
194 << preStepPoint->GetPosition().z() << ","
195 << preStepPoint->GetGlobalTime() << ","
196 << preStepPoint->GetMomentum().x() << ","
197 << preStepPoint->GetMomentum().y() << ","
198 << preStepPoint->GetMomentum().z() << ","
199 << preStepPoint->GetTotalEnergy() << ","
200 << postStepPoint->GetPosition().x() << ","
201 << postStepPoint->GetPosition().y() << ","
202 << postStepPoint->GetPosition().z() << ","
203 << postStepPoint->GetGlobalTime() << ","
204 << postStepPoint->GetMomentum().x() << ","
205 << postStepPoint->GetMomentum().y() << ","
206 << postStepPoint->GetMomentum().z() << ","
207 << postStepPoint->GetTotalEnergy() << std::endl;
208}

Member Data Documentation

◆ m_absorbers

std::vector<float> m_absorbers
protected

The absorbers defined at given radii where tracks across them will be destroyed.

Definition at line 69 of file SteppingAction.h.

◆ m_maxNumberSteps

int m_maxNumberSteps
protected

The maximum number of steps before the track transportation is stopped and the track is killed.

Definition at line 66 of file SteppingAction.h.

◆ m_storeTrajectories

bool m_storeTrajectories = false
protected

if true, check if the track has attached trajectory info and append step information if necessary

Definition at line 68 of file SteppingAction.h.

◆ m_writeSimSteps

bool m_writeSimSteps {false}
private

Flag for writing out the simulation steps.

Definition at line 77 of file SteppingAction.h.

77{false};
The documentation for this class was generated from the following files: