Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
B03PhysicsList.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 //
28 //
29 //
30 // $Id: B03PhysicsList.cc 75089 2013-10-25 23:25:21Z dwright $
31 //
32 
33 #include "globals.hh"
34 #include <iomanip>
35 
36 #include "B03PhysicsList.hh"
37 
38 #include "G4ParticleDefinition.hh"
39 #include "G4ParticleWithCuts.hh"
40 #include "G4ProcessManager.hh"
41 #include "G4ProcessVector.hh"
42 #include "G4ParticleTypes.hh"
43 #include "G4ParticleTable.hh"
44 #include "G4BosonConstructor.hh"
45 #include "G4LeptonConstructor.hh"
46 #include "G4MesonConstructor.hh"
47 #include "G4BaryonConstructor.hh"
48 #include "G4IonConstructor.hh"
50 #include "G4Material.hh"
51 #include "G4MaterialTable.hh"
52 #include "G4SystemOfUnits.hh"
53 
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
55 
58 fBiasWorldName(parallelname)
59 {
60  fParaWorldName.clear();
61  SetVerboseLevel(1);
62 }
63 
64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
65 
67 {
68  fParaWorldName.clear();
69 }
70 
71 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
72 
74 {
75  // In this method, static member functions should be called
76  // for all particles which you want to use.
77  // This ensures that objects of these particle types will be
78  // created in the program.
79 
86 }
87 
88 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
89 
91 {
92  // Construct all bosons
93  G4BosonConstructor pConstructor;
94  pConstructor.ConstructParticle();
95 }
96 
97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
98 
100 {
101  // Construct all leptons
102  G4LeptonConstructor pConstructor;
103  pConstructor.ConstructParticle();
104 }
105 
106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
107 
109 {
110  // Construct all mesons
111  G4MesonConstructor pConstructor;
112  pConstructor.ConstructParticle();
113 }
114 
115 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
116 
118 {
119  // Construct all barions
120  G4BaryonConstructor pConstructor;
121  pConstructor.ConstructParticle();
122 }
123 
124 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
125 
127 {
128  // Construct light ions
129  G4IonConstructor pConstructor;
130  pConstructor.ConstructParticle();
131 }
132 
133 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
134 
136 {
137  // Construct resonaces and quarks
138  G4ShortLivedConstructor pConstructor;
139  pConstructor.ConstructParticle();
140 }
141 
142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
143 
145 {
149  ConstructEM();
151  ConstructHad();
153 }
154 
155 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
156 
157 #include "G4ComptonScattering.hh"
158 #include "G4GammaConversion.hh"
159 #include "G4PhotoElectricEffect.hh"
160 
161 #include "G4eMultipleScattering.hh"
162 #include "G4MuMultipleScattering.hh"
163 #include "G4hMultipleScattering.hh"
164 
165 #include "G4eIonisation.hh"
166 #include "G4eBremsstrahlung.hh"
167 #include "G4eplusAnnihilation.hh"
168 
169 #include "G4MuIonisation.hh"
170 #include "G4MuBremsstrahlung.hh"
171 #include "G4MuPairProduction.hh"
172 
173 #include "G4hIonisation.hh"
174 
176 {
178  particleIterator->reset();
179  while( (*particleIterator)() ){
180  G4ParticleDefinition* particle = particleIterator->value();
181  G4ProcessManager* pmanager = particle->GetProcessManager();
182  G4String particleName = particle->GetParticleName();
183 
184  if (particleName == "gamma") {
185  // gamma
186  // Construct processes for gamma
187  pmanager->AddDiscreteProcess(new G4GammaConversion());
188  pmanager->AddDiscreteProcess(new G4ComptonScattering());
189  pmanager->AddDiscreteProcess(new G4PhotoElectricEffect());
190 
191  } else if (particleName == "e-") {
192  //electron
193  // Construct processes for electron
194  pmanager->AddProcess(new G4eMultipleScattering(),-1,1,1);
195  pmanager->AddProcess(new G4eIonisation(),-1,2,2);
196  pmanager->AddProcess(new G4eBremsstrahlung(),-1,-1,3);
197 
198  } else if (particleName == "e+") {
199  //positron
200  // Construct processes for positron
201  pmanager->AddProcess(new G4eMultipleScattering(),-1,1,1);
202 
203  pmanager->AddProcess(new G4eIonisation(),-1,2,2);
204  pmanager->AddProcess(new G4eBremsstrahlung(),-1,-1,3);
205  pmanager->AddProcess(new G4eplusAnnihilation(),0,-1,4);
206 
207  } else if( particleName == "mu+" ||
208  particleName == "mu-" ) {
209  //muon
210  // Construct processes for muon+
211  pmanager->AddProcess(new G4MuMultipleScattering(),-1,1,1);
212  pmanager->AddProcess(new G4MuIonisation(),-1,2,2);
213  pmanager->AddProcess(new G4MuBremsstrahlung(),-1,-1,3);
214  pmanager->AddProcess(new G4MuPairProduction(),-1,-1,4);
215 
216  } else if( particleName == "GenericIon" ) {
217  pmanager->AddProcess(new G4hMultipleScattering(),-1,1,1);
218  pmanager->AddProcess(new G4hIonisation(),-1,2,2);
219  } else {
220  if ((particle->GetPDGCharge() != 0.0) &&
221  (particle->GetParticleName() != "chargedgeantino")&&
222  (!particle->IsShortLived()) ) {
223  // all others charged particles except geantino
224  pmanager->AddProcess(new G4hMultipleScattering(),-1,1,1);
225  pmanager->AddProcess(new G4hIonisation(),-1,2,2);
226  }
227  }
228  }
229 }
230 
231 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
232 
233 // Hadron Processes
234 
235 #include "G4HadronElasticProcess.hh"
236 #include "G4HadronFissionProcess.hh"
237 #include "G4HadronCaptureProcess.hh"
238 
264 
265 // Low-energy Models
266 
267 #include "G4HadronElastic.hh"
268 #include "G4NeutronRadCapture.hh"
269 #include "G4LFission.hh"
270 
271 // -- generator models
272 #include "G4TheoFSGenerator.hh"
273 #include "G4ExcitationHandler.hh"
274 #include "G4CompetitiveFission.hh"
276 #include "G4Fancy3DNucleus.hh"
277 #include "G4CascadeInterface.hh"
278 #include "G4StringModel.hh"
279 #include "G4PreCompoundModel.hh"
280 #include "G4FTFModel.hh"
281 #include "G4QGSMFragmentation.hh"
283 #include "G4ExcitedStringDecay.hh"
284 #include "G4QMDReaction.hh"
286 
287 // Cross sections
288 #include "G4IonsShenCrossSection.hh"
289 #include "G4TripathiCrossSection.hh"
291 
292 //
293 // ConstructHad()
294 //
295 // Makes discrete physics processes for the hadrons
296 // The processes are: Elastic scattering, Inelastic scattering,
297 // Fission (for neutron only), and Capture (neutron).
298 //
299 
301 {
302  // this will be the model class for high energies
303  G4TheoFSGenerator* theTheoModel = new G4TheoFSGenerator;
304  G4TheoFSGenerator* antiBHighEnergyModel = new G4TheoFSGenerator;
305 
306  // Evaporation logic
307  G4ExcitationHandler* theHandler = new G4ExcitationHandler;
308  theHandler->SetMinEForMultiFrag(3*MeV);
309 
310  // Pre equilibrium stage
311  G4PreCompoundModel* thePreEquilib = new G4PreCompoundModel(theHandler);
312 
313  // a no-cascade generator-precompound interaface
314  G4GeneratorPrecompoundInterface* theCascade =
316  theCascade->SetDeExcitation(thePreEquilib);
317 
318  // Bertini cascade
319  G4CascadeInterface* bertini = new G4CascadeInterface;
320  bertini->SetMaxEnergy(22*MeV);
321 
322  // here come the high energy parts
323  G4VPartonStringModel* theStringModel;
324  theStringModel = new G4FTFModel;
325  theTheoModel->SetTransport(theCascade);
326  theTheoModel->SetHighEnergyGenerator(theStringModel);
327  theTheoModel->SetMinEnergy(19*GeV);
328  theTheoModel->SetMaxEnergy(100*TeV);
329 
330  G4VLongitudinalStringDecay* theFragmentation = new G4QGSMFragmentation;
331  G4ExcitedStringDecay* theStringDecay =
332  new G4ExcitedStringDecay(theFragmentation);
333  theStringModel->SetFragmentationModel(theStringDecay);
334 
335  // high energy model for anti-baryons
336  antiBHighEnergyModel = new G4TheoFSGenerator("ANTI-FTFP");
337  G4FTFModel* antiBStringModel = new G4FTFModel;
338  G4ExcitedStringDecay* stringDecay =
340  antiBStringModel->SetFragmentationModel(stringDecay);
341 
342  G4GeneratorPrecompoundInterface* antiBCascade =
344  G4PreCompoundModel* preEquilib =
346  antiBCascade->SetDeExcitation(preEquilib);
347 
348  antiBHighEnergyModel->SetTransport(antiBCascade);
349  antiBHighEnergyModel->SetHighEnergyGenerator(antiBStringModel);
350  antiBHighEnergyModel->SetMinEnergy(0.0);
351  antiBHighEnergyModel->SetMaxEnergy(20*TeV);
352 
353  // Light ion models
355  binaryCascade->SetMinEnergy(0.0);
356  binaryCascade->SetMaxEnergy(110*MeV);
357 
358  G4QMDReaction* qmd = new G4QMDReaction;
359  qmd->SetMinEnergy(100*MeV);
360  qmd->SetMaxEnergy(10*GeV);
361 
365 
366  // Elastic process
367  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
368  G4HadronElastic* theElasticModel = new G4HadronElastic;
369  theElasticProcess->RegisterMe(theElasticModel);
370 
372  particleIterator->reset();
373  while ((*particleIterator)()) {
374  G4ParticleDefinition* particle = particleIterator->value();
375  G4ProcessManager* pmanager = particle->GetProcessManager();
376  G4String particleName = particle->GetParticleName();
377 
378  if (particleName == "pi+") {
379  pmanager->AddDiscreteProcess(theElasticProcess);
380  G4PionPlusInelasticProcess* theInelasticProcess =
381  new G4PionPlusInelasticProcess("inelastic");
382  theInelasticProcess->RegisterMe(bertini);
383  theInelasticProcess->RegisterMe(theTheoModel);
384  pmanager->AddDiscreteProcess(theInelasticProcess);
385  } else if (particleName == "pi-") {
386  pmanager->AddDiscreteProcess(theElasticProcess);
387  G4PionMinusInelasticProcess* theInelasticProcess =
388  new G4PionMinusInelasticProcess("inelastic");
389  theInelasticProcess->RegisterMe(bertini);
390  theInelasticProcess->RegisterMe(theTheoModel);
391  pmanager->AddDiscreteProcess(theInelasticProcess);
392  } else if (particleName == "kaon+") {
393  pmanager->AddDiscreteProcess(theElasticProcess);
394  G4KaonPlusInelasticProcess* theInelasticProcess =
395  new G4KaonPlusInelasticProcess("inelastic");
396  theInelasticProcess->RegisterMe(bertini);
397  theInelasticProcess->RegisterMe(theTheoModel);
398  pmanager->AddDiscreteProcess(theInelasticProcess);
399  }
400  else if (particleName == "kaon0S") {
401  pmanager->AddDiscreteProcess(theElasticProcess);
402  G4KaonZeroSInelasticProcess* theInelasticProcess =
403  new G4KaonZeroSInelasticProcess("inelastic");
404  theInelasticProcess->RegisterMe(bertini);
405  theInelasticProcess->RegisterMe(theTheoModel);
406  pmanager->AddDiscreteProcess(theInelasticProcess);
407  }
408  else if (particleName == "kaon0L") {
409  pmanager->AddDiscreteProcess(theElasticProcess);
410  G4KaonZeroLInelasticProcess* theInelasticProcess =
411  new G4KaonZeroLInelasticProcess("inelastic");
412  theInelasticProcess->RegisterMe(bertini);
413  theInelasticProcess->RegisterMe(theTheoModel);
414  pmanager->AddDiscreteProcess(theInelasticProcess);
415  }
416  else if (particleName == "kaon-") {
417  pmanager->AddDiscreteProcess(theElasticProcess);
418  G4KaonMinusInelasticProcess* theInelasticProcess =
419  new G4KaonMinusInelasticProcess("inelastic");
420  theInelasticProcess->RegisterMe(bertini);
421  theInelasticProcess->RegisterMe(theTheoModel);
422  pmanager->AddDiscreteProcess(theInelasticProcess);
423  }
424  else if (particleName == "proton") {
425  pmanager->AddDiscreteProcess(theElasticProcess);
426  G4ProtonInelasticProcess* theInelasticProcess =
427  new G4ProtonInelasticProcess("inelastic");
428  theInelasticProcess->RegisterMe(bertini);
429  theInelasticProcess->RegisterMe(theTheoModel);
430  pmanager->AddDiscreteProcess(theInelasticProcess);
431  }
432  else if (particleName == "anti_proton") {
433  pmanager->AddDiscreteProcess(theElasticProcess);
434  G4AntiProtonInelasticProcess* theInelasticProcess =
435  new G4AntiProtonInelasticProcess("inelastic");
436  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
437  pmanager->AddDiscreteProcess(theInelasticProcess);
438 
439  } else if (particleName == "neutron") {
440  // elastic scattering
441  pmanager->AddDiscreteProcess(theElasticProcess);
442 
443  // inelastic scattering
444  G4NeutronInelasticProcess* theInelasticProcess =
445  new G4NeutronInelasticProcess("inelastic");
446  theInelasticProcess->RegisterMe(bertini);
447  theInelasticProcess->RegisterMe(theTheoModel);
448  pmanager->AddDiscreteProcess(theInelasticProcess);
449 
450  // fission
451  G4HadronFissionProcess* theFissionProcess = new G4HadronFissionProcess;
452  G4LFission* theFissionModel = new G4LFission;
453  theFissionProcess->RegisterMe(theFissionModel);
454  pmanager->AddDiscreteProcess(theFissionProcess);
455 
456  // capture
457  G4HadronCaptureProcess* theCaptureProcess = new G4HadronCaptureProcess;
458  G4NeutronRadCapture* theCaptureModel = new G4NeutronRadCapture;
459  theCaptureProcess->RegisterMe(theCaptureModel);
460  pmanager->AddDiscreteProcess(theCaptureProcess);
461 
462  } else if (particleName == "anti_neutron") {
463  pmanager->AddDiscreteProcess(theElasticProcess);
464  G4AntiNeutronInelasticProcess* theInelasticProcess =
465  new G4AntiNeutronInelasticProcess("inelastic");
466  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
467  pmanager->AddDiscreteProcess(theInelasticProcess);
468 
469  } else if (particleName == "lambda") {
470  pmanager->AddDiscreteProcess(theElasticProcess);
471  G4LambdaInelasticProcess* theInelasticProcess =
472  new G4LambdaInelasticProcess("inelastic");
473  theInelasticProcess->RegisterMe(bertini);
474  theInelasticProcess->RegisterMe(theTheoModel);
475  pmanager->AddDiscreteProcess(theInelasticProcess);
476  }
477  else if (particleName == "anti_lambda") {
478  pmanager->AddDiscreteProcess(theElasticProcess);
479  G4AntiLambdaInelasticProcess* theInelasticProcess =
480  new G4AntiLambdaInelasticProcess("inelastic");
481  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
482  pmanager->AddDiscreteProcess(theInelasticProcess);
483  }
484  else if (particleName == "sigma+") {
485  pmanager->AddDiscreteProcess(theElasticProcess);
486  G4SigmaPlusInelasticProcess* theInelasticProcess =
487  new G4SigmaPlusInelasticProcess("inelastic");
488  theInelasticProcess->RegisterMe(bertini);
489  theInelasticProcess->RegisterMe(theTheoModel);
490  pmanager->AddDiscreteProcess(theInelasticProcess);
491  }
492  else if (particleName == "sigma-") {
493  pmanager->AddDiscreteProcess(theElasticProcess);
494  G4SigmaMinusInelasticProcess* theInelasticProcess =
495  new G4SigmaMinusInelasticProcess("inelastic");
496  theInelasticProcess->RegisterMe(bertini);
497  theInelasticProcess->RegisterMe(theTheoModel);
498  pmanager->AddDiscreteProcess(theInelasticProcess);
499  }
500  else if (particleName == "anti_sigma+") {
501  pmanager->AddDiscreteProcess(theElasticProcess);
502  G4AntiSigmaPlusInelasticProcess* theInelasticProcess =
503  new G4AntiSigmaPlusInelasticProcess("inelastic");
504  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
505  pmanager->AddDiscreteProcess(theInelasticProcess);
506  }
507  else if (particleName == "anti_sigma-") {
508  pmanager->AddDiscreteProcess(theElasticProcess);
509  G4AntiSigmaMinusInelasticProcess* theInelasticProcess =
510  new G4AntiSigmaMinusInelasticProcess("inelastic");
511  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
512  pmanager->AddDiscreteProcess(theInelasticProcess);
513  }
514  else if (particleName == "xi0") {
515  pmanager->AddDiscreteProcess(theElasticProcess);
516  G4XiZeroInelasticProcess* theInelasticProcess =
517  new G4XiZeroInelasticProcess("inelastic");
518  theInelasticProcess->RegisterMe(bertini);
519  theInelasticProcess->RegisterMe(theTheoModel);
520  pmanager->AddDiscreteProcess(theInelasticProcess);
521  }
522  else if (particleName == "xi-") {
523  pmanager->AddDiscreteProcess(theElasticProcess);
524  G4XiMinusInelasticProcess* theInelasticProcess =
525  new G4XiMinusInelasticProcess("inelastic");
526  theInelasticProcess->RegisterMe(bertini);
527  theInelasticProcess->RegisterMe(theTheoModel);
528  pmanager->AddDiscreteProcess(theInelasticProcess);
529  }
530  else if (particleName == "anti_xi0") {
531  pmanager->AddDiscreteProcess(theElasticProcess);
532  G4AntiXiZeroInelasticProcess* theInelasticProcess =
533  new G4AntiXiZeroInelasticProcess("inelastic");
534  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
535  pmanager->AddDiscreteProcess(theInelasticProcess);
536  }
537  else if (particleName == "anti_xi-") {
538  pmanager->AddDiscreteProcess(theElasticProcess);
539  G4AntiXiMinusInelasticProcess* theInelasticProcess =
540  new G4AntiXiMinusInelasticProcess("inelastic");
541  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
542  pmanager->AddDiscreteProcess(theInelasticProcess);
543  }
544  else if (particleName == "deuteron") {
545  pmanager->AddDiscreteProcess(theElasticProcess);
546  G4DeuteronInelasticProcess* theInelasticProcess =
547  new G4DeuteronInelasticProcess("inelastic");
548  theInelasticProcess->RegisterMe(binaryCascade);
549  theInelasticProcess->RegisterMe(qmd);
550  theInelasticProcess->AddDataSet(shenXS);
551  theInelasticProcess->AddDataSet(tripXS);
552  theInelasticProcess->AddDataSet(tripLightXS);
553  pmanager->AddDiscreteProcess(theInelasticProcess);
554  }
555  else if (particleName == "triton") {
556  pmanager->AddDiscreteProcess(theElasticProcess);
557  G4TritonInelasticProcess* theInelasticProcess =
558  new G4TritonInelasticProcess("inelastic");
559  theInelasticProcess->RegisterMe(binaryCascade);
560  theInelasticProcess->RegisterMe(qmd);
561  theInelasticProcess->AddDataSet(shenXS);
562  theInelasticProcess->AddDataSet(tripXS);
563  theInelasticProcess->AddDataSet(tripLightXS);
564  pmanager->AddDiscreteProcess(theInelasticProcess);
565  }
566  else if (particleName == "alpha") {
567  pmanager->AddDiscreteProcess(theElasticProcess);
568  G4AlphaInelasticProcess* theInelasticProcess =
569  new G4AlphaInelasticProcess("inelastic");
570  theInelasticProcess->RegisterMe(binaryCascade);
571  theInelasticProcess->RegisterMe(qmd);
572  theInelasticProcess->AddDataSet(shenXS);
573  theInelasticProcess->AddDataSet(tripXS);
574  theInelasticProcess->AddDataSet(tripLightXS);
575  pmanager->AddDiscreteProcess(theInelasticProcess);
576 
577  } else if (particleName == "omega-") {
578  pmanager->AddDiscreteProcess(theElasticProcess);
579  G4OmegaMinusInelasticProcess* theInelasticProcess =
580  new G4OmegaMinusInelasticProcess("inelastic");
581  theInelasticProcess->RegisterMe(bertini);
582  theInelasticProcess->RegisterMe(theTheoModel);
583  pmanager->AddDiscreteProcess(theInelasticProcess);
584 
585  } else if (particleName == "anti_omega-") {
586  pmanager->AddDiscreteProcess(theElasticProcess);
587  G4AntiOmegaMinusInelasticProcess* theInelasticProcess =
588  new G4AntiOmegaMinusInelasticProcess("inelastic");
589  theInelasticProcess->RegisterMe(antiBHighEnergyModel);
590  pmanager->AddDiscreteProcess(theInelasticProcess);
591  }
592  }
593 }
594 
595 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
596 
598 {;}
599 
600 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
601 
602 #include "G4Decay.hh"
604 {
605  G4Decay* theDecayProcess = new G4Decay();
607  particleIterator->reset();
608  while( (*particleIterator)() ){
609  G4ParticleDefinition* particle = particleIterator->value();
610  G4ProcessManager* pmanager = particle->GetProcessManager();
611  if (theDecayProcess->IsApplicable(*particle)) {
612  pmanager ->AddProcess(theDecayProcess);
613  pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
614  pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
615  }
616  }
617 }
618 
619 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
620 
622 {
623  if (verboseLevel >0)
624  {
625  G4cout << "B03PhysicsList::SetCuts:";
626  G4cout << "CutLength : " << defaultCutValue/mm << " (mm)" << G4endl;
627  }
628  // "G4VUserPhysicsList::SetCutsWithDefault" method sets
629  // the default cut value for all particle types
631 }
632 
633 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
634 
635 #include "G4ParallelWorldProcess.hh"
637 
638  G4int npw = fParaWorldName.size();
639  for ( G4int i = 0; i < npw; i++){
640  G4String procName = "ParaWorldProc_"+fParaWorldName[i];
641  G4ParallelWorldProcess* theParallelWorldProcess
642  = new G4ParallelWorldProcess(procName);
643  theParallelWorldProcess->SetParallelWorld(fParaWorldName[i]);
644 
646  particleIterator->reset();
647  while( (*particleIterator)() ){
648  G4ParticleDefinition* particle = particleIterator->value();
649  G4ProcessManager* pmanager = particle->GetProcessManager();
650  pmanager->AddProcess(theParallelWorldProcess);
651  if(theParallelWorldProcess->IsAtRestRequired(particle))
652  {pmanager->SetProcessOrdering(theParallelWorldProcess, idxAtRest, 9900);}
653  pmanager->SetProcessOrderingToSecond(theParallelWorldProcess, idxAlongStep);
654  pmanager->SetProcessOrdering(theParallelWorldProcess, idxPostStep, 9900);
655  }
656  }
657 
658 }
659 
660 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
661 
662 #include "G4ImportanceProcess.hh"
663 #include "G4IStore.hh"
665 
666  G4cout << " Preparing Importance Sampling with GhostWorld "
667  << fBiasWorldName << G4endl;
668 
670  G4GeometrySampler fGeomSampler(fBiasWorldName,"neutron");
671  fGeomSampler.SetParallel(true); // parallelworld
672  // fGeomSampler.SetWorld(iStore->GetParallelWorldVolumePointer());
673  // fGeomSampler->PrepareImportanceSampling(G4IStore::
674  // GetInstance(fBiasWorldName), 0);
675  static G4bool first = true;
676  if(first) {
677  fGeomSampler.PrepareImportanceSampling(iStore, 0);
678 
679  fGeomSampler.Configure();
680  G4cout << " GeomSampler Configured!!! " << G4endl;
681  first = false;
682  }
683 
684 #ifdef G4MULTITHREADED
685  if(!G4Threading::IsMasterThread()) fGeomSampler.AddProcess();
686 #else
687  G4cout << " Running in singlethreaded mode!!! " << G4endl;
688 #endif
689 
690 }
691 
692 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
virtual void ConstructParticle()
void ConstructAllLeptons()
static void ConstructParticle()
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
static constexpr double MeV
Definition: G4SIunits.hh:214
void AddBiasingProcess()
G4int AddDiscreteProcess(G4VProcess *aProcess, G4int ord=ordDefault)
static constexpr double mm
Definition: G4SIunits.hh:115
virtual void PrepareImportanceSampling(G4VIStore *istore, const G4VImportanceAlgorithm *ialg)
void SetParallelWorld(G4String parallelWorldName)
void RegisterMe(G4HadronicInteraction *a)
#define G4endl
Definition: G4ios.hh:61
static void ConstructParticle()
G4bool IsAtRestRequired(G4ParticleDefinition *)
std::vector< G4String > fParaWorldName
static void ConstructParticle()
const G4String & GetParticleName() const
G4double GetPDGCharge() const
G4bool IsMasterThread()
Definition: G4Threading.cc:130
virtual void ConstructProcess()
B03PhysicsList(G4String)
virtual G4bool IsApplicable(const G4ParticleDefinition &) override
Definition: G4Decay.cc:89
void SetMinEForMultiFrag(G4double anE)
virtual void ConstructGeneral()
void SetProcessOrdering(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt, G4int ordDoIt=ordDefault)
void SetProcessOrderingToSecond(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt)
static constexpr double TeV
Definition: G4SIunits.hh:218
void ConstructAllBaryons()
G4ParticleTable::G4PTblDicIterator * GetParticleIterator() const
bool G4bool
Definition: G4Types.hh:79
G4int AddProcess(G4VProcess *aProcess, G4int ordAtRestDoIt=ordInActive, G4int ordAlongSteptDoIt=ordInActive, G4int ordPostStepDoIt=ordInActive)
void SetDeExcitation(G4VPreCompoundModel *ptr)
void SetTransport(G4VIntraNuclearTransportModel *const value)
virtual void ConstructEM()
virtual void ConstructHad()
void AddScoringProcess()
static void ConstructParticle()
Definition of the B03PhysicsList class.
static void ConstructParticle()
G4String fBiasWorldName
void SetParallel(G4bool paraflag)
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
virtual void AddProcess()
void SetHighEnergyGenerator(G4VHighEnergyGenerator *const value)
int G4int
Definition: G4Types.hh:78
virtual ~B03PhysicsList()
G4ProcessManager * GetProcessManager() const
virtual void ConstructLeptHad()
static G4IStore * GetInstance()
Definition: G4IStore.cc:236
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
G4GLOB_DLL std::ostream G4cout
void ConstructAllBosons()
void ConstructAllMesons()
void SetFragmentationModel(G4VStringFragmentation *aModel)
virtual void Configure()
void SetVerboseLevel(G4int value)
virtual void SetCuts()
static constexpr double GeV
Definition: G4SIunits.hh:217
void ConstructAllShortLiveds()