// $Id$
// Category: physics
//
+// Author: I. Hrivnacova
+//
+// Class TG4PhysicsManager
+// -----------------------
// See the class description in the header file.
#include "TG4PhysicsManager.h"
+#include "TG4ModularPhysicsList.h"
#include "TG4ParticlesManager.h"
#include "TG4G3PhysicsManager.h"
+#include "TG4PhysicsConstructorGeneral.h"
#include "TG4PhysicsConstructorEM.h"
-#include "TG4PhysicsConstructorOptical.h"
+#include "TG4PhysicsConstructorMuon.h"
#include "TG4PhysicsConstructorHadron.h"
+#include "TG4PhysicsConstructorIon.h"
+#include "TG4PhysicsConstructorOptical.h"
#include "TG4PhysicsConstructorSpecialCuts.h"
#include "TG4PhysicsConstructorSpecialControls.h"
#include "TG4GeometryServices.h"
#include "TG4G3Cut.h"
#include "TG4G3Control.h"
+#include "TG4G3Units.h"
#include "TG4Limits.h"
#include "AliDecayer.h"
#include <G4ParticleDefinition.hh>
+#include <G4OpBoundaryProcess.hh>
#include <G4VProcess.hh>
-#include <G4VModularPhysicsList.hh>
#include <G3MedTable.hh>
#include <TDatabasePDG.h>
TG4PhysicsManager* TG4PhysicsManager::fgInstance = 0;
//_____________________________________________________________________________
-TG4PhysicsManager::TG4PhysicsManager(G4VModularPhysicsList* physicsList)
- : fPhysicsList(physicsList),
+TG4PhysicsManager::TG4PhysicsManager(TG4ModularPhysicsList* physicsList)
+ : TG4Verbose("physicsManager"),
+ fMessenger(this),
+ fPhysicsList(physicsList),
fDecayer(0),
fSetEMPhysics(true),
- fSetOpticalPhysics(false),
+ fSetMuonPhysics(true),
fSetHadronPhysics(false),
+ fSetOpticalPhysics(false),
fSetSpecialCutsPhysics(false),
fSetSpecialControlsPhysics(false)
fG3PhysicsManager = new TG4G3PhysicsManager();
// fill process name map
- FillProcessMap();
+ // FillProcessMap();
}
//_____________________________________________________________________________
-TG4PhysicsManager::TG4PhysicsManager(){
+TG4PhysicsManager::TG4PhysicsManager()
+ : TG4Verbose("physicsManager"),
+ fMessenger(this) {
//
- delete fDecayer;
- delete fParticlesManager;
- delete fG3PhysicsManager;
}
//_____________________________________________________________________________
-TG4PhysicsManager::TG4PhysicsManager(const TG4PhysicsManager& right) {
+TG4PhysicsManager::TG4PhysicsManager(const TG4PhysicsManager& right)
+ : TG4Verbose("physicsManager"),
+ fMessenger(this) {
//
TG4Globals::Exception(
"Attempt to copy TG4PhysicsManager singleton.");
//_____________________________________________________________________________
TG4PhysicsManager::~TG4PhysicsManager() {
//
+ delete fDecayer;
+ delete fParticlesManager;
+ delete fG3PhysicsManager;
}
// operators
//_____________________________________________________________________________
void TG4PhysicsManager::FillProcessMap()
{
-// Fills fProcessMap.
+// Fills fProcessMCMap.
// The default G4 process names are used in the map.
+// Not used - the map is filled in physics constructors
+// only with processes that are really built.
// ---
// multiple scattering
- fProcessMap.Add("msc", kPMultipleScattering);
- fProcessMap.Add("Imsc", kPMultipleScattering);
+ fProcessMCMap.Add("msc", kPMultipleScattering);
+ fProcessMCMap.Add("Imsc", kPMultipleScattering);
// continuous energy loss
// !! including delta rays
- fProcessMap.Add("eIoni", kPEnergyLoss);
- fProcessMap.Add("IeIoni", kPEnergyLoss);
- fProcessMap.Add("LowEnergyIoni", kPEnergyLoss);
- fProcessMap.Add("hIoni", kPEnergyLoss);
- fProcessMap.Add("IhIoni", kPEnergyLoss);
- fProcessMap.Add("hLowEIoni", kPEnergyLoss);
- fProcessMap.Add("MuIoni", kPEnergyLoss);
- fProcessMap.Add("IMuIonisation", kPEnergyLoss);
- fProcessMap.Add("ionIoni", kPEnergyLoss);
- fProcessMap.Add("ionLowEIoni", kPEnergyLoss);
- fProcessMap.Add("PAIonisation", kPEnergyLoss);
+ fProcessMCMap.Add("eIoni", kPEnergyLoss);
+ fProcessMCMap.Add("IeIoni", kPEnergyLoss);
+ fProcessMCMap.Add("LowEnergyIoni", kPEnergyLoss);
+ fProcessMCMap.Add("hIoni", kPEnergyLoss);
+ fProcessMCMap.Add("IhIoni", kPEnergyLoss);
+ fProcessMCMap.Add("hLowEIoni", kPEnergyLoss);
+ fProcessMCMap.Add("MuIoni", kPEnergyLoss);
+ fProcessMCMap.Add("IMuIonisation", kPEnergyLoss);
+ fProcessMCMap.Add("ionIoni", kPEnergyLoss);
+ fProcessMCMap.Add("ionLowEIoni", kPEnergyLoss);
+ fProcessMCMap.Add("PAIonisation", kPEnergyLoss);
// bending in mag. field
// kPMagneticFieldL
// particle decay
- fProcessMap.Add("Decay", kPDecay);
+ fProcessMCMap.Add("Decay", kPDecay);
// photon pair production or
// muon direct pair production
- fProcessMap.Add("conv", kPPair);
- fProcessMap.Add("LowEnConversion", kPPair);
- fProcessMap.Add("MuPairProd", kPPair);
- fProcessMap.Add("IMuPairProduction", kPPair);
+ fProcessMCMap.Add("conv", kPPair);
+ fProcessMCMap.Add("LowEnConversion", kPPair);
+ fProcessMCMap.Add("MuPairProd", kPPair);
+ fProcessMCMap.Add("IMuPairProduction", kPPair);
// Compton scattering
- fProcessMap.Add("compt", kPCompton);
- fProcessMap.Add("LowEnCompton", kPCompton);
- fProcessMap.Add("polarCompt", kPCompton);
+ fProcessMCMap.Add("compt", kPCompton);
+ fProcessMCMap.Add("LowEnCompton", kPCompton);
+ fProcessMCMap.Add("polarCompt", kPCompton);
// photoelectric effect
- fProcessMap.Add("phot", kPPhotoelectric);
- fProcessMap.Add("LowEnPhotoElec", kPPhotoelectric);
+ fProcessMCMap.Add("phot", kPPhotoelectric);
+ fProcessMCMap.Add("LowEnPhotoElec", kPPhotoelectric);
// bremsstrahlung
- fProcessMap.Add("eBrem", kPBrem);
- fProcessMap.Add("IeBrem", kPBrem);
- fProcessMap.Add("MuBrems", kPBrem);
- fProcessMap.Add("IMuBremsstrahlung", kPBrem);
- fProcessMap.Add("LowEnBrem", kPBrem);
+ fProcessMCMap.Add("eBrem", kPBrem);
+ fProcessMCMap.Add("IeBrem", kPBrem);
+ fProcessMCMap.Add("MuBrems", kPBrem);
+ fProcessMCMap.Add("IMuBremsstrahlung", kPBrem);
+ fProcessMCMap.Add("LowEnBrem", kPBrem);
// delta-ray production
// kPDeltaRay
// has to be distinguished from kPEnergyLoss on flight
// positron annihilation
- fProcessMap.Add("annihil", kPAnnihilation);
- fProcessMap.Add("Iannihil", kPAnnihilation);
+ fProcessMCMap.Add("annihil", kPAnnihilation);
+ fProcessMCMap.Add("Iannihil", kPAnnihilation);
// hadronic interaction
// kPHadronic
// kPNuclearFission
// nuclear absorption
- fProcessMap.Add("PionMinusAbsorptionAtRest", kPNuclearAbsorption);
- fProcessMap.Add("PiMinusAbsorptionAtRest", kPNuclearAbsorption);
- fProcessMap.Add("KaonMinusAbsorption", kPNuclearAbsorption);
- fProcessMap.Add("KaonMinusAbsorptionAtRest", kPNuclearAbsorption);
+ fProcessMCMap.Add("PionMinusAbsorptionAtRest", kPNuclearAbsorption);
+ fProcessMCMap.Add("PiMinusAbsorptionAtRest", kPNuclearAbsorption);
+ fProcessMCMap.Add("KaonMinusAbsorption", kPNuclearAbsorption);
+ fProcessMCMap.Add("KaonMinusAbsorptionAtRest", kPNuclearAbsorption);
// antiproton annihilation
- fProcessMap.Add("AntiProtonAnnihilationAtRest", kPPbarAnnihilation);
- // fProcessMap.Add("AntiNeutronAnnihilationAtRest", not defined);
+ fProcessMCMap.Add("AntiProtonAnnihilationAtRest", kPPbarAnnihilation);
+ // fProcessMCMap.Add("AntiNeutronAnnihilationAtRest", not defined);
// neutron capture
- fProcessMap.Add("NeutronCaptureAtRest", kPNCapture);
- // fProcessMap.Add("LCapture", hadron capture not defined);
+ fProcessMCMap.Add("NeutronCaptureAtRest", kPNCapture);
+ // fProcessMCMap.Add("LCapture", hadron capture not defined);
// hadronic elastic incoherent scattering
- fProcessMap.Add("LElastic", kPHElastic);
+ fProcessMCMap.Add("LElastic", kPHElastic);
// hadronic inelastic scattering
- fProcessMap.Add("inelastic", kPHInhelastic);
+ fProcessMCMap.Add("inelastic", kPHInhelastic);
// muon nuclear interaction
- fProcessMap.Add("MuNucl", kPMuonNuclear);
+ fProcessMCMap.Add("MuNucl", kPMuonNuclear);
// exceeded time of flight cut
// kPTOFlimit
// kPPhotoFission
// Rayleigh scattering
- fProcessMap.Add("Rayleigh Scattering", kPRayleigh);
+ fProcessMCMap.Add("Rayleigh Scattering", kPRayleigh);
// no mechanism is active, usually at the entrance of a new volume
- fProcessMap.Add("Transportation", kPNull);
+ fProcessMCMap.Add("Transportation", kPNull);
// particle has fallen below energy threshold and tracking stops
// kPStop
// Cerenkov photon absorption
- fProcessMap.Add("Absorption", kPLightAbsorption);
+ fProcessMCMap.Add("Absorption", kPLightAbsorption);
// Cerenkov photon reflection/refraction
// kPLightScattering, kPLightReflection, kPLightRefraction
// has to be inquired from the G4OpBoundary process
// synchrotron radiation
- fProcessMap.Add("SynchrotronRadiation", kPSynchrotron);
+ fProcessMCMap.Add("SynchrotronRadiation", kPSynchrotron);
}
//_____________________________________________________________________________
}
// get/create user limits
- G4UserLimits* limits = medium->GetLimits();
- TG4Limits* tg4Limits;
- if (limits) {
- tg4Limits = dynamic_cast<TG4Limits*> (limits);
- if (!tg4Limits)
- G4Exception("TG4PhysicsManager::GstparCut: Wrong limits type.");
- }
- else {
- tg4Limits = new TG4Limits();
- medium->SetLimits(tg4Limits);
-
- // add verbose
- G4cout << "TG4PhysicsManager::GstparCut: new TG4Limits() for medium "
- << itmed << " has been created." << G4endl;
+ TG4Limits* limits
+ = TG4GeometryServices::Instance()->GetLimits(medium->GetLimits());
+
+ if (!limits) {
+ limits = new TG4Limits(*fG3PhysicsManager->GetCutVector(),
+ *fG3PhysicsManager->GetControlVector());
+ medium->SetLimits(limits);
+
+ if (VerboseLevel() > 1) {
+ G4cout << "TG4PhysicsManager::GstparCut: new TG4Limits() for medium "
+ << itmed << " has been created." << G4endl;
+ }
}
+
+ // add units
+ if (par == kTOFMAX) parval *= TG4G3Units::Time();
+ else parval *= TG4G3Units::Energy();
+
// set parameter
- tg4Limits->SetG3Cut(par, parval*GeV);
+ limits->SetG3Cut(par, parval);
}
//_____________________________________________________________________________
void TG4PhysicsManager::GstparControl(G4int itmed, TG4G3Control par,
- G4double parval)
+ TG4G3ControlValue parval)
{
// Sets special tracking medium parameter.
// It is applied to all logical volumes that use the specified
}
// get/create user limits
- G4UserLimits* limits = medium->GetLimits();
- TG4Limits* tg4Limits;
- if (limits) {
- tg4Limits = dynamic_cast<TG4Limits*> (limits);
- if (!tg4Limits)
- G4Exception("TG4PhysicsManager::GstparControl: Wrong limits type.");
- }
- else {
- tg4Limits = new TG4Limits();
- medium->SetLimits(tg4Limits);
-
- // add verbose
- G4cout << "TG4PhysicsManager::GstparControl: new TG4Limits() for medium"
- << itmed << " has been created." << G4endl;
- }
+ TG4Limits* limits
+ = TG4GeometryServices::Instance()->GetLimits(medium->GetLimits());
+
+ if (!limits) {
+ limits = new TG4Limits(*fG3PhysicsManager->GetCutVector(),
+ *fG3PhysicsManager->GetControlVector());
+ medium->SetLimits(limits);
+
+ if (VerboseLevel() > 1) {
+ G4cout << "TG4PhysicsManager::GstparControl: new TG4Limits() for medium"
+ << itmed << " has been created." << G4endl;
+ }
+ }
+
// set parameter
- tg4Limits->SetG3Control(par, parval);
+ limits->SetG3Control(par, parval);
}
// public methods
//_____________________________________________________________________________
-void TG4PhysicsManager::BuildPhysics()
-{
-// Empty function - not needed in G4.
-// (Physics is built within /run/initialize.)
-// ---
-
- TG4Globals::Warning(
- "TG4PhysicsManager::BuildPhysics: is empty function in G4 MC.");
-}
-
-//_____________________________________________________________________________
-void TG4PhysicsManager::Gstpar(Int_t itmed, const char *param,
- Float_t parval)
+void TG4PhysicsManager::Gstpar(Int_t itmed, const char *param, Float_t parval)
{
// Passes the tracking medium parameter to TG4Limits.
// The tracking medium parameter is set only in case
}
else {
TG4G3Control control;
- if (fG3PhysicsManager->CheckControlWithTheVector(name, parval, control)) {
- GstparControl(itmed, control, parval);
+ TG4G3ControlValue controlValue;
+ if (fG3PhysicsManager
+ ->CheckControlWithTheVector(name, parval, control, controlValue)) {
+ GstparControl(itmed, control, controlValue);
fG3PhysicsManager->Lock();
}
else if (cut==kNoG3Cuts && control==kNoG3Controls) {
// and registeres them in the modular physics list.
// ---
+ // general physics
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorGeneral(VerboseLevel()));
+
// electromagnetic physics
if (fSetEMPhysics)
- fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorEM());
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorEM(VerboseLevel()));
- // optical physics
- if (fSetOpticalPhysics)
- fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorOptical());
+ // muon physics
+ if (fSetMuonPhysics && fSetEMPhysics)
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorMuon(VerboseLevel()));
// hadron physics
- if (fSetHadronPhysics)
- fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorHadron());
+ if (fSetEMPhysics || fSetHadronPhysics) {
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorIon(
+ VerboseLevel(), fSetEMPhysics, fSetHadronPhysics));
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorHadron(
+ VerboseLevel(), fSetEMPhysics, fSetHadronPhysics));
+ }
+ // optical physics
+ if (fSetOpticalPhysics)
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorOptical(VerboseLevel()));
+
+ // special processes
if (fSetSpecialCutsPhysics)
- fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorSpecialCuts());
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorSpecialCuts(VerboseLevel()));
if (fSetSpecialControlsPhysics)
- fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorSpecialControls());
+ fPhysicsList->RegisterPhysics(
+ new TG4PhysicsConstructorSpecialControls(VerboseLevel()));
+
+ // warn about not allowed combinations
+ if (fSetMuonPhysics && !fSetEMPhysics) {
+ G4String text = "TG4PhysicsManager::CreatePhysicsConstructors:\n";
+ text = text + " Muon physics cannot be constructed without EM physics.\n";
+ text = text + " SetMuon control was ignored.";
+ TG4Globals::Warning(text);
+ }
// all created physics constructors are deleted
- // in the G4VModularPhysicsList destructor
+ // in the TG4ModularPhysicsList destructor
}
//_____________________________________________________________________________
// ---
fG3PhysicsManager->CheckLock();
- TG4G3Cut g3Cut = fG3PhysicsManager->GetG3Cut(cutName);
- if (g3Cut != kNoG3Cuts)
- fG3PhysicsManager->SetCut(g3Cut, cutValue);
- else {
+ TG4G3Cut g3Cut = TG4G3CutVector::GetCut(cutName);
+
+ if (g3Cut == kNoG3Cuts) {
G4String text = "TG4PhysicsManager::SetCut:\n";
text = text + " Parameter " + cutName;
text = text + " is not implemented.";
TG4Globals::Warning(text);
+ return;
}
+
+ // add units
+ if (g3Cut == kTOFMAX) cutValue *= TG4G3Units::Time();
+ else cutValue *= TG4G3Units::Energy();
+
+ fG3PhysicsManager->SetCut(g3Cut, cutValue);
}
//_____________________________________________________________________________
-void TG4PhysicsManager::SetProcess(const char* controlName, Int_t controlValue)
+void TG4PhysicsManager::SetProcess(const char* controlName, Int_t value)
{
// Sets the specified process control.
// ---
fG3PhysicsManager->CheckLock();
- TG4G3Control control = fG3PhysicsManager->GetG3Control(controlName);
- if (control != kNoG3Controls)
+ TG4G3Control control = TG4G3ControlVector::GetControl(controlName);
+
+ if (control != kNoG3Controls) {
+ TG4G3ControlValue controlValue
+ = TG4G3ControlVector::GetControlValue(value, control);
fG3PhysicsManager->SetProcess(control, controlValue);
+ }
else {
G4String text = "TG4PhysicsManager::SetProcess:\n";
text = text + " Parameter " + controlName;
void TG4PhysicsManager::SetProcessActivation()
{
// (In)Activates built processes according
-// to the setup in fControlVector.
+// to the setup in TG4G3PhysicsManager::fControlVector.
// ---
- if (fPhysicsList) {
- // temporarily excluded
- // fPhysicsList->SetProcessActivation();
- }
- else {
- G4String text = "TG4PhysicsManager::SetProcessActivation:\n";
- text = text + " There is no physics list set.";
- TG4Globals::Exception(text);
- }
+ fPhysicsList->SetProcessActivation();
}
if (!process) return kPNoProcess;
- G4String name = process->GetProcessName();
- G4int code = fProcessMap.GetSecond(name);
+ return fProcessMCMap.GetMCProcess(process);
+}
+
+//_____________________________________________________________________________
+AliMCProcess TG4PhysicsManager::GetOpBoundaryStatus(const G4VProcess* process)
+{
+// Returns the AliMCProcess code according to the OpBoundary process
+// status.
+// ---
+
+ if (!process) return kPNoProcess;
+
+#ifdef TGEANT4_DEBUG
+ G4OpBoundaryProcess* opBoundary
+ = dynamic_cast<G4OpBoundaryProcess*>(process);
+
+ if (!opBoundary)
+ TG4Globals::Exception(
+ "TG4PhysicsManager::GetOpBoundaryStatus: Wrong process type.");
+ return kPNoProcess;
+ }
- if (code == 0) return kPNoProcess;
+ return opBoundary;
+#else
+ G4OpBoundaryProcess* opBoundary = (G4OpBoundaryProcess*)process;
+#endif
+
+ switch (opBoundary->GetStatus()) {
+ // reflection
+ case FresnelReflection:
+ case TotalInternalReflection:
+ case LambertianReflection:
+ case LobeReflection:
+ case SpikeReflection:
+ case BackScattering:
+ return kPLightReflection;
+ ;;
+
+ // refraction
+ case FresnelRefraction:
+ return kPLightRefraction;
+ ;;
+
+ // absorption
+ case Absorption:
+ case Detection:
+ return kPLightAbsorption;
+ ;;
+ }
- return (AliMCProcess)code;
+ // should not happen
+ return kPNoProcess;
}