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