Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4CrossSectionDataStore.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//
27// -------------------------------------------------------------------
28//
29// GEANT4 Class file
30//
31//
32// File name: G4CrossSectionDataStore
33//
34// Modifications:
35// 23.01.2009 V.Ivanchenko add destruction of data sets
36// 29.04.2010 G.Folger modifictaions for integer A & Z
37// 14.03.2011 V.Ivanchenko fixed DumpPhysicsTable
38// 15.08.2011 G.Folger, V.Ivanchenko, T.Koi, D.Wright redesign the class
39// 07.03.2013 M.Maire cosmetic in DumpPhysicsTable
40//
41//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
42//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
43
45#include "G4SystemOfUnits.hh"
46#include "G4UnitsTable.hh"
47#include "Randomize.hh"
48#include "G4Nucleus.hh"
49
50#include "G4DynamicParticle.hh"
51#include "G4Isotope.hh"
52#include "G4Element.hh"
53#include "G4Material.hh"
54#include "G4NistManager.hh"
55#include <algorithm>
56
57
58//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
59
61 nDataSetList(0), verboseLevel(0),fastPathFlags(),fastPathParams(),
62 counters(),fastPathCache()
63{
65 currentMaterial = elmMaterial = nullptr;
66 currentElement = nullptr; //ALB 14-Aug-2012 Coverity fix.
67 matParticle = elmParticle = nullptr;
68 matKinEnergy = elmKinEnergy = matCrossSection = elmCrossSection = 0.0;
69}
70
71//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
72
74{}
75
76
77//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
80 const G4Material* mat , G4bool requiresSlowPath)
81{
82 //The fast-path algorithm
83 // requiresSlowPath == true => Use slow path independently of other conditions
84 //A. Dotti: modifications to this algorithm following the studies of P. Ruth and R. Fowler
85 // on speeding up the cross-sections calculations. Their algorithm is called in the
86 // following "fast-path" while the normal approach to calcualte cross-section is
87 // referred to as "slow-path".
88 //Some details on the fast-path algorithm:
89 //The idea is to use a cached approximation of the material cross-section.
90 //Starting points:
91 //1- We need the material cross-section for navigation purposes: e.g. to calculate the PIL.
92 //2- If this interaction occurs at the end of a step we need to select the G4Element on which
93 // nucleus the interaction is actually happening. This is done calculating the element cross-section
94 // throwing a random number and selecting the appropriate nucleus (see SampleZandA function)
95 //3- To calculate the material cross section needed for 1- we use the G4Element cross sections.
96 // Since the material cross-section is simply the weighted sum of the element cross-sections.
97 //4- The slow path algorithm here accomplishes two use cases (not very good design IMHO): it
98 // calculates the material cross-section and it updates the xsecelem array that contains the
99 // relative cross-sections used by SampleZandA to select the element on which the interaction
100 // occurs.
101 //The idea of the fast-path algorithm is to replace the request for 1- with a faster calculation
102 //of the material cross-section from an approximation contained in cache.
103 //There two complications to take into account:
104 //A- If I use a fast parametrization for material cross-section, I still need to do the full
105 // calculations if an interaction occurs, because I need to calculate the element cross-sections.
106 // Since the function that updates xsecelem is the same (this one) I need to be sure that
107 // if I call this method for SampleAandZ the xsecelem is updated.
108 //B- It exists the possibility to be even fast the the fast-path algorithm: this happens when
109 // to select the element of the interaction via SampleZandI I call again this method exactly
110 // with the same conditions as when I called this method to calculate the material cross-section.
111 // In such a case xsecelem is updated and we do not need to do much more. This happens when
112 // for example a neutron undergoes an interaction at the end of the step.
113 //Dealing with B- complicates a bit the algorithm.
114 //In summary:
115 // If no fast-path algo is available (or user does not want that), go with the old plain algorithm
116 // If a fast-path algo is avilable, use it whenever possible.
117 //
118 // In general we expect user to selectively decide for which processes, materials and particle combinations
119 // we want to use the fast-path. If this is activated we also expect that the cross-section fast-path
120 // cache is created during the run initialization via calls to this method.
121 //
122 //fastPathFlags contains control flags for the fast-path algorithm:
123 // .prevCalcUsedFastPath == true => Previous call to GetCrossSection used the fast-path
124 // it is used in the decision to assess if xsecelem is correctly set-up
125 // .useFastPathIfAvailable == true => User requested the use of fast-path algorithm
126 // .initializationPhase == true => If true we are in Geant4 Init phase before the event-loop
127
128 //Check user-request, does he want fast-path? if not
129 // OR
130 // we want fast-path and we are in initialization phase?
131 if ( !fastPathFlags.useFastPathIfAvailable
132 || (fastPathFlags.useFastPathIfAvailable&&fastPathFlags.initializationPhase) ) {
133 //Traditional algorithm is requested
134 requiresSlowPath=true;
135 }
136
137 //Logging for performance calculations and counter, active only in FPDEBUG mode
138 counters.MethodCalled();
139 //Measure number of cycles
141
142 //This is the cache entry of the fast-path cross-section parametrization
144 //Did user request fast-path in first place and are we not in the initialization phase
145 if ( fastPathFlags.useFastPathIfAvailable && !fastPathFlags.initializationPhase ) {
146 //Important: if it is in initialization phase we should NOT use fast path: we are going to build it
147 //G4FastPathHadronicCrossSection::G4CrossSectionDataStore_Key searchkey = {part->GetParticleDefinition(),mat};
148 entry = fastPathCache[{part->GetParticleDefinition(),mat}];
149 }
150
151 //Super-fast-path: are we calling again this method for exactly the same conditions
152 //of the triplet {particle,material,energy}?
153 if(mat == currentMaterial && part->GetDefinition() == matParticle
154 && part->GetKineticEnergy() == matKinEnergy)
155 {
156 G4FastPathHadronicCrossSection::logInvocationTriedOneLine(entry);
157 //If there is no user-request for the fast-path in first place?
158 //It means we built the xsecelem for sure, let's return immediately
159 if ( !fastPathFlags.useFastPathIfAvailable ) {
160 return matCrossSection;
161 } else {
162 //Check that the last time we called this method we used the slow
163 //path: we need the data-member xsecelem to be setup correctly for the current
164 //interaction. This is ensured only if: we will do the slow path right now or we
165 //did it exactly for the same conditions of the last call.
166 if ( !fastPathFlags.prevCalcUsedFastPath && ! requiresSlowPath ) {
167 counters.HitOneLine();
168 G4FastPathHadronicCrossSection::logInvocationOneLine(entry);
169 //Good everything is setup correctly, exit!
170 return matCrossSection;
171 } else {
172 //We need to follow the slow-path because
173 //xsecelem is not calculated correctly
174 requiresSlowPath = true;
175 }
176 }
177 }
178
179 //Ok, now check if we have cached for this {particle,material,energy} the cross-section
180 //in this case let's return immediately, if we are not forced to take the slow path
181 //(e.g. as before if the xsecelem is not up-to-date we need to take the slow-path).
182 //Note that this is not equivalent to the previous ultra-fast check: we now have a map here
183 //So we can have for example a different particle.
184 if ( entry != nullptr && entry->energy == part->GetKineticEnergy() ) {
185 G4FastPathHadronicCrossSection::logHit(entry);
186 if ( !requiresSlowPath ) {
187 return entry->crossSection;
188 }
189 }
190
191 currentMaterial = mat;
192 matParticle = part->GetDefinition();
193 matKinEnergy = part->GetKineticEnergy();
194 matCrossSection = 0;
195
196 //Now check if the cache entry has a fast-path cross-section calculation available
198 if ( entry != nullptr && ! requiresSlowPath ) {
199 fast_entry = entry->fastPath;
200 assert(fast_entry!=nullptr && !fastPathFlags.initializationPhase);
201 }
202
203 //Each fast-path cross-section has a minimum value of validity, if energy is below
204 //that skip fast-path algorithm
205 if ( fast_entry != nullptr && part->GetKineticEnergy() < fast_entry->min_cutoff )
206 {
207 assert(requiresSlowPath==false);
208 requiresSlowPath = true;
209 }
210
211 //Ready to use the fast-path calculation
212 if ( !requiresSlowPath && fast_entry != nullptr ) {
213 counters.FastPath();
214 //Retrieve cross-section from fast-path cache
215 matCrossSection = fast_entry->GetCrossSection(part->GetKineticEnergy());
216 fastPathFlags.prevCalcUsedFastPath=true;
217 } else {
218 counters.SlowPath();
219 //Remember that we are now doing the full calculation: xsecelem will
220 //be made valid
221 fastPathFlags.prevCalcUsedFastPath=false;
222
223 G4int nElements = mat->GetNumberOfElements();
224 const G4double* nAtomsPerVolume = mat->GetVecNbOfAtomsPerVolume();
225
226 if(G4int(xsecelm.size()) < nElements) { xsecelm.resize(nElements); }
227
228 for(G4int i=0; i<nElements; ++i) {
229 matCrossSection += nAtomsPerVolume[i] *
230 GetCrossSection(part, (*mat->GetElementVector())[i], mat);
231 xsecelm[i] = matCrossSection;
232 }
233 }
234 //Stop measurement of cpu cycles
236
237 if ( entry != nullptr ) {
238 entry->energy = part->GetKineticEnergy();
239 entry->crossSection = matCrossSection;
240 }
241 //Some logging of timing
242 G4FastPathHadronicCrossSection::logTiming(entry,fast_entry,timing);
243 return matCrossSection;
244}
245
246//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
247
248void
250{
251 const G4FastPathHadronicCrossSection::cycleCountEntry* entry = fastPathCache[{pd,mat}];
252 if ( entry != nullptr ) {
253 if ( entry->fastPath != nullptr ) {
254 os<<*entry->fastPath;
255 } else {
256 os<<"#Cache entry for {"<<(pd!=nullptr?pd->GetParticleName():"UNDEFINED")<<",";
257 os<<(mat!=nullptr?mat->GetName():"UNDEFINED")<<"} found, but no fast path defined";
258 }
259 } else {
260 os<<"#Cache entry for {"<<(pd!=nullptr?pd->GetParticleName():"UNDEFINED")<<",";
261 os<<(mat!=nullptr?mat->GetName():"UNDEFINED")<<"} not found.";
262 }
263}
264
265//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
266
269 const G4Material* mat)
270{
271 if(mat == currentMaterial && part->GetDefinition() == matParticle
272 && part->GetKineticEnergy() == matKinEnergy) {
273 return matCrossSection;
274 }
275 currentMaterial = mat;
276 matParticle = part->GetDefinition();
277 matKinEnergy = part->GetKineticEnergy();
278 matCrossSection = 0.0;
279
280 size_t nElements = mat->GetNumberOfElements();
281 const G4double* nAtomsPerVolume = mat->GetVecNbOfAtomsPerVolume();
282
283 if(xsecelm.size() < nElements) { xsecelm.resize(nElements); }
284
285 for(size_t i=0; i<nElements; ++i) {
286 matCrossSection += nAtomsPerVolume[i] *
287 GetCrossSection(part, mat->GetElement(i), mat);
288 xsecelm[i] = matCrossSection;
289 }
290 return matCrossSection;
291}
292
293//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
294
297 const G4Element* elm,
298 const G4Material* mat)
299{
300 if(mat == elmMaterial && elm == currentElement &&
301 part->GetDefinition() == elmParticle &&
302 part->GetKineticEnergy() == elmKinEnergy)
303 { return elmCrossSection; }
304
305 elmMaterial = mat;
306 currentElement = elm;
307 elmParticle = part->GetDefinition();
308 elmKinEnergy = part->GetKineticEnergy();
309 elmCrossSection = 0.0;
310
311 G4int i = nDataSetList-1;
312 G4int Z = elm->GetZasInt();
313 if (elm->GetNaturalAbundanceFlag() &&
314 dataSetList[i]->IsElementApplicable(part, Z, mat)) {
315
316 // element wise cross section
317 elmCrossSection = dataSetList[i]->GetElementCrossSection(part, Z, mat);
318
319 //G4cout << "Element wise " << elmParticle->GetParticleName()
320 // << " xsec(barn)= " << elmCrossSection/barn
321 // << " E(MeV)= " << elmKinEnergy/MeV
322 // << " Z= " << Z << " AbundFlag= " << elm->GetNaturalAbandancesFlag()
323 // <<G4endl;
324
325 } else {
326 // isotope wise cross section
327 size_t nIso = elm->GetNumberOfIsotopes();
328
329 // user-defined isotope abundances
330 const G4double* abundVector = elm->GetRelativeAbundanceVector();
331
332 for (size_t j=0; j<nIso; ++j) {
333 if(abundVector[j] > 0.0) {
334 const G4Isotope* iso = elm->GetIsotope(j);
335 elmCrossSection += abundVector[j]*
336 GetIsoCrossSection(part, Z, iso->GetN(), iso, elm, mat, i);
337 //G4cout << "Isotope wise " << elmParticle->GetParticleName()
338 // << " xsec(barn)= " << elmCrossSection/barn
339 // << " E(MeV)= " << elmKinEnergy/MeV
340 // << " Z= " << Z << " A= " << iso->GetN() << " j= " << j << G4endl;
341 }
342 }
343 }
344 //G4cout << " E(MeV)= " << elmKinEnergy/MeV
345 // << "xsec(barn)= " << elmCrossSection/barn <<G4endl;
346 return elmCrossSection;
347}
348
349//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
350
352G4CrossSectionDataStore::GetIsoCrossSection(const G4DynamicParticle* part,
353 G4int Z, G4int A,
354 const G4Isotope* iso,
355 const G4Element* elm,
356 const G4Material* mat,
357 G4int idx)
358{
359 // this methods is called after the check that dataSetList[idx]
360 // depend on isotopes, so for this DataSet only isotopes are checked
361
362 // isotope-wise cross section does exist
363 if(dataSetList[idx]->IsIsoApplicable(part, Z, A, elm, mat) ) {
364 return dataSetList[idx]->GetIsoCrossSection(part, Z, A, iso, elm, mat);
365
366 } else {
367 // seach for other dataSet
368 for (G4int j = nDataSetList-1; j >= 0; --j) {
369 if (dataSetList[j]->IsElementApplicable(part, Z, mat)) {
370 return dataSetList[j]->GetElementCrossSection(part, Z, mat);
371 } else if (dataSetList[j]->IsIsoApplicable(part, Z, A, elm, mat)) {
372 return dataSetList[j]->GetIsoCrossSection(part, Z, A, iso, elm, mat);
373 }
374 }
375 }
377 ed << "No isotope cross section found for "
378 << part->GetDefinition()->GetParticleName()
379 << " off Element " << elm->GetName()
380 << " in " << mat->GetName() << " Z= " << Z << " A= " << A
381 << " E(MeV)= " << part->GetKineticEnergy()/MeV << G4endl;
382 G4Exception("G4CrossSectionDataStore::GetIsoCrossSection", "had001",
383 FatalException, ed);
384 return 0.0;
385}
386
387//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
388
391 G4int Z, G4int A,
392 const G4Isotope* iso,
393 const G4Element* elm,
394 const G4Material* mat)
395{
396 for (G4int i = nDataSetList-1; i >= 0; --i) {
397 if (dataSetList[i]->IsIsoApplicable(part, Z, A, elm, mat) ) {
398 return dataSetList[i]->GetIsoCrossSection(part, Z, A, iso, elm, mat);
399 }
400 }
402 ed << "No isotope cross section found for "
403 << part->GetDefinition()->GetParticleName()
404 << " off Element " << elm->GetName()
405 << " in " << mat->GetName() << " Z= " << Z << " A= " << A
406 << " E(MeV)= " << part->GetKineticEnergy()/MeV << G4endl;
407 G4Exception("G4CrossSectionDataStore::GetCrossSection", "had001",
408 FatalException, ed);
409 return 0.0;
410}
411
412//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
413
414const G4Element*
416 const G4Material* mat,
417 G4Nucleus& target)
418{
419 size_t nElements = mat->GetNumberOfElements();
420 const G4Element* anElement = mat->GetElement(0);
421
422 // select element from a compound
423 if(1 < nElements) {
424 G4double cross = matCrossSection*G4UniformRand();
425 for(size_t i=0; i<nElements; ++i) {
426 if(cross <= xsecelm[i]) {
427 anElement = mat->GetElement(i);
428 break;
429 }
430 }
431 }
432
433 G4int Z = anElement->GetZasInt();
434 const G4Isotope* iso = nullptr;
435
436 G4int i = nDataSetList-1;
437 if (dataSetList[i]->IsElementApplicable(part, Z, mat)) {
438
439 //----------------------------------------------------------------
440 // element-wise cross section
441 // isotope cross section is not computed
442 //----------------------------------------------------------------
443 size_t nIso = anElement->GetNumberOfIsotopes();
444 iso = anElement->GetIsotope(0);
445
446 // more than 1 isotope
447 if(1 < nIso) {
448 iso = dataSetList[i]->SelectIsotope(anElement,
449 part->GetKineticEnergy(),
450 part->GetLogKineticEnergy());
451 }
452 } else {
453
454 //----------------------------------------------------------------
455 // isotope-wise cross section
456 // isotope cross section is computed
457 //----------------------------------------------------------------
458 size_t nIso = anElement->GetNumberOfIsotopes();
459 iso = anElement->GetIsotope(0);
460
461 // more than 1 isotope
462 if(1 < nIso) {
463 const G4double* abundVector = anElement->GetRelativeAbundanceVector();
464 if(xseciso.size() < nIso) { xseciso.resize(nIso); }
465
466 G4double cross = 0.0;
467 size_t j;
468 for (j = 0; j<nIso; ++j) {
469 G4double xsec = 0.0;
470 if(abundVector[j] > 0.0) {
471 iso = anElement->GetIsotope(j);
472 xsec = abundVector[j]*
473 GetIsoCrossSection(part, Z, iso->GetN(), iso, anElement, mat, i);
474 }
475 cross += xsec;
476 xseciso[j] = cross;
477 }
478 cross *= G4UniformRand();
479 for (j = 0; j<nIso; ++j) {
480 if(cross <= xseciso[j]) {
481 iso = anElement->GetIsotope(j);
482 break;
483 }
484 }
485 }
486 }
487 target.SetIsotope(iso);
488 return anElement;
489}
490
491//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
492
493void
495{
496 if (nDataSetList == 0) {
498 ed << "No cross section is registered for "
499 << aParticleType.GetParticleName() << G4endl;
500 G4Exception("G4CrossSectionDataStore::BuildPhysicsTable", "had001",
501 FatalException, ed);
502 return;
503 }
504 for (G4int i=0; i<nDataSetList; ++i) {
505 dataSetList[i]->BuildPhysicsTable(aParticleType);
506 }
507 //A.Dotti: if fast-path has been requested we can now create the surrogate
508 // model for fast path.
509 if ( fastPathFlags.useFastPathIfAvailable ) {
510 fastPathFlags.initializationPhase = true;
511 using my_value_type=G4FastPathHadronicCrossSection::G4CrossSectionDataStore_Requests::value_type;
512 //Loop on all requests, if particle matches create the corresponding fsat-path
513 std::for_each( requests.begin() , requests.end() ,
514 [&aParticleType,this](const my_value_type& req) {
515 if ( aParticleType == *req.part_mat.first ) {
516 G4FastPathHadronicCrossSection::cycleCountEntry* entry =
517 new G4FastPathHadronicCrossSection::cycleCountEntry(aParticleType.GetParticleName(),req.part_mat.second);
518 entry->fastPath =
519 new G4FastPathHadronicCrossSection::fastPathEntry(&aParticleType,req.part_mat.second,req.min_cutoff);
520 entry->fastPath->Initialize(this);
521 fastPathCache[req.part_mat] = entry;
522 }
523 }
524 );
525 fastPathFlags.initializationPhase = false;
526 }
527}
528
529//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
530
532 const G4Material* mat, G4double min_cutoff)
533{
534 assert(pdef!=nullptr&&mat!=nullptr);
536 if ( requests.insert( { key , min_cutoff } ).second ) {
538 ed << "Attempting to request FastPath for couple: <"
539 << pdef->GetParticleName() << ", " <<mat->GetName()
540 << "> but combination already exists" << G4endl;
541 G4Exception("G4CrossSectionDataStore::ActivateFastPath", "had001",
542 FatalException, ed);
543 }
544}
545
546//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
547
548void
550{
551 // Print out all cross section data sets used and the energies at
552 // which they apply
553
554 if (nDataSetList == 0) {
555 G4cout << "WARNING - G4CrossSectionDataStore::DumpPhysicsTable: "
556 << " no data sets registered" << G4endl;
557 return;
558 }
559
560 for (G4int i = nDataSetList-1; i >= 0; --i) {
561 G4double e1 = dataSetList[i]->GetMinKinEnergy();
562 G4double e2 = dataSetList[i]->GetMaxKinEnergy();
563 G4cout
564 << " Cr_sctns: " << std::setw(25) << dataSetList[i]->GetName() << ": "
565 << G4BestUnit(e1, "Energy")
566 << " ---> "
567 << G4BestUnit(e2, "Energy") << "\n";
568 if (dataSetList[i]->GetName() == "G4CrossSectionPairGG") {
569 dataSetList[i]->DumpPhysicsTable(aParticleType);
570 }
571 }
572}
573
574//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
575#include <typeinfo>
577 std::ofstream& outFile) const
578{
579 // Write cross section data set info to html physics list
580 // documentation page
581
582 G4double ehi = 0;
583 G4double elo = 0;
584 G4String physListName(std::getenv("G4PhysListName"));
585 for (G4int i = nDataSetList-1; i > 0; i--) {
586 elo = dataSetList[i]->GetMinKinEnergy()/GeV;
587 ehi = dataSetList[i]->GetMaxKinEnergy()/GeV;
588 outFile << " <li><b><a href=\"" << physListName << "_"
589 << dataSetList[i]->GetName() << ".html\"> "
590 << dataSetList[i]->GetName() << "</a> from "
591 << elo << " GeV to " << ehi << " GeV </b></li>\n";
592 //G4cerr << i << ": XS for " << pD.GetParticleName() << " : " << dataSetList[i]->GetName()
593 // << " typeid : " << typeid(dataSetList[i]).name()<< G4endl;
594 PrintCrossSectionHtml(dataSetList[i]);
595 }
596
597 G4double defaultHi = dataSetList[0]->GetMaxKinEnergy()/GeV;
598 if (ehi < defaultHi) {
599 outFile << " <li><b><a href=\"" << dataSetList[0]->GetName() << ".html\"> "
600 << dataSetList[0]->GetName() << "</a> from "
601 << ehi << " GeV to " << defaultHi << " GeV </b></li>\n";
602 PrintCrossSectionHtml(dataSetList[0]);
603 }
604}
605
606//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
607
609{
610 G4String dirName(std::getenv("G4PhysListDocDir"));
611 G4String physListName(std::getenv("G4PhysListName"));
612
613 G4String pathName = dirName + "/" + physListName + "_" + HtmlFileName(cs->GetName());
614 std::ofstream outCS;
615 outCS.open(pathName);
616 outCS << "<html>\n";
617 outCS << "<head>\n";
618 outCS << "<title>Description of " << cs->GetName()
619 << "</title>\n";
620 outCS << "</head>\n";
621 outCS << "<body>\n";
622
623 cs->CrossSectionDescription(outCS);
624
625 outCS << "</body>\n";
626 outCS << "</html>\n";
627
628}
629//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
630
631G4String G4CrossSectionDataStore::HtmlFileName(const G4String & in) const
632{
633 G4String str(in);
634 // replace blanks by _ C++11 version:
635 std::transform(str.begin(), str.end(), str.begin(), [](char ch) {
636 return ch == ' ' ? '_' : ch;
637 });
638 str=str + ".html";
639 return str;
640}
641
642//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
643
645{
646 if(p->ForAllAtomsAndEnergies()) {
647 dataSetList.clear();
648 nDataSetList = 0;
649 }
650 dataSetList.push_back(p);
651 ++nDataSetList;
652}
653
654
655
656//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
657
659{
660 if(p->ForAllAtomsAndEnergies()) {
661 dataSetList.clear();
662 dataSetList.push_back(p);
663 nDataSetList = 1;
664 } else {
665 if ( i > dataSetList.size() ) i = dataSetList.size();
666 std::vector< G4VCrossSectionDataSet* >::iterator it = dataSetList.end() - i;
667 dataSetList.insert(it , p);
668 ++nDataSetList;
669 }
670}
671//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
double A(double temperature)
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
#define G4BestUnit(a, b)
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
#define G4UniformRand()
Definition: Randomize.hh:52
void DumpHtml(const G4ParticleDefinition &, std::ofstream &) const
void BuildPhysicsTable(const G4ParticleDefinition &)
void AddDataSet(G4VCrossSectionDataSet *)
void PrintCrossSectionHtml(const G4VCrossSectionDataSet *cs) const
void DumpPhysicsTable(const G4ParticleDefinition &)
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Material *)
void DumpFastPath(const G4ParticleDefinition *, const G4Material *, std::ostream &os)
void ActivateFastPath(const G4ParticleDefinition *, const G4Material *, G4double)
G4double GetCrossSection(const G4DynamicParticle *, const G4Material *)
const G4Element * SampleZandA(const G4DynamicParticle *, const G4Material *, G4Nucleus &target)
G4double GetLogKineticEnergy() const
const G4ParticleDefinition * GetParticleDefinition() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:166
G4bool GetNaturalAbundanceFlag() const
Definition: G4Element.hh:261
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:169
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:158
const G4String & GetName() const
Definition: G4Element.hh:126
G4int GetZasInt() const
Definition: G4Element.hh:131
G4int GetN() const
Definition: G4Isotope.hh:93
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:200
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
const G4String & GetName() const
Definition: G4Material.hh:175
static G4NistManager * Instance()
void SetIsotope(const G4Isotope *iso)
Definition: G4Nucleus.hh:122
const G4String & GetParticleName() const
const G4String & GetName() const
virtual void CrossSectionDescription(std::ostream &) const
std::pair< const G4ParticleDefinition *, const G4Material * > G4CrossSectionDataStore_Key