Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4RToEConvForGamma.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 //
27 // $Id: G4RToEConvForGamma.cc 70745 2013-06-05 10:54:00Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------
31 // GEANT 4 class implementation file/ History:
32 // 5 Oct. 2002, H.Kuirashige : Structure created based on object model
33 // --------------------------------------------------------------
34 
35 #include "G4RToEConvForGamma.hh"
36 #include "G4ParticleDefinition.hh"
37 #include "G4ParticleTable.hh"
38 #include "G4Material.hh"
39 #include "G4PhysicsLogVector.hh"
40 
41 #include "G4ios.hh"
42 #include "G4SystemOfUnits.hh"
43 
46  Z(-1),
47  s200keV(0.), s1keV(0.),
48  tmin(0.), tlow(0.),
49  smin(0.), slow(0.),
50  cmin(0.), clow(0.), chigh(0.)
51 {
53  if (theParticle ==0) {
54 #ifdef G4VERBOSE
55  if (GetVerboseLevel()>0) {
56  G4cout << " G4RToEConvForGamma::G4RToEConvForGamma() ";
57  G4cout << " Gamma is not defined !!" << G4endl;
58  }
59 #endif
60  }
61 }
62 
64 {
65 }
66 
67 
68 // ***********************************************************************
69 // ******************* BuildAbsorptionLengthVector ***********************
70 // ***********************************************************************
72  const G4Material* aMaterial,
73  G4RangeVector* absorptionLengthVector )
74 {
75  // fill the absorption length vector for this material
76  // absorption length is defined here as
77  //
78  // absorption length = 5./ macroscopic absorption cross section
79  //
80  const G4CrossSectionTable* aCrossSectionTable = (G4CrossSectionTable*)(theLossTable);
81  const G4ElementVector* elementVector = aMaterial->GetElementVector();
82  const G4double* atomicNumDensityVector = aMaterial->GetAtomicNumDensityVector();
83 
84  // fill absorption length vector
85  G4int NumEl = aMaterial->GetNumberOfElements();
86  G4double absorptionLengthMax = 0.0;
87  for (size_t ibin=0; ibin<size_t(TotBin); ibin++) {
88  G4double SIGMA = 0. ;
89  for (size_t iel=0; iel<size_t(NumEl); iel++) {
90  G4int IndEl = (*elementVector)[iel]->GetIndex();
91  SIGMA += atomicNumDensityVector[iel]*
92  (*((*aCrossSectionTable)[IndEl]))[ibin];
93  }
94  // absorption length=5./SIGMA
95  absorptionLengthVector->PutValue(ibin, 5./SIGMA);
96  if (absorptionLengthMax < 5./SIGMA ) absorptionLengthMax = 5./SIGMA;
97  }
98 }
99 
100 
101 
102 // ***********************************************************************
103 // ********************** ComputeCrossSection ****************************
104 // ***********************************************************************
106  G4double KineticEnergy)
107 {
108  // Compute the "absorption" cross section of the photon "absorption"
109  // cross section means here the sum of the cross sections of the
110  // pair production, Compton scattering and photoelectric processes
111  const G4double t1keV = 1.*keV;
112  const G4double t200keV = 200.*keV;
113  const G4double t100MeV = 100.*MeV;
114 
115  // compute Z dependent quantities in the case of a new AtomicNumber
116  if(std::abs(AtomicNumber-Z)>0.1) {
117  Z = AtomicNumber;
118  G4double Zsquare = Z*Z;
119  G4double Zlog = std::log(Z);
120  G4double Zlogsquare = Zlog*Zlog;
121 
122  s200keV = (0.2651-0.1501*Zlog+0.02283*Zlogsquare)*Zsquare;
123  tmin = (0.552+218.5/Z+557.17/Zsquare)*MeV;
124  smin = (0.01239+0.005585*Zlog-0.000923*Zlogsquare)*std::exp(1.5*Zlog);
125  cmin=std::log(s200keV/smin)/(std::log(tmin/t200keV)*std::log(tmin/t200keV));
126  tlow = 0.2*std::exp(-7.355/std::sqrt(Z))*MeV;
127  slow = s200keV*std::exp(0.042*Z*std::log(t200keV/tlow)*std::log(t200keV/tlow));
128  s1keV = 300.*Zsquare;
129  clow =std::log(s1keV/slow)/std::log(tlow/t1keV);
130 
131  chigh=(7.55e-5-0.0542e-5*Z)*Zsquare*Z/std::log(t100MeV/tmin);
132  }
133 
134  // calculate the cross section (using an approximate empirical formula)
135  G4double xs;
136  if ( KineticEnergy<tlow ) {
137  if(KineticEnergy<t1keV) xs = slow*std::exp(clow*std::log(tlow/t1keV));
138  else xs = slow*std::exp(clow*std::log(tlow/KineticEnergy));
139 
140  } else if ( KineticEnergy<t200keV ) {
141  xs = s200keV
142  * std::exp(0.042*Z*std::log(t200keV/KineticEnergy)*std::log(t200keV/KineticEnergy));
143 
144  } else if( KineticEnergy<tmin ){
145  xs = smin
146  * std::exp(cmin*std::log(tmin/KineticEnergy)*std::log(tmin/KineticEnergy));
147 
148  } else {
149  xs = smin + chigh*std::log(KineticEnergy/tmin);
150 
151  }
152  return xs * barn;
153 }
154 
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static const G4double tlow
static G4ParticleTable * GetParticleTable()
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double keV
Definition: G4SIunits.hh:216
#define G4endl
Definition: G4ios.hh:61
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:191
Float_t Z
double G4double
Definition: G4Types.hh:76
G4double ComputeCrossSection(G4double AtomicNumber, G4double KineticEnergy)
std::vector< G4Element * > G4ElementVector
int G4int
Definition: G4Types.hh:78
static constexpr double barn
Definition: G4SIunits.hh:105
void BuildAbsorptionLengthVector(const G4Material *aMaterial, G4RangeVector *rangeVector)
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * theParticle
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:217
void PutValue(size_t index, G4double theValue)
size_t GetNumberOfElements() const
Definition: G4Material.hh:187