Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4ParticleGunMessenger Class Reference

#include <G4ParticleGunMessenger.hh>

+ Inheritance diagram for G4ParticleGunMessenger:

Public Member Functions

 G4ParticleGunMessenger (G4ParticleGun *fPtclGun)
 
 ~G4ParticleGunMessenger ()
 
void SetNewValue (G4UIcommand *command, G4String newValues)
 
G4String GetCurrentValue (G4UIcommand *command)
 
- Public Member Functions inherited from G4UImessenger
 G4UImessenger ()
 
 G4UImessenger (const G4String &path, const G4String &dsc)
 
virtual ~G4UImessenger ()
 
virtual G4String GetCurrentValue (G4UIcommand *command)
 
virtual void SetNewValue (G4UIcommand *command, G4String newValue)
 
G4bool operator== (const G4UImessenger &messenger) const
 

Additional Inherited Members

- Protected Member Functions inherited from G4UImessenger
G4String ItoS (G4int i)
 
G4String DtoS (G4double a)
 
G4String BtoS (G4bool b)
 
G4int StoI (G4String s)
 
G4double StoD (G4String s)
 
G4bool StoB (G4String s)
 
void AddUIcommand (G4UIcommand *newCommand)
 
void CreateDirectory (const G4String &path, const G4String &dsc)
 
template<typename T >
T * CreateCommand (const G4String &cname, const G4String &dsc)
 
- Protected Attributes inherited from G4UImessenger
G4UIdirectorybaseDir
 
G4String baseDirName
 

Detailed Description

Definition at line 54 of file G4ParticleGunMessenger.hh.

Constructor & Destructor Documentation

◆ G4ParticleGunMessenger()

G4ParticleGunMessenger::G4ParticleGunMessenger ( G4ParticleGun fPtclGun)

Definition at line 47 of file G4ParticleGunMessenger.cc.

