BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSOutputROOTEventModel.hh
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
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 BDSOUTPUTROOTEVENTMODEL_H
20#define BDSOUTPUTROOTEVENTMODEL_H
21
22#include "BDSOutputROOTEventCollimatorInfo.hh"
23
24#include <map>
25#include <string>
26#include <vector>
27
28#include "Rtypes.h"
29#include "TObject.h"
30#include "TVector3.h"
31#include "TRotation.h"
32
33#ifndef __ROOTBUILD__
34#include "G4RotationMatrix.hh"
35#include "G4String.hh"
36#include "G4Transform3D.hh"
37#include "G4Types.hh"
38#endif
39
46class BDSOutputROOTEventModel: public TObject
47{
48public:
53
55 void Flush();
56
58 int findNearestElement(const TVector3& point) const;
59
60#ifndef __ROOTBUILD__
62 explicit BDSOutputROOTEventModel(G4bool storeCollimatorInfoIn);
63
65 TRotation ConvertToROOT(const G4RotationMatrix* rm) const;
66 TRotation ConvertToROOT(const G4RotationMatrix& rm) const;
67 TVector3 ConvertToROOT(const G4ThreeVector& v) const;
68
70 virtual void Fill(const std::vector<G4int>& collimatorIndicesIn = {},
71 const std::map<G4String, G4int>& collimatorIndicesByNameIn = {},
72 const std::vector<BDSOutputROOTEventCollimatorInfo>& collimatorInfoIn = {},
73 const std::vector<G4String>& collimatorBranchNamesIn = {},
74 const std::map<G4String, G4Transform3D>* scorerMeshPlacements = nullptr,
75 const std::map<short int, G4String>* materialIDToNameUnique = nullptr,
76 G4bool storeTrajectory = false);
77#endif
78
79 int n;
80 std::vector<std::string> componentName;
81 std::vector<std::string> placementName;
82 std::vector<std::string> componentType;
83 std::vector<float> length;
84 std::vector<float> angle;
85 std::vector<TVector3> staPos;
86 std::vector<TVector3> midPos;
87 std::vector<TVector3> endPos;
88 std::vector<TRotation> staRot;
89 std::vector<TRotation> midRot;
90 std::vector<TRotation> endRot;
91 std::vector<TVector3> staRefPos;
92 std::vector<TVector3> midRefPos;
93 std::vector<TVector3> endRefPos;
94 std::vector<TRotation> staRefRot;
95 std::vector<TRotation> midRefRot;
96 std::vector<TRotation> endRefRot;
97 std::vector<float> tilt;
98 std::vector<float> offsetX;
99 std::vector<float> offsetY;
100 std::vector<float> staS;
101 std::vector<float> midS;
102 std::vector<float> endS;
103 std::vector<std::string> beamPipeType;
104 std::vector<double> beamPipeAper1;
105 std::vector<double> beamPipeAper2;
106 std::vector<double> beamPipeAper3;
107 std::vector<double> beamPipeAper4;
108 std::vector<std::string> material;
109 std::vector<float> k1;
110 std::vector<float> k2;
111 std::vector<float> k3;
112 std::vector<float> k4;
113 std::vector<float> k5;
114 std::vector<float> k6;
115 std::vector<float> k7;
116 std::vector<float> k8;
117 std::vector<float> k9;
118 std::vector<float> k10;
119 std::vector<float> k11;
120 std::vector<float> k12;
121 std::vector<float> k1s;
122 std::vector<float> k2s;
123 std::vector<float> k3s;
124 std::vector<float> k4s;
125 std::vector<float> k5s;
126 std::vector<float> k6s;
127 std::vector<float> k7s;
128 std::vector<float> k8s;
129 std::vector<float> k9s;
130 std::vector<float> k10s;
131 std::vector<float> k11s;
132 std::vector<float> k12s;
133 std::vector<float> ks;
134 std::vector<float> hkick;
135 std::vector<float> vkick;
136 std::vector<float> bField;
137 std::vector<float> eField;
138 std::vector<float> e1;
139 std::vector<float> e2;
140 std::vector<float> hgap;
141 std::vector<float> fint;
142 std::vector<float> fintx;
143 std::vector<float> fintk2;
144 std::vector<float> fintxk2;
145
148
151 std::vector<int> collimatorIndices;
152
154 std::map<std::string, int> collimatorIndicesByName;
155
157 std::vector<BDSOutputROOTEventCollimatorInfo> collimatorInfo;
158
160 std::vector<std::string> collimatorBranchNamesUnique;
161
162 std::map<std::string, TVector3> scoringMeshTranslation;
163 std::map<std::string, TRotation> scoringMeshRotation;
164 std::vector<std::string> scoringMeshName;
165
166 std::map<short int, std::string> materialIDToName;
167 std::map<std::string, short int> materialNameToID;
168
169 std::vector<std::string> samplerNamesUnique;
170 std::vector<double> samplerSPosition;
171 std::vector<std::string> samplerCNamesUnique;
172 std::vector<std::string> samplerSNamesUnique;
173 std::map<std::string, double> samplerCRadius;
174 std::map<std::string, double> samplerSRadius;
175
176 ClassDef(BDSOutputROOTEventModel, 6);
177};
178
179#endif
Information stored per model representing accelerator.
std::vector< BDSOutputROOTEventCollimatorInfo > collimatorInfo
Collimator information explicitly.
bool storeCollimatorInfo
Whether optional collimator information was stored.
std::vector< float > vkick
Vertical fractional momentum kick.
std::vector< float > eField
E field in V/m.
std::vector< std::string > collimatorBranchNamesUnique
Vector of all collimator branch names in event tree used to load data.
int nCollimators
Number of collimators in beam line.
virtual ~BDSOutputROOTEventModel()
Destructor.
int findNearestElement(const TVector3 &point) const
Find element index closest to point.
std::vector< float > bField
B field in T.
std::vector< float > hkick
Horizontal fractional momentum kick.
std::vector< int > collimatorIndices
std::vector< std::string > material
Material associated with element if any.
virtual void Fill(const std::vector< G4int > &collimatorIndicesIn={}, const std::map< G4String, G4int > &collimatorIndicesByNameIn={}, const std::vector< BDSOutputROOTEventCollimatorInfo > &collimatorInfoIn={}, const std::vector< G4String > &collimatorBranchNamesIn={}, const std::map< G4String, G4Transform3D > *scorerMeshPlacements=nullptr, const std::map< short int, G4String > *materialIDToNameUnique=nullptr, G4bool storeTrajectory=false)
Fill root output.
void Flush()
Initialise all members.
BDSOutputROOTEventModel()
Default constructor.
std::map< std::string, int > collimatorIndicesByName
Similar cache but by name of collimator as built by BDSIM.
std::vector< float > ks
Solenoid strength.
TRotation ConvertToROOT(const G4RotationMatrix *rm) const
Utility function.