Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/utils/src/G4LossTableManager.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  1 //
  2 // ********************************************************************
  3 // * License and Disclaimer                                           *
  4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.                             *
 10 // *                                                                  *
 11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                                                  *
 18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // ********************************************************************
 25 //
 26 // -------------------------------------------------------------------
 27 //
 28 // GEANT4 Class file
 29 //
 30 //
 31 // File name:     G4LossTableManager
 32 //
 33 // Author:        Vladimir Ivanchenko
 34 //
 35 // Creation date: 03.01.2002
 36 //
 37 // Modifications: by V.Ivanchenko
 38 //
 39 //
 40 // Class Description:
 41 //
 42 // -------------------------------------------------------------------
 43 //
 44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 46 
 47 #include "G4LossTableManager.hh"
 48 #include "G4SystemOfUnits.hh"
 49 
 50 #include "G4VMultipleScattering.hh"
 51 #include "G4VEmProcess.hh"
 52 
 53 #include "G4EmParameters.hh"
 54 #include "G4EmSaturation.hh"
 55 #include "G4EmConfigurator.hh"
 56 #include "G4ElectronIonPair.hh"
 57 #include "G4NIELCalculator.hh"
 58 #include "G4EmCorrections.hh"
 59 #include "G4LossTableBuilder.hh"
 60 #include "G4VAtomDeexcitation.hh"
 61 #include "G4VSubCutProducer.hh"
 62 
 63 #include "G4PhysicsTable.hh"
 64 #include "G4ParticleDefinition.hh"
 65 #include "G4MaterialCutsCouple.hh"
 66 #include "G4ProcessManager.hh"
 67 #include "G4Electron.hh"
 68 #include "G4Proton.hh"
 69 #include "G4ProductionCutsTable.hh"
 70 #include "G4PhysicsTableHelper.hh"
 71 #include "G4EmTableType.hh"
 72 #include "G4Region.hh"
 73 #include "G4PhysicalConstants.hh"
 74 
 75 #include "G4Gamma.hh"
 76 #include "G4Positron.hh"
 77 #include "G4OpticalPhoton.hh"
 78 #include "G4Neutron.hh"
 79 #include "G4MuonPlus.hh"
 80 #include "G4MuonMinus.hh"
 81 #include "G4GenericIon.hh"
 82 
 83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 84 
 85 static std::once_flag applyOnce;
 86 G4ThreadLocal G4LossTableManager* G4LossTableManager::instance = nullptr;
 87 
 88 G4LossTableManager* G4LossTableManager::Instance()
 89 {
 90   if(nullptr == instance) {
 91     static G4ThreadLocalSingleton<G4LossTableManager> inst;
 92     instance = inst.Instance();
 93   }
 94   return instance;
 95 }
 96 
 97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 98 
 99 G4LossTableManager::~G4LossTableManager()
100 {
101   for (auto const & p : loss_vector) { delete p; }
102   for (auto const & p : msc_vector) { delete p; }
103   for (auto const & p : emp_vector) { delete p; }
104   for (auto const & p : p_vector) { delete p; }
105 
106   std::size_t mod = mod_vector.size();
107   std::size_t fmod = fmod_vector.size();
108   for (std::size_t a=0; a<mod; ++a) {
109     if( nullptr != mod_vector[a] ) { 
110       for (std::size_t b=0; b<fmod; ++b) {
111         if((G4VEmModel*)(fmod_vector[b]) == mod_vector[a]) {
112           fmod_vector[b] = nullptr;
113         }
114       }
115       delete mod_vector[a]; 
116       mod_vector[a] = nullptr;
117     }
118   }
119   for (auto const & p : fmod_vector) { delete p; }
120 
121   Clear();
122   delete tableBuilder;
123   delete emCorrections;
124   delete emConfigurator;
125   delete emElectronIonPair;
126   delete nielCalculator;
127   delete atomDeexcitation;
128   delete subcutProducer;
129 }
130 
131 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
132 
133 G4LossTableManager::G4LossTableManager()
134 {
135   theParameters = G4EmParameters::Instance();
136   theElectron = G4Electron::Electron();
137   
138   // only one thread is the master 
139   std::call_once(applyOnce, [this]() { isMaster = true; });
140   verbose = isMaster ? theParameters->Verbose() : theParameters->WorkerVerbose();
141 
142   tableBuilder = new G4LossTableBuilder(isMaster);
143   emCorrections = new G4EmCorrections(verbose);
144 
145   std::size_t n = 70;
146   loss_vector.reserve(n);
147   part_vector.reserve(n);
148   base_part_vector.reserve(n);
149   dedx_vector.reserve(n);
150   range_vector.reserve(n);
151   inv_range_vector.reserve(n);
152   tables_are_built.reserve(n);
153   isActive.reserve(n);
154   msc_vector.reserve(10);
155   emp_vector.reserve(16);
156   mod_vector.reserve(150);
157   fmod_vector.reserve(60);
158 }
159 
160 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
161 
162 void G4LossTableManager::Clear()
163 {
164   all_tables_are_built = false;
165   currentLoss = nullptr;
166   currentParticle = nullptr;
167   if(n_loss) {
168     dedx_vector.clear();
169     range_vector.clear();
170     inv_range_vector.clear();
171     loss_map.clear();
172     loss_vector.clear();
173     part_vector.clear();
174     base_part_vector.clear();
175     tables_are_built.clear();
176     isActive.clear();
177     n_loss = 0;
178   }
179 }
180 
181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
182 
183 void G4LossTableManager::Register(G4VEnergyLossProcess* p)
184 {
185   if (nullptr == p) { return; }
186   for (G4int i=0; i<n_loss; ++i) {
187     if(loss_vector[i] == p) { return; }
188   }
189   if(verbose > 1) {
190     G4cout << "G4LossTableManager::Register G4VEnergyLossProcess : " 
191            << p->GetProcessName() << "  idx= " << n_loss << G4endl;
192   }
193   ++n_loss;
194   loss_vector.push_back(p);
195   part_vector.push_back(nullptr);
196   base_part_vector.push_back(nullptr);
197   dedx_vector.push_back(nullptr);
198   range_vector.push_back(nullptr);
199   inv_range_vector.push_back(nullptr);
200   tables_are_built.push_back(false);
201   isActive.push_back(true);
202   all_tables_are_built = false;
203 }
204 
205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
206 
207 void G4LossTableManager::ResetParameters()
208 {
209   // initialisation once per run
210   if (!resetParam) { return; }
211   resetParam = false;
212   startInitialisation = true;
213   verbose = theParameters->Verbose();
214   if(!isMaster) {
215     verbose = theParameters->WorkerVerbose();
216   } else {
217     if(verbose > 0) { theParameters->Dump(); }
218   }
219 
220   tableBuilder->InitialiseBaseMaterials(); 
221   if (nullptr != nielCalculator) { nielCalculator->Initialise(); } 
222 
223   emCorrections->SetVerbose(verbose); 
224   if(nullptr != emConfigurator) { emConfigurator->SetVerbose(verbose); };
225   if(nullptr != emElectronIonPair) { emElectronIonPair->SetVerbose(verbose); };
226   if(nullptr != atomDeexcitation) {
227     atomDeexcitation->SetVerboseLevel(verbose);
228     atomDeexcitation->InitialiseAtomicDeexcitation();
229   }
230   if (1 < verbose) {
231     G4cout << "====== G4LossTableManager::ResetParameters " 
232            << " Nloss=" << loss_vector.size()
233      << " run=" << run << " master=" << isMaster
234      << G4endl;
235   }
236 }
237 
238 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
239 
240 void G4LossTableManager::DeRegister(G4VEnergyLossProcess* p)
241 {
242   if (nullptr == p) { return; }
243   for (G4int i=0; i<n_loss; ++i) {
244     if(loss_vector[i] == p) { 
245       loss_vector[i] = nullptr;
246       break;
247     }
248   }
249 }
250 
251 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
252 
253 void G4LossTableManager::Register(G4VMultipleScattering* p)
254 {
255   if (nullptr == p) { return; }
256   std::size_t n = msc_vector.size();
257   for (std::size_t i=0; i<n; ++i) {
258     if(msc_vector[i] == p) { return; }
259   }
260   if(verbose > 1) {
261     G4cout << "G4LossTableManager::Register G4VMultipleScattering : " 
262            << p->GetProcessName() << "  idx= " << msc_vector.size() << G4endl;
263   }
264   msc_vector.push_back(p);
265 }
266 
267 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
268 
269 void G4LossTableManager::DeRegister(G4VMultipleScattering* p)
270 {
271   if (nullptr == p) { return; }
272   std::size_t msc = msc_vector.size();
273   for (std::size_t i=0; i<msc; ++i) {
274     if(msc_vector[i] == p) { 
275       msc_vector[i] = nullptr;
276       break;
277     }
278   }
279 }
280 
281 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
282 
283 void G4LossTableManager::Register(G4VEmProcess* p)
284 {
285   if (nullptr == p) { return; }
286   std::size_t n = emp_vector.size();
287   for (std::size_t i=0; i<n; ++i) {
288     if(emp_vector[i] == p) { return; }
289   }
290   if(verbose > 1) {
291     G4cout << "G4LossTableManager::Register G4VEmProcess : " 
292            << p->GetProcessName() << "  idx= " << emp_vector.size() << G4endl;
293   }
294   emp_vector.push_back(p);
295 }
296 
297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
298 
299 void G4LossTableManager::DeRegister(G4VEmProcess* p)
300 {
301   if (nullptr == p) { return; }
302   std::size_t emp = emp_vector.size();
303   for (std::size_t i=0; i<emp; ++i) {
304     if(emp_vector[i] == p) { 
305       emp_vector[i] = nullptr; 
306       break;
307     }
308   }
309 }
310 
311 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
312 
313 void G4LossTableManager::Register(G4VProcess* p)
314 {
315   if (nullptr == p) { return; }
316   std::size_t n = p_vector.size();
317   for (std::size_t i=0; i<n; ++i) {
318     if(p_vector[i] == p) { return; }
319   }
320   if(verbose > 1) {
321     G4cout << "G4LossTableManager::Register G4VProcess : " 
322            << p->GetProcessName() << "  idx= " << p_vector.size() << G4endl;
323   }
324   p_vector.push_back(p);
325 }
326 
327 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
328 
329 void G4LossTableManager::DeRegister(G4VProcess* p)
330 {
331   if (nullptr == p) { return; }
332   std::size_t emp = p_vector.size();
333   for (std::size_t i=0; i<emp; ++i) {
334     if(p_vector[i] == p) { 
335       p_vector[i] = nullptr;
336       break;
337     }
338   }
339 }
340 
341 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
342 
343 void G4LossTableManager::Register(G4VEmModel* p)
344 {
345   mod_vector.push_back(p);
346   if(verbose > 1) {
347     G4cout << "G4LossTableManager::Register G4VEmModel : " 
348            << p->GetName() << "  " << p << "  " << mod_vector.size() << G4endl;
349   }
350 }
351 
352 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
353 
354 void G4LossTableManager::DeRegister(G4VEmModel* p)
355 {
356   //G4cout << "G4LossTableManager::DeRegister G4VEmModel : " << p << G4endl;
357   std::size_t n = mod_vector.size();
358   for (std::size_t i=0; i<n; ++i) {
359     if(mod_vector[i] == p) { 
360       mod_vector[i] = nullptr; 
361       break;
362     }
363   }
364 }
365 
366 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
367 
368 void G4LossTableManager::Register(G4VEmFluctuationModel* p)
369 {
370   fmod_vector.push_back(p);
371   if(verbose > 1) {
372     G4cout << "G4LossTableManager::Register G4VEmFluctuationModel : " 
373            << p->GetName() << "  " << fmod_vector.size() << G4endl;
374   }
375 }
376 
377 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
378 
379 void G4LossTableManager::DeRegister(G4VEmFluctuationModel* p)
380 {
381   std::size_t n = fmod_vector.size();
382   for (std::size_t i=0; i<n; ++i) {
383     if(fmod_vector[i] == p) { fmod_vector[i] = nullptr; }
384   }
385 }
386 
387 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
388 
389 void G4LossTableManager::RegisterExtraParticle(
390      const G4ParticleDefinition* part,
391      G4VEnergyLossProcess* p)
392 { 
393   if (nullptr == p || nullptr == part) { return; }
394   for (G4int i=0; i<n_loss; ++i) {
395     if(loss_vector[i] == p) { return; }
396   }
397   if(verbose > 1) {
398     G4cout << "G4LossTableManager::RegisterExtraParticle "
399            << part->GetParticleName() << "  G4VEnergyLossProcess : " 
400            << p->GetProcessName() << "  idx= " << n_loss << G4endl;
401   }
402   ++n_loss;
403   loss_vector.push_back(p);
404   part_vector.push_back(part);
405   base_part_vector.push_back(p->BaseParticle());
406   dedx_vector.push_back(nullptr);
407   range_vector.push_back(nullptr);
408   inv_range_vector.push_back(nullptr);
409   tables_are_built.push_back(false);
410   all_tables_are_built = false;
411 }
412 
413 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
414 
415 G4VEnergyLossProcess* 
416 G4LossTableManager::GetEnergyLossProcess(const G4ParticleDefinition *aParticle)
417 {
418   if(aParticle != currentParticle) {
419     currentParticle = aParticle;
420     std::map<PD,G4VEnergyLossProcess*,std::less<PD> >::const_iterator pos;
421     if ((pos = loss_map.find(aParticle)) != loss_map.end()) {
422       currentLoss = (*pos).second;
423     } else {
424       currentLoss = nullptr;
425       if(0.0 != aParticle->GetPDGCharge() && 
426    (pos = loss_map.find(theGenericIon)) != loss_map.end()) {
427   currentLoss = (*pos).second;
428       }
429     }
430   }
431   return currentLoss;
432 }
433 
434 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
435 
436 void 
437 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
438                                         G4VEnergyLossProcess* p)
439 {
440   if (1 < verbose) {
441     G4cout << "G4LossTableManager::PreparePhysicsTable for " 
442            << particle->GetParticleName() 
443            << " and " << p->GetProcessName() << " run= " << run 
444            << "   loss_vector " << loss_vector.size()
445      << " run=" << run << " master=" << isMaster
446      << G4endl;
447   }
448 
449   // start initialisation for the first run
450   if( -1 == run ) {
451     if (nullptr != emConfigurator) { emConfigurator->PrepareModels(particle, p); }
452 
453     // initialise particles for given process
454     for (G4int j=0; j<n_loss; ++j) {
455       if (p == loss_vector[j] && nullptr == part_vector[j]) { 
456         part_vector[j] = particle;
457         if (particle->GetParticleName() == "GenericIon") {
458           theGenericIon = particle;
459         }
460       }
461     }
462   }
463   ResetParameters();
464 }
465 
466 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
467 
468 void 
469 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
470                                         G4VEmProcess* p)
471 {
472   if (1 < verbose) {
473     G4cout << "G4LossTableManager::PreparePhysicsTable for " 
474            << particle->GetParticleName() 
475            << " and " << p->GetProcessName()
476      << " run=" << run << " master=" << isMaster
477      << G4endl;
478   }
479 
480   // start initialisation for the first run
481   if( -1 == run ) {
482     if (nullptr != emConfigurator) { emConfigurator->PrepareModels(particle, p); }
483   }
484 
485   ResetParameters();
486 }
487 
488 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
489 
490 void 
491 G4LossTableManager::PreparePhysicsTable(const G4ParticleDefinition* particle,
492                                         G4VMultipleScattering* p)
493 {
494   if (1 < verbose) {
495     G4cout << "G4LossTableManager::PreparePhysicsTable for " 
496            << particle->GetParticleName() 
497            << " and " << p->GetProcessName()
498      << " run=" << run << " master=" << isMaster
499      << G4endl;
500   }
501 
502   // start initialisation for the first run
503   if ( -1 == run ) {
504     if (nullptr != emConfigurator) { emConfigurator->PrepareModels(particle, p); }
505   } 
506   
507   ResetParameters();
508 }
509 
510 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
511 
512 void 
513 G4LossTableManager::BuildPhysicsTable(const G4ParticleDefinition*)
514 {
515   if(-1 == run && startInitialisation) {
516     if (nullptr != emConfigurator) { emConfigurator->Clear(); }
517   }
518   if (startInitialisation) { resetParam = true; }
519 }
520 
521 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
522 
523 void G4LossTableManager::LocalPhysicsTables(
524      const G4ParticleDefinition* aParticle,
525      G4VEnergyLossProcess* p)
526 {
527   if (1 < verbose) {
528     G4cout << "### G4LossTableManager::LocalPhysicsTable() for "
529            << aParticle->GetParticleName()
530            << " and process " << p->GetProcessName()
531            << G4endl;
532   }
533 
534   if(-1 == run && startInitialisation) {
535     if (nullptr != emConfigurator) { emConfigurator->Clear(); }
536     firstParticle = aParticle; 
537   }
538 
539   if (startInitialisation) {
540     ++run;
541     if (1 < verbose) {
542       G4cout << "===== G4LossTableManager::LocalPhysicsTable() for run "
543              << run << " =====" << G4endl;
544     }
545     currentParticle = nullptr;
546     startInitialisation = false;
547     resetParam = true;
548     for (G4int i=0; i<n_loss; ++i) {
549       if (nullptr != loss_vector[i]) {
550         tables_are_built[i] = false;
551       } else {
552         tables_are_built[i] = true;
553         part_vector[i] = nullptr;
554       }
555     }
556   }
557 
558   all_tables_are_built= true;
559   for (G4int i=0; i<n_loss; ++i) {
560     if(p == loss_vector[i]) {
561       tables_are_built[i] = true;
562       isActive[i] = true;
563       part_vector[i] = p->Particle(); 
564       base_part_vector[i] = p->BaseParticle(); 
565       dedx_vector[i] = p->DEDXTable();
566       range_vector[i] = p->RangeTableForLoss();
567       inv_range_vector[i] = p->InverseRangeTable();
568       if(0 == run && p->IsIonisationProcess()) {
569         loss_map[part_vector[i]] = p;
570       }
571 
572       if(1 < verbose) { 
573         G4cout << i <<".   "<< p->GetProcessName(); 
574         if(part_vector[i]) {
575           G4cout << "  for "  << part_vector[i]->GetParticleName();
576         }
577         G4cout << "  active= " << isActive[i]
578                << "  table= " << tables_are_built[i]
579                << "  isIonisation= " << p->IsIonisationProcess()
580                << G4endl;
581       }
582       break;
583     } else if(!tables_are_built[i]) {
584       all_tables_are_built = false;
585     }
586   }
587 
588   if(1 < verbose) {
589     G4cout << "### G4LossTableManager::LocalPhysicsTable end"
590            << G4endl;
591   }
592   if(all_tables_are_built) { 
593     if(1 < verbose) {
594       G4cout << "%%%%% All dEdx and Range tables for worker are ready for run "
595              << run << " %%%%%" << G4endl;
596     }
597   }
598 }
599 
600 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
601 
602 void G4LossTableManager::BuildPhysicsTable(
603      const G4ParticleDefinition* aParticle,
604      G4VEnergyLossProcess* p)
605 {
606   if(1 < verbose) {
607     G4cout << "### G4LossTableManager::BuildPhysicsTable() for "
608            << aParticle->GetParticleName()
609            << " and process " << p->GetProcessName() << G4endl;
610   }
611   // clear configurator
612   if(-1 == run && startInitialisation) {
613     if( nullptr != emConfigurator) { emConfigurator->Clear(); }
614     firstParticle = aParticle; 
615   }
616   if(startInitialisation) {
617     ++run;
618     resetParam = true;
619     startInitialisation = false;
620     if(1 < verbose) {
621       G4cout << "===== G4LossTableManager::BuildPhysicsTable() for run "
622              << run << " ===== " << atomDeexcitation << G4endl;
623     }
624     currentParticle = nullptr;
625     all_tables_are_built = false;
626 
627     for (G4int i=0; i<n_loss; ++i) {
628       G4VEnergyLossProcess* el = loss_vector[i];
629 
630       if(nullptr != el) {
631         isActive[i] = true;
632         part_vector[i] = el->Particle(); 
633         base_part_vector[i] = el->BaseParticle(); 
634         tables_are_built[i] = false;  
635         if(1 < verbose) { 
636           G4cout << i <<".   "<< el->GetProcessName();
637           if(el->Particle()) {
638             G4cout << "  for "  << el->Particle()->GetParticleName();
639           }
640           G4cout << "  active= " << isActive[i]
641                  << "  table= " << tables_are_built[i]
642                  << "  isIonisation= " << el->IsIonisationProcess();
643           if(base_part_vector[i]) { 
644             G4cout << "  base particle " 
645                    << base_part_vector[i]->GetParticleName();
646           }
647           G4cout << G4endl;
648         }
649       } else {
650         tables_are_built[i] = true;
651         part_vector[i] = nullptr;
652         isActive[i] = false;
653       }
654     }
655   }
656 
657   if (all_tables_are_built) { 
658     theParameters->SetIsPrintedFlag(true);
659     return; 
660   }
661 
662   // Build tables for given particle
663   all_tables_are_built = true;
664 
665   for(G4int i=0; i<n_loss; ++i) {
666     if(p == loss_vector[i] && !tables_are_built[i] && nullptr == base_part_vector[i]) {
667       const G4ParticleDefinition* curr_part = part_vector[i];
668       if(1 < verbose) {
669         G4cout << "### Build Table for " << p->GetProcessName()
670                << " and " << curr_part->GetParticleName() 
671                << "  " << tables_are_built[i] << "  " << base_part_vector[i] 
672                << G4endl;
673       }
674       G4VEnergyLossProcess* curr_proc = BuildTables(curr_part);
675       if(curr_proc) { 
676         CopyTables(curr_part, curr_proc); 
677         if(p == curr_proc && 0 == run && p->IsIonisationProcess()) { 
678           loss_map[aParticle] = p; 
679           //G4cout << "G4LossTableManager::BuildPhysicsTable: " 
680           //     << aParticle->GetParticleName()
681           //         << " added to map " << p <<  G4endl;
682         }
683       }
684     }
685     if ( !tables_are_built[i] ) { all_tables_are_built = false; }
686   }
687   if(1 < verbose) {
688     G4cout << "### G4LossTableManager::BuildPhysicsTable end: "
689            << "all_tables_are_built= " << all_tables_are_built << " "
690            << aParticle->GetParticleName() << " proc: " << p << G4endl;
691   }
692 }
693 
694 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
695 
696 void G4LossTableManager::CopyTables(const G4ParticleDefinition* part,
697                                     G4VEnergyLossProcess* base_proc)
698 {
699   for (G4int j=0; j<n_loss; ++j) {
700 
701     G4VEnergyLossProcess* proc = loss_vector[j];
702 
703     if (!tables_are_built[j] && part == base_part_vector[j]) {
704       tables_are_built[j] = true;
705       // for base particle approach only ionisation table should be used
706       proc->SetDEDXTable(base_proc->IonisationTable(),fRestricted);
707       proc->SetDEDXTable(base_proc->DEDXunRestrictedTable(),fTotal);
708       proc->SetCSDARangeTable(base_proc->CSDARangeTable());
709       proc->SetRangeTableForLoss(base_proc->RangeTableForLoss());
710       proc->SetInverseRangeTable(base_proc->InverseRangeTable());
711       proc->SetLambdaTable(base_proc->LambdaTable());
712       if(proc->IsIonisationProcess()) { 
713         range_vector[j] = base_proc->RangeTableForLoss();
714         inv_range_vector[j] = base_proc->InverseRangeTable();
715         loss_map[part_vector[j]] = proc; 
716         //G4cout << "G4LossTableManager::CopyTable " 
717         //       << part_vector[j]->GetParticleName()
718         //       << " added to map " << proc << G4endl;
719       }
720       if (1 < verbose) {
721          G4cout << "   CopyTables for " << proc->GetProcessName()
722                 << " for " << part_vector[j]->GetParticleName()
723                 << " base_part= " << part->GetParticleName()
724                 << " tables are assigned"
725                 << G4endl;
726       }
727     }
728   }
729 }
730 
731 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
732 
733 G4VEnergyLossProcess* G4LossTableManager::BuildTables(
734                       const G4ParticleDefinition* aParticle)
735 {
736   if(1 < verbose) {
737     G4cout << "   G4LossTableManager::BuildTables(part) for "
738            << aParticle->GetParticleName() << G4endl;
739   }
740 
741   std::vector<G4PhysicsTable*> t_list;  
742   std::vector<G4VEnergyLossProcess*> loss_list;
743   std::vector<G4bool> build_flags;
744   G4VEnergyLossProcess* em = nullptr;
745   G4VEnergyLossProcess* p = nullptr;
746   G4int iem = 0;
747   G4PhysicsTable* dedx = nullptr;
748   G4int i;
749 
750   G4ProcessVector* pvec = 
751     aParticle->GetProcessManager()->GetProcessList();
752   G4int nvec = (G4int)pvec->size();
753 
754   for (i=0; i<n_loss; ++i) {
755     p = loss_vector[i];
756     if (nullptr != p) {
757       G4bool yes = (aParticle == part_vector[i]);
758 
759       // possible case of process sharing between particle/anti-particle
760       if(!yes) {
761         auto ptr = static_cast<G4VProcess*>(p);
762         for(G4int j=0; j<nvec; ++j) {
763           //G4cout << "j= " << j << " " << (*pvec)[j] << " " << ptr << G4endl;
764           if(ptr == (*pvec)[j]) {
765             yes = true;
766             break;
767           }
768         }
769       }      
770       // process belong to this particle
771       if(yes && isActive[i]) {
772         if (p->IsIonisationProcess() || !em) {
773           em = p;
774           iem= i;
775         }
776         // tables may be shared between particle/anti-particle
777         G4bool val = false;
778         if (!tables_are_built[i]) {
779           val = true;
780           dedx = p->BuildDEDXTable(fRestricted);
781           //G4cout << "===Build DEDX table for " << p->GetProcessName()
782           // << " idx= " << i << " dedx:" << dedx << " " << dedx->length() << G4endl;
783           p->SetDEDXTable(dedx,fRestricted);
784           tables_are_built[i] = true;
785         } else {
786           dedx = p->DEDXTable();
787         }
788         t_list.push_back(dedx);
789         loss_list.push_back(p);
790         build_flags.push_back(val);
791       }
792     }
793   }
794 
795   G4int n_dedx = (G4int)t_list.size();
796   if (0 == n_dedx || !em) {
797     G4cout << "G4LossTableManager WARNING: no DEDX processes for " 
798            << aParticle->GetParticleName() << G4endl;
799     return nullptr;
800   }
801   G4int nSubRegions = em->NumberOfSubCutoffRegions();
802 
803   if (1 < verbose) {
804     G4cout << "     Start to build the sum of " << n_dedx << " processes"
805            << " iem= " << iem << " em= " << em->GetProcessName()
806            << " buildCSDARange= " << theParameters->BuildCSDARange()
807            << " nSubRegions= " << nSubRegions;
808     if(subcutProducer) { 
809       G4cout << " SubCutProducer " << subcutProducer->GetName(); 
810     }
811     G4cout << G4endl;
812   }
813   // do not build tables if producer class is defined
814   if(subcutProducer) { nSubRegions = 0; }
815 
816   dedx = em->DEDXTable();
817   em->SetDEDXTable(dedx, fIsIonisation);
818 
819   if (1 < n_dedx) {
820     dedx = nullptr;
821     dedx = G4PhysicsTableHelper::PreparePhysicsTable(dedx);
822     tableBuilder->BuildDEDXTable(dedx, t_list);
823     em->SetDEDXTable(dedx, fRestricted);
824   }
825 
826   dedx_vector[iem] = dedx;
827 
828   G4PhysicsTable* range = em->RangeTableForLoss();
829   if(!range) range  = G4PhysicsTableHelper::PreparePhysicsTable(range);
830   range_vector[iem] = range;
831 
832   G4PhysicsTable* invrange = em->InverseRangeTable();
833   if(!invrange) invrange = G4PhysicsTableHelper::PreparePhysicsTable(invrange);
834   inv_range_vector[iem]  = invrange;
835 
836   tableBuilder->BuildRangeTable(dedx, range);
837   tableBuilder->BuildInverseRangeTable(range, invrange);
838 
839   em->SetRangeTableForLoss(range);
840   em->SetInverseRangeTable(invrange);
841 
842   std::vector<G4PhysicsTable*> listCSDA;
843 
844   for (i=0; i<n_dedx; ++i) {
845     p = loss_list[i];
846     if(build_flags[i]) {
847       p->SetLambdaTable(p->BuildLambdaTable(fRestricted));
848     }
849     if(theParameters->BuildCSDARange()) { 
850       dedx = p->BuildDEDXTable(fTotal);
851       p->SetDEDXTable(dedx,fTotal);
852       listCSDA.push_back(dedx); 
853     }     
854   }
855 
856   if(theParameters->BuildCSDARange()) {
857     G4PhysicsTable* dedxCSDA = em->DEDXunRestrictedTable();
858     if (1 < n_dedx) {
859       dedxCSDA = G4PhysicsTableHelper::PreparePhysicsTable(dedxCSDA);
860       tableBuilder->BuildDEDXTable(dedxCSDA, listCSDA);
861       em->SetDEDXTable(dedxCSDA,fTotal);
862     }
863     G4PhysicsTable* rCSDA = em->CSDARangeTable();
864     if(!rCSDA) { rCSDA = G4PhysicsTableHelper::PreparePhysicsTable(rCSDA); }
865     tableBuilder->BuildRangeTable(dedxCSDA, rCSDA);
866     em->SetCSDARangeTable(rCSDA);
867   }
868 
869   if (1 < verbose) {
870     G4cout << "G4LossTableManager::BuildTables: Tables are built for "
871            << aParticle->GetParticleName()
872            << "; ionisation process: " << em->GetProcessName()
873            << "  " << em
874            << G4endl;
875   }
876   return em;
877 }
878 
879 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
880 
881 void G4LossTableManager::ParticleHaveNoLoss(
882      const G4ParticleDefinition* aParticle)
883 {
884   G4ExceptionDescription ed;
885   ed << "Energy loss process not found for " << aParticle->GetParticleName() 
886      << " !";
887   G4Exception("G4LossTableManager::ParticleHaveNoLoss", "em0001",
888               FatalException, ed);
889 }
890 
891 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
892 
893 void G4LossTableManager::SetVerbose(G4int val)
894 {
895   verbose = val;
896 }
897 
898 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
899 
900 const std::vector<G4VEnergyLossProcess*>& 
901 G4LossTableManager::GetEnergyLossProcessVector()
902 {
903   return loss_vector;
904 }
905 
906 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
907 
908 const std::vector<G4VEmProcess*>& G4LossTableManager::GetEmProcessVector()
909 {
910   return emp_vector;
911 }
912 
913 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
914 
915 const std::vector<G4VMultipleScattering*>& 
916 G4LossTableManager::GetMultipleScatteringVector()
917 {
918   return msc_vector;
919 }
920 
921 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
922 
923 G4EmSaturation* G4LossTableManager::EmSaturation()
924 {
925   return theParameters->GetEmSaturation();
926 }
927 
928 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
929 
930 G4EmConfigurator* G4LossTableManager::EmConfigurator()
931 {
932   if(!emConfigurator) {
933     emConfigurator = new G4EmConfigurator(verbose); 
934   }
935   return emConfigurator;
936 }
937 
938 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
939 
940 G4ElectronIonPair* G4LossTableManager::ElectronIonPair()
941 {
942   if(!emElectronIonPair) { 
943     emElectronIonPair = new G4ElectronIonPair(verbose);
944   }
945   return emElectronIonPair;
946 }
947 
948 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
949 
950 void G4LossTableManager::SetNIELCalculator(G4NIELCalculator* ptr)
951 {
952   if(nullptr != ptr && ptr != nielCalculator) {
953     delete nielCalculator;
954     nielCalculator = ptr;
955   }
956 }
957 
958 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
959 
960 G4NIELCalculator* G4LossTableManager::NIELCalculator()
961 {
962   if(!nielCalculator) { 
963     nielCalculator = new G4NIELCalculator(nullptr, verbose); 
964   }
965   return nielCalculator;
966 }
967 
968 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
969  
970 void G4LossTableManager::SetAtomDeexcitation(G4VAtomDeexcitation* p)
971 {
972   if(atomDeexcitation != p) {
973     delete atomDeexcitation;
974     atomDeexcitation = p;
975   }
976 }
977 
978 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
979 
980 void G4LossTableManager::SetSubCutProducer(G4VSubCutProducer* p) 
981 {
982   if(subcutProducer != p) {
983     delete subcutProducer;
984     subcutProducer = p;
985   }
986 }
987 
988 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
989 
990 void G4LossTableManager::PrintEWarning(G4String tit, G4double /*val*/)
991 {
992   G4String ss = "G4LossTableManager::" + tit; 
993   G4ExceptionDescription ed;
994   /*
995   ed << "Parameter is out of range: " << val 
996      << " it will have no effect!\n" << " ## " 
997      << " nbins= " << nbinsLambda 
998      << " nbinsPerDecade= " << nbinsPerDecade 
999      << " Emin(keV)= " << minKinEnergy/keV 
1000      << " Emax(GeV)= " << maxKinEnergy/GeV;
1001   */
1002   G4Exception(ss, "em0044", JustWarning, ed);
1003 }
1004 
1005 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1006 
1007 void G4LossTableManager::DumpHtml()
1008 {
1009   // Automatic generation of html documentation page for physics lists
1010   // List processes and models for the most important
1011   // particles in descending order of importance
1012   // NB. for model names with length > 18 characters the .rst file needs
1013   //  to be edited by hand. Or modify G4EmModelManager::DumpModelList
1014 
1015   char* dirName = std::getenv("G4PhysListDocDir");
1016   char* physList = std::getenv("G4PhysListName");
1017   if (dirName && physList) {
1018     G4String physListName = G4String(physList);
1019     G4String pathName = G4String(dirName) + "/" + physListName + ".rst";
1020 
1021     std::ofstream outFile;
1022     outFile.open(pathName);
1023    
1024     outFile << physListName << G4endl;
1025     outFile << std::string(physListName.length(), '=') << G4endl;
1026 
1027     std::vector<G4ParticleDefinition*> particles {
1028         G4Gamma::Gamma(),
1029         G4Electron::Electron(),
1030         G4Positron::Positron(),
1031         G4Proton::ProtonDefinition(),
1032         G4MuonPlus::MuonPlusDefinition(),
1033         G4MuonMinus::MuonMinusDefinition(),
1034       };
1035    
1036     std::vector<G4VEmProcess*> emproc_vector = GetEmProcessVector();
1037     std::vector<G4VEnergyLossProcess*> enloss_vector = 
1038       GetEnergyLossProcessVector();
1039     std::vector<G4VMultipleScattering*> mscat_vector =
1040       GetMultipleScatteringVector();
1041     
1042     for (auto theParticle : particles) {
1043       outFile << G4endl << "**" << theParticle->GetParticleName()
1044               << "**" << G4endl << G4endl << " .. code-block:: none" << G4endl;
1045 
1046       G4ProcessManager* pm = theParticle->GetProcessManager();
1047       G4ProcessVector*  pv = pm->GetProcessList();
1048       G4int plen = pm->GetProcessListLength();
1049 
1050       for (auto emproc : emproc_vector) {
1051         for (G4int i = 0; i < plen; ++i) {
1052           G4VProcess* proc = (*pv)[i];
1053           if (proc == emproc) {
1054             outFile << G4endl;
1055             proc->ProcessDescription(outFile);
1056             break;
1057           }
1058         }
1059       }
1060 
1061       for (auto mscproc : mscat_vector) {
1062         for (G4int i = 0; i < plen; ++i) {
1063           G4VProcess* proc = (*pv)[i];
1064           if (proc == mscproc) {
1065             outFile << G4endl;
1066             proc->ProcessDescription(outFile);
1067             break;
1068           }
1069         }
1070       }
1071 
1072       for (auto enlossproc : enloss_vector) {
1073         for (G4int i = 0; i < plen; ++i) {
1074           G4VProcess* proc = (*pv)[i];
1075           if (proc == enlossproc) {
1076             outFile << G4endl;
1077             proc->ProcessDescription(outFile);
1078             break;
1079           }
1080         }
1081       }
1082     }
1083     outFile.close();
1084   }
1085 }
1086 
1087 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1088 
1089