| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 /// \file XSHistoManager.cc 26 /// \file XSHistoManager.cc
27 /// \brief Create a set of profiles for XS st 27 /// \brief Create a set of profiles for XS study.
28 // 28 //
29 // Adapted from hadronic/Hadr00/src/HistoMana 29 // Adapted from hadronic/Hadr00/src/HistoManager.cc
30 // Author: G.Hugo, 06 January 2023 30 // Author: G.Hugo, 06 January 2023
31 // 31 //
32 // ******************************************* 32 // ***************************************************************************
33 // 33 //
34 // XSHistoManager 34 // XSHistoManager
35 // 35 //
36 /// Create a set of profiles for XS study. 36 /// Create a set of profiles for XS study.
37 /// 37 ///
38 /// All profiles are G4H1. << 38 /// All profiles are G4H1.
39 /// They are created and filled via G4VAnalys 39 /// They are created and filled via G4VAnalysisManager.
40 /// 40 ///
41 /// The profiles can be dumped to all usual f << 41 /// The profiles can be dumped to all usual formats, including ROOT
42 /// (via G4VAnalysisManager). 42 /// (via G4VAnalysisManager).
43 /// They are also dumped in a format compatib 43 /// They are also dumped in a format compatible with Flair
44 /// (via tools::histo::flair). 44 /// (via tools::histo::flair).
45 // 45 //
46 // ******************************************* 46 // ***************************************************************************
47 47
48 #include "XSHistoManager.hh" 48 #include "XSHistoManager.hh"
49 49
50 #include "G4Element.hh" << 50 #include "G4ios.hh"
51 #include "G4HadronicProcessStore.hh" << 51
52 #include "G4Material.hh" <<
53 #include "G4NistManager.hh" <<
54 #include "G4ParticleDefinition.hh" 52 #include "G4ParticleDefinition.hh"
55 #include "G4ParticleTable.hh" 53 #include "G4ParticleTable.hh"
56 #include "G4ios.hh" <<
57 54
58 // #include "G4AnalysisManager.hh" << 55 #include "G4NistManager.hh"
59 #include "tools_histo_flair.hh" << 56 #include "G4Element.hh"
>> 57 #include "G4Material.hh"
>> 58
>> 59 #include "G4HadronicProcessStore.hh"
60 60
>> 61 //#include "G4AnalysisManager.hh"
61 #include "G4RootAnalysisManager.hh" 62 #include "G4RootAnalysisManager.hh"
62 63
>> 64 #include "tools_histo_flair.hh"
>> 65
>> 66
63 //....oooOO0OOooo........oooOO0OOooo........oo 67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
64 68
65 XSHistoManager::XSHistoManager() << 69 XSHistoManager::XSHistoManager() :
66 : fMessenger(new XSHistoManagerMessenger(thi << 70 fMessenger(new XSHistoManagerMessenger(this)),
67 fOutputFileName("all_XS"), << 71 fOutputFileName("all_XS"),
68 fRootOutputFileName("all_XS.root"), << 72 fRootOutputFileName("all_XS.root"),
69 fFlairOutputFileName("all_XS.hist"), << 73 fFlairOutputFileName("all_XS.hist"),
70 fParticle(nullptr), << 74 fParticle(nullptr),
71 fElement(nullptr), << 75 fElement(nullptr),
72 fMaterial(nullptr), << 76 fMaterial(nullptr),
73 fNumBins(10000), << 77 fNumBins(10000),
74 fMinKineticEnergy(1. * keV), << 78 fMinKineticEnergy(1.*keV),
75 fMaxKineticEnergy(10. * TeV), << 79 fMaxKineticEnergy(10.*TeV),
76 fFunctionName("none"), << 80 fFunctionName("none"),
77 fBinSchemeName("log"), << 81 fBinSchemeName("log"),
78 fRootEnergyUnit("MeV"), << 82 fRootEnergyUnit("MeV"),
79 fAnalysisManager(G4RootAnalysisManager::In << 83 fAnalysisManager(G4RootAnalysisManager::Instance()),
80 fElasticXSIndex(0), << 84 fElasticXSIndex(0),
81 fInelasticXSIndex(1), << 85 fInelasticXSIndex(1),
82 fCaptureXSIndex(2), << 86 fCaptureXSIndex(2),
83 fFissionXSIndex(3), << 87 fFissionXSIndex(3),
84 fChargeExchangeXSIndex(4), << 88 fChargeExchangeXSIndex(4),
85 fTotalXSIndex(5), << 89 fTotalXSIndex(5),
86 fElasticPerVolumeXSIndex(6), << 90 fElasticPerVolumeXSIndex(6),
87 fInelasticPerVolumeXSIndex(7) << 91 fInelasticPerVolumeXSIndex(7)
88 { 92 {
89 // G4NistManager::Instance()->ListMaterials( << 93 //G4NistManager::Instance()->ListMaterials("all");
90 } 94 }
91 95
>> 96
92 // ******************************************* 97 // ***************************************************************************
93 // Set output files names: 2 formats supported 98 // Set output files names: 2 formats supported, ROOT and Flair.
94 // ******************************************* 99 // ***************************************************************************
95 void XSHistoManager::SetOutputFileName(const G << 100 void XSHistoManager::SetOutputFileName(const G4String& outputFileName) {
96 { <<
97 fOutputFileName = outputFileName; 101 fOutputFileName = outputFileName;
98 fRootOutputFileName = outputFileName + ".roo 102 fRootOutputFileName = outputFileName + ".root";
99 fFlairOutputFileName = outputFileName + ".hi 103 fFlairOutputFileName = outputFileName + ".hist";
100 } 104 }
101 105
>> 106
102 // ******************************************* 107 // ***************************************************************************
103 // Set the particle considered for XS study. 108 // Set the particle considered for XS study.
104 // ******************************************* 109 // ***************************************************************************
105 void XSHistoManager::SetParticle(const G4Strin << 110 void XSHistoManager::SetParticle(const G4String& particleName) {
106 { <<
107 fParticle = G4ParticleTable::GetParticleTabl 111 fParticle = G4ParticleTable::GetParticleTable()->FindParticle(particleName);
108 } 112 }
109 113
>> 114
110 // ******************************************* 115 // ***************************************************************************
111 // Set the target element considered for XS st 116 // Set the target element considered for XS study.
112 // ******************************************* 117 // ***************************************************************************
113 void XSHistoManager::SetElement(const G4String << 118 void XSHistoManager::SetElement(const G4String& elementName) {
114 { <<
115 fElement = G4NistManager::Instance()->FindOr 119 fElement = G4NistManager::Instance()->FindOrBuildElement(elementName);
116 // Also needs to set material! 120 // Also needs to set material!
117 SetMaterial(elementName); 121 SetMaterial(elementName);
118 } 122 }
119 123
>> 124
120 // ******************************************* 125 // ***************************************************************************
121 // Set the target material considered for XS s 126 // Set the target material considered for XS study.
122 // ******************************************* 127 // ***************************************************************************
123 void XSHistoManager::SetMaterial(const G4Strin << 128 void XSHistoManager::SetMaterial(const G4String& materialName) {
124 { << 129
125 // Check that material is not set already. 130 // Check that material is not set already.
126 if (fMaterial) { 131 if (fMaterial) {
127 G4ExceptionDescription msg; 132 G4ExceptionDescription msg;
128 msg << "Please use UI command /allXS/eleme 133 msg << "Please use UI command /allXS/elementName"
129 << " OR UI command /allXS/nonElementar 134 << " OR UI command /allXS/nonElementaryMaterialName,"
130 << " BUT NOT BOTH!" << G4endl; << 135 << " BUT NOT BOTH!"
131 G4Exception("XSHistoManager::SetMaterial", << 136 << G4endl;
132 FatalException, msg); << 137 G4Exception("XSHistoManager::SetMaterial",
>> 138 "A target material is already defined.",
>> 139 FatalException,
>> 140 msg);
133 } 141 }
134 142
135 fMaterial = G4NistManager::Instance()->FindO 143 fMaterial = G4NistManager::Instance()->FindOrBuildMaterial("G4_" + materialName);
136 } 144 }
137 145
>> 146
138 // ******************************************* 147 // ***************************************************************************
139 // Open output file + create all profiles cons 148 // Open output file + create all profiles considered for XS study.
140 // All profiles are G4H1, created via G4VAnaly 149 // All profiles are G4H1, created via G4VAnalysisManager.
141 // ******************************************* 150 // ***************************************************************************
142 void XSHistoManager::Book() << 151 void XSHistoManager::Book() {
143 { << 152
144 // Check all XSHistoManager data is set prop 153 // Check all XSHistoManager data is set properly.
145 CheckInput(); 154 CheckInput();
146 155
147 // Open file. 156 // Open file.
148 if (!fAnalysisManager->OpenFile(fRootOutputF 157 if (!fAnalysisManager->OpenFile(fRootOutputFileName)) {
>> 158
149 G4ExceptionDescription msg; 159 G4ExceptionDescription msg;
150 msg << "Booking profiles: cannot open file << 160 msg << "Booking profiles: cannot open file " << fRootOutputFileName
151 G4Exception("XSHistoManager::Book", "Canno << 161 << G4endl;
>> 162 G4Exception("XSHistoManager::Book",
>> 163 "Cannot open file",
>> 164 FatalException,
>> 165 msg);
152 } 166 }
153 G4cout << "### XSHistoManager::Book: Success << 167 G4cout << "### XSHistoManager::Book: Successfully opened file "
154 << " for dumping profiles." << G4endl << 168 << fRootOutputFileName
>> 169 << " for dumping profiles."
>> 170 << G4endl;
155 171
156 // Create all G4H1, and keep track of each h 172 // Create all G4H1, and keep track of each histo index in fXSProfileIndex.
157 const G4int elasticXSProfileIndex = << 173 const G4int elasticXSProfileIndex = fAnalysisManager->CreateH1("ElasticXS",
158 fAnalysisManager->CreateH1("ElasticXS", "E << 174 "Elastic XS",
159 fMaxKineticEner << 175 fNumBins,
160 fXSProfileIndex.insert(std::make_pair(fElast << 176 fMinKineticEnergy,
161 << 177 fMaxKineticEnergy,
162 const G4int inelasticXSProfileIndex = << 178 fRootEnergyUnit,
163 fAnalysisManager->CreateH1("InelasticXS", << 179 fFunctionName,
164 fMaxKineticEner << 180 fBinSchemeName);
165 fXSProfileIndex.insert(std::make_pair(fInela << 181 fXSProfileIndex.insert(std::make_pair(fElasticXSIndex,
166 << 182 elasticXSProfileIndex));
167 const G4int captureXSProfileIndex = << 183
168 fAnalysisManager->CreateH1("CaptureXS", "C << 184 const G4int inelasticXSProfileIndex = fAnalysisManager->CreateH1("InelasticXS",
169 fMaxKineticEner << 185 "Inelastic XS",
170 fXSProfileIndex.insert(std::make_pair(fCaptu << 186 fNumBins,
171 << 187 fMinKineticEnergy,
172 const G4int fissionXSProfileIndex = << 188 fMaxKineticEnergy,
173 fAnalysisManager->CreateH1("FissionXS", "F << 189 fRootEnergyUnit,
174 fMaxKineticEner << 190 fFunctionName,
175 fXSProfileIndex.insert(std::make_pair(fFissi << 191 fBinSchemeName);
176 << 192 fXSProfileIndex.insert(std::make_pair(fInelasticXSIndex,
177 const G4int chargeExchangeXSProfileIndex = f << 193 inelasticXSProfileIndex));
178 "ChargeExchangeXS", "Charge exchange XS", << 194
179 fRootEnergyUnit, fFunctionName, fBinScheme << 195 const G4int captureXSProfileIndex = fAnalysisManager->CreateH1("CaptureXS",
180 fXSProfileIndex.insert(std::make_pair(fCharg << 196 "Capture XS",
181 << 197 fNumBins,
182 const G4int totalXSProfileIndex = << 198 fMinKineticEnergy,
183 fAnalysisManager->CreateH1("TotalXS", "Tot << 199 fMaxKineticEnergy,
184 fMaxKineticEner << 200 fRootEnergyUnit,
185 fXSProfileIndex.insert(std::make_pair(fTotal << 201 fFunctionName,
186 << 202 fBinSchemeName);
187 const G4int elasticPerVolumeXSProfileIndex = << 203 fXSProfileIndex.insert(std::make_pair(fCaptureXSIndex,
188 "ElasticPerVolumeXS", "Elastic XS per volu << 204 captureXSProfileIndex));
189 fRootEnergyUnit, fFunctionName, fBinScheme << 205
190 fXSProfileIndex.insert(std::make_pair(fElast << 206 const G4int fissionXSProfileIndex = fAnalysisManager->CreateH1("FissionXS",
191 << 207 "Fission XS",
192 const G4int inelasticPerVolumeXSProfileIndex << 208 fNumBins,
193 "InelasticPerVolumeXS", "Inelastic XS per << 209 fMinKineticEnergy,
194 fMaxKineticEnergy, fRootEnergyUnit, fFunct << 210 fMaxKineticEnergy,
195 fXSProfileIndex.insert( << 211 fRootEnergyUnit,
196 std::make_pair(fInelasticPerVolumeXSIndex, << 212 fFunctionName,
>> 213 fBinSchemeName);
>> 214 fXSProfileIndex.insert(std::make_pair(fFissionXSIndex,
>> 215 fissionXSProfileIndex));
>> 216
>> 217 const G4int chargeExchangeXSProfileIndex = fAnalysisManager->CreateH1("ChargeExchangeXS",
>> 218 "Charge exchange XS",
>> 219 fNumBins,
>> 220 fMinKineticEnergy,
>> 221 fMaxKineticEnergy,
>> 222 fRootEnergyUnit,
>> 223 fFunctionName,
>> 224 fBinSchemeName);
>> 225 fXSProfileIndex.insert(std::make_pair(fChargeExchangeXSIndex,
>> 226 chargeExchangeXSProfileIndex));
>> 227
>> 228 const G4int totalXSProfileIndex = fAnalysisManager->CreateH1("TotalXS",
>> 229 "Total XS",
>> 230 fNumBins,
>> 231 fMinKineticEnergy,
>> 232 fMaxKineticEnergy,
>> 233 fRootEnergyUnit,
>> 234 fFunctionName,
>> 235 fBinSchemeName);
>> 236 fXSProfileIndex.insert(std::make_pair(fTotalXSIndex,
>> 237 totalXSProfileIndex));
>> 238
>> 239 const G4int elasticPerVolumeXSProfileIndex = fAnalysisManager->CreateH1("ElasticPerVolumeXS",
>> 240 "Elastic XS per volume",
>> 241 fNumBins,
>> 242 fMinKineticEnergy,
>> 243 fMaxKineticEnergy,
>> 244 fRootEnergyUnit,
>> 245 fFunctionName,
>> 246 fBinSchemeName);
>> 247 fXSProfileIndex.insert(std::make_pair(fElasticPerVolumeXSIndex,
>> 248 elasticPerVolumeXSProfileIndex));
>> 249
>> 250 const G4int inelasticPerVolumeXSProfileIndex =
>> 251 fAnalysisManager->CreateH1("InelasticPerVolumeXS",
>> 252 "Inelastic XS per volume",
>> 253 fNumBins,
>> 254 fMinKineticEnergy,
>> 255 fMaxKineticEnergy,
>> 256 fRootEnergyUnit,
>> 257 fFunctionName,
>> 258 fBinSchemeName);
>> 259 fXSProfileIndex.insert(std::make_pair(fInelasticPerVolumeXSIndex,
>> 260 inelasticPerVolumeXSProfileIndex));
197 } 261 }
198 262
>> 263
199 // ******************************************* 264 // ***************************************************************************
200 // Fill all plots, then dump them into relevan 265 // Fill all plots, then dump them into relevant formats.
201 // ******************************************* 266 // ***************************************************************************
202 void XSHistoManager::EndOfRun() << 267 void XSHistoManager::EndOfRun() {
203 { << 268 G4cout << "### XSHistoManager::EndOfRun: Compute & fill XS for "
204 G4cout << "### XSHistoManager::EndOfRun: Com << 269 << fParticle->GetParticleName()
205 << " in " << (fElement ? fElement->Ge << 270 << " in " << (fElement ? fElement->GetName() : fMaterial->GetName())
206 << 271 << G4endl;
207 G4HadronicProcessStore* const store = G4Hadr << 272
208 << 273 G4HadronicProcessStore* const store = G4HadronicProcessStore::Instance();
209 // Fill XS profiles. << 274
210 const G4double logMinKineticEnergy = std::lo << 275 // Fill XS profiles.
211 const G4double logMaxKineticEnergy = std::lo << 276 const G4double logMinKineticEnergy = std::log10(fMinKineticEnergy);
212 const G4double deltaLogKineticEnergy = (logM << 277 const G4double logMaxKineticEnergy = std::log10(fMaxKineticEnergy);
213 << 278 const G4double deltaLogKineticEnergy =
214 G4double logKineticEnergy = logMinKineticEne << 279 (logMaxKineticEnergy - logMinKineticEnergy) / fNumBins;
215 << 280
216 // Loop on all kinetic energies of interest. << 281 G4double logKineticEnergy = logMinKineticEnergy - deltaLogKineticEnergy/2.;
217 for (G4int binIndex = 0; binIndex < fNumBins << 282
218 logKineticEnergy += deltaLogKineticEnergy; << 283 // Loop on all kinetic energies of interest.
219 const G4double kineticEnergy = std::pow(10 << 284 for (G4int binIndex = 0; binIndex < fNumBins; ++binIndex) {
220 << 285
221 G4double totalXS = 0.; << 286 logKineticEnergy += deltaLogKineticEnergy;
222 if (fElement) { << 287 const G4double kineticEnergy = std::pow(10., logKineticEnergy) * MeV;
223 // ELASTIC (ELEMENTARY MATERIAL) << 288
224 const G4double elasticXS = << 289 G4double totalXS = 0.;
225 store->GetElasticCrossSectionPerAtom(f << 290 if (fElement) {
226 fAnalysisManager->FillH1(fXSProfileIndex << 291 // ELASTIC (ELEMENTARY MATERIAL)
227 totalXS += elasticXS; << 292 const G4double elasticXS = store->GetElasticCrossSectionPerAtom(fParticle,
228 << 293 kineticEnergy,
229 // INELASTIC (ELEMENTARY MATERIAL) << 294 fElement,
230 const G4double inelasticXS = << 295 fMaterial);
231 store->GetInelasticCrossSectionPerAtom << 296 fAnalysisManager->FillH1(fXSProfileIndex[fElasticXSIndex],
232 fAnalysisManager->FillH1(fXSProfileIndex << 297 kineticEnergy,
233 inelasticXS / b << 298 elasticXS/barn);
234 totalXS += inelasticXS; << 299 totalXS += elasticXS;
235 << 300
236 if (fParticle == G4Neutron::Definition() << 301 // INELASTIC (ELEMENTARY MATERIAL)
237 // NEUTRON CAPTURE (ELEMENTARY MATERIA << 302 const G4double inelasticXS = store->GetInelasticCrossSectionPerAtom(fParticle,
238 const G4double captureXS = << 303 kineticEnergy,
239 store->GetCaptureCrossSectionPerAtom << 304 fElement,
240 fAnalysisManager->FillH1(fXSProfileInd << 305 fMaterial);
241 totalXS += captureXS; << 306 fAnalysisManager->FillH1(fXSProfileIndex[fInelasticXSIndex],
242 << 307 kineticEnergy,
243 // FISSION (ELEMENTARY MATERIAL) << 308 inelasticXS/barn);
244 const G4double fissionXS = << 309 totalXS += inelasticXS;
245 store->GetFissionCrossSectionPerAtom << 310
246 totalXS += fissionXS; << 311 if (fParticle == G4Neutron::Definition()) {
247 fAnalysisManager->FillH1(fXSProfileInd << 312 // NEUTRON CAPTURE (ELEMENTARY MATERIAL)
248 } << 313 const G4double captureXS = store->GetCaptureCrossSectionPerAtom(fParticle,
249 << 314 kineticEnergy,
250 // CHARGE EXCHANGE (ELEMENTARY MATERIAL) << 315 fElement,
251 const G4double chargeExchangeXS = << 316 fMaterial);
252 store->GetChargeExchangeCrossSectionPe << 317 fAnalysisManager->FillH1(fXSProfileIndex[fCaptureXSIndex],
253 fAnalysisManager->FillH1(fXSProfileIndex << 318 kineticEnergy,
254 chargeExchangeX << 319 captureXS/barn);
255 totalXS += chargeExchangeXS; << 320 totalXS += captureXS;
256 << 321
257 // TOTAL (ELEMENTARY MATERIAL) << 322 // FISSION (ELEMENTARY MATERIAL)
258 fAnalysisManager->FillH1(fXSProfileIndex << 323 const G4double fissionXS = store->GetFissionCrossSectionPerAtom(fParticle,
259 } << 324 kineticEnergy,
260 << 325 fElement,
261 if (fMaterial) { << 326 fMaterial);
262 const G4double materialSurfacicDensity = << 327 totalXS += fissionXS;
263 (fMaterial ? fMaterial->GetDensity() / << 328 fAnalysisManager->FillH1(fXSProfileIndex[fFissionXSIndex],
264 << 329 kineticEnergy,
265 // ELASTIC << 330 fissionXS/barn);
266 const G4double elasticPerVolumeXS = << 331 }
267 store->GetElasticCrossSectionPerVolume << 332
268 fAnalysisManager->FillH1(fXSProfileIndex << 333 // CHARGE EXCHANGE (ELEMENTARY MATERIAL)
269 elasticPerVolum << 334 const G4double chargeExchangeXS =
270 << 335 store->GetChargeExchangeCrossSectionPerAtom(fParticle,
271 // INELASTIC << 336 kineticEnergy,
272 const G4double inelasticPerVolumeXS = << 337 fElement,
273 store->GetInelasticCrossSectionPerVolu << 338 fMaterial);
274 fAnalysisManager->FillH1(fXSProfileIndex << 339 fAnalysisManager->FillH1(fXSProfileIndex[fChargeExchangeXSIndex],
275 inelasticPerVol << 340 kineticEnergy,
276 } << 341 chargeExchangeXS/barn);
277 } << 342 totalXS += chargeExchangeXS;
278 << 343
279 // DUMP G4H1 PLOTS INTO ROOT FILE << 344 // TOTAL (ELEMENTARY MATERIAL)
280 DumpAllG4H1IntoRootFile(); << 345 fAnalysisManager->FillH1(fXSProfileIndex[fTotalXSIndex],
281 << 346 kineticEnergy,
282 // DUMP G4H1 PLOTS INTO FLAIR FILE << 347 totalXS/barn);
283 DumpAllG4H1IntoFlairFile(); << 348 }
284 << 349
285 // Close and clear fAnalysisManager. << 350 if (fMaterial) {
286 fAnalysisManager->CloseFile(); << 351 const G4double materialSurfacicDensity = (fMaterial ?
287 fAnalysisManager->Clear(); << 352 fMaterial->GetDensity() / (g/cm2)
>> 353 : 1.);
>> 354
>> 355 // ELASTIC
>> 356 const G4double elasticPerVolumeXS =
>> 357 store->GetElasticCrossSectionPerVolume(fParticle,
>> 358 kineticEnergy,
>> 359 fMaterial);
>> 360 fAnalysisManager->FillH1(fXSProfileIndex[fElasticPerVolumeXSIndex],
>> 361 kineticEnergy,
>> 362 elasticPerVolumeXS/materialSurfacicDensity);
>> 363
>> 364 // INELASTIC
>> 365 const G4double inelasticPerVolumeXS =
>> 366 store->GetInelasticCrossSectionPerVolume(fParticle,
>> 367 kineticEnergy,
>> 368 fMaterial);
>> 369 fAnalysisManager->FillH1(fXSProfileIndex[fInelasticPerVolumeXSIndex],
>> 370 kineticEnergy,
>> 371 inelasticPerVolumeXS/materialSurfacicDensity);
>> 372 }
>> 373
>> 374 }
>> 375
>> 376
>> 377 // DUMP G4H1 PLOTS INTO ROOT FILE
>> 378 DumpAllG4H1IntoRootFile();
>> 379
>> 380 // DUMP G4H1 PLOTS INTO FLAIR FILE
>> 381 DumpAllG4H1IntoFlairFile();
>> 382
>> 383
>> 384 // Close and clear fAnalysisManager.
>> 385 fAnalysisManager->CloseFile();
>> 386 fAnalysisManager->Clear();
288 } 387 }
289 388
>> 389
290 // ******************************************* 390 // ***************************************************************************
291 // Checks that particle and material are set << 391 // Checks that particle and material are set
292 // (all others have relevant default values). 392 // (all others have relevant default values).
293 // ******************************************* 393 // ***************************************************************************
294 void XSHistoManager::CheckInput() << 394 void XSHistoManager::CheckInput() {
295 { << 395
296 if (!fParticle) { 396 if (!fParticle) {
297 G4ExceptionDescription msg; 397 G4ExceptionDescription msg;
298 msg << "Please add a particle to study XS: << 398 msg << "Please add a particle to study XS: UI command /allXS/particleName"
299 G4Exception("XSHistoManager::CheckInput()" << 399 << G4endl;
300 FatalException, msg); << 400 G4Exception("XSHistoManager::CheckInput()",
>> 401 "Print XS: no input particle defined.",
>> 402 FatalException,
>> 403 msg);
301 } 404 }
302 405
303 if (!fMaterial) { 406 if (!fMaterial) {
304 G4ExceptionDescription msg; 407 G4ExceptionDescription msg;
305 msg << "Please add a material to study XS: 408 msg << "Please add a material to study XS:"
306 << " UI command /allXS/elementName for 409 << " UI command /allXS/elementName for an elementary material,"
307 << " or UI command /allXS/nonElementar 410 << " or UI command /allXS/nonElementaryMaterialName for a compound/mixture material."
308 << G4endl; 411 << G4endl;
309 G4Exception("XSHistoManager::CheckInput()" << 412 G4Exception("XSHistoManager::CheckInput()",
310 FatalException, msg); << 413 "Print XS: no target material defined.",
>> 414 FatalException,
>> 415 msg);
311 } 416 }
312 } 417 }
313 418
>> 419
314 // ******************************************* 420 // ***************************************************************************
315 // DUMP G4H1 PLOTS INTO ROOT FILE (via G4VAnal 421 // DUMP G4H1 PLOTS INTO ROOT FILE (via G4VAnalysisManager).
316 // ******************************************* 422 // ***************************************************************************
317 void XSHistoManager::DumpAllG4H1IntoRootFile() << 423 void XSHistoManager::DumpAllG4H1IntoRootFile() const {
318 { << 424
319 if (!fAnalysisManager->Write()) { << 425 if (!fAnalysisManager->Write()) {
320 G4ExceptionDescription message; << 426 G4ExceptionDescription message;
321 message << "Could not write ROOT file."; << 427 message << "Could not write ROOT file.";
322 G4Exception("XSHistoManager::EndOfRun()", << 428 G4Exception("XSHistoManager::EndOfRun()",
323 } << 429 "I/O Error",
324 G4cout << "### All profiles saved to " << fR << 430 FatalException,
>> 431 message);
>> 432 }
>> 433 G4cout << "### All profiles saved to " << fRootOutputFileName << G4endl;
325 } 434 }
326 435
>> 436
327 // ******************************************* 437 // ***************************************************************************
328 // DUMP G4H1 PLOTS INTO FLAIR FILE (via tools: 438 // DUMP G4H1 PLOTS INTO FLAIR FILE (via tools::histo::flair).
329 // ******************************************* 439 // ***************************************************************************
330 void XSHistoManager::DumpAllG4H1IntoFlairFile( << 440 void XSHistoManager::DumpAllG4H1IntoFlairFile() const {
331 { << 441
332 std::ofstream output; 442 std::ofstream output;
333 output.open(fFlairOutputFileName, std::ios_b 443 output.open(fFlairOutputFileName, std::ios_base::out);
334 auto const rootAnalysisManager = dynamic_cas 444 auto const rootAnalysisManager = dynamic_cast<G4RootAnalysisManager*>(fAnalysisManager);
335 445
336 G4int indexInOutputFile = 1; 446 G4int indexInOutputFile = 1;
337 for (G4int xsIndex = fElasticXSIndex; xsInde 447 for (G4int xsIndex = fElasticXSIndex; xsIndex <= fInelasticPerVolumeXSIndex; ++xsIndex) {
>> 448
338 const G4int histoIndex = fXSProfileIndex.a 449 const G4int histoIndex = fXSProfileIndex.at(xsIndex);
339 const G4String& histoName = fAnalysisManag 450 const G4String& histoName = fAnalysisManager->GetH1Name(histoIndex);
340 const auto& histo = rootAnalysisManager->G 451 const auto& histo = rootAnalysisManager->GetH1(histoIndex);
341 << 452
342 tools::histo::flair::dumpG4H1ProfileInFlai << 453 tools::histo::flair::dumpG4H1ProfileInFlairFormat(output,
>> 454 indexInOutputFile,
>> 455 histoName,
>> 456 histo,
343 457 tools::histo::flair::Abscissa::KineticEnergy,
344 458 fBinSchemeName);
345 ++indexInOutputFile; 459 ++indexInOutputFile;
346 } 460 }
347 output.close(); 461 output.close();
348 G4cout << "### All profiles saved to " << fF 462 G4cout << "### All profiles saved to " << fFlairOutputFileName << G4endl;
349 } 463 }
350 464
351 //....oooOO0OOooo........oooOO0OOooo........oo 465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
352 466