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