#include <G4EmModelManager.hh>
Public Member Functions | |
G4EmModelManager () | |
~G4EmModelManager () | |
void | Clear () |
const G4DataVector * | Initialise (const G4ParticleDefinition *, const G4ParticleDefinition *, G4double, G4int) |
void | FillDEDXVector (G4PhysicsVector *, const G4MaterialCutsCouple *, G4EmTableType t=fRestricted) |
void | FillLambdaVector (G4PhysicsVector *, const G4MaterialCutsCouple *, G4bool startFromNull=true, G4EmTableType t=fRestricted) |
G4VEmModel * | GetModel (G4int, G4bool ver=false) |
void | AddEmModel (G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *) |
void | UpdateEmModel (const G4String &, G4double, G4double) |
void | DumpModelList (G4int verb) |
G4VEmModel * | SelectModel (G4double &energy, size_t &index) |
const G4DataVector * | Cuts () const |
const G4DataVector * | SubCutoff () const |
void | SetFluoFlag (G4bool val) |
G4int | NumberOfModels () const |
Definition at line 142 of file G4EmModelManager.hh.
G4EmModelManager::G4EmModelManager | ( | ) |
Definition at line 124 of file G4EmModelManager.cc.
00124 : 00125 nEmModels(0), 00126 nRegions(0), 00127 particle(0), 00128 verboseLevel(0) 00129 { 00130 maxSubCutInRange = 0.7*mm; 00131 models.reserve(4); 00132 flucModels.reserve(4); 00133 regions.reserve(4); 00134 orderOfModels.reserve(4); 00135 isUsed.reserve(4); 00136 severalModels = true; 00137 fluoFlag = false; 00138 currRegionModel = 0; 00139 currModel = 0; 00140 theCuts = 0; 00141 theSubCuts = 0; 00142 }
G4EmModelManager::~G4EmModelManager | ( | ) |
Definition at line 146 of file G4EmModelManager.cc.
References Clear().
00147 { 00148 verboseLevel = 0; // no verbosity at destruction 00149 Clear(); 00150 delete theSubCuts; 00151 }
void G4EmModelManager::AddEmModel | ( | G4int | , | |
G4VEmModel * | , | |||
G4VEmFluctuationModel * | , | |||
const G4Region * | ||||
) |
Definition at line 171 of file G4EmModelManager.cc.
References G4VEmModel::DefineForRegion(), G4cout, and G4endl.
Referenced by G4VMultipleScattering::AddEmModel(), G4VEnergyLossProcess::AddEmModel(), G4VEmProcess::AddEmModel(), and G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel().
00173 { 00174 if(!p) { 00175 G4cout << "G4EmModelManager::AddEmModel WARNING: no model defined." 00176 << G4endl; 00177 return; 00178 } 00179 models.push_back(p); 00180 flucModels.push_back(fm); 00181 regions.push_back(r); 00182 orderOfModels.push_back(num); 00183 isUsed.push_back(0); 00184 p->DefineForRegion(r); 00185 ++nEmModels; 00186 }
void G4EmModelManager::Clear | ( | ) |
Definition at line 155 of file G4EmModelManager.cc.
References G4cout, G4endl, and CLHEP::detail::n.
Referenced by Initialise(), and ~G4EmModelManager().
00156 { 00157 if(1 < verboseLevel) { 00158 G4cout << "G4EmModelManager::Clear()" << G4endl; 00159 } 00160 size_t n = setOfRegionModels.size(); 00161 if(n > 0) { 00162 for(size_t i=0; i<n; ++i) { 00163 delete setOfRegionModels[i]; 00164 setOfRegionModels[i] = 0; 00165 } 00166 } 00167 }
const G4DataVector * G4EmModelManager::Cuts | ( | ) | const [inline] |
void G4EmModelManager::DumpModelList | ( | G4int | verb | ) |
Definition at line 663 of file G4EmModelManager.cc.
References G4InuclParticleNames::an, G4BestUnit, G4cout, G4endl, G4Region::GetName(), CLHEP::detail::n, and G4InuclParticleNames::nn.
Referenced by G4VMultipleScattering::BuildPhysicsTable(), G4VMultipleScattering::PrintInfoDefinition(), G4VEnergyLossProcess::PrintInfoDefinition(), and G4VEmProcess::PrintInfoDefinition().
00664 { 00665 if(verb == 0) { return; } 00666 for(G4int i=0; i<nRegions; ++i) { 00667 G4RegionModels* r = setOfRegionModels[i]; 00668 const G4Region* reg = r->Region(); 00669 G4int n = r->NumberOfModels(); 00670 if(n > 0) { 00671 G4cout << " ===== EM models for the G4Region " << reg->GetName() 00672 << " ======" << G4endl;; 00673 for(G4int j=0; j<n; ++j) { 00674 G4VEmModel* model = models[r->ModelIndex(j)]; 00675 G4double emin = 00676 std::max(r->LowEdgeEnergy(j),model->LowEnergyActivationLimit()); 00677 G4double emax = 00678 std::min(r->LowEdgeEnergy(j+1),model->HighEnergyActivationLimit()); 00679 G4cout << std::setw(20); 00680 G4cout << model->GetName() << " : Emin= " 00681 << std::setw(8) << G4BestUnit(emin,"Energy") 00682 << " Emax= " 00683 << std::setw(8) << G4BestUnit(emax,"Energy"); 00684 G4PhysicsTable* table = model->GetCrossSectionTable(); 00685 if(table) { 00686 size_t kk = table->size(); 00687 for(size_t k=0; k<kk; ++k) { 00688 G4PhysicsVector* v = (*table)[k]; 00689 if(v) { 00690 G4int nn = v->GetVectorLength() - 1; 00691 G4cout << " Table with " << nn << " bins Emin= " 00692 << std::setw(6) << G4BestUnit(v->Energy(0),"Energy") 00693 << " Emax= " 00694 << std::setw(6) << G4BestUnit(v->Energy(nn),"Energy"); 00695 break; 00696 } 00697 } 00698 } 00699 G4VEmAngularDistribution* an = model->GetAngularDistribution(); 00700 if(an) { G4cout << " " << an->GetName(); } 00701 if(fluoFlag && model->DeexcitationFlag()) { G4cout << " FluoActive"; } 00702 G4cout << G4endl; 00703 } 00704 } 00705 if(1 == nEmModels) { break; } 00706 } 00707 }
void G4EmModelManager::FillDEDXVector | ( | G4PhysicsVector * | , | |
const G4MaterialCutsCouple * | , | |||
G4EmTableType | t = fRestricted | |||
) |
Definition at line 510 of file G4EmModelManager.cc.
References DBL_MAX, G4PhysicsVector::Energy(), fSubRestricted, fTotal, G4cout, G4endl, G4MaterialCutsCouple::GetIndex(), G4MaterialCutsCouple::GetMaterial(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), G4PhysicsVector::GetVectorLength(), G4InuclParticleNames::k0, G4RegionModels::LowEdgeEnergy(), G4RegionModels::ModelIndex(), G4RegionModels::NumberOfModels(), and G4PhysicsVector::PutValue().
Referenced by G4VEnergyLossProcess::BuildDEDXTable().
00513 { 00514 size_t i = couple->GetIndex(); 00515 G4double cut = (*theCuts)[i]; 00516 G4double emin = 0.0; 00517 00518 if(fTotal == tType) { cut = DBL_MAX; } 00519 else if(fSubRestricted == tType) { 00520 emin = cut; 00521 if(theSubCuts) { emin = (*theSubCuts)[i]; } 00522 } 00523 00524 if(1 < verboseLevel) { 00525 G4cout << "G4EmModelManager::FillDEDXVector() for " 00526 << couple->GetMaterial()->GetName() 00527 << " cut(MeV)= " << cut 00528 << " emin(MeV)= " << emin 00529 << " Type " << tType 00530 << " for " << particle->GetParticleName() 00531 << G4endl; 00532 } 00533 00534 G4int reg = 0; 00535 if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; } 00536 const G4RegionModels* regModels = setOfRegionModels[reg]; 00537 G4int nmod = regModels->NumberOfModels(); 00538 00539 // Calculate energy losses vector 00540 00541 //G4cout << "nmod= " << nmod << G4endl; 00542 size_t totBinsLoss = aVector->GetVectorLength(); 00543 G4double del = 0.0; 00544 G4int k0 = 0; 00545 00546 for(size_t j=0; j<totBinsLoss; ++j) { 00547 00548 G4double e = aVector->Energy(j); 00549 00550 // Choose a model of energy losses 00551 G4int k = 0; 00552 if (nmod > 1) { 00553 k = nmod; 00554 do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k)); 00555 //G4cout << "k= " << k << G4endl; 00556 if(k > 0 && k != k0) { 00557 k0 = k; 00558 G4double elow = regModels->LowEdgeEnergy(k); 00559 G4double dedx1 = ComputeDEDX(models[regModels->ModelIndex(k-1)], 00560 couple,elow,cut,emin); 00561 G4double dedx2 = ComputeDEDX(models[regModels->ModelIndex(k)], 00562 couple,elow,cut,emin); 00563 del = 0.0; 00564 if(dedx2 > 0.0) { del = (dedx1/dedx2 - 1.0)*elow; } 00565 //G4cout << "elow= " << elow 00566 // << " dedx1= " << dedx1 << " dedx2= " << dedx2 << G4endl; 00567 } 00568 } 00569 G4double dedx = 00570 ComputeDEDX(models[regModels->ModelIndex(k)],couple,e,cut,emin); 00571 dedx *= (1.0 + del/e); 00572 00573 if(2 < verboseLevel) { 00574 G4cout << "Material= " << couple->GetMaterial()->GetName() 00575 << " E(MeV)= " << e/MeV 00576 << " dEdx(MeV/mm)= " << dedx*mm/MeV 00577 << " del= " << del*mm/MeV<< " k= " << k 00578 << " modelIdx= " << regModels->ModelIndex(k) 00579 << G4endl; 00580 } 00581 if(dedx < 0.0) { dedx = 0.0; } 00582 aVector->PutValue(j, dedx); 00583 } 00584 }
void G4EmModelManager::FillLambdaVector | ( | G4PhysicsVector * | , | |
const G4MaterialCutsCouple * | , | |||
G4bool | startFromNull = true , |
|||
G4EmTableType | t = fRestricted | |||
) |
Definition at line 588 of file G4EmModelManager.cc.
References G4VEmModel::CrossSection(), DBL_MAX, G4PhysicsVector::Energy(), fIsCrossSectionPrim, fSubRestricted, G4cout, G4endl, G4MaterialCutsCouple::GetIndex(), G4MaterialCutsCouple::GetMaterial(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), G4PhysicsVector::GetVectorLength(), G4InuclParticleNames::k0, G4RegionModels::LowEdgeEnergy(), G4RegionModels::ModelIndex(), G4RegionModels::NumberOfModels(), and G4PhysicsVector::PutValue().
Referenced by G4VEnergyLossProcess::BuildLambdaTable().
00592 { 00593 size_t i = couple->GetIndex(); 00594 G4double cut = (*theCuts)[i]; 00595 G4double tmax = DBL_MAX; 00596 if (fSubRestricted == tType) { 00597 tmax = cut; 00598 if(theSubCuts) { cut = (*theSubCuts)[i]; } 00599 } 00600 00601 G4int reg = 0; 00602 if(nRegions > 1 && nEmModels > 1) { reg = idxOfRegionModels[i]; } 00603 const G4RegionModels* regModels = setOfRegionModels[reg]; 00604 G4int nmod = regModels->NumberOfModels(); 00605 if(1 < verboseLevel) { 00606 G4cout << "G4EmModelManager::FillLambdaVector() for " 00607 << particle->GetParticleName() 00608 << " in " << couple->GetMaterial()->GetName() 00609 << " Ecut(MeV)= " << cut 00610 << " Emax(MeV)= " << tmax 00611 << " Type " << tType 00612 << " nmod= " << nmod 00613 << G4endl; 00614 } 00615 00616 // Calculate lambda vector 00617 size_t totBinsLambda = aVector->GetVectorLength(); 00618 G4double del = 0.0; 00619 G4int k0 = 0; 00620 G4int k = 0; 00621 G4VEmModel* mod = models[regModels->ModelIndex(0)]; 00622 for(size_t j=0; j<totBinsLambda; ++j) { 00623 00624 G4double e = aVector->Energy(j); 00625 00626 // Choose a model 00627 if (nmod > 1) { 00628 k = nmod; 00629 do {--k;} while (k>0 && e <= regModels->LowEdgeEnergy(k)); 00630 if(k > 0 && k != k0) { 00631 k0 = k; 00632 G4double elow = regModels->LowEdgeEnergy(k); 00633 G4VEmModel* mod1 = models[regModels->ModelIndex(k-1)]; 00634 G4double xs1 = mod1->CrossSection(couple,particle,elow,cut,tmax); 00635 mod = models[regModels->ModelIndex(k)]; 00636 G4double xs2 = mod->CrossSection(couple,particle,elow,cut,tmax); 00637 del = 0.0; 00638 if(xs2 > 0.0) { del = (xs1/xs2 - 1.0)*elow; } 00639 //G4cout << "New model k=" << k << " E(MeV)= " << e/MeV 00640 // << " Elow(MeV)= " << elow/MeV << " del= " << del << G4endl; 00641 } 00642 } 00643 G4double cross = mod->CrossSection(couple,particle,e,cut,tmax); 00644 cross *= (1.0 + del/e); 00645 if(fIsCrossSectionPrim == tType) { cross *= e; } 00646 00647 if(j==0 && startFromNull) { cross = 0.0; } 00648 00649 if(2 < verboseLevel) { 00650 G4cout << "FillLambdaVector: " << j << ". e(MeV)= " << e/MeV 00651 << " cross(1/mm)= " << cross*mm 00652 << " del= " << del*mm << " k= " << k 00653 << " modelIdx= " << regModels->ModelIndex(k) 00654 << G4endl; 00655 } 00656 if(cross < 0.0) { cross = 0.0; } 00657 aVector->PutValue(j, cross); 00658 } 00659 }
G4VEmModel * G4EmModelManager::GetModel | ( | G4int | , | |
G4bool | ver = false | |||
) |
Definition at line 209 of file G4EmModelManager.cc.
References G4cout, G4endl, and G4ParticleDefinition::GetParticleName().
Referenced by G4VMultipleScattering::GetModelByIndex(), G4VEnergyLossProcess::GetModelByIndex(), G4VEmProcess::GetModelByIndex(), G4VMultipleScattering::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), G4VEmProcess::PreparePhysicsTable(), G4VMultipleScattering::SetIonisation(), G4VMultipleScattering::StartTracking(), and G4VMultipleScattering::StorePhysicsTable().
00210 { 00211 G4VEmModel* model = 0; 00212 if(i >= 0 && i < nEmModels) { model = models[i]; } 00213 else if(verboseLevel > 0 && ver) { 00214 G4cout << "G4EmModelManager::GetModel WARNING: " 00215 << "index " << i << " is wrong Nmodels= " 00216 << nEmModels; 00217 if(particle) G4cout << " for " << particle->GetParticleName(); 00218 G4cout<< G4endl; 00219 } 00220 return model; 00221 }
const G4DataVector * G4EmModelManager::Initialise | ( | const G4ParticleDefinition * | , | |
const G4ParticleDefinition * | , | |||
G4double | , | |||
G4int | ||||
) |
Definition at line 226 of file G4EmModelManager.cc.
References Clear(), G4ProductionCutsTable::ConvertRangeToEnergy(), DBL_MAX, FatalException, G4cout, G4endl, G4Exception(), G4Gamma::Gamma(), G4ProductionCutsTable::GetEnergyCutsVector(), G4RegionStore::GetInstance(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), G4ProductionCuts::GetProductionCut(), G4ProductionCutsTable::GetProductionCutsTable(), G4RegionStore::GetRegion(), G4ProductionCutsTable::GetTableSize(), G4VEmModel::Initialise(), CLHEP::detail::n, G4InuclParticleNames::nn, and G4Positron::Positron().
Referenced by G4AdjointBremsstrahlungModel::AdjointCrossSection(), G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecond(), G4VMultipleScattering::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), and G4VEmProcess::PreparePhysicsTable().
00230 { 00231 verboseLevel = val; 00232 G4String partname = p->GetParticleName(); 00233 if(1 < verboseLevel) { 00234 G4cout << "G4EmModelManager::Initialise() for " 00235 << partname << G4endl; 00236 } 00237 // Are models defined? 00238 if(nEmModels < 1) { 00239 G4ExceptionDescription ed; 00240 ed << "No models found out for " << p->GetParticleName() 00241 << " !" << G4endl; 00242 G4Exception("G4EmModelManager::Initialise","em0002", 00243 FatalException, ed); 00244 } 00245 00246 particle = p; 00247 Clear(); // needed if run is not first 00248 00249 G4RegionStore* regionStore = G4RegionStore::GetInstance(); 00250 const G4Region* world = 00251 regionStore->GetRegion("DefaultRegionForTheWorld", false); 00252 00253 // Identify the list of regions with different set of models 00254 nRegions = 1; 00255 std::vector<const G4Region*> setr; 00256 setr.push_back(world); 00257 G4bool isWorld = false; 00258 00259 for (G4int ii=0; ii<nEmModels; ++ii) { 00260 const G4Region* r = regions[ii]; 00261 if ( r == 0 || r == world) { 00262 isWorld = true; 00263 regions[ii] = world; 00264 } else { 00265 G4bool newRegion = true; 00266 if (nRegions>1) { 00267 for (G4int j=1; j<nRegions; ++j) { 00268 if ( r == setr[j] ) newRegion = false; 00269 } 00270 } 00271 if (newRegion) { 00272 setr.push_back(r); 00273 nRegions++; 00274 } 00275 } 00276 } 00277 // Are models defined? 00278 if(!isWorld) { 00279 G4ExceptionDescription ed; 00280 ed << "No models defined for the World volume for " << p->GetParticleName() 00281 << " !" << G4endl; 00282 G4Exception("G4EmModelManager::Initialise","em0002", 00283 FatalException, ed); 00284 } 00285 00286 G4ProductionCutsTable* theCoupleTable= 00287 G4ProductionCutsTable::GetProductionCutsTable(); 00288 size_t numOfCouples = theCoupleTable->GetTableSize(); 00289 00290 // prepare vectors, shortcut for the case of only 1 model 00291 if(nRegions > 1 && nEmModels > 1) { 00292 if(numOfCouples > idxOfRegionModels.size()) { 00293 idxOfRegionModels.resize(numOfCouples); 00294 } 00295 } 00296 size_t nr = 1; 00297 if(nEmModels > 1) { nr = nRegions; } 00298 if(nr > setOfRegionModels.size()) { setOfRegionModels.resize(nr); } 00299 00300 std::vector<G4int> modelAtRegion(nEmModels); 00301 std::vector<G4int> modelOrd(nEmModels); 00302 G4DataVector eLow(nEmModels+1); 00303 G4DataVector eHigh(nEmModels); 00304 00305 // Order models for regions 00306 for (G4int reg=0; reg<nRegions; ++reg) { 00307 const G4Region* region = setr[reg]; 00308 G4int n = 0; 00309 00310 for (G4int ii=0; ii<nEmModels; ++ii) { 00311 00312 G4VEmModel* model = models[ii]; 00313 if ( region == regions[ii] ) { 00314 00315 G4double tmin = model->LowEnergyLimit(); 00316 G4double tmax = model->HighEnergyLimit(); 00317 G4int ord = orderOfModels[ii]; 00318 G4bool push = true; 00319 G4bool insert = false; 00320 G4int idx = n; 00321 00322 if(1 < verboseLevel) { 00323 G4cout << "Model #" << ii 00324 << " <" << model->GetName() << "> for region <"; 00325 if (region) G4cout << region->GetName(); 00326 G4cout << "> " 00327 << " tmin(MeV)= " << tmin/MeV 00328 << "; tmax(MeV)= " << tmax/MeV 00329 << "; order= " << ord 00330 << G4endl; 00331 } 00332 00333 if(n > 0) { 00334 00335 // extend energy range to previous models 00336 tmin = std::min(tmin, eHigh[n-1]); 00337 tmax = std::max(tmax, eLow[0]); 00338 //G4cout << "tmin= " << tmin << " tmax= " 00339 // << tmax << " ord= " << ord <<G4endl; 00340 // empty energy range 00341 if( tmax - tmin <= eV) push = false; 00342 // low-energy model 00343 else if (tmax == eLow[0]) { 00344 push = false; 00345 insert = true; 00346 idx = 0; 00347 // resolve intersections 00348 } else if(tmin < eHigh[n-1]) { 00349 // compare order 00350 for(G4int k=0; k<n; ++k) { 00351 // new model has lower application 00352 if(ord >= modelOrd[k]) { 00353 if(tmin < eHigh[k] && tmin >= eLow[k]) tmin = eHigh[k]; 00354 if(tmax <= eHigh[k] && tmax > eLow[k]) tmax = eLow[k]; 00355 if(tmax > eHigh[k] && tmin < eLow[k]) { 00356 if(tmax - eHigh[k] > eLow[k] - tmin) tmin = eHigh[k]; 00357 else tmax = eLow[k]; 00358 } 00359 if( tmax - tmin <= eV) { 00360 push = false; 00361 break; 00362 } 00363 } 00364 } 00365 //G4cout << "tmin= " << tmin << " tmax= " 00366 // << tmax << " push= " << push << " idx= " << idx <<G4endl; 00367 if(push) { 00368 if (tmax == eLow[0]) { 00369 push = false; 00370 insert = true; 00371 idx = 0; 00372 // continue resolve intersections 00373 } else if(tmin < eHigh[n-1]) { 00374 // last energy interval 00375 if(tmin > eLow[n-1] && tmax >= eHigh[n-1]) { 00376 eHigh[n-1] = tmin; 00377 // first energy interval 00378 } else if(tmin <= eLow[0] && tmax < eHigh[0]) { 00379 eLow[0] = tmax; 00380 push = false; 00381 insert = true; 00382 idx = 0; 00383 } else { 00384 // find energy interval to replace 00385 for(G4int k=0; k<n; ++k) { 00386 if(tmin <= eLow[k] && tmax >= eHigh[k]) { 00387 push = false; 00388 modelAtRegion[k] = ii; 00389 modelOrd[k] = ord; 00390 isUsed[ii] = 1; 00391 } 00392 } 00393 } 00394 } 00395 } 00396 } 00397 } 00398 if(insert) { 00399 for(G4int k=n-1; k>=idx; --k) { 00400 modelAtRegion[k+1] = modelAtRegion[k]; 00401 modelOrd[k+1] = modelOrd[k]; 00402 eLow[k+1] = eLow[k]; 00403 eHigh[k+1] = eHigh[k]; 00404 } 00405 } 00406 //G4cout << "push= " << push << " insert= " << insert 00407 //<< " idx= " << idx <<G4endl; 00408 if (push || insert) { 00409 ++n; 00410 modelAtRegion[idx] = ii; 00411 modelOrd[idx] = ord; 00412 eLow[idx] = tmin; 00413 eHigh[idx] = tmax; 00414 isUsed[ii] = 1; 00415 } 00416 } 00417 } 00418 eLow[0] = 0.0; 00419 eLow[n] = eHigh[n-1]; 00420 00421 if(1 < verboseLevel) { 00422 G4cout << "New G4RegionModels set with " << n << " models for region <"; 00423 if (region) G4cout << region->GetName(); 00424 G4cout << "> Elow(MeV)= "; 00425 for(G4int ii=0; ii<=n; ++ii) {G4cout << eLow[ii]/MeV << " ";} 00426 G4cout << G4endl; 00427 } 00428 G4RegionModels* rm = new G4RegionModels(n, modelAtRegion, eLow, region); 00429 setOfRegionModels[reg] = rm; 00430 if(1 == nEmModels) { break; } 00431 } 00432 00433 currRegionModel = setOfRegionModels[0]; 00434 00435 // Access to materials and build cuts 00436 size_t idx = 1; 00437 if(secondaryParticle) { 00438 if( secondaryParticle == G4Gamma::Gamma() ) { idx = 0; } 00439 else if( secondaryParticle == G4Positron::Positron()) { idx = 2; } 00440 } 00441 00442 //theCuts = theCoupleTable->GetEnergyCutsVector(idx); 00443 theCuts = static_cast<const G4DataVector*>(theCoupleTable->GetEnergyCutsVector(idx)); 00444 00445 if(minSubRange < 1.0) { 00446 if( !theSubCuts ) { theSubCuts = new G4DataVector(); } 00447 theSubCuts->resize(numOfCouples,DBL_MAX); 00448 } 00449 for(size_t i=0; i<numOfCouples; ++i) { 00450 00451 const G4MaterialCutsCouple* couple = 00452 theCoupleTable->GetMaterialCutsCouple(i); 00453 const G4Material* material = couple->GetMaterial(); 00454 const G4ProductionCuts* pcuts = couple->GetProductionCuts(); 00455 00456 G4int reg = 0; 00457 if(nRegions > 1 && nEmModels > 1) { 00458 reg = nRegions; 00459 do {--reg;} while (reg>0 && pcuts != (setr[reg]->GetProductionCuts())); 00460 idxOfRegionModels[i] = reg; 00461 } 00462 if(1 < verboseLevel) { 00463 G4cout << "G4EmModelManager::Initialise() for " 00464 << material->GetName() 00465 << " indexOfCouple= " << i 00466 << " indexOfRegion= " << reg 00467 << G4endl; 00468 } 00469 00470 G4double cut = (*theCuts)[i]; 00471 if(secondaryParticle) { 00472 00473 // compute subcut 00474 if( cut < DBL_MAX && minSubRange < 1.0) { 00475 G4double subcut = minSubRange*cut; 00476 G4double rcut = std::min(minSubRange*pcuts->GetProductionCut(idx), 00477 maxSubCutInRange); 00478 G4double tcutmax = 00479 theCoupleTable->ConvertRangeToEnergy(secondaryParticle,material,rcut); 00480 if(tcutmax < subcut) { subcut = tcutmax; } 00481 (*theSubCuts)[i] = subcut; 00482 } 00483 } 00484 } 00485 00486 // initialize models 00487 G4int nn = 0; 00488 severalModels = true; 00489 for(G4int jj=0; jj<nEmModels; ++jj) { 00490 if(1 == isUsed[jj]) { 00491 ++nn; 00492 currModel = models[jj]; 00493 currModel->Initialise(particle, *theCuts); 00494 if(flucModels[jj]) { flucModels[jj]->InitialiseMe(particle); } 00495 } 00496 } 00497 if(1 == nn) { severalModels = false; } 00498 00499 if(1 < verboseLevel) { 00500 G4cout << "G4EmModelManager for " << particle->GetParticleName() 00501 << " is initialised; nRegions= " << nRegions 00502 << G4endl; 00503 } 00504 00505 return theCuts; 00506 }
G4int G4EmModelManager::NumberOfModels | ( | ) | const [inline] |
Definition at line 264 of file G4EmModelManager.hh.
Referenced by G4VEnergyLossProcess::NumberOfModels(), G4VMultipleScattering::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), G4VEmProcess::PreparePhysicsTable(), and G4VMultipleScattering::StorePhysicsTable().
G4VEmModel * G4EmModelManager::SelectModel | ( | G4double & | energy, | |
size_t & | index | |||
) | [inline] |
Definition at line 229 of file G4EmModelManager.hh.
References G4RegionModels::SelectIndex().
Referenced by G4VMultipleScattering::SelectModel(), G4VEnergyLossProcess::SelectModel(), G4VEmProcess::SelectModel(), G4VEnergyLossProcess::SelectModelForMaterial(), and G4VEmProcess::SelectModelForMaterial().
00231 { 00232 if(severalModels) { 00233 if(nRegions > 1) { 00234 currRegionModel = setOfRegionModels[idxOfRegionModels[index]]; 00235 } 00236 currModel = models[currRegionModel->SelectIndex(kinEnergy)]; 00237 } 00238 return currModel; 00239 }
void G4EmModelManager::SetFluoFlag | ( | G4bool | val | ) | [inline] |
Definition at line 257 of file G4EmModelManager.hh.
Referenced by G4VEmProcess::PreparePhysicsTable().
const G4DataVector * G4EmModelManager::SubCutoff | ( | ) | const [inline] |
Definition at line 250 of file G4EmModelManager.hh.
Referenced by G4VEnergyLossProcess::PreparePhysicsTable().
Definition at line 190 of file G4EmModelManager.cc.
References G4cout, and G4endl.
Referenced by G4VEnergyLossProcess::UpdateEmModel(), and G4VEmProcess::UpdateEmModel().
00192 { 00193 if (nEmModels > 0) { 00194 for(G4int i=0; i<nEmModels; ++i) { 00195 if(nam == models[i]->GetName()) { 00196 models[i]->SetLowEnergyLimit(emin); 00197 models[i]->SetHighEnergyLimit(emax); 00198 break; 00199 } 00200 } 00201 } 00202 G4cout << "G4EmModelManager::UpdateEmModel WARNING: no model <" 00203 << nam << "> is found out" 00204 << G4endl; 00205 }