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