Geant4_10
G4CascadeParamMessenger.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: G4CascadeParamMessenger.cc 72016 2013-07-03 16:24:15Z mkelsey $
27 // Define simple UI commands as alternative to environment variables
28 //
29 // 20130304 M. Kelsey -- Add flag to collect and display cascade structure
30 // 20130621 M. Kelsey -- Add flag for CHECK_ECONS, replacing #ifdef's; add
31 // flag to use three-body momentum parametrizations
32 // 20130703 M. Kelsey -- Add flag for USE_PHASESPACE
33 
35 #include "G4CascadeParameters.hh"
36 #include "G4UIcmdWithABool.hh"
37 #include "G4UIcmdWithADouble.hh"
38 #include "G4UIcmdWithAString.hh"
39 #include "G4UIcmdWithAnInteger.hh"
41 #include "G4UIcommand.hh"
42 #include "G4UIcommandTree.hh"
43 #include "G4UIdirectory.hh"
44 #include "G4UImanager.hh"
45 
46 
47 // Constructor and destructor
48 
50  : G4UImessenger(), theParams(params), cmdDir(0), localCmdDir(false) {
51  // NOTE: Put under same top-level tree as EM
52  CreateDirectory("/process/had/cascade/","Bertini-esque cascade parameters");
53 
54  verboseCmd = CreateCommand<G4UIcmdWithAnInteger>("verbose",
55  "Enable information messages");
56  balanceCmd = CreateCommand<G4UIcmdWithABool>("checkBalance",
57  "Enable internal conservation checking");
58  reportCmd = CreateCommand<G4UIcmdWithoutParameter>("report",
59  "Dump all non-default parameter settings");
60  usePreCoCmd = CreateCommand<G4UIcmdWithABool>("usePreCompound",
61  "Use PreCompoundModel for nuclear de-excitation");
62  doCoalCmd = CreateCommand<G4UIcmdWithABool>("doCoalescence",
63  "Apply final-state nucleon clustering");
64  historyCmd = CreateCommand<G4UIcmdWithABool>("showHistory",
65  "Collect and report full structure of cascade");
66  use3BodyCmd = CreateCommand<G4UIcmdWithABool>("use3BodyMom",
67  "Use three-body momentum parametrizations");
68  usePSCmd = CreateCommand<G4UIcmdWithABool>("usePhaseSpace",
69  "Use Kopylov N-body momentum generator");
70  randomFileCmd = CreateCommand<G4UIcmdWithAString>("randomFile",
71  "Save random-engine to file at each interaction");
72  nucUseBestCmd = CreateCommand<G4UIcmdWithABool>("useBestNuclearModel",
73  "Use all physical-units for nuclear structure");
74  nucRad2parCmd = CreateCommand<G4UIcmdWithADouble>("useTwoParamNuclearRadius",
75  "Use R = C1*cbrt(A) + C2/cbrt(A)");
76  nucRadScaleCmd = CreateCommand<G4UIcmdWithADouble>("nuclearRadiusScale",
77  "Set length scale for nuclear model");
78  nucRadSmallCmd = CreateCommand<G4UIcmdWithADouble>("smallNucleusRadius",
79  "Set radius of A<4 nuclei");
80  nucRadAlphaCmd = CreateCommand<G4UIcmdWithADouble>("alphaRadiusScale",
81  "Fraction of small-radius for He-4");
82  nucRadTrailingCmd = CreateCommand<G4UIcmdWithADouble>("shadowningRadius",
83  "Effective nucleon radius for trailing effect");
84  nucFermiScaleCmd = CreateCommand<G4UIcmdWithADouble>("fermiScale",
85  "Scale factor for fermi momentum");
86  nucXsecScaleCmd = CreateCommand<G4UIcmdWithADouble>("crossSectionScale",
87  "Scale fator for total cross-sections");
88  nucGammaQDCmd = CreateCommand<G4UIcmdWithADouble>("gammaQuasiDeutScale",
89  "Scale factor for gamma-quasideutron cross-sections");
90  coalDPmax2Cmd = CreateCommand<G4UIcmdWithADouble>("cluster2DPmax",
91  "Maximum momentum for p-n clusters");
92  coalDPmax3Cmd = CreateCommand<G4UIcmdWithADouble>("cluster3DPmax",
93  "Maximum momentum for ppn/pnn clusters");
94  coalDPmax4Cmd = CreateCommand<G4UIcmdWithADouble>("cluster4DPmax",
95  "Maximum momentum for alpha clusters");
96 }
97 
99  delete verboseCmd;
100  delete balanceCmd;
101  delete reportCmd;
102  delete usePreCoCmd;
103  delete doCoalCmd;
104  delete historyCmd;
105  delete use3BodyCmd;
106  delete usePSCmd;
107  delete randomFileCmd;
108  delete nucUseBestCmd;
109  delete nucRad2parCmd;
110  delete nucRadScaleCmd;
111  delete nucRadSmallCmd;
112  delete nucRadAlphaCmd;
113  delete nucRadTrailingCmd;
114  delete nucFermiScaleCmd;
115  delete nucXsecScaleCmd;
116  delete nucGammaQDCmd;
117  delete coalDPmax2Cmd;
118  delete coalDPmax3Cmd;
119  delete coalDPmax4Cmd;
120  if (localCmdDir) delete cmdDir;
121 }
122 
123 
124 // Create or reuse existing UIdirectory path
125 
127  const char* desc) {
129  if (!UIman) return;
130 
131  // Directory path must be absolute, prepend "/" if ncessary
132  G4String fullPath = path;
133  if (fullPath(0) != '/') fullPath.prepend("/");
134  if (fullPath(fullPath.length()-1) != '/') fullPath.append("/");
135 
136  // See if input path has already been registered
137  G4UIcommand* foundPath = UIman->GetTree()->FindPath(fullPath);
138  if (foundPath) cmdDir = dynamic_cast<G4UIdirectory*>(foundPath);
139 
140  if (!cmdDir) { // Create local deletable directory
141  localCmdDir = true;
142  cmdDir = new G4UIdirectory(fullPath.c_str());
143  cmdDir->SetGuidance(desc);
144  }
145 }
146 
147 
148 // Use command argument (literal string) to set envvar maps in container
149 
151  if (cmd == reportCmd) theParams->DumpConfig(G4cout);
152 
153  if (cmd == verboseCmd)
154  theParams->G4CASCADE_VERBOSE = strdup(arg.c_str());
155 
156  if (cmd == balanceCmd)
157  theParams->G4CASCADE_CHECK_ECONS = StoB(arg) ? strdup(arg.c_str()) : 0;
158 
159  if (cmd == usePreCoCmd)
160  theParams->G4CASCADE_USE_PRECOMPOUND = StoB(arg) ? strdup(arg.c_str()) : 0;
161 
162  if (cmd == doCoalCmd)
163  theParams->G4CASCADE_DO_COALESCENCE = StoB(arg) ? strdup(arg.c_str()) : 0;
164 
165  if (cmd == historyCmd)
166  theParams->G4CASCADE_SHOW_HISTORY = StoB(arg) ? strdup(arg.c_str()) : 0;
167 
168  if (cmd == use3BodyCmd)
169  theParams->G4CASCADE_USE_3BODYMOM = StoB(arg) ? strdup(arg.c_str()) : 0;
170 
171  if (cmd == usePSCmd)
172  theParams->G4CASCADE_USE_PHASESPACE = StoB(arg) ? strdup(arg.c_str()) : 0;
173 
174  if (cmd == randomFileCmd)
175  theParams->G4CASCADE_RANDOM_FILE = arg.empty() ? 0 : strdup(arg.c_str());
176 
177  if (cmd == nucUseBestCmd)
178  theParams->G4NUCMODEL_USE_BEST = StoB(arg) ? strdup(arg.c_str()) : 0;
179 
180  if (cmd == nucRad2parCmd)
181  theParams->G4NUCMODEL_RAD_2PAR = strdup(arg.c_str());
182 
183  if (cmd == nucRadScaleCmd)
184  theParams->G4NUCMODEL_RAD_SCALE = strdup(arg.c_str());
185 
186  if (cmd == nucRadSmallCmd)
187  theParams->G4NUCMODEL_RAD_SMALL = strdup(arg.c_str());
188 
189  if (cmd == nucRadAlphaCmd)
190  theParams->G4NUCMODEL_RAD_ALPHA = strdup(arg.c_str());
191 
192  if (cmd == nucRadTrailingCmd)
193  theParams->G4NUCMODEL_RAD_TRAILING = strdup(arg.c_str());
194 
195  if (cmd == nucFermiScaleCmd)
196  theParams->G4NUCMODEL_FERMI_SCALE = strdup(arg.c_str());
197 
198  if (cmd == nucXsecScaleCmd)
199  theParams->G4NUCMODEL_XSEC_SCALE = strdup(arg.c_str());
200 
201  if (cmd == nucGammaQDCmd)
202  theParams->G4NUCMODEL_GAMMAQD = strdup(arg.c_str());
203 
204  if (cmd == coalDPmax2Cmd)
205  theParams->DPMAX_2CLUSTER = strdup(arg.c_str());
206 
207  if (cmd == coalDPmax3Cmd)
208  theParams->DPMAX_3CLUSTER = strdup(arg.c_str());
209 
210  if (cmd == coalDPmax4Cmd)
211  theParams->DPMAX_4CLUSTER = strdup(arg.c_str());
212 
213  theParams->Initialize(); // Update numerical values from settings
214 }
char cmd[1024]
Definition: tracer.cxx:25
G4UIcommand * FindPath(const char *commandPath) const
static G4UImanager * GetUIpointer()
Definition: G4UImanager.cc:58
G4bool StoB(G4String s)
G4String & prepend(const char *)
G4GLOB_DLL std::ostream G4cout
void SetGuidance(const char *aGuidance)
Definition: G4UIcommand.hh:161
virtual void SetNewValue(G4UIcommand *command, G4String newValue)
G4UIcommandTree * GetTree() const
Definition: G4UImanager.hh:206
G4String & append(const G4String &)
G4CascadeParamMessenger(G4CascadeParameters *params)
void CreateDirectory(const char *path, const char *desc)