48 :fParticleGun(fPtclGun),fShootIon(false)
49{
50 particleTable = G4ParticleTable::GetParticleTable();
51
52 gunDirectory = new G4UIdirectory("/gun/");
53 gunDirectory->SetGuidance("Particle Gun control commands.");
54
55 listCmd = new G4UIcmdWithoutParameter("/gun/List",this);
56 listCmd->SetGuidance("List available particles.");
57 listCmd->SetGuidance(" Invoke G4ParticleTable.");
58
59 particleCmd = new G4UIcmdWithAString("/gun/particle",this);
60 particleCmd->SetGuidance("Set particle to be generated.");
61 particleCmd->SetGuidance(" (geantino is default)");
62 particleCmd->SetGuidance(" (ion can be specified for shooting ions)");
63 particleCmd->SetParameterName("particleName",true);
64 particleCmd->SetDefaultValue("geantino");
65 G4String candidateList;
66 G4int nPtcl = particleTable->entries();
67 for(G4int i=0;i<nPtcl;i++)
68 {
69 G4ParticleDefinition* pd = particleTable->GetParticle(i);
70 if( !(pd->IsShortLived()) || pd->GetDecayTable() )
71 {
72 candidateList += pd->GetParticleName();
73 candidateList += " ";
74 }
75 }
76 candidateList += "ion ";
77 particleCmd->SetCandidates(candidateList);
78
79 directionCmd = new G4UIcmdWith3Vector("/gun/direction",this);
80 directionCmd->SetGuidance("Set momentum direction.");
81 directionCmd->SetGuidance("Direction needs not to be a unit vector.");
82 directionCmd->SetParameterName("ex","ey","ez",true,true);
83 directionCmd->SetRange("ex != 0 || ey != 0 || ez != 0");
84
85 energyCmd = new G4UIcmdWithADoubleAndUnit("/gun/energy",this);
86 energyCmd->SetGuidance("Set kinetic energy.");
87 energyCmd->SetParameterName("Energy",true,true);
88 energyCmd->SetDefaultUnit("GeV");
89 //energyCmd->SetUnitCategory("Energy");
90 //energyCmd->SetUnitCandidates("eV keV MeV GeV TeV");
91
92 momCmd = new G4UIcmdWith3VectorAndUnit("/gun/momentum",this);
93 momCmd->SetGuidance("Set momentum. This command is equivalent to two commands /gun/direction and /gun/momentumAmp");
94 momCmd->SetParameterName("px","py","pz",true,true);
95 momCmd->SetRange("px != 0 || py != 0 || pz != 0");
96 momCmd->SetDefaultUnit("GeV");
97
98 momAmpCmd = new G4UIcmdWithADoubleAndUnit("/gun/momentumAmp",this);
99 momAmpCmd->SetGuidance("Set absolute value of momentum.");
100 momAmpCmd->SetGuidance("Direction should be set by /gun/direction command.");
101 momAmpCmd->SetGuidance("This command should be used alternatively with /gun/energy.");
102 momAmpCmd->SetParameterName("Momentum",true,true);
103 momAmpCmd->SetDefaultUnit("GeV");
104
105 positionCmd = new G4UIcmdWith3VectorAndUnit("/gun/position",this);
106 positionCmd->SetGuidance("Set starting position of the particle.");
107 positionCmd->SetParameterName("X","Y","Z",true,true);
108 positionCmd->SetDefaultUnit("cm");
109 //positionCmd->SetUnitCategory("Length");
110 //positionCmd->SetUnitCandidates("microm mm cm m km");
111
112 timeCmd = new G4UIcmdWithADoubleAndUnit("/gun/time",this);
113 timeCmd->SetGuidance("Set initial time of the particle.");
114 timeCmd->SetParameterName("t0",true,true);
115 timeCmd->SetDefaultUnit("ns");
116 //timeCmd->SetUnitCategory("Time");
117 //timeCmd->SetUnitCandidates("ns ms s");
118
119 polCmd = new G4UIcmdWith3Vector("/gun/polarization",this);
120 polCmd->SetGuidance("Set polarization.");
121 polCmd->SetParameterName("Px","Py","Pz",true,true);
122 polCmd->SetRange("Px>=-1.&&Px<=1.&&Py>=-1.&&Py<=1.&&Pz>=-1.&&Pz<=1.");
123
124 numberCmd = new G4UIcmdWithAnInteger("/gun/number",this);
125 numberCmd->SetGuidance("Set number of particles to be generated.");
126 numberCmd->SetParameterName("N",true,true);
127 numberCmd->SetRange("N>0");
128
129 ionCmd = new G4UIcommand("/gun/ion",this);
130 ionCmd->SetGuidance("Set properties of ion to be generated.");
131 ionCmd->SetGuidance("[usage] /gun/ion Z A Q");
132 ionCmd->SetGuidance(" Z:(int) AtomicNumber");
133 ionCmd->SetGuidance(" A:(int) AtomicMass");
134 ionCmd->SetGuidance(" Q:(int) Charge of Ion (in unit of e)");
135 ionCmd->SetGuidance(" E:(double) Excitation energy (in keV)");
136
137 G4UIparameter* param;
138 param = new G4UIparameter("Z",'i',false);
139 param->SetDefaultValue("1");
140 ionCmd->SetParameter(param);
141 param = new G4UIparameter("A",'i',false);
142 param->SetDefaultValue("1");
143 ionCmd->SetParameter(param);
144 param = new G4UIparameter("Q",'i',true);
145 param->SetDefaultValue("0");
146 ionCmd->SetParameter(param);
147 param = new G4UIparameter("E",'d',true);
148 param->SetDefaultValue("0.0");
149 ionCmd->SetParameter(param);
150
151 // set initial value to G4ParticleGun
153 fParticleGun->SetParticleMomentumDirection( G4ThreeVector(1.0,0.0,0.0) );
154 fParticleGun->SetParticleEnergy( 1.0*GeV );
155 fParticleGun->SetParticlePosition(G4ThreeVector(0.0*cm, 0.0*cm, 0.0*cm));
156 fParticleGun->SetParticleTime( 0.0*ns );
157}
CLHEP::Hep3Vector G4ThreeVector
int G4int
Definition: G4Types.hh:66
static G4Geantino * Geantino()
Definition: G4Geantino.cc:87
G4DecayTable * GetDecayTable() const
const G4String & GetParticleName() const
void SetParticleMomentumDirection(G4ParticleMomentum aMomentumDirection)
void SetParticleDefinition(G4ParticleDefinition *aParticleDefinition)
void SetParticleEnergy(G4double aKineticEnergy)
G4int entries() const
static G4ParticleTable * GetParticleTable()
G4ParticleDefinition * GetParticle(G4int index)
void SetDefaultUnit(const char *defUnit)
void SetParameterName(const char *theNameX, const char *theNameY, const char *theNameZ, G4bool omittable, G4bool currentAsDefault=false)
void SetParameterName(const char *theNameX, const char *theNameY, const char *theNameZ, G4bool omittable, G4bool currentAsDefault=false)
void SetDefaultUnit(const char *defUnit)
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
void SetCandidates(const char *candidateList)
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
void SetDefaultValue(const char *defVal)
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
void SetParameter(G4UIparameter *const newParameter)
Definition: G4UIcommand.hh:147
void SetGuidance(const char *aGuidance)
Definition: G4UIcommand.hh:156
void SetRange(const char *rs)
Definition: G4UIcommand.hh:120
void SetDefaultValue(const char *theDefaultValue)
void SetParticleTime(G4double aTime)
void SetParticlePosition(G4ThreeVector aPosition)
#define ns
Definition: xmlparse.cc:597

◆ ~G4ParticleGunMessenger()

G4ParticleGunMessenger::~G4ParticleGunMessenger ( )

Definition at line 159 of file G4ParticleGunMessenger.cc.

160{
161 delete listCmd;
162 delete particleCmd;
163 delete directionCmd;
164 delete energyCmd;
165 delete momCmd;
166 delete momAmpCmd;
167 delete positionCmd;
168 delete timeCmd;
169 delete polCmd;
170 delete numberCmd;
171 delete ionCmd;
172 delete gunDirectory;
173}

Member Function Documentation

◆ GetCurrentValue()

G4String G4ParticleGunMessenger::GetCurrentValue ( G4UIcommand command)
virtual

