Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4DNACPA100ElasticModel.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 // CPA100 elastic model class for electrons
27 //
28 // Based on the work of M. Terrissol and M. C. Bordage
29 //
30 // Users are requested to cite the following papers:
31 // - M. Terrissol, A. Baudre, Radiat. Prot. Dosim. 31 (1990) 175-177
32 // - M.C. Bordage, J. Bordes, S. Edel, M. Terrissol, X. Franceries,
33 // M. Bardies, N. Lampe, S. Incerti, Phys. Med. 32 (2016) 1833-1840
34 //
35 // Authors of this class:
36 // M.C. Bordage, M. Terrissol, S. Edel, J. Bordes, S. Incerti
37 //
38 // 15.01.2014: creation
39 //
40 
42 #include "G4PhysicalConstants.hh"
43 #include "G4SystemOfUnits.hh"
45 
46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 
48 using namespace std;
49 
50 // #define CPA100_VERBOSE
51 
52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
53 
55  const G4String& nam)
56 :G4VEmModel(nam),isInitialised(false)
57 {
58 
59  //killBelowEnergy = 11. * eV; // Default
60  //killBelowEnergy = 10.481 * eV;
61  //lowEnergyLimit = 11 * eV;
62  //highEnergyLimit = 255955 * eV;
64  SetHighEnergyLimit(255955 * eV);
65 
66  verboseLevel= 0;
67  // Verbosity scale:
68  // 0 = nothing
69  // 1 = warning for energy non-conservation
70  // 2 = details of energy budget
71  // 3 = calculation of cross sections, file openings, sampling of atoms
72  // 4 = entering in methods
73 
74 #ifdef UEHARA_VERBOSE
75  if( verboseLevel>0 )
76  {
77  G4cout << "CPA100 Elastic model is constructed " << G4endl
78  << "Energy range: "
79  << lowEnergyLimit / eV << " eV - "
80  << highEnergyLimit / keV << " keV"
81  << G4endl;
82  }
83 #endif
84 
87 
88  // Selection of stationary mode
89 
90  statCode = false;
91 }
92 
93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
94 
96 {
97  // For total cross section
98 
99  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
100  for (pos = tableData.begin(); pos != tableData.end(); ++pos)
101  {
102  G4DNACrossSectionDataSet* table = pos->second;
103  delete table;
104  }
105 
106  // For final state
107 
108  eVecm.clear();
109 
110 }
111 
112 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
113 
115  particle,
116  const G4DataVector& /*cuts*/)
117 {
118 
119 #ifdef UEHARA_VERBOSE
120  if (verboseLevel > 3)
121  G4cout << "Calling G4DNACPA100ElasticModel::Initialise()" << G4endl;
122 #endif
123 
124  if(particle->GetParticleName() != "e-")
125  {
126  G4Exception("*** WARNING: the G4DNACPA100ElasticModel is "
127  "not intented to be used with another particle than the electron",
128  "",FatalException,"") ;
129  }
130 
131  // Energy limits
132 
133  if (LowEnergyLimit() < 11.*eV)
134  {
135  G4cout << "G4DNACPA100ElasticModel: low energy limit increased from " <<
136  LowEnergyLimit()/eV << " eV to " << 11 << " eV" << G4endl;
137  SetLowEnergyLimit(11.*eV);
138  }
139 
140  if (HighEnergyLimit() > 255955.*eV)
141  {
142  G4cout << "G4DNACPA100ElasticModel: high energy limit decreased from " <<
143  HighEnergyLimit()/keV << " keV to " << 255.955 << " keV"
144  << G4endl;
145  SetHighEnergyLimit(255955.*eV);
146  }
147 
148  // Reading of data files
149 
150  G4double scaleFactor = 1e-20*m*m;
151 
152  G4String fileElectron("dna/sigma_elastic_e_cpa100");
153 
156 
157  // *** ELECTRON
158 
159  // For total cross section
160 
161  electron = electronDef->GetParticleName();
162 
163  tableFile[electron] = fileElectron;
164 
165  G4DNACrossSectionDataSet* tableE =
167  eV,scaleFactor );
168 
169  /*
170  G4DNACrossSectionDataSet* tableE =
171  new G4DNACrossSectionDataSet(new G4DNACPA100LogLogInterpolation,
172  eV,scaleFactor );
173  */
174 
175  tableE->LoadData(fileElectron);
176 
177  tableData[electron] = tableE;
178 
179  // For final state
180 
181  char *path = getenv("G4LEDATA");
182 
183  if (!path)
184  {
185  G4Exception("G4DNACPA100ElasticModel::Initialise","em0006",
186  FatalException,"G4LEDATA environment variable not set.");
187  return;
188  }
189 
190  std::ostringstream eFullFileName;
191 
192  eFullFileName << path
193  << "/dna/sigmadiff_cumulated_elastic_e_cpa100.dat";
194  std::ifstream eDiffCrossSection(eFullFileName.str().c_str());
195 
196  if (!eDiffCrossSection)
197  G4Exception("G4DNACPA100ElasticModel::Initialise","em0003",
199  "Missing data file:/dna/sigmadiff_cumulated_elastic_e_cpa100.dat");
200 
201  // March 25th, 2014 - Vaclav Stepan, Sebastien Incerti
202  // Added clear for MT
203 
204  eTdummyVec.clear();
205  eVecm.clear();
206  eDiffCrossSectionData.clear();
207 
208  //
209 
210  eTdummyVec.push_back(0.);
211 
212  while(!eDiffCrossSection.eof())
213  {
214  G4double tDummy;
215  G4double eDummy;
216  eDiffCrossSection>>tDummy>>eDummy;
217 
218  // SI : mandatory eVecm initialization
219 
220  if (tDummy != eTdummyVec.back())
221  {
222  eTdummyVec.push_back(tDummy);
223  eVecm[tDummy].push_back(0.);
224  }
225 
226  eDiffCrossSection>>eDiffCrossSectionData[tDummy][eDummy];
227 
228  if (eDummy != eVecm[tDummy].back()) eVecm[tDummy].push_back(eDummy);
229 
230  }
231 
232  // End final state
233 
234 #ifdef UEHARA_VERBOSE
235  if (verboseLevel > 2)
236  G4cout << "Loaded cross section files for CPA100 Elastic model" << G4endl;
237 #endif
238 
239 #ifdef UEHARA_VERBOSE
240  if( verboseLevel>0 )
241  {
242  G4cout << "CPA100 Elastic model is initialized " << G4endl
243  << "Energy range: "
244  << LowEnergyLimit() / eV << " eV - "
245  << HighEnergyLimit() / keV << " keV"
246  << G4endl;
247  }
248 #endif
249 
250  // Initialize water density pointer
253 
254  if (isInitialised) { return; }
256  isInitialised = true;
257 
258 }
259 
260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
261 
263 (const G4Material* material,
264  const G4ParticleDefinition* p,
265  G4double ekin,
266  G4double,
267  G4double)
268 {
269 
270 #ifdef UEHARA_VERBOSE
271  if (verboseLevel > 3)
272  G4cout <<
273  "Calling CrossSectionPerVolume() of G4DNACPA100ElasticModel" << G4endl;
274 #endif
275 
276  // Calculate total cross section for model
277 
278  G4double sigma=0;
279 
280  G4double waterDensity = (*fpMolWaterDensity)[material->GetIndex()];
281 
282  if(waterDensity!= 0.0)
283  {
284  const G4String& particleName = p->GetParticleName();
285 
286  if (ekin < HighEnergyLimit() && ekin >= LowEnergyLimit())
287  {
288  //SI : XS must not be zero otherwise sampling of secondaries
289  // method ignored
290 
291  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
292  pos = tableData.find(particleName);
293 
294  if (pos != tableData.end())
295  {
296  G4DNACrossSectionDataSet* table = pos->second;
297  if (table != 0)
298  {
299  sigma = table->FindValue(ekin);
300 
301  //
302  //Dump in non-MT mode
303  //
304  /*
305  G4double minEnergy = 10.481 * eV;
306  G4double maxEnergy = 255955. * eV;
307  G4int nEnergySteps = 1000;
308  G4double energy(minEnergy);
309  G4double
310  stpEnergy(std::pow(maxEnergy/energy,
311  1./static_cast<G4double>(nEnergySteps-1)));
312  G4int step(nEnergySteps);
313  system ("rm -rf elastic-cpa100.out");
314  FILE* myFile=fopen("elastic-cpa100.out","a");
315  while (step>0)
316  {
317  step--;
318  fprintf (myFile,"%16.9le %16.9le\n",
319  energy/eV,
320  table->FindValue(energy)/(1e-20*m*m));
321  energy*=stpEnergy;
322  }
323  fclose (myFile);
324  abort();
325  */
326  //
327  // end of dump
328  //
329 
330 
331  }
332  }
333  else
334  {
335  G4Exception("G4DNACPA100ElasticModel::ComputeCrossSectionPerVolume",
336  "em0002",
337  FatalException,"Model not applicable to particle type.");
338  }
339 
340  }
341 
342 #ifdef UEHARA_VERBOSE
343  if (verboseLevel > 2)
344  {
345  G4cout << "__________________________________" << G4endl;
346  G4cout << "G4DNACPA100ElasticModel - XS INFO START" << G4endl;
347  G4cout << "Kinetic energy(eV)=" << ekin/eV << " particle : " << particleName << G4endl;
348  G4cout << "Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl;
349  G4cout << "Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl;
350  // G4cout << " - Cross section per water molecule (cm^-1)="
351  // << sigma*material->GetAtomicNumDensityVector()[1]/(1./cm) << G4endl;
352  G4cout << "G4DNACPA100ElasticModel - XS INFO END" << G4endl;
353  }
354 #endif
355 
356  }
357 
358  return sigma*waterDensity;
359 }
360 
361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
362 
363 void G4DNACPA100ElasticModel::SampleSecondaries(std::vector<G4DynamicParticle*>* /*fvect*/,
364  const G4MaterialCutsCouple* /*couple*/,
365  const G4DynamicParticle* aDynamicElectron,
366  G4double,
367  G4double)
368 {
369 #ifdef UEHARA_VERBOSE
370  if (verboseLevel > 3)
371  G4cout << "Calling SampleSecondaries() of G4DNACPA100ElasticModel" << G4endl;
372 #endif
373 
374  G4double electronEnergy0 = aDynamicElectron->GetKineticEnergy();
375 
376  {
377  G4double cosTheta = RandomizeCosTheta(electronEnergy0);
378  G4double phi = 2. * pi * G4UniformRand();
379 
380  G4ThreeVector zVers = aDynamicElectron->GetMomentumDirection();
381 
382  //G4ThreeVector xVers = zVers.orthogonal();
383  //G4ThreeVector yVers = zVers.cross(xVers);
384  //G4double xDir = std::sqrt(1. - cosTheta*cosTheta);
385  //G4double yDir = xDir;
386  //xDir *= std::cos(phi);
387  //yDir *= std::sin(phi);
388 
389  // Computation of scattering angles (from Subroutine DIRAN in CPA100)
390 
391  G4double CT1, ST1, CF1, SF1, CT2, ST2, CF2, SF2;
392  G4double sinTheta = std::sqrt (1-cosTheta*cosTheta);
393 
394  CT1=0;
395  ST1=0;
396  CF1=0;
397  SF1=0;
398  CT2=0;
399  ST2=0;
400  CF2=0;
401  SF2=0;
402 
403  CT1 = zVers.z();
404  ST1=std::sqrt(1.-CT1*CT1);
405 
406  if (ST1!=0) CF1 = zVers.x()/ST1; else CF1 = std::cos(2. * pi * G4UniformRand());
407  if (ST1!=0) SF1 = zVers.y()/ST1; else SF1 = std::sqrt(1.-CF1*CF1);
408 
409  G4double A3, A4, A5, A2, A1;
410  A3=0;
411  A4=0;
412  A5=0;
413  A2=0;
414  A1=0;
415 
416  A3 = sinTheta*std::cos(phi);
417  A4 = A3*CT1 + ST1*cosTheta;
418  A5 = sinTheta * std::sin(phi);
419  A2 = A4 * SF1 + A5 * CF1;
420  A1 = A4 * CF1 - A5 * SF1;
421 
422  CT2 = CT1*cosTheta - ST1*A3;
423  ST2 = std::sqrt(1.-CT2*CT2);
424 
425  if (ST2==0) ST2=1E-6;
426  CF2 = A1/ST2;
427  SF2 = A2/ST2;
428 
429  /*
430  G4cout << "CT1=" << CT1 << G4endl;
431  G4cout << "ST1=" << ST1 << G4endl;
432  G4cout << "CF1=" << CF1 << G4endl;
433  G4cout << "SF1=" << SF1 << G4endl;
434  G4cout << "cosTheta=" << cosTheta << G4endl;
435  G4cout << "sinTheta=" << sinTheta << G4endl;
436  G4cout << "cosPhi=" << std::cos(phi) << G4endl;
437  G4cout << "sinPhi=" << std::sin(phi) << G4endl;
438  G4cout << "CT2=" << CT2 << G4endl;
439  G4cout << "ST2=" << ST2 << G4endl;
440  G4cout << "CF2=" << CF2 << G4endl;
441  G4cout << "SF2=" << SF2 << G4endl;
442  */
443 
444  G4ThreeVector zPrimeVers(ST2*CF2,ST2*SF2,CT2);
445 
446  //
447 
449 
450  //fParticleChangeForGamma->SetProposedKineticEnergy(electronEnergy0);
451 
452  if (!statCode)
453 
455  (electronEnergy0-1.214E-4*(1.-cosTheta)*electronEnergy0);
456 
457  else fParticleChangeForGamma->SetProposedKineticEnergy(electronEnergy0);
458 
459  fParticleChangeForGamma->ProposeLocalEnergyDeposit(1.214E-4*(1.-cosTheta)*electronEnergy0);
460 
461  }
462 
463 }
464 
465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
466 
468  (G4ParticleDefinition * particleDefinition, G4double k, G4double integrDiff)
469 {
470 
471  G4double theta = 0.;
472  G4double valueT1 = 0;
473  G4double valueT2 = 0;
474  G4double valueE21 = 0;
475  G4double valueE22 = 0;
476  G4double valueE12 = 0;
477  G4double valueE11 = 0;
478  G4double xs11 = 0;
479  G4double xs12 = 0;
480  G4double xs21 = 0;
481  G4double xs22 = 0;
482 
483  if (particleDefinition == G4Electron::ElectronDefinition())
484  {
485  std::vector<G4double>::iterator t2 = std::upper_bound(eTdummyVec.begin(),eTdummyVec.end(), k);
486  std::vector<G4double>::iterator t1 = t2-1;
487 
488  std::vector<G4double>::iterator e12 = std::upper_bound(eVecm[(*t1)].begin(),eVecm[(*t1)].end(),
489  integrDiff);
490  std::vector<G4double>::iterator e11 = e12-1;
491 
492  std::vector<G4double>::iterator e22 = std::upper_bound(eVecm[(*t2)].begin(),eVecm[(*t2)].end(),
493  integrDiff);
494  std::vector<G4double>::iterator e21 = e22-1;
495 
496  valueT1 =*t1;
497  valueT2 =*t2;
498  valueE21 =*e21;
499  valueE22 =*e22;
500  valueE12 =*e12;
501  valueE11 =*e11;
502 
503  xs11 = eDiffCrossSectionData[valueT1][valueE11];
504  xs12 = eDiffCrossSectionData[valueT1][valueE12];
505  xs21 = eDiffCrossSectionData[valueT2][valueE21];
506  xs22 = eDiffCrossSectionData[valueT2][valueE22];
507 
508  //FOR CPA100
509  //if(k==valueT1) xs22 = eDiffCrossSectionData[valueT1][valueE12];
510 
511  }
512 
513  if (xs11==0 && xs12==0 && xs21==0 && xs22==0) return (0.);
514 
515  // FOR CPA100
516 
517 
518  theta = QuadInterpolator ( valueE11, valueE12,
519  valueE21, valueE22,
520  xs11, xs12,
521  xs21, xs22,
522  valueT1, valueT2,
523  k, integrDiff );
524 
525  return theta;
526 
527 
528  //FOR CPA100
529  //return xs22;
530 }
531 
532 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
533 
535  G4double e2,
536  G4double e,
537  G4double xs1,
538  G4double xs2)
539 {
540  G4double d1 = std::log(xs1);
541  G4double d2 = std::log(xs2);
542  G4double value = std::exp(d1 + (d2 - d1)*(e - e1)/ (e2 - e1));
543  return value;
544 }
545 
546 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
547 
549  G4double e2,
550  G4double e,
551  G4double xs1,
552  G4double xs2)
553 {
554  G4double d1 = xs1;
555  G4double d2 = xs2;
556  G4double value = (d1 + (d2 - d1)*(e - e1)/ (e2 - e1));
557  return value;
558 }
559 
560 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
561 
563  G4double e2,
564  G4double e,
565  G4double xs1,
566  G4double xs2)
567 {
568  G4double a = (std::log10(xs2)-std::log10(xs1)) / (std::log10(e2)-std::log10(e1));
569  G4double b = std::log10(xs2) - a*std::log10(e2);
570  G4double sigma = a*std::log10(e) + b;
571  G4double value = (std::pow(10.,sigma));
572  return value;
573 }
574 
575 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
576 
577 
579  G4double e21, G4double e22,
580  G4double xs11, G4double xs12,
581  G4double xs21, G4double xs22,
583  G4double t, G4double e)
584 {
585  // Log-Log
586 /*
587  G4double interpolatedvalue1 = LogLogInterpolate(e11, e12, e, xs11, xs12);
588  G4double interpolatedvalue2 = LogLogInterpolate(e21, e22, e, xs21, xs22);
589  G4double value = LogLogInterpolate(t1, t2, t, interpolatedvalue1, interpolatedvalue2);
590 
591 
592  // Lin-Log
593  G4double interpolatedvalue1 = LinLogInterpolate(e11, e12, e, xs11, xs12);
594  G4double interpolatedvalue2 = LinLogInterpolate(e21, e22, e, xs21, xs22);
595  G4double value = LinLogInterpolate(t1, t2, t, interpolatedvalue1, interpolatedvalue2);
596 */
597 
598  // Lin-Lin
599  G4double interpolatedvalue1 = LinLinInterpolate(e11, e12, e, xs11, xs12);
600  G4double interpolatedvalue2 = LinLinInterpolate(e21, e22, e, xs21, xs22);
601  G4double value = LinLinInterpolate(t1, t2, t, interpolatedvalue1, interpolatedvalue2);
602 
603  return value;
604 }
605 
606 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
607 
609 {
610 
611  G4double integrdiff=0; // PROBABILITY between 0 and 1.
612  G4double uniformRand=G4UniformRand();
613  integrdiff = uniformRand;
614 
615  G4double cosTheta=0.;
616 
617  // 1 - COS THETA is read from the data file
618  cosTheta = 1 - Theta(G4Electron::ElectronDefinition(),k/eV,integrdiff);
619 
620  //
621  //
622  //Dump
623  //
624  //G4cout << "theta=" << theta << G4endl;
625  //G4cout << "cos theta=" << std::cos(theta*pi/180) << G4endl;
626  //G4cout << "sin theta=" << std::sin(theta*pi/180) << G4endl;
627  //G4cout << "acos(cos theta)=" << std::acos(cosTheta) << G4endl;
628  //G4cout << "cos theta="<< cosTheta << G4endl;
629  //G4cout << "1 - cos theta="<< 1. - cosTheta << G4endl;
630  //G4cout << "sin theta=" << std::sqrt(1-cosTheta*cosTheta) << G4endl;
631  //
632  /*
633  G4double minProb = 0; // we scan probability between 0 and one
634  G4double maxProb = 1;
635  G4int nProbSteps = 100;
636  G4double prob(minProb);
637  G4double stepProb((maxProb-minProb)/static_cast<G4double>(nProbSteps));
638  G4int step(nProbSteps);
639  system ("rm -rf elastic-cumul-cpa100-100keV.out");
640  FILE* myFile=fopen("elastic-cumul-cpa100-100keV.out","a");
641  while (step>=0)
642  {
643  step--;
644  fprintf (myFile,"%16.9le %16.9le\n",
645  prob,
646  Theta(G4Electron::ElectronDefinition(),100000,prob)); // SELECT NRJ IN eV !!!
647  prob=prob+stepProb;
648  }
649  fclose (myFile);
650  abort();
651  */
652  //
653  // end of dump
654  //
655 
656  return cosTheta;
657 }
virtual G4bool LoadData(const G4String &argFileName)
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
static const G4double pos
TTree * t1
Definition: plottest35.C:26
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:707
G4ParticleChangeForGamma * fParticleChangeForGamma
static constexpr double keV
Definition: G4SIunits.hh:216
size_t GetIndex() const
Definition: G4Material.hh:262
#define G4endl
Definition: G4ios.hh:61
const G4ThreeVector & GetMomentumDirection() const
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
double z() const
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:608
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:616
const std::vector< double > * GetNumMolPerVolTableFor(const G4Material *) const
Retrieve a table of molecular densities (number of molecules per unit volume) in the G4 unit system f...
G4double RandomizeCosTheta(G4double k)
static constexpr double m
Definition: G4SIunits.hh:129
double G4double
Definition: G4Types.hh:76
G4double LinLogInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)
static const G4double d2
const XML_Char int const XML_Char * value
Definition: expat.h:331
static G4DNAMolecularMaterial * Instance()
static const G4double d1
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
#define G4UniformRand()
Definition: Randomize.hh:53
virtual G4double CrossSectionPerVolume(const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
static constexpr double eV
Definition: G4SIunits.hh:215
G4double QuadInterpolator(G4double e11, G4double e12, G4double e21, G4double e22, G4double x11, G4double x12, G4double x21, G4double x22, G4double t1, G4double t2, G4double t, G4double e)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
Hep3Vector unit() const
TTree * t2
Definition: plottest35.C:36
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:714
std::vector< G4double > eTdummyVec
G4DNACPA100ElasticModel(const G4ParticleDefinition *p=0, const G4String &nam="DNACPA100ElasticModel")
const std::vector< G4double > * fpMolWaterDensity
G4double LinLinInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)
G4double GetKineticEnergy() const
virtual G4double FindValue(G4double e, G4int componentId=0) const
static constexpr double cm
Definition: G4SIunits.hh:119
G4GLOB_DLL std::ostream G4cout
double x() const
G4double Theta(G4ParticleDefinition *aParticleDefinition, G4double k, G4double integrDiff)
static constexpr double pi
Definition: G4SIunits.hh:75
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:131
double y() const
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:609
G4double LogLogInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)