| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 #ifndef PRIMARYGENERATORACTION_HH 26 #ifndef PRIMARYGENERATORACTION_HH
27 #define PRIMARYGENERATORACTION_HH 27 #define PRIMARYGENERATORACTION_HH
28 28
29 #include "G4String.hh" 29 #include "G4String.hh"
30 #include "G4VUserPrimaryGeneratorAction.hh" 30 #include "G4VUserPrimaryGeneratorAction.hh"
31 #include "CLHEP/Units/SystemOfUnits.h" 31 #include "CLHEP/Units/SystemOfUnits.h"
32 32
33 class G4ParticleGun; 33 class G4ParticleGun;
34 class G4Event; 34 class G4Event;
35 class G4Box; 35 class G4Box;
36 class PrimaryGeneratorMessenger; 36 class PrimaryGeneratorMessenger;
37 37
38 #ifdef WITHROOT 38 #ifdef WITHROOT
39 #include <RtypesCore.h> 39 #include <RtypesCore.h>
40 class TFile; 40 class TFile;
41 class TTreeReader; 41 class TTreeReader;
42 template <typename T> class TTreeReaderValue; 42 template <typename T> class TTreeReaderValue;
43 #endif 43 #endif
44 44
45 /** 45 /**
46 * @brief Primary generator 46 * @brief Primary generator
47 * 47 *
48 * Primary generator action. 48 * Primary generator action.
49 * 49 *
50 * By default particle gun is used. It can be 50 * By default particle gun is used. It can be controlled with standard UI
51 * commands (/gun/) and with additional ones i 51 * commands (/gun/) and with additional ones introduced by the messenger.
52 * "/HGCalTestbeam/generator/momentumSpread <V 52 * "/HGCalTestbeam/generator/momentumSpread <VALUE>" to change constant
53 * particle energy to Gaussian distribution wi 53 * particle energy to Gaussian distribution with sigma expressed in units of the
54 * initial energy (e.g. value 0.05 means sigma 54 * initial energy (e.g. value 0.05 means sigma of 0.05 * E).
55 * By default it equals to 0 and constant ener 55 * By default it equals to 0 and constant energy value is used.
56 * "/HGCalTestbeam/generator/beamSpread <none/ 56 * "/HGCalTestbeam/generator/beamSpread <none/Gaussian/flat>" to define type of
57 * beam position spread. By default none is us 57 * beam position spread. By default none is used.
58 * "/HGCalTestbeam/generator/beamSpreadX <SIZE 58 * "/HGCalTestbeam/generator/beamSpreadX <SIZE>" to define size of beam spread
59 * along x axis. It is sigma of a Gaussian dis 59 * along x axis. It is sigma of a Gaussian distribution, or half-width of a
60 * flat distribution. 60 * flat distribution.
61 * "/HGCalTestbeam/generator/beamSpreadY <SIZE 61 * "/HGCalTestbeam/generator/beamSpreadY <SIZE>" to define size of beam spread
62 * along y axis. It is sigma of a Gaussian dis 62 * along y axis. It is sigma of a Gaussian distribution, or half-width of a
63 * flat distribution. 63 * flat distribution.
64 * "/HGCalTestbeam/generator/fBeamZ0 <POSITION 64 * "/HGCalTestbeam/generator/fBeamZ0 <POSITION>" to define beam position along z
65 * axis. By default edge of the world volume i 65 * axis. By default edge of the world volume is used.
66 * 66 *
67 * If installation was done with ROOT package 67 * If installation was done with ROOT package (CMake was able to locate it),
68 * an additional option of input read from the 68 * an additional option of input read from the ROOT file is enabled.
69 * It can be activated with "/HGCalTestbeam/ge 69 * It can be activated with "/HGCalTestbeam/generator/fReadInputFile true".
70 * "/HGCalTestbeam/generator/fPathInputFile <F 70 * "/HGCalTestbeam/generator/fPathInputFile <FILE>" sets the path to the input
71 * file. 71 * file.
72 * "/HGCalTestbeam/generator/startFromEvent <N 72 * "/HGCalTestbeam/generator/startFromEvent <N>" allows to start simulation from
73 * Nth event. 73 * Nth event.
74 * Please note that in current implementation 74 * Please note that in current implementation input from file needs to be
75 * executed in a non-multithreaded mode (or wi 75 * executed in a non-multithreaded mode (or with 1 thread).
76 * 76 *
77 */ 77 */
78 78
79 class PrimaryGeneratorAction : public G4VUserP 79 class PrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction {
80 public: 80 public:
81 PrimaryGeneratorAction(); 81 PrimaryGeneratorAction();
82 virtual ~PrimaryGeneratorAction(); 82 virtual ~PrimaryGeneratorAction();
83 83
84 virtual void GeneratePrimaries(G4Event *); 84 virtual void GeneratePrimaries(G4Event *);
85 85
86 const G4ParticleGun *GetParticleGun() const 86 const G4ParticleGun *GetParticleGun() const { return fParticleGun; }
87 87
88 #ifdef WITHROOT 88 #ifdef WITHROOT
89 /// Open input file with list of particles 89 /// Open input file with list of particles
90 void OpenInput(); 90 void OpenInput();
91 91
92 /// Set flag indicating that particles shoul 92 /// Set flag indicating that particles should be read from file
93 inline void SetIfUseInputFiles(G4bool aUseIn 93 inline void SetIfUseInputFiles(G4bool aUseInputFiles) {
94 fReadInputFile = aUseInputFiles; 94 fReadInputFile = aUseInputFiles;
95 }; 95 };
96 /// Set the flag indicating that particles s 96 /// Set the flag indicating that particles should be read from file
97 inline G4bool GetIfUseInputFiles() { return 97 inline G4bool GetIfUseInputFiles() { return fReadInputFile; }
98 /// Set the path to the input file 98 /// Set the path to the input file
99 inline void SetInputFiles(G4String aInputFil 99 inline void SetInputFiles(G4String aInputFiles) {
100 fPathInputFile = aInputFiles; 100 fPathInputFile = aInputFiles;
101 }; 101 };
102 /// Get the path to the input file 102 /// Get the path to the input file
103 inline G4String GetInputFiles() const { retu 103 inline G4String GetInputFiles() const { return fPathInputFile; }
104 /// Set ID of the first event to be read for 104 /// Set ID of the first event to be read for the simulation
105 inline void SetStartFromEvent(G4int aStartFr 105 inline void SetStartFromEvent(G4int aStartFromEvent) {
106 fStartFromEvent = aStartFromEvent; 106 fStartFromEvent = aStartFromEvent;
107 }; 107 };
108 /// Get ID of the first event to be read for 108 /// Get ID of the first event to be read for the simulation
109 inline G4int GetStartFromEvent() const { ret 109 inline G4int GetStartFromEvent() const { return fStartFromEvent; }
110 #endif 110 #endif
111 /// Set sigma of the Gaussian distribution f 111 /// Set sigma of the Gaussian distribution for the momentum spread
112 /// @param[in] aMomentumSpread sigma of Gaus 112 /// @param[in] aMomentumSpread sigma of Gaussian distribution expressed in
113 /// units of initial energy (e.g. 0.05 means 113 /// units of initial energy (e.g. 0.05 means sigma = 0.05 * E)
114 inline void SetMomentumSpread(G4double aMome 114 inline void SetMomentumSpread(G4double aMomentumSpread) {
115 fMomentumGaussianSpread = aMomentumSpread; 115 fMomentumGaussianSpread = aMomentumSpread;
116 }; 116 };
117 /// Get sigma of the Gaussian distribution f 117 /// Get sigma of the Gaussian distribution for the momentum spread
118 inline G4double GetMomentumSpread() const { 118 inline G4double GetMomentumSpread() const { return fMomentumGaussianSpread; }
119 /// Set type of beam position spread 119 /// Set type of beam position spread
120 /// @param[in] aType Type of beam position s 120 /// @param[in] aType Type of beam position spread: "none", "Gaussian" or
121 /// "flat". By default "none" is used. 121 /// "flat". By default "none" is used.
122 inline void SetBeamSpreadType(G4String aType 122 inline void SetBeamSpreadType(G4String aType) {
123 if (aType == "none") 123 if (aType == "none")
124 fBeamType = eNone; 124 fBeamType = eNone;
125 if (aType == "Gaussian") 125 if (aType == "Gaussian")
126 fBeamType = eGaussian; 126 fBeamType = eGaussian;
127 if (aType == "flat") 127 if (aType == "flat")
128 fBeamType = eFlat; 128 fBeamType = eFlat;
129 } 129 }
130 /// Get type of beam position spread 130 /// Get type of beam position spread
131 inline G4String GetBeamSpreadType() const { 131 inline G4String GetBeamSpreadType() const {
132 switch (fBeamType) { 132 switch (fBeamType) {
133 case eNone: 133 case eNone:
134 return "none"; 134 return "none";
135 break; 135 break;
136 case eGaussian: 136 case eGaussian:
137 return "Gaussian"; 137 return "Gaussian";
138 break; 138 break;
139 case eFlat: 139 case eFlat:
140 return "flat"; 140 return "flat";
141 break; 141 break;
142 } 142 }
143 return ""; 143 return "";
144 } 144 }
145 /// Set size of beam position spread along X 145 /// Set size of beam position spread along X axis
146 /// @param[in] aBeamSpreadX Size of beam pos 146 /// @param[in] aBeamSpreadX Size of beam position spread. Sigma for Gaussian
147 /// distribution or half-width of flat distr 147 /// distribution or half-width of flat distribution.
148 inline void SetBeamSpreadX(G4double aBeamSpr 148 inline void SetBeamSpreadX(G4double aBeamSpreadX) {
149 fSigmaBeamX = aBeamSpreadX; 149 fSigmaBeamX = aBeamSpreadX;
150 } 150 }
151 /// Get size of beam position spread along X 151 /// Get size of beam position spread along X axis
152 inline G4double GetBeamSpreadX() const { ret 152 inline G4double GetBeamSpreadX() const { return fSigmaBeamX; }
153 /// Set size of beam position spread along Y 153 /// Set size of beam position spread along Y axis
154 /// @param[in] aBeamSpreadY Size of beam pos 154 /// @param[in] aBeamSpreadY Size of beam position spread. Sigma for Gaussian
155 /// distribution or half-width of flat distr 155 /// distribution or half-width of flat distribution.
156 inline void SetBeamSpreadY(G4double aBeamSpr 156 inline void SetBeamSpreadY(G4double aBeamSpreadY) {
157 fSigmaBeamY = aBeamSpreadY; 157 fSigmaBeamY = aBeamSpreadY;
158 } 158 }
159 /// Get size of beam position spread along Y 159 /// Get size of beam position spread along Y axis
160 inline G4double GetBeamSpreadY() const { ret 160 inline G4double GetBeamSpreadY() const { return fSigmaBeamY; }
161 /// Set initial beam position along Z axis 161 /// Set initial beam position along Z axis
162 /// By default edge of world volume is used 162 /// By default edge of world volume is used
163 inline void SetBeamZ0(G4double aBeamZ0) { fB 163 inline void SetBeamZ0(G4double aBeamZ0) { fBeamZ0 = aBeamZ0; }
164 /// Get initial beam position along Z axis 164 /// Get initial beam position along Z axis
165 inline G4double GetBeamZ0() const { return f 165 inline G4double GetBeamZ0() const { return fBeamZ0; }
166 166
167 private: 167 private:
168 /// Pointer to the particle gun 168 /// Pointer to the particle gun
169 G4ParticleGun *fParticleGun; 169 G4ParticleGun *fParticleGun;
170 /// Pointer to the world volume for initial 170 /// Pointer to the world volume for initial beam position
171 G4Box *fEnvelopeBox; 171 G4Box *fEnvelopeBox;
172 /// Pointer to the messenger with custom UI 172 /// Pointer to the messenger with custom UI commands
173 PrimaryGeneratorMessenger *fMessenger; 173 PrimaryGeneratorMessenger *fMessenger;
174 /// enum describing the beam position spread 174 /// enum describing the beam position spread in transverse plane
175 enum eBeamType { eNone, eGaussian, eFlat }; 175 enum eBeamType { eNone, eGaussian, eFlat };
176 /// Type of beam position spread in transver 176 /// Type of beam position spread in transverse dimension
177 eBeamType fBeamType = eBeamType::eNone; 177 eBeamType fBeamType = eBeamType::eNone;
178 /// Size of beam position spread along X axi 178 /// Size of beam position spread along X axis
179 /// Sigma for Gaussian, and half-width for f 179 /// Sigma for Gaussian, and half-width for flat distribution
180 G4double fSigmaBeamX = 0; 180 G4double fSigmaBeamX = 0;
181 /// Size of beam position spread along Y axi 181 /// Size of beam position spread along Y axis
182 /// Sigma for Gaussian, and half-width for f 182 /// Sigma for Gaussian, and half-width for flat distribution
183 G4double fSigmaBeamY = 0; 183 G4double fSigmaBeamY = 0;
184 /// Initial beam position along Z axis 184 /// Initial beam position along Z axis
185 G4double fBeamZ0 = -999 * CLHEP::m; 185 G4double fBeamZ0 = -999 * CLHEP::m;
186 /// Sigma of Gaussian momentum spread 186 /// Sigma of Gaussian momentum spread
187 G4double fMomentumGaussianSpread = 0; 187 G4double fMomentumGaussianSpread = 0;
188 188
189 #ifdef WITHROOT 189 #ifdef WITHROOT
190 /// Flag indicating if primaries should be r 190 /// Flag indicating if primaries should be read from file instead of using
191 /// the particle gun 191 /// the particle gun
192 G4bool fReadInputFile = false; 192 G4bool fReadInputFile = false;
193 /// Path to the input file 193 /// Path to the input file
194 G4String fPathInputFile = ""; 194 G4String fPathInputFile = "";
195 /// ID of the first event in the file to be 195 /// ID of the first event in the file to be used in this simulation
196 G4int fStartFromEvent = 0; 196 G4int fStartFromEvent = 0;
197 /// Counter of event 197 /// Counter of event
198 G4int fEventCounter = -1; 198 G4int fEventCounter = -1;
199 /// Pointer to the input file 199 /// Pointer to the input file
200 TFile *fInputFile = nullptr; 200 TFile *fInputFile = nullptr;
201 /// Pointer to the tree containing particles 201 /// Pointer to the tree containing particles
202 TTreeReader *fHgcalReader = nullptr; 202 TTreeReader *fHgcalReader = nullptr;
203 /// Reader of event ID 203 /// Reader of event ID
204 TTreeReaderValue<Float_t> *fHgcalEventId; 204 TTreeReaderValue<Float_t> *fHgcalEventId;
205 /// Reader of particle PDG 205 /// Reader of particle PDG
206 TTreeReaderValue<Float_t> *fHgcalPdgId; 206 TTreeReaderValue<Float_t> *fHgcalPdgId;
207 /// Reader of particle X position (in mm) 207 /// Reader of particle X position (in mm)
208 TTreeReaderValue<Float_t> *fHgcalPosX; 208 TTreeReaderValue<Float_t> *fHgcalPosX;
209 /// Reader of particle Y position (in mm) 209 /// Reader of particle Y position (in mm)
210 TTreeReaderValue<Float_t> *fHgcalPosY; 210 TTreeReaderValue<Float_t> *fHgcalPosY;
211 /// Reader of particle Z position (in mm) 211 /// Reader of particle Z position (in mm)
212 // TTreeReaderValue<Float_t> *fHgcalPosZ; << 212 TTreeReaderValue<Float_t> *fHgcalPosZ;
213 /// Reader of particle X momentum (in MeV) 213 /// Reader of particle X momentum (in MeV)
214 TTreeReaderValue<Float_t> *fHgcalMomX; 214 TTreeReaderValue<Float_t> *fHgcalMomX;
215 /// Reader of particle Y momentum (in MeV) 215 /// Reader of particle Y momentum (in MeV)
216 TTreeReaderValue<Float_t> *fHgcalMomY; 216 TTreeReaderValue<Float_t> *fHgcalMomY;
217 /// Reader of particle Z momentum (in MeV) 217 /// Reader of particle Z momentum (in MeV)
218 TTreeReaderValue<Float_t> *fHgcalMomZ; 218 TTreeReaderValue<Float_t> *fHgcalMomZ;
219 #endif 219 #endif
220 }; 220 };
221 221
222 //....oooOO0OOooo........oooOO0OOooo........oo 222 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
223 223
224 #endif 224 #endif
225 225