Geant4  10.02.p01
G4CascadeData.icc
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: G4CascadeData.icc 71890 2013-06-27 22:14:13Z mkelsey $
27 //
28 // 20100803 M. Kelsey -- Move implementations to this .icc file. Use name
29 // string to report output
30 // 20110718 M. Kelsey -- Add inelastic cross-section sum to deal with
31 // suppressing elastic scattering off free nucleons (hydrogen)
32 // 20110719 M. Kelsey -- Add ctor argument for two-body initial state
33 // 20110725 M. Kelsey -- Save initial state as data member
34 // 20110728 M. Kelsey -- Fix Coverity #22954, recursive #include.
35 // 20110923 M. Kelsey -- Add optional ostream& argument to print() fns,
36 // drop all "inline" keywords
37 // 20120608 M. Kelsey -- Fix variable-name "shadowing" compiler warnings.
38 // 20130627 M. Kelsey -- Use new function to print particle name strings.
39 
40 #ifndef G4_CASCADE_DATA_ICC
41 #define G4_CASCADE_DATA_ICC
42 
43 #include "G4InuclParticleNames.hh"
44 #include <iostream>
45 #include <iomanip>
46 
47 
48 // Fill cumulative cross-section arrays from input data
49 template <int NE,int N2,int N3,int N4,int N5,int N6,int N7,int N8,int N9>
50 void G4CascadeData<NE,N2,N3,N4,N5,N6,N7,N8,N9>::initialize() {
51  // Initialize index offsets for cross-section array (can't do globally)
52  index[0] = 0; index[1] = N02; index[2] = N23; index[3] = N24;
53  index[4] = N25; index[5] = N26; index[6] = N27; index[7] = N28;
54  index[8] = N29;
55 
56  // Initialize multiplicity array
57  for (G4int im = 0; im < NM; im++) {
58  G4int start = index[im];
59  G4int stop = index[im+1];
60  for (G4int k = 0; k < NE; k++) {
61  multiplicities[im][k] = 0.0;
62  for (G4int i = start; i < stop; i++) {
63  multiplicities[im][k] += crossSections[i][k];
64  }
65  }
66  }
67 
68  // Initialize total cross section arrays
69  for (G4int k = 0; k < NE; k++) {
70  sum[k] = 0.0;
71  for (G4int im = 0; im < NM; im++) {
72  sum[k] += multiplicities[im][k];
73  }
74  }
75 
76  // Identify elastic scattering channel and subtract from inclusive
77  G4int i2b = 0;
78  for (i2b=index[0]; i2b<index[1]; i2b++) {
79  if (x2bfs[i2b][0]*x2bfs[i2b][1] == initialState) break; // Found it
80  }
81 
82  for (G4int k = 0; k < NE; k++) {
83  if (i2b<index[1]) inelastic[k] = tot[k] - crossSections[i2b][k];
84  else inelastic[k] = tot[k]; // FIXME: No elastic channel in table!
85  }
86 }
87 
88 
89 // Dump complete final state tables, all multiplicities
90 template <int NE,int N2,int N3,int N4,int N5,int N6,int N7,int N8,int N9>
91 void G4CascadeData<NE,N2,N3,N4,N5,N6,N7,N8,N9>::print(std::ostream& os) const {
92  os << "\n " << name << " Total cross section:" << G4endl;
93  printXsec(tot, os);
94  os << "\n Summed cross section:" << G4endl;
95  printXsec(sum, os);
96  os << "\n Inelastic cross section:" << G4endl;
97  printXsec(inelastic, os);
98  os << "\n Individual channel cross sections" << G4endl;
99 
100  for (int im=2; im<NM+2; im++) print(im, os);
101  return;
102 }
103 
104 // Dump tables for specified multiplicity
105 template <int NE,int N2,int N3,int N4,int N5,int N6,int N7,int N8,int N9>
106 void G4CascadeData<NE,N2,N3,N4,N5,N6,N7,N8,N9>::
107 print(G4int mult, std::ostream& os) const {
108  if (mult < 0) { // Old interface used mult == -1 for all
109  print(os);
110  return;
111  }
112 
113  G4int im = mult-2; // Convert multiplicity to array index
114 
115  G4int start = index[im];
116  G4int stop = index[im+1];
117  os << "\n Mulitplicity " << mult << " (indices " << start << " to "
118  << stop-1 << ") summed cross section:" << G4endl;
119 
120  printXsec(multiplicities[im], os);
121 
122  for (G4int i=start; i<stop; i++) {
123  G4int ichan=i-start;
124  os << "\n final state x" << mult << "bfs[" << ichan << "] : ";
125  for (G4int fsi=0; fsi<mult; fsi++) {
126  switch (mult) {
127  case 2: os << " " << G4InuclParticleNames::nameShort(x2bfs[ichan][fsi]);
128  break;
129  case 3: os << " " << G4InuclParticleNames::nameShort(x3bfs[ichan][fsi]);
130  break;
131  case 4: os << " " << G4InuclParticleNames::nameShort(x4bfs[ichan][fsi]);
132  break;
133  case 5: os << " " << G4InuclParticleNames::nameShort(x5bfs[ichan][fsi]);
134  break;
135  case 6: os << " " << G4InuclParticleNames::nameShort(x6bfs[ichan][fsi]);
136  break;
137  case 7: os << " " << G4InuclParticleNames::nameShort(x7bfs[ichan][fsi]);
138  break;
139  case 8: os << " " << G4InuclParticleNames::nameShort(x8bfs[ichan][fsi]);
140  break;
141  case 9: os << " " << G4InuclParticleNames::nameShort(x9bfs[ichan][fsi]);
142  break;
143  default: ;
144  }
145  }
146  os << " -- cross section [" << i << "]:" << G4endl;
147  printXsec(crossSections[i], os);
148  }
149 }
150 
151 // Dump individual cross-section table, two lines of 12 values
152 template <int NE,int N2,int N3,int N4,int N5,int N6,int N7,int N8,int N9>
153 void G4CascadeData<NE,N2,N3,N4,N5,N6,N7,N8,N9>::
154 printXsec(const G4double (&xsec)[NE], std::ostream& os) const {
155  for (G4int k=0; k<NE; k++) {
156  os << " " << std::setw(6) << xsec[k];
157  if ((k+1)%10 == 0) os << G4endl;
158  }
159  os << G4endl;
160 }
161 
162 #endif /* G4_CASCADE_DATA_ICC */