Geant4  10.02.p01
G4EnergyLossMessenger.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 // $Id: G4EnergyLossMessenger.cc 90095 2015-05-13 12:10:25Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 // File name: G4EnergyLossMessenger
33 //
34 // Author: Michel Maire
35 //
36 // Creation date: 17-03-2011 (original version of 22-06-2000)
37 //
38 // Modifications:
39 // 10-01-06 SetStepLimits -> SetStepFunction (V.Ivanchenko)
40 // 10-01-06 PreciseRange -> CSDARange (V.Ivanchenko)
41 // 10-05-06 Add command MscStepLimit (V.Ivanchenko)
42 // 10-10-06 Add DEDXBinning command (V.Ivanchenko)
43 // 07-02-07 Add MscLateralDisplacement command (V.Ivanchenko)
44 // 12-02-07 Add SetSkin, SetLinearLossLimit (V.Ivanchenko)
45 // 15-03-07 Send a message "/run/physicsModified" if reinitialisation
46 // is needed after the command (V.Ivanchenko)
47 // 16-03-07 modify /process/eLoss/minsubsec command (V.Ivanchenko)
48 // 18-05-07 add /process/msc directory and commands (V.Ivanchenko)
49 // 11-03-08 add /process/em directory and commands (V.Ivanchenko)
50 //
51 // -------------------------------------------------------------------
52 //
53 
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
56 
57 #include "G4EnergyLossMessenger.hh"
58 
59 #include "G4UIdirectory.hh"
60 #include "G4UIcommand.hh"
61 #include "G4UIparameter.hh"
62 #include "G4UIcmdWithABool.hh"
63 #include "G4UIcmdWithAnInteger.hh"
64 #include "G4UIcmdWithADouble.hh"
66 #include "G4UIcmdWithAString.hh"
67 #include "G4EmProcessOptions.hh"
68 #include "G4UImanager.hh"
69 #include "G4MscStepLimitType.hh"
70 #include "G4EmProcessOptions.hh"
71 
72 #include <sstream>
73 
74 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
75 
77 {
78 
79  opt = 0;
80  SubSecCmd = new G4UIcmdWithABool("/process/eLoss/subsec",this);
81  SubSecCmd->SetGuidance("Switch true/false the subcutoff generation.");
82  SubSecCmd->SetParameterName("choice",true);
85 
86  StepFuncCmd = new G4UIcommand("/process/eLoss/StepFunction",this);
87  StepFuncCmd->SetGuidance("Set the energy loss step limitation parameters.");
88  StepFuncCmd->SetGuidance(" dRoverR : max Range variation per step");
89  StepFuncCmd->SetGuidance(" finalRange: range for final step");
90 
91  G4UIparameter* dRoverRPrm = new G4UIparameter("dRoverR",'d',false);
92  dRoverRPrm->SetGuidance("max Range variation per step (fractional number)");
93  dRoverRPrm->SetParameterRange("dRoverR>0. && dRoverR<=1.");
94  StepFuncCmd->SetParameter(dRoverRPrm);
95 
96  G4UIparameter* finalRangePrm = new G4UIparameter("finalRange",'d',false);
97  finalRangePrm->SetGuidance("range for final step");
98  finalRangePrm->SetParameterRange("finalRange>0.");
99  StepFuncCmd->SetParameter(finalRangePrm);
100 
101  G4UIparameter* unitPrm = new G4UIparameter("unit",'s',true);
102  unitPrm->SetGuidance("unit of finalRange");
103  unitPrm->SetDefaultValue("mm");
104  G4String unitCandidates =
106  unitPrm->SetParameterCandidates(unitCandidates);
107 
108  StepFuncCmd->SetParameter(unitPrm);
110 
111  IntegCmd = new G4UIcmdWithABool("/process/eLoss/integral",this);
112  IntegCmd->SetGuidance("Switch true/false the integral option");
113  IntegCmd->SetParameterName("integ",true);
114  IntegCmd->SetDefaultValue(true);
116 
117  deexCmd = new G4UIcommand("/process/em/deexcitation",this);
118  deexCmd->SetGuidance("Set deexcitation flags per G4Region.");
119  deexCmd->SetGuidance(" regName : G4Region name");
120  deexCmd->SetGuidance(" flagFluo : Fluorescence");
121  deexCmd->SetGuidance(" flagAuger : Auger");
122  deexCmd->SetGuidance(" flagPIXE : PIXE");
123 
124  G4UIparameter* regName = new G4UIparameter("regName",'s',false);
125  deexCmd->SetParameter(regName);
126 
127  G4UIparameter* flagFluo = new G4UIparameter("flagFluo",'s',false);
128  deexCmd->SetParameter(flagFluo);
129 
130  G4UIparameter* flagAuger = new G4UIparameter("flagAuger",'s',false);
131  deexCmd->SetParameter(flagAuger);
132 
133  G4UIparameter* flagPIXE = new G4UIparameter("flagPIXE",'s',false);
134  deexCmd->SetParameter(flagPIXE);
135 
136  bfCmd = new G4UIcommand("/process/em/setBiasingFactor",this);
137  bfCmd->SetGuidance("Set factor for the process cross section.");
138  bfCmd->SetGuidance(" procName : process name");
139  bfCmd->SetGuidance(" procFact : factor");
140  bfCmd->SetGuidance(" flagFact : flag to change weight");
141 
142  G4UIparameter* procName = new G4UIparameter("procName",'s',false);
143  bfCmd->SetParameter(procName);
144 
145  G4UIparameter* procFact = new G4UIparameter("procFact",'d',false);
146  bfCmd->SetParameter(procFact);
147 
148  G4UIparameter* flagFact = new G4UIparameter("flagFact",'s',false);
149  bfCmd->SetParameter(flagFact);
151 
152  fiCmd = new G4UIcommand("/process/em/setForcedInteraction",this);
153  fiCmd->SetGuidance("Set factor for the process cross section.");
154  fiCmd->SetGuidance(" procNam : process name");
155  fiCmd->SetGuidance(" regNam : region name");
156  fiCmd->SetGuidance(" tlength : fixed target length");
157  fiCmd->SetGuidance(" tflag : flag to change weight");
158 
159  G4UIparameter* procNam = new G4UIparameter("procNam",'s',false);
160  fiCmd->SetParameter(procNam);
161 
162  G4UIparameter* regNam = new G4UIparameter("regNam",'s',false);
163  fiCmd->SetParameter(regNam);
164 
165  G4UIparameter* tlength = new G4UIparameter("tlength",'d',false);
166  fiCmd->SetParameter(tlength);
167 
168  G4UIparameter* unitT = new G4UIparameter("unitT",'s',true);
169  fiCmd->SetParameter(unitT);
170  unitT->SetGuidance("unit of tlength");
171 
172  G4UIparameter* flagT = new G4UIparameter("tflag",'s',true);
173  fiCmd->SetParameter(flagT);
175 
176  brCmd = new G4UIcommand("/process/em/setSecBiasing",this);
177  brCmd->SetGuidance("Set bremsstrahlung or delta-e- splitting/Russian roullette per region.");
178  brCmd->SetGuidance(" bProcNam : process name");
179  brCmd->SetGuidance(" bRegNam : region name");
180  brCmd->SetGuidance(" bFactor : number of splitted gamma or probability of Russian roulette");
181  brCmd->SetGuidance(" bEnergy : max energy of a secondary for this biasing method");
182 
183  G4UIparameter* bProcNam = new G4UIparameter("bProcNam",'s',false);
184  brCmd->SetParameter(bProcNam);
185 
186  G4UIparameter* bRegNam = new G4UIparameter("bRegNam",'s',false);
187  brCmd->SetParameter(bRegNam);
188 
189  G4UIparameter* bFactor = new G4UIparameter("bFactor",'d',false);
190  brCmd->SetParameter(bFactor);
191 
192  G4UIparameter* bEnergy = new G4UIparameter("bEnergy",'d',false);
193  brCmd->SetParameter(bEnergy);
194 
195  G4UIparameter* bUnit = new G4UIparameter("bUnit",'s',true);
196  brCmd->SetParameter(bUnit);
197  brCmd->SetGuidance("unit of energy");
198 
200 }
201 
202 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
203 
205 {
206  delete opt;
207  delete SubSecCmd;
208  delete StepFuncCmd;
209  delete IntegCmd;
210  delete bfCmd;
211  delete fiCmd;
212  delete brCmd;
213 }
214 
215 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
216 
218 {
219  if(!opt) { opt = new G4EmProcessOptions(); }
220 
221  if(command == SubSecCmd) {
223  G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified");
224  } else if (command == StepFuncCmd) {
225  G4double v1,v2;
226  G4String unt;
227  std::istringstream is(newValue);
228  is >> v1 >> v2 >> unt;
229  v2 *= G4UIcommand::ValueOf(unt);
230  opt->SetStepFunction(v1,v2);
231  } else if (command == deexCmd) {
232  G4String s1 (""), s2(""), s3(""), s4("");
233  G4bool b2(false), b3(false), b4(false);
234  std::istringstream is(newValue);
235  is >> s1 >> s2 >> s3 >> s4;
236  if(s2 == "true") { b2 = true; }
237  if(s3 == "true") { b3 = true; }
238  if(s4 == "true") { b4 = true; }
240  G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified");
241  } else if (command == IntegCmd) {
243  } else if (command == bfCmd) {
244  G4double v1(1.0);
245  G4String s0(""),s1("");
246  std::istringstream is(newValue);
247  is >> s0 >> v1 >> s1;
248  G4bool yes = false;
249  if(s1 == "true") { yes = true; }
250  opt->SetProcessBiasingFactor(s0,v1,yes);
251  G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified");
252  } else if (command == fiCmd) {
253  G4double v1(0.0);
254  G4String s1(""),s2(""),s3(""),unt("mm");
255  std::istringstream is(newValue);
256  is >> s1 >> s2 >> v1 >> unt >> s3;
257  G4bool yes = false;
258  if(s3 == "true") { yes = true; }
259  v1 *= G4UIcommand::ValueOf(unt);
260  opt->ActivateForcedInteraction(s1,v1,s2,yes);
261  G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified");
262  } else if (command == brCmd) {
263  G4double fb(1.0),en(1.e+30);
264  G4String s1(""),s2(""),unt("MeV");
265  std::istringstream is(newValue);
266  is >> s1 >> s2 >> fb >> en >> unt;
267  en *= G4UIcommand::ValueOf(unt);
268  if (s1=="phot"||s1=="compt"||s1=="conv")
270  else opt->ActivateSecondaryBiasing(s1,s2,fb,en);
271  G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified");
272  }
273 }
274 
275 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
void SetParameter(G4UIparameter *const newParameter)
Definition: G4UIcommand.hh:152
G4UIcmdWithABool * IntegCmd
void SetParameterRange(const char *theRange)
void SetParameterCandidates(const char *theString)
void SetNewValue(G4UIcommand *, G4String)
void SetDefaultValue(const char *theDefaultValue)
void SetStepFunction(G4double v1, G4double v2)
void ActivateSecondaryBiasingForGamma(const G4String &name, const G4String &region, G4double factor, G4double energyLimit)
static G4bool GetNewBoolValue(const char *paramString)
void SetDefaultValue(G4bool defVal)
static G4UImanager * GetUIpointer()
Definition: G4UImanager.cc:58
void SetParameterName(const char *theName, G4bool omittable, G4bool currentAsDefault=false)
static G4String UnitsList(const char *unitCategory)
Definition: G4UIcommand.cc:320
static const G4double b3
bool G4bool
Definition: G4Types.hh:79
void SetGuidance(const char *aGuidance)
Definition: G4UIcommand.hh:161
static const G4double b2
void AvailableForStates(G4ApplicationState s1)
Definition: G4UIcommand.cc:239
void SetDeexcitationActiveRegion(const G4String &rname="", G4bool valDeexcitation=true, G4bool valAuger=true, G4bool valPIXE=true)
static G4double ValueOf(const char *unitName)
Definition: G4UIcommand.cc:308
G4EmProcessOptions * opt
G4UIcmdWithABool * SubSecCmd
void ActivateSecondaryBiasing(const G4String &name, const G4String &region, G4double factor, G4double energyLimit)
void ActivateForcedInteraction(const G4String &name, G4double length=0.0, const G4String &region="", G4bool flag=true)
void SetIntegral(G4bool val)
double G4double
Definition: G4Types.hh:76
void SetGuidance(const char *theGuidance)
void SetSubCutoff(G4bool val, const G4Region *r=0)
static const G4double b4
static G4String CategoryOf(const char *unitName)
Definition: G4UIcommand.cc:315
void SetProcessBiasingFactor(const G4String &name, G4double val, G4bool flag=true)
G4int ApplyCommand(const char *aCommand)
Definition: G4UImanager.cc:446