Geant4  10.00.p01
G4RadioactiveDecaymessenger.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 //
29 #include "G4NuclearLevelStore.hh"
30 
31 #include <sstream>
32 
34 //
36 (G4RadioactiveDecay* theRadioactiveDecayContainer1)
37 :theRadioactiveDecayContainer(theRadioactiveDecayContainer1)
38 {
39  //
40  //
41  // main directory for control of the RDM
42  //
43  //
44  grdmDirectory = new G4UIdirectory("/grdm/");
45  grdmDirectory->SetGuidance("Controls for the Radioactive Decay Module.");
46  //
47  //
48  // Command to define the limits on nucleus the RDM will treat.
49  //
50  nucleuslimitsCmd = new
51  G4UIcmdWithNucleusLimits("/grdm/nucleusLimits",this);
52  nucleuslimitsCmd->SetGuidance
53  ("Set the atomic weight and number limits for the RDM.");
54  nucleuslimitsCmd->SetParameterName("aMin","aMax","zMin","zMax",true);
55  //
56 
57  //
58  // The next command contols whether the decay will be treated analoguely or
59  // with variance reduction
60  //
61  analoguemcCmd = new G4UIcmdWithABool ("/grdm/analogueMC",this);
62  analoguemcCmd->SetGuidance("false: variance reduction method; true: analogue method");
63  analoguemcCmd->SetParameterName("AnalogueMC",true);
64  analoguemcCmd->SetDefaultValue(true);
65  //
66  // The next command contols whether beta decay will be treated faithfully or
67  // in fast mode
68  //
69  fbetaCmd = new G4UIcmdWithABool ("/grdm/fBeta",this);
70  fbetaCmd->SetGuidance("false: use 3-body decay, true: use histogram method");
71  fbetaCmd->SetParameterName("fBeta",true);
72  fbetaCmd->SetDefaultValue(false);
73 
74  //
75  //
76  // Command to selete a logical volume for RDM.
77  //
78  avolumeCmd = new
79  G4UIcmdWithAString("/grdm/selectVolume",this);
80  avolumeCmd->SetGuidance
81  ("Suppply a logical volumes name to add it to the RDM apply list");
82  avolumeCmd->SetParameterName("aVolume",false);
83  //
84  //
85  //
86  // Command to de-selete a logical volume for RDM.
87  //
88  deavolumeCmd = new
89  G4UIcmdWithAString("/grdm/deselectVolume",this);
90  deavolumeCmd->SetGuidance
91  ("Suppply a logical volumes name to remove it from the RDM apply list");
92  deavolumeCmd->SetParameterName("aVolume",false);
93  //
94  //
95  // Command to selete all logical volumes for RDM.
96  //
97  allvolumesCmd = new
98  G4UIcmdWithoutParameter("/grdm/allVolumes",this);
99  allvolumesCmd->SetGuidance
100  (" apply RDM to all logical volumes. No parameter required.");
101  // allvolumeCmd->SetParameterName("AddAVolume",true);
102 
103  //
104  // Command to de-selete a logical volume for RDM.
105  //
106  deallvolumesCmd = new
107  G4UIcmdWithoutParameter("/grdm/noVolumes",this);
108  deallvolumesCmd->SetGuidance
109  (" RDM is not applied to any logical volumes");
110 
111  // deallvolumesCmd->SetParameterName("RemoveAVolume",true);
112  //
113  // The next command contols whether the branching ratio biasing will be applied or not
114  //
115  brbiasCmd = new G4UIcmdWithABool ("/grdm/BRbias",this);
116  brbiasCmd->SetGuidance("false: no biasing; true: all branches are treated as equal");
117  brbiasCmd->SetParameterName("BRBias",true);
118  brbiasCmd->SetDefaultValue(true);
119  //
120  // Command contols whether ICM will be applied or not
121  //
122  icmCmd = new G4UIcmdWithABool ("/grdm/applyICM",this);
123  icmCmd->SetGuidance("True: ICM is applied; false: no");
124  icmCmd->SetParameterName("applyICM",true);
125  icmCmd->SetDefaultValue(true);
126  //icmCmd->AvailableForStates(G4State_PreInit);
127  //
128  // Command contols whether ARM will be applied or not
129  //
130  armCmd = new G4UIcmdWithABool ("/grdm/applyARM",this);
131  armCmd->SetGuidance("True: ARM is applied; false: no");
132  armCmd->SetParameterName("applyARM",true);
133  armCmd->SetDefaultValue(true);
134  //armCmd->AvailableForStates(G4State_PreInit);
135  //
136  // Command to set the h-l thresold for isomer production
137  //
138  hlthCmd = new G4UIcmdWithADoubleAndUnit("/grdm/hlThreshold",this);
139  hlthCmd->SetGuidance("Set the h-l threshold for isomer production");
140  hlthCmd->SetParameterName("hlThreshold",false);
141  // hlthCmd->SetRange("hlThreshold>0.");
142  hlthCmd->SetUnitCategory("Time");
143  // hlthCmd->AvailableForStates(G4State_PreInit);
144  //
145  // Command to define the incident particle source time profile.
146  //
147  sourcetimeprofileCmd = new
148  G4UIcmdWithAString("/grdm/sourceTimeProfile",this);
149  sourcetimeprofileCmd->SetGuidance
150  ("Supply the name of the ascii file containing the source particle time profile");
151  sourcetimeprofileCmd->SetParameterName("STimeProfile",true);
152  sourcetimeprofileCmd->SetDefaultValue("source.data");
153  //
154  //
155  // Command to define the incident particle source time profile.
156  //
157  decaybiasprofileCmd = new
158  G4UIcmdWithAString("/grdm/decayBiasProfile",this);
159  decaybiasprofileCmd->SetGuidance
160  ("Supply the name of the ascii file containing the decay bias time profile");
161  decaybiasprofileCmd->SetParameterName("DBiasProfile",true);
162  decaybiasprofileCmd->SetDefaultValue("bias.data");
163 
164  //
165  // Command to set the directional bias (collimation) vector
166  //
167  colldirCmd = new G4UIcmdWith3Vector("/grdm/decayDirection",this);
168  colldirCmd->SetGuidance("Supply the direction vector for decay products");
169  colldirCmd->SetParameterName("X","Y","Z",false);
170 
171  //
172  // Command to set the directional bias (collimation) half angle ("cone")
173  //
174  collangleCmd = new G4UIcmdWithADoubleAndUnit("/grdm/decayHalfAngle",this);
175  collangleCmd->SetGuidance
176  ("Supply maximum angle from direction vector for decay products");
177  collangleCmd->SetParameterName("halfAngle",false);
178  collangleCmd->SetUnitCategory("Angle");
179 
180  //
181  // This command setup the nuclei spliting parameter
182  //
183  splitnucleiCmd = new G4UIcmdWithAnInteger("/grdm/splitNuclei",this);
184  splitnucleiCmd->SetGuidance("Set number of spliting for the isotopes.");
185  splitnucleiCmd->SetParameterName("NSplit",true);
186  splitnucleiCmd->SetDefaultValue(1);
187  splitnucleiCmd->SetRange("NSplit>=1");
188 
189  //
190  // This command setup the verbose level of radioactive decay
191  //
192  verboseCmd = new G4UIcmdWithAnInteger("/grdm/verbose",this);
193  verboseCmd->SetGuidance("Set verbose level: 0, 1, 2 or 3");
194  verboseCmd->SetParameterName("VerboseLevel",true);
195  verboseCmd->SetDefaultValue(1);
196  verboseCmd->SetRange("VerboseLevel>=0");
197 
198  //
199  //This commansd allows the user to define its own decay datafile for
200  // a given isotope
201  //
202  userDecayDataCmd = new G4UIcommand("/grdm/setRadioactiveDecayFile",this);
203  G4UIparameter* Z_para= new G4UIparameter("Z_isotope",'i',true);
204  Z_para->SetParameterRange("Z_isotope > 0");
205  Z_para->SetGuidance("Z: Charge number of isotope");
206 
207 
208  G4UIparameter* A_para= new G4UIparameter("A_isotope",'i',true);
209  A_para->SetParameterRange("A_isotope > 1");
210  A_para->SetGuidance("A: mass number of isotope");
211 
212  G4UIparameter* FileName_para= new G4UIparameter("file_name",'s',true);
213  FileName_para->SetGuidance("Name of the user data file");
214  userDecayDataCmd->SetParameter(Z_para);
215  userDecayDataCmd->SetParameter(A_para);
216  userDecayDataCmd->SetParameter(FileName_para);
217 
218  //
219  //This commands allows the user to define its own evaporation data file for
220  // a given isotope
221  //
222  userEvaporationDataCmd = new G4UIcommand("/grdm/setPhotoEvaporationFile",this);
223  userEvaporationDataCmd->SetParameter(Z_para);
224  userEvaporationDataCmd->SetParameter(A_para);
225  userEvaporationDataCmd->SetParameter(FileName_para);
226 
227 
228 }
230 //
232 {
233  delete grdmDirectory;
234  delete nucleuslimitsCmd;
235  delete sourcetimeprofileCmd;
236  delete decaybiasprofileCmd;
237  delete analoguemcCmd;
238  delete fbetaCmd;
239  delete brbiasCmd;
240  delete splitnucleiCmd;
241  delete verboseCmd;
242  delete avolumeCmd;
243  delete deavolumeCmd;
244  delete allvolumesCmd;
245  delete deallvolumesCmd;
246  delete icmCmd;
247  delete armCmd;
248  delete hlthCmd;
249  delete userDecayDataCmd;
250  delete userEvaporationDataCmd;
251  delete colldirCmd;
252  delete collangleCmd;
253 
254 }
256 //
258 {
260  SetNucleusLimits(nucleuslimitsCmd->GetNewNucleusLimitsValue(newValues));}
261  else if (command==analoguemcCmd) {theRadioactiveDecayContainer->
262  SetAnalogueMonteCarlo(analoguemcCmd->GetNewBoolValue(newValues));}
263  else if (command==fbetaCmd) {theRadioactiveDecayContainer->
264  SetFBeta(fbetaCmd->GetNewBoolValue(newValues));}
265  else if (command==avolumeCmd) {theRadioactiveDecayContainer->
266  SelectAVolume(newValues);}
267  else if (command==deavolumeCmd) {theRadioactiveDecayContainer->
268  DeselectAVolume(newValues);}
269  else if (command==allvolumesCmd) {theRadioactiveDecayContainer->
270  SelectAllVolumes();}
271  else if (command==deallvolumesCmd) {theRadioactiveDecayContainer->
272  DeselectAllVolumes();}
273  else if (command==brbiasCmd) {theRadioactiveDecayContainer->
274  SetBRBias(brbiasCmd->GetNewBoolValue(newValues));}
276  SetSourceTimeProfile(newValues);}
278  SetDecayBias(newValues);}
279  else if (command==splitnucleiCmd) {theRadioactiveDecayContainer->
280  SetSplitNuclei(splitnucleiCmd->GetNewIntValue(newValues));}
281  else if (command==verboseCmd) {theRadioactiveDecayContainer->
282  SetVerboseLevel(verboseCmd->GetNewIntValue(newValues));}
283  else if (command==icmCmd ) {theRadioactiveDecayContainer->
284  SetICM(icmCmd->GetNewBoolValue(newValues));}
285  else if (command==armCmd ) {theRadioactiveDecayContainer->
286  SetARM(armCmd->GetNewBoolValue(newValues));}
287  else if (command==hlthCmd ) {theRadioactiveDecayContainer->
288  SetHLThreshold(hlthCmd->GetNewDoubleValue(newValues));}
289 
290  else if (command ==userDecayDataCmd){
291  G4int Z,A;
292  G4String file_name;
293  const char* nv = (const char*)newValues;
294  std::istringstream is(nv);
295  is >> Z>>A>>file_name;
297  }
298  else if (command ==userEvaporationDataCmd){
299  G4int Z,A;
300  G4String file_name;
301  const char* nv = (const char*)newValues;
302  std::istringstream is(nv);
303  is >> Z>>A>>file_name;
305  }
306  else if (command==colldirCmd) {theRadioactiveDecayContainer->
307  SetDecayDirection(colldirCmd->GetNew3VectorValue(newValues));}
308  else if (command==collangleCmd) {theRadioactiveDecayContainer->
309  SetDecayHalfAngle(collangleCmd->GetNewDoubleValue(newValues));}
310 }
311 
312 
313 
314 
315 
316 
static G4int GetNewIntValue(const char *paramString)
void SetParameterRange(const char *theRange)
G4UIcmdWithoutParameter * deallvolumesCmd
static G4NuclearLevelStore * GetInstance()
static G4double GetNewDoubleValue(const char *paramString)
int G4int
Definition: G4Types.hh:78
static G4bool GetNewBoolValue(const char *paramString)
void AddUserEvaporationDataFile(G4int Z, G4int A, const G4String &filename)
G4UIcmdWithoutParameter * allvolumesCmd
static G4ThreeVector GetNew3VectorValue(const char *paramString)
G4UIcmdWithADoubleAndUnit * collangleCmd
static const G4double A[nN]
void SetNewValue(G4UIcommand *command, G4String newValues)
G4NucleusLimits GetNewNucleusLimitsValue(G4String paramString)
G4RadioactiveDecay * theRadioactiveDecayContainer
G4RadioactiveDecaymessenger(G4RadioactiveDecay *theRadioactiveDecayContainer)
void SetGuidance(const char *theGuidance)
G4UIcmdWithNucleusLimits * nucleuslimitsCmd
void AddUserDecayDataFile(G4int Z, G4int A, G4String filename)
G4UIcmdWithADoubleAndUnit * hlthCmd