Reimplemented from G4UImessenger.

Definition at line 210 of file G4ParticleGunMessenger.cc.

211{
212 G4String cv;
213
214 if( command==directionCmd )
215 { cv = directionCmd->ConvertToString(fParticleGun->GetParticleMomentumDirection()); }
216 else if( command==particleCmd )
217 { cv = fParticleGun->GetParticleDefinition()->GetParticleName(); }
218 else if( command==energyCmd )
219 {
220 G4double ene = fParticleGun->GetParticleEnergy();
221 if(ene == 0.)
222 { G4cerr << " G4ParticleGun: was defined in terms of momentum." << G4endl; }
223 else
224 { cv = energyCmd->ConvertToString(ene,"GeV"); }
225 }
226 else if( command==momCmd || command==momAmpCmd )
227 {
228 G4double mom = fParticleGun->GetParticleMomentum();
229 if(mom == 0.)
230 { G4cerr << " G4ParticleGun: was defined in terms of kinetic energy." << G4endl; }
231 else
232 {
233 if( command==momCmd )
234 { cv = momCmd->ConvertToString(mom*(fParticleGun->GetParticleMomentumDirection()),"GeV"); }
235 else
236 { cv = momAmpCmd->ConvertToString(mom,"GeV"); }
237 }
238 }
239 else if( command==positionCmd )
240 { cv = positionCmd->ConvertToString(fParticleGun->GetParticlePosition(),"cm"); }
241 else if( command==timeCmd )
242 { cv = timeCmd->ConvertToString(fParticleGun->GetParticleTime(),"ns"); }
243 else if( command==polCmd )
244 { cv = polCmd->ConvertToString(fParticleGun->GetParticlePolarization()); }
245 else if( command==numberCmd )
246 { cv = numberCmd->ConvertToString(fParticleGun->GetNumberOfParticles()); }
247 else if( command==ionCmd )
248 {
249 if (fShootIon) {
250 cv = ItoS(fAtomicNumber) + " " + ItoS(fAtomicMass) + " ";
251 cv += ItoS(fIonCharge);
252 } else {
253 cv = "";
254 }
255 }
256 return cv;
257}
double G4double
Definition: G4Types.hh:64
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cerr
G4ThreeVector GetParticlePolarization() const
G4ParticleMomentum GetParticleMomentumDirection() const
G4double GetParticleMomentum() const
G4ParticleDefinition * GetParticleDefinition() const
G4int GetNumberOfParticles() const
G4double GetParticleEnergy() const
static G4String ConvertToString(G4bool boolVal)
Definition: G4UIcommand.cc:349
G4String ItoS(G4int i)
G4ThreeVector GetParticlePosition()

◆ SetNewValue()

void G4ParticleGunMessenger::SetNewValue ( G4UIcommand command,
G4String  newValues 
)
virtual

Reimplemented from G4UImessenger.

Definition at line 175 of file G4ParticleGunMessenger.cc.

176{
177 if( command==listCmd )
178 { particleTable->DumpTable(); }
179 else if( command==particleCmd )
180 {
181 if (newValues =="ion") {
182 fShootIon = true;
183 } else {
184 fShootIon = false;
185 G4ParticleDefinition* pd = particleTable->FindParticle(newValues);
186 if(pd != 0)
187 { fParticleGun->SetParticleDefinition( pd ); }
188 }
189 }
190 else if( command==directionCmd )
191 { fParticleGun->SetParticleMomentumDirection(directionCmd->GetNew3VectorValue(newValues)); }
192 else if( command==energyCmd )
193 { fParticleGun->SetParticleEnergy(energyCmd->GetNewDoubleValue(newValues)); }
194 else if( command==momCmd )
195 { fParticleGun->SetParticleMomentum(momCmd->GetNew3VectorValue(newValues)); }
196 else if( command==momAmpCmd )
197 { fParticleGun->SetParticleMomentum(momAmpCmd->GetNewDoubleValue(newValues)); }
198 else if( command==positionCmd )
199 { fParticleGun->SetParticlePosition(positionCmd->GetNew3VectorValue(newValues)); }
200 else if( command==timeCmd )
201 { fParticleGun->SetParticleTime(timeCmd->GetNewDoubleValue(newValues)); }
202 else if( command==polCmd )
203 { fParticleGun->SetParticlePolarization(polCmd->GetNew3VectorValue(newValues)); }
204 else if( command==numberCmd )
205 { fParticleGun->SetNumberOfParticles(numberCmd->GetNewIntValue(newValues)); }
206 else if( command==ionCmd )
207 { IonCommand(newValues); }
208}
void SetNumberOfParticles(G4int i)
void SetParticlePolarization(G4ThreeVector aVal)
void SetParticleMomentum(G4double aMomentum)
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
void DumpTable(const G4String &particle_name="ALL")
static G4ThreeVector GetNew3VectorValue(const char *paramString)
static G4ThreeVector GetNew3VectorValue(const char *paramString)
static G4double GetNewDoubleValue(const char *paramString)
static G4int GetNewIntValue(const char *paramString)

The documentation for this class was generated from the following files: