00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 #include "BaBar.hh"
00024 #include <math.h>
00025 #include "BgsPhysicsList.hh"
00026 #include <iomanip>
00027
00028
00029
00030
00031
00032 #include "BgsGenocide.hh"
00033 #include "BgsGentleGenocide.hh"
00034
00035
00036
00037
00038 #include "G4ParticleDefinition.hh"
00039 #include "G4ParticleWithCuts.hh"
00040 #include "G4ProcessManager.hh"
00041 #include "G4ProcessVector.hh"
00042 #include "G4ParticleTypes.hh"
00043 #include "G4ParticleTable.hh"
00044 #include "G4FastSimulationManagerProcess.hh"
00045 #include "G4ComptonScattering.hh"
00046 #include "G4GammaConversion.hh"
00047 #include "G4PhotoElectricEffect.hh"
00048 #include "G4VMultipleScattering.hh"
00049 #include "G4MultipleScattering.hh"
00050 #include "G4eIonisation.hh"
00051 #include "G4eBremsstrahlung.hh"
00052 #include "G4eplusAnnihilation.hh"
00053 #include "G4MuIonisation.hh"
00054 #include "G4MuBremsstrahlung.hh"
00055 #include "G4MuPairProduction.hh"
00056 #include "G4hIonisation.hh"
00057 #include "G4Decay.hh"
00058 #include "G4Material.hh"
00059 #include "G4ProductionCuts.hh"
00060 #include "G4ProductionCutsTable.hh"
00061 #include "G4MaterialCutsCouple.hh"
00062
00063
00064
00065 #include "G4GammaNuclearReaction.hh"
00066 #include "G4ElectroNuclearReaction.hh"
00067 #include "G4TheoFSGenerator.hh"
00068 #include "G4GeneratorPrecompoundInterface.hh"
00069 #include "G4QGSModel.hh"
00070 #include "G4GammaParticipants.hh"
00071 #include "G4QGSMFragmentation.hh"
00072 #include "G4ExcitedStringDecay.hh"
00073 #include "G4PhotoNuclearProcess.hh"
00074 #include "G4ElectronNuclearProcess.hh"
00075 #include "G4PositronNuclearProcess.hh"
00076
00077
00078
00079 #include "G4HadronElasticProcess.hh"
00080 #include "G4HadronCaptureProcess.hh"
00081 #include "G4HadronFissionProcess.hh"
00082 #include "G4PionPlusInelasticProcess.hh"
00083 #include "G4PionMinusInelasticProcess.hh"
00084 #include "G4KaonPlusInelasticProcess.hh"
00085 #include "G4KaonMinusInelasticProcess.hh"
00086 #include "G4KaonZeroLInelasticProcess.hh"
00087 #include "G4KaonZeroSInelasticProcess.hh"
00088 #include "G4ProtonInelasticProcess.hh"
00089 #include "G4AntiProtonInelasticProcess.hh"
00090 #include "G4NeutronInelasticProcess.hh"
00091 #include "G4AntiNeutronInelasticProcess.hh"
00092 #include "G4LambdaInelasticProcess.hh"
00093 #include "G4AntiLambdaInelasticProcess.hh"
00094 #include "G4DeuteronInelasticProcess.hh"
00095 #include "G4TritonInelasticProcess.hh"
00096 #include "G4AlphaInelasticProcess.hh"
00097 #include "G4ShortLivedConstructor.hh"
00098
00099
00100
00101
00102 #include "G4LElastic.hh"
00103 #include "G4LEAntiProtonInelastic.hh"
00104 #include "G4LEAntiNeutronInelastic.hh"
00105 #include "G4LEAntiLambdaInelastic.hh"
00106 #include "G4LEDeuteronInelastic.hh"
00107 #include "G4LETritonInelastic.hh"
00108 #include "G4LEAlphaInelastic.hh"
00109
00110 #include "G4NeutronHPElastic.hh"
00111 #include "G4NeutronHPElasticData.hh"
00112 #include "G4NeutronHPInelastic.hh"
00113 #include "G4NeutronHPInelasticData.hh"
00114 #include "G4NeutronHPFission.hh"
00115 #include "G4NeutronHPFissionData.hh"
00116 #include "G4NeutronHPCapture.hh"
00117 #include "G4NeutronHPCaptureData.hh"
00118 #include "G4LFission.hh"
00119 #include "G4LCapture.hh"
00120
00121 #include "G4CascadeInterface.hh"
00122 #include "CLHEP/Random/RanecuEngine.h"
00123
00124
00125
00126 #include "G4PiNuclearCrossSection.hh"
00127 using std::ostream;
00128
00129
00130
00131 BgsPhysicsList::BgsPhysicsList() : G4VUserPhysicsList(),
00132
00133
00134 first(false)
00135
00136
00137 {
00138 SetVerboseLevel(2);
00139 bertini_model = new G4CascadeInterface;
00140 }
00141
00142
00143 BgsPhysicsList::~BgsPhysicsList()
00144 {
00145
00146 }
00147
00148 void
00149 BgsPhysicsList::ConstructParticle()
00150 {
00151
00152
00153
00154
00155
00156
00157 ConstructBosons();
00158 ConstructLeptons();
00159 ConstructMesons();
00160 ConstructBaryons();
00161 ConstructIons();
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192 G4ShortLivedConstructor shortLived;
00193 shortLived.ConstructParticle();
00194 }
00195
00196 void
00197 BgsPhysicsList::ConstructBosons()
00198 {
00199
00200 G4Geantino::GeantinoDefinition();
00201 G4ChargedGeantino::ChargedGeantinoDefinition();
00202
00203
00204 G4Gamma::GammaDefinition();
00205
00206
00207 G4OpticalPhoton::OpticalPhotonDefinition();
00208 }
00209
00210 void
00211 BgsPhysicsList::ConstructLeptons()
00212 {
00213
00214 G4Electron::ElectronDefinition();
00215 G4Positron::PositronDefinition();
00216 G4MuonPlus::MuonPlusDefinition();
00217 G4MuonMinus::MuonMinusDefinition();
00218 G4TauPlus::TauPlusDefinition();
00219 G4TauMinus::TauMinusDefinition();
00220
00221 G4NeutrinoE::NeutrinoEDefinition();
00222 G4AntiNeutrinoE::AntiNeutrinoEDefinition();
00223 G4NeutrinoMu::NeutrinoMuDefinition();
00224 G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
00225 G4NeutrinoTau::NeutrinoTauDefinition();
00226 G4AntiNeutrinoTau::AntiNeutrinoTauDefinition();
00227 }
00228
00229 void
00230 BgsPhysicsList::ConstructMesons()
00231 {
00232
00233 G4PionPlus ::PionPlusDefinition();
00234 G4PionMinus ::PionMinusDefinition();
00235 G4PionZero ::PionZeroDefinition();
00236 G4Eta ::EtaDefinition();
00237 G4EtaPrime ::EtaPrimeDefinition();
00238
00239 G4KaonPlus ::KaonPlusDefinition();
00240 G4KaonMinus ::KaonMinusDefinition();
00241 G4KaonZero ::KaonZeroDefinition();
00242 G4AntiKaonZero ::AntiKaonZeroDefinition();
00243 G4KaonZeroLong ::KaonZeroLongDefinition();
00244 G4KaonZeroShort::KaonZeroShortDefinition();
00245 }
00246
00247 void
00248 BgsPhysicsList::ConstructBaryons()
00249 {
00250
00251 G4Proton ::ProtonDefinition();
00252 G4AntiProton ::AntiProtonDefinition();
00253 G4Neutron ::NeutronDefinition();
00254 G4AntiNeutron ::AntiNeutronDefinition();
00255 G4Lambda ::LambdaDefinition();
00256 G4AntiLambda ::AntiLambdaDefinition();
00257 G4SigmaPlus ::SigmaPlusDefinition();
00258 G4SigmaZero ::SigmaZeroDefinition();
00259 G4SigmaMinus ::SigmaMinusDefinition();
00260 G4AntiSigmaPlus ::AntiSigmaPlusDefinition();
00261 G4AntiSigmaZero ::AntiSigmaZeroDefinition();
00262 G4AntiSigmaMinus::AntiSigmaMinusDefinition();
00263 G4XiZero ::XiZeroDefinition();
00264 G4XiMinus ::XiMinusDefinition();
00265 G4AntiXiZero ::AntiXiZeroDefinition();
00266 G4AntiXiMinus ::AntiXiMinusDefinition();
00267 G4OmegaMinus ::OmegaMinusDefinition();
00268 G4AntiOmegaMinus::AntiOmegaMinusDefinition();
00269 G4Deuteron ::DeuteronDefinition();
00270 G4Triton ::TritonDefinition();
00271 G4He3 ::He3Definition();
00272 G4Alpha ::AlphaDefinition();
00273 }
00274
00275 void
00276 BgsPhysicsList::ConstructIons()
00277 {
00278 G4GenericIon::GenericIonDefinition();
00279 }
00280
00281 void
00282 BgsPhysicsList::ConstructProcess()
00283 {
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309 AddTransportation();
00310 std::cout << "+^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^+ " <<
00311 std::endl;
00312 std::cout << "+v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v+ " <<
00313 std::endl;
00314 std::cout << " +--- ---+ " <<
00315 std::endl;
00316 std::cout << " +--- G4Transportation ---+ " <<
00317 std::endl;
00318 std::cout << " +--- USED !! ---+ " <<
00319 std::endl;
00320 std::cout << " +--- ---+ " <<
00321 std::endl;
00322 std::cout << "+^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^+ " <<
00323 std::endl;
00324 std::cout << "+v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v^v+ " <<
00325 std::endl;
00326
00327
00328 ConstructEM();
00329 ConstructLeptHad();
00330 ConstructHad();
00331 ConstructGeneral();
00332 ConstructNeutrinoGenocide();
00333 ConstructIonFix();
00334 ConstructNeutronFix();
00335 }
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374
00375
00376 void
00377 BgsPhysicsList::ConstructEM()
00378 {
00379 theParticleIterator->reset();
00380 while( (*theParticleIterator)() ){
00381 G4ParticleDefinition* particle = theParticleIterator->value();
00382 G4ProcessManager* pmanager = particle->GetProcessManager();
00383 G4String particleName = particle->GetParticleName();
00384
00385 if (particleName == "gamma") {
00386
00387
00388
00389 pmanager->AddDiscreteProcess( new G4PhotoElectricEffect());
00390 pmanager->AddDiscreteProcess( new G4ComptonScattering());
00391 pmanager->AddDiscreteProcess( new G4GammaConversion());
00392
00393 } else if (particleName == "e-") {
00394
00395
00396
00397 G4MultipleScattering* ms = new G4MultipleScattering();
00398
00399 ms->SetRangeFactor(0.02);
00400
00401
00402
00403 G4eIonisation *ionizationProcess = new G4eIonisation();
00404
00405
00406 pmanager->AddProcess( ms, -1, 1,1);
00407 pmanager->AddProcess( ionizationProcess, -1, 2,2);
00408 pmanager->AddProcess( new G4eBremsstrahlung(), -1,-1,3);
00409
00410 } else if (particleName == "e+") {
00411
00412
00413
00414
00415 G4MultipleScattering* ms = new G4MultipleScattering();
00416
00417 ms->SetRangeFactor(0.02);
00418
00419
00420
00421 G4eIonisation *ionizationProcess = new G4eIonisation();
00422
00423
00424 pmanager->AddProcess( ms, -1, 1,1);
00425 pmanager->AddProcess( ionizationProcess, -1, 2,2);
00426 pmanager->AddProcess( new G4eBremsstrahlung(), -1,-1,3);
00427 pmanager->AddProcess( new G4eplusAnnihilation(), 0,-1,4);
00428
00429 } else if( particleName == "mu+" ||
00430 particleName == "mu-" ) {
00431
00432
00433
00434
00435 G4MultipleScattering* ms = new G4MultipleScattering();
00436
00437 ms->SetRangeFactor(0.02);
00438
00439
00440
00441 G4MuIonisation *ionizationProcess = new G4MuIonisation();
00442
00443
00444 pmanager->AddProcess( ms, -1, 1,1);
00445 pmanager->AddProcess( ionizationProcess, -1, 2,2);
00446 pmanager->AddProcess( new G4MuBremsstrahlung(), -1,-1,3);
00447 pmanager->AddProcess( new G4MuPairProduction(), -1,-1,4);
00448
00449 } else if ((!particle->IsShortLived()) &&
00450 (particle->GetPDGCharge() != 0.0) &&
00451 (particle->GetParticleName() != "chargedgeantino")) {
00452
00453 G4int AP=1;
00454 if (particle->GetParticleName() == "GenericIon") {
00455
00456 std::cerr <<"*********************************************************************" <<std::endl;
00457 std::cerr << "*** Disabling G4MultipleScattering process for particle " <<particle->GetParticleName() << std::endl;
00458 std::cerr <<"*********************************************************************" <<std::endl;
00459 } else {
00460
00461 G4MultipleScattering* ms = new G4MultipleScattering();
00462
00463 ms->SetRangeFactor(0.02);
00464
00465
00466 pmanager->AddProcess( ms, -1,AP,AP);
00467 AP++;
00468 }
00469 G4hIonisation *ionizationProcess = new G4hIonisation();
00470
00471
00472 pmanager->AddProcess( ionizationProcess, -1,AP,AP);
00473 }
00474 }
00475 }
00476
00477
00478
00479
00480
00481
00482
00483
00484
00485
00486
00487
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499
00500
00501
00502
00503
00504
00505
00506
00507
00508
00509
00510
00511
00512
00513
00514
00515
00516 void
00517 BgsPhysicsList::SetStatusEmProcess()
00518 {
00519 }
00520
00521
00522
00523
00524
00525
00526
00527
00528
00529
00530
00531
00532
00533
00534
00535
00536
00537
00538 void
00539 BgsPhysicsList::ConstructGeneral()
00540 {
00541
00542 G4Decay* theDecayProcess = new G4Decay();
00543 theParticleIterator->reset();
00544 while( (*theParticleIterator)() ){
00545 G4ParticleDefinition* particle = theParticleIterator->value();
00546 G4ProcessManager* pmanager = particle->GetProcessManager();
00547 if (theDecayProcess->IsApplicable(*particle)) {
00548 pmanager ->AddProcess( theDecayProcess );
00549
00550 pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
00551 pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
00552
00553
00554 }
00555 }
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567
00568
00569
00570
00571
00572
00573
00574
00575
00576
00577
00578
00579
00580
00581
00582
00583
00584
00585
00586
00587
00588
00589
00590
00591
00592
00593
00594
00595
00596
00597
00598
00599
00600
00601
00602
00603
00604
00605
00606 }
00607
00608
00609
00610 void
00611 BgsPhysicsList::ConstructLeptHad()
00612 {
00613
00614
00615
00616
00617
00618 G4GammaNuclearReaction* lowEGammaModel = new G4GammaNuclearReaction();
00619 lowEGammaModel->SetMaxEnergy(3.5*GeV);
00620 G4PhotoNuclearProcess* thePhotoNuclearProcess = new
00621 G4PhotoNuclearProcess();
00622 thePhotoNuclearProcess->RegisterMe(lowEGammaModel);
00623
00624
00625
00626
00627
00628
00629
00630
00631
00632
00633
00634
00635
00636
00637 G4TheoFSGenerator* highEGammaModel = new G4TheoFSGenerator();
00638 G4GeneratorPrecompoundInterface* preComModel =
00639 new G4GeneratorPrecompoundInterface();
00640 highEGammaModel->SetTransport(preComModel);
00641
00642 G4QGSModel<G4GammaParticipants>* theStringModel =
00643 new G4QGSModel<G4GammaParticipants>;
00644 G4QGSMFragmentation* fragModel = new G4QGSMFragmentation();
00645 G4ExcitedStringDecay* stringDecay =
00646 new G4ExcitedStringDecay(fragModel);
00647 theStringModel->SetFragmentationModel(stringDecay);
00648
00649 highEGammaModel->SetHighEnergyGenerator(theStringModel);
00650 highEGammaModel->SetMinEnergy(3.*GeV);
00651 highEGammaModel->SetMaxEnergy(20.*GeV);
00652
00653 thePhotoNuclearProcess->RegisterMe(highEGammaModel);
00654
00655 G4ProcessManager* gamMan = G4Gamma::Gamma()->GetProcessManager();
00656
00657 gamMan->AddDiscreteProcess(thePhotoNuclearProcess);
00658
00659
00660
00661 G4ElectroNuclearReaction* theElectronReaction =
00662 new G4ElectroNuclearReaction();
00663 G4ElectronNuclearProcess* theElectronNuclearProcess =
00664 new G4ElectronNuclearProcess();
00665 theElectronNuclearProcess->RegisterMe(theElectronReaction);
00666
00667 G4ProcessManager* electronMan =
00668 G4Electron::Electron()->GetProcessManager();
00669
00670 electronMan->AddProcess(theElectronNuclearProcess, -1, -1, 4);
00671
00672
00673
00674
00675
00676
00677
00678
00679
00680
00681
00682
00683
00684
00685
00686
00687
00688
00689 G4PositronNuclearProcess* thePositronNuclearProcess =
00690 new G4PositronNuclearProcess();
00691 thePositronNuclearProcess->RegisterMe(theElectronReaction);
00692
00693 G4ProcessManager* positronMan =
00694 G4Positron::Positron()->GetProcessManager();
00695 positronMan->AddProcess(thePositronNuclearProcess, -1, -1, 5);
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707 }
00708
00709 void
00710 BgsPhysicsList::ConstructHad()
00711 {
00712
00713
00714
00715
00716 G4PiNuclearCrossSection* piNucCS = new G4PiNuclearCrossSection();
00717
00718 theParticleIterator->reset();
00719 while( (*theParticleIterator)() ){
00720 G4ParticleDefinition* particle = theParticleIterator->value();
00721 G4ProcessManager* pmanager = particle->GetProcessManager();
00722 G4String particleName = particle->GetParticleName();
00723
00724
00725
00726 if (particleName == "pi+") {
00727
00728
00729 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00730 G4LElastic *model = new G4LElastic();
00731 process->RegisterMe(model);
00732 pmanager->AddDiscreteProcess(process);
00733
00734
00735
00736 G4PionPlusInelasticProcess *inel_process = new
00737 G4PionPlusInelasticProcess();
00738 inel_process->AddDataSet(piNucCS);
00739 inel_process->RegisterMe(bertini_model);
00740 pmanager->AddDiscreteProcess(inel_process);
00741
00742
00743
00744
00745
00746 } else if (particleName == "pi-") {
00747
00748
00749 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00750 G4LElastic *model = new G4LElastic();
00751 process->RegisterMe(model);
00752 pmanager->AddDiscreteProcess(process);
00753
00754
00755
00756 G4PionMinusInelasticProcess *inel_process = new
00757 G4PionMinusInelasticProcess();
00758 inel_process->AddDataSet(piNucCS);
00759 inel_process->RegisterMe(bertini_model);
00760 pmanager->AddDiscreteProcess(inel_process);
00761
00762
00763
00764
00765
00766 } else if (particleName == "kaon+") {
00767
00768
00769 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00770 G4LElastic *model = new G4LElastic();
00771 process->RegisterMe(model);
00772 pmanager->AddDiscreteProcess(process);
00773
00774
00775
00776 G4KaonPlusInelasticProcess* inel_process = new
00777 G4KaonPlusInelasticProcess();
00778 inel_process->RegisterMe(bertini_model);
00779 pmanager->AddDiscreteProcess(inel_process);
00780
00781
00782
00783
00784
00785 } else if (particleName == "kaon-") {
00786
00787
00788 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00789 G4LElastic *model = new G4LElastic();
00790 process->RegisterMe(model);
00791 pmanager->AddDiscreteProcess(process);
00792
00793
00794
00795 G4KaonMinusInelasticProcess* inel_process = new
00796 G4KaonMinusInelasticProcess();
00797 inel_process->RegisterMe(bertini_model);
00798 pmanager->AddDiscreteProcess(inel_process);
00799
00800
00801
00802
00803
00804 } else if (particleName == "kaon0L") {
00805
00806
00807 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00808 G4LElastic *model = new G4LElastic();
00809 process->RegisterMe(model);
00810 pmanager->AddDiscreteProcess(process);
00811
00812
00813
00814 G4KaonZeroLInelasticProcess* inel_process = new
00815 G4KaonZeroLInelasticProcess();
00816 inel_process->RegisterMe(bertini_model);
00817 pmanager->AddDiscreteProcess(inel_process);
00818
00819
00820
00821
00822
00823 } else if (particleName == "kaon0S") {
00824
00825
00826 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00827 G4LElastic *model = new G4LElastic();
00828 process->RegisterMe(model);
00829 pmanager->AddDiscreteProcess(process);
00830
00831
00832
00833 G4KaonZeroSInelasticProcess* inel_process =
00834 new G4KaonZeroSInelasticProcess();
00835 inel_process->RegisterMe(bertini_model);
00836 pmanager->AddDiscreteProcess(inel_process);
00837
00838
00839
00840
00841
00842 } else if (particleName == "proton") {
00843
00844
00845 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00846 G4LElastic *model = new G4LElastic();
00847 process->RegisterMe(model);
00848 pmanager->AddDiscreteProcess(process);
00849
00850
00851
00852 G4ProtonInelasticProcess *inel_process = new
00853 G4ProtonInelasticProcess();
00854 inel_process->RegisterMe(bertini_model);
00855 pmanager->AddDiscreteProcess(inel_process);
00856
00857
00858
00859
00860
00861 } else if (particleName == "anti_proton") {
00862
00863
00864 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00865 G4LElastic *model = new G4LElastic();
00866 process->RegisterMe(model);
00867 pmanager->AddDiscreteProcess(process);
00868
00869
00870
00871 G4AntiProtonInelasticProcess *inel_process =
00872 new G4AntiProtonInelasticProcess();
00873 G4LEAntiProtonInelastic *inel_model = new G4LEAntiProtonInelastic();
00874 inel_process->RegisterMe(inel_model);
00875 pmanager->AddDiscreteProcess(inel_process);
00876
00877
00878
00879
00880
00881 } else if (particleName == "neutron") {
00882
00883
00884 if(1){
00885 G4cout << "High precision neutron models chosen" << G4endl;
00886
00887 putenv("G4NEUTRONHPDATA=/afs/ihep.ac.cn/bes3/offline/sw/packages/geant4/4.9.0/slc4_ia32_gcc346/geant4.9.0.p01/data/G4NDL3.11/");
00888
00889
00890 G4HadronElasticProcess* el_process = new G4HadronElasticProcess();
00891
00892
00893 G4NeutronHPElastic* hpel_model = new G4NeutronHPElastic();
00894 G4NeutronHPElasticData* el_data = new G4NeutronHPElasticData();
00895 el_process->AddDataSet(el_data);
00896 el_process->RegisterMe(hpel_model);
00897
00898
00899 G4LElastic* el_model = new G4LElastic();
00900 el_model->SetMinEnergy(19.9*MeV);
00901 el_process->RegisterMe(el_model);
00902
00903 pmanager->AddDiscreteProcess(el_process);
00904
00905
00906
00907
00908 G4NeutronInelasticProcess* inel_process =
00909 new G4NeutronInelasticProcess();
00910
00911
00912 G4NeutronHPInelastic* hpinel_model = new G4NeutronHPInelastic();
00913 G4NeutronHPInelasticData* hpinel_data = new
00914 G4NeutronHPInelasticData();
00915 inel_process->AddDataSet(hpinel_data);
00916 inel_process->RegisterMe(hpinel_model);
00917
00918
00919 G4CascadeInterface* neutron_bertini = new G4CascadeInterface;
00920 neutron_bertini->SetMinEnergy(19.9*MeV);
00921 inel_process->RegisterMe(neutron_bertini);
00922
00923 pmanager->AddDiscreteProcess(inel_process);
00924
00925
00926
00927
00928 G4HadronCaptureProcess* cap_process = new G4HadronCaptureProcess();
00929
00930
00931 G4NeutronHPCapture* hpcap_model = new G4NeutronHPCapture();
00932 G4NeutronHPCaptureData* hpcap_data = new G4NeutronHPCaptureData();
00933 cap_process->AddDataSet(hpcap_data);
00934 cap_process->RegisterMe(hpcap_model);
00935
00936
00937
00938 G4LCapture* cap_model = new G4LCapture();
00939 cap_model->SetMinEnergy(19.9*MeV);
00940 cap_process->RegisterMe(cap_model);
00941
00942 pmanager->AddDiscreteProcess(cap_process);
00943
00944
00945
00946
00947
00948 G4HadronFissionProcess* fis_process = new G4HadronFissionProcess();
00949
00950
00951 G4NeutronHPFission* hpfis_model = new G4NeutronHPFission();
00952 G4NeutronHPFissionData* hpfis_data = new G4NeutronHPFissionData();
00953 fis_process->AddDataSet(hpfis_data);
00954 fis_process->RegisterMe(hpfis_model);
00955
00956
00957
00958 G4LFission* fis_model = new G4LFission();
00959 fis_model->SetMinEnergy(19.9*MeV);
00960 fis_process->RegisterMe(fis_model);
00961
00962 pmanager->AddDiscreteProcess(fis_process);
00963
00964
00965
00966
00967 } else {
00968
00969
00970 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00971 G4LElastic *model = new G4LElastic();
00972 process->RegisterMe(model);
00973 pmanager->AddDiscreteProcess(process);
00974
00975
00976 G4NeutronInelasticProcess *inel_process =
00977 new G4NeutronInelasticProcess();
00978 inel_process->RegisterMe(bertini_model);
00979 pmanager->AddDiscreteProcess(inel_process);
00980
00981
00982 }
00983
00984
00985
00986 } else if (particleName == "anti_neutron") {
00987
00988
00989 G4HadronElasticProcess *process = new G4HadronElasticProcess();
00990 G4LElastic *model = new G4LElastic();
00991 process->RegisterMe(model);
00992 pmanager->AddDiscreteProcess(process);
00993
00994
00995
00996 G4AntiNeutronInelasticProcess *inel_process =
00997 new G4AntiNeutronInelasticProcess();
00998 G4LEAntiNeutronInelastic *inel_model = new G4LEAntiNeutronInelastic();
00999 inel_process->RegisterMe(inel_model);
01000 pmanager->AddDiscreteProcess(inel_process);
01001
01002
01003
01004
01005
01006 } else if (particleName == "lambda") {
01007
01008
01009 G4HadronElasticProcess *process = new G4HadronElasticProcess();
01010 G4LElastic *model = new G4LElastic();
01011 process->RegisterMe(model);
01012 pmanager->AddDiscreteProcess(process);
01013
01014
01015
01016 G4LambdaInelasticProcess* inel_process = new
01017 G4LambdaInelasticProcess();
01018 inel_process->RegisterMe(bertini_model);
01019 pmanager->AddDiscreteProcess(inel_process);
01020
01021
01022
01023
01024
01025 } else if (particleName == "anti_lambda") {
01026
01027
01028 G4HadronElasticProcess *process = new G4HadronElasticProcess();
01029 G4LElastic *model = new G4LElastic();
01030 process->RegisterMe(model);
01031 pmanager->AddDiscreteProcess(process);
01032
01033
01034
01035 G4AntiLambdaInelasticProcess *inel_process =
01036 new G4AntiLambdaInelasticProcess();
01037 G4LEAntiLambdaInelastic *inel_model = new G4LEAntiLambdaInelastic();
01038 inel_process->RegisterMe(inel_model);
01039 pmanager->AddDiscreteProcess(inel_process);
01040
01041
01042
01043
01044
01045 } else if (particleName == "deuteron") {
01046
01047
01048 G4HadronElasticProcess *process = new G4HadronElasticProcess();
01049 G4LElastic *model = new G4LElastic();
01050 process->RegisterMe(model);
01051 pmanager->AddDiscreteProcess(process);
01052
01053
01054
01055 G4DeuteronInelasticProcess *inel_process =
01056 new G4DeuteronInelasticProcess();
01057 G4LEDeuteronInelastic *inel_model = new G4LEDeuteronInelastic();
01058 inel_process->RegisterMe(inel_model);
01059 pmanager->AddDiscreteProcess(inel_process);
01060
01061
01062
01063
01064
01065 } else if (particleName == "triton") {
01066
01067
01068 G4HadronElasticProcess *process = new G4HadronElasticProcess();
01069 G4LElastic *model = new G4LElastic();
01070 process->RegisterMe(model);
01071 pmanager->AddDiscreteProcess(process);
01072
01073
01074
01075 G4TritonInelasticProcess *inel_process =
01076 new G4TritonInelasticProcess();
01077 G4LETritonInelastic *inel_model = new G4LETritonInelastic();
01078 inel_process->RegisterMe(inel_model);
01079 pmanager->AddDiscreteProcess(inel_process);
01080
01081
01082
01083
01084
01085 } else if (particleName == "alpha") {
01086
01087
01088 G4HadronElasticProcess *process = new G4HadronElasticProcess();
01089 G4LElastic *model = new G4LElastic();
01090 process->RegisterMe(model);
01091 pmanager->AddDiscreteProcess(process);
01092
01093
01094
01095 G4AlphaInelasticProcess *inel_process = new G4AlphaInelasticProcess();
01096 G4LEAlphaInelastic *inel_model = new G4LEAlphaInelastic();
01097 inel_process->RegisterMe(inel_model);
01098 pmanager->AddDiscreteProcess(inel_process);
01099
01100
01101 }
01102 }
01103 }
01104
01105
01106
01107
01108
01109 void BgsPhysicsList::ConstructNeutrinoGenocide()
01110 {
01111 BgsGenocide *genocide = new BgsGenocide();
01112
01113 theParticleIterator->reset();
01114 while( (*theParticleIterator)() ){
01115 G4ParticleDefinition* particle = theParticleIterator->value();
01116 G4ProcessManager* pmanager = particle->GetProcessManager();
01117 G4String particleName = particle->GetParticleName();
01118
01119 if (particleName == "nu_e" ||
01120 particleName == "nu_mu" ||
01121 particleName == "nu_tau" ||
01122 particleName == "anti_nu_e" ||
01123 particleName == "anti_nu_mu" ||
01124 particleName == "anti_nu_tau" ||
01125
01126
01127
01128 particleName == "kaon0" ||
01129 particleName == "anti_kaon0" ) {
01130
01131 pmanager->AddProcess(genocide, -1, -1, 5);
01132
01133
01134 }
01135 }
01136 }
01137
01138
01139
01140
01141
01142 void BgsPhysicsList::ConstructIonFix()
01143 {
01144 BgsGenocide *genocide = new
01145
01146
01147 BgsGentleGenocide(0.01*MeV,
01148 60);
01149
01150 theParticleIterator->reset();
01151 while( (*theParticleIterator)() ){
01152 G4ParticleDefinition* particle = theParticleIterator->value();
01153 G4ProcessManager* pmanager = particle->GetProcessManager();
01154 G4String particleName = particle->GetParticleName();
01155
01156 if ( particleName == "triton" ||
01157 particleName == "alpha" ||
01158 particleName == "proton" ||
01159 particleName == "deuteron" ) {
01160
01161 pmanager->AddProcess(genocide, -1, -1, 5);
01162
01163
01164 }
01165 }
01166 }
01167
01168
01169
01170
01171
01172 void BgsPhysicsList::ConstructNeutronFix()
01173 {
01174 BgsGenocide *genocide = new
01175 BgsGentleGenocide(0.01*MeV,0);
01176
01177 theParticleIterator->reset();
01178 while( (*theParticleIterator)() ){
01179 G4ParticleDefinition* particle = theParticleIterator->value();
01180 G4ProcessManager* pmanager = particle->GetProcessManager();
01181 G4String particleName = particle->GetParticleName();
01182
01183 if ( particleName == "neutron" ) {
01184
01185 pmanager->AddProcess(genocide, -1, -1, 1);
01186
01187
01188 }
01189 }
01190 }
01191
01192
01193
01194
01195
01196 void BgsPhysicsList::SetCuts()
01197 {
01198
01199
01200 SetDefaultCutValue(0.7*mm);
01201 SetCutsWithDefault();
01202
01203
01204
01205
01206
01207 }