Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4LossTableManager Class Reference

#include <G4LossTableManager.hh>

Public Member Functions

 ~G4LossTableManager ()
 
void Clear ()
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VEmProcess *p, G4bool theMaster)
 
void PreparePhysicsTable (const G4ParticleDefinition *aParticle, G4VMultipleScattering *p, G4bool theMaster)
 
void BuildPhysicsTable (const G4ParticleDefinition *aParticle)
 
void BuildPhysicsTable (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
void LocalPhysicsTables (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
G4double GetDEDX (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetSubDEDX (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetRange (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetCSDARange (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetRangeFromRestricteDEDX (const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double GetEnergy (const G4ParticleDefinition *aParticle, G4double range, const G4MaterialCutsCouple *couple)
 
G4double GetDEDXDispersion (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double &length)
 
void Register (G4VEnergyLossProcess *p)
 
void DeRegister (G4VEnergyLossProcess *p)
 
void Register (G4VMultipleScattering *p)
 
void DeRegister (G4VMultipleScattering *p)
 
void Register (G4VEmProcess *p)
 
void DeRegister (G4VEmProcess *p)
 
void Register (G4VEmModel *p)
 
void DeRegister (G4VEmModel *p)
 
void Register (G4VEmFluctuationModel *p)
 
void DeRegister (G4VEmFluctuationModel *p)
 
void RegisterExtraParticle (const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
 
void SetVerbose (G4int val)
 
void SetAtomDeexcitation (G4VAtomDeexcitation *)
 
void SetSubCutProducer (G4VSubCutProducer *)
 
G4bool IsMaster () const
 
G4VEnergyLossProcessGetEnergyLossProcess (const G4ParticleDefinition *)
 
const std::vector
< G4VEnergyLossProcess * > & 
GetEnergyLossProcessVector ()
 
const std::vector
< G4VEmProcess * > & 
GetEmProcessVector ()
 
const std::vector
< G4VMultipleScattering * > & 
GetMultipleScatteringVector ()
 
G4EmSaturationEmSaturation ()
 
G4EmConfiguratorEmConfigurator ()
 
G4ElectronIonPairElectronIonPair ()
 
G4EmCorrectionsEmCorrections ()
 
G4VAtomDeexcitationAtomDeexcitation ()
 
G4VSubCutProducerSubCutProducer ()
 
G4LossTableBuilderGetTableBuilder ()
 

Static Public Member Functions

static G4LossTableManagerInstance ()
 

Friends

class G4ThreadLocalSingleton< G4LossTableManager >
 

Detailed Description

Definition at line 100 of file G4LossTableManager.hh.

Constructor & Destructor Documentation

G4LossTableManager::~G4LossTableManager ( )

Definition at line 123 of file G4LossTableManager.cc.

124 {
125  //G4cout << "### G4LossTableManager::~G4LossTableManager() "<< this << G4endl;
126  for (G4int i=0; i<n_loss; ++i) {
127  //G4cout << "### eloss #" << i << G4endl;
128  if( loss_vector[i] ) {
129  delete loss_vector[i];
130  }
131  }
132  size_t msc = msc_vector.size();
133  for (size_t j=0; j<msc; ++j) {
134  if( msc_vector[j] ) { delete msc_vector[j]; }
135  }
136  size_t emp = emp_vector.size();
137  for (size_t k=0; k<emp; ++k) {
138  if( emp_vector[k] ) { delete emp_vector[k]; }
139  }
140  size_t mod = mod_vector.size();
141  size_t fmod = fmod_vector.size();
142  //G4cout << " Nmod" << mod << " Nfluc= " << fmod << G4endl;
143  for (size_t a=0; a<mod; ++a) {
144  //G4cout << "Delete model #" << a << " " << mod_vector[a] << G4endl;
145  if( nullptr != mod_vector[a] ) {
146  for (size_t b=0; b<fmod; ++b) {
147  if((G4VEmModel*)(fmod_vector[b]) == mod_vector[a]) {
148  fmod_vector[b] = nullptr;
149  }
150  }
151  delete mod_vector[a];
152  mod_vector[a] = nullptr;
153  }
154  }
155  for (size_t b=0; b<fmod; ++b) {
156  if( fmod_vector[b] ) { delete fmod_vector[b]; }
157  }
158  Clear();
159  delete tableBuilder;
160  delete emCorrections;
161  delete emConfigurator;
162  delete emElectronIonPair;
163  delete atomDeexcitation;
164  delete subcutProducer;
165 }
int G4int
Definition: G4Types.hh:78

Here is the call graph for this function:

Member Function Documentation

G4VAtomDeexcitation * G4LossTableManager::AtomDeexcitation ( )
inline

Definition at line 444 of file G4LossTableManager.hh.

445 {
446  return atomDeexcitation;
447 }

Here is the caller graph for this function:

void G4LossTableManager::BuildPhysicsTable ( const G4ParticleDefinition aParticle)

Definition at line 499 of file G4LossTableManager.cc.

500 {
501  if(-1 == run && startInitialisation) {
502  if(emConfigurator) { emConfigurator->Clear(); }
503  }
504 }

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 588 of file G4LossTableManager.cc.

591 {
592  if(1 < verbose) {
593  G4cout << "### G4LossTableManager::BuildPhysicsTable() for "
594  << aParticle->GetParticleName()
595  << " and process " << p->GetProcessName() << G4endl;
596  }
597  // clear configurator
598  if(-1 == run && startInitialisation) {
599  if(emConfigurator) { emConfigurator->Clear(); }
600  firstParticle = aParticle;
601  }
602  if(startInitialisation) {
603  ++run;
604  if(1 < verbose) {
605  G4cout << "===== G4LossTableManager::BuildPhysicsTable() for run "
606  << run << " ===== " << atomDeexcitation << G4endl;
607  }
608  currentParticle = nullptr;
609  all_tables_are_built= true;
610  }
611 
612  // initialisation before any table is built
613  if ( startInitialisation && aParticle == firstParticle ) {
614 
615  startInitialisation = false;
616  if(1 < verbose) {
617  G4cout << "### G4LossTableManager start initilisation for first particle "
618  << firstParticle->GetParticleName()
619  << G4endl;
620  }
621  for (G4int i=0; i<n_loss; ++i) {
622  G4VEnergyLossProcess* el = loss_vector[i];
623 
624  if(el) {
625  isActive[i] = true;
626  base_part_vector[i] = el->BaseParticle();
627  tables_are_built[i] = false;
628  all_tables_are_built= false;
629  if(!isActive[i]) {
630  el->SetIonisation(false);
631  tables_are_built[i] = true;
632  }
633 
634  if(1 < verbose) {
635  G4cout << i <<". "<< el->GetProcessName();
636  if(el->Particle()) {
637  G4cout << " for " << el->Particle()->GetParticleName();
638  }
639  G4cout << " active= " << isActive[i]
640  << " table= " << tables_are_built[i]
641  << " isIonisation= " << el->IsIonisationProcess();
642  if(base_part_vector[i]) {
643  G4cout << " base particle "
644  << base_part_vector[i]->GetParticleName();
645  }
646  G4cout << G4endl;
647  }
648  } else {
649  tables_are_built[i] = true;
650  part_vector[i] = nullptr;
651  isActive[i] = false;
652  }
653  }
654  }
655 
656  if (all_tables_are_built) { return; }
657 
658  // Build tables for given particle
659  all_tables_are_built = true;
660 
661  for(G4int i=0; i<n_loss; ++i) {
662  if(p == loss_vector[i] && !tables_are_built[i] && !base_part_vector[i]) {
663  const G4ParticleDefinition* curr_part = part_vector[i];
664  if(1 < verbose) {
665  G4cout << "### Build Table for " << p->GetProcessName()
666  << " and " << curr_part->GetParticleName()
667  << " " << tables_are_built[i] << " " << base_part_vector[i] << G4endl;
668  }
669  G4VEnergyLossProcess* curr_proc = BuildTables(curr_part);
670  if(curr_proc) {
671  CopyTables(curr_part, curr_proc);
672  if(p == curr_proc && 0 == run && p->IsIonisationProcess()) {
673  loss_map[aParticle] = p;
674  //G4cout << "G4LossTableManager::BuildPhysicsTable: " << aParticle->GetParticleName()
675  // << " added to map " << p << G4endl;
676  }
677  }
678  }
679  if ( !tables_are_built[i] ) { all_tables_are_built = false; }
680  }
681  if(1 < verbose) {
682  G4cout << "### G4LossTableManager::BuildPhysicsTable end: "
683  << "all_tables_are_built= " << all_tables_are_built << " "
684  << aParticle->GetParticleName() << " proc: " << p << G4endl;
685  }
686  if(all_tables_are_built) {
687  if(1 < verbose) {
688  G4cout << "%%%%% All dEdx and Range tables are built for master run= "
689  << run << " %%%%%" << G4endl;
690  }
691  }
692 }
void SetIonisation(G4bool val)
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * BaseParticle() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
const G4ParticleDefinition * Particle() const
#define G4endl
Definition: G4ios.hh:61
G4bool IsIonisationProcess() const

Here is the call graph for this function:

void G4LossTableManager::Clear ( )

Definition at line 198 of file G4LossTableManager.cc.

199 {
200  all_tables_are_built = false;
201  currentLoss = nullptr;
202  currentParticle = nullptr;
203  if(n_loss)
204  {
205  dedx_vector.clear();
206  range_vector.clear();
207  inv_range_vector.clear();
208  loss_map.clear();
209  loss_vector.clear();
210  part_vector.clear();
211  base_part_vector.clear();
212  tables_are_built.clear();
213  isActive.clear();
214  n_loss = 0;
215  }
216 }

Here is the caller graph for this function:

void G4LossTableManager::DeRegister ( G4VEnergyLossProcess p)

Definition at line 263 of file G4LossTableManager.cc.

264 {
265  if(!p) { return; }
266  for (G4int i=0; i<n_loss; ++i) {
267  if(loss_vector[i] == p) { loss_vector[i] = nullptr; }
268  }
269 }
int G4int
Definition: G4Types.hh:78

Here is the caller graph for this function:

void G4LossTableManager::DeRegister ( G4VMultipleScattering p)

Definition at line 289 of file G4LossTableManager.cc.

290 {
291  if(!p) { return; }
292  size_t msc = msc_vector.size();
293  for (size_t i=0; i<msc; ++i) {
294  if(msc_vector[i] == p) { msc_vector[i] = nullptr; }
295  }
296 }
void G4LossTableManager::DeRegister ( G4VEmProcess p)

Definition at line 316 of file G4LossTableManager.cc.

317 {
318  if(!p) { return; }
319  size_t emp = emp_vector.size();
320  for (size_t i=0; i<emp; ++i) {
321  if(emp_vector[i] == p) { emp_vector[i] = nullptr; }
322  }
323 }
void G4LossTableManager::DeRegister ( G4VEmModel p)

Definition at line 338 of file G4LossTableManager.cc.

339 {
340  //G4cout << "G4LossTableManager::DeRegister G4VEmModel : " << p << G4endl;
341  size_t n = mod_vector.size();
342  for (size_t i=0; i<n; ++i) {
343  if(mod_vector[i] == p) {
344  mod_vector[i] = nullptr;
345  break;
346  }
347  }
348 }
void G4LossTableManager::DeRegister ( G4VEmFluctuationModel p)

Definition at line 363 of file G4LossTableManager.cc.

364 {
365  size_t n = fmod_vector.size();
366  for (size_t i=0; i<n; ++i) {
367  if(fmod_vector[i] == p) { fmod_vector[i] = nullptr; }
368  }
369 }
G4ElectronIonPair * G4LossTableManager::ElectronIonPair ( )

Definition at line 980 of file G4LossTableManager.cc.

981 {
982  if(!emElectronIonPair) {
983  emElectronIonPair = new G4ElectronIonPair(verbose);
984  }
985  return emElectronIonPair;
986 }
G4EmConfigurator * G4LossTableManager::EmConfigurator ( )

Definition at line 972 of file G4LossTableManager.cc.

973 {
974  if(!emConfigurator) { emConfigurator = new G4EmConfigurator(verbose); }
975  return emConfigurator;
976 }
G4EmCorrections * G4LossTableManager::EmCorrections ( )
inline

Definition at line 437 of file G4LossTableManager.hh.

438 {
439  return emCorrections;
440 }

Here is the caller graph for this function:

G4EmSaturation * G4LossTableManager::EmSaturation ( )

Definition at line 965 of file G4LossTableManager.cc.

966 {
967  return theParameters->GetEmSaturation();
968 }
G4EmSaturation * GetEmSaturation()

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 361 of file G4LossTableManager.hh.

364 {
365  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
366  G4double x = DBL_MAX;
367  if(currentLoss) { x = currentLoss->GetCSDARange(kineticEnergy, couple); }
368  return x;
369 }
G4double GetCSDARange(G4double &kineticEnergy, const G4MaterialCutsCouple *)
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 335 of file G4LossTableManager.hh.

338 {
339  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
340  G4double x = 0.0;
341  if(currentLoss) { x = currentLoss->GetDEDX(kineticEnergy, couple); }
342  return x;
343 }
G4double GetDEDX(G4double &kineticEnergy, const G4MaterialCutsCouple *)
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 416 of file G4LossTableManager.hh.

420 {
421  const G4ParticleDefinition* aParticle = dp->GetParticleDefinition();
422  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
423  G4double x = 0.0;
424  if(currentLoss) { currentLoss->GetDEDXDispersion(couple, dp, length); }
425  return x;
426 }
const G4ParticleDefinition * GetParticleDefinition() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
G4double GetDEDXDispersion(const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double length)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

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

Definition at line 950 of file G4LossTableManager.cc.

951 {
952  return emp_vector;
953 }
G4double G4LossTableManager::GetEnergy ( const G4ParticleDefinition aParticle,
G4double  range,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 403 of file G4LossTableManager.hh.

406 {
407  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
408  G4double x = 0;
409  if(currentLoss) { x = currentLoss->GetKineticEnergy(range, couple); }
410  return x;
411 }
G4double GetKineticEnergy(G4double &range, const G4MaterialCutsCouple *)
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4VEnergyLossProcess * G4LossTableManager::GetEnergyLossProcess ( const G4ParticleDefinition aParticle)
inline

Definition at line 313 of file G4LossTableManager.hh.

314 {
315  //G4cout << "G4LossTableManager::GetEnergyLossProcess: "
316  //<< aParticle << " " << currentParticle << " " << currentLoss << G4endl;
317  if(aParticle != currentParticle) {
318  currentParticle = aParticle;
319  std::map<PD,G4VEnergyLossProcess*,std::less<PD> >::const_iterator pos;
320  if ((pos = loss_map.find(aParticle)) != loss_map.end()) {
321  currentLoss = (*pos).second;
322  } else {
323  currentLoss = 0;
324  if ((pos = loss_map.find(theGenericIon)) != loss_map.end()) {
325  currentLoss = (*pos).second;
326  }
327  }
328  }
329  return currentLoss;
330 }
static const G4double pos

Here is the caller graph for this function:

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

Definition at line 943 of file G4LossTableManager.cc.

944 {
945  return loss_vector;
946 }

Here is the caller graph for this function:

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

Definition at line 958 of file G4LossTableManager.cc.

959 {
960  return msc_vector;
961 }
G4double G4LossTableManager::GetRange ( const G4ParticleDefinition aParticle,
G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 388 of file G4LossTableManager.hh.

391 {
392  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
393  G4double x = DBL_MAX;
394  if(currentLoss) {
395  x = currentLoss->GetRange(kineticEnergy, couple);
396  }
397  return x;
398 }
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
G4double GetRange(G4double &kineticEnergy, const G4MaterialCutsCouple *)

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 374 of file G4LossTableManager.hh.

378 {
379  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
380  G4double x = DBL_MAX;
381  if(currentLoss) { x = currentLoss->GetRangeForLoss(kineticEnergy, couple); }
382  return x;
383 }
G4double GetRangeForLoss(G4double &kineticEnergy, const G4MaterialCutsCouple *)
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 348 of file G4LossTableManager.hh.

351 {
352  if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
353  G4double x = 0.0;
354  if(currentLoss) { x = currentLoss->GetDEDXForSubsec(kineticEnergy, couple); }
355  return x;
356 }
G4double GetDEDXForSubsec(G4double &kineticEnergy, const G4MaterialCutsCouple *)
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4LossTableBuilder * G4LossTableManager::GetTableBuilder ( )
inline

Definition at line 458 of file G4LossTableManager.hh.

459 {
460  return tableBuilder;
461 }

Here is the caller graph for this function:

G4LossTableManager * G4LossTableManager::Instance ( void  )
static

Definition at line 112 of file G4LossTableManager.cc.

113 {
114  if(!instance) {
116  instance = inst.Instance();
117  }
118  return instance;
119 }

Here is the call graph for this function:

G4bool G4LossTableManager::IsMaster ( ) const
inline

Definition at line 430 of file G4LossTableManager.hh.

431 {
432  return isMaster;
433 }

Here is the caller graph for this function:

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

Definition at line 508 of file G4LossTableManager.cc.

511 {
512  if(1 < verbose) {
513  G4cout << "### G4LossTableManager::LocalPhysicsTable() for "
514  << aParticle->GetParticleName()
515  << " and process " << p->GetProcessName()
516  << G4endl;
517  }
518 
519  if(-1 == run && startInitialisation) {
520  if(emConfigurator) { emConfigurator->Clear(); }
521  firstParticle = aParticle;
522  }
523 
524  if(startInitialisation) {
525  ++run;
526  if(1 < verbose) {
527  G4cout << "===== G4LossTableManager::LocalPhysicsTable() for run "
528  << run << " =====" << G4endl;
529  }
530  currentParticle = nullptr;
531  startInitialisation = false;
532  for (G4int i=0; i<n_loss; ++i) {
533  if(loss_vector[i]) {
534  tables_are_built[i] = false;
535  } else {
536  tables_are_built[i] = true;
537  part_vector[i] = nullptr;
538  }
539  }
540  }
541 
542  all_tables_are_built= true;
543  for (G4int i=0; i<n_loss; ++i) {
544  if(p == loss_vector[i]) {
545  tables_are_built[i] = true;
546  isActive[i] = true;
547  part_vector[i] = p->Particle();
548  base_part_vector[i] = p->BaseParticle();
549  dedx_vector[i] = p->DEDXTable();
550  range_vector[i] = p->RangeTableForLoss();
551  inv_range_vector[i] = p->InverseRangeTable();
552  if(0 == run && p->IsIonisationProcess()) {
553  loss_map[part_vector[i]] = p;
554  //G4cout << "G4LossTableManager::LocalPhysicsTable " << part_vector[i]->GetParticleName()
555  // << " added to map " << p << G4endl;
556  }
557 
558  if(1 < verbose) {
559  G4cout << i <<". "<< p->GetProcessName();
560  if(part_vector[i]) {
561  G4cout << " for " << part_vector[i]->GetParticleName();
562  }
563  G4cout << " active= " << isActive[i]
564  << " table= " << tables_are_built[i]
565  << " isIonisation= " << p->IsIonisationProcess()
566  << G4endl;
567  }
568  break;
569  } else if(!tables_are_built[i]) {
570  all_tables_are_built = false;
571  }
572  }
573 
574  if(1 < verbose) {
575  G4cout << "### G4LossTableManager::LocalPhysicsTable end"
576  << G4endl;
577  }
578  if(all_tables_are_built) {
579  if(1 < verbose) {
580  G4cout << "%%%%% All dEdx and Range tables for worker are ready for run "
581  << run << " %%%%%" << G4endl;
582  }
583  }
584 }
G4PhysicsTable * RangeTableForLoss() const
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
G4PhysicsTable * DEDXTable() const
const G4ParticleDefinition * BaseParticle() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4PhysicsTable * InverseRangeTable() const
const G4ParticleDefinition * Particle() const
#define G4endl
Definition: G4ios.hh:61
G4bool IsIonisationProcess() const

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 400 of file G4LossTableManager.cc.

403 {
404  if (1 < verbose) {
405  G4cout << "G4LossTableManager::PreparePhysicsTable for "
406  << particle->GetParticleName()
407  << " and " << p->GetProcessName() << " run= " << run
408  << " loss_vector " << loss_vector.size() << G4endl;
409  }
410 
411  isMaster = theMaster;
412 
413  if(!startInitialisation) {
414  ResetParameters();
415  if (1 < verbose) {
416  G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
417  << G4endl;
418  }
419  }
420 
421  // start initialisation for the first run
422  if( -1 == run ) {
423  if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
424 
425  // initialise particles for given process
426  for (G4int j=0; j<n_loss; ++j) {
427  if (p == loss_vector[j] && !part_vector[j]) {
428  part_vector[j] = particle;
429  if(particle->GetParticleName() == "GenericIon") {
430  theGenericIon = particle;
431  }
432  }
433  }
434  }
435  startInitialisation = true;
436 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void PrepareModels(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 441 of file G4LossTableManager.cc.

443 {
444  if (1 < verbose) {
445  G4cout << "G4LossTableManager::PreparePhysicsTable for "
446  << particle->GetParticleName()
447  << " and " << p->GetProcessName() << G4endl;
448  }
449  isMaster = theMaster;
450 
451  if(!startInitialisation) {
452  ResetParameters();
453  if (1 < verbose) {
454  G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
455  << G4endl;
456  }
457  }
458 
459  // start initialisation for the first run
460  if( -1 == run ) {
461  if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
462  }
463  startInitialisation = true;
464 }
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void PrepareModels(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

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

Definition at line 469 of file G4LossTableManager.cc.

472 {
473  if (1 < verbose) {
474  G4cout << "G4LossTableManager::PreparePhysicsTable for "
475  << particle->GetParticleName()
476  << " and " << p->GetProcessName() << G4endl;
477  }
478 
479  isMaster = theMaster;
480 
481  if(!startInitialisation) {
482  ResetParameters();
483  if (1 < verbose) {
484  G4cout << "====== G4LossTableManager::PreparePhysicsTable start ====="
485  << G4endl;
486  }
487  }
488 
489  // start initialisation for the first run
490  if( -1 == run ) {
491  if(emConfigurator) { emConfigurator->PrepareModels(particle, p); }
492  }
493  startInitialisation = true;
494 }
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void PrepareModels(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4LossTableManager::Register ( G4VEnergyLossProcess p)

Definition at line 220 of file G4LossTableManager.cc.

221 {
222  if(!p) { return; }
223  for (G4int i=0; i<n_loss; ++i) {
224  if(loss_vector[i] == p) { return; }
225  }
226  if(verbose > 1) {
227  G4cout << "G4LossTableManager::Register G4VEnergyLossProcess : "
228  << p->GetProcessName() << " idx= " << n_loss << G4endl;
229  }
230  ++n_loss;
231  loss_vector.push_back(p);
232  part_vector.push_back(nullptr);
233  base_part_vector.push_back(nullptr);
234  dedx_vector.push_back(nullptr);
235  range_vector.push_back(nullptr);
236  inv_range_vector.push_back(nullptr);
237  tables_are_built.push_back(false);
238  isActive.push_back(true);
239  all_tables_are_built = false;
240 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4LossTableManager::Register ( G4VMultipleScattering p)

Definition at line 273 of file G4LossTableManager.cc.

274 {
275  if(!p) { return; }
276  G4int n = msc_vector.size();
277  for (G4int i=0; i<n; ++i) {
278  if(msc_vector[i] == p) { return; }
279  }
280  if(verbose > 1) {
281  G4cout << "G4LossTableManager::Register G4VMultipleScattering : "
282  << p->GetProcessName() << " idx= " << msc_vector.size() << G4endl;
283  }
284  msc_vector.push_back(p);
285 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4LossTableManager::Register ( G4VEmProcess p)

Definition at line 300 of file G4LossTableManager.cc.

301 {
302  if(!p) { return; }
303  G4int n = emp_vector.size();
304  for (G4int i=0; i<n; ++i) {
305  if(emp_vector[i] == p) { return; }
306  }
307  if(verbose > 1) {
308  G4cout << "G4LossTableManager::Register G4VEmProcess : "
309  << p->GetProcessName() << " idx= " << emp_vector.size() << G4endl;
310  }
311  emp_vector.push_back(p);
312 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4LossTableManager::Register ( G4VEmModel p)

Definition at line 327 of file G4LossTableManager.cc.

328 {
329  mod_vector.push_back(p);
330  if(verbose > 1) {
331  G4cout << "G4LossTableManager::Register G4VEmModel : "
332  << p->GetName() << " " << p << " " << mod_vector.size() << G4endl;
333  }
334 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
const G4String & GetName() const
Definition: G4VEmModel.hh:794

Here is the call graph for this function:

void G4LossTableManager::Register ( G4VEmFluctuationModel p)

Definition at line 352 of file G4LossTableManager.cc.

353 {
354  fmod_vector.push_back(p);
355  if(verbose > 1) {
356  G4cout << "G4LossTableManager::Register G4VEmFluctuationModel : "
357  << p->GetName() << " " << fmod_vector.size() << G4endl;
358  }
359 }
const G4String & GetName() const
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

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

Definition at line 373 of file G4LossTableManager.cc.

376 {
377  if(!p || !part) { return; }
378  for (G4int i=0; i<n_loss; ++i) {
379  if(loss_vector[i] == p) { return; }
380  }
381  if(verbose > 1) {
382  G4cout << "G4LossTableManager::RegisterExtraParticle "
383  << part->GetParticleName() << " G4VEnergyLossProcess : "
384  << p->GetProcessName() << " idx= " << n_loss << G4endl;
385  }
386  ++n_loss;
387  loss_vector.push_back(p);
388  part_vector.push_back(part);
389  base_part_vector.push_back(p->BaseParticle());
390  dedx_vector.push_back(nullptr);
391  range_vector.push_back(nullptr);
392  inv_range_vector.push_back(nullptr);
393  tables_are_built.push_back(false);
394  all_tables_are_built = false;
395 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * BaseParticle() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4LossTableManager::SetAtomDeexcitation ( G4VAtomDeexcitation p)

Definition at line 990 of file G4LossTableManager.cc.

991 {
992  if(atomDeexcitation != p) {
993  delete atomDeexcitation;
994  atomDeexcitation = p;
995  }
996 }
const char * p
Definition: xmltok.h:285

Here is the caller graph for this function:

void G4LossTableManager::SetSubCutProducer ( G4VSubCutProducer p)

Definition at line 1000 of file G4LossTableManager.cc.

1001 {
1002  if(subcutProducer != p) {
1003  delete subcutProducer;
1004  subcutProducer = p;
1005  }
1006 }
const char * p
Definition: xmltok.h:285
void G4LossTableManager::SetVerbose ( G4int  val)

Definition at line 935 of file G4LossTableManager.cc.

936 {
937  verbose = val;
938 }
G4VSubCutProducer * G4LossTableManager::SubCutProducer ( )
inline

Definition at line 451 of file G4LossTableManager.hh.

452 {
453  return subcutProducer;
454 }

Here is the caller graph for this function:

Friends And Related Function Documentation

Definition at line 103 of file G4LossTableManager.hh.


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