Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Updated list of magnetic field steppers and EM physlists for Geant4 10.4 #22073

Merged
merged 5 commits into from Feb 5, 2018
Merged

Updated list of magnetic field steppers and EM physlists for Geant4 10.4 #22073

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
7f337a1
use 10.4 steppers
civanch Jan 31, 2018
8475a90
update EM physics classes
civanch Jan 31, 2018
09180f7
correct configuration of protons
civanch Feb 1, 2018
961b06f
correct configuration of protons
civanch Feb 1, 2018
bf10d06
change LOgInfo
civanch Feb 2, 2018
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
3 changes: 1 addition & 2 deletions SimG4Core/Application/src/ParametrisedEMPhysics.cc
View file
Open in desktop
Expand Up @@ -221,9 +221,8 @@ void ParametrisedEMPhysics::ConstructProcess() {
}
// enable fluorescence
bool fluo = theParSet.getParameter<bool>("FlagFluo");
if(fluo) {
if(fluo && !G4LossTableManager::Instance()->AtomDeexcitation()) {
G4VAtomDeexcitation* de = new G4UAtomicDeexcitation();
G4LossTableManager::Instance()->SetAtomDeexcitation(de);
de->SetFluo(true);
}
}
14 changes: 9 additions & 5 deletions SimG4Core/MagneticField/src/FieldStepper.cc
View file
Open in desktop
Expand Up @@ -6,8 +6,9 @@
#include "G4SimpleRunge.hh"
#include "G4SimpleHeum.hh"
#include "G4CashKarpRKF45.hh"
#include "G4RKG3_Stepper.hh"
#include "G4ExplicitEuler.hh"
#include "G4DormandPrince745.hh"
#include "G4BogackiShampine45.hh"
#include "G4TsitourasRK45.hh"
#include "G4ImplicitEuler.hh"
#include "G4HelixExplicitEuler.hh"
#include "G4HelixImplicitEuler.hh"
Expand Down Expand Up @@ -38,8 +39,9 @@ G4MagIntegratorStepper * FieldStepper::select(const std::string & ss)
else if (ss == "G4SimpleRunge") theStepper = new G4SimpleRunge(theEquation);
else if (ss == "G4SimpleHeum") theStepper = new G4SimpleHeum(theEquation);
else if (ss == "G4CashKarpRKF45") theStepper = new G4CashKarpRKF45(theEquation);
else if (ss == "G4RKG3_Stepper") theStepper = new G4RKG3_Stepper(theEquation);
else if (ss == "G4ExplicitEuler") theStepper = new G4ExplicitEuler(theEquation);
else if (ss == "G4DormandPrince745") theStepper = new G4DormandPrince745(theEquation);
else if (ss == "G4BogackiShampine45") theStepper = new G4BogackiShampine45(theEquation);
else if (ss == "G4TsitourasRK45") theStepper = new G4TsitourasRK45(theEquation);
else if (ss == "G4ImplicitEuler") theStepper = new G4ImplicitEuler(theEquation);
else if (ss == "G4HelixExplicitEuler") theStepper = new G4HelixExplicitEuler(theEquation);
else if (ss == "G4HelixImplicitEuler") theStepper = new G4HelixImplicitEuler(theEquation);
Expand All @@ -48,8 +50,10 @@ G4MagIntegratorStepper * FieldStepper::select(const std::string & ss)
else
{
edm::LogWarning("SimG4CoreMagneticField")
<< " FieldStepper invalid choice, defaulting to G4ClassicalRK4 ";
<< " FieldStepper <" << ss << "> is not known, defaulting to G4ClassicalRK4 ";
theStepper = new G4ClassicalRK4(theEquation);
}
edm::LogInfo("SimG4CoreMagneticField")
<< "### FieldStepper: <" << ss << ">";
return theStepper;
}
43 changes: 25 additions & 18 deletions SimG4Core/PhysicsLists/src/CMSEmStandardPhysics.cc
View file
Open in desktop
Expand Up @@ -78,6 +78,7 @@ CMSEmStandardPhysics::CMSEmStandardPhysics(G4int ver) :
param->SetDefaults();
param->SetVerbose(verbose);
param->SetApplyCuts(true);
param->SetStepFunction(0.8, 1*CLHEP::mm);
param->SetMscRangeFactor(0.2);
param->SetMscStepLimitType(fMinimal);
SetPhysicsType(bElectromagnetic);
Expand Down Expand Up @@ -154,12 +155,18 @@ void CMSEmStandardPhysics::ConstructProcess() {
G4hPairProduction* pp = nullptr;

// muon & hadron multiple scattering
G4MuMultipleScattering* mumsc= nullptr;
G4hMultipleScattering* pimsc= nullptr;
G4MuMultipleScattering* mumsc = nullptr;
G4hMultipleScattering* pimsc = nullptr;
G4hMultipleScattering* kmsc = nullptr;
G4hMultipleScattering* pmsc = nullptr;
G4hMultipleScattering* hmsc = nullptr;

// muon and hadron single scattering
G4CoulombScattering* muss = nullptr;
G4CoulombScattering* piss = nullptr;
G4CoulombScattering* kss = nullptr;
G4CoulombScattering* pss = nullptr;

// high energy limit for e+- scattering models and bremsstrahlung
G4double highEnergyLimit = 100*MeV;

Expand All @@ -177,16 +184,14 @@ void CMSEmStandardPhysics::ConstructProcess() {
} else if (particleName == "e-") {

G4eIonisation* eioni = new G4eIonisation();
eioni->SetStepFunction(0.8, 1.0*mm);

G4eMultipleScattering* msc = new G4eMultipleScattering;
msc->SetStepLimitType(fMinimal);
UrbanMscModel93* msc1 = new UrbanMscModel93();
G4WentzelVIModel* msc2 = new G4WentzelVIModel();
msc1->SetHighEnergyLimit(highEnergyLimit);
msc2->SetLowEnergyLimit(highEnergyLimit);
msc->AddEmModel(0, msc1);
msc->AddEmModel(0, msc2);
msc->SetEmModel(msc1);
msc->SetEmModel(msc2);

G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
G4CoulombScattering* ss = new G4CoulombScattering();
Expand All @@ -203,16 +208,14 @@ void CMSEmStandardPhysics::ConstructProcess() {
} else if (particleName == "e+") {

G4eIonisation* eioni = new G4eIonisation();
eioni->SetStepFunction(0.8, 1.0*mm);

G4eMultipleScattering* msc = new G4eMultipleScattering;
msc->SetStepLimitType(fMinimal);
UrbanMscModel93* msc1 = new UrbanMscModel93();
G4WentzelVIModel* msc2 = new G4WentzelVIModel();
msc1->SetHighEnergyLimit(highEnergyLimit);
msc2->SetLowEnergyLimit(highEnergyLimit);
msc->AddEmModel(0, msc1);
msc->AddEmModel(0, msc2);
msc->SetEmModel(msc1);
msc->SetEmModel(msc2);

G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
G4CoulombScattering* ss = new G4CoulombScattering();
Expand All @@ -234,13 +237,14 @@ void CMSEmStandardPhysics::ConstructProcess() {
mub = new G4MuBremsstrahlung();
mup = new G4MuPairProduction();
mumsc = new G4MuMultipleScattering();
mumsc->AddEmModel(0, new G4WentzelVIModel());
mumsc->SetEmModel(new G4WentzelVIModel());
muss = new G4CoulombScattering();
}
ph->RegisterProcess(mumsc, particle);
ph->RegisterProcess(new G4MuIonisation(), particle);
ph->RegisterProcess(mub, particle);
ph->RegisterProcess(mup, particle);
ph->RegisterProcess(new G4CoulombScattering(), particle);
ph->RegisterProcess(muss, particle);

} else if (particleName == "alpha" ||
particleName == "He3" ) {
Expand All @@ -263,13 +267,14 @@ void CMSEmStandardPhysics::ConstructProcess() {
pib = new G4hBremsstrahlung();
pip = new G4hPairProduction();
pimsc = new G4hMultipleScattering();
pimsc->AddEmModel(0, new G4WentzelVIModel());
pimsc->SetEmModel(new G4WentzelVIModel());
piss = new G4CoulombScattering();
}
ph->RegisterProcess(pimsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(pib, particle);
ph->RegisterProcess(pip, particle);
ph->RegisterProcess(new G4CoulombScattering(), particle);
ph->RegisterProcess(piss, particle);

} else if (particleName == "kaon+" ||
particleName == "kaon-" ) {
Expand All @@ -278,13 +283,14 @@ void CMSEmStandardPhysics::ConstructProcess() {
kb = new G4hBremsstrahlung();
kp = new G4hPairProduction();
kmsc = new G4hMultipleScattering();
kmsc->AddEmModel(0, new G4WentzelVIModel());
kmsc->SetEmModel(new G4WentzelVIModel());
kss = new G4CoulombScattering();
}
ph->RegisterProcess(kmsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(kb, particle);
ph->RegisterProcess(kp, particle);
ph->RegisterProcess(new G4CoulombScattering(), particle);
ph->RegisterProcess(kss, particle);

} else if (particleName == "proton" ||
particleName == "anti_proton") {
Expand All @@ -293,14 +299,15 @@ void CMSEmStandardPhysics::ConstructProcess() {
pb = new G4hBremsstrahlung();
pp = new G4hPairProduction();
pmsc = new G4hMultipleScattering();
pmsc->AddEmModel(0, new G4WentzelVIModel());
}
pmsc->SetEmModel(new G4WentzelVIModel());
pss = new G4CoulombScattering();

ph->RegisterProcess(pmsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(pb, particle);
ph->RegisterProcess(pp, particle);
ph->RegisterProcess(new G4CoulombScattering(), particle);
ph->RegisterProcess(pss, particle);

} else if (particleName == "B+" ||
particleName == "B-" ||
Expand Down
36 changes: 20 additions & 16 deletions SimG4Core/PhysicsLists/src/CMSEmStandardPhysics95.cc
View file
Open in desktop
Expand Up @@ -77,6 +77,7 @@ CMSEmStandardPhysics95::CMSEmStandardPhysics95(const G4String& name, G4int ver,
param->SetDefaults();
param->SetVerbose(verbose);
param->SetApplyCuts(true);
param->SetStepFunction(0.8, 1*CLHEP::mm);
param->SetMscRangeFactor(0.2);
param->SetMscStepLimitType(fMinimal);
SetPhysicsType(bElectromagnetic);
Expand Down Expand Up @@ -132,11 +133,12 @@ void CMSEmStandardPhysics95::ConstructParticle() {

void CMSEmStandardPhysics95::ConstructProcess()
{
// This EM builder takes default models of Geant4 10 EMV.
// This EM builder takes default models of Geant4 EMV.
// Multiple scattering by Urban for all particles

G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();

// muon & hadron bremsstrahlung and pair production
G4MuBremsstrahlung* mub = nullptr;
G4MuPairProduction* mup = nullptr;
G4hBremsstrahlung* pib = nullptr;
Expand All @@ -146,7 +148,12 @@ void CMSEmStandardPhysics95::ConstructProcess()
G4hBremsstrahlung* pb = nullptr;
G4hPairProduction* pp = nullptr;

G4hMultipleScattering* hmsc = nullptr;
// muon & hadron multiple scattering
G4MuMultipleScattering* mumsc = nullptr;
G4hMultipleScattering* pimsc = nullptr;
G4hMultipleScattering* kmsc = nullptr;
G4hMultipleScattering* pmsc = nullptr;
G4hMultipleScattering* hmsc = nullptr;

// Add standard EM Processes
G4ParticleTable* table = G4ParticleTable::GetParticleTable();
Expand All @@ -166,10 +173,7 @@ void CMSEmStandardPhysics95::ConstructProcess()
} else if (particleName == "e-") {

G4eIonisation* eioni = new G4eIonisation();
eioni->SetStepFunction(0.8, 1.0*mm);
G4eMultipleScattering* msc = new G4eMultipleScattering();
msc->SetStepLimitType(fMinimal);

G4eBremsstrahlung* ebrem = new G4eBremsstrahlung();

ph->RegisterProcess(msc, particle);
Expand All @@ -179,10 +183,7 @@ void CMSEmStandardPhysics95::ConstructProcess()
} else if (particleName == "e+") {

G4eIonisation* eioni = new G4eIonisation();
eioni->SetStepFunction(0.8, 1.0*mm);
G4eMultipleScattering* msc = new G4eMultipleScattering();
msc->SetStepLimitType(fMinimal);

G4eBremsstrahlung* ebrem = new G4eBremsstrahlung();

ph->RegisterProcess(msc, particle);
Expand All @@ -191,14 +192,14 @@ void CMSEmStandardPhysics95::ConstructProcess()
ph->RegisterProcess(new G4eplusAnnihilation(), particle);

} else if (particleName == "mu+" ||
particleName == "mu-" ) {
particleName == "mu-" ) {

if(nullptr == mub) {
mub = new G4MuBremsstrahlung();
mup = new G4MuPairProduction();
mumsc = new G4MuMultipleScattering();
}
G4MuMultipleScattering* msc = new G4MuMultipleScattering();
ph->RegisterProcess(msc, particle);
ph->RegisterProcess(mumsc, particle);
ph->RegisterProcess(new G4MuIonisation(), particle);
ph->RegisterProcess(mub, particle);
ph->RegisterProcess(mup, particle);
Expand All @@ -223,8 +224,9 @@ void CMSEmStandardPhysics95::ConstructProcess()
if(nullptr == pib) {
pib = new G4hBremsstrahlung();
pip = new G4hPairProduction();
pimsc = new G4hMultipleScattering();
}
ph->RegisterProcess(new G4hMultipleScattering(), particle);
ph->RegisterProcess(pimsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(pib, particle);
ph->RegisterProcess(pip, particle);
Expand All @@ -234,21 +236,23 @@ void CMSEmStandardPhysics95::ConstructProcess()

if(nullptr == kb) {
kb = new G4hBremsstrahlung();
kp = new G4hPairProduction();
kp = new G4hPairProduction();
kmsc = new G4hMultipleScattering();
}
ph->RegisterProcess(new G4hMultipleScattering(), particle);
ph->RegisterProcess(kmsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(kb, particle);
ph->RegisterProcess(kp, particle);

} else if (particleName == "proton" ||
} else if (particleName == "proton" ||
particleName == "anti_proton") {

if(nullptr == pb) {
pb = new G4hBremsstrahlung();
pp = new G4hPairProduction();
}
ph->RegisterProcess(new G4hMultipleScattering(), particle);
pmsc = new G4hMultipleScattering();
ph->RegisterProcess(pmsc, particle);
ph->RegisterProcess(new G4hIonisation(), particle);
ph->RegisterProcess(pb, particle);
ph->RegisterProcess(pp, particle);
Expand Down