Geant4  10.02.p03
G4EmModelManager Class Reference

#include <G4EmModelManager.hh>

Collaboration diagram for G4EmModelManager:

Public Member Functions

 G4EmModelManager ()
 
 ~G4EmModelManager ()
 
void Clear ()
 
const G4DataVectorInitialise (const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
 
void FillDEDXVector (G4PhysicsVector *, const G4MaterialCutsCouple *, G4EmTableType t=fRestricted)
 
void FillLambdaVector (G4PhysicsVector *, const G4MaterialCutsCouple *, G4bool startFromNull=true, G4EmTableType t=fRestricted)
 
G4VEmModelGetModel (G4int, G4bool ver=false)
 
void AddEmModel (G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
 
void UpdateEmModel (const G4String &, G4double, G4double)
 
void DumpModelList (G4int verb)
 
G4VEmModelSelectModel (G4double &energy, size_t &index)
 
const G4DataVectorCuts () const
 
const G4DataVectorSubCutoff () const
 
void SetFluoFlag (G4bool val)
 
G4int NumberOfModels () const
 

Private Member Functions

G4double ComputeDEDX (G4VEmModel *model, const G4MaterialCutsCouple *, G4double kinEnergy, G4double cutEnergy, G4double minEnergy)
 
 G4EmModelManager (G4EmModelManager &)
 
G4EmModelManageroperator= (const G4EmModelManager &right)
 

Private Attributes

const G4DataVectortheCuts
 
G4DataVectortheCutsNew
 
G4DataVectortheSubCuts
 
std::vector< G4VEmModel * > models
 
std::vector< G4VEmFluctuationModel * > flucModels
 
std::vector< const G4Region * > regions
 
std::vector< G4intorderOfModels
 
std::vector< G4intisUsed
 
G4int nEmModels
 
G4int nRegions
 
std::vector< G4intidxOfRegionModels
 
std::vector< G4RegionModels * > setOfRegionModels
 
G4double maxSubCutInRange
 
const G4ParticleDefinitionparticle
 
G4int verboseLevel
 
G4bool severalModels
 
G4bool fluoFlag
 
G4RegionModelscurrRegionModel
 
G4VEmModelcurrModel
 

Detailed Description

Definition at line 143 of file G4EmModelManager.hh.

Constructor & Destructor Documentation

◆ G4EmModelManager() [1/2]

G4EmModelManager::G4EmModelManager ( )

Definition at line 123 of file G4EmModelManager.cc.

123  :
124  nEmModels(0),
125  nRegions(0),
126  particle(0),
127  verboseLevel(0)
128 {
129  maxSubCutInRange = 0.7*mm;
130  models.reserve(4);
131  flucModels.reserve(4);
132  regions.reserve(4);
133  orderOfModels.reserve(4);
134  isUsed.reserve(4);
135  severalModels = true;
136  fluoFlag = false;
137  currRegionModel = nullptr;
138  currModel = nullptr;
139  theCuts = nullptr;
140  theCutsNew = nullptr;
141  theSubCuts = nullptr;
142 }
const G4DataVector * theCuts
std::vector< G4int > orderOfModels
std::vector< G4int > isUsed
const G4ParticleDefinition * particle
G4DataVector * theSubCuts
std::vector< G4VEmModel * > models
std::vector< const G4Region * > regions
G4RegionModels * currRegionModel
G4VEmModel * currModel
std::vector< G4VEmFluctuationModel * > flucModels
static const double mm
Definition: G4SIunits.hh:114
G4DataVector * theCutsNew

◆ ~G4EmModelManager()

G4EmModelManager::~G4EmModelManager ( )

Definition at line 146 of file G4EmModelManager.cc.

147 {
148  verboseLevel = 0; // no verbosity at destruction
149  Clear();
150  delete theCutsNew;
151  delete theSubCuts;
152 }
G4DataVector * theSubCuts
G4DataVector * theCutsNew
Here is the call graph for this function:

◆ G4EmModelManager() [2/2]

G4EmModelManager::G4EmModelManager ( G4EmModelManager )
private

Member Function Documentation

◆ AddEmModel()

void G4EmModelManager::AddEmModel ( G4int  num,
G4VEmModel p,
G4VEmFluctuationModel fm,
const G4Region r 
)

Definition at line 172 of file G4EmModelManager.cc.

174 {
175  if(!p) {
176  G4cout << "G4EmModelManager::AddEmModel WARNING: no model defined."
177  << G4endl;
178  return;
179  }
180  models.push_back(p);
181  flucModels.push_back(fm);
182  regions.push_back(r);
183  orderOfModels.push_back(num);
184  isUsed.push_back(0);
185  p->DefineForRegion(r);
186  ++nEmModels;
187 }
virtual void DefineForRegion(const G4Region *)
Definition: G4VEmModel.cc:332
std::vector< G4int > orderOfModels
G4GLOB_DLL std::ostream G4cout
std::vector< G4int > isUsed
std::vector< G4VEmModel * > models
std::vector< const G4Region * > regions
#define G4endl
Definition: G4ios.hh:61
std::vector< G4VEmFluctuationModel * > flucModels
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Clear()

void G4EmModelManager::Clear ( )

Definition at line 156 of file G4EmModelManager.cc.

157 {
158  if(1 < verboseLevel) {
159  G4cout << "G4EmModelManager::Clear()" << G4endl;
160  }
161  size_t n = setOfRegionModels.size();
162  if(n > 0) {
163  for(size_t i=0; i<n; ++i) {
164  delete setOfRegionModels[i];
165  setOfRegionModels[i] = 0;
166  }
167  }
168 }
Char_t n[5]
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
std::vector< G4RegionModels * > setOfRegionModels
Here is the caller graph for this function:

◆ ComputeDEDX()

G4double G4EmModelManager::ComputeDEDX ( G4VEmModel model,
const G4MaterialCutsCouple couple,
G4double  kinEnergy,
G4double  cutEnergy,
G4double  minEnergy 
)
inlineprivate

Definition at line 275 of file G4EmModelManager.hh.

280 {
281  G4double dedx = 0.0;
282  if(model && cut > emin) {
283  dedx = model->ComputeDEDX(couple,particle,e,cut);
284  if(emin > 0.0) {dedx -= model->ComputeDEDX(couple,particle,e,emin);}
285  }
286  return dedx;
287 }
const G4ParticleDefinition * particle
G4double ComputeDEDX(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.hh:490
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Cuts()

const G4DataVector * G4EmModelManager::Cuts ( ) const
inline

Definition at line 246 of file G4EmModelManager.hh.

247 {
248  return theCuts;
249 }
const G4DataVector * theCuts

◆ DumpModelList()

void G4EmModelManager::DumpModelList ( G4int  verb)

Definition at line 720 of file G4EmModelManager.cc.

721 {
722  if(verb == 0) { return; }
723  for(G4int i=0; i<nRegions; ++i) {
725  const G4Region* reg = r->Region();
726  G4int n = r->NumberOfModels();
727  if(n > 0) {
728  G4cout << " ===== EM models for the G4Region " << reg->GetName()
729  << " ======" << G4endl;;
730  for(G4int j=0; j<n; ++j) {
731  G4VEmModel* model = models[r->ModelIndex(j)];
732  G4double emin =
734  G4double emax =
736  G4cout << std::setw(20);
737  G4cout << model->GetName() << " : Emin= "
738  << std::setw(8) << G4BestUnit(emin,"Energy")
739  << " Emax= "
740  << std::setw(8) << G4BestUnit(emax,"Energy");
741  G4PhysicsTable* table = model->GetCrossSectionTable();
742  if(table) {
743  size_t kk = table->size();
744  for(size_t k=0; k<kk; ++k) {
745  G4PhysicsVector* v = (*table)[k];
746  if(v) {
747  G4int nn = v->GetVectorLength() - 1;
748  G4cout << " Table with " << nn << " bins Emin= "
749  << std::setw(6) << G4BestUnit(v->Energy(0),"Energy")
750  << " Emax= "
751  << std::setw(6) << G4BestUnit(v->Energy(nn),"Energy");
752  break;
753  }
754  }
755  }
757  if(an) { G4cout << " " << an->GetName(); }
758  if(fluoFlag && model->DeexcitationFlag()) {
759  G4cout << " FluoActive";
760  }
761  G4cout << G4endl;
762  }
763  }
764  if(1 == nEmModels) { break; }
765  }
766  if(theCutsNew) {
767  G4cout << " ===== Limit on energy threshold has been applied "
768  << G4endl;
769  }
770 }
const G4String & GetName() const
const G4String & GetName() const
Definition: G4VEmModel.hh:795
G4double LowEdgeEnergy(G4int n) const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:617
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:826
G4int ModelIndex(G4int n) const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
int G4int
Definition: G4Types.hh:78
G4int NumberOfModels() const
static const G4double reg
Char_t n[5]
G4GLOB_DLL std::ostream G4cout
const G4Region * Region() const
G4bool DeexcitationFlag() const
Definition: G4VEmModel.hh:683
G4double LowEnergyActivationLimit() const
Definition: G4VEmModel.hh:655
size_t GetVectorLength() const
const G4String & GetName() const
static const G4double emax
std::vector< G4VEmModel * > models
#define G4endl
Definition: G4ios.hh:61
G4double Energy(size_t index) const
double G4double
Definition: G4Types.hh:76
G4double HighEnergyActivationLimit() const
Definition: G4VEmModel.hh:648
G4DataVector * theCutsNew
std::vector< G4RegionModels * > setOfRegionModels
Here is the call graph for this function:
Here is the caller graph for this function:

◆ FillDEDXVector()

void G4EmModelManager::FillDEDXVector ( G4PhysicsVector aVector,
const G4MaterialCutsCouple couple,
G4EmTableType  t = fRestricted 
)

Definition at line 563 of file G4EmModelManager.cc.

566 {
567  size_t i = couple->GetIndex();
568  G4double cut = (*theCuts)[i];
569  G4double emin = 0.0;
570 
571  if(fTotal == tType) { cut = DBL_MAX; }
572  else if(fSubRestricted == tType) {
573  emin = cut;
574  if(theSubCuts) { emin = (*theSubCuts)[i]; }
575  }
576 
577  if(1 < verboseLevel) {
578  G4cout << "G4EmModelManager::FillDEDXVector() for "
579  << couple->GetMaterial()->GetName()
580  << " cut(MeV)= " << cut
581  << " emin(MeV)= " << emin
582  << " Type " << tType
583  << " for " << particle->GetParticleName()
584  << G4endl;
585  }
586 
587  G4int reg = 0;
588  if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; }
589  const G4RegionModels* regModels = setOfRegionModels[reg];
590  G4int nmod = regModels->NumberOfModels();
591 
592  // Calculate energy losses vector
593 
594  //G4cout << "nmod= " << nmod << G4endl;
595  size_t totBinsLoss = aVector->GetVectorLength();
596  G4double del = 0.0;
597  G4int k0 = 0;
598 
599  for(size_t j=0; j<totBinsLoss; ++j) {
600 
601  G4double e = aVector->Energy(j);
602 
603  // Choose a model of energy losses
604  G4int k = 0;
605  if (nmod > 1) {
606  k = nmod;
607  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
608  do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k));
609  //G4cout << "k= " << k << G4endl;
610  if(k > 0 && k != k0) {
611  k0 = k;
612  G4double elow = regModels->LowEdgeEnergy(k);
613  G4double dedx1 = ComputeDEDX(models[regModels->ModelIndex(k-1)],
614  couple,elow,cut,emin);
615  G4double dedx2 = ComputeDEDX(models[regModels->ModelIndex(k)],
616  couple,elow,cut,emin);
617  del = 0.0;
618  if(dedx2 > 0.0) { del = (dedx1/dedx2 - 1.0)*elow; }
619  //G4cout << "elow= " << elow
620  // << " dedx1= " << dedx1 << " dedx2= " << dedx2 << G4endl;
621  }
622  }
623  G4double dedx =
624  ComputeDEDX(models[regModels->ModelIndex(k)],couple,e,cut,emin);
625  dedx *= (1.0 + del/e);
626 
627  if(2 < verboseLevel) {
628  G4cout << "Material= " << couple->GetMaterial()->GetName()
629  << " E(MeV)= " << e/MeV
630  << " dEdx(MeV/mm)= " << dedx*mm/MeV
631  << " del= " << del*mm/MeV<< " k= " << k
632  << " modelIdx= " << regModels->ModelIndex(k)
633  << G4endl;
634  }
635  if(dedx < 0.0) { dedx = 0.0; }
636  aVector->PutValue(j, dedx);
637  }
638 }
static const double MeV
Definition: G4SIunits.hh:211
std::vector< G4int > idxOfRegionModels
const G4Material * GetMaterial() const
G4double LowEdgeEnergy(G4int n) const
G4int ModelIndex(G4int n) const
int G4int
Definition: G4Types.hh:78
G4int NumberOfModels() const
static const G4double reg
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * particle
void PutValue(size_t index, G4double theValue)
size_t GetVectorLength() const
G4DataVector * theSubCuts
std::vector< G4VEmModel * > models
G4double ComputeDEDX(G4VEmModel *model, const G4MaterialCutsCouple *, G4double kinEnergy, G4double cutEnergy, G4double minEnergy)
#define G4endl
Definition: G4ios.hh:61
G4double Energy(size_t index) const
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Material.hh:178
#define DBL_MAX
Definition: templates.hh:83
static const double mm
Definition: G4SIunits.hh:114
std::vector< G4RegionModels * > setOfRegionModels
Here is the call graph for this function:
Here is the caller graph for this function:

◆ FillLambdaVector()

void G4EmModelManager::FillLambdaVector ( G4PhysicsVector aVector,
const G4MaterialCutsCouple couple,
G4bool  startFromNull = true,
G4EmTableType  t = fRestricted 
)

Definition at line 642 of file G4EmModelManager.cc.

646 {
647  size_t i = couple->GetIndex();
648  G4double cut = (*theCuts)[i];
649  G4double tmax = DBL_MAX;
650  if (fSubRestricted == tType) {
651  tmax = cut;
652  if(theSubCuts) { cut = (*theSubCuts)[i]; }
653  }
654 
655  G4int reg = 0;
656  if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; }
657  const G4RegionModels* regModels = setOfRegionModels[reg];
658  G4int nmod = regModels->NumberOfModels();
659  if(1 < verboseLevel) {
660  G4cout << "G4EmModelManager::FillLambdaVector() for "
662  << " in " << couple->GetMaterial()->GetName()
663  << " Emin(MeV)= " << aVector->Energy(0)
664  << " Emax(MeV)= " << aVector->GetMaxEnergy()
665  << " cut= " << cut
666  << " Type " << tType
667  << " nmod= " << nmod
668  << " theSubCuts " << theSubCuts
669  << G4endl;
670  }
671 
672  // Calculate lambda vector
673  size_t totBinsLambda = aVector->GetVectorLength();
674  G4double del = 0.0;
675  G4int k0 = 0;
676  G4int k = 0;
677  G4VEmModel* mod = models[regModels->ModelIndex(0)];
678  for(size_t j=0; j<totBinsLambda; ++j) {
679 
680  G4double e = aVector->Energy(j);
681 
682  // Choose a model
683  if (nmod > 1) {
684  k = nmod;
685  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
686  do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k));
687  if(k > 0 && k != k0) {
688  k0 = k;
689  G4double elow = regModels->LowEdgeEnergy(k);
690  G4VEmModel* mod1 = models[regModels->ModelIndex(k-1)];
691  G4double xs1 = mod1->CrossSection(couple,particle,elow,cut,tmax);
692  mod = models[regModels->ModelIndex(k)];
693  G4double xs2 = mod->CrossSection(couple,particle,elow,cut,tmax);
694  del = 0.0;
695  if(xs2 > 0.0) { del = (xs1/xs2 - 1.0)*elow; }
696  //G4cout << "New model k=" << k << " E(MeV)= " << e/MeV
697  // << " Elow(MeV)= " << elow/MeV << " del= " << del << G4endl;
698  }
699  }
700  G4double cross = mod->CrossSection(couple,particle,e,cut,tmax);
701  cross *= (1.0 + del/e);
702  if(fIsCrossSectionPrim == tType) { cross *= e; }
703 
704  if(j==0 && startFromNull) { cross = 0.0; }
705 
706  if(2 < verboseLevel) {
707  G4cout << "FillLambdaVector: " << j << ". e(MeV)= " << e/MeV
708  << " cross(1/mm)= " << cross*mm
709  << " del= " << del*mm << " k= " << k
710  << " modelIdx= " << regModels->ModelIndex(k)
711  << G4endl;
712  }
713  if(cross < 0.0) { cross = 0.0; }
714  aVector->PutValue(j, cross);
715  }
716 }
static const double MeV
Definition: G4SIunits.hh:211
std::vector< G4int > idxOfRegionModels
const G4Material * GetMaterial() const
G4double LowEdgeEnergy(G4int n) const
G4int ModelIndex(G4int n) const
int G4int
Definition: G4Types.hh:78
G4int NumberOfModels() const
static const G4double reg
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
G4double CrossSection(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:501
const G4ParticleDefinition * particle
void PutValue(size_t index, G4double theValue)
size_t GetVectorLength() const
G4DataVector * theSubCuts
std::vector< G4VEmModel * > models
G4double GetMaxEnergy() const
#define G4endl
Definition: G4ios.hh:61
G4double Energy(size_t index) const
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Material.hh:178
#define DBL_MAX
Definition: templates.hh:83
static const double mm
Definition: G4SIunits.hh:114
std::vector< G4RegionModels * > setOfRegionModels
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetModel()

G4VEmModel * G4EmModelManager::GetModel ( G4int  i,
G4bool  ver = false 
)

Definition at line 210 of file G4EmModelManager.cc.

211 {
212  G4VEmModel* model = nullptr;
213  if(i >= 0 && i < nEmModels) { model = models[i]; }
214  else if(verboseLevel > 0 && ver) {
215  G4cout << "G4EmModelManager::GetModel WARNING: "
216  << "index " << i << " is wrong Nmodels= "
217  << nEmModels;
218  if(particle) G4cout << " for " << particle->GetParticleName();
219  G4cout<< G4endl;
220  }
221  return model;
222 }
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * particle
std::vector< G4VEmModel * > models
#define G4endl
Definition: G4ios.hh:61
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Initialise()

const G4DataVector * G4EmModelManager::Initialise ( const G4ParticleDefinition part,
const G4ParticleDefinition secPart,
G4double  minSubRange,
G4int  verb 
)

Definition at line 227 of file G4EmModelManager.cc.

231 {
232  verboseLevel = val;
233  G4String partname = p->GetParticleName();
234  // if(partname == "proton") verboseLevel = 2;
235  //else verboseLevel = 0;
236  if(1 < verboseLevel) {
237  G4cout << "G4EmModelManager::Initialise() for "
238  << partname << G4endl;
239  }
240  // Are models defined?
241  if(nEmModels < 1) {
243  ed << "No models found out for " << p->GetParticleName()
244  << " !";
245  G4Exception("G4EmModelManager::Initialise","em0002",
246  FatalException, ed);
247  }
248 
249  particle = p;
250  Clear(); // needed if run is not first
251 
252  G4RegionStore* regionStore = G4RegionStore::GetInstance();
253  const G4Region* world =
254  regionStore->GetRegion("DefaultRegionForTheWorld", false);
255 
256  // Identify the list of regions with different set of models
257  nRegions = 1;
258  std::vector<const G4Region*> setr;
259  setr.push_back(world);
260  G4bool isWorld = false;
261 
262  for (G4int ii=0; ii<nEmModels; ++ii) {
263  const G4Region* r = regions[ii];
264  if ( r == 0 || r == world) {
265  isWorld = true;
266  regions[ii] = world;
267  } else {
268  G4bool newRegion = true;
269  if (nRegions>1) {
270  for (G4int j=1; j<nRegions; ++j) {
271  if ( r == setr[j] ) { newRegion = false; }
272  }
273  }
274  if (newRegion) {
275  setr.push_back(r);
276  nRegions++;
277  }
278  }
279  }
280  // Are models defined?
281  if(!isWorld) {
283  ed << "No models defined for the World volume for "
284  << p->GetParticleName() << " !";
285  G4Exception("G4EmModelManager::Initialise","em0002",
286  FatalException, ed);
287  }
288 
289  G4ProductionCutsTable* theCoupleTable=
291  size_t numOfCouples = theCoupleTable->GetTableSize();
292 
293  // prepare vectors, shortcut for the case of only 1 model
294  // or only one region
295  if(nRegions > 1 && nEmModels > 1) {
296  idxOfRegionModels.resize(numOfCouples,0);
297  setOfRegionModels.resize((size_t)nRegions,0);
298  } else {
299  idxOfRegionModels.resize(1,0);
300  setOfRegionModels.resize(1,0);
301  }
302 
303  std::vector<G4int> modelAtRegion(nEmModels);
304  std::vector<G4int> modelOrd(nEmModels);
305  G4DataVector eLow(nEmModels+1);
306  G4DataVector eHigh(nEmModels);
307 
308  if(1 < verboseLevel) {
309  G4cout << " Nregions= " << nRegions
310  << " Nmodels= " << nEmModels << G4endl;
311  }
312 
313  // Order models for regions
314  for (G4int reg=0; reg<nRegions; ++reg) {
315  const G4Region* region = setr[reg];
316  G4int n = 0;
317 
318  for (G4int ii=0; ii<nEmModels; ++ii) {
319 
320  G4VEmModel* model = models[ii];
321  if ( region == regions[ii] ) {
322 
323  G4double tmin = model->LowEnergyLimit();
324  G4double tmax = model->HighEnergyLimit();
325  G4int ord = orderOfModels[ii];
326  G4bool push = true;
327  G4bool insert = false;
328  G4int idx = n;
329 
330  if(1 < verboseLevel) {
331  G4cout << "Model #" << ii
332  << " <" << model->GetName() << "> for region <";
333  if (region) G4cout << region->GetName();
334  G4cout << "> "
335  << " tmin(MeV)= " << tmin/MeV
336  << "; tmax(MeV)= " << tmax/MeV
337  << "; order= " << ord
338  << G4endl;
339  }
340 
341  if(n > 0) {
342 
343  // extend energy range to previous models
344  tmin = std::min(tmin, eHigh[n-1]);
345  tmax = std::max(tmax, eLow[0]);
346  //G4cout << "tmin= " << tmin << " tmax= "
347  // << tmax << " ord= " << ord <<G4endl;
348  // empty energy range
349  if( tmax - tmin <= eV) push = false;
350  // low-energy model
351  else if (tmax == eLow[0]) {
352  push = false;
353  insert = true;
354  idx = 0;
355  // resolve intersections
356  } else if(tmin < eHigh[n-1]) {
357  // compare order
358  for(G4int k=0; k<n; ++k) {
359  // new model has lower application
360  if(ord >= modelOrd[k]) {
361  if(tmin < eHigh[k] && tmin >= eLow[k]) tmin = eHigh[k];
362  if(tmax <= eHigh[k] && tmax > eLow[k]) tmax = eLow[k];
363  if(tmax > eHigh[k] && tmin < eLow[k]) {
364  if(tmax - eHigh[k] > eLow[k] - tmin) tmin = eHigh[k];
365  else tmax = eLow[k];
366  }
367  if( tmax - tmin <= eV) {
368  push = false;
369  break;
370  }
371  }
372  }
373  //G4cout << "tmin= " << tmin << " tmax= "
374  // << tmax << " push= " << push << " idx= " << idx <<G4endl;
375  if(push) {
376  if (tmax == eLow[0]) {
377  push = false;
378  insert = true;
379  idx = 0;
380  // continue resolve intersections
381  } else if(tmin < eHigh[n-1]) {
382  // last energy interval
383  if(tmin > eLow[n-1] && tmax >= eHigh[n-1]) {
384  eHigh[n-1] = tmin;
385  // first energy interval
386  } else if(tmin <= eLow[0] && tmax < eHigh[0]) {
387  eLow[0] = tmax;
388  push = false;
389  insert = true;
390  idx = 0;
391  } else {
392  // find energy interval to replace
393  for(G4int k=0; k<n; ++k) {
394  if(tmin <= eLow[k] && tmax >= eHigh[k]) {
395  push = false;
396  modelAtRegion[k] = ii;
397  modelOrd[k] = ord;
398  isUsed[ii] = 1;
399  }
400  }
401  }
402  }
403  }
404  }
405  }
406  if(insert) {
407  for(G4int k=n-1; k>=idx; --k) {
408  modelAtRegion[k+1] = modelAtRegion[k];
409  modelOrd[k+1] = modelOrd[k];
410  eLow[k+1] = eLow[k];
411  eHigh[k+1] = eHigh[k];
412  }
413  }
414  //G4cout << "push= " << push << " insert= " << insert
415  //<< " idx= " << idx <<G4endl;
416  if (push || insert) {
417  ++n;
418  modelAtRegion[idx] = ii;
419  modelOrd[idx] = ord;
420  eLow[idx] = tmin;
421  eHigh[idx] = tmax;
422  isUsed[ii] = 1;
423  }
424  }
425  }
426  eLow[0] = 0.0;
427  eLow[n] = eHigh[n-1];
428 
429  if(1 < verboseLevel) {
430  G4cout << "New G4RegionModels set with " << n << " models for region <";
431  if (region) { G4cout << region->GetName(); }
432  G4cout << "> Elow(MeV)= ";
433  for(G4int ii=0; ii<=n; ++ii) {G4cout << eLow[ii]/MeV << " ";}
434  G4cout << G4endl;
435  }
436  G4RegionModels* rm = new G4RegionModels(n, modelAtRegion, eLow, region);
437  setOfRegionModels[reg] = rm;
438  // shortcut
439  if(1 == nEmModels) { break; }
440  }
441 
443  currModel = models[0];
444 
445  // Access to materials and build cuts
446  size_t idx = 1;
447  if(secondaryParticle) {
448  if( secondaryParticle == G4Gamma::Gamma() ) { idx = 0; }
449  else if( secondaryParticle == G4Electron::Electron()) { idx = 1; }
450  else if( secondaryParticle == G4Positron::Positron()) { idx = 2; }
451  else { idx = 3; }
452  }
453 
454  theCuts =
455  static_cast<const G4DataVector*>(theCoupleTable->GetEnergyCutsVector(idx));
456 
457  // for the second run the check on cuts should be repeated
458  if(theCutsNew) { *theCutsNew = *theCuts; }
459 
460  if(minSubRange < 1.0) {
461  if( !theSubCuts ) { theSubCuts = new G4DataVector(); }
462  theSubCuts->resize(numOfCouples,DBL_MAX);
463  }
464 
465  // G4cout << "========Start define cuts" << G4endl;
466  // define cut values
467  for(size_t i=0; i<numOfCouples; ++i) {
468 
469  const G4MaterialCutsCouple* couple =
470  theCoupleTable->GetMaterialCutsCouple(i);
471  const G4Material* material = couple->GetMaterial();
472  const G4ProductionCuts* pcuts = couple->GetProductionCuts();
473 
474  G4int reg = 0;
475  if(nRegions > 1 && nEmModels > 1) {
476  reg = nRegions;
477  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
478  do {--reg;} while (reg>0 && pcuts != (setr[reg]->GetProductionCuts()));
479  idxOfRegionModels[i] = reg;
480  }
481  if(1 < verboseLevel) {
482  G4cout << "G4EmModelManager::Initialise() for "
483  << material->GetName()
484  << " indexOfCouple= " << i
485  << " indexOfRegion= " << reg
486  << G4endl;
487  }
488 
489  G4double cut = (*theCuts)[i];
490  if(secondaryParticle) {
491 
492  // compute subcut
493  if( cut < DBL_MAX && minSubRange < 1.0) {
494  G4double subcut = minSubRange*cut;
495  G4double rcut = std::min(minSubRange*pcuts->GetProductionCut(idx),
497  G4double tcutmax =
498  theCoupleTable->ConvertRangeToEnergy(secondaryParticle,
499  material,rcut);
500  if(tcutmax < subcut) { subcut = tcutmax; }
501  (*theSubCuts)[i] = subcut;
502  }
503 
504  // note that idxOfRegionModels[] not always filled
505  G4int inn = 0;
506  G4int nnm = 1;
507  if(nRegions > 1 && nEmModels > 1) {
508  inn = idxOfRegionModels[i];
509  }
510  // check cuts and introduce upper limits
511  //G4cout << "idx= " << i << " cut(keV)= " << cut/keV << G4endl;
514 
515  //G4cout << "idx= " << i << " Nmod= " << nnm << G4endl;
516 
517  for(G4int jj=0; jj<nnm; ++jj) {
518  //G4cout << "jj= " << jj << " modidx= "
519  // << currRegionModel->ModelIndex(jj) << G4endl;
521  G4double cutlim = currModel->MinEnergyCut(particle,couple);
522  if(cutlim > cut) {
523  if(!theCutsNew) { theCutsNew = new G4DataVector(*theCuts); }
524  (*theCutsNew)[i] = cutlim;
525  /*
526  G4cout << "### " << partname << " energy loss model in "
527  << material->GetName()
528  << " Cut was changed from " << cut/keV << " keV to "
529  << cutlim/keV << " keV " << " due to "
530  << currModel->GetName() << G4endl;
531  */
532  }
533  }
534  }
535  }
536  if(theCutsNew) { theCuts = theCutsNew; }
537 
538  // initialize models
539  G4int nn = 0;
540  severalModels = true;
541  for(G4int jj=0; jj<nEmModels; ++jj) {
542  if(1 == isUsed[jj]) {
543  ++nn;
544  currModel = models[jj];
546  if(flucModels[jj]) { flucModels[jj]->InitialiseMe(particle); }
547  }
548  }
549  if(1 == nn) { severalModels = false; }
550 
551  if(1 < verboseLevel) {
552  G4cout << "G4EmModelManager for " << partname
553  << " is initialised; nRegions= " << nRegions
554  << " severalModels: " << severalModels
555  << G4endl;
556  }
557 
558  return theCuts;
559 }
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:641
const G4String & GetName() const
Definition: G4VEmModel.hh:795
static const double MeV
Definition: G4SIunits.hh:211
std::vector< G4int > idxOfRegionModels
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
const G4Material * GetMaterial() const
G4ProductionCuts * GetProductionCuts() const
const G4DataVector * theCuts
G4int ModelIndex(G4int n) const
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
virtual G4double MinEnergyCut(const G4ParticleDefinition *, const G4MaterialCutsCouple *)
Definition: G4VEmModel.cc:387
int G4int
Definition: G4Types.hh:78
G4int NumberOfModels() const
std::vector< G4int > orderOfModels
static const G4double reg
G4double GetProductionCut(G4int index) const
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0
static G4RegionStore * GetInstance()
string material
Definition: eplot.py:19
Char_t n[5]
G4GLOB_DLL std::ostream G4cout
std::vector< G4int > isUsed
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:634
bool G4bool
Definition: G4Types.hh:79
const G4ParticleDefinition * particle
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
const G4String & GetName() const
G4DataVector * theSubCuts
std::vector< G4VEmModel * > models
static G4ProductionCutsTable * GetProductionCutsTable()
static const double eV
Definition: G4SIunits.hh:212
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
static G4Positron * Positron()
Definition: G4Positron.cc:94
std::vector< const G4Region * > regions
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
G4RegionModels * currRegionModel
G4double ConvertRangeToEnergy(const G4ParticleDefinition *particle, const G4Material *material, G4double range)
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
G4VEmModel * currModel
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Material.hh:178
#define DBL_MAX
Definition: templates.hh:83
std::vector< G4VEmFluctuationModel * > flucModels
G4DataVector * theCutsNew
std::vector< G4RegionModels * > setOfRegionModels
Here is the call graph for this function:
Here is the caller graph for this function:

◆ NumberOfModels()

G4int G4EmModelManager::NumberOfModels ( ) const
inline

Definition at line 267 of file G4EmModelManager.hh.

268 {
269  return nEmModels;
270 }
Here is the caller graph for this function:

◆ operator=()

G4EmModelManager& G4EmModelManager::operator= ( const G4EmModelManager right)
private

◆ SelectModel()

G4VEmModel * G4EmModelManager::SelectModel ( G4double energy,
size_t &  index 
)
inline

Definition at line 232 of file G4EmModelManager.hh.

234 {
235  if(severalModels) {
236  if(nRegions > 1) {
238  }
240  }
241  return currModel;
242 }
std::vector< G4int > idxOfRegionModels
Int_t index
G4int SelectIndex(G4double e) const
std::vector< G4VEmModel * > models
G4RegionModels * currRegionModel
G4VEmModel * currModel
std::vector< G4RegionModels * > setOfRegionModels
Here is the caller graph for this function:

◆ SetFluoFlag()

void G4EmModelManager::SetFluoFlag ( G4bool  val)
inline

Definition at line 260 of file G4EmModelManager.hh.

261 {
262  fluoFlag = val;
263 }
Here is the caller graph for this function:

◆ SubCutoff()

const G4DataVector * G4EmModelManager::SubCutoff ( ) const
inline

Definition at line 253 of file G4EmModelManager.hh.

254 {
255  return theSubCuts;
256 }
G4DataVector * theSubCuts
Here is the caller graph for this function:

◆ UpdateEmModel()

void G4EmModelManager::UpdateEmModel ( const G4String nam,
G4double  emin,
G4double  emax 
)

Definition at line 191 of file G4EmModelManager.cc.

193 {
194  if (nEmModels > 0) {
195  for(G4int i=0; i<nEmModels; ++i) {
196  if(nam == models[i]->GetName()) {
197  models[i]->SetLowEnergyLimit(emin);
198  models[i]->SetHighEnergyLimit(emax);
199  break;
200  }
201  }
202  }
203  G4cout << "G4EmModelManager::UpdateEmModel WARNING: no model <"
204  << nam << "> is found out"
205  << G4endl;
206 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
static const G4double emax
std::vector< G4VEmModel * > models
#define G4endl
Definition: G4ios.hh:61
Here is the caller graph for this function:

Member Data Documentation

◆ currModel

G4VEmModel* G4EmModelManager::currModel
private

Definition at line 226 of file G4EmModelManager.hh.

◆ currRegionModel

G4RegionModels* G4EmModelManager::currRegionModel
private

Definition at line 225 of file G4EmModelManager.hh.

◆ flucModels

std::vector<G4VEmFluctuationModel*> G4EmModelManager::flucModels
private

Definition at line 205 of file G4EmModelManager.hh.

◆ fluoFlag

G4bool G4EmModelManager::fluoFlag
private

Definition at line 222 of file G4EmModelManager.hh.

◆ idxOfRegionModels

std::vector<G4int> G4EmModelManager::idxOfRegionModels
private

Definition at line 213 of file G4EmModelManager.hh.

◆ isUsed

std::vector<G4int> G4EmModelManager::isUsed
private

Definition at line 208 of file G4EmModelManager.hh.

◆ maxSubCutInRange

G4double G4EmModelManager::maxSubCutInRange
private

Definition at line 216 of file G4EmModelManager.hh.

◆ models

std::vector<G4VEmModel*> G4EmModelManager::models
private

Definition at line 204 of file G4EmModelManager.hh.

◆ nEmModels

G4int G4EmModelManager::nEmModels
private

Definition at line 210 of file G4EmModelManager.hh.

◆ nRegions

G4int G4EmModelManager::nRegions
private

Definition at line 211 of file G4EmModelManager.hh.

◆ orderOfModels

std::vector<G4int> G4EmModelManager::orderOfModels
private

Definition at line 207 of file G4EmModelManager.hh.

◆ particle

const G4ParticleDefinition* G4EmModelManager::particle
private

Definition at line 218 of file G4EmModelManager.hh.

◆ regions

std::vector<const G4Region*> G4EmModelManager::regions
private

Definition at line 206 of file G4EmModelManager.hh.

◆ setOfRegionModels

std::vector<G4RegionModels*> G4EmModelManager::setOfRegionModels
private

Definition at line 214 of file G4EmModelManager.hh.

◆ severalModels

G4bool G4EmModelManager::severalModels
private

Definition at line 221 of file G4EmModelManager.hh.

◆ theCuts

const G4DataVector* G4EmModelManager::theCuts
private

Definition at line 200 of file G4EmModelManager.hh.

◆ theCutsNew

G4DataVector* G4EmModelManager::theCutsNew
private

Definition at line 201 of file G4EmModelManager.hh.

◆ theSubCuts

G4DataVector* G4EmModelManager::theSubCuts
private

Definition at line 202 of file G4EmModelManager.hh.

◆ verboseLevel

G4int G4EmModelManager::verboseLevel
private

Definition at line 220 of file G4EmModelManager.hh.


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