BgsPhysicsList Class Reference

#include <BgsPhysicsList.hh>

List of all members.

Public Member Functions

 BgsPhysicsList ()
virtual ~BgsPhysicsList ()
void SetCuts ()
void SetStatusEmProcess ()

Protected Member Functions

void ConstructParticle ()
void ConstructProcess ()
void ConstructBosons ()
void ConstructLeptons ()
void ConstructMesons ()
void ConstructBaryons ()
void ConstructIons ()
void ConstructGeneral ()
void ConstructEM ()
void ConstructLeptHad ()
void ConstructHad ()
void ConstructNeutrinoGenocide ()
void ConstructIonFix ()
void ConstructNeutronFix ()

Private Attributes

G4bool first
G4CascadeInterface * bertini_model


Detailed Description

Definition at line 55 of file BgsPhysicsList.hh.


Constructor & Destructor Documentation

BgsPhysicsList::BgsPhysicsList (  ) 

Definition at line 131 of file BgsPhysicsList.cc.

References bertini_model.

00131                                : G4VUserPhysicsList(),//bes
00132 //    control(theControl),
00133 //    first(true),
00134     first(false) //bes
00135 //    physicsRegistrar(pr),
00136 //    theLooperDeath(0)                 // looperdeath process
00137 {
00138   SetVerboseLevel(2);
00139   bertini_model = new G4CascadeInterface;
00140 }

BgsPhysicsList::~BgsPhysicsList (  )  [virtual]

Definition at line 143 of file BgsPhysicsList.cc.

00144 {
00145 //  delete theLooperDeath;
00146 }


Member Function Documentation

void BgsPhysicsList::ConstructBaryons (  )  [protected]

Definition at line 248 of file BgsPhysicsList.cc.

Referenced by ConstructParticle().

00249 {
00250   //  baryons
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 }

void BgsPhysicsList::ConstructBosons (  )  [protected]

Definition at line 197 of file BgsPhysicsList.cc.

Referenced by ConstructParticle().

00198 {
00199   // pseudo-particles
00200   G4Geantino::GeantinoDefinition();
00201   G4ChargedGeantino::ChargedGeantinoDefinition();
00202 
00203   // gamma
00204   G4Gamma::GammaDefinition();
00205 
00206   // optical photon
00207   G4OpticalPhoton::OpticalPhotonDefinition();
00208 }

void BgsPhysicsList::ConstructEM (  )  [protected]

Definition at line 377 of file BgsPhysicsList.cc.

Referenced by ConstructProcess().

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       // gamma
00387       // Construct processes for gamma
00388 
00389       pmanager->AddDiscreteProcess( new G4PhotoElectricEffect());
00390       pmanager->AddDiscreteProcess( new G4ComptonScattering());
00391       pmanager->AddDiscreteProcess( new G4GammaConversion());
00392 
00393     } else if (particleName == "e-") {
00394       //electron
00395       // Construct processes for electron
00396       //change for geant4.9.0.p01
00397       G4MultipleScattering* ms = new G4MultipleScattering();
00398       //ms->MscStepLimitation(false,0.02);
00399       ms->SetRangeFactor(0.02);
00400 
00401       //      ms->SetLateralDisplasmentFlag(false);
00402       
00403       G4eIonisation *ionizationProcess = new G4eIonisation();
00404       //      ionizationProcess->SetLinearLossLimit(1.0);
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       //positron
00412       // Construct processes for positron
00413 
00414       //change for geant4.9.p01
00415       G4MultipleScattering* ms = new G4MultipleScattering();
00416       //ms->MscStepLimitation(false,0.02);
00417       ms->SetRangeFactor(0.02);
00418 
00419       //      ms->SetLateralDisplasmentFlag(false);
00420       
00421       G4eIonisation *ionizationProcess = new G4eIonisation();
00422       //      ionizationProcess->SetLinearLossLimit(1.0);
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       //muon  
00432       // Construct processes for muon+
00433 
00434       //change for geant4.9.0.p01
00435       G4MultipleScattering* ms = new G4MultipleScattering();
00436       //ms->MscStepLimitation(false,0.02);
00437       ms->SetRangeFactor(0.02);
00438 
00439      //      ms->SetLateralDisplasmentFlag(false);
00440       
00441       G4MuIonisation *ionizationProcess = new G4MuIonisation();
00442       //      ionizationProcess->SetLinearLossLimit(1.0);
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       // all others charged particles except geantino
00453       G4int AP=1;
00454       if (particle->GetParticleName() == "GenericIon") {
00455         //ostream& o = ErrMsg(warning);
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   //change for Geant4.9.0.p01
00461         G4MultipleScattering* ms = new G4MultipleScattering();
00462   //ms->MscStepLimitation(false,0.02);
00463   ms->SetRangeFactor(0.02);
00464 
00465         //              ms->SetLateralDisplasmentFlag(false);
00466         pmanager->AddProcess( ms, -1,AP,AP);
00467         AP++;
00468       }
00469       G4hIonisation *ionizationProcess = new G4hIonisation();
00470       //      ionizationProcess->SetLinearLossLimit(1.0); 
00471       
00472       pmanager->AddProcess( ionizationProcess, -1,AP,AP);
00473     }
00474   }
00475 }

void BgsPhysicsList::ConstructGeneral (  )  [protected]

Definition at line 539 of file BgsPhysicsList.cc.

Referenced by ConstructProcess().

00540 {
00541   // Add Decay Process
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       // set ordering for PostStepDoIt and AtRestDoIt
00550       pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
00551       pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
00552 //      physicsRegistrar->
00553 //      Register( theDecayProcess, pmanager, "dcay", GVertex::decay );
00554     }
00555   }
00556 /*
00557   if (control->GetLooperCut()>0) {
00558 
00559     // Set special process to kill loopers
00560     theLooperDeath = new BgsLooperDeath( control->GetLooperCut()*MeV );
00561     theParticleIterator->reset();
00562     while( (*theParticleIterator)() ){
00563       G4ParticleDefinition* particle = theParticleIterator->value();
00564       if (theLooperDeath->IsApplicable(*particle)) { 
00565         G4ProcessManager* pmanager = particle->GetProcessManager();
00566         pmanager->AddProcess(theLooperDeath, -1, -1, 5);
00567         physicsRegistrar->
00568           Register( theLooperDeath, pmanager, "loop", GVertex::looperDeath
00569 );
00570       }
00571     }
00572     ErrMsg(warning) << "Loopers with pt < " << control->GetLooperCut() 
00573                     << " MeV will be killed" << endmsg;
00574   } 
00575 
00576   if (control->GetMaxNumberOfSteps()>0) {
00577 
00578     //
00579     // Set special process to kill runaway particles
00580     // Only needed if dE/dx is turned off!
00581     //
00582     // Do not abuse!
00583     //
00584     theStepDeath = new BgsChargedStepDeath( control->GetMaxNumberOfSteps()
00585 );
00586     
00587     theParticleIterator->reset();
00588     while( (*theParticleIterator)() ){
00589       G4ParticleDefinition* particle = theParticleIterator->value();
00590       if (theStepDeath->IsApplicable(*particle)) { 
00591         G4ProcessManager* pmanager = particle->GetProcessManager();
00592         pmanager->AddProcess(theStepDeath, -1, -1, 5);
00593         physicsRegistrar->
00594           Register( theStepDeath, pmanager, "maxStep", GVertex::runAway );
00595       }
00596     }
00597     ErrMsg(warning) 
00598       << "\n  Charged particles will be killed if they take more than " 
00599       << control->GetMaxNumberOfSteps() << " steps.\n"
00600       << "  If you do not understand this message, you should be very
00601 concerned.\n" 
00602       << "  If this message appears in production, you should be very
00603 upset." << endmsg;
00604   }
00605  */ 
00606 }

void BgsPhysicsList::ConstructHad (  )  [protected]

Definition at line 710 of file BgsPhysicsList.cc.

References bertini_model.

Referenced by ConstructProcess().

00711 {
00712   // One process handles hadronic elastic processes for all hadrons.
00713   // However hadronic inelastic processes are unique to each hadron.
00714 
00715   // For pi+ and pi- only, substitute pi-Nuclear cross sections 
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     // * pi+ *
00725     // *******
00726     if (particleName == "pi+") {
00727 
00728       // Elastic process
00729       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00730       G4LElastic *model =  new G4LElastic();
00731       process->RegisterMe(model);
00732       pmanager->AddDiscreteProcess(process);
00733 
00734       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00742 //                               GVertex::hadronInelastic);
00743     // *******
00744     // * pi- *
00745     // *******
00746     } else if (particleName == "pi-") {
00747 
00748       // Elastic process
00749       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00750       G4LElastic *model =  new G4LElastic();
00751       process->RegisterMe(model);
00752       pmanager->AddDiscreteProcess(process);
00753 
00754       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00762 //                               GVertex::hadronInelastic);
00763     // *******
00764     // * K+  *
00765     // *******
00766     } else if (particleName == "kaon+") {
00767 
00768       // Elastic process
00769       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00770       G4LElastic *model =  new G4LElastic();
00771       process->RegisterMe(model);
00772       pmanager->AddDiscreteProcess(process);
00773 
00774       // Inelastic process
00775 
00776       G4KaonPlusInelasticProcess* inel_process = new
00777 G4KaonPlusInelasticProcess();
00778       inel_process->RegisterMe(bertini_model);
00779       pmanager->AddDiscreteProcess(inel_process);
00780 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00781 //                               GVertex::hadronInelastic);
00782     // *******
00783     // * K-  *
00784     // *******
00785     } else if (particleName == "kaon-") {
00786 
00787       // Elastic process
00788       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00789       G4LElastic *model =  new G4LElastic();
00790       process->RegisterMe(model);
00791       pmanager->AddDiscreteProcess(process);
00792 
00793       // Inelastic process
00794 
00795       G4KaonMinusInelasticProcess* inel_process = new
00796 G4KaonMinusInelasticProcess();
00797       inel_process->RegisterMe(bertini_model);
00798       pmanager->AddDiscreteProcess(inel_process);
00799 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00800 //                               GVertex::hadronInelastic);
00801     // *******
00802     // * K0L *
00803     // *******
00804     } else if (particleName == "kaon0L") {
00805 
00806       // Elastic process
00807       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00808       G4LElastic *model =  new G4LElastic();
00809       process->RegisterMe(model);
00810       pmanager->AddDiscreteProcess(process);
00811 
00812       // Inelastic process
00813 
00814       G4KaonZeroLInelasticProcess* inel_process = new
00815 G4KaonZeroLInelasticProcess();
00816       inel_process->RegisterMe(bertini_model);
00817       pmanager->AddDiscreteProcess(inel_process);
00818 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00819 //                               GVertex::hadronInelastic);
00820     // *******
00821     // * K0S *
00822     // *******
00823     } else if (particleName == "kaon0S") {
00824 
00825       // Elastic process
00826       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00827       G4LElastic *model =  new G4LElastic();
00828       process->RegisterMe(model);
00829       pmanager->AddDiscreteProcess(process);
00830 
00831       // Inelastic process
00832 
00833       G4KaonZeroSInelasticProcess* inel_process = 
00834       new G4KaonZeroSInelasticProcess();
00835       inel_process->RegisterMe(bertini_model);
00836       pmanager->AddDiscreteProcess(inel_process);
00837 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00838 //                               GVertex::hadronInelastic);
00839     // *******
00840     // * p   *
00841     // *******
00842     } else if (particleName == "proton") {
00843 
00844       // Elastic process
00845       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00846       G4LElastic *model =  new G4LElastic();
00847       process->RegisterMe(model);
00848       pmanager->AddDiscreteProcess(process);
00849   
00850       // Inelastic process
00851 
00852       G4ProtonInelasticProcess *inel_process = new
00853 G4ProtonInelasticProcess();
00854       inel_process->RegisterMe(bertini_model);
00855       pmanager->AddDiscreteProcess(inel_process);
00856 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00857 //                               GVertex::hadronInelastic);
00858     // *********
00859     // * p-bar *
00860     // *********
00861     } else if (particleName == "anti_proton") {
00862 
00863       // Elastic process
00864       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00865       G4LElastic *model =  new G4LElastic();
00866       process->RegisterMe(model);
00867       pmanager->AddDiscreteProcess(process);
00868 
00869       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
00877 //                               GVertex::hadronInelastic);
00878     // *******
00879     // * n   *
00880     // *******
00881     } else if (particleName == "neutron") {
00882 
00883       //if (control->UseHPNeutrons()) {
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         // Elastic process
00890         G4HadronElasticProcess* el_process = new G4HadronElasticProcess();
00891 
00892         // High precision model and data below 20 MeV
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         // LEP model above 20 MeV
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         //physicsRegistrar->Register(el_process,pmanager,"hade",
00905         //                         GVertex::hadronElastic);
00906 
00907         // Inelastic process
00908         G4NeutronInelasticProcess* inel_process = 
00909                                        new G4NeutronInelasticProcess();
00910 
00911         // High precision model and data below 20 MeV
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         // Bertini model above 20 MeV
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         //physicsRegistrar->Register(inel_process,pmanager,"hadi",
00925         //                         GVertex::hadronInelastic);
00926 
00927         // Capture process
00928         G4HadronCaptureProcess* cap_process = new G4HadronCaptureProcess();
00929 
00930         // High precision model and data below 20 MeV
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         // LEP model above 20 MeV - default cross sections are used here
00937         // hence no need to explicitly invoke AddDataSet method
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         // Note: need to update GVertex to include hadronCapture
00944 //       physicsRegistrar->Register(cap_process,pmanager,"hadi",
00945 //                                 GVertex::hadronInelastic);
00946       
00947         // Fission process
00948         G4HadronFissionProcess* fis_process = new G4HadronFissionProcess();
00949 
00950         // High precision model and data below 20 MeV
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         // LEP model above 20 MeV - default cross sections are used here
00957         // hence no need to explicitly invoke AddDataSet method
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         // Note: need to update GVertex to include hadronFission
00964 //        physicsRegistrar->Register(fis_process,pmanager,"hadi",
00965 //                                 GVertex::hadronInelastic);
00966       
00967       } else {
00968 
00969         // Elastic process
00970         G4HadronElasticProcess *process = new G4HadronElasticProcess();
00971         G4LElastic *model =  new G4LElastic();
00972         process->RegisterMe(model);
00973         pmanager->AddDiscreteProcess(process);
00974 
00975         // Inelastic process
00976         G4NeutronInelasticProcess *inel_process = 
00977                                      new G4NeutronInelasticProcess();
00978         inel_process->RegisterMe(bertini_model);
00979         pmanager->AddDiscreteProcess(inel_process);
00980 //        physicsRegistrar->Register(inel_process,pmanager,"hadi",
00981 //                                 GVertex::hadronInelastic);
00982       }
00983     // *********
00984     // * n-bar *
00985     // *********
00986     } else if (particleName == "anti_neutron") {
00987 
00988       // Elastic process
00989       G4HadronElasticProcess *process = new G4HadronElasticProcess();
00990       G4LElastic *model =  new G4LElastic();
00991       process->RegisterMe(model);
00992       pmanager->AddDiscreteProcess(process);
00993 
00994       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01002 //                               GVertex::hadronInelastic);
01003     // **********
01004     // * lambda *
01005     // **********
01006     } else if (particleName == "lambda") {
01007 
01008       // Elastic process
01009       G4HadronElasticProcess *process = new G4HadronElasticProcess();
01010       G4LElastic *model =  new G4LElastic();
01011       process->RegisterMe(model);
01012       pmanager->AddDiscreteProcess(process);
01013 
01014       // Inelastic process
01015 
01016       G4LambdaInelasticProcess* inel_process = new
01017 G4LambdaInelasticProcess();
01018       inel_process->RegisterMe(bertini_model);
01019       pmanager->AddDiscreteProcess(inel_process);
01020 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01021 //                               GVertex::hadronInelastic);
01022     // **************
01023     // * lambda-bar *
01024     // **************
01025     } else if (particleName == "anti_lambda") {
01026 
01027       // Elastic process
01028       G4HadronElasticProcess *process = new G4HadronElasticProcess();
01029       G4LElastic *model =  new G4LElastic();
01030       process->RegisterMe(model);
01031       pmanager->AddDiscreteProcess(process);
01032 
01033       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01041 //                               GVertex::hadronInelastic);
01042     // **************
01043     // * deuteron   *
01044     // **************
01045     } else if (particleName == "deuteron") {
01046 
01047       // Elastic process
01048       G4HadronElasticProcess *process = new G4HadronElasticProcess();
01049       G4LElastic *model =  new G4LElastic();
01050       process->RegisterMe(model);
01051       pmanager->AddDiscreteProcess(process);
01052 
01053       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01061 //                               GVertex::hadronInelastic);
01062     // **************
01063     // * triton     *
01064     // **************
01065     } else if (particleName == "triton") {
01066 
01067       // Elastic process
01068       G4HadronElasticProcess *process = new G4HadronElasticProcess();
01069       G4LElastic *model =  new G4LElastic();
01070       process->RegisterMe(model);
01071       pmanager->AddDiscreteProcess(process);
01072 
01073       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01081 //                               GVertex::hadronInelastic);
01082     // **************
01083     // * alpha      *
01084     // **************
01085     } else if (particleName == "alpha") {
01086 
01087       // Elastic process
01088       G4HadronElasticProcess *process = new G4HadronElasticProcess();
01089       G4LElastic *model =  new G4LElastic();
01090       process->RegisterMe(model);
01091       pmanager->AddDiscreteProcess(process);
01092 
01093       // Inelastic process
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 //      physicsRegistrar->Register(inel_process,pmanager,"hadi",
01100 //                               GVertex::hadronInelastic);
01101     }
01102   }  // while 
01103 }

void BgsPhysicsList::ConstructIonFix (  )  [protected]

Definition at line 1142 of file BgsPhysicsList.cc.

Referenced by ConstructProcess().

01143 {
01144   BgsGenocide *genocide = new
01145 //BgsGentleGenocide(control->GetIonEnergyCut()*MeV,
01146 //                                              60);
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 //        physicsRegistrar->Register( genocide, pmanager, "ionFix", 
01163 //                                  GVertex::minimumEnergy );
01164     }
01165   }
01166 }

void BgsPhysicsList::ConstructIons (  )  [protected]

Definition at line 276 of file BgsPhysicsList.cc.

Referenced by ConstructParticle().

00277 {  
00278   G4GenericIon::GenericIonDefinition();
00279 }

void BgsPhysicsList::ConstructLeptHad (  )  [protected]

Definition at line 611 of file BgsPhysicsList.cc.

References rb::Gamma().

Referenced by ConstructProcess().

00612 {
00613   //
00614   // Gamma-nuclear process
00615   //
00616 
00617   // low energy part
00618   G4GammaNuclearReaction* lowEGammaModel = new G4GammaNuclearReaction();
00619   lowEGammaModel->SetMaxEnergy(3.5*GeV);
00620   G4PhotoNuclearProcess* thePhotoNuclearProcess = new
00621 G4PhotoNuclearProcess();
00622   thePhotoNuclearProcess->RegisterMe(lowEGammaModel);
00623 
00624   // bias the cross section
00625   //
00626 /*  double thePhotoNuclearBias = control->GetPhotoNuclearBias();
00627   if (thePhotoNuclearBias != 1.) {
00628     thePhotoNuclearProcess->BiasCrossSectionByFactor(thePhotoNuclearBias);
00629 
00630   // print out a warning if biasing
00631   //
00632     ErrMsg(warning) << "*** Biasing the photo-nuclear process by factor "
00633                     << thePhotoNuclearBias << endmsg; 
00634   }
00635 */
00636   // high energy part
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   // Electron-nuclear process
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   // bias the cross section
00673   //
00674 /*
00675   G4double theElectroNuclearBias = control->GetElectroNuclearBias();
00676   if (theElectroNuclearBias != 1.) {
00677 
00678 theElectronNuclearProcess->BiasCrossSectionByFactor(theElectroNuclearBias);
00679 
00680   // print out a warning if biasing
00681   //
00682     ErrMsg(warning) << "*** Biasing the electron-nuclear process by factor "
00683                     << theElectroNuclearBias << endmsg; 
00684   }
00685 */
00686   //
00687   // Positron-nuclear process
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   // bias the cross section
00698   //
00699 /*
00700   if (theElectroNuclearBias != 1.) { 
00701 
00702 thePositronNuclearProcess->BiasCrossSectionByFactor(theElectroNuclearBias);
00703     ErrMsg(warning) << "*** Biasing the positron-nuclear process by factor "
00704                     << theElectroNuclearBias << endmsg; 
00705   }
00706   */
00707 }

void BgsPhysicsList::ConstructLeptons (  )  [protected]

Definition at line 211 of file BgsPhysicsList.cc.

Referenced by ConstructParticle().

00212 {
00213   // leptons
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 }

void BgsPhysicsList::ConstructMesons (  )  [protected]

Definition at line 230 of file BgsPhysicsList.cc.

Referenced by ConstructParticle().

00231 {
00232   //  mesons
00233   G4PionPlus     ::PionPlusDefinition();
00234   G4PionMinus    ::PionMinusDefinition();
00235   G4PionZero     ::PionZeroDefinition();
00236   G4Eta          ::EtaDefinition();
00237   G4EtaPrime     ::EtaPrimeDefinition();
00238   //  G4RhoZero      ::RhoZeroDefinition();
00239   G4KaonPlus     ::KaonPlusDefinition();
00240   G4KaonMinus    ::KaonMinusDefinition();
00241   G4KaonZero     ::KaonZeroDefinition();
00242   G4AntiKaonZero ::AntiKaonZeroDefinition();
00243   G4KaonZeroLong ::KaonZeroLongDefinition();
00244   G4KaonZeroShort::KaonZeroShortDefinition();
00245 }

void BgsPhysicsList::ConstructNeutrinoGenocide (  )  [protected]

Definition at line 1109 of file BgsPhysicsList.cc.

Referenced by ConstructProcess().

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         // temporary fix for K0, K0bar until CHIPS photonuclear model isfixed
01127 
01128         particleName == "kaon0"  ||
01129         particleName == "anti_kaon0"    ) {
01130         
01131         pmanager->AddProcess(genocide, -1, -1, 5);
01132 //        physicsRegistrar->Register( genocide, pmanager, "neutrinoDeath", 
01133 //                                  GVertex::neutrino );
01134     }
01135   }
01136 }

void BgsPhysicsList::ConstructNeutronFix (  )  [protected]

Definition at line 1172 of file BgsPhysicsList.cc.

Referenced by ConstructProcess().

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 //        physicsRegistrar->Register( genocide, pmanager, "neutronFix", 
01187 //                                  GVertex::minimumEnergy );
01188     }
01189   }
01190 }

void BgsPhysicsList::ConstructParticle (  )  [protected]

Definition at line 149 of file BgsPhysicsList.cc.

References ConstructBaryons(), ConstructBosons(), ConstructIons(), ConstructLeptons(), and ConstructMesons().

00150 {
00151 
00152   // In this method, static member functions should be called
00153   // for all particles which you want to use.
00154   // This ensures that objects of these particle types will be
00155   // created in the program.
00156 
00157   ConstructBosons();
00158   ConstructLeptons();
00159   ConstructMesons();
00160   ConstructBaryons();
00161   ConstructIons();
00162 /*
00163   if (first) {
00164     first = false;
00165 
00166     //
00167     // Check to make sure our particles have valid PDG encodings
00168     //
00169     BgsPDGEncode encode(false);
00170 
00171     G4ParticleTable *table = G4ParticleTable::GetParticleTable();
00172 
00173     G4int nProb = 0;
00174     G4int n = table->entries();
00175     while(n--) {
00176       G4ParticleDefinition *part = table->GetParticle(n);
00177       if (encode.pdt(part) == 0) {
00178         nProb++;
00179         std::cerr << "The Geant4 particle \""
00180                       << part->GetParticleName()
00181                       << "\" is not recognized by BgsPDGEncode" << std::endl;
00182       }
00183     }
00184   
00185     if(nProb > 0) std::cerr << "One or more PDG encoding errors" <<
00186 std::endl;
00187   }
00188 */
00189   // Add short lived particles for high energy models, 
00190   // but don't check PDG codes - they are not propagated in Bogus anyway
00191    
00192   G4ShortLivedConstructor shortLived;
00193   shortLived.ConstructParticle();
00194 }

void BgsPhysicsList::ConstructProcess (  )  [protected]

Definition at line 282 of file BgsPhysicsList.cc.

References ConstructEM(), ConstructGeneral(), ConstructHad(), ConstructIonFix(), ConstructLeptHad(), ConstructNeutrinoGenocide(), and ConstructNeutronFix().

00283 {
00284   //if (control->UseBgsTran()) {
00285 /*  if(0) {
00286     AddBgsTransportation(control->GetMaxTrackingStepSize(),
00287                          control->GetMaxVacStepSize(), 
00288                          control->GetVacName());
00289     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^+ " << 
00290             std::endl;
00291     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+ " <<
00292             std::endl;
00293     std::cout << " +---                                                  ---+ " <<
00294             std::endl;
00295     std::cout << " +---                 BgsTransportation                ---+ " <<
00296       std::endl;
00297     std::cout << " +---                      USED !!                     ---+ " <<
00298       std::endl;
00299     std::cout << " +---                                                  ---+ " <<
00300       std::endl;
00301     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^+ " <<
00302       std::endl;
00303     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+ " <<
00304       std::endl;
00305     
00306   }
00307   else {
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   //  AddParameterisation();
00328   ConstructEM();
00329   ConstructLeptHad();
00330   ConstructHad();
00331   ConstructGeneral();
00332   ConstructNeutrinoGenocide();
00333   ConstructIonFix();
00334   ConstructNeutronFix();
00335 }

void BgsPhysicsList::SetCuts (  ) 

Definition at line 1196 of file BgsPhysicsList.cc.

01197 {
01198   // Set default cuts, all volumes 
01199 
01200   SetDefaultCutValue(0.7*mm);
01201   SetCutsWithDefault();
01202 
01203   // Enable print out of cuts after tables are built
01204   // This is now done in BgsRunAction
01205   //
01206   //    if (verboseLevel > 1) DumpCutValuesTable(); 
01207 }

void BgsPhysicsList::SetStatusEmProcess (  ) 

Definition at line 517 of file BgsPhysicsList.cc.

00518 {
00519 }


Member Data Documentation

G4CascadeInterface* BgsPhysicsList::bertini_model [private]

Definition at line 126 of file BgsPhysicsList.hh.

Referenced by BgsPhysicsList(), and ConstructHad().

G4bool BgsPhysicsList::first [private]

Definition at line 117 of file BgsPhysicsList.hh.


Generated on Tue Nov 29 23:17:54 2016 for BOSS_7.0.2 by  doxygen 1.4.7