00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081 #include "G4LossTableManager.hh"
00082 #include "G4SystemOfUnits.hh"
00083 #include "G4EnergyLossMessenger.hh"
00084 #include "G4PhysicsTable.hh"
00085 #include "G4ParticleDefinition.hh"
00086 #include "G4MaterialCutsCouple.hh"
00087 #include "G4ProcessManager.hh"
00088 #include "G4Electron.hh"
00089 #include "G4Proton.hh"
00090 #include "G4VMultipleScattering.hh"
00091 #include "G4VEmProcess.hh"
00092 #include "G4ProductionCutsTable.hh"
00093 #include "G4PhysicsTableHelper.hh"
00094 #include "G4EmCorrections.hh"
00095 #include "G4EmSaturation.hh"
00096 #include "G4EmConfigurator.hh"
00097 #include "G4ElectronIonPair.hh"
00098 #include "G4EmTableType.hh"
00099 #include "G4LossTableBuilder.hh"
00100 #include "G4VAtomDeexcitation.hh"
00101 #include "G4Region.hh"
00102
00103 G4LossTableManager* G4LossTableManager::theInstance = 0;
00104
00105
00106
00107 G4LossTableManager* G4LossTableManager::Instance()
00108 {
00109 if(0 == theInstance) {
00110 static G4LossTableManager manager;
00111 theInstance = &manager;
00112 }
00113 return theInstance;
00114 }
00115
00116
00117
00118 G4LossTableManager::~G4LossTableManager()
00119 {
00120 for (G4int i=0; i<n_loss; ++i) {
00121 if( loss_vector[i] ) { delete loss_vector[i]; }
00122 }
00123 size_t msc = msc_vector.size();
00124 for (size_t j=0; j<msc; ++j) {
00125 if( msc_vector[j] ) { delete msc_vector[j]; }
00126 }
00127 size_t emp = emp_vector.size();
00128 for (size_t k=0; k<emp; ++k) {
00129 if( emp_vector[k] ) { delete emp_vector[k]; }
00130 }
00131 size_t mod = mod_vector.size();
00132 for (size_t a=0; a<mod; ++a) {
00133 if( mod_vector[a] ) { delete mod_vector[a]; }
00134 }
00135 size_t fmod = fmod_vector.size();
00136 for (size_t b=0; b<fmod; ++b) {
00137 if( fmod_vector[b] ) { delete fmod_vector[b]; }
00138 }
00139 Clear();
00140 delete theMessenger;
00141 delete tableBuilder;
00142 delete emCorrections;
00143 delete emSaturation;
00144 delete emConfigurator;
00145 delete emElectronIonPair;
00146 delete atomDeexcitation;
00147 }
00148
00149
00150
00151 G4LossTableManager::G4LossTableManager()
00152 {
00153 n_loss = 0;
00154 run = 0;
00155 startInitialisation = false;
00156 all_tables_are_built = false;
00157 currentLoss = 0;
00158 currentParticle = 0;
00159 firstParticle = 0;
00160 lossFluctuationFlag = true;
00161 subCutoffFlag = false;
00162 rndmStepFlag = false;
00163 minSubRange = 0.0;
00164 maxRangeVariation = 1.0;
00165 maxFinalStep = 0.0;
00166 minKinEnergy = 0.1*keV;
00167 maxKinEnergy = 10.0*TeV;
00168 nbinsLambda = 77;
00169 nbinsPerDecade = 7;
00170 maxKinEnergyForMuons = 10.*TeV;
00171 integral = true;
00172 integralActive = false;
00173 buildCSDARange = false;
00174 minEnergyActive = false;
00175 maxEnergyActive = false;
00176 maxEnergyForMuonsActive = false;
00177 stepFunctionActive = false;
00178 flagLPM = true;
00179 splineFlag = true;
00180 bremsTh = DBL_MAX;
00181 factorForAngleLimit = 1.0;
00182 verbose = 1;
00183 theMessenger = new G4EnergyLossMessenger();
00184 theElectron = G4Electron::Electron();
00185 tableBuilder = new G4LossTableBuilder();
00186 emCorrections= new G4EmCorrections();
00187 emSaturation = new G4EmSaturation();
00188 emConfigurator = new G4EmConfigurator(verbose);
00189 emElectronIonPair = new G4ElectronIonPair();
00190 tableBuilder->SetSplineFlag(splineFlag);
00191 atomDeexcitation = 0;
00192 }
00193
00194
00195
00196 void G4LossTableManager::Clear()
00197 {
00198 all_tables_are_built = false;
00199 currentLoss = 0;
00200 currentParticle = 0;
00201 if(n_loss)
00202 {
00203 dedx_vector.clear();
00204 range_vector.clear();
00205 inv_range_vector.clear();
00206 loss_map.clear();
00207 loss_vector.clear();
00208 part_vector.clear();
00209 base_part_vector.clear();
00210 tables_are_built.clear();
00211 isActive.clear();
00212 n_loss = 0;
00213 }
00214 }
00215
00216
00217
00218 void G4LossTableManager::Register(G4VEnergyLossProcess* p)
00219 {
00220 if(!p) { return; }
00221 for (G4int i=0; i<n_loss; ++i) {
00222 if(loss_vector[i] == p) { return; }
00223 }
00224 if(verbose > 1) {
00225 G4cout << "G4LossTableManager::Register G4VEnergyLossProcess : "
00226 << p->GetProcessName() << " idx= " << n_loss << G4endl;
00227 }
00228 ++n_loss;
00229 loss_vector.push_back(p);
00230 part_vector.push_back(0);
00231 base_part_vector.push_back(0);
00232 dedx_vector.push_back(0);
00233 range_vector.push_back(0);
00234 inv_range_vector.push_back(0);
00235 tables_are_built.push_back(false);
00236 isActive.push_back(true);
00237 all_tables_are_built = false;
00238 if(!lossFluctuationFlag) { p->SetLossFluctuations(false); }
00239 if(subCutoffFlag) { p->ActivateSubCutoff(true); }
00240 if(rndmStepFlag) { p->SetRandomStep(true); }
00241 if(stepFunctionActive) { p->SetStepFunction(maxRangeVariation, maxFinalStep); }
00242 if(integralActive) { p->SetIntegral(integral); }
00243 if(minEnergyActive) { p->SetMinKinEnergy(minKinEnergy); }
00244 if(maxEnergyActive) { p->SetMaxKinEnergy(maxKinEnergy); }
00245 }
00246
00247
00248
00249 void G4LossTableManager::DeRegister(G4VEnergyLossProcess* p)
00250 {
00251 if(!p) { return; }
00252 for (G4int i=0; i<n_loss; ++i) {
00253 if(loss_vector[i] == p) { loss_vector[i] = 0; }
00254 }
00255 }
00256
00257
00258
00259 void G4LossTableManager::Register(G4VMultipleScattering* p)
00260 {
00261 if(!p) { return; }
00262 G4int n = msc_vector.size();
00263 for (G4int i=0; i<n; ++i) {
00264 if(msc_vector[i] == p) { return; }
00265 }
00266 if(verbose > 1) {
00267 G4cout << "G4LossTableManager::Register G4VMultipleScattering : "
00268 << p->GetProcessName() << " idx= " << msc_vector.size() << G4endl;
00269 }
00270 msc_vector.push_back(p);
00271 }
00272
00273
00274
00275 void G4LossTableManager::DeRegister(G4VMultipleScattering* p)
00276 {
00277 if(!p) { return; }
00278 size_t msc = msc_vector.size();
00279 for (size_t i=0; i<msc; ++i) {
00280 if(msc_vector[i] == p) { msc_vector[i] = 0; }
00281 }
00282 }
00283
00284
00285
00286 void G4LossTableManager::Register(G4VEmProcess* p)
00287 {
00288 if(!p) { return; }
00289 G4int n = emp_vector.size();
00290 for (G4int i=0; i<n; ++i) {
00291 if(emp_vector[i] == p) { return; }
00292 }
00293 if(verbose > 1) {
00294 G4cout << "G4LossTableManager::Register G4VEmProcess : "
00295 << p->GetProcessName() << " idx= " << emp_vector.size() << G4endl;
00296 }
00297 emp_vector.push_back(p);
00298 }
00299
00300
00301
00302 void G4LossTableManager::DeRegister(G4VEmProcess* p)
00303 {
00304 if(!p) { return; }
00305 size_t emp = emp_vector.size();
00306 for (size_t i=0; i<emp; ++i) {
00307 if(emp_vector[i] == p) { emp_vector[i] = 0; }
00308 }
00309 }
00310
00311
00312
00313 void G4LossTableManager::Register(G4VEmModel* p)
00314 {
00315 mod_vector.push_back(p);
00316 if(verbose > 1) {
00317 G4cout << "G4LossTableManager::Register G4VEmModel : "
00318 << p->GetName() << G4endl;
00319 }
00320 }
00321
00322
00323
00324 void G4LossTableManager::DeRegister(G4VEmModel* p)
00325 {
00326 size_t n = mod_vector.size();
00327 for (size_t i=0; i<n; ++i) {
00328 if(mod_vector[i] == p) { mod_vector[i] = 0; }
00329 }
00330 }
00331
00332
00333
00334 void G4LossTableManager::Register(G4VEmFluctuationModel* p)
00335 {
00336 fmod_vector.push_back(p);
00337 if(verbose > 1) {
00338 G4cout << "G4LossTableManager::Register G4VEmFluctuationModel : "
00339 << p->GetName() << G4endl;
00340 }
00341 }
00342
00343
00344
00345 void G4LossTableManager::DeRegister(G4VEmFluctuationModel* p)
00346 {
00347 size_t n = fmod_vector.size();
00348 for (size_t i=0; i<n; ++i) {
00349 if(fmod_vector[i] == p) { fmod_vector[i] = 0; }
00350 }
00351 }
00352
00353
00354
00355 void G4LossTableManager::RegisterIon(const G4ParticleDefinition* ion,
00356 G4VEnergyLossProcess* p)
00357 {
00358 loss_map[ion] = p;
00359 }
00360
00361
00362
00363 void G4LossTableManager::RegisterExtraParticle(
00364 const G4ParticleDefinition* part,
00365 G4VEnergyLossProcess* p)
00366 {
00367 if(!p || !part) { return; }
00368 for (G4int i=0; i<n_loss; ++i) {
00369 if(loss_vector[i] == p) { return; }
00370 }
00371 if(verbose > 1) {
00372 G4cout << "G4LossTableManager::RegisterExtraParticle "
00373 << part->GetParticleName() << " G4VEnergyLossProcess : "
00374 << p->GetProcessName() << " idx= " << n_loss << G4endl;
00375 }
00376 ++n_loss;
00377 loss_vector.push_back(p);
00378 part_vector.push_back(part);
00379 base_part_vector.push_back(p->BaseParticle());
00380 dedx_vector.push_back(0);
00381 range_vector.push_back(0);
00382 inv_range_vector.push_back(0);
00383 tables_are_built.push_back(false);
00384 all_tables_are_built = false;
00385 }
00386
00387
00388
00389 void
00390 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
00391 G4VEnergyLossProcess* p)
00392 {
00393 if (1 < verbose) {
00394 G4cout << "G4LossTableManager::PreparePhysicsTable for "
00395 << particle->GetParticleName()
00396 << " and " << p->GetProcessName() << " run= " << run
00397 << " loss_vector " << loss_vector.size() << G4endl;
00398 }
00399 if(!startInitialisation) { tableBuilder->SetInitialisationFlag(false); }
00400
00401
00402 startInitialisation = true;
00403
00404 if( 0 == run ) {
00405 emConfigurator->PrepareModels(particle, p);
00406
00407
00408 for (G4int j=0; j<n_loss; ++j) {
00409 if (p == loss_vector[j]) {
00410 if (!part_vector[j]) { part_vector[j] = particle; }
00411 }
00412 }
00413 }
00414 }
00415
00416
00417
00418 void
00419 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
00420 G4VEmProcess* p)
00421 {
00422 if (1 < verbose) {
00423 G4cout << "G4LossTableManager::PreparePhysicsTable for "
00424 << particle->GetParticleName()
00425 << " and " << p->GetProcessName() << G4endl;
00426 }
00427 if(!startInitialisation) { tableBuilder->SetInitialisationFlag(false); }
00428
00429
00430 if( 0 == run ) {
00431 emConfigurator->PrepareModels(particle, p);
00432 }
00433 startInitialisation = true;
00434 }
00435
00436
00437
00438 void
00439 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
00440 G4VMultipleScattering* p)
00441 {
00442 if (1 < verbose) {
00443 G4cout << "G4LossTableManager::PreparePhysicsTable for "
00444 << particle->GetParticleName()
00445 << " and " << p->GetProcessName() << G4endl;
00446 }
00447 if(!startInitialisation) { tableBuilder->SetInitialisationFlag(false); }
00448
00449
00450 if( 0 == run ) {
00451 emConfigurator->PrepareModels(particle, p);
00452 }
00453 }
00454
00455
00456
00457 void
00458 G4LossTableManager::BuildPhysicsTable(const G4ParticleDefinition*)
00459 {
00460 if(0 == run && startInitialisation) {
00461 emConfigurator->Clear();
00462 }
00463 }
00464
00465
00466
00467 void G4LossTableManager::BuildPhysicsTable(
00468 const G4ParticleDefinition* aParticle,
00469 G4VEnergyLossProcess* p)
00470 {
00471 if(1 < verbose) {
00472 G4cout << "### G4LossTableManager::BuildDEDXTable() is requested for "
00473 << aParticle->GetParticleName()
00474 << " and process " << p->GetProcessName() << G4endl;
00475 }
00476
00477 if(0 == run && startInitialisation) {
00478 emConfigurator->Clear();
00479 firstParticle = aParticle;
00480 }
00481 if(startInitialisation && atomDeexcitation) {
00482 atomDeexcitation->InitialiseAtomicDeexcitation();
00483 }
00484 startInitialisation = false;
00485
00486
00487 if ( aParticle == firstParticle ) {
00488 all_tables_are_built = true;
00489
00490 if(1 < verbose) {
00491 G4cout << "### G4LossTableManager start initilisation for first particle "
00492 << firstParticle->GetParticleName()
00493 << G4endl;
00494 }
00495 for (G4int i=0; i<n_loss; ++i) {
00496 G4VEnergyLossProcess* el = loss_vector[i];
00497
00498 if(el) {
00499 const G4ProcessManager* pm = el->GetProcessManager();
00500 isActive[i] = false;
00501 if(pm) { isActive[i] = pm->GetProcessActivation(el); }
00502 if(0 == run) { base_part_vector[i] = el->BaseParticle(); }
00503 tables_are_built[i] = false;
00504 all_tables_are_built= false;
00505 if(!isActive[i]) { el->SetIonisation(false); }
00506
00507 if(1 < verbose) {
00508 G4cout << i <<". "<< el->GetProcessName();
00509 if(el->Particle()) {
00510 G4cout << " for " << el->Particle()->GetParticleName();
00511 }
00512 G4cout << " active= " << isActive[i]
00513 << " table= " << tables_are_built[i]
00514 << " isIonisation= " << el->IsIonisationProcess();
00515 if(base_part_vector[i]) {
00516 G4cout << " base particle " << base_part_vector[i]->GetParticleName();
00517 }
00518 G4cout << G4endl;
00519 }
00520 } else {
00521 tables_are_built[i] = true;
00522 part_vector[i] = 0;
00523 }
00524 }
00525 ++run;
00526 currentParticle = 0;
00527 }
00528
00529
00530 SetParameters(aParticle, p);
00531
00532 if (all_tables_are_built) { return; }
00533
00534
00535 all_tables_are_built = true;
00536
00537 for(G4int i=0; i<n_loss; ++i) {
00538 if(p == loss_vector[i] && !tables_are_built[i] && !base_part_vector[i]) {
00539 const G4ParticleDefinition* curr_part = part_vector[i];
00540 if(1 < verbose) {
00541 G4cout << "### BuildPhysicsTable for " << p->GetProcessName()
00542 << " and " << curr_part->GetParticleName()
00543 << " start BuildTable " << G4endl;
00544 }
00545 G4VEnergyLossProcess* curr_proc = BuildTables(curr_part);
00546 if(curr_proc) { CopyTables(curr_part, curr_proc); }
00547 }
00548 if ( !tables_are_built[i] ) { all_tables_are_built = false; }
00549 }
00550
00551 if(1 < verbose) {
00552 G4cout << "### G4LossTableManager::BuildDEDXTable end: "
00553 << "all_tables_are_built= " << all_tables_are_built
00554 << G4endl;
00555
00556 if(all_tables_are_built) {
00557 G4cout << "### All dEdx and Range tables are built #####" << G4endl;
00558 }
00559 }
00560 }
00561
00562
00563
00564 void G4LossTableManager::CopyTables(const G4ParticleDefinition* part,
00565 G4VEnergyLossProcess* base_proc)
00566 {
00567 for (G4int j=0; j<n_loss; ++j) {
00568
00569 G4VEnergyLossProcess* proc = loss_vector[j];
00570
00571
00572
00573 if (!tables_are_built[j] && part == base_part_vector[j]) {
00574 tables_are_built[j] = true;
00575 proc->SetDEDXTable(base_proc->DEDXTable(),fRestricted);
00576 proc->SetDEDXTable(base_proc->DEDXTableForSubsec(),fSubRestricted);
00577 proc->SetDEDXTable(base_proc->DEDXunRestrictedTable(),fTotal);
00578 proc->SetCSDARangeTable(base_proc->CSDARangeTable());
00579 proc->SetRangeTableForLoss(base_proc->RangeTableForLoss());
00580 proc->SetInverseRangeTable(base_proc->InverseRangeTable());
00581 proc->SetLambdaTable(base_proc->LambdaTable());
00582 proc->SetSubLambdaTable(base_proc->SubLambdaTable());
00583 proc->SetIonisation(base_proc->IsIonisationProcess());
00584 loss_map[part_vector[j]] = proc;
00585 if (1 < verbose) {
00586 G4cout << "For " << proc->GetProcessName()
00587 << " for " << part_vector[j]->GetParticleName()
00588 << " base_part= " << part->GetParticleName()
00589 << " tables are assigned "
00590 << G4endl;
00591 }
00592 }
00593
00594 if (theElectron == part && theElectron == proc->SecondaryParticle() ) {
00595 proc->SetSecondaryRangeTable(base_proc->RangeTableForLoss());
00596 }
00597 }
00598 }
00599
00600
00601
00602 G4VEnergyLossProcess* G4LossTableManager::BuildTables(
00603 const G4ParticleDefinition* aParticle)
00604 {
00605 if(1 < verbose) {
00606 G4cout << "G4LossTableManager::BuildTables() for "
00607 << aParticle->GetParticleName() << G4endl;
00608 }
00609
00610 std::vector<G4PhysicsTable*> t_list;
00611 std::vector<G4VEnergyLossProcess*> loss_list;
00612 loss_list.clear();
00613 G4VEnergyLossProcess* em = 0;
00614 G4VEnergyLossProcess* p = 0;
00615 G4int iem = 0;
00616 G4PhysicsTable* dedx = 0;
00617 G4int i;
00618
00619 for (i=0; i<n_loss; ++i) {
00620 p = loss_vector[i];
00621 if (p && aParticle == part_vector[i] && !tables_are_built[i]) {
00622 if ((p->IsIonisationProcess() && isActive[i]) ||
00623 !em || (em && !isActive[iem]) ) {
00624 em = p;
00625 iem= i;
00626 }
00627 dedx = p->BuildDEDXTable(fRestricted);
00628
00629
00630 p->SetDEDXTable(dedx,fRestricted);
00631 t_list.push_back(dedx);
00632 loss_list.push_back(p);
00633 tables_are_built[i] = true;
00634 }
00635 }
00636
00637 G4int n_dedx = t_list.size();
00638 if (0 == n_dedx || !em) {
00639 G4cout << "G4LossTableManager WARNING: no DEDX processes for "
00640 << aParticle->GetParticleName() << G4endl;
00641 return 0;
00642 }
00643 G4int nSubRegions = em->NumberOfSubCutoffRegions();
00644
00645 if (1 < verbose) {
00646 G4cout << "G4LossTableManager::BuildTables() start to build range tables"
00647 << " and the sum of " << n_dedx << " processes"
00648 << " iem= " << iem << " em= " << em->GetProcessName()
00649 << " buildCSDARange= " << buildCSDARange
00650 << " nSubRegions= " << nSubRegions
00651 << G4endl;
00652 }
00653
00654 dedx = em->IonisationTable();
00655 if (1 < n_dedx) {
00656 em->SetDEDXTable(dedx, fIsIonisation);
00657 dedx = 0;
00658 dedx = G4PhysicsTableHelper::PreparePhysicsTable(dedx);
00659 tableBuilder->BuildDEDXTable(dedx, t_list);
00660 em->SetDEDXTable(dedx, fRestricted);
00661 }
00662 dedx_vector[iem] = dedx;
00663
00664 G4PhysicsTable* range = em->RangeTableForLoss();
00665 if(!range) range = G4PhysicsTableHelper::PreparePhysicsTable(range);
00666 range_vector[iem] = range;
00667
00668 G4PhysicsTable* invrange = em->InverseRangeTable();
00669 if(!invrange) invrange = G4PhysicsTableHelper::PreparePhysicsTable(invrange);
00670 inv_range_vector[iem] = invrange;
00671
00672 G4bool flag = em->IsIonisationProcess();
00673 tableBuilder->BuildRangeTable(dedx, range, flag);
00674 tableBuilder->BuildInverseRangeTable(range, invrange, flag);
00675
00676
00677
00678 em->SetRangeTableForLoss(range);
00679 em->SetInverseRangeTable(invrange);
00680
00681
00682
00683 std::vector<G4PhysicsTable*> listSub;
00684 std::vector<G4PhysicsTable*> listCSDA;
00685
00686 for (i=0; i<n_dedx; ++i) {
00687 p = loss_list[i];
00688 p->SetIonisation(false);
00689 p->SetLambdaTable(p->BuildLambdaTable(fRestricted));
00690 if (0 < nSubRegions) {
00691 dedx = p->BuildDEDXTable(fSubRestricted);
00692 p->SetDEDXTable(dedx,fSubRestricted);
00693 listSub.push_back(dedx);
00694 p->SetSubLambdaTable(p->BuildLambdaTable(fSubRestricted));
00695 if(p != em) em->AddCollaborativeProcess(p);
00696 }
00697 if(buildCSDARange) {
00698 dedx = p->BuildDEDXTable(fTotal);
00699 p->SetDEDXTable(dedx,fTotal);
00700 listCSDA.push_back(dedx);
00701 }
00702 }
00703
00704 if (0 < nSubRegions) {
00705 G4PhysicsTable* dedxSub = em->IonisationTableForSubsec();
00706 if (1 < listSub.size()) {
00707 em->SetDEDXTable(dedxSub, fIsSubIonisation);
00708 dedxSub = 0;
00709 dedxSub = G4PhysicsTableHelper::PreparePhysicsTable(dedxSub);
00710 tableBuilder->BuildDEDXTable(dedxSub, listSub);
00711 em->SetDEDXTable(dedxSub, fSubRestricted);
00712 }
00713 }
00714 if(buildCSDARange) {
00715 G4PhysicsTable* dedxCSDA = em->DEDXunRestrictedTable();
00716 if (1 < n_dedx) {
00717 dedxCSDA = 0;
00718 dedxCSDA = G4PhysicsTableHelper::PreparePhysicsTable(dedxCSDA);
00719 tableBuilder->BuildDEDXTable(dedxCSDA, listCSDA);
00720 em->SetDEDXTable(dedxCSDA,fTotal);
00721 }
00722 G4PhysicsTable* rCSDA = em->CSDARangeTable();
00723 if(!rCSDA) rCSDA = G4PhysicsTableHelper::PreparePhysicsTable(rCSDA);
00724 tableBuilder->BuildRangeTable(dedxCSDA, rCSDA, flag);
00725 em->SetCSDARangeTable(rCSDA);
00726 }
00727
00728 em->SetIonisation(true);
00729 loss_map[aParticle] = em;
00730
00731 if (1 < verbose) {
00732 G4cout << "G4LossTableManager::BuildTables: Tables are built for "
00733 << aParticle->GetParticleName()
00734 << "; ionisation process: " << em->GetProcessName()
00735 << G4endl;
00736 }
00737 return em;
00738 }
00739
00740
00741
00742 G4EnergyLossMessenger* G4LossTableManager::GetMessenger()
00743 {
00744 return theMessenger;
00745 }
00746
00747
00748
00749 void G4LossTableManager::ParticleHaveNoLoss(
00750 const G4ParticleDefinition* aParticle)
00751 {
00752 G4ExceptionDescription ed;
00753 ed << "Energy loss process not found for " << aParticle->GetParticleName()
00754 << " !" << G4endl;
00755 G4Exception("G4LossTableManager::ParticleHaveNoLoss", "em0001",
00756 FatalException, ed);
00757
00758 }
00759
00760
00761
00762 G4bool G4LossTableManager::BuildCSDARange() const
00763 {
00764 return buildCSDARange;
00765 }
00766
00767
00768
00769 void G4LossTableManager::SetLossFluctuations(G4bool val)
00770 {
00771 lossFluctuationFlag = val;
00772 for(G4int i=0; i<n_loss; ++i) {
00773 if(loss_vector[i]) { loss_vector[i]->SetLossFluctuations(val); }
00774 }
00775 }
00776
00777
00778
00779 void G4LossTableManager::SetSubCutoff(G4bool val, const G4Region* r)
00780 {
00781 subCutoffFlag = val;
00782 for(G4int i=0; i<n_loss; ++i) {
00783 if(loss_vector[i]) { loss_vector[i]->ActivateSubCutoff(val, r); }
00784 }
00785 }
00786
00787
00788
00789 void G4LossTableManager::SetIntegral(G4bool val)
00790 {
00791 integral = val;
00792 integralActive = true;
00793 for(G4int i=0; i<n_loss; ++i) {
00794 if(loss_vector[i]) { loss_vector[i]->SetIntegral(val); }
00795 }
00796 size_t emp = emp_vector.size();
00797 for (size_t k=0; k<emp; ++k) {
00798 if(emp_vector[k]) { emp_vector[k]->SetIntegral(val); }
00799 }
00800 }
00801
00802
00803
00804 void G4LossTableManager::SetMinSubRange(G4double val)
00805 {
00806 minSubRange = val;
00807 for(G4int i=0; i<n_loss; ++i) {
00808 if(loss_vector[i]) { loss_vector[i]->SetMinSubRange(val); }
00809 }
00810 }
00811
00812
00813
00814 void G4LossTableManager::SetRandomStep(G4bool val)
00815 {
00816 rndmStepFlag = val;
00817 for(G4int i=0; i<n_loss; ++i) {
00818 if(loss_vector[i]) { loss_vector[i]->SetRandomStep(val); }
00819 }
00820 }
00821
00822
00823
00824 void G4LossTableManager::SetMinEnergy(G4double val)
00825 {
00826 minEnergyActive = true;
00827 minKinEnergy = val;
00828 for(G4int i=0; i<n_loss; ++i) {
00829 if(loss_vector[i]) { loss_vector[i]->SetMinKinEnergy(val); }
00830 }
00831 size_t emp = emp_vector.size();
00832 for (size_t k=0; k<emp; ++k) {
00833 if(emp_vector[k]) { emp_vector[k]->SetMinKinEnergy(val); }
00834 }
00835 }
00836
00837
00838
00839 void G4LossTableManager::SetMaxEnergy(G4double val)
00840 {
00841 maxEnergyActive = true;
00842 maxKinEnergy = val;
00843 for(G4int i=0; i<n_loss; ++i) {
00844 if(loss_vector[i]) { loss_vector[i]->SetMaxKinEnergy(val); }
00845 }
00846 size_t emp = emp_vector.size();
00847 for (size_t k=0; k<emp; ++k) {
00848 if(emp_vector[k]) { emp_vector[k]->SetMaxKinEnergy(val); }
00849 }
00850 }
00851
00852
00853
00854 void G4LossTableManager::SetMaxEnergyForCSDARange(G4double val)
00855 {
00856 for(G4int i=0; i<n_loss; ++i) {
00857 if(loss_vector[i]) { loss_vector[i]->SetMaxKinEnergyForCSDARange(val); }
00858 }
00859 }
00860
00861
00862
00863 void G4LossTableManager::SetMaxEnergyForMuons(G4double val)
00864 {
00865 maxEnergyForMuonsActive = true;
00866 maxKinEnergyForMuons = val;
00867 }
00868
00869
00870
00871 void G4LossTableManager::SetDEDXBinning(G4int val)
00872 {
00873 for(G4int i=0; i<n_loss; ++i) {
00874 if(loss_vector[i]) { loss_vector[i]->SetDEDXBinning(val); }
00875 }
00876 }
00877
00878
00879
00880 void G4LossTableManager::SetDEDXBinningForCSDARange(G4int val)
00881 {
00882 for(G4int i=0; i<n_loss; ++i) {
00883 if(loss_vector[i]) { loss_vector[i]->SetDEDXBinningForCSDARange(val); }
00884 }
00885 }
00886
00887
00888
00889 void G4LossTableManager::SetLambdaBinning(G4int val)
00890 {
00891 G4int n = val/G4int(std::log10(maxKinEnergy/minKinEnergy) + 0.5);
00892 if(n < 5) {
00893 G4cout << "G4LossTableManager::SetLambdaBinning WARNING "
00894 << "too small number of bins " << val << " ignored"
00895 << G4endl;
00896 return;
00897 }
00898 nbinsLambda = val;
00899 nbinsPerDecade = n;
00900 size_t emp = emp_vector.size();
00901 for (size_t k=0; k<emp; ++k) {
00902 if(emp_vector[k]) { emp_vector[k]->SetLambdaBinning(val); }
00903 }
00904 }
00905
00906
00907
00908 G4int G4LossTableManager::GetNumberOfBinsPerDecade() const
00909 {
00910 return nbinsPerDecade;
00911 }
00912
00913
00914
00915 void G4LossTableManager::SetVerbose(G4int val)
00916 {
00917 verbose = val;
00918 for(G4int i=0; i<n_loss; ++i) {
00919 if(loss_vector[i]) { loss_vector[i]->SetVerboseLevel(val); }
00920 }
00921 size_t msc = msc_vector.size();
00922 for (size_t j=0; j<msc; ++j) {
00923 if(msc_vector[j]) { msc_vector[j]->SetVerboseLevel(val); }
00924 }
00925 size_t emp = emp_vector.size();
00926 for (size_t k=0; k<emp; ++k) {
00927 if(emp_vector[k]) { emp_vector[k]->SetVerboseLevel(val); }
00928 }
00929 emConfigurator->SetVerbose(val);
00930
00931
00932 emSaturation->SetVerbose(val);
00933 emElectronIonPair->SetVerbose(val);
00934 if(atomDeexcitation) { atomDeexcitation->SetVerboseLevel(val); }
00935 }
00936
00937
00938
00939 void G4LossTableManager::SetStepFunction(G4double v1, G4double v2)
00940 {
00941 stepFunctionActive = true;
00942 maxRangeVariation = v1;
00943 maxFinalStep = v2;
00944 for(G4int i=0; i<n_loss; ++i) {
00945 if(loss_vector[i]) { loss_vector[i]->SetStepFunction(v1, v2); }
00946 }
00947 }
00948
00949
00950
00951 void G4LossTableManager::SetLinearLossLimit(G4double val)
00952 {
00953 for(G4int i=0; i<n_loss; ++i) {
00954 if(loss_vector[i]) { loss_vector[i]->SetLinearLossLimit(val); }
00955 }
00956 }
00957
00958
00959
00960 void G4LossTableManager::SetBuildCSDARange(G4bool val)
00961 {
00962 buildCSDARange = val;
00963 }
00964
00965
00966
00967 void
00968 G4LossTableManager::SetParameters(const G4ParticleDefinition* aParticle,
00969 G4VEnergyLossProcess* p)
00970 {
00971 if(stepFunctionActive) { p->SetStepFunction(maxRangeVariation, maxFinalStep); }
00972 if(integralActive) { p->SetIntegral(integral); }
00973 if(minEnergyActive) { p->SetMinKinEnergy(minKinEnergy); }
00974 if(maxEnergyActive) { p->SetMaxKinEnergy(maxKinEnergy); }
00975 p->SetVerboseLevel(verbose);
00976 if(maxEnergyForMuonsActive) {
00977 G4double dm = std::abs(aParticle->GetPDGMass() - 105.7*MeV);
00978 if(dm < 5.*MeV) { p->SetMaxKinEnergy(maxKinEnergyForMuons); }
00979 }
00980 }
00981
00982
00983
00984 const std::vector<G4VEnergyLossProcess*>&
00985 G4LossTableManager::GetEnergyLossProcessVector()
00986 {
00987 return loss_vector;
00988 }
00989
00990
00991
00992 const std::vector<G4VEmProcess*>& G4LossTableManager::GetEmProcessVector()
00993 {
00994 return emp_vector;
00995 }
00996
00997
00998
00999 const std::vector<G4VMultipleScattering*>&
01000 G4LossTableManager::GetMultipleScatteringVector()
01001 {
01002 return msc_vector;
01003 }
01004
01005
01006
01007 void G4LossTableManager::SetLPMFlag(G4bool val)
01008 {
01009 flagLPM = val;
01010 }
01011
01012
01013
01014 G4bool G4LossTableManager::LPMFlag() const
01015 {
01016 return flagLPM;
01017 }
01018
01019
01020
01021 void G4LossTableManager::SetSplineFlag(G4bool val)
01022 {
01023 splineFlag = val;
01024 tableBuilder->SetSplineFlag(splineFlag);
01025 }
01026
01027
01028
01029 G4bool G4LossTableManager::SplineFlag() const
01030 {
01031 return splineFlag;
01032 }
01033
01034
01035
01036 void G4LossTableManager::SetBremsstrahlungTh(G4double val)
01037 {
01038 bremsTh = val;
01039 }
01040
01041
01042
01043 G4double G4LossTableManager::BremsstrahlungTh() const
01044 {
01045 return bremsTh;
01046 }
01047
01048
01049
01050 void G4LossTableManager::SetFactorForAngleLimit(G4double val)
01051 {
01052 if(val > 0.0) { factorForAngleLimit = val; }
01053 }
01054
01055
01056
01057 G4double G4LossTableManager::FactorForAngleLimit() const
01058 {
01059 return factorForAngleLimit;
01060 }
01061
01062
01063
01064 G4double G4LossTableManager::MinKinEnergy() const
01065 {
01066 return minKinEnergy;
01067 }
01068
01069
01070
01071 G4double G4LossTableManager::MaxKinEnergy() const
01072 {
01073 return maxKinEnergy;
01074 }
01075
01076
01077
01078 G4EmCorrections* G4LossTableManager::EmCorrections()
01079 {
01080 return emCorrections;
01081 }
01082
01083
01084
01085 G4EmSaturation* G4LossTableManager::EmSaturation()
01086 {
01087 return emSaturation;
01088 }
01089
01090
01091
01092 G4EmConfigurator* G4LossTableManager::EmConfigurator()
01093 {
01094 return emConfigurator;
01095 }
01096
01097
01098
01099 G4ElectronIonPair* G4LossTableManager::ElectronIonPair()
01100 {
01101 return emElectronIonPair;
01102 }
01103
01104
01105
01106 G4VAtomDeexcitation* G4LossTableManager::AtomDeexcitation()
01107 {
01108 return atomDeexcitation;
01109 }
01110
01111
01112
01113 G4LossTableBuilder* G4LossTableManager::GetTableBuilder()
01114 {
01115 return tableBuilder;
01116 }
01117
01118
01119
01120 void G4LossTableManager::SetAtomDeexcitation(G4VAtomDeexcitation* p)
01121 {
01122 atomDeexcitation = p;
01123 }
01124
01125
01126
01127 G4VEnergyLossProcess*
01128 G4LossTableManager::GetEnergyLossProcess(const G4ParticleDefinition *aParticle)
01129 {
01130 if(aParticle != currentParticle) {
01131 currentParticle = aParticle;
01132 std::map<PD,G4VEnergyLossProcess*,std::less<PD> >::const_iterator pos;
01133 if ((pos = loss_map.find(aParticle)) != loss_map.end()) {
01134 currentLoss = (*pos).second;
01135 } else {
01136 currentLoss = 0;
01137
01138 }
01139 }
01140 return currentLoss;
01141 }
01142
01143
01144
01145 G4double G4LossTableManager::GetDEDX(const G4ParticleDefinition *aParticle,
01146 G4double kineticEnergy,
01147 const G4MaterialCutsCouple *couple)
01148 {
01149 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01150 G4double x = 0.0;
01151 if(currentLoss) { x = currentLoss->GetDEDX(kineticEnergy, couple); }
01152 else { x = G4EnergyLossTables::GetDEDX(currentParticle,
01153 kineticEnergy,couple,false); }
01154 return x;
01155 }
01156
01157
01158
01159 G4double G4LossTableManager::GetSubDEDX(const G4ParticleDefinition *aParticle,
01160 G4double kineticEnergy,
01161 const G4MaterialCutsCouple *couple)
01162 {
01163 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01164 G4double x = 0.0;
01165 if(currentLoss) { x = currentLoss->GetDEDXForSubsec(kineticEnergy, couple); }
01166 return x;
01167 }
01168
01169
01170
01171 G4double G4LossTableManager::GetCSDARange(const G4ParticleDefinition *aParticle,
01172 G4double kineticEnergy,
01173 const G4MaterialCutsCouple *couple)
01174 {
01175 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01176 G4double x = DBL_MAX;
01177 if(currentLoss) { x = currentLoss->GetCSDARange(kineticEnergy, couple); }
01178 return x;
01179 }
01180
01181
01182
01183 G4double
01184 G4LossTableManager::GetRangeFromRestricteDEDX(const G4ParticleDefinition *aParticle,
01185 G4double kineticEnergy,
01186 const G4MaterialCutsCouple *couple)
01187 {
01188 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01189 G4double x;
01190 if(currentLoss) { x = currentLoss->GetRangeForLoss(kineticEnergy, couple); }
01191 else { x = G4EnergyLossTables::GetRange(currentParticle,kineticEnergy,
01192 couple,false); }
01193 return x;
01194 }
01195
01196
01197
01198 G4double G4LossTableManager::GetRange(const G4ParticleDefinition *aParticle,
01199 G4double kineticEnergy,
01200 const G4MaterialCutsCouple *couple)
01201 {
01202 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01203 G4double x;
01204 if(currentLoss) { x = currentLoss->GetRange(kineticEnergy, couple); }
01205 else { x = G4EnergyLossTables::GetRange(currentParticle,kineticEnergy,
01206 couple,false); }
01207 return x;
01208 }
01209
01210
01211
01212 G4double G4LossTableManager::GetEnergy(const G4ParticleDefinition *aParticle,
01213 G4double range,
01214 const G4MaterialCutsCouple *couple)
01215 {
01216 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01217 G4double x;
01218 if(currentLoss) { x = currentLoss->GetKineticEnergy(range, couple); }
01219 else { x = G4EnergyLossTables::GetPreciseEnergyFromRange(currentParticle,range,
01220 couple,false); }
01221 return x;
01222 }
01223
01224
01225
01226 G4double G4LossTableManager::GetDEDXDispersion(const G4MaterialCutsCouple *couple,
01227 const G4DynamicParticle* dp,
01228 G4double& length)
01229 {
01230 const G4ParticleDefinition* aParticle = dp->GetParticleDefinition();
01231 if(aParticle != currentParticle) { GetEnergyLossProcess(aParticle); }
01232 G4double x = 0.0;
01233 if(currentLoss) { currentLoss->GetDEDXDispersion(couple, dp, length); }
01234 return x;
01235 }
01236
01237