BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSPhysicsCutsAndLimits.cc
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#include "BDSDebug.hh"
20#include "BDSPhysicsCutsAndLimits.hh"
21#include "BDSProcessUserSpecialCutsPDG.hh"
22
23#include "globals.hh"
24#include "G4Gamma.hh"
25#include "G4Electron.hh"
26#include "G4Positron.hh"
27#include "G4Proton.hh"
28#include "G4StepLimiter.hh"
29#include "G4Types.hh"
30#include "G4UserSpecialCuts.hh"
31#include "G4Version.hh"
32
33#include <set>
34
35BDSPhysicsCutsAndLimits::BDSPhysicsCutsAndLimits(const std::set<G4int>& pdgsToExcludeFromCuts):
36 G4VPhysicsConstructor("BDSPhysicsCutsAndLimits"),
37 useParticleExclusionFromCuts(!pdgsToExcludeFromCuts.empty())
38{
39 stepLimiter = new G4StepLimiter();
40 specialCuts = new G4UserSpecialCuts();
41 bdsSpecialCuts = new BDSProcessUserSpecialCutsPDG(pdgsToExcludeFromCuts);
42 if (useParticleExclusionFromCuts)
43 {
44 G4cout << __METHOD_NAME__ << "Excluding the following particle IDs from minimumKineticEnergy, minimumRange, maximumTrackLength cuts\n{ ";
45 for (const auto& pdgID : pdgsToExcludeFromCuts)
46 {G4cout << pdgID << " ";}
47 G4cout << "}" << G4endl;
48 }
49}
50
51BDSPhysicsCutsAndLimits::~BDSPhysicsCutsAndLimits()
52{
53 delete stepLimiter;
54 delete specialCuts;
55 delete bdsSpecialCuts;
56}
57
59{
60 G4Gamma::Gamma();
61 G4Electron::Electron();
62 G4Positron::Positron();
63 G4Proton::Proton();
64}
65
67{
68 if (Activated())
69 {return;}
70
71 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
72
73 auto cutsProcess = useParticleExclusionFromCuts ? bdsSpecialCuts : specialCuts;
74
75#if G4VERSION_NUMBER > 1029
76 auto aParticleIterator = GetParticleIterator();
77#endif
78 aParticleIterator->reset();
79 while( (*aParticleIterator)())
80 {
81 G4ParticleDefinition* particle = aParticleIterator->value();
82
83 // Flag as applying range production cuts. These only ever
84 // apply in geant4 to gamma, e+- and proton
85 if((particle->GetParticleName()=="gamma")||
86 (particle->GetParticleName()=="e-")||
87 (particle->GetParticleName()=="e+")||
88 (particle->GetParticleName()=="proton"))
89 {particle->SetApplyCutsFlag(true);}
90
91 // apply general cuts processes to all particles
92 ph->RegisterProcess(stepLimiter, particle); // this is for MaxAllowedStep
93 ph->RegisterProcess(cutsProcess, particle); // this is for all other limits
94 }
95
97}
98
G4StepLimiter * stepLimiter
Step limit process for MaxAllowedStep.
G4UserSpecialCuts * specialCuts
Process for all other limits.
virtual void ConstructProcess()
Construct and attach step limiter and cuts processes.
virtual void ConstructParticle()
Construct gamma, e+- and proton - the minimum this sets limits for.
Apply regular cuts but not to a specified set of PDG IDs.
G4bool Activated() const
Get whether this instance has been activated.
Definition: BDSSingleUse.hh:37
void SetActivated()
Flag this instance as activated for later querying.
Definition: BDSSingleUse.hh:40