Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/management/src/G4HadronicProcessStore.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 ]

  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 //
 29 // GEANT4 Class file
 30 //
 31 //
 32 // File name:     G4HadronicProcessStore
 33 //
 34 // Author:        Vladimir Ivanchenko
 35 //
 36 // Creation date: 09.05.2008
 37 //
 38 // Modifications:
 39 // 23.01.2009 V.Ivanchenko add destruction of processes
 40 // 12.05.2020 A.Ribon introduced general verbose level in hadronics
 41 //
 42 // Class Description:
 43 // Singleton to store hadronic processes, to provide access to processes
 44 // and to printout information about processes
 45 //
 46 // -------------------------------------------------------------------
 47 //
 48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 50 
 51 #include "G4HadronicProcessStore.hh"
 52 #include "G4SystemOfUnits.hh"
 53 #include "G4UnitsTable.hh"
 54 #include "G4Element.hh"
 55 #include "G4ProcessManager.hh"
 56 #include "G4Electron.hh"
 57 #include "G4Proton.hh"
 58 #include "G4ParticleTable.hh"
 59 #include "G4HadronicInteractionRegistry.hh"
 60 #include "G4CrossSectionDataSetRegistry.hh"
 61 #include "G4HadronicEPTestMessenger.hh"
 62 #include "G4HadronicParameters.hh"
 63 #include "G4HadronicProcessType.hh"
 64 #include <algorithm>
 65 
 66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 67 
 68 G4HadronicProcessStore* G4HadronicProcessStore::Instance()
 69 {
 70   static thread_local auto* _instance = new G4HadronicProcessStore{};
 71   return _instance;
 72 }
 73 
 74 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 75 
 76 G4HadronicProcessStore::~G4HadronicProcessStore()
 77 {
 78   Clean();
 79   delete theEPTestMessenger;
 80 }
 81 
 82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 83 
 84 void G4HadronicProcessStore::Clean()
 85 {
 86   for (auto const& itr : process)
 87     delete itr;
 88   process.clear();
 89 
 90   for (auto const& itr : extraProcess)
 91     delete itr;
 92   extraProcess.clear();
 93 
 94   for (auto const& it : ep_map)
 95     delete it.second;
 96 
 97   m_map.clear();
 98   p_map.clear();
 99 
100   n_extra = 0;
101   n_proc = 0;
102 }
103 
104 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
105 
106 G4HadronicProcessStore::G4HadronicProcessStore()
107 {
108   theGenericIon = 
109     G4ParticleTable::GetParticleTable()->FindParticle("GenericIon");
110   param = G4HadronicParameters::Instance();
111   theEPTestMessenger = new G4HadronicEPTestMessenger(this);
112 }
113 
114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
115 
116 G4double G4HadronicProcessStore::GetCrossSectionPerAtom(
117                                  const G4ParticleDefinition* part,
118                                  G4double energy,
119                                  const G4VProcess* proc,
120                                  const G4Element*  element,
121          const G4Material* material)
122 {
123   G4double cross = 0.;    
124   G4int subType = proc->GetProcessSubType();      
125   if (subType == fHadronElastic)   
126     cross = GetElasticCrossSectionPerAtom(part,energy,element,material);
127   else if (subType == fHadronInelastic)   
128     cross = GetInelasticCrossSectionPerAtom(part,energy,element,material);
129   else if (subType == fCapture)   
130     cross = GetCaptureCrossSectionPerAtom(part,energy,element,material);      
131   else if (subType == fFission)   
132     cross = GetFissionCrossSectionPerAtom(part,energy,element,material); 
133   else if (subType == fChargeExchange)   
134     cross = GetChargeExchangeCrossSectionPerAtom(part,energy,element,material);
135   return cross;
136 }      
137 
138 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
139 
140 G4double G4HadronicProcessStore::GetCrossSectionPerVolume(
141                                  const G4ParticleDefinition* part,
142                                  G4double energy,
143                                  const G4VProcess* proc,
144                                  const G4Material* material)
145 {
146   G4double cross = 0.;    
147   G4int subType = proc->GetProcessSubType();      
148   if (subType == fHadronElastic)   
149     cross = GetElasticCrossSectionPerVolume(part,energy,material);
150   else if (subType == fHadronInelastic)   
151     cross = GetInelasticCrossSectionPerVolume(part,energy,material);
152   else if (subType == fCapture)   
153     cross = GetCaptureCrossSectionPerVolume(part,energy,material);      
154   else if (subType == fFission)   
155     cross = GetFissionCrossSectionPerVolume(part,energy,material); 
156   else if (subType == fChargeExchange)   
157     cross = GetChargeExchangeCrossSectionPerVolume(part,energy,material);
158   return cross;
159 }      
160 
161 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
162 
163 G4double G4HadronicProcessStore::GetElasticCrossSectionPerVolume(
164     const G4ParticleDefinition *aParticle,
165     G4double kineticEnergy,
166     const G4Material *material)
167 {
168   G4double cross = 0.0;
169   const G4ElementVector* theElementVector = material->GetElementVector();
170   const G4double* theAtomNumDensityVector = 
171     material->GetVecNbOfAtomsPerVolume();
172   size_t nelm = material->GetNumberOfElements();
173   for (size_t i=0; i<nelm; ++i) {
174     const G4Element* elm = (*theElementVector)[i];
175     cross += theAtomNumDensityVector[i]*
176       GetElasticCrossSectionPerAtom(aParticle,kineticEnergy,elm,material);
177   }
178   return cross;
179 }
180 
181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
182 
183 G4double G4HadronicProcessStore::GetElasticCrossSectionPerAtom(
184     const G4ParticleDefinition *aParticle,
185     G4double kineticEnergy,
186     const G4Element *anElement, const G4Material* mat)
187 {
188   G4HadronicProcess* hp = FindProcess(aParticle, fHadronElastic);
189   G4double cross = 0.0;
190   localDP.SetKineticEnergy(kineticEnergy);
191   if(hp) {
192     cross = hp->GetElementCrossSection(&localDP,anElement,mat);
193   }
194   return cross;
195 }
196 
197 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
198 
199 G4double G4HadronicProcessStore::GetElasticCrossSectionPerIsotope(
200     const G4ParticleDefinition*,
201     G4double,
202     G4int, G4int)
203 {
204   return 0.0;
205 }
206 
207 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
208 
209 G4double G4HadronicProcessStore::GetInelasticCrossSectionPerVolume(
210     const G4ParticleDefinition *aParticle,
211     G4double kineticEnergy,
212     const G4Material *material)
213 {
214   G4double cross = 0.0;
215   const G4ElementVector* theElementVector = material->GetElementVector();
216   const G4double* theAtomNumDensityVector = 
217     material->GetVecNbOfAtomsPerVolume();
218   size_t nelm = material->GetNumberOfElements();
219   for (size_t i=0; i<nelm; ++i) {
220     const G4Element* elm = (*theElementVector)[i];
221     cross += theAtomNumDensityVector[i]*
222       GetInelasticCrossSectionPerAtom(aParticle,kineticEnergy,elm,material);
223   }
224   return cross;
225 }
226 
227 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
228 
229 G4double G4HadronicProcessStore::GetInelasticCrossSectionPerAtom(
230     const G4ParticleDefinition *aParticle,
231     G4double kineticEnergy,
232     const G4Element *anElement, const G4Material* mat)
233 {
234   G4HadronicProcess* hp = FindProcess(aParticle, fHadronInelastic);
235   localDP.SetKineticEnergy(kineticEnergy);
236   G4double cross = 0.0;
237   if(hp) { 
238     cross = hp->GetElementCrossSection(&localDP,anElement,mat);
239   }
240   return cross;
241 }
242 
243 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
244 
245 G4double G4HadronicProcessStore::GetInelasticCrossSectionPerIsotope(
246     const G4ParticleDefinition *,
247     G4double,
248     G4int, G4int)
249 {
250   return 0.0;
251 }
252 
253 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
254 
255 G4double G4HadronicProcessStore::GetCaptureCrossSectionPerVolume(
256     const G4ParticleDefinition *aParticle,
257     G4double kineticEnergy,
258     const G4Material *material)
259 {
260   G4double cross = 0.0;
261   const G4ElementVector* theElementVector = material->GetElementVector();
262   const G4double* theAtomNumDensityVector = 
263     material->GetVecNbOfAtomsPerVolume();
264   size_t nelm = material->GetNumberOfElements();
265   for (size_t i=0; i<nelm; ++i) {
266     const G4Element* elm = (*theElementVector)[i];
267     cross += theAtomNumDensityVector[i]*
268       GetCaptureCrossSectionPerAtom(aParticle,kineticEnergy,elm,material);
269   }
270   return cross;
271 }
272 
273 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
274 
275 G4double G4HadronicProcessStore::GetCaptureCrossSectionPerAtom(
276     const G4ParticleDefinition *aParticle,
277     G4double kineticEnergy,
278     const G4Element *anElement, const G4Material* mat)
279 {
280   G4HadronicProcess* hp = FindProcess(aParticle, fCapture);
281   localDP.SetKineticEnergy(kineticEnergy);
282   G4double cross = 0.0;
283   if(hp) {
284     cross = hp->GetElementCrossSection(&localDP,anElement,mat);
285   }
286   return cross;
287 }
288 
289 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
290 
291 G4double G4HadronicProcessStore::GetCaptureCrossSectionPerIsotope(
292     const G4ParticleDefinition *,
293     G4double,
294     G4int, G4int)
295 {
296   return 0.0;
297 }
298 
299 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
300 
301 G4double G4HadronicProcessStore::GetFissionCrossSectionPerVolume(
302     const G4ParticleDefinition *aParticle,
303     G4double kineticEnergy,
304     const G4Material *material)
305 {
306   G4double cross = 0.0;
307   const G4ElementVector* theElementVector = material->GetElementVector();
308   const G4double* theAtomNumDensityVector = 
309     material->GetVecNbOfAtomsPerVolume();
310   size_t nelm = material->GetNumberOfElements();
311   for (size_t i=0; i<nelm; i++) {
312     const G4Element* elm = (*theElementVector)[i];
313     cross += theAtomNumDensityVector[i]*
314       GetFissionCrossSectionPerAtom(aParticle,kineticEnergy,elm,material);
315   }
316   return cross;
317 }
318 
319 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
320 
321 G4double G4HadronicProcessStore::GetFissionCrossSectionPerAtom(
322     const G4ParticleDefinition *aParticle,
323     G4double kineticEnergy,
324     const G4Element *anElement, const G4Material* mat)
325 {
326   G4HadronicProcess* hp = FindProcess(aParticle, fFission);
327   localDP.SetKineticEnergy(kineticEnergy);
328   G4double cross = 0.0;
329   if(hp) {
330     cross = hp->GetElementCrossSection(&localDP,anElement,mat);
331   }
332   return cross;
333 }
334 
335 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
336 
337 G4double G4HadronicProcessStore::GetFissionCrossSectionPerIsotope(
338     const G4ParticleDefinition *,
339     G4double,
340     G4int, G4int)
341 {
342   return 0.0;
343 }
344 
345 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
346 
347 G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerVolume(
348     const G4ParticleDefinition *aParticle,
349     G4double kineticEnergy,
350     const G4Material *material)
351 {
352   G4double cross = 0.0;
353   const G4ElementVector* theElementVector = material->GetElementVector();
354   const G4double* theAtomNumDensityVector = 
355     material->GetVecNbOfAtomsPerVolume();
356   size_t nelm = material->GetNumberOfElements();
357   for (size_t i=0; i<nelm; ++i) {
358     const G4Element* elm = (*theElementVector)[i];
359     cross += theAtomNumDensityVector[i]*
360     GetChargeExchangeCrossSectionPerAtom(aParticle,kineticEnergy,elm,material);
361   }
362   return cross;
363 }
364 
365 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
366 
367 G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerAtom(
368     const G4ParticleDefinition *aParticle,
369     G4double kineticEnergy,
370     const G4Element *anElement, const G4Material* mat)
371 {
372   G4HadronicProcess* hp = FindProcess(aParticle, fChargeExchange);
373   localDP.SetKineticEnergy(kineticEnergy);
374   G4double cross = 0.0;
375   if(hp) {
376     cross = hp->GetElementCrossSection(&localDP,anElement,mat);
377   }
378   return cross;
379 }
380 
381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
382 
383 G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerIsotope(
384     const G4ParticleDefinition *,
385     G4double,
386     G4int, G4int)
387 {
388   return 0.0;
389 }
390 
391 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
392 
393 void G4HadronicProcessStore::Register(G4HadronicProcess* proc) 
394 { 
395   for(G4int i=0; i<n_proc; ++i) {
396     if(process[i] == proc) { return; }
397   }
398   if(1 < param->GetVerboseLevel()) {
399     G4cout << "G4HadronicProcessStore::Register hadronic " << n_proc
400      << "  " << proc->GetProcessName() << G4endl;
401   }
402   ++n_proc;
403   process.push_back(proc);
404 }
405 
406 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
407 
408 void G4HadronicProcessStore::RegisterParticle(G4HadronicProcess* proc, 
409                 const G4ParticleDefinition* part) 
410 { 
411   G4int i=0;
412   for(; i<n_proc; ++i) {if(process[i] == proc) break;}
413   G4int j=0;
414   for(; j<n_part; ++j) {if(particle[j] == part) break;}
415 
416   if(1 < param->GetVerboseLevel()) {
417     G4cout << "G4HadronicProcessStore::RegisterParticle " 
418      << part->GetParticleName()
419      << " for  " << proc->GetProcessName() << G4endl;
420   }
421   if(j == n_part) {
422     ++n_part;
423     particle.push_back(part);
424     wasPrinted.push_back(0);
425   }
426   
427   // the pair should be added?
428   if(i < n_proc) {
429     std::multimap<PD,HP,std::less<PD> >::iterator it;
430     for(it=p_map.lower_bound(part); it!=p_map.upper_bound(part); ++it) {
431       if(it->first == part) {
432   HP process2 = (it->second);
433   if(proc == process2) { return; }
434       }
435     }
436   }
437   
438   p_map.insert(std::multimap<PD,HP>::value_type(part,proc));
439 }
440 
441 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
442 
443 void G4HadronicProcessStore::RegisterInteraction(G4HadronicProcess* proc,
444              G4HadronicInteraction* mod)
445 {
446   G4int i=0;
447   for(; i<n_proc; ++i) {if(process[i] == proc) { break; }}
448   G4int k=0;
449   for(; k<n_model; ++k) {if(model[k] == mod) { break; }}
450    
451   m_map.insert(std::multimap<HP,HI>::value_type(proc,mod));
452     
453   if(k == n_model) {
454     ++n_model;
455     model.push_back(mod);
456     modelName.push_back(mod->GetModelName());
457   }
458 }
459 
460 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
461 
462 void G4HadronicProcessStore::DeRegister(G4HadronicProcess* proc)
463 {
464   for(G4int i=0; i<n_proc; ++i) {
465     if(process[i] == proc) {
466       process[i] = nullptr;
467       DeRegisterExtraProcess((G4VProcess*)proc);      
468       return;
469     }
470   }
471 } 
472 
473 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
474 
475 void G4HadronicProcessStore::RegisterExtraProcess(G4VProcess* proc)
476 {
477   for(G4int i=0; i<n_extra; ++i) {
478     if(extraProcess[i] == proc) { return; }
479   }
480   G4HadronicProcess* hproc = static_cast<G4HadronicProcess*>(proc);
481   if(hproc) {
482     for(G4int i=0; i<n_proc; ++i) {
483       if(process[i] == hproc) { return; }
484     }
485   }
486   if(1 < param->GetVerboseLevel()) {
487     G4cout << "Extra Process: " << n_extra 
488      << "  " <<  proc->GetProcessName() << G4endl;
489   }
490   ++n_extra;
491   extraProcess.push_back(proc);
492 }
493 
494 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
495 
496 void G4HadronicProcessStore::RegisterParticleForExtraProcess(
497                              G4VProcess* proc,
498            const G4ParticleDefinition* part)
499 {
500   G4int i=0;
501   for(; i<n_extra; ++i) { if(extraProcess[i] == proc) { break; } }
502   G4int j=0;
503   for(; j<n_part; ++j) { if(particle[j] == part) { break; } }
504 
505   if(j == n_part) {
506     ++n_part;
507     particle.push_back(part);
508     wasPrinted.push_back(0);
509   }
510   
511   // the pair should be added?
512   if(i < n_extra) {
513     std::multimap<PD,G4VProcess*,std::less<PD> >::iterator it;
514     for(it=ep_map.lower_bound(part); it!=ep_map.upper_bound(part); ++it) {
515       if(it->first == part) {
516   G4VProcess* process2 = (it->second);
517   if(proc == process2) { return; }
518       }
519     }
520   }
521   
522   ep_map.insert(std::multimap<PD,G4VProcess*>::value_type(part,proc));
523 } 
524 
525 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
526 
527 void G4HadronicProcessStore::DeRegisterExtraProcess(G4VProcess* proc)
528 {
529   for(G4int i=0; i<n_extra; ++i) {
530     if(extraProcess[i] == proc) {
531       extraProcess[i] = nullptr;
532       if(1 < param->GetVerboseLevel()) {
533   G4cout << "Extra Process: " << i << "  " 
534          <<proc->GetProcessName()<< " is deregisted " << G4endl;
535       }
536       return;
537     }
538   }
539 }
540 
541 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
542 
543 void G4HadronicProcessStore::SetBuildXSTable(G4bool val)
544 {
545   buildXSTable = val;
546 }
547 
548 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
549 
550 G4bool G4HadronicProcessStore::GetBuildXSTable() const
551 {
552   return buildXSTable;
553 }
554 
555 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
556 
557 void G4HadronicProcessStore::PrintInfo(const G4ParticleDefinition* part) 
558 {
559   // Trigger particle/process/model printout only when last particle is 
560   // registered
561   if(buildTableStart && part == particle[n_part - 1]) {
562     buildTableStart = false;
563     Dump(param->GetVerboseLevel());
564     if (!(param->GetPhysListDocDir()).empty()) DumpHtml();
565     G4HadronicInteractionRegistry::Instance()->InitialiseModels();
566   }
567 }
568 
569 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
570 
571 void G4HadronicProcessStore::DumpHtml()
572 {
573   // Automatic generation of html documentation page for physics lists
574   // List processes, models and cross sections for the most important
575   // particles in descending order of importance
576 
577   const G4String& dir = param->GetPhysListDocDir();
578   const G4String& pl = param->GetPhysListName();
579   if (!dir.empty() && !pl.empty()) {
580 
581     // Open output file with path name
582     G4String pathName = dir + "/" + pl + ".html";
583     std::ofstream outFile;
584     outFile.open(pathName);
585 
586     // Write physics list summary file
587     outFile << "<html>\n";
588     outFile << "<head>\n";
589     outFile << "<title>Physics List Summary</title>\n";
590     outFile << "</head>\n";
591     outFile << "<body>\n";
592     outFile << "<h2> Summary of Hadronic Processes, Models and Cross Sections"
593       << " for Physics List " << pl << "</h2>\n";
594     outFile << "<ul>\n";
595 
596     PrintHtml(G4Proton::Proton(), outFile);
597     PrintHtml(G4Neutron::Neutron(), outFile);
598     PrintHtml(G4PionPlus::PionPlus(), outFile); 
599     PrintHtml(G4PionMinus::PionMinus(), outFile);
600     PrintHtml(G4Gamma::Gamma(), outFile);
601     PrintHtml(G4Electron::Electron(), outFile);
602     //    PrintHtml(G4MuonMinus::MuonMinus(), outFile);
603     PrintHtml(G4Positron::Positron(), outFile);
604     PrintHtml(G4KaonPlus::KaonPlus(), outFile);
605     PrintHtml(G4KaonMinus::KaonMinus(), outFile);
606     PrintHtml(G4Lambda::Lambda(), outFile);
607     PrintHtml(G4Alpha::Alpha(), outFile);
608     PrintHtml(G4GenericIon::GenericIon(), outFile);
609 
610     outFile << "</ul>\n";
611     outFile << "</body>\n";
612     outFile << "</html>\n";
613     outFile.close();
614   }
615 }
616 
617 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
618 
619 void G4HadronicProcessStore::PrintHtml(const G4ParticleDefinition* theParticle,
620                                        std::ofstream& outFile)
621 {
622   // Automatic generation of html documentation page for physics lists
623   // List processes for the most important particles in descending order
624   // of importance
625  
626   outFile << "<br> <li><h2><font color=\" ff0000 \">" 
627           << theParticle->GetParticleName() << "</font></h2></li>\n";
628 
629   typedef std::multimap<PD,HP,std::less<PD> > PDHPmap;
630   typedef std::multimap<HP,HI,std::less<HP> > HPHImap;
631 
632   std::pair<PDHPmap::iterator, PDHPmap::iterator> itpart =
633                         p_map.equal_range(theParticle);
634 
635   const G4String& pl = param->GetPhysListName();
636 
637   // Loop over processes assigned to particle
638   G4HadronicProcess* theProcess;
639   for (PDHPmap::iterator it = itpart.first; it != itpart.second; ++it) {
640     theProcess = (*it).second;
641     outFile << "<br> &nbsp;&nbsp; <b><font color=\" 0000ff \">process : "
642             << theProcess->GetProcessName() << "</font></b>\n";
643     outFile << "<ul>\n";
644     outFile << "  <li>";
645     theProcess->ProcessDescription(outFile);
646     outFile << "  <li><b><font color=\" 00AA00 \">models : </font></b>\n";
647     // Loop over models assigned to process
648     std::pair<HPHImap::iterator, HPHImap::iterator> itmod =
649                         m_map.equal_range(theProcess);
650 
651     outFile << "    <ul>\n";
652 
653     for (HPHImap::iterator jt = itmod.first; jt != itmod.second; ++jt) {
654       outFile << "    <li><b><a href=\"" << pl << "_" 
655             << HtmlFileName((*jt).second->GetModelName()) << "\"> "
656               << (*jt).second->GetModelName() << "</a>" 
657               << " from " << (*jt).second->GetMinEnergy()/GeV
658               << " GeV to " << (*jt).second->GetMaxEnergy()/GeV
659               << " GeV </b></li>\n";
660 
661       // Print ModelDescription, ignore that we overwrite files n-times.
662       PrintModelHtml((*jt).second);
663 
664     }
665     outFile << "    </ul>\n";
666     outFile << "  </li>\n";
667 
668     // List cross sections assigned to process
669     outFile << "  <li><b><font color=\" 00AA00 \">cross sections : </font></b>\n";
670     outFile << "    <ul>\n";
671     theProcess->GetCrossSectionDataStore()->DumpHtml(*theParticle, outFile);
672     //        << " \n";
673     outFile << "    </ul>\n";
674 
675     outFile << "  </li>\n";
676     outFile << "</ul>\n";
677 
678   }
679 
680   // Loop over extra (G4VProcess) processes
681   std::multimap<PD,G4VProcess*,std::less<PD> >::iterator itp;
682   for (itp=ep_map.lower_bound(theParticle); itp!=ep_map.upper_bound(theParticle); ++itp) {
683     if (itp->first == theParticle) {
684       G4VProcess* proc = (itp->second);
685       outFile << "<br> &nbsp;&nbsp; <b><font color=\" 0000ff \">process : "
686               << proc->GetProcessName() << "</font></b>\n";
687       outFile << "<ul>\n";
688       outFile << "  <li>";
689       proc->ProcessDescription(outFile);
690       outFile << "  </li>\n";
691       outFile << "</ul>\n";
692     }
693   }
694 
695 } // PrintHtml for particle
696 
697 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
698 
699 void 
700 G4HadronicProcessStore::PrintModelHtml(const G4HadronicInteraction * mod) const
701 {
702   const G4String& dir = param->GetPhysListDocDir();
703   const G4String& pl = param->GetPhysListName();
704   G4String pathName = dir + "/" + pl + "_" + HtmlFileName(mod->GetModelName());
705   std::ofstream outModel;
706   outModel.open(pathName);
707   outModel << "<html>\n";
708   outModel << "<head>\n";
709   outModel << "<title>Description of " << mod->GetModelName() 
710      << "</title>\n";
711   outModel << "</head>\n";
712   outModel << "<body>\n";
713 
714   mod->ModelDescription(outModel);
715 
716   outModel << "</body>\n";
717   outModel << "</html>\n";
718 }
719 
720 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
721 
722 G4String G4HadronicProcessStore::HtmlFileName(const G4String & in) const
723 {
724   G4String str(in);
725 
726   // replace blanks:
727   std::transform(str.begin(), str.end(), str.begin(), [](char ch) {
728       return ch == ' ' ? '_' : ch;
729     });
730   str=str + ".html";    
731   return str;
732 }
733 
734 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
735 
736 void G4HadronicProcessStore::Dump(G4int verb)
737 {
738   G4int level = std::max(param->GetVerboseLevel(), verb);
739   if (0 == level) return;
740   
741  G4cout 
742    << "\n====================================================================\n"
743    << std::setw(60) << "HADRONIC PROCESSES SUMMARY (verbose level " 
744    << level << ")" << G4endl;
745   
746  for (G4int i=0; i<n_part; ++i) {
747     PD part = particle[i];
748     G4String pname = part->GetParticleName();
749     G4bool yes = false;
750 
751     if (level == 1 && (pname == "proton" || 
752            pname == "neutron" ||
753                        pname == "deuteron" ||
754                        pname == "triton" ||
755                        pname == "He3" ||
756                        pname == "alpha" ||
757            pname == "pi+" ||
758            pname == "pi-" ||
759                        pname == "gamma" ||
760                        pname == "e+" ||
761                        pname == "e-" ||
762                        pname == "nu_e" ||
763                        pname == "anti_nu_e" ||
764                        pname == "nu_mu" ||
765                        pname == "anti_nu_mu" ||
766                        pname == "mu+" ||
767                        pname == "mu-" ||
768            pname == "kaon+" ||
769            pname == "kaon-" ||
770            pname == "lambda" ||
771            pname == "anti_lambda" ||
772            pname == "sigma-" ||
773            pname == "D-" ||
774            pname == "B-" ||
775            pname == "GenericIon" ||
776            pname == "hypertriton" ||
777            pname == "anti_neutron" ||
778            pname == "anti_proton" ||
779                        pname == "anti_deuteron" ||
780                        pname == "anti_triton" ||
781                        pname == "anti_He3" ||
782                        pname == "anti_alpha" ||
783                        pname == "anti_hypertriton")) yes = true;
784     if (level > 1) yes = true;     
785     if (yes) {
786       // main processes
787       std::multimap<PD,HP,std::less<PD> >::iterator it;
788 
789       for (it=p_map.lower_bound(part); it!=p_map.upper_bound(part); ++it) {
790   if (it->first == part) {
791     HP proc = (it->second);
792     G4int j=0;
793     for (; j<n_proc; ++j) {
794       if (process[j] == proc) { Print(j, i); }
795     }
796   }
797       }
798       
799       // extra processes
800       std::multimap<PD,G4VProcess*,std::less<PD> >::iterator itp;
801       for(itp=ep_map.lower_bound(part); itp!=ep_map.upper_bound(part); ++itp) {
802   if(itp->first == part) {
803     G4VProcess* proc = (itp->second);
804     if (wasPrinted[i] == 0) {
805             G4cout << "-------------------------------------------------------------------------\n"
806        << std::setw(50) << "Hadronic Processes for " 
807        << part->GetParticleName() << "\n";    
808       wasPrinted[i] = 1;
809     }
810     G4cout << "  Process: " << proc->GetProcessName() << G4endl;
811   }
812       }
813     }
814   }
815 }
816 
817 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
818 
819 void G4HadronicProcessStore::Print(G4int idxProc, G4int idxPart)
820 {
821   G4HadronicProcess* proc = process[idxProc];
822   const G4ParticleDefinition* part = particle[idxPart];
823   if(part == nullptr || proc == nullptr) { return; }
824   if (wasPrinted[idxPart] == 0) {
825     G4cout << "-----------------------------------------------------------------------\n"  
826            << std::setw(50) << "Hadronic Processes for " 
827      << part->GetParticleName() << "\n";
828     wasPrinted[idxPart] = 1;     
829   }
830   
831   G4cout << "  Process: " << proc->GetProcessName();
832 
833   // Append the string "/n" (i.e. "per nucleon") on the kinetic energy of ions.
834   G4String stringEnergyPerNucleon = "";
835   if (part == G4GenericIon::Definition() || 
836       std::abs( part->GetBaryonNumber() ) > 1) {
837     stringEnergyPerNucleon = "/n";
838   }
839   // print cross section factor
840   if(param->ApplyFactorXS()) {
841     G4int pdg = part->GetPDGEncoding();
842     G4int subType = proc->GetProcessSubType();
843     G4double fact = 1.0;
844     if(subType == fHadronInelastic) {
845       if(pdg == 2212 || pdg == 2112) {
846         fact = param->XSFactorNucleonInelastic();
847       } else if(std::abs(pdg) == 211) {
848         fact = param->XSFactorPionInelastic();
849       } else {
850         fact = param->XSFactorHadronInelastic();
851       }
852     } else if(subType == fHadronElastic) {
853       if(pdg == 2212 || pdg == 2112) {
854         fact = param->XSFactorNucleonElastic();
855       } else if(std::abs(pdg) == 211) {
856         fact = param->XSFactorPionElastic();
857       } else {
858         fact = param->XSFactorHadronElastic();
859       }
860     }
861     if(std::abs(fact - 1.0) > 1.e-6) {
862       G4cout << "        XSfactor= " << fact; 
863     }
864   }
865 
866   HI hi = 0;
867   std::multimap<HP,HI,std::less<HP> >::iterator ih;
868   for(ih=m_map.lower_bound(proc); ih!=m_map.upper_bound(proc); ++ih) {
869     if(ih->first == proc) {
870       hi = ih->second;
871       G4int i=0;
872       for(; i<n_model; ++i) {
873   if(model[i] == hi) { break; }
874       }
875       G4cout << "\n        Model: " << std::setw(25) << modelName[i] << ": "  
876        << G4BestUnit(hi->GetMinEnergy(), "Energy")
877        << stringEnergyPerNucleon << " ---> " 
878        << G4BestUnit(hi->GetMaxEnergy(), "Energy") 
879        << stringEnergyPerNucleon;
880     }
881   }
882   G4cout << G4endl;
883   
884   G4CrossSectionDataStore* csds = proc->GetCrossSectionDataStore();
885   csds->DumpPhysicsTable(*part);
886 }
887 
888 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
889 
890 void G4HadronicProcessStore::SetVerbose(G4int val)
891 // this code is obsolete - not optimal change verbose in each thread
892 {
893   G4int i;
894   for(i=0; i<n_proc; ++i) {
895     if(process[i]) { process[i]->SetVerboseLevel(val); }
896   }
897   for(i=0; i<n_model; ++i) {
898     if(model[i]) { model[i]->SetVerboseLevel(val); }
899   }
900 }
901 
902 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
903 
904 G4int G4HadronicProcessStore::GetVerbose()
905 {
906   return param->GetVerboseLevel();
907 }
908 
909 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
910 
911 G4HadronicProcess* G4HadronicProcessStore::FindProcess(
912    const G4ParticleDefinition* part, G4HadronicProcessType subType)
913 {
914   bool isNew = false;
915   G4HadronicProcess* hp = nullptr;
916   localDP.SetDefinition(part);
917 
918   if(part != currentParticle) {
919     const G4ParticleDefinition* p = part;
920     if(p->GetBaryonNumber() > 4 && p->GetParticleType() == "nucleus") {
921       p = theGenericIon;
922     }
923     if(p !=  currentParticle) { 
924       isNew = true;
925       currentParticle = p;
926     }
927   }
928   if(!isNew) { 
929     if(!currentProcess) {
930       isNew = true;
931     } else if(subType == currentProcess->GetProcessSubType()) {
932       hp = currentProcess;
933     } else {
934       isNew = true;
935     }
936   }
937   if(isNew) {
938     std::multimap<PD,HP,std::less<PD> >::iterator it;
939     for(it=p_map.lower_bound(currentParticle); 
940   it!=p_map.upper_bound(currentParticle); ++it) {
941       if(it->first == currentParticle && 
942    subType == (it->second)->GetProcessSubType()) {
943   hp = it->second;
944   break;
945       }
946     }  
947     currentProcess = hp;
948   }
949   return hp;
950 }
951 
952 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
953 
954 void G4HadronicProcessStore::SetEpReportLevel(G4int level)
955 {
956   G4cout << " Setting energy/momentum report level to " << level 
957          << " for " << process.size() << " hadronic processes " << G4endl;
958   for (auto& theProcess : process) {
959     theProcess->SetEpReportLevel(level);
960   }
961 }
962 
963 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
964 
965 void G4HadronicProcessStore::SetProcessAbsLevel(G4double abslevel)
966 {
967   G4cout << " Setting absolute energy/momentum test level to " << abslevel 
968    << G4endl;
969   for (auto& theProcess : process) {
970     G4double rellevel = theProcess->GetEnergyMomentumCheckLevels().first;
971     theProcess->SetEnergyMomentumCheckLevels(rellevel, abslevel);
972   }
973 }
974 
975 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
976 
977 void G4HadronicProcessStore::SetProcessRelLevel(G4double rellevel)
978 {
979   G4cout << " Setting relative energy/momentum test level to " << rellevel 
980    << G4endl;
981   for (auto& theProcess : process) {
982     G4double abslevel = theProcess->GetEnergyMomentumCheckLevels().second;
983     theProcess->SetEnergyMomentumCheckLevels(rellevel, abslevel);
984   }
985 }
986 
987 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
988