// $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 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
+//_____________________________________________________________________________
TG4PhysicsManager&
TG4PhysicsManager::operator=(const TG4PhysicsManager& right)
{
// private methods
+//_____________________________________________________________________________
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("MuBrem", 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);
}
+//_____________________________________________________________________________
+void TG4PhysicsManager::GstparCut(G4int itmed, TG4G3Cut par, G4double parval)
+{
+// Sets special tracking medium parameter.
+// It is applied to all logical volumes that use the specified
+// tracking medium.
+// ---
-// public methods
+ // get medium from table
+ G3MedTableEntry* medium = G3Med.get(itmed);
+ if (!medium) {
+ G4String text = "TG4PhysicsManager::GstparCut: \n";
+ text = text + " Medium not found.";
+ G4Exception(text);
+ }
+
+ // get/create user limits
+ 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
+ limits->SetG3Cut(par, parval);
+}
-void TG4PhysicsManager::BuildPhysics()
+
+//_____________________________________________________________________________
+void TG4PhysicsManager::GstparControl(G4int itmed, TG4G3Control par,
+ TG4G3ControlValue parval)
{
-// Empty function - not needed in G4.
-// (Physics is built within /run/initialize.)
+// Sets special tracking medium parameter.
+// It is applied to all logical volumes that use the specified
+// tracking medium.
// ---
- TG4Globals::Warning(
- "TG4PhysicsManager::BuildPhysics: is empty function in G4 MC.");
-}
+ // get medium from table
+ G3MedTableEntry* medium = G3Med.get(itmed);
+ if (!medium) {
+ G4String text = "TG4PhysicsManager::GstparControl: \n";
+ text = text + " Medium not found.";
+ G4Exception(text);
+ }
+
+ // get/create user limits
+ 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
+ limits->SetG3Control(par, parval);
+}
+// public methods
+
+//_____________________________________________________________________________
+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
+// its value is different from the "global" physics setup.
+// (If: CheckCut/ControlWithG3Defaults is used checking
+// is performed with respect to G3 default values.)
+// When any cut/control parameter is set in limits
+// the physics manager is locked and the physics setup
+// cannot be changed.
+// Applying this TG4Limits to particles and physical
+// processes is still in development.
+//
+// Geant3 desription:
+// ==================
+// To change the value of cut or mechanism "CHPAR"
+// to a new value PARVAL for tracking medium ITMED
+// The data structure JTMED contains the standard tracking
+// parameters (CUTS and flags to control the physics processes) which
+// are used by default for all tracking media. It is possible to
+// redefine individually with GSTPAR any of these parameters for a
+// given tracking medium.
+// ITMED tracking medium number
+// CHPAR is a character string (variable name)
+// PARVAL must be given as a floating point.
+// ---
+
+ G4String name = TG4GeometryServices::Instance()->CutName(param);
+ TG4G3Cut cut;
+ if (fG3PhysicsManager->CheckCutWithTheVector(name, parval, cut)) {
+ GstparCut(itmed, cut, parval);
+ fG3PhysicsManager->Lock();
+ }
+ else {
+ TG4G3Control control;
+ TG4G3ControlValue controlValue;
+ if (fG3PhysicsManager
+ ->CheckControlWithTheVector(name, parval, control, controlValue)) {
+ GstparControl(itmed, control, controlValue);
+ fG3PhysicsManager->Lock();
+ }
+ else if (cut==kNoG3Cuts && control==kNoG3Controls) {
+ G4String text = "TG4PhysicsManager::Gstpar:";
+ text = text + name;
+ text = text + " parameter is not yet implemented.";
+ TG4Globals::Warning(text);
+ }
+ }
+}
+
+//_____________________________________________________________________________
void TG4PhysicsManager::CreatePhysicsConstructors()
{
// Creates the selected physics constructors
// 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
}
+//_____________________________________________________________________________
void TG4PhysicsManager::SetCut(const char* cutName, Float_t cutValue)
{
// Sets the specified cut.
// ---
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;
}
}
+//_____________________________________________________________________________
Float_t TG4PhysicsManager::Xsec(char* ch, Float_t p1, Int_t i1, Int_t i2)
{
// Not yet implemented -> gives exception.
return 0.;
}
+//_____________________________________________________________________________
Int_t TG4PhysicsManager::IdFromPDG(Int_t pdgID) const
{
// G4 does not use the integer particle identifiers
return pdgID;
}
+//_____________________________________________________________________________
Int_t TG4PhysicsManager::PDGFromId(Int_t mcID) const
{
// G4 does not use integer particle identifiers
return mcID;
}
+//_____________________________________________________________________________
void TG4PhysicsManager::DefineParticles()
{
// ======
fParticlesManager->MapParticles();
}
-
+//_____________________________________________________________________________
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();
}
+//_____________________________________________________________________________
AliMCProcess TG4PhysicsManager::GetMCProcess(const G4VProcess* process)
{
// Returns the AliMCProcess code of the specified G4 process.
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;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff