BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBunchSixTrackLink.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 "BDSBunchSixTrackLink.hh"
20#include "BDSDebug.hh"
21#include "BDSException.hh"
22#include "BDSIonDefinition.hh"
23#include "BDSParticleDefinition.hh"
24
25#include "globals.hh"
26#include "G4IonTable.hh"
27#include "G4ParticleDefinition.hh"
28#include "G4ParticleTable.hh"
29#include "G4String.hh"
30
31#include <vector>
32
33BDSBunchSixTrackLink::BDSBunchSixTrackLink():
34 BDSBunch("sixtracklink"),
35 currentIndex(0),
36 currentExternalParticleID(0),
37 currentExternalParentID(0),
38 currentParticleDefinition(nullptr),
39 size(0)
40{;}
41
42BDSBunchSixTrackLink::~BDSBunchSixTrackLink()
43{;}
44
46{
47 if (currentIndex >= size)
48 {
49 G4cout << __METHOD_NAME__ << "looping to start of bunch" << G4endl;
50 currentIndex = 0;
51 }
52
53 G4int ci = currentIndex;
54 currentIndex++;
55
56 auto particle = particles[ci];
57 currentParticleDefinition = particle->particleDefinition;
59 //UpdateGeant4ParticleDefinition(particleDefinition->PDGID()); // TBC
61
62 currentExternalParticleID = particle->externalParticleID;
63 currentExternalParentID = particle->externalParentID;
64
65 return particle->coords;
66}
67
69 const BDSParticleCoordsFull& coordsIn,
70 int externalParticleID,
71 int externalParentID)
72{
73 particles.emplace_back(new BDSParticleExternal(particleDefinitionIn, coordsIn, externalParticleID, externalParentID));
74 size = (G4int)particles.size();
75}
76
78{
79 currentIndex = 0;
80 size = 0;
81 for (auto p : particles)
82 {delete p;}
83 particles.clear();
84}
85
87{
88 G4ParticleDefinition* newParticleDefinition = nullptr;
90 {
91 G4ParticleTable::G4PTblEncodingDictionary* encoding = G4ParticleTable::fEncodingDictionary;
92 auto search = encoding->find(pdgID);
93 if (search != encoding->end())
94 {newParticleDefinition = search->second;}
95 else
96 {throw BDSException(__METHOD_NAME__,"PDG ID \"" + std::to_string(pdgID) + "not found in particle table");}
97 }
98 else
99 {
100 G4IonTable* ionTable = G4ParticleTable::GetParticleTable()->GetIonTable();
102 newParticleDefinition = ionTable->GetIon(ionDefinition->Z(),
103 ionDefinition->A(),
104 ionDefinition->ExcitationEnergy());
105 }
106 particleDefinition->UpdateG4ParticleDefinition(newParticleDefinition);
107 // Note we don't need to take care of electrons here. These are automatically
108 // allocated by Geant4 when it converts the primary vertex to a dynamic particle
109 // (in the process of constructing a track from it) (done in G4PrimaryTransformer)
110 // this relies on the charge being set correctly - Geant4 detects this isn't the same
111 // as Z and adds electrons accordingly.
112#if G4VERSION_NUMBER > 1049
113 // in the case of ions the particle definition is only available now
114 // fix the looping thresholds now it's available
115 BDS::FixGeant105ThreshholdsForParticle(newParticleDefinition);
116#endif
117}
118
120{
121 if (!currentParticleDefinition->IsAnIon())
122 {return;}
123
124 G4IonTable* ionTable = G4ParticleTable::GetParticleTable()->GetIonTable();
125 BDSIonDefinition* ionDefinition = currentParticleDefinition->IonDefinition();
126 G4ParticleDefinition* ionParticleDef = ionTable->GetIon(ionDefinition->Z(),
127 ionDefinition->A(),
128 ionDefinition->ExcitationEnergy());
129 currentParticleDefinition->UpdateG4ParticleDefinition(ionParticleDef);
130 // Note we don't need to take care of electrons here. These are automatically
131 // allocated by Geant4 when it converts the primary vertex to a dynamic particle
132 // (in the process of constructing a track from it) (done in G4PrimaryTransformer)
133 // this relies on the charge being set correctly - Geant4 detects this isn't the same
134 // as Z and adds electrons accordingly.
135#if G4VERSION_NUMBER > 1049
136 // in the case of ions the particle definition is only available now
137 // fix the looping thresholds now it's available
139#endif
140}
The base class for bunch distribution generators.
Definition: BDSBunch.hh:47
G4bool particleDefinitionHasBeenUpdated
Definition: BDSBunch.hh:203
BDSParticleDefinition * particleDefinition
Particle definition for bunch - this class owns it.
Definition: BDSBunch.hh:199
General exception with possible name of object and message.
Definition: BDSException.hh:35
Class to parse an ion particle definition.
G4int A() const
Accessor.
G4double ExcitationEnergy() const
Accessor.
G4int Z() const
Accessor.
A set of particle coordinates including energy and weight.
Wrapper for particle definition.
BDSIonDefinition * IonDefinition() const
Accessor.
G4bool IsAnIon() const
Accessor.
void UpdateG4ParticleDefinition(G4ParticleDefinition *particleIn)
A set of particle coordinates from an external interface.
void FixGeant105ThreshholdsForParticle(const G4ParticleDefinition *particleDefinition)