BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSOutputROOTEventAperture.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 "BDSOutputROOTEventAperture.hh"
20#include "BDSOutputROOTParticleData.hh"
21
22#ifndef __ROOTBUILD__
23#include "CLHEP/Units/SystemOfUnits.h"
24#include "BDSHitApertureImpact.hh"
25#endif
26
28
29BDSOutputROOTParticleData* BDSOutputROOTEventAperture::particleTable = nullptr;
30
31BDSOutputROOTEventAperture::BDSOutputROOTEventAperture()
32{
33 FlushLocal();
34}
35
36BDSOutputROOTEventAperture::~BDSOutputROOTEventAperture()
37{;}
38
39#ifndef __ROOTBUILD__
40void BDSOutputROOTEventAperture::Fill(const BDSHitApertureImpact* hit,
41 G4bool isPrimaryFirstImpact)
42{
43 n++;
44 energy.push_back((float)hit->totalEnergy / CLHEP::GeV);
45 S.push_back((double)hit->S / CLHEP::m);
46 weight.push_back((float)hit->weight);
47 isPrimary.push_back(hit->parentID == 0);
48 firstPrimaryImpact.push_back(isPrimaryFirstImpact);
49 partID.push_back((int)hit->partID);
50 turn.push_back((int)hit->turnsTaken);
51 x.push_back((float)hit->x / CLHEP::m);
52 y.push_back((float)hit->y / CLHEP::m);
53 xp.push_back((float)hit->xp / CLHEP::rad);
54 yp.push_back((float)hit->yp / CLHEP::rad);
55 T.push_back((float)hit->globalTime / CLHEP::ns);
56 kineticEnergy.push_back((float)hit->preStepKineticEnergy / CLHEP::GeV);
57 G4int pid = hit->partID;
58
59 G4bool isIonL = false; // 'L' for local variable
60 if (particleTable)
61 {isIonL = particleTable->IsIon(pid);}
62 if (particleTable && isIonL) // avoid nested ifs with duplicated setting of variables by doing it this way
63 {
64 auto& ionInfo = particleTable->GetIonInfo(pid);
65 isIon.push_back(true);
66 ionA.push_back((int)ionInfo.a);
67 ionZ.push_back((int)ionInfo.z);
68 nElectrons.push_back(0);
69 }
70 else
71 {
72 G4bool isIonAgain = hit->nElectrons > 0;
73 isIon.push_back(isIonAgain);
74 ionA.push_back(0);
75 ionZ.push_back(0);
76 nElectrons.push_back(hit->nElectrons);
77 }
78
79 trackID.push_back((int)hit->trackID);
80 parentID.push_back((int)hit->parentID);
81 modelID.push_back((int)hit->beamlineIndex);
82}
83
84#endif
85
86void BDSOutputROOTEventAperture::Fill(const BDSOutputROOTEventAperture* other)
87{
88 if (!other)
89 {return;}
90 n = other->n;
91 energy = other->energy;
92 S = other->S;
93 weight = other->weight;
94 isPrimary = other->isPrimary;
96 partID = other->partID;
97 turn = other->turn;
98 x = other->x;
99 y = other->y;
100 xp = other->xp;
101 yp = other->yp;
102 T = other->T;
104 isIon = other->isIon;
105 ionA = other->ionA;
106 ionZ = other->ionZ;
107 nElectrons = other->nElectrons;
108 trackID = other->trackID;
109 parentID = other->parentID;
110 modelID = other->modelID;
111}
112
114{
115 n = 0;
116 energy.clear();
117 S.clear();
118 weight.clear();
119 isPrimary.clear();
120 firstPrimaryImpact.clear();
121 partID.clear();
122 turn.clear();
123 x.clear();
124 y.clear();
125 xp.clear();
126 yp.clear();
127 T.clear();
128 kineticEnergy.clear();
129 isIon.clear();
130 ionA.clear();
131 ionZ.clear();
132 nElectrons.clear();
133 trackID.clear();
134 parentID.clear();
135 modelID.clear();
136}
Snapshot of information for particle passing through a collimator.
Data stored for energy deposition hits per event.
std::vector< float > kineticEnergy
Kinetic energy in GeV at pre step point.
std::vector< int > parentID
ParentID that created the deposit.
std::vector< int > turn
Turn number.
std::vector< float > energy
Total energy of particle.
std::vector< int > trackID
TrackID that created the deposit.
std::vector< bool > firstPrimaryImpact
Whether the first time the primary is passing through.
std::vector< int > modelID
Geometry model index.
std::vector< float > weight
Weight associated with loss.
std::vector< float > T
Global time (time since beginning of event).
std::vector< double > S
Global curvilinear S coordinate.
std::vector< int > partID
ParticleID that create the deposit.
Geant4 particle data for particles used in simulation.
bool IsIon(const int &pdgID) const
This doesn't count a proton (even with electrons) as an ion.