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