Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4PhysListFactory.cc
이 파일의 문서화 페이지로 가기
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: G4PhysListFactory.cc 107319 2017-11-08 16:29:22Z gcosmo $
27 //
28 //---------------------------------------------------------------------------
29 //
30 // ClassName: G4PhysListFactory
31 //
32 // Author: 21 April 2008 V. Ivanchenko
33 //
34 // Modified:
35 //
36 // 2014.08.05 K.L.Genser used provision for Hadronic Physics Variant M in
37 // Shielding for ShieldingM
38 //
39 //----------------------------------------------------------------------------
40 //
41 
42 #include "G4PhysListFactory.hh"
43 #include "FTFP_BERT.hh"
44 #include "FTFP_BERT_HP.hh"
45 #include "FTFP_BERT_TRV.hh"
46 #include "FTFP_BERT_ATL.hh"
47 #include "FTFQGSP_BERT.hh"
48 #include "FTFP_INCLXX.hh"
49 #include "FTFP_INCLXX_HP.hh"
50 #include "FTF_BIC.hh"
51 #include "LBE.hh"
52 #include "QBBC.hh"
53 #include "QGSP_BERT.hh"
54 #include "QGSP_BERT_HP.hh"
55 #include "QGSP_BIC.hh"
56 #include "QGSP_BIC_HP.hh"
57 #include "QGSP_BIC_AllHP.hh"
58 #include "QGSP_FTFP_BERT.hh"
59 #include "QGS_BIC.hh"
60 #include "QGSP_INCLXX.hh"
61 #include "QGSP_INCLXX_HP.hh"
62 #include "Shielding.hh"
63 #include "ShieldingLEND.hh"
64 #include "NuBeam.hh"
65 
66 #include "G4EmStandardPhysics.hh"
71 #include "G4EmStandardPhysicsGS.hh"
72 #include "G4EmStandardPhysicsSS.hh"
73 #include "G4EmLivermorePhysics.hh"
74 #include "G4EmPenelopePhysics.hh"
76 #include "G4UImessenger.hh"
77 
79  : defName("FTFP_BERT"),verbose(1),theMessenger(nullptr)
80 {
81  nlists_hadr = 23;
82  G4String ss[23] = {
83  "FTFP_BERT","FTFP_BERT_TRV","FTFP_BERT_ATL","FTFP_BERT_HP","FTFQGSP_BERT",
84  "FTFP_INCLXX","FTFP_INCLXX_HP","FTF_BIC", "LBE","QBBC",
85  "QGSP_BERT","QGSP_BERT_HP","QGSP_BIC","QGSP_BIC_HP","QGSP_BIC_AllHP",
86  "QGSP_FTFP_BERT","QGSP_INCLXX","QGSP_INCLXX_HP","QGS_BIC",
87  "Shielding","ShieldingLEND","ShieldingM","NuBeam"};
88  for(size_t i=0; i<nlists_hadr; ++i) {
89  listnames_hadr.push_back(ss[i]);
90  }
91 
92  nlists_em = 9;
93  G4String s1[9] = {"","_EMV","_EMX","_EMY","_EMZ","_LIV","_PEN","__GS","__SS"};
94  for(size_t i=0; i<nlists_em; ++i) {
95  listnames_em.push_back(s1[i]);
96  }
97 }
98 
100 {
101  delete theMessenger;
102 }
103 
106 {
107  // instantiate PhysList by environment variable "PHYSLIST"
108  G4String name = "";
109  char* path = getenv("PHYSLIST");
110  if (path) {
111  name = G4String(path);
112  } else {
113  name = defName;
114  G4cout << "### G4PhysListFactory WARNING: "
115  << " environment variable PHYSLIST is not defined"
116  << G4endl
117  << " Default Physics Lists " << name
118  << " is instantiated"
119  << G4endl;
120  }
121  return GetReferencePhysList(name);
122 }
123 
126 {
127  // analysis on the string
128  size_t n = name.size();
129 
130  // last characters in the string
131  size_t em_opt = 0;
132  G4String em_name = "";
133 
134  // check EM options
135  if(n > 4) {
136  em_name = name.substr(n - 4, 4);
137  for(size_t i=1; i<nlists_em; ++i) {
138  if(listnames_em[i] == em_name) {
139  em_opt = i;
140  n -= 4;
141  break;
142  }
143  }
144  if(0 == em_opt) { em_name = ""; }
145  }
146 
147  // hadronic pHysics List
148  G4String had_name = name.substr(0, n);
149 
150  if(0 < verbose) {
151  G4cout << "G4PhysListFactory::GetReferencePhysList <" << had_name
152  << em_name << "> EMoption= " << em_opt << G4endl;
153  }
154  G4VModularPhysicsList* p = nullptr;
155  if(had_name == "FTFP_BERT") {p = new FTFP_BERT(verbose);}
156  else if(had_name == "FTFP_BERT_HP") {p = new FTFP_BERT_HP(verbose);}
157  else if(had_name == "FTFP_BERT_TRV") {p = new FTFP_BERT_TRV(verbose);}
158  else if(had_name == "FTFP_BERT_ATL") {p = new FTFP_BERT_ATL(verbose);}
159  else if(had_name == "FTFQGSP_BERT") {p = new FTFQGSP_BERT(verbose);}
160  else if(had_name == "FTFP_INCLXX") {p = new FTFP_INCLXX(verbose);}
161  else if(had_name == "FTFP_INCLXX_HP") {p = new FTFP_INCLXX_HP(verbose);}
162  else if(had_name == "FTF_BIC") {p = new FTF_BIC(verbose);}
163  else if(had_name == "LBE") {p = new LBE();}
164  else if(had_name == "QBBC") {p = new QBBC(verbose);}
165  else if(had_name == "QGSP_BERT") {p = new QGSP_BERT(verbose);}
166  else if(had_name == "QGSP_BERT_HP") {p = new QGSP_BERT_HP(verbose);}
167  else if(had_name == "QGSP_BIC") {p = new QGSP_BIC(verbose);}
168  else if(had_name == "QGSP_BIC_HP") {p = new QGSP_BIC_HP(verbose);}
169  else if(had_name == "QGSP_BIC_AllHP") {p = new QGSP_BIC_AllHP(verbose);}
170  else if(had_name == "QGSP_FTFP_BERT") {p = new QGSP_FTFP_BERT(verbose);}
171  else if(had_name == "QGSP_INCLXX") {p = new QGSP_INCLXX(verbose);}
172  else if(had_name == "QGSP_INCLXX_HP") {p = new QGSP_INCLXX_HP(verbose);}
173  else if(had_name == "QGS_BIC") {p = new QGS_BIC(verbose);}
174  else if(had_name == "Shielding") {p = new Shielding(verbose);}
175  else if(had_name == "ShieldingLEND") {p = new ShieldingLEND(verbose);}
176  else if(had_name == "ShieldingM") {p = new Shielding(verbose,"HP","M");}
177  else if(had_name == "NuBeam") {p = new NuBeam(verbose);}
178  else {
179  G4cout << "### G4PhysListFactory WARNING: "
180  << "PhysicsList " << had_name << " is not known"
181  << G4endl;
182  }
183  if(p) {
184  G4cout << "<<< Reference Physics List " << had_name
185  << em_name << " is built" << G4endl;
186  G4int ver = p->GetVerboseLevel();
187  p->SetVerboseLevel(0);
188  if(0 < em_opt && had_name != "LBE") {
189  if(1 == em_opt) {
191  } else if(2 == em_opt) {
193  } else if(3 == em_opt) {
195  } else if(4 == em_opt) {
197  } else if(5 == em_opt) {
199  } else if(6 == em_opt) {
201  } else if(7 == em_opt) {
203  } else if(8 == em_opt) {
205  }
206  }
207  p->SetVerboseLevel(ver);
209  }
210  G4cout << G4endl;
211  return p;
212 }
213 
215 {
216  G4bool res = false;
217  size_t n = name.size();
218  if(n > 4) {
219  G4String em_name = name.substr(n - 4, 4);
220  for(size_t i=1; i<nlists_em; ++i) {
221  if(listnames_em[i] == em_name) {
222  n -= 4;
223  break;
224  }
225  }
226  }
227  G4String had_name = name.substr(0, n);
228  for(size_t i=0; i<nlists_hadr; ++i) {
229  if(had_name == listnames_hadr[i]) {
230  res = true;
231  break;
232  }
233  }
234  return res;
235 }
236 
237 const std::vector<G4String>&
239 {
240  return listnames_hadr;
241 }
242 
243 const std::vector<G4String>&
245 {
246  return listnames_em;
247 }
248 
TINCLXXPhysicsListHelper< G4VModularPhysicsList, true, true > FTFP_INCLXX_HP
TNuBeam< G4VModularPhysicsList > NuBeam
Definition: NuBeam.hh:64
const XML_Char * name
Definition: expat.h:151
G4UImessenger * theMessenger
TQGS_BIC< G4VModularPhysicsList > QGS_BIC
Definition: QGS_BIC.hh:63
TINCLXXPhysicsListHelper< G4VModularPhysicsList, false, false > QGSP_INCLXX
#define G4endl
Definition: G4ios.hh:61
const char * p
Definition: xmltok.h:285
TLBE< G4VModularPhysicsList > LBE
Definition: LBE.hh:123
const std::vector< G4String > & AvailablePhysListsEM() const
G4bool IsReferencePhysList(const G4String &)
G4VModularPhysicsList * GetReferencePhysList(const G4String &)
Definition: QBBC.hh:44
const std::vector< G4String > & AvailablePhysLists() const
TFTFQGSP_BERT< G4VModularPhysicsList > FTFQGSP_BERT
Definition: FTFQGSP_BERT.hh:64
TQGSP_BIC_AllHP< G4VModularPhysicsList > QGSP_BIC_AllHP
bool G4bool
Definition: G4Types.hh:79
TINCLXXPhysicsListHelper< G4VModularPhysicsList, false, true > FTFP_INCLXX
Definition: FTFP_INCLXX.hh:41
TINCLXXPhysicsListHelper< G4VModularPhysicsList, true, false > QGSP_INCLXX_HP
TFTF_BIC< G4VModularPhysicsList > FTF_BIC
Definition: FTF_BIC.hh:62
TFTFP_BERT_HP< G4VModularPhysicsList > FTFP_BERT_HP
Definition: FTFP_BERT_HP.hh:63
TFTFP_BERT< G4VModularPhysicsList > FTFP_BERT
Definition: FTFP_BERT.hh:63
TQGSP_BIC_HP< G4VModularPhysicsList > QGSP_BIC_HP
Definition: QGSP_BIC_HP.hh:63
void ReplacePhysics(G4VPhysicsConstructor *)
int G4int
Definition: G4Types.hh:78
TShieldingLEND< G4VModularPhysicsList > ShieldingLEND
TQGSP_BIC< G4VModularPhysicsList > QGSP_BIC
Definition: QGSP_BIC.hh:62
TQGSP_BERT< G4VModularPhysicsList > QGSP_BERT
Definition: QGSP_BERT.hh:63
TQGSP_BERT_HP< G4VModularPhysicsList > QGSP_BERT_HP
Definition: QGSP_BERT_HP.hh:62
std::vector< G4String > listnames_hadr
G4VModularPhysicsList * ReferencePhysList()
G4GLOB_DLL std::ostream G4cout
TFTFP_BERT_TRV< G4VModularPhysicsList > FTFP_BERT_TRV
Char_t n[5]
TShielding< G4VModularPhysicsList > Shielding
Definition: Shielding.hh:68
std::vector< G4String > listnames_em
TQGSP_FTFP_BERT< G4VModularPhysicsList > QGSP_FTFP_BERT
TFTFP_BERT_ATL< G4VModularPhysicsList > FTFP_BERT_ATL
void SetVerboseLevel(G4int value)