Geant4-11
Public Member Functions | Static Public Member Functions | Private Types | Private Member Functions | Private Attributes | Static Private Attributes | Friends
G4LossTableManager Class Reference

#include <G4LossTableManager.hh>

Public Member Functions

G4VAtomDeexcitationAtomDeexcitation ()
 
void BuildPhysicsTable (const G4ParticleDefinition *aParticle)
 
void BuildPhysicsTable (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
void DeRegister (G4VEmFluctuationModel *p)
 
void DeRegister (G4VEmModel *p)
 
void DeRegister (G4VEmProcess *p)
 
void DeRegister (G4VEnergyLossProcess *p)
 
void DeRegister (G4VMultipleScattering *p)
 
void DeRegister (G4VProcess *p)
 
void DumpHtml ()
 
G4ElectronIonPairElectronIonPair ()
 
G4EmConfiguratorEmConfigurator ()
 
G4EmCorrectionsEmCorrections ()
 
G4EmSaturationEmSaturation ()
 
 G4LossTableManager (G4LossTableManager &)=delete
 
G4double GetCSDARange (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetDEDX (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetDEDXDispersion (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double &length)
 
G4VEmProcessGetElectronGeneralProcess ()
 
const std::vector< G4VEmProcess * > & GetEmProcessVector ()
 
G4double GetEnergy (const G4ParticleDefinition *aParticle, G4double range, const G4MaterialCutsCouple *couple)
 
G4VEnergyLossProcessGetEnergyLossProcess (const G4ParticleDefinition *)
 
const std::vector< G4VEnergyLossProcess * > & GetEnergyLossProcessVector ()
 
G4VEmProcessGetGammaGeneralProcess ()
 
const std::vector< G4VMultipleScattering * > & GetMultipleScatteringVector ()
 
G4VEmProcessGetPositronGeneralProcess ()
 
G4double GetRange (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetRangeFromRestricteDEDX (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4LossTableBuilderGetTableBuilder ()
 
G4bool IsMaster () const
 
void LocalPhysicsTables (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
G4NIELCalculatorNIELCalculator ()
 
G4LossTableManageroperator= (const G4LossTableManager &right)=delete
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VEmProcess *p, G4bool theMaster)
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VMultipleScattering *p, G4bool theMaster)
 
void Register (G4VEmFluctuationModel *p)
 
void Register (G4VEmModel *p)
 
void Register (G4VEmProcess *p)
 
void Register (G4VEnergyLossProcess *p)
 
void Register (G4VMultipleScattering *p)
 
void Register (G4VProcess *p)
 
void RegisterExtraParticle (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
void ResetParameters ()
 
void SetAtomDeexcitation (G4VAtomDeexcitation *)
 
void SetElectronGeneralProcess (G4VEmProcess *)
 
void SetGammaGeneralProcess (G4VEmProcess *)
 
void SetNIELCalculator (G4NIELCalculator *)
 
void SetPositronGeneralProcess (G4VEmProcess *)
 
void SetSubCutProducer (G4VSubCutProducer *)
 
void SetVerbose (G4int val)
 
G4VSubCutProducerSubCutProducer ()
 
 ~G4LossTableManager ()
 

Static Public Member Functions

static G4LossTableManagerInstance ()
 

Private Types

typedef const G4ParticleDefinitionPD
 

Private Member Functions

G4VEnergyLossProcessBuildTables (const G4ParticleDefinition *aParticle)
 
void Clear ()
 
void CopyDEDXTables ()
 
void CopyTables (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *)
 
 G4LossTableManager ()
 
void ParticleHaveNoLoss (const G4ParticleDefinition *aParticle)
 
void PrintEWarning (G4String, G4double)
 

Private Attributes

G4bool all_tables_are_built
 
G4VAtomDeexcitationatomDeexcitation
 
std::vector< PDbase_part_vector
 
G4VEnergyLossProcesscurrentLoss
 
PD currentParticle
 
std::vector< G4PhysicsTable * > dedx_vector
 
G4VEmProcesseGeneral
 
G4EmConfiguratoremConfigurator
 
G4EmCorrectionsemCorrections
 
G4ElectronIonPairemElectronIonPair
 
std::vector< G4VEmProcess * > emp_vector
 
PD firstParticle
 
std::vector< G4VEmFluctuationModel * > fmod_vector
 
G4VEmProcessgGeneral
 
std::vector< G4PhysicsTable * > inv_range_vector
 
std::vector< G4boolisActive
 
G4bool isMaster
 
std::map< PD, G4VEnergyLossProcess *, std::less< PD > > loss_map
 
std::vector< G4VEnergyLossProcess * > loss_vector
 
std::vector< G4VEmModel * > mod_vector
 
std::vector< G4VMultipleScattering * > msc_vector
 
G4int n_loss
 
G4NIELCalculatornielCalculator
 
std::vector< G4VProcess * > p_vector
 
std::vector< PDpart_vector
 
G4VEmProcesspGeneral
 
std::vector< G4PhysicsTable * > range_vector
 
G4int run
 
G4bool startInitialisation
 
G4VSubCutProducersubcutProducer
 
G4LossTableBuildertableBuilder
 
std::vector< G4booltables_are_built
 
PD theElectron
 
PD theGenericIon
 
G4EmParameterstheParameters
 
G4int verbose
 

Static Private Attributes

static G4ThreadLocal G4LossTableManagerinstance = nullptr
 

Friends

class G4ThreadLocalSingleton< G4LossTableManager >
 

Detailed Description

Definition at line 77 of file G4LossTableManager.hh.

Member Typedef Documentation

◆ PD

Definition at line 256 of file G4LossTableManager.hh.

Constructor & Destructor Documentation

◆ ~G4LossTableManager()

G4LossTableManager::~G4LossTableManager ( )

Definition at line 98 of file G4LossTableManager.cc.

99{
100 for (G4int i=0; i<n_loss; ++i) {
101 delete loss_vector[i];
102 }
103 size_t msc = msc_vector.size();
104 for (size_t j=0; j<msc; ++j) {
105 delete msc_vector[j];
106 }
107 size_t emp = emp_vector.size();
108 for (size_t k=0; k<emp; ++k) {
109 delete emp_vector[k];
110 }
111 emp = p_vector.size();
112 for (size_t k=0; k<emp; ++k) {
113 delete p_vector[k];
114 }
115 size_t mod = mod_vector.size();
116 size_t fmod = fmod_vector.size();
117 //G4cout << " Nmod" << mod << " Nfluc= " << fmod << G4endl;
118 for (size_t a=0; a<mod; ++a) {
119 //G4cout << "Delete model #" << a << " " << mod_vector[a] << G4endl;
120 if( nullptr != mod_vector[a] ) {
121 for (size_t b=0; b<fmod; ++b) {
122 if((G4VEmModel*)(fmod_vector[b]) == mod_vector[a]) {
123 fmod_vector[b] = nullptr;
124 }
125 }
126 delete mod_vector[a];
127 mod_vector[a] = nullptr;
128 }
129 }
130 for (size_t b=0; b<fmod; ++b) {
131 delete fmod_vector[b];
132 }
133 Clear();
134 delete tableBuilder;
135 delete emCorrections;
136 delete emConfigurator;
137 delete emElectronIonPair;
138 delete nielCalculator;
139 delete atomDeexcitation;
140 delete subcutProducer;
141}
int G4int
Definition: G4Types.hh:85
G4ElectronIonPair * emElectronIonPair
std::vector< G4VEmProcess * > emp_vector
G4VSubCutProducer * subcutProducer
std::vector< G4VMultipleScattering * > msc_vector
G4EmConfigurator * emConfigurator
std::vector< G4VEmModel * > mod_vector
std::vector< G4VEnergyLossProcess * > loss_vector
std::vector< G4VProcess * > p_vector
G4LossTableBuilder * tableBuilder
G4EmCorrections * emCorrections
std::vector< G4VEmFluctuationModel * > fmod_vector
G4NIELCalculator * nielCalculator
G4VAtomDeexcitation * atomDeexcitation

References atomDeexcitation, Clear(), emConfigurator, emCorrections, emElectronIonPair, emp_vector, fmod_vector, loss_vector, mod_vector, msc_vector, n_loss, nielCalculator, p_vector, subcutProducer, and tableBuilder.

◆ G4LossTableManager() [1/2]

G4LossTableManager::G4LossTableManager ( G4LossTableManager )
delete

◆ G4LossTableManager() [2/2]

G4LossTableManager::G4LossTableManager ( )
private

Definition at line 145 of file G4LossTableManager.cc.

146{
148 n_loss = 0;
149 run = -1;
150 startInitialisation = false;
151 all_tables_are_built = false;
152 currentLoss = nullptr;
153 currentParticle = nullptr;
154 firstParticle = nullptr;
155 isMaster = true;
158 theGenericIon= nullptr;
161 isMaster = false;
162 }
165 emConfigurator = nullptr;
166 emElectronIonPair = nullptr;
167 atomDeexcitation = nullptr;
168 subcutProducer = nullptr;
169 nielCalculator = nullptr;
170 gGeneral = nullptr;
171 eGeneral = pGeneral = nullptr;
172}
static G4Electron * Electron()
Definition: G4Electron.cc:93
static G4EmParameters * Instance()
G4int Verbose() const
G4int WorkerVerbose() const
G4VEnergyLossProcess * currentLoss
G4EmParameters * theParameters
G4bool IsWorkerThread()
Definition: G4Threading.cc:123
Definition: run.py:1

References all_tables_are_built, atomDeexcitation, currentLoss, currentParticle, eGeneral, G4Electron::Electron(), emConfigurator, emCorrections, emElectronIonPair, firstParticle, gGeneral, G4EmParameters::Instance(), isMaster, G4Threading::IsWorkerThread(), n_loss, nielCalculator, pGeneral, startInitialisation, subcutProducer, tableBuilder, theElectron, theGenericIon, theParameters, G4EmParameters::Verbose(), verbose, and G4EmParameters::WorkerVerbose().

Member Function Documentation

◆ AtomDeexcitation()

G4VAtomDeexcitation * G4LossTableManager::AtomDeexcitation ( )
inline

◆ BuildPhysicsTable() [1/2]

void G4LossTableManager::BuildPhysicsTable ( const G4ParticleDefinition aParticle)

◆ BuildPhysicsTable() [2/2]

void G4LossTableManager::BuildPhysicsTable ( const G4ParticleDefinition aParticle,
G4VEnergyLossProcess p 
)

Definition at line 627 of file G4LossTableManager.cc.

630{
631 if(1 < verbose) {
632 G4cout << "### G4LossTableManager::BuildPhysicsTable() for "
633 << aParticle->GetParticleName()
634 << " and process " << p->GetProcessName() << G4endl;
635 }
636 // clear configurator
637 if(-1 == run && startInitialisation) {
639 firstParticle = aParticle;
640 }
642 ++run;
643 if(1 < verbose) {
644 G4cout << "===== G4LossTableManager::BuildPhysicsTable() for run "
645 << run << " ===== " << atomDeexcitation << G4endl;
646 }
647 currentParticle = nullptr;
649 }
650
651 // initialisation before any table is built
652 if ( startInitialisation && aParticle == firstParticle ) {
653
654 startInitialisation = false;
655 if(1 < verbose) {
656 G4cout << "### G4LossTableManager start initialisation for first particle "
658 << G4endl;
659 }
660
662
663 for (G4int i=0; i<n_loss; ++i) {
665
666 if(el) {
667 isActive[i] = true;
668 base_part_vector[i] = el->BaseParticle();
669 tables_are_built[i] = false;
671 if(!isActive[i]) {
672 el->SetIonisation(false);
673 tables_are_built[i] = true;
674 }
675
676 if(1 < verbose) {
677 G4cout << i <<". "<< el->GetProcessName();
678 if(el->Particle()) {
679 G4cout << " for " << el->Particle()->GetParticleName();
680 }
681 G4cout << " active= " << isActive[i]
682 << " table= " << tables_are_built[i]
683 << " isIonisation= " << el->IsIonisationProcess();
684 if(base_part_vector[i]) {
685 G4cout << " base particle "
686 << base_part_vector[i]->GetParticleName();
687 }
688 G4cout << G4endl;
689 }
690 } else {
691 tables_are_built[i] = true;
692 part_vector[i] = nullptr;
693 isActive[i] = false;
694 }
695 }
696 }
697
700 return;
701 }
702
703 // Build tables for given particle
705
706 for(G4int i=0; i<n_loss; ++i) {
707 if(p == loss_vector[i] && !tables_are_built[i] && !base_part_vector[i]) {
708 const G4ParticleDefinition* curr_part = part_vector[i];
709 if(1 < verbose) {
710 G4cout << "### Build Table for " << p->GetProcessName()
711 << " and " << curr_part->GetParticleName()
712 << " " << tables_are_built[i] << " " << base_part_vector[i]
713 << G4endl;
714 }
715 G4VEnergyLossProcess* curr_proc = BuildTables(curr_part);
716 if(curr_proc) {
717 CopyTables(curr_part, curr_proc);
718 if(p == curr_proc && 0 == run && p->IsIonisationProcess()) {
719 loss_map[aParticle] = p;
720 //G4cout << "G4LossTableManager::BuildPhysicsTable: "
721 // << aParticle->GetParticleName()
722 // << " added to map " << p << G4endl;
723 }
724 }
725 }
726 if ( !tables_are_built[i] ) { all_tables_are_built = false; }
727 }
728 if(1 < verbose) {
729 G4cout << "### G4LossTableManager::BuildPhysicsTable end: "
730 << "all_tables_are_built= " << all_tables_are_built << " "
731 << aParticle->GetParticleName() << " proc: " << p << G4endl;
732 }
734 if(1 < verbose) {
735 G4cout << "%%%%% All dEdx and Range tables are built for master run= "
736 << run << " %%%%%" << G4endl;
737 }
738 }
739}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
void SetIsPrintedFlag(G4bool val)
std::vector< PD > base_part_vector
std::vector< G4bool > isActive
void CopyTables(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *)
std::map< PD, G4VEnergyLossProcess *, std::less< PD > > loss_map
std::vector< G4bool > tables_are_built
std::vector< PD > part_vector
G4VEnergyLossProcess * BuildTables(const G4ParticleDefinition *aParticle)
const G4String & GetParticleName() const
const G4ParticleDefinition * BaseParticle() const
const G4ParticleDefinition * Particle() const
G4bool IsIonisationProcess() const
void SetIonisation(G4bool val)
const G4String & GetProcessName() const
Definition: G4VProcess.hh:382

References all_tables_are_built, atomDeexcitation, base_part_vector, G4VEnergyLossProcess::BaseParticle(), BuildTables(), G4EmConfigurator::Clear(), CopyTables(), currentParticle, emConfigurator, firstParticle, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), G4NIELCalculator::Initialise(), isActive, G4VEnergyLossProcess::IsIonisationProcess(), loss_map, loss_vector, n_loss, nielCalculator, part_vector, G4VEnergyLossProcess::Particle(), run, G4VEnergyLossProcess::SetIonisation(), G4EmParameters::SetIsPrintedFlag(), startInitialisation, tables_are_built, theParameters, and verbose.

◆ BuildTables()

G4VEnergyLossProcess * G4LossTableManager::BuildTables ( const G4ParticleDefinition aParticle)
private

Definition at line 784 of file G4LossTableManager.cc.

786{
787 if(1 < verbose) {
788 G4cout << "G4LossTableManager::BuildTables() for "
789 << aParticle->GetParticleName() << G4endl;
790 }
791
792 std::vector<G4PhysicsTable*> t_list;
793 std::vector<G4VEnergyLossProcess*> loss_list;
794 std::vector<G4bool> build_flags;
795 G4VEnergyLossProcess* em = nullptr;
796 G4VEnergyLossProcess* p = nullptr;
797 G4int iem = 0;
798 G4PhysicsTable* dedx = nullptr;
799 G4int i;
800
801 G4ProcessVector* pvec =
802 aParticle->GetProcessManager()->GetProcessList();
803 G4int nvec = pvec->size();
804
805 for (i=0; i<n_loss; ++i) {
806 p = loss_vector[i];
807 if (p) {
808 G4bool yes = (aParticle == part_vector[i]);
809
810 // possible case of process sharing between particle/anti-particle
811 if(!yes) {
812 G4VProcess* ptr = static_cast<G4VProcess*>(p);
813 for(G4int j=0; j<nvec; ++j) {
814 //G4cout << "j= " << j << " " << (*pvec)[j] << " " << ptr << G4endl;
815 if(ptr == (*pvec)[j]) {
816 yes = true;
817 break;
818 }
819 }
820 }
821 // process belong to this particle
822 if(yes && isActive[i]) {
823 if (p->IsIonisationProcess() || !em) {
824 em = p;
825 iem= i;
826 }
827 // tables may be shared between particle/anti-particle
828 G4bool val = false;
829 if (!tables_are_built[i]) {
830 val = true;
831 dedx = p->BuildDEDXTable(fRestricted);
832 //G4cout << "Build DEDX table for " << p->GetProcessName()
833 // << " idx= " << i << dedx << " " << dedx->length() << G4endl;
834 p->SetDEDXTable(dedx,fRestricted);
835 tables_are_built[i] = true;
836 } else {
837 dedx = p->DEDXTable();
838 }
839 t_list.push_back(dedx);
840 loss_list.push_back(p);
841 build_flags.push_back(val);
842 }
843 }
844 }
845
846 G4int n_dedx = t_list.size();
847 if (0 == n_dedx || !em) {
848 G4cout << "G4LossTableManager WARNING: no DEDX processes for "
849 << aParticle->GetParticleName() << G4endl;
850 return 0;
851 }
852 G4int nSubRegions = em->NumberOfSubCutoffRegions();
853
854 if (1 < verbose) {
855 G4cout << "G4LossTableManager::BuildTables() start to build range tables"
856 << " and the sum of " << n_dedx << " processes"
857 << " iem= " << iem << " em= " << em->GetProcessName()
858 << " buildCSDARange= " << theParameters->BuildCSDARange()
859 << " nSubRegions= " << nSubRegions;
860 if(subcutProducer) {
861 G4cout << " SubCutProducer " << subcutProducer->GetName();
862 }
863 G4cout << G4endl;
864 }
865 // do not build tables if producer class is defined
866 if(subcutProducer) { nSubRegions = 0; }
867
868 dedx = em->DEDXTable();
869 em->SetIonisation(true);
870 em->SetDEDXTable(dedx, fIsIonisation);
871
872 if (1 < n_dedx) {
873 dedx = 0;
875 tableBuilder->BuildDEDXTable(dedx, t_list);
876 em->SetDEDXTable(dedx, fRestricted);
877 }
878
879 /*
880 if(2==run && "e-" == aParticle->GetParticleName()) {
881 G4cout << "G4LossTableManager::BuildTables for e- " << dedx << G4endl;
882 G4cout << (*dedx) << G4endl;
883 G4cout << "%%%%% Instance ID= " << (*dedx)[0]->GetInstanceID() << G4endl;
884 G4cout << "%%%%% LastValue= " << (*dedx)[0]->GetLastValue() << G4endl;
885 G4cout << "%%%%% 1.2 " << (*(dedx))[0]->Value(1.2) << G4endl;
886 }
887 */
888 dedx_vector[iem] = dedx;
889
890 G4PhysicsTable* range = em->RangeTableForLoss();
891 if(!range) range = G4PhysicsTableHelper::PreparePhysicsTable(range);
892 range_vector[iem] = range;
893
894 G4PhysicsTable* invrange = em->InverseRangeTable();
895 if(!invrange) invrange = G4PhysicsTableHelper::PreparePhysicsTable(invrange);
896 inv_range_vector[iem] = invrange;
897
898 tableBuilder->BuildRangeTable(dedx, range);
899 tableBuilder->BuildInverseRangeTable(range, invrange);
900
901 // if(1<verbose) G4cout << *dedx << G4endl;
902
903 em->SetRangeTableForLoss(range);
904 em->SetInverseRangeTable(invrange);
905
906 // if(1<verbose) G4cout << *range << G4endl;
907
908 std::vector<G4PhysicsTable*> listSub;
909 std::vector<G4PhysicsTable*> listCSDA;
910
911 for (i=0; i<n_dedx; ++i) {
912 p = loss_list[i];
913 if(p != em) { p->SetIonisation(false); }
914 if(build_flags[i]) {
916 }
918 dedx = p->BuildDEDXTable(fTotal);
919 p->SetDEDXTable(dedx,fTotal);
920 listCSDA.push_back(dedx);
921 }
922 }
923
925 G4PhysicsTable* dedxCSDA = em->DEDXunRestrictedTable();
926 if (1 < n_dedx) {
927 dedxCSDA = nullptr;
929 tableBuilder->BuildDEDXTable(dedxCSDA, listCSDA);
930 em->SetDEDXTable(dedxCSDA,fTotal);
931 }
932 G4PhysicsTable* rCSDA = em->CSDARangeTable();
933 if(!rCSDA) { rCSDA = G4PhysicsTableHelper::PreparePhysicsTable(rCSDA); }
934 tableBuilder->BuildRangeTable(dedxCSDA, rCSDA);
935 em->SetCSDARangeTable(rCSDA);
936 }
937
938 if (1 < verbose) {
939 G4cout << "G4LossTableManager::BuildTables: Tables are built for "
940 << aParticle->GetParticleName()
941 << "; ionisation process: " << em->GetProcessName()
942 << " " << em
943 << G4endl;
944 }
945 return em;
946}
@ fTotal
@ fRestricted
@ fIsIonisation
bool G4bool
Definition: G4Types.hh:86
G4bool BuildCSDARange() const
void BuildRangeTable(const G4PhysicsTable *dedxTable, G4PhysicsTable *rangeTable)
void BuildDEDXTable(G4PhysicsTable *dedxTable, const std::vector< G4PhysicsTable * > &)
void BuildInverseRangeTable(const G4PhysicsTable *rangeTable, G4PhysicsTable *invRangeTable)
std::vector< G4PhysicsTable * > range_vector
std::vector< G4PhysicsTable * > dedx_vector
std::vector< G4PhysicsTable * > inv_range_vector
G4ProcessManager * GetProcessManager() const
static G4PhysicsTable * PreparePhysicsTable(G4PhysicsTable *physTable)
void push_back(G4PhysicsVector *)
G4ProcessVector * GetProcessList() const
std::size_t size() const
G4PhysicsTable * RangeTableForLoss() const
G4PhysicsTable * InverseRangeTable() const
G4PhysicsTable * CSDARangeTable() const
void SetRangeTableForLoss(G4PhysicsTable *p)
G4int NumberOfSubCutoffRegions() const
G4PhysicsTable * BuildDEDXTable(G4EmTableType tType=fRestricted)
void SetInverseRangeTable(G4PhysicsTable *p)
void SetDEDXTable(G4PhysicsTable *p, G4EmTableType tType)
G4PhysicsTable * BuildLambdaTable(G4EmTableType tType=fRestricted)
void SetLambdaTable(G4PhysicsTable *p)
void SetCSDARangeTable(G4PhysicsTable *pRange)
G4PhysicsTable * DEDXunRestrictedTable() const
G4PhysicsTable * DEDXTable() const
const G4String & GetName() const

References G4EmParameters::BuildCSDARange(), G4VEnergyLossProcess::BuildDEDXTable(), G4LossTableBuilder::BuildDEDXTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4VEnergyLossProcess::BuildLambdaTable(), G4LossTableBuilder::BuildRangeTable(), G4VEnergyLossProcess::CSDARangeTable(), dedx_vector, G4VEnergyLossProcess::DEDXTable(), G4VEnergyLossProcess::DEDXunRestrictedTable(), fIsIonisation, fRestricted, fTotal, G4cout, G4endl, G4VSubCutProducer::GetName(), G4ParticleDefinition::GetParticleName(), G4ProcessManager::GetProcessList(), G4ParticleDefinition::GetProcessManager(), G4VProcess::GetProcessName(), inv_range_vector, G4VEnergyLossProcess::InverseRangeTable(), isActive, G4VEnergyLossProcess::IsIonisationProcess(), loss_vector, n_loss, G4VEnergyLossProcess::NumberOfSubCutoffRegions(), part_vector, G4PhysicsTableHelper::PreparePhysicsTable(), G4PhysicsTable::push_back(), range_vector, G4VEnergyLossProcess::RangeTableForLoss(), G4VEnergyLossProcess::SetCSDARangeTable(), G4VEnergyLossProcess::SetDEDXTable(), G4VEnergyLossProcess::SetInverseRangeTable(), G4VEnergyLossProcess::SetIonisation(), G4VEnergyLossProcess::SetLambdaTable(), G4VEnergyLossProcess::SetRangeTableForLoss(), G4ProcessVector::size(), subcutProducer, tableBuilder, tables_are_built, theParameters, and verbose.

Referenced by BuildPhysicsTable().

◆ Clear()

void G4LossTableManager::Clear ( )
private

Definition at line 176 of file G4LossTableManager.cc.

177{
178 all_tables_are_built = false;
179 currentLoss = nullptr;
180 currentParticle = nullptr;
181 if(n_loss) {
182 dedx_vector.clear();
183 range_vector.clear();
184 inv_range_vector.clear();
185 loss_map.clear();
186 loss_vector.clear();
187 part_vector.clear();
188 base_part_vector.clear();
189 tables_are_built.clear();
190 isActive.clear();
191 n_loss = 0;
192 }
193}

References all_tables_are_built, base_part_vector, currentLoss, currentParticle, dedx_vector, inv_range_vector, isActive, loss_map, loss_vector, n_loss, part_vector, range_vector, and tables_are_built.

Referenced by ~G4LossTableManager().

◆ CopyDEDXTables()

void G4LossTableManager::CopyDEDXTables ( )
private

◆ CopyTables()

void G4LossTableManager::CopyTables ( const G4ParticleDefinition aParticle,
G4VEnergyLossProcess base_proc 
)
private

Definition at line 743 of file G4LossTableManager.cc.

745{
746 for (G4int j=0; j<n_loss; ++j) {
747
749
750 if (!tables_are_built[j] && part == base_part_vector[j]) {
751 tables_are_built[j] = true;
752 proc->SetDEDXTable(base_proc->IonisationTable(),fRestricted);
753 proc->SetDEDXTable(base_proc->DEDXunRestrictedTable(),fTotal);
754 proc->SetCSDARangeTable(base_proc->CSDARangeTable());
755 proc->SetRangeTableForLoss(base_proc->RangeTableForLoss());
756 proc->SetInverseRangeTable(base_proc->InverseRangeTable());
757 proc->SetLambdaTable(base_proc->LambdaTable());
758 proc->SetIonisation(base_proc->IsIonisationProcess());
759 if(proc->IsIonisationProcess()) {
760 range_vector[j] = base_proc->RangeTableForLoss();
761 inv_range_vector[j] = base_proc->InverseRangeTable();
762 loss_map[part_vector[j]] = proc;
763 //G4cout << "G4LossTableManager::CopyTable "
764 // << part_vector[j]->GetParticleName()
765 // << " added to map " << proc << G4endl;
766 }
767 if (1 < verbose) {
768 G4cout << "For " << proc->GetProcessName()
769 << " for " << part_vector[j]->GetParticleName()
770 << " base_part= " << part->GetParticleName()
771 << " tables are assigned"
772 << G4endl;
773 }
774 }
775
776 if (theElectron == part && theElectron == proc->SecondaryParticle() ) {
777 proc->SetSecondaryRangeTable(base_proc->RangeTableForLoss());
778 }
779 }
780}
void SetSecondaryRangeTable(G4PhysicsTable *p)
G4PhysicsTable * IonisationTable() const
G4PhysicsTable * LambdaTable() const
const G4ParticleDefinition * SecondaryParticle() const

References base_part_vector, G4VEnergyLossProcess::CSDARangeTable(), G4VEnergyLossProcess::DEDXunRestrictedTable(), fRestricted, fTotal, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), inv_range_vector, G4VEnergyLossProcess::InverseRangeTable(), G4VEnergyLossProcess::IonisationTable(), G4VEnergyLossProcess::IsIonisationProcess(), G4VEnergyLossProcess::LambdaTable(), loss_map, loss_vector, n_loss, part_vector, range_vector, G4VEnergyLossProcess::RangeTableForLoss(), G4VEnergyLossProcess::SecondaryParticle(), G4VEnergyLossProcess::SetCSDARangeTable(), G4VEnergyLossProcess::SetDEDXTable(), G4VEnergyLossProcess::SetInverseRangeTable(), G4VEnergyLossProcess::SetIonisation(), G4VEnergyLossProcess::SetLambdaTable(), G4VEnergyLossProcess::SetRangeTableForLoss(), G4VEnergyLossProcess::SetSecondaryRangeTable(), tables_are_built, theElectron, and verbose.

Referenced by BuildPhysicsTable().

◆ DeRegister() [1/6]

void G4LossTableManager::DeRegister ( G4VEmFluctuationModel p)

Definition at line 380 of file G4LossTableManager.cc.

381{
382 size_t n = fmod_vector.size();
383 for (size_t i=0; i<n; ++i) {
384 if(fmod_vector[i] == p) { fmod_vector[i] = nullptr; }
385 }
386}

References fmod_vector, and CLHEP::detail::n.

◆ DeRegister() [2/6]

void G4LossTableManager::DeRegister ( G4VEmModel p)

Definition at line 355 of file G4LossTableManager.cc.

356{
357 //G4cout << "G4LossTableManager::DeRegister G4VEmModel : " << p << G4endl;
358 size_t n = mod_vector.size();
359 for (size_t i=0; i<n; ++i) {
360 if(mod_vector[i] == p) {
361 mod_vector[i] = nullptr;
362 break;
363 }
364 }
365}

References mod_vector, and CLHEP::detail::n.

◆ DeRegister() [3/6]

void G4LossTableManager::DeRegister ( G4VEmProcess p)

Definition at line 300 of file G4LossTableManager.cc.

301{
302 if(!p) { return; }
303 size_t emp = emp_vector.size();
304 for (size_t i=0; i<emp; ++i) {
305 if(emp_vector[i] == p) {
306 emp_vector[i] = nullptr;
307 break;
308 }
309 }
310}

References emp_vector.

◆ DeRegister() [4/6]

void G4LossTableManager::DeRegister ( G4VEnergyLossProcess p)

◆ DeRegister() [5/6]

void G4LossTableManager::DeRegister ( G4VMultipleScattering p)

Definition at line 270 of file G4LossTableManager.cc.

271{
272 if(!p) { return; }
273 size_t msc = msc_vector.size();
274 for (size_t i=0; i<msc; ++i) {
275 if(msc_vector[i] == p) {
276 msc_vector[i] = nullptr;
277 break;
278 }
279 }
280}

References msc_vector.

◆ DeRegister() [6/6]

void G4LossTableManager::DeRegister ( G4VProcess p)

Definition at line 330 of file G4LossTableManager.cc.

331{
332 if(!p) { return; }
333 size_t emp = p_vector.size();
334 for (size_t i=0; i<emp; ++i) {
335 if(p_vector[i] == p) {
336 p_vector[i] = nullptr;
337 break;
338 }
339 }
340}

References p_vector.

◆ DumpHtml()

void G4LossTableManager::DumpHtml ( )

Definition at line 1076 of file G4LossTableManager.cc.

1077{
1078 // Automatic generation of html documentation page for physics lists
1079 // List processes and models for the most important
1080 // particles in descending order of importance
1081 // NB. for model names with length > 18 characters the .rst file needs
1082 // to be edited by hand. Or modify G4EmModelManager::DumpModelList
1083
1084 char* dirName = std::getenv("G4PhysListDocDir");
1085 char* physList = std::getenv("G4PhysListName");
1086 if (dirName && physList) {
1087 G4String physListName = G4String(physList);
1088 G4String pathName = G4String(dirName) + "/" + physListName + ".rst";
1089
1090 std::ofstream outFile;
1091 outFile.open(pathName);
1092
1093 outFile << physListName << G4endl;
1094 outFile << std::string(physListName.length(), '=') << G4endl;
1095
1096 std::vector<G4ParticleDefinition*> particles {
1103 };
1104
1105 std::vector<G4VEmProcess*> emproc_vector = GetEmProcessVector();
1106 std::vector<G4VEnergyLossProcess*> enloss_vector =
1108 std::vector<G4VMultipleScattering*> mscat_vector =
1110
1111 for (auto theParticle : particles) {
1112 outFile << G4endl << "**" << theParticle->GetParticleName()
1113 << "**" << G4endl << G4endl << " .. code-block:: none" << G4endl;
1114
1115 G4ProcessManager* pm = theParticle->GetProcessManager();
1116 G4ProcessVector* pv = pm->GetProcessList();
1117 G4int plen = pm->GetProcessListLength();
1118
1119 for (auto emproc : emproc_vector) {
1120 for (G4int i = 0; i < plen; ++i) {
1121 G4VProcess* proc = (*pv)[i];
1122 if (proc == emproc) {
1123 outFile << G4endl;
1124 proc->ProcessDescription(outFile);
1125 break;
1126 }
1127 }
1128 }
1129
1130 for (auto mscproc : mscat_vector) {
1131 for (G4int i = 0; i < plen; ++i) {
1132 G4VProcess* proc = (*pv)[i];
1133 if (proc == mscproc) {
1134 outFile << G4endl;
1135 proc->ProcessDescription(outFile);
1136 break;
1137 }
1138 }
1139 }
1140
1141 for (auto enlossproc : enloss_vector) {
1142 for (G4int i = 0; i < plen; ++i) {
1143 G4VProcess* proc = (*pv)[i];
1144 if (proc == enlossproc) {
1145 outFile << G4endl;
1146 proc->ProcessDescription(outFile);
1147 break;
1148 }
1149 }
1150 }
1151 }
1152 outFile.close();
1153 }
1154}
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
const std::vector< G4VEmProcess * > & GetEmProcessVector()
const std::vector< G4VMultipleScattering * > & GetMultipleScatteringVector()
const std::vector< G4VEnergyLossProcess * > & GetEnergyLossProcessVector()
static G4MuonMinus * MuonMinusDefinition()
Definition: G4MuonMinus.cc:94
static G4MuonPlus * MuonPlusDefinition()
Definition: G4MuonPlus.cc:93
static G4Positron * Positron()
Definition: G4Positron.cc:93
G4int GetProcessListLength() const
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:87
virtual void ProcessDescription(std::ostream &outfile) const
Definition: G4VProcess.cc:175

References G4Electron::Electron(), G4endl, G4Gamma::Gamma(), GetEmProcessVector(), GetEnergyLossProcessVector(), GetMultipleScatteringVector(), G4ProcessManager::GetProcessList(), G4ProcessManager::GetProcessListLength(), G4MuonMinus::MuonMinusDefinition(), G4MuonPlus::MuonPlusDefinition(), G4Positron::Positron(), G4VProcess::ProcessDescription(), and G4Proton::ProtonDefinition().

◆ ElectronIonPair()

G4ElectronIonPair * G4LossTableManager::ElectronIonPair ( )

Definition at line 1009 of file G4LossTableManager.cc.

References emElectronIonPair, and verbose.

◆ EmConfigurator()

G4EmConfigurator * G4LossTableManager::EmConfigurator ( )

◆ EmCorrections()

G4EmCorrections * G4LossTableManager::EmCorrections ( )
inline

◆ EmSaturation()

G4EmSaturation * G4LossTableManager::EmSaturation ( )

Definition at line 992 of file G4LossTableManager.cc.

993{
995}
G4EmSaturation * GetEmSaturation()

References G4EmParameters::GetEmSaturation(), and theParameters.

Referenced by G4OpticalPhysics::ConstructProcess().

◆ GetCSDARange()

G4double G4LossTableManager::GetCSDARange ( const G4ParticleDefinition aParticle,
G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 318 of file G4LossTableManager.hh.

321{
322 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
323 return currentLoss ? currentLoss->GetCSDARange(kineticEnergy, couple) : DBL_MAX;
324}
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
G4double GetCSDARange(G4double kineticEnergy, const G4MaterialCutsCouple *)
#define DBL_MAX
Definition: templates.hh:62

References currentLoss, currentParticle, DBL_MAX, G4VEnergyLossProcess::GetCSDARange(), and GetEnergyLossProcess().

Referenced by G4EmCalculator::GetCSDARange().

◆ GetDEDX()

G4double G4LossTableManager::GetDEDX ( const G4ParticleDefinition aParticle,
G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

◆ GetDEDXDispersion()

G4double G4LossTableManager::GetDEDXDispersion ( const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double length 
)
inline

Definition at line 363 of file G4LossTableManager.hh.

367{
368 const G4ParticleDefinition* aParticle = dp->GetParticleDefinition();
369 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
370 return currentLoss ? currentLoss->GetDEDXDispersion(couple, dp, length) : 0.0;
371}
const G4ParticleDefinition * GetParticleDefinition() const
G4double GetDEDXDispersion(const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double length)

References currentLoss, currentParticle, G4VEnergyLossProcess::GetDEDXDispersion(), GetEnergyLossProcess(), and G4DynamicParticle::GetParticleDefinition().

◆ GetElectronGeneralProcess()

G4VEmProcess * G4LossTableManager::GetElectronGeneralProcess ( )
inline

Definition at line 431 of file G4LossTableManager.hh.

432{
433 return eGeneral;
434}

References eGeneral.

Referenced by G4EmExtraPhysics::ConstructGammaElectroNuclear().

◆ GetEmProcessVector()

const std::vector< G4VEmProcess * > & G4LossTableManager::GetEmProcessVector ( )

Definition at line 977 of file G4LossTableManager.cc.

978{
979 return emp_vector;
980}

References emp_vector.

Referenced by DumpHtml(), and G4EmCalculator::FindDiscreteProcess().

◆ GetEnergy()

G4double G4LossTableManager::GetEnergy ( const G4ParticleDefinition aParticle,
G4double  range,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 352 of file G4LossTableManager.hh.

355{
356 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
357 return currentLoss ? currentLoss->GetKineticEnergy(range, couple) : 0.0;
358}
G4double GetKineticEnergy(G4double range, const G4MaterialCutsCouple *)

References currentLoss, currentParticle, GetEnergyLossProcess(), and G4VEnergyLossProcess::GetKineticEnergy().

Referenced by G4EmCalculator::GetKinEnergy(), and G4EnergyLossTables::GetPreciseEnergyFromRange().

◆ GetEnergyLossProcess()

G4VEnergyLossProcess * G4LossTableManager::GetEnergyLossProcess ( const G4ParticleDefinition aParticle)

Definition at line 417 of file G4LossTableManager.cc.

418{
419 //G4cout << "G4LossTableManager::GetEnergyLossProcess: "
420 //<< aParticle << " " << currentParticle << " " << currentLoss << G4endl;
421 if(aParticle != currentParticle) {
422 currentParticle = aParticle;
423 std::map<PD,G4VEnergyLossProcess*,std::less<PD> >::const_iterator pos;
424 if ((pos = loss_map.find(aParticle)) != loss_map.end()) {
425 currentLoss = (*pos).second;
426 } else {
427 currentLoss = nullptr;
428 if ((pos = loss_map.find(theGenericIon)) != loss_map.end()) {
429 currentLoss = (*pos).second;
430 }
431 }
432 }
433 return currentLoss;
434}
static const G4double pos

References currentLoss, currentParticle, loss_map, pos, and theGenericIon.

Referenced by G4EmBiasingManager::ApplyRangeCut(), G4EmCalculator::FindEnergyLossProcess(), GetCSDARange(), GetDEDX(), GetDEDXDispersion(), GetEnergy(), GetRange(), GetRangeFromRestricteDEDX(), and G4VMultipleScattering::StartTracking().

◆ GetEnergyLossProcessVector()

const std::vector< G4VEnergyLossProcess * > & G4LossTableManager::GetEnergyLossProcessVector ( )

◆ GetGammaGeneralProcess()

G4VEmProcess * G4LossTableManager::GetGammaGeneralProcess ( )
inline

◆ GetMultipleScatteringVector()

const std::vector< G4VMultipleScattering * > & G4LossTableManager::GetMultipleScatteringVector ( )

Definition at line 985 of file G4LossTableManager.cc.

986{
987 return msc_vector;
988}

References msc_vector.

Referenced by DumpHtml(), and G4EmCalculator::FindMscProcess().

◆ GetPositronGeneralProcess()

G4VEmProcess * G4LossTableManager::GetPositronGeneralProcess ( )
inline

Definition at line 445 of file G4LossTableManager.hh.

446{
447 return pGeneral;
448}

Referenced by G4EmExtraPhysics::ConstructGammaElectroNuclear().

◆ GetRange()

G4double G4LossTableManager::GetRange ( const G4ParticleDefinition aParticle,
G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

◆ GetRangeFromRestricteDEDX()

G4double G4LossTableManager::GetRangeFromRestricteDEDX ( const G4ParticleDefinition aParticle,
G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 329 of file G4LossTableManager.hh.

333{
334 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
335 return currentLoss ? currentLoss->GetRange(kineticEnergy, couple) : DBL_MAX;
336}

References currentLoss, currentParticle, DBL_MAX, GetEnergyLossProcess(), and G4VEnergyLossProcess::GetRange().

Referenced by G4EmCalculator::GetRangeFromRestricteDEDX().

◆ GetTableBuilder()

G4LossTableBuilder * G4LossTableManager::GetTableBuilder ( )
inline

◆ Instance()

G4LossTableManager * G4LossTableManager::Instance ( void  )
static

Definition at line 87 of file G4LossTableManager.cc.

88{
89 if(!instance) {
91 instance = inst.Instance();
92 }
93 return instance;
94}
static G4ThreadLocal G4LossTableManager * instance

References G4ThreadLocalSingleton< T >::Instance(), and instance.

Referenced by G4EmModelActivator::ActivateEmOptions(), G4EmModelActivator::ActivateMicroElec(), G4EmModelActivator::ActivatePAI(), G4EmDNAPhysicsActivator::AddElectronModels0(), G4EmDNAPhysicsActivator::AddElectronModels2(), G4EmDNAPhysicsActivator::AddElectronModels4(), G4EmDNAPhysicsActivator::AddElectronModels4a(), G4EmDNAPhysicsActivator::AddElectronModels6(), G4EmDNAPhysicsActivator::AddElectronModels6a(), G4EmDNAPhysicsActivator::AddElectronModels7(), G4EmDNAPhysicsActivator::AddGenericIonModels0(), G4EmDNAPhysicsActivator::AddHeliumModels0(), G4EmDNAPhysicsActivator::AddProtonModels0(), G4ITStepProcessor::ApplyProductionCut(), G4SteppingManager::ApplyProductionCut(), G4EmBiasingManager::ApplyRangeCut(), G4BertiniElectroNuclearBuilder::Build(), G4GammaGeneralProcess::BuildPhysicsTable(), G4EmCalculator::ComputeDEDXForCutInRange(), G4EmCalculator::ComputeElectronicDEDX(), G4EmExtraPhysics::ConstructGammaElectroNuclear(), LBE::ConstructGeneral(), G4EmDNAPhysics_option1::ConstructProcess(), G4EmDNAPhysics_option2::ConstructProcess(), G4EmDNAPhysics_option3::ConstructProcess(), G4EmDNAPhysics_option4::ConstructProcess(), G4EmDNAPhysics_option5::ConstructProcess(), G4EmDNAPhysics_option6::ConstructProcess(), G4EmDNAPhysics_option7::ConstructProcess(), G4EmDNAPhysics_option8::ConstructProcess(), G4EmDNAPhysics_stationary_option2::ConstructProcess(), G4EmDNAPhysics_stationary_option4::ConstructProcess(), G4EmDNAPhysics_stationary_option6::ConstructProcess(), G4EmExtraPhysics::ConstructProcess(), G4RadioactiveDecayPhysics::ConstructProcess(), G4EmDNAPhysicsActivator::ConstructProcess(), G4EmStandardPhysics::ConstructProcess(), G4EmStandardPhysics_option1::ConstructProcess(), G4EmStandardPhysics_option2::ConstructProcess(), G4EmStandardPhysics_option3::ConstructProcess(), G4EmStandardPhysics_option4::ConstructProcess(), G4EmStandardPhysicsSS::ConstructProcess(), G4OpticalPhysics::ConstructProcess(), G4ECDecay::DecayIt(), G4ITDecay::DecayIt(), G4AnnihiToMuPair::G4AnnihiToMuPair(), G4AtimaEnergyLossModel::G4AtimaEnergyLossModel(), G4BetheBlochModel::G4BetheBlochModel(), G4BraggModel::G4BraggModel(), G4EmCalculator::G4EmCalculator(), G4GammaConversionToMuons::G4GammaConversionToMuons(), G4GoudsmitSaundersonMscModel::G4GoudsmitSaundersonMscModel(), G4ionIonisation::G4ionIonisation(), G4IonParametrisedLossModel::G4IonParametrisedLossModel(), G4LindhardSorensenIonModel::G4LindhardSorensenIonModel(), G4MuBetheBlochModel::G4MuBetheBlochModel(), G4NIELCalculator::G4NIELCalculator(), G4SynchrotronRadiation::G4SynchrotronRadiation(), G4UAtomicDeexcitation::G4UAtomicDeexcitation(), G4UrbanAdjointMscModel::G4UrbanAdjointMscModel(), G4UrbanMscModel::G4UrbanMscModel(), G4UserSpecialCuts::G4UserSpecialCuts(), G4VEmFluctuationModel::G4VEmFluctuationModel(), G4VEmModel::G4VEmModel(), G4VEmProcess::G4VEmProcess(), G4VEnergyLossProcess::G4VEnergyLossProcess(), G4VMultipleScattering::G4VMultipleScattering(), G4VTransitionRadiation::G4VTransitionRadiation(), G4EnergyLossTables::GetDEDX(), G4VMscModel::GetParticleChangeForMSC(), G4EnergyLossTables::GetPreciseDEDX(), G4EnergyLossTables::GetPreciseEnergyFromRange(), G4EnergyLossTables::GetPreciseRangeFromEnergy(), G4EnergyLossTables::GetRange(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermorePhotoElectricModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4PenelopeComptonModel::Initialise(), G4PenelopeIonisationModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4BraggIonModel::Initialise(), G4KleinNishinaModel::Initialise(), G4PEEffectFluoModel::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNACPA100IonisationModel::Initialise(), G4DNAEmfietzoglouIonisationModel::Initialise(), G4DNARelativisticIonisationModel::Initialise(), G4GammaGeneralProcess::InitialiseProcess(), PhysicsList::PhysicsList(), G4Cerenkov::PostStepGetPhysicalInteractionLength(), G4EmBuilder::PrepareEMPhysics(), G4GammaGeneralProcess::PreparePhysicsTable(), G4ContinuousGainOfEnergy::SetDynamicMassCharge(), and G4EmSaturation::VisibleEnergyDeposition().

◆ IsMaster()

G4bool G4LossTableManager::IsMaster ( ) const
inline

◆ LocalPhysicsTables()

void G4LossTableManager::LocalPhysicsTables ( const G4ParticleDefinition aParticle,
G4VEnergyLossProcess p 
)

Definition at line 547 of file G4LossTableManager.cc.

550{
551 if(1 < verbose) {
552 G4cout << "### G4LossTableManager::LocalPhysicsTable() for "
553 << aParticle->GetParticleName()
554 << " and process " << p->GetProcessName()
555 << G4endl;
556 }
557
558 if(-1 == run && startInitialisation) {
560 firstParticle = aParticle;
561 }
562
564 ++run;
565 if(1 < verbose) {
566 G4cout << "===== G4LossTableManager::LocalPhysicsTable() for run "
567 << run << " =====" << G4endl;
568 }
569 currentParticle = nullptr;
570 startInitialisation = false;
571 for (G4int i=0; i<n_loss; ++i) {
572 if(loss_vector[i]) {
573 tables_are_built[i] = false;
574 } else {
575 tables_are_built[i] = true;
576 part_vector[i] = nullptr;
577 }
578 }
579 }
580
582 for (G4int i=0; i<n_loss; ++i) {
583 if(p == loss_vector[i]) {
584 tables_are_built[i] = true;
585 isActive[i] = true;
586 part_vector[i] = p->Particle();
588 dedx_vector[i] = p->DEDXTable();
591 if(0 == run && p->IsIonisationProcess()) {
592 loss_map[part_vector[i]] = p;
593 //G4cout << "G4LossTableManager::LocalPhysicsTable " << part_vector[i]->GetParticleName()
594 // << " added to map " << p << G4endl;
595 }
596
597 if(1 < verbose) {
598 G4cout << i <<". "<< p->GetProcessName();
599 if(part_vector[i]) {
600 G4cout << " for " << part_vector[i]->GetParticleName();
601 }
602 G4cout << " active= " << isActive[i]
603 << " table= " << tables_are_built[i]
604 << " isIonisation= " << p->IsIonisationProcess()
605 << G4endl;
606 }
607 break;
608 } else if(!tables_are_built[i]) {
609 all_tables_are_built = false;
610 }
611 }
612
613 if(1 < verbose) {
614 G4cout << "### G4LossTableManager::LocalPhysicsTable end"
615 << G4endl;
616 }
618 if(1 < verbose) {
619 G4cout << "%%%%% All dEdx and Range tables for worker are ready for run "
620 << run << " %%%%%" << G4endl;
621 }
622 }
623}

References all_tables_are_built, base_part_vector, G4VEnergyLossProcess::BaseParticle(), G4EmConfigurator::Clear(), currentParticle, dedx_vector, G4VEnergyLossProcess::DEDXTable(), emConfigurator, firstParticle, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), inv_range_vector, G4VEnergyLossProcess::InverseRangeTable(), isActive, G4VEnergyLossProcess::IsIonisationProcess(), loss_map, loss_vector, n_loss, part_vector, G4VEnergyLossProcess::Particle(), range_vector, G4VEnergyLossProcess::RangeTableForLoss(), run, startInitialisation, tables_are_built, and verbose.

Referenced by G4VEnergyLossProcess::BuildPhysicsTable().

◆ NIELCalculator()

G4NIELCalculator * G4LossTableManager::NIELCalculator ( )

Definition at line 1029 of file G4LossTableManager.cc.

1030{
1031 if(!nielCalculator) {
1032 nielCalculator = new G4NIELCalculator(nullptr, verbose);
1033 }
1034 return nielCalculator;
1035}

References nielCalculator, and verbose.

◆ operator=()

G4LossTableManager & G4LossTableManager::operator= ( const G4LossTableManager right)
delete

◆ ParticleHaveNoLoss()

void G4LossTableManager::ParticleHaveNoLoss ( const G4ParticleDefinition aParticle)
private

Definition at line 950 of file G4LossTableManager.cc.

952{
954 ed << "Energy loss process not found for " << aParticle->GetParticleName()
955 << " !";
956 G4Exception("G4LossTableManager::ParticleHaveNoLoss", "em0001",
957 FatalException, ed);
958}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40

References FatalException, G4Exception(), and G4ParticleDefinition::GetParticleName().

◆ PreparePhysicsTable() [1/3]

void G4LossTableManager::PreparePhysicsTable ( const G4ParticleDefinition aParticle,
G4VEmProcess p,
G4bool  theMaster 
)

Definition at line 480 of file G4LossTableManager.cc.

482{
483 if (1 < verbose) {
484 G4cout << "G4LossTableManager::PreparePhysicsTable for "
485 << particle->GetParticleName()
486 << " and " << p->GetProcessName() << G4endl;
487 }
488 isMaster = theMaster;
489
490 if(!startInitialisation) {
492 if (1 < verbose) {
493 G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
494 << G4endl;
495 }
496 }
497
498 // start initialisation for the first run
499 if( -1 == run ) {
500 if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
501 }
502 startInitialisation = true;
503}
void PrepareModels(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)

References emConfigurator, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), isMaster, G4EmConfigurator::PrepareModels(), ResetParameters(), startInitialisation, and verbose.

◆ PreparePhysicsTable() [2/3]

void G4LossTableManager::PreparePhysicsTable ( const G4ParticleDefinition aParticle,
G4VEnergyLossProcess p,
G4bool  theMaster 
)

Definition at line 439 of file G4LossTableManager.cc.

442{
443 if (1 < verbose) {
444 G4cout << "G4LossTableManager::PreparePhysicsTable for "
445 << particle->GetParticleName()
446 << " and " << p->GetProcessName() << " run= " << run
447 << " loss_vector " << loss_vector.size() << G4endl;
448 }
449
450 isMaster = theMaster;
451
452 if(!startInitialisation) {
454 if (1 < verbose) {
455 G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
456 << G4endl;
457 }
458 }
459
460 // start initialisation for the first run
461 if( -1 == run ) {
462 if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
463
464 // initialise particles for given process
465 for (G4int j=0; j<n_loss; ++j) {
466 if (p == loss_vector[j] && !part_vector[j]) {
467 part_vector[j] = particle;
468 if(particle->GetParticleName() == "GenericIon") {
469 theGenericIon = particle;
470 }
471 }
472 }
473 }
474 startInitialisation = true;
475}

References emConfigurator, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), isMaster, loss_vector, n_loss, part_vector, G4EmConfigurator::PrepareModels(), ResetParameters(), startInitialisation, theGenericIon, and verbose.

Referenced by G4VEmProcess::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), and G4VMultipleScattering::PreparePhysicsTable().

◆ PreparePhysicsTable() [3/3]

void G4LossTableManager::PreparePhysicsTable ( const G4ParticleDefinition aParticle,
G4VMultipleScattering p,
G4bool  theMaster 
)

Definition at line 508 of file G4LossTableManager.cc.

511{
512 if (1 < verbose) {
513 G4cout << "G4LossTableManager::PreparePhysicsTable for "
514 << particle->GetParticleName()
515 << " and " << p->GetProcessName() << G4endl;
516 }
517
518 isMaster = theMaster;
519
520 if(!startInitialisation) {
522 if (1 < verbose) {
523 G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
524 << G4endl;
525 }
526 }
527
528 // start initialisation for the first run
529 if( -1 == run ) {
530 if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
531 }
532 startInitialisation = true;
533}

References emConfigurator, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), isMaster, G4EmConfigurator::PrepareModels(), ResetParameters(), startInitialisation, and verbose.

◆ PrintEWarning()

void G4LossTableManager::PrintEWarning ( G4String  tit,
G4double   
)
private

Definition at line 1059 of file G4LossTableManager.cc.

1060{
1061 G4String ss = "G4LossTableManager::" + tit;
1063 /*
1064 ed << "Parameter is out of range: " << val
1065 << " it will have no effect!\n" << " ## "
1066 << " nbins= " << nbinsLambda
1067 << " nbinsPerDecade= " << nbinsPerDecade
1068 << " Emin(keV)= " << minKinEnergy/keV
1069 << " Emax(GeV)= " << maxKinEnergy/GeV;
1070 */
1071 G4Exception(ss, "em0044", JustWarning, ed);
1072}
@ JustWarning

References G4Exception(), and JustWarning.

◆ Register() [1/6]

void G4LossTableManager::Register ( G4VEmFluctuationModel p)

Definition at line 369 of file G4LossTableManager.cc.

370{
371 fmod_vector.push_back(p);
372 if(verbose > 1) {
373 G4cout << "G4LossTableManager::Register G4VEmFluctuationModel : "
374 << p->GetName() << " " << fmod_vector.size() << G4endl;
375 }
376}
const G4String & GetName() const

References fmod_vector, G4cout, G4endl, G4VEmFluctuationModel::GetName(), and verbose.

◆ Register() [2/6]

void G4LossTableManager::Register ( G4VEmModel p)

Definition at line 344 of file G4LossTableManager.cc.

345{
346 mod_vector.push_back(p);
347 if(verbose > 1) {
348 G4cout << "G4LossTableManager::Register G4VEmModel : "
349 << p->GetName() << " " << p << " " << mod_vector.size() << G4endl;
350 }
351}
const G4String & GetName() const
Definition: G4VEmModel.hh:837

References G4cout, G4endl, G4VEmModel::GetName(), mod_vector, and verbose.

◆ Register() [3/6]

void G4LossTableManager::Register ( G4VEmProcess p)

Definition at line 284 of file G4LossTableManager.cc.

285{
286 if(!p) { return; }
287 G4int n = emp_vector.size();
288 for (G4int i=0; i<n; ++i) {
289 if(emp_vector[i] == p) { return; }
290 }
291 if(verbose > 1) {
292 G4cout << "G4LossTableManager::Register G4VEmProcess : "
293 << p->GetProcessName() << " idx= " << emp_vector.size() << G4endl;
294 }
295 emp_vector.push_back(p);
296}

References emp_vector, G4cout, G4endl, G4VProcess::GetProcessName(), CLHEP::detail::n, and verbose.

◆ Register() [4/6]

void G4LossTableManager::Register ( G4VEnergyLossProcess p)

Definition at line 197 of file G4LossTableManager.cc.

198{
199 if(!p) { return; }
200 for (G4int i=0; i<n_loss; ++i) {
201 if(loss_vector[i] == p) { return; }
202 }
203 if(verbose > 1) {
204 G4cout << "G4LossTableManager::Register G4VEnergyLossProcess : "
205 << p->GetProcessName() << " idx= " << n_loss << G4endl;
206 }
207 ++n_loss;
208 loss_vector.push_back(p);
209 part_vector.push_back(nullptr);
210 base_part_vector.push_back(nullptr);
211 dedx_vector.push_back(nullptr);
212 range_vector.push_back(nullptr);
213 inv_range_vector.push_back(nullptr);
214 tables_are_built.push_back(false);
215 isActive.push_back(true);
216 all_tables_are_built = false;
217}

References all_tables_are_built, base_part_vector, dedx_vector, G4cout, G4endl, G4VProcess::GetProcessName(), inv_range_vector, isActive, loss_vector, n_loss, part_vector, range_vector, tables_are_built, and verbose.

Referenced by G4AnnihiToMuPair::G4AnnihiToMuPair(), G4GammaConversionToMuons::G4GammaConversionToMuons(), G4SynchrotronRadiation::G4SynchrotronRadiation(), G4VEmFluctuationModel::G4VEmFluctuationModel(), G4VEmModel::G4VEmModel(), G4VEmProcess::G4VEmProcess(), G4VEnergyLossProcess::G4VEnergyLossProcess(), G4VMultipleScattering::G4VMultipleScattering(), and G4VTransitionRadiation::G4VTransitionRadiation().

◆ Register() [5/6]

void G4LossTableManager::Register ( G4VMultipleScattering p)

Definition at line 254 of file G4LossTableManager.cc.

255{
256 if(!p) { return; }
257 G4int n = msc_vector.size();
258 for (G4int i=0; i<n; ++i) {
259 if(msc_vector[i] == p) { return; }
260 }
261 if(verbose > 1) {
262 G4cout << "G4LossTableManager::Register G4VMultipleScattering : "
263 << p->GetProcessName() << " idx= " << msc_vector.size() << G4endl;
264 }
265 msc_vector.push_back(p);
266}

References G4cout, G4endl, G4VProcess::GetProcessName(), msc_vector, CLHEP::detail::n, and verbose.

◆ Register() [6/6]

void G4LossTableManager::Register ( G4VProcess p)

Definition at line 314 of file G4LossTableManager.cc.

315{
316 if(!p) { return; }
317 G4int n = p_vector.size();
318 for (G4int i=0; i<n; ++i) {
319 if(p_vector[i] == p) { return; }
320 }
321 if(verbose > 1) {
322 G4cout << "G4LossTableManager::Register G4VProcess : "
323 << p->GetProcessName() << " idx= " << p_vector.size() << G4endl;
324 }
325 p_vector.push_back(p);
326}

References G4cout, G4endl, G4VProcess::GetProcessName(), CLHEP::detail::n, p_vector, and verbose.

◆ RegisterExtraParticle()

void G4LossTableManager::RegisterExtraParticle ( const G4ParticleDefinition aParticle,
G4VEnergyLossProcess p 
)

Definition at line 390 of file G4LossTableManager.cc.

393{
394 if(!p || !part) { return; }
395 for (G4int i=0; i<n_loss; ++i) {
396 if(loss_vector[i] == p) { return; }
397 }
398 if(verbose > 1) {
399 G4cout << "G4LossTableManager::RegisterExtraParticle "
400 << part->GetParticleName() << " G4VEnergyLossProcess : "
401 << p->GetProcessName() << " idx= " << n_loss << G4endl;
402 }
403 ++n_loss;
404 loss_vector.push_back(p);
405 part_vector.push_back(part);
406 base_part_vector.push_back(p->BaseParticle());
407 dedx_vector.push_back(nullptr);
408 range_vector.push_back(nullptr);
409 inv_range_vector.push_back(nullptr);
410 tables_are_built.push_back(false);
411 all_tables_are_built = false;
412}

References all_tables_are_built, base_part_vector, G4VEnergyLossProcess::BaseParticle(), dedx_vector, G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), inv_range_vector, loss_vector, n_loss, part_vector, range_vector, tables_are_built, and verbose.

Referenced by G4VEnergyLossProcess::PreparePhysicsTable().

◆ ResetParameters()

void G4LossTableManager::ResetParameters ( )

Definition at line 221 of file G4LossTableManager.cc.

222{
224 if(!isMaster) {
226 } else {
227 if(verbose > 0) { theParameters->Dump(); }
228 }
231 if(nullptr != emConfigurator) { emConfigurator->SetVerbose(verbose); };
233 if(nullptr != atomDeexcitation) {
236 }
237}
void SetVerbose(G4int value)
void SetVerbose(G4int verb)
void SetInitialisationFlag(G4bool flag)

References atomDeexcitation, G4EmParameters::Dump(), emConfigurator, emCorrections, emElectronIonPair, G4VAtomDeexcitation::InitialiseAtomicDeexcitation(), isMaster, G4LossTableBuilder::SetInitialisationFlag(), G4EmConfigurator::SetVerbose(), G4EmCorrections::SetVerbose(), G4ElectronIonPair::SetVerbose(), G4VAtomDeexcitation::SetVerboseLevel(), tableBuilder, theParameters, G4EmParameters::Verbose(), verbose, and G4EmParameters::WorkerVerbose().

Referenced by G4RadioactiveDecayPhysics::ConstructProcess(), and PreparePhysicsTable().

◆ SetAtomDeexcitation()

void G4LossTableManager::SetAtomDeexcitation ( G4VAtomDeexcitation p)

◆ SetElectronGeneralProcess()

void G4LossTableManager::SetElectronGeneralProcess ( G4VEmProcess ptr)
inline

Definition at line 424 of file G4LossTableManager.hh.

425{
426 eGeneral = ptr;
427}

References eGeneral.

◆ SetGammaGeneralProcess()

void G4LossTableManager::SetGammaGeneralProcess ( G4VEmProcess ptr)
inline

◆ SetNIELCalculator()

void G4LossTableManager::SetNIELCalculator ( G4NIELCalculator ptr)

Definition at line 1019 of file G4LossTableManager.cc.

1020{
1021 if(ptr && ptr != nielCalculator) {
1022 delete nielCalculator;
1023 nielCalculator = ptr;
1024 }
1025}

References nielCalculator.

Referenced by G4NIELCalculator::G4NIELCalculator().

◆ SetPositronGeneralProcess()

void G4LossTableManager::SetPositronGeneralProcess ( G4VEmProcess ptr)
inline

Definition at line 438 of file G4LossTableManager.hh.

439{
440 pGeneral = ptr;
441}

References pGeneral.

◆ SetSubCutProducer()

void G4LossTableManager::SetSubCutProducer ( G4VSubCutProducer p)

Definition at line 1049 of file G4LossTableManager.cc.

1050{
1051 if(subcutProducer != p) {
1052 delete subcutProducer;
1053 subcutProducer = p;
1054 }
1055}

References subcutProducer.

◆ SetVerbose()

void G4LossTableManager::SetVerbose ( G4int  val)

Definition at line 962 of file G4LossTableManager.cc.

963{
964 verbose = val;
965}

References verbose.

◆ SubCutProducer()

G4VSubCutProducer * G4LossTableManager::SubCutProducer ( )
inline

Definition at line 396 of file G4LossTableManager.hh.

397{
398 return subcutProducer;
399}

References subcutProducer.

Referenced by G4VEnergyLossProcess::PreparePhysicsTable().

Friends And Related Function Documentation

◆ G4ThreadLocalSingleton< G4LossTableManager >

Definition at line 445 of file G4LossTableManager.hh.

Field Documentation

◆ all_tables_are_built

G4bool G4LossTableManager::all_tables_are_built
private

◆ atomDeexcitation

G4VAtomDeexcitation* G4LossTableManager::atomDeexcitation
private

◆ base_part_vector

std::vector<PD> G4LossTableManager::base_part_vector
private

◆ currentLoss

G4VEnergyLossProcess* G4LossTableManager::currentLoss
private

◆ currentParticle

PD G4LossTableManager::currentParticle
private

◆ dedx_vector

std::vector<G4PhysicsTable*> G4LossTableManager::dedx_vector
private

◆ eGeneral

G4VEmProcess* G4LossTableManager::eGeneral
private

◆ emConfigurator

G4EmConfigurator* G4LossTableManager::emConfigurator
private

◆ emCorrections

G4EmCorrections* G4LossTableManager::emCorrections
private

◆ emElectronIonPair

G4ElectronIonPair* G4LossTableManager::emElectronIonPair
private

◆ emp_vector

std::vector<G4VEmProcess*> G4LossTableManager::emp_vector
private

◆ firstParticle

PD G4LossTableManager::firstParticle
private

◆ fmod_vector

std::vector<G4VEmFluctuationModel*> G4LossTableManager::fmod_vector
private

Definition at line 271 of file G4LossTableManager.hh.

Referenced by DeRegister(), Register(), and ~G4LossTableManager().

◆ gGeneral

G4VEmProcess* G4LossTableManager::gGeneral
private

◆ instance

G4ThreadLocal G4LossTableManager * G4LossTableManager::instance = nullptr
staticprivate

Definition at line 254 of file G4LossTableManager.hh.

Referenced by Instance().

◆ inv_range_vector

std::vector<G4PhysicsTable*> G4LossTableManager::inv_range_vector
private

◆ isActive

std::vector<G4bool> G4LossTableManager::isActive
private

◆ isMaster

G4bool G4LossTableManager::isMaster
private

◆ loss_map

std::map<PD,G4VEnergyLossProcess*,std::less<PD> > G4LossTableManager::loss_map
private

◆ loss_vector

std::vector<G4VEnergyLossProcess*> G4LossTableManager::loss_vector
private

◆ mod_vector

std::vector<G4VEmModel*> G4LossTableManager::mod_vector
private

Definition at line 270 of file G4LossTableManager.hh.

Referenced by DeRegister(), Register(), and ~G4LossTableManager().

◆ msc_vector

std::vector<G4VMultipleScattering*> G4LossTableManager::msc_vector
private

◆ n_loss

G4int G4LossTableManager::n_loss
private

◆ nielCalculator

G4NIELCalculator* G4LossTableManager::nielCalculator
private

◆ p_vector

std::vector<G4VProcess*> G4LossTableManager::p_vector
private

Definition at line 272 of file G4LossTableManager.hh.

Referenced by DeRegister(), Register(), and ~G4LossTableManager().

◆ part_vector

std::vector<PD> G4LossTableManager::part_vector
private

◆ pGeneral

G4VEmProcess* G4LossTableManager::pGeneral
private

Definition at line 292 of file G4LossTableManager.hh.

Referenced by G4LossTableManager(), and SetPositronGeneralProcess().

◆ range_vector

std::vector<G4PhysicsTable*> G4LossTableManager::range_vector
private

◆ run

G4int G4LossTableManager::run
private

Definition at line 296 of file G4LossTableManager.hh.

Referenced by BuildPhysicsTable(), and LocalPhysicsTables().

◆ startInitialisation

G4bool G4LossTableManager::startInitialisation
private

◆ subcutProducer

G4VSubCutProducer* G4LossTableManager::subcutProducer
private

◆ tableBuilder

G4LossTableBuilder* G4LossTableManager::tableBuilder
private

◆ tables_are_built

std::vector<G4bool> G4LossTableManager::tables_are_built
private

◆ theElectron

PD G4LossTableManager::theElectron
private

Definition at line 277 of file G4LossTableManager.hh.

Referenced by CopyTables(), and G4LossTableManager().

◆ theGenericIon

PD G4LossTableManager::theGenericIon
private

◆ theParameters

G4EmParameters* G4LossTableManager::theParameters
private

◆ verbose

G4int G4LossTableManager::verbose
private

The documentation for this class was generated from the following files: