Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4OpticalPhysicsMessenger.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 //----------------------------------------------------------------------------
28 //
29 // ClassName: G4OpticalPhysicsMessenger
30 //
31 // Author: P.Gumplinger 30.09.2009 //
32 //
33 // Modified: P.Gumplinger 29.09.2011
34 // (based on code from I. Hrivnacova)
35 //
36 //----------------------------------------------------------------------------
37 //
38 
40 #include "G4OpticalPhysics.hh"
41 
42 #include "G4UIcommand.hh"
43 #include "G4UIdirectory.hh"
44 
45 #include "G4UIcommand.hh"
46 #include "G4UIdirectory.hh"
47 #include "G4UIcmdWithABool.hh"
48 #include "G4UIcmdWithAString.hh"
49 #include "G4UIcmdWithADouble.hh"
50 #include "G4UIcmdWithAnInteger.hh"
52 #include "G4UIparameter.hh"
53 
54 #include "G4ParticleTable.hh"
55 #include "G4ProcessManager.hh"
56 #include "G4ParticleDefinition.hh"
57 
58 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
59 
61  G4OpticalPhysics* opticalPhysics)
62  : G4UImessenger(),
63  fOpticalPhysics(opticalPhysics),
64  fSelectedProcessIndex(kNoProcess),
65  fActivateProcessCmd(0),
66  fSetOpProcessVerboseCmd(0),
67  fSetCerenkovMaxPhotonsCmd(0),
68  fSetCerenkovMaxBetaChangeCmd(0),
69  fSetCerenkovStackPhotonsCmd(0),
70  fSetScintillationYieldFactorCmd(0),
71  fSetScintillationByParticleTypeCmd(0),
72  fSetScintillationTrackInfoCmd(0),
73  fSetScintillationStackPhotonsCmd(0),
74  fSetWLSTimeProfileCmd(0),
75  fSetTrackSecondariesFirstCmd(0),
76  fSetFiniteRiseTimeCmd(0),
77  fSetInvokeSDCmd(0)
78 {
79  G4bool toBeBroadcasted = false;
80  fDir = new G4UIdirectory("/process/optical/defaults/",toBeBroadcasted);
81  fDir->SetGuidance("Commands related to the optical physics simulation engine.");
82  fDir2 = new G4UIdirectory("/process/optical/",toBeBroadcasted);
83  fDir2->SetGuidance("Commands related to the optical physics simulation engine.");
84 
85  fActivateProcessCmd= new G4UIcommand("/process/optical/processActivation", this);
86  fActivateProcessCmd->SetGuidance("Activate/deactivate the specified optical process");
87  G4UIparameter* par = new G4UIparameter("proc_name",'s',false);
88  G4String candidates;
89  for ( G4int i=0; i<kNoProcess; i++ ) {
90  candidates += G4OpticalProcessName(i);
91  candidates += G4String(" ");
92  }
93  par->SetParameterCandidates(candidates);
94  par->SetGuidance("the process name");
95  fActivateProcessCmd->SetParameter(par);
96  par = new G4UIparameter("flag",'b',true);
97  par->SetDefaultValue(true);
98  par->SetGuidance("activation flag");
99  fActivateProcessCmd->SetParameter(par);
100  fActivateProcessCmd->AvailableForStates(G4State_PreInit);
101 
102 
103  fSetOpProcessVerboseCmd = new G4UIcmdWithAnInteger("/process/optical/verbose", this);
104  fSetOpProcessVerboseCmd->SetGuidance("Set default verbosity level for optical processes");
105  fSetOpProcessVerboseCmd->SetParameterName("ver", true);
106  fSetOpProcessVerboseCmd->SetDefaultValue(1);
107  fSetOpProcessVerboseCmd->SetRange("ver>=0");
108  fSetOpProcessVerboseCmd->AvailableForStates(G4State_PreInit);
109 
110 
111  fSetTrackSecondariesFirstCmd = new G4UIcommand("/process/optical/setTrackSecondariesFirst", this);
112  fSetTrackSecondariesFirstCmd->SetGuidance("Activate/deactivate tracking of secondaries before finishing their parent track");
113  par = new G4UIparameter("proc_name",'s',false);
114  par->SetParameterCandidates(candidates);
115  fSetTrackSecondariesFirstCmd->SetParameter(par);
116  par = new G4UIparameter("flag",'b',false);
117  par->SetDefaultValue(true);
118  fSetTrackSecondariesFirstCmd->SetParameter(par);
119  fSetTrackSecondariesFirstCmd->AvailableForStates(G4State_PreInit);
120 
121  //This are repetition of process specific ui commands needed by ATLICE and VMC
122  //(UI messages needed to exist in PreInit state before actual processes are instantiated)
123  fSetCerenkovMaxPhotonsCmd = new G4UIcmdWithAnInteger("/process/optical/defaults/cerenkov/setMaxPhotons", this);
124  fSetCerenkovMaxPhotonsCmd->SetGuidance("Set default maximum number of photons per step");
125  fSetCerenkovMaxPhotonsCmd->SetGuidance("Note this command is used to set the default value,");
126  fSetCerenkovMaxPhotonsCmd->SetGuidance("if process is not active command will not have effect.");
127  fSetCerenkovMaxPhotonsCmd->SetParameterName("CerenkovMaxPhotons", false);
128  fSetCerenkovMaxPhotonsCmd->SetRange("CerenkovMaxPhotons>=0");
129  fSetCerenkovMaxPhotonsCmd->AvailableForStates(G4State_PreInit);
130 
131  fSetCerenkovMaxBetaChangeCmd = new G4UIcmdWithADouble("/process/optical/defaults/cerenkov/setMaxBetaChange", this);
132  fSetCerenkovMaxBetaChangeCmd->SetGuidance("Set default maximum change of beta of parent particle per step");
133  fSetCerenkovMaxBetaChangeCmd->SetGuidance("Note this command is used to set the default value,");
134  fSetCerenkovMaxBetaChangeCmd->SetGuidance("if process is not active command will not have effect.");
135  fSetCerenkovMaxBetaChangeCmd->SetParameterName("CerenkovMaxBetaChange", false);
136  fSetCerenkovMaxBetaChangeCmd->SetRange("CerenkovMaxBetaChange>=0");
137  fSetCerenkovMaxBetaChangeCmd->AvailableForStates(G4State_PreInit);
138 
139  fSetCerenkovStackPhotonsCmd = new G4UIcmdWithABool("/process/optical/defaults/cerenkov/setStackPhotons", this);
140  fSetCerenkovStackPhotonsCmd->SetGuidance("Set default whether or not to stack secondary Cerenkov photons");
141  fSetCerenkovStackPhotonsCmd->SetGuidance("Note this command is used to set the default value,");
142  fSetCerenkovStackPhotonsCmd->SetGuidance("if process is not active command will not have effect.");
143  fSetCerenkovStackPhotonsCmd->SetParameterName("CerenkovStackPhotons", true);
144  fSetCerenkovStackPhotonsCmd->AvailableForStates(G4State_PreInit);
145 
146  fSetScintillationYieldFactorCmd = new G4UIcmdWithADouble("/process/optical/defaults/scintillation/setYieldFactor", this);
147  fSetScintillationYieldFactorCmd->SetGuidance("Set scintillation yield factor");
148  fSetScintillationYieldFactorCmd->SetGuidance("Note this command is used to set the default value,");
149  fSetScintillationYieldFactorCmd->SetGuidance("if process is not active command will not have effect.");
150  fSetScintillationYieldFactorCmd->SetParameterName("ScintillationYieldFactor", false);
151  fSetScintillationYieldFactorCmd->SetRange("ScintillationYieldFactor>=0");
152  fSetScintillationYieldFactorCmd->AvailableForStates(G4State_PreInit, G4State_Idle);
153 
154  fSetScintillationByParticleTypeCmd = new G4UIcmdWithABool("/process/optical/defaults/scintillation/setByParticleType", this);
155  fSetScintillationByParticleTypeCmd->SetGuidance("Activate/Inactivate scintillation process by particle type");
156  fSetScintillationByParticleTypeCmd->SetGuidance("Note this command is used to set the default value,");
157  fSetScintillationByParticleTypeCmd->SetGuidance("if process is not active command will not have effect.");
158  fSetScintillationByParticleTypeCmd->SetParameterName("ScintillationByParticleTypeActivation", false);
159  fSetScintillationByParticleTypeCmd->AvailableForStates(G4State_PreInit);
160 
161  fSetScintillationTrackInfoCmd = new G4UIcmdWithABool("/process/optical/defaults/scintillation/setTrackInfo", this);
162  fSetScintillationTrackInfoCmd->SetGuidance("Activate/Inactivate scintillation TrackInformation");
163  fSetScintillationTrackInfoCmd->SetParameterName("ScintillationTrackInfo", false);
164  fSetScintillationTrackInfoCmd->AvailableForStates(G4State_PreInit);
165 
166  fSetFiniteRiseTimeCmd = new G4UIcmdWithABool("/process/optical/defaults/scintillation/setFiniteRiseTime", this);
167  fSetFiniteRiseTimeCmd->SetGuidance("Set option of a finite rise-time for G4Scintillation");
168  fSetFiniteRiseTimeCmd->SetGuidance("If set, the G4Scintillation process expects the user to have set the constant material property FAST/SLOWSCINTILLATIONRISETIME");
169  fSetFiniteRiseTimeCmd->SetGuidance("Note this command is used to set the default value,");
170  fSetFiniteRiseTimeCmd->SetGuidance("if process is not active command will not have effect.");
171  fSetFiniteRiseTimeCmd->SetParameterName("FiniteRiseTime", false);
172  fSetFiniteRiseTimeCmd->AvailableForStates(G4State_PreInit);
173 
174  fSetScintillationStackPhotonsCmd = new G4UIcmdWithABool("/process/optical/defaults/scintillation/setStackPhotons", this);
175  fSetScintillationStackPhotonsCmd->SetGuidance("Set default whether or not to stack secondary Scintillation photons");
176  fSetScintillationStackPhotonsCmd->SetGuidance("Note this command is used to set the default value,");
177  fSetScintillationStackPhotonsCmd->SetGuidance("if process is not active command will not have effect.");
178  fSetScintillationStackPhotonsCmd->SetParameterName("ScintillationStackPhotons", true);
179  fSetScintillationStackPhotonsCmd->AvailableForStates(G4State_PreInit);
180 
181  fSetWLSTimeProfileCmd = new G4UIcmdWithAString("/process/optical/defaults/wls/setTimeProfile", this);
182  fSetWLSTimeProfileCmd->SetGuidance("Set the WLS time profile (delta or exponential)");
183  fSetWLSTimeProfileCmd->SetParameterName("WLSTimeProfile", false);
184  fSetWLSTimeProfileCmd->SetCandidates("delta exponential");
185  fSetWLSTimeProfileCmd->AvailableForStates(G4State_PreInit);
186 
187  fSetInvokeSDCmd = new G4UIcmdWithABool("/process/optical/defaults/boundary/setInvokeSD", this);
188  fSetInvokeSDCmd->SetGuidance("Set option for calling InvokeSD in G4OpBoundaryProcess");
189  fSetInvokeSDCmd->SetParameterName("InvokeSD", false);
190  fSetInvokeSDCmd->AvailableForStates(G4State_PreInit);
191 }
192 
194 {
195 // Destructor
196 
197  delete fDir;
198  delete fDir2;
199  delete fActivateProcessCmd;
200  delete fSetOpProcessVerboseCmd;
201  delete fSetCerenkovMaxPhotonsCmd;
202  delete fSetCerenkovMaxBetaChangeCmd;
203  delete fSetCerenkovStackPhotonsCmd;
204  delete fSetScintillationYieldFactorCmd;
205  delete fSetScintillationByParticleTypeCmd;
206  delete fSetScintillationTrackInfoCmd;
207  delete fSetScintillationStackPhotonsCmd;
208  delete fSetWLSTimeProfileCmd;
209  delete fSetTrackSecondariesFirstCmd;
210  delete fSetFiniteRiseTimeCmd;
211  delete fSetInvokeSDCmd;
212 }
213 
214 #include <iostream>
216  G4String newValue)
217 {
219  if (command == fActivateProcessCmd) {
220  std::istringstream is(newValue.data());
221  G4String pn;
222  G4String flag;
223  is >> pn >> flag;
224  if ( pn == "Cerenkov" ) {
225  fSelectedProcessIndex = kCerenkov;
226  } else if ( pn == "Scintillation" ) {
227  fSelectedProcessIndex = kScintillation;
228  } else if ( pn == "OpAbsorption" ) {
229  fSelectedProcessIndex = kAbsorption;
230  } else if ( pn == "OpRayleigh" ) {
231  fSelectedProcessIndex = kRayleigh;
232  } else if ( pn == "OpMieHG" ) {
233  fSelectedProcessIndex = kMieHG;
234  } else if ( pn == "OpBoundary" ) {
235  fSelectedProcessIndex = kBoundary;
236  } else if ( pn == "OpWLS" ) {
237  fSelectedProcessIndex = kWLS;
238  } else {
240  msg << "Not allowed process name: "<<pn<<" (UI: "<<newValue<<")";
241  G4Exception("G4OpticalPhysicsMessenger::SetNewValue(...)","Optical001",FatalException,msg);
242  }
244  fOpticalPhysics->Configure(fSelectedProcessIndex,value);
245  }
246  else if (command == fSetTrackSecondariesFirstCmd )
247  {
248  std::istringstream is(newValue.data());
249  G4String pn;
250  G4String flag;
251  is >> pn >> flag;
252  if ( pn == "Cerenkov" ) {
253  fSelectedProcessIndex = kCerenkov;
254  } else if ( pn == "Scintillation" ) {
255  fSelectedProcessIndex = kScintillation;
256  } else if ( pn == "OpAbsorption" ) {
257  fSelectedProcessIndex = kAbsorption;
258  } else if ( pn == "OpRayleigh" ) {
259  fSelectedProcessIndex = kRayleigh;
260  } else if ( pn == "OpMieHG" ) {
261  fSelectedProcessIndex = kMieHG;
262  } else if ( pn == "OpBoundary" ) {
263  fSelectedProcessIndex = kBoundary;
264  } else if ( pn == "OpWLS" ) {
265  fSelectedProcessIndex = kWLS;
266  } else {
268  msg << "Not allowed process name: "<<pn<<" (UI: "<<newValue<<")";
269  G4Exception("G4OpticalPhysicsMessenger::SetNewValue(...)","Optical001",FatalException,msg);
270  }
272  fOpticalPhysics->SetTrackSecondariesFirst(fSelectedProcessIndex,value);
273  }
274  else if (command == fSetOpProcessVerboseCmd) {
275  fOpticalPhysics->SetVerboseLevel(fSetOpProcessVerboseCmd->GetNewIntValue(newValue));
276  }
277  else if (command == fSetCerenkovMaxPhotonsCmd) {
278  fOpticalPhysics
280  fSetCerenkovMaxPhotonsCmd->GetNewIntValue(newValue));
281  }
282  else if (command == fSetCerenkovMaxBetaChangeCmd) {
283  fOpticalPhysics
285  fSetCerenkovMaxBetaChangeCmd->GetNewDoubleValue(newValue));
286  }
287  else if (command == fSetCerenkovStackPhotonsCmd) {
288  fOpticalPhysics
290  fSetCerenkovStackPhotonsCmd->GetNewBoolValue(newValue));
291  }
292  else if (command == fSetScintillationYieldFactorCmd) {
293  fOpticalPhysics
295  fSetScintillationYieldFactorCmd->GetNewDoubleValue(newValue));
296  }
297  else if (command == fSetScintillationByParticleTypeCmd) {
298  fOpticalPhysics
300  fSetScintillationByParticleTypeCmd->GetNewBoolValue(newValue));
301  }
302  else if (command == fSetScintillationTrackInfoCmd) {
303  fOpticalPhysics
305  fSetScintillationTrackInfoCmd->GetNewBoolValue(newValue));
306  }
307  else if (command == fSetFiniteRiseTimeCmd) {
308  fOpticalPhysics
310  fSetFiniteRiseTimeCmd->GetNewBoolValue(newValue));
311  }
312  else if (command == fSetScintillationStackPhotonsCmd) {
313  fOpticalPhysics
315  fSetScintillationStackPhotonsCmd->GetNewBoolValue(newValue));
316  }
317  else if (command == fSetWLSTimeProfileCmd) {
318  fOpticalPhysics->SetWLSTimeProfile(newValue);
319  }
320  else if (command == fSetInvokeSDCmd) {
321  fOpticalPhysics
322  ->SetInvokeSD(fSetInvokeSDCmd->GetNewBoolValue(newValue));
323  }
324 }
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
virtual void SetNewValue(G4UIcommand *, G4String)
void SetParameter(G4UIparameter *const newParameter)
Definition: G4UIcommand.hh:152
Number of processes, no selected process.
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
static G4int GetNewIntValue(const char *paramString)
void SetMaxBetaChangePerStep(G4double)
void SetParameterCandidates(const char *theString)
void SetInvokeSD(G4bool)
void SetDefaultValue(const char *theDefaultValue)
Scintillation process index.
Mie scattering process index.
void SetFiniteRiseTime(G4bool)
Absorption process index.
int G4int
Definition: G4Types.hh:78
static G4bool GetNewBoolValue(const char *paramString)
void SetScintillationByParticleType(G4bool)
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
void SetMaxNumPhotonsPerStep(G4int)
const XML_Char int const XML_Char * value
Definition: expat.h:331
static G4bool ConvertToBool(const char *st)
Definition: G4UIcommand.cc:437
bool G4bool
Definition: G4Types.hh:79
void Configure(G4OpticalProcessIndex, G4bool)
static G4double GetNewDoubleValue(const char *paramString)
void SetRange(const char *rs)
Definition: G4UIcommand.hh:125
void SetGuidance(const char *aGuidance)
Definition: G4UIcommand.hh:161
void AvailableForStates(G4ApplicationState s1)
Definition: G4UIcommand.cc:240
void SetVerboseLevel(G4int value)
Wave Length Shifting process index.
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
const char * data() const
G4String G4OpticalProcessName(G4int)
Return the name for a given optical process index.
Boundary process index.
void SetScintillationYieldFactor(G4double)
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
void SetScintillationStackPhotons(G4bool)
Cerenkov process index.
void SetScintillationTrackInfo(G4bool)
void SetWLSTimeProfile(G4String)
void SetTrackSecondariesFirst(G4OpticalProcessIndex, G4bool)
void SetCandidates(const char *candidateList)
void SetDefaultValue(G4int defVal)
Rayleigh scattering process index.
void SetGuidance(const char *theGuidance)
void SetCerenkovStackPhotons(G4bool)
G4OpticalPhysicsMessenger(G4OpticalPhysics *)