Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4EmModelManager Class Reference

#include <G4EmModelManager.hh>

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)
 
void AddEmModel (G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
 
void UpdateEmModel (const G4String &model_name, G4double emin, G4double emax)
 
G4VEmModelGetModel (G4int idx, G4bool ver=false)
 
G4VEmModelGetRegionModel (G4int idx, size_t index_couple)
 
G4int NumberOfRegionModels (size_t index_couple) const
 
void DumpModelList (std::ostream &out, G4int verb)
 
G4VEmModelSelectModel (G4double &energy, size_t &index)
 
const G4DataVectorCuts () const
 
const G4DataVectorSubCutoff () const
 
void SetFluoFlag (G4bool val)
 
G4int NumberOfModels () const
 
 G4EmModelManager (G4EmModelManager &)=delete
 
G4EmModelManageroperator= (const G4EmModelManager &right)=delete
 

Detailed Description

Definition at line 141 of file G4EmModelManager.hh.

Constructor & Destructor Documentation

◆ G4EmModelManager() [1/2]

G4EmModelManager::G4EmModelManager ( )

Definition at line 98 of file G4EmModelManager.cc.

98 :
99 nEmModels(0),
100 nRegions(0),
101 particle(0),
102 verboseLevel(0)
103{
104 maxSubCutInRange = 0.7*mm;
105 models.reserve(4);
106 flucModels.reserve(4);
107 regions.reserve(4);
108 orderOfModels.reserve(4);
109 isUsed.reserve(4);
110 severalModels = true;
111 fluoFlag = false;
112 currRegionModel = nullptr;
113 currModel = nullptr;
114 theCuts = nullptr;
115 theCutsNew = nullptr;
116 theSubCuts = nullptr;
117}

◆ ~G4EmModelManager()

G4EmModelManager::~G4EmModelManager ( )

Definition at line 121 of file G4EmModelManager.cc.

122{
123 verboseLevel = 0; // no verbosity at destruction
124 Clear();
125 delete theCutsNew;
126 delete theSubCuts;
127}

◆ G4EmModelManager() [2/2]

G4EmModelManager::G4EmModelManager ( G4EmModelManager )
delete

Member Function Documentation

◆ AddEmModel()

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

Definition at line 147 of file G4EmModelManager.cc.

149{
150 if(!p) {
151 G4cout << "G4EmModelManager::AddEmModel WARNING: no model defined."
152 << G4endl;
153 return;
154 }
155 models.push_back(p);
156 flucModels.push_back(fm);
157 regions.push_back(r);
158 orderOfModels.push_back(num);
159 isUsed.push_back(0);
160 p->DefineForRegion(r);
161 ++nEmModels;
162}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
virtual void DefineForRegion(const G4Region *)
Definition: G4VEmModel.cc:378

Referenced by G4VMultipleScattering::AddEmModel(), G4VEmProcess::AddEmModel(), G4VEnergyLossProcess::AddEmModel(), and G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel().

◆ Clear()

void G4EmModelManager::Clear ( )

Definition at line 131 of file G4EmModelManager.cc.

132{
133 if(1 < verboseLevel) {
134 G4cout << "G4EmModelManager::Clear()" << G4endl;
135 }
136 size_t n = setOfRegionModels.size();
137 if(n > 0) {
138 for(size_t i=0; i<n; ++i) {
139 delete setOfRegionModels[i];
140 setOfRegionModels[i] = nullptr;
141 }
142 }
143}

Referenced by Initialise(), and ~G4EmModelManager().

◆ Cuts()

const G4DataVector * G4EmModelManager::Cuts ( ) const
inline

Definition at line 254 of file G4EmModelManager.hh.

255{
256 return theCuts;
257}

◆ DumpModelList()

void G4EmModelManager::DumpModelList ( std::ostream &  out,
G4int  verb 
)

Definition at line 773 of file G4EmModelManager.cc.

774{
775 if(verb == 0) { return; }
776 for(G4int i=0; i<nRegions; ++i) {
777 G4RegionModels* r = setOfRegionModels[i];
778 const G4Region* reg = r->Region();
779 G4int n = r->NumberOfModels();
780 if(n > 0) {
781 out << " ===== EM models for the G4Region " << reg->GetName()
782 << " ======" << G4endl;
783 for(G4int j=0; j<n; ++j) {
784 G4VEmModel* model = models[r->ModelIndex(j)];
785 G4double emin =
786 std::max(r->LowEdgeEnergy(j),model->LowEnergyActivationLimit());
787 G4double emax =
788 std::min(r->LowEdgeEnergy(j+1),model->HighEnergyActivationLimit());
789 if(emax > emin) {
790 out << std::setw(20);
791 out << model->GetName() << " : Emin="
792 << std::setw(5) << G4BestUnit(emin,"Energy")
793 << " Emax="
794 << std::setw(5) << G4BestUnit(emax,"Energy");
795 G4PhysicsTable* table = model->GetCrossSectionTable();
796 if(table) {
797 size_t kk = table->size();
798 for(size_t k=0; k<kk; ++k) {
799 const G4PhysicsVector* v = (*table)[k];
800 if(v) {
801 G4int nn = v->GetVectorLength() - 1;
802 out << " Nbins=" << nn << " "
803 << std::setw(3) << G4BestUnit(v->Energy(0),"Energy")
804 << " - "
805 << std::setw(3) << G4BestUnit(v->Energy(nn),"Energy");
806 break;
807 }
808 }
809 }
811 if(an) { out << " " << an->GetName(); }
812 if(fluoFlag && model->DeexcitationFlag()) {
813 out << " Fluo";
814 }
815 out << G4endl;
816 G4VMscModel* msc = dynamic_cast<G4VMscModel*>(model);
817 if(msc != nullptr) msc->DumpParameters(out);
818 }
819 }
820 }
821 if(1 == nEmModels) { break; }
822 }
823 if(theCutsNew) {
824 out << " ===== Limit on energy threshold has been applied " << G4endl;
825 }
826}
#define G4BestUnit(a, b)
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
G4double Energy(std::size_t index) const
std::size_t GetVectorLength() const
const G4String & GetName() const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:611
G4double HighEnergyActivationLimit() const
Definition: G4VEmModel.hh:659
G4bool DeexcitationFlag() const
Definition: G4VEmModel.hh:694
const G4String & GetName() const
Definition: G4VEmModel.hh:827
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:860
G4double LowEnergyActivationLimit() const
Definition: G4VEmModel.hh:666
void DumpParameters(std::ostream &out) const
Definition: G4VMscModel.cc:160

◆ FillDEDXVector()

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

Definition at line 616 of file G4EmModelManager.cc.

619{
620 size_t i = couple->GetIndex();
621 G4double cut = (*theCuts)[i];
622 G4double emin = 0.0;
623
624 if(fTotal == tType) { cut = DBL_MAX; }
625 else if(fSubRestricted == tType) {
626 emin = cut;
627 if(theSubCuts) { emin = (*theSubCuts)[i]; }
628 }
629
630 if(1 < verboseLevel) {
631 G4cout << "G4EmModelManager::FillDEDXVector() for "
632 << couple->GetMaterial()->GetName()
633 << " cut(MeV)= " << cut
634 << " emin(MeV)= " << emin
635 << " Type " << tType
636 << " for " << particle->GetParticleName()
637 << G4endl;
638 }
639
640 G4int reg = 0;
641 if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; }
642 const G4RegionModels* regModels = setOfRegionModels[reg];
643 G4int nmod = regModels->NumberOfModels();
644
645 // Calculate energy losses vector
646
647 //G4cout << "nmod= " << nmod << G4endl;
648 size_t totBinsLoss = aVector->GetVectorLength();
649 G4double del = 0.0;
650 G4int k0 = 0;
651
652 for(size_t j=0; j<totBinsLoss; ++j) {
653
654 G4double e = aVector->Energy(j);
655
656 // Choose a model of energy losses
657 G4int k = 0;
658 if (nmod > 1) {
659 k = nmod;
660 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
661 do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k));
662 //G4cout << "k= " << k << G4endl;
663 if(k > 0 && k != k0) {
664 k0 = k;
665 G4double elow = regModels->LowEdgeEnergy(k);
666 G4double dedx1 = ComputeDEDX(models[regModels->ModelIndex(k-1)],
667 couple,elow,cut,emin);
668 G4double dedx2 = ComputeDEDX(models[regModels->ModelIndex(k)],
669 couple,elow,cut,emin);
670 del = 0.0;
671 if(dedx2 > 0.0) { del = (dedx1/dedx2 - 1.0)*elow; }
672 //G4cout << "elow= " << elow
673 // << " dedx1= " << dedx1 << " dedx2= " << dedx2 << G4endl;
674 }
675 }
676 G4double dedx =
677 ComputeDEDX(models[regModels->ModelIndex(k)],couple,e,cut,emin);
678 dedx *= (1.0 + del/e);
679
680 if(2 < verboseLevel) {
681 G4cout << "Material= " << couple->GetMaterial()->GetName()
682 << " E(MeV)= " << e/MeV
683 << " dEdx(MeV/mm)= " << dedx*mm/MeV
684 << " del= " << del*mm/MeV<< " k= " << k
685 << " modelIdx= " << regModels->ModelIndex(k)
686 << G4endl;
687 }
688 if(dedx < 0.0) { dedx = 0.0; }
689 aVector->PutValue(j, dedx);
690 }
691}
@ fSubRestricted
@ fTotal
const G4Material * GetMaterial() const
const G4String & GetName() const
Definition: G4Material.hh:175
const G4String & GetParticleName() const
void PutValue(std::size_t index, G4double theValue)
#define DBL_MAX
Definition: templates.hh:62

Referenced by G4VEnergyLossProcess::BuildDEDXTable().

◆ FillLambdaVector()

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

Definition at line 695 of file G4EmModelManager.cc.

699{
700 size_t i = couple->GetIndex();
701 G4double cut = (*theCuts)[i];
702 G4double tmax = DBL_MAX;
703 if (fSubRestricted == tType) {
704 tmax = cut;
705 if(theSubCuts) { cut = (*theSubCuts)[i]; }
706 }
707
708 G4int reg = 0;
709 if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; }
710 const G4RegionModels* regModels = setOfRegionModels[reg];
711 G4int nmod = regModels->NumberOfModels();
712 if(1 < verboseLevel) {
713 G4cout << "G4EmModelManager::FillLambdaVector() for "
714 << particle->GetParticleName()
715 << " in " << couple->GetMaterial()->GetName()
716 << " Emin(MeV)= " << aVector->Energy(0)
717 << " Emax(MeV)= " << aVector->GetMaxEnergy()
718 << " cut= " << cut
719 << " Type " << tType
720 << " nmod= " << nmod
721 << " theSubCuts " << theSubCuts
722 << G4endl;
723 }
724
725 // Calculate lambda vector
726 size_t totBinsLambda = aVector->GetVectorLength();
727 G4double del = 0.0;
728 G4int k0 = 0;
729 G4int k = 0;
730 G4VEmModel* mod = models[regModels->ModelIndex(0)];
731 for(size_t j=0; j<totBinsLambda; ++j) {
732
733 G4double e = aVector->Energy(j);
734
735 // Choose a model
736 if (nmod > 1) {
737 k = nmod;
738 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
739 do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k));
740 if(k > 0 && k != k0) {
741 k0 = k;
742 G4double elow = regModels->LowEdgeEnergy(k);
743 G4VEmModel* mod1 = models[regModels->ModelIndex(k-1)];
744 G4double xs1 = mod1->CrossSection(couple,particle,elow,cut,tmax);
745 mod = models[regModels->ModelIndex(k)];
746 G4double xs2 = mod->CrossSection(couple,particle,elow,cut,tmax);
747 del = 0.0;
748 if(xs2 > 0.0) { del = (xs1/xs2 - 1.0)*elow; }
749 //G4cout << "New model k=" << k << " E(MeV)= " << e/MeV
750 // << " Elow(MeV)= " << elow/MeV << " del= " << del << G4endl;
751 }
752 }
753 G4double cross = mod->CrossSection(couple,particle,e,cut,tmax);
754 cross *= (1.0 + del/e);
755 if(fIsCrossSectionPrim == tType) { cross *= e; }
756
757 if(j==0 && startFromNull) { cross = 0.0; }
758
759 if(2 < verboseLevel) {
760 G4cout << "FillLambdaVector: " << j << ". e(MeV)= " << e/MeV
761 << " cross(1/mm)= " << cross*mm
762 << " del= " << del*mm << " k= " << k
763 << " modelIdx= " << regModels->ModelIndex(k)
764 << G4endl;
765 }
766 cross = std::max(cross, 0.0);
767 aVector->PutValue(j, cross);
768 }
769}
@ fIsCrossSectionPrim
G4double GetMaxEnergy() const
G4double CrossSection(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:529

Referenced by G4VEnergyLossProcess::BuildLambdaTable().

◆ GetModel()

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

Definition at line 185 of file G4EmModelManager.cc.

186{
187 G4VEmModel* model = nullptr;
188 if(i < nEmModels) { model = models[i]; }
189 else if(verboseLevel > 0 && ver) {
190 G4cout << "G4EmModelManager::GetModel WARNING: "
191 << "index " << i << " is wrong Nmodels= "
192 << nEmModels;
193 if(particle) { G4cout << " for " << particle->GetParticleName(); }
194 G4cout<< G4endl;
195 }
196 return model;
197}

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

◆ GetRegionModel()

G4VEmModel * G4EmModelManager::GetRegionModel ( G4int  idx,
size_t  index_couple 
)

Definition at line 201 of file G4EmModelManager.cc.

202{
203 G4RegionModels* rm = setOfRegionModels[idxOfRegionModels[idx]];
204 G4VEmModel* mod = models[rm->ModelIndex(k)];
205 return mod;
206}

Referenced by G4VEmProcess::GetRegionModel().

◆ Initialise()

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

Definition at line 219 of file G4EmModelManager.cc.

223{
224 verboseLevel = val;
225 G4String partname = p->GetParticleName();
226 if(1 < verboseLevel) {
227 G4cout << "G4EmModelManager::Initialise() for "
228 << partname << " Nmodels= " << nEmModels << G4endl;
229 }
230 // Are models defined?
231 if(nEmModels < 1) {
233 ed << "No models found out for " << p->GetParticleName()
234 << " !";
235 G4Exception("G4EmModelManager::Initialise","em0002",
236 FatalException, ed);
237 }
238
239 particle = p;
240 Clear(); // needed if run is not first
241
243 const G4Region* world =
244 regionStore->GetRegion("DefaultRegionForTheWorld", false);
245
246 // Identify the list of regions with different set of models
247 nRegions = 1;
248 std::vector<const G4Region*> setr;
249 setr.push_back(world);
250 G4bool isWorld = false;
251
252 for (G4int ii=0; ii<nEmModels; ++ii) {
253 const G4Region* r = regions[ii];
254 if ( r == 0 || r == world) {
255 isWorld = true;
256 regions[ii] = world;
257 } else {
258 G4bool newRegion = true;
259 if (nRegions>1) {
260 for (G4int j=1; j<nRegions; ++j) {
261 if ( r == setr[j] ) { newRegion = false; }
262 }
263 }
264 if (newRegion) {
265 setr.push_back(r);
266 nRegions++;
267 }
268 }
269 }
270 // Are models defined?
271 if(!isWorld) {
273 ed << "No models defined for the World volume for "
274 << p->GetParticleName() << " !";
275 G4Exception("G4EmModelManager::Initialise","em0002",
276 FatalException, ed);
277 }
278
279 G4ProductionCutsTable* theCoupleTable=
281 size_t numOfCouples = theCoupleTable->GetTableSize();
282
283 // prepare vectors, shortcut for the case of only 1 model
284 // or only one region
285 if(nRegions > 1 && nEmModels > 1) {
286 idxOfRegionModels.resize(numOfCouples,0);
287 setOfRegionModels.resize((size_t)nRegions,0);
288 } else {
289 idxOfRegionModels.resize(1,0);
290 setOfRegionModels.resize(1,0);
291 }
292
293 std::vector<G4int> modelAtRegion(nEmModels);
294 std::vector<G4int> modelOrd(nEmModels);
295 G4DataVector eLow(nEmModels+1);
296 G4DataVector eHigh(nEmModels);
297
298 if(1 < verboseLevel) {
299 G4cout << " Nregions= " << nRegions
300 << " Nmodels= " << nEmModels << G4endl;
301 }
302
303 // Order models for regions
304 for (G4int reg=0; reg<nRegions; ++reg) {
305 const G4Region* region = setr[reg];
306 G4int n = 0;
307
308 for (G4int ii=0; ii<nEmModels; ++ii) {
309
310 G4VEmModel* model = models[ii];
311 if ( region == regions[ii] ) {
312
313 G4double tmin = model->LowEnergyLimit();
314 G4double tmax = model->HighEnergyLimit();
315 G4int ord = orderOfModels[ii];
316 G4bool push = true;
317 G4bool insert = false;
318 G4int idx = n;
319
320 if(1 < verboseLevel) {
321 G4cout << "Model #" << ii
322 << " <" << model->GetName() << "> for region <";
323 if (region) G4cout << region->GetName();
324 G4cout << "> "
325 << " tmin(MeV)= " << tmin/MeV
326 << "; tmax(MeV)= " << tmax/MeV
327 << "; order= " << ord
328 << "; tminAct= " << model->LowEnergyActivationLimit()/MeV
329 << "; tmaxAct= " << model->HighEnergyActivationLimit()/MeV
330 << G4endl;
331 }
332
333 static const G4double limitdelta = 0.01*eV;
334 if(n > 0) {
335
336 // extend energy range to previous models
337 tmin = std::min(tmin, eHigh[n-1]);
338 tmax = std::max(tmax, eLow[0]);
339 //G4cout << "tmin= " << tmin << " tmax= "
340 // << tmax << " ord= " << ord <<G4endl;
341 // empty energy range
342 if( tmax - tmin <= limitdelta) { push = false; }
343 // low-energy model
344 else if (tmax == eLow[0]) {
345 push = false;
346 insert = true;
347 idx = 0;
348 // resolve intersections
349 } else if(tmin < eHigh[n-1]) {
350 // compare order
351 for(G4int k=0; k<n; ++k) {
352 // new model has higher order parameter,
353 // so, its application area may be reduced
354 // to avoid intersections
355 if(ord >= modelOrd[k]) {
356 if(tmin < eHigh[k] && tmin >= eLow[k]) { tmin = eHigh[k]; }
357 if(tmax <= eHigh[k] && tmax > eLow[k]) { tmax = eLow[k]; }
358 if(tmax > eHigh[k] && tmin < eLow[k]) {
359 if(tmax - eHigh[k] > eLow[k] - tmin) { tmin = eHigh[k]; }
360 else { tmax = eLow[k]; }
361 }
362 if( tmax - tmin <= limitdelta) {
363 push = false;
364 break;
365 }
366 }
367 }
368 // this model has lower order parameter than possible
369 // other models, with which there may be intersections
370 // so, appliction area of such models may be reduced
371
372 // insert below the first model
373 if (tmax <= eLow[0]) {
374 push = false;
375 insert = true;
376 idx = 0;
377 // resolve intersections
378 } else if(tmin < eHigh[n-1]) {
379 // last energy interval
380 if(tmin > eLow[n-1] && tmax >= eHigh[n-1]) {
381 eHigh[n-1] = tmin;
382 // first energy interval
383 } else if(tmin <= eLow[0] && tmax < eHigh[0]) {
384 eLow[0] = tmax;
385 push = false;
386 insert = true;
387 idx = 0;
388 // loop over all models
389 } else {
390 for(G4int k=n-1; k>=0; --k) {
391 if(tmin <= eLow[k] && tmax >= eHigh[k]) {
392 // full overlap exclude previous model
393 isUsed[modelAtRegion[k]] = 0;
394 idx = k;
395 if(k < n-1) {
396 // shift upper models and change index
397 for(G4int kk=k; kk<n-1; ++kk) {
398 modelAtRegion[kk] = modelAtRegion[kk+1];
399 modelOrd[kk] = modelOrd[kk+1];
400 eLow[kk] = eLow[kk+1];
401 eHigh[kk] = eHigh[kk+1];
402 }
403 ++k;
404 }
405 --n;
406 } else {
407 // partially reduce previous model area
408 if(tmin <= eLow[k] && tmax > eLow[k]) {
409 eLow[k] = tmax;
410 idx = k;
411 insert = true;
412 push = false;
413 } else if(tmin < eHigh[k] && tmax >= eHigh[k]) {
414 eHigh[k] = tmin;
415 idx = k + 1;
416 if(idx < n) {
417 insert = true;
418 push = false;
419 }
420 } else if(tmin > eLow[k] && tmax < eHigh[k]) {
421 if(eHigh[k] - tmax > tmin - eLow[k]) {
422 eLow[k] = tmax;
423 idx = k;
424 insert = true;
425 push = false;
426 } else {
427 eHigh[k] = tmin;
428 idx = k + 1;
429 if(idx < n) {
430 insert = true;
431 push = false;
432 }
433 }
434 }
435 }
436 }
437 }
438 }
439 }
440 }
441 // provide space for the new model
442 if(insert) {
443 for(G4int k=n-1; k>=idx; --k) {
444 modelAtRegion[k+1] = modelAtRegion[k];
445 modelOrd[k+1] = modelOrd[k];
446 eLow[k+1] = eLow[k];
447 eHigh[k+1] = eHigh[k];
448 }
449 }
450 //G4cout << "push= " << push << " insert= " << insert
451 // << " idx= " << idx <<G4endl;
452 // the model is added
453 if (push || insert) {
454 ++n;
455 modelAtRegion[idx] = ii;
456 modelOrd[idx] = ord;
457 eLow[idx] = tmin;
458 eHigh[idx] = tmax;
459 isUsed[ii] = 1;
460 }
461 // exclude models with zero energy range
462 for(G4int k=n-1; k>=0; --k) {
463 if(eHigh[k] - eLow[k] <= limitdelta) {
464 isUsed[modelAtRegion[k]] = 0;
465 if(k < n-1) {
466 for(G4int kk=k; kk<n-1; ++kk) {
467 modelAtRegion[kk] = modelAtRegion[kk+1];
468 modelOrd[kk] = modelOrd[kk+1];
469 eLow[kk] = eLow[kk+1];
470 eHigh[kk] = eHigh[kk+1];
471 }
472 }
473 --n;
474 }
475 }
476 }
477 }
478 eLow[0] = 0.0;
479 eLow[n] = eHigh[n-1];
480
481 if(1 < verboseLevel) {
482 G4cout << "### New G4RegionModels set with " << n
483 << " models for region <";
484 if (region) { G4cout << region->GetName(); }
485 G4cout << "> Elow(MeV)= ";
486 for(G4int iii=0; iii<=n; ++iii) {G4cout << eLow[iii]/MeV << " ";}
487 G4cout << G4endl;
488 }
489 G4RegionModels* rm = new G4RegionModels(n, modelAtRegion, eLow, region);
490 setOfRegionModels[reg] = rm;
491 // shortcut
492 if(1 == nEmModels) { break; }
493 }
494
495 currRegionModel = setOfRegionModels[0];
496 currModel = models[0];
497
498 // Access to materials and build cuts
499 size_t idx = 1;
500 if(secondaryParticle) {
501 if( secondaryParticle == G4Gamma::Gamma() ) { idx = 0; }
502 else if( secondaryParticle == G4Electron::Electron()) { idx = 1; }
503 else if( secondaryParticle == G4Positron::Positron()) { idx = 2; }
504 else { idx = 3; }
505 }
506
507 theCuts =
508 static_cast<const G4DataVector*>(theCoupleTable->GetEnergyCutsVector(idx));
509
510 // for the second run the check on cuts should be repeated
511 if(theCutsNew) { *theCutsNew = *theCuts; }
512
513 if(minSubRange < 1.0) {
514 if( !theSubCuts ) { theSubCuts = new G4DataVector(); }
515 theSubCuts->resize(numOfCouples,DBL_MAX);
516 }
517
518 // G4cout << "========Start define cuts" << G4endl;
519 // define cut values
520 for(size_t i=0; i<numOfCouples; ++i) {
521
522 const G4MaterialCutsCouple* couple =
523 theCoupleTable->GetMaterialCutsCouple(i);
524 const G4Material* material = couple->GetMaterial();
525 const G4ProductionCuts* pcuts = couple->GetProductionCuts();
526
527 G4int reg = 0;
528 if(nRegions > 1 && nEmModels > 1) {
529 reg = nRegions;
530 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
531 do {--reg;} while (reg>0 && pcuts != (setr[reg]->GetProductionCuts()));
532 idxOfRegionModels[i] = reg;
533 }
534 if(1 < verboseLevel) {
535 G4cout << "G4EmModelManager::Initialise() for "
536 << material->GetName()
537 << " indexOfCouple= " << i
538 << " indexOfRegion= " << reg
539 << G4endl;
540 }
541
542 G4double cut = (*theCuts)[i];
543 if(secondaryParticle) {
544
545 // compute subcut
546 if( cut < DBL_MAX && minSubRange < 1.0) {
547 G4double subcut = minSubRange*cut;
548 G4double rcut = std::min(minSubRange*pcuts->GetProductionCut(idx),
549 maxSubCutInRange);
550 G4double tcutmax =
551 theCoupleTable->ConvertRangeToEnergy(secondaryParticle,
552 material,rcut);
553 if(tcutmax < subcut) { subcut = tcutmax; }
554 (*theSubCuts)[i] = subcut;
555 }
556
557 // note that idxOfRegionModels[] not always filled
558 G4int inn = 0;
559 G4int nnm = 1;
560 if(nRegions > 1 && nEmModels > 1) {
561 inn = idxOfRegionModels[i];
562 }
563 // check cuts and introduce upper limits
564 //G4cout << "idx= " << i << " cut(keV)= " << cut/keV << G4endl;
565 currRegionModel = setOfRegionModels[inn];
566 nnm = currRegionModel->NumberOfModels();
567
568 //G4cout << "idx= " << i << " Nmod= " << nnm << G4endl;
569
570 for(G4int jj=0; jj<nnm; ++jj) {
571 //G4cout << "jj= " << jj << " modidx= "
572 // << currRegionModel->ModelIndex(jj) << G4endl;
573 currModel = models[currRegionModel->ModelIndex(jj)];
574 G4double cutlim = currModel->MinEnergyCut(particle,couple);
575 if(cutlim > cut) {
576 if(!theCutsNew) { theCutsNew = new G4DataVector(*theCuts); }
577 (*theCutsNew)[i] = cutlim;
578 /*
579 G4cout << "### " << partname << " energy loss model in "
580 << material->GetName()
581 << " Cut was changed from " << cut/keV << " keV to "
582 << cutlim/keV << " keV " << " due to "
583 << currModel->GetName() << G4endl;
584 */
585 }
586 }
587 }
588 }
589 if(theCutsNew) { theCuts = theCutsNew; }
590
591 // initialize models
592 G4int nn = 0;
593 severalModels = true;
594 for(G4int jj=0; jj<nEmModels; ++jj) {
595 if(1 == isUsed[jj]) {
596 ++nn;
597 currModel = models[jj];
598 currModel->Initialise(particle, *theCuts);
599 if(flucModels[jj]) { flucModels[jj]->InitialiseMe(particle); }
600 }
601 }
602 if(1 == nn) { severalModels = false; }
603
604 if(1 < verboseLevel) {
605 G4cout << "G4EmModelManager for " << partname
606 << " is initialised; nRegions= " << nRegions
607 << " severalModels: " << severalModels
608 << G4endl;
609 }
610
611 return theCuts;
612}
@ 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
bool G4bool
Definition: G4Types.hh:86
static G4Electron * Electron()
Definition: G4Electron.cc:93
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
G4ProductionCuts * GetProductionCuts() const
static G4Positron * Positron()
Definition: G4Positron.cc:93
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
const std::vector< G4double > * GetEnergyCutsVector(std::size_t pcIdx) const
static G4ProductionCutsTable * GetProductionCutsTable()
G4double ConvertRangeToEnergy(const G4ParticleDefinition *particle, const G4Material *material, G4double range)
G4double GetProductionCut(G4int index) const
static G4RegionStore * GetInstance()
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:652
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:645
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0
virtual G4double MinEnergyCut(const G4ParticleDefinition *, const G4MaterialCutsCouple *)
Definition: G4VEmModel.cc:433

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

◆ NumberOfModels()

◆ NumberOfRegionModels()

G4int G4EmModelManager::NumberOfRegionModels ( size_t  index_couple) const

Definition at line 210 of file G4EmModelManager.cc.

211{
212 G4RegionModels* rm = setOfRegionModels[idxOfRegionModels[idx]];
213 return rm->NumberOfModels();
214}

Referenced by G4VEmProcess::GetNumberOfRegionModels().

◆ operator=()

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

◆ SelectModel()

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

Definition at line 240 of file G4EmModelManager.hh.

242{
243 if(severalModels) {
244 if(nRegions > 1) {
245 currRegionModel = setOfRegionModels[idxOfRegionModels[index]];
246 }
247 currModel = models[currRegionModel->SelectIndex(kinEnergy)];
248 }
249 return currModel;
250}

Referenced by G4VEnergyLossProcess::SelectModel(), G4VMultipleScattering::SelectModel(), G4VEmProcess::SelectModel(), G4VEnergyLossProcess::SelectModelForMaterial(), and G4VEmProcess::SelectModelForMaterial().

◆ SetFluoFlag()

void G4EmModelManager::SetFluoFlag ( G4bool  val)
inline

Definition at line 268 of file G4EmModelManager.hh.

269{
270 fluoFlag = val;
271}

Referenced by G4VEmProcess::PreparePhysicsTable().

◆ SubCutoff()

const G4DataVector * G4EmModelManager::SubCutoff ( ) const
inline

Definition at line 261 of file G4EmModelManager.hh.

262{
263 return theSubCuts;
264}

Referenced by G4VEnergyLossProcess::PreparePhysicsTable().

◆ UpdateEmModel()

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

Definition at line 166 of file G4EmModelManager.cc.

168{
169 if (nEmModels > 0) {
170 for(G4int i=0; i<nEmModels; ++i) {
171 if(nam == models[i]->GetName()) {
172 models[i]->SetLowEnergyLimit(emin);
173 models[i]->SetHighEnergyLimit(emax);
174 break;
175 }
176 }
177 }
178 G4cout << "G4EmModelManager::UpdateEmModel WARNING: no model <"
179 << nam << "> is found out"
180 << G4endl;
181}

Referenced by G4VEmProcess::UpdateEmModel(), and G4VEnergyLossProcess::UpdateEmModel().


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