| 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 // $Id: G4EnergyLossForExtrapolator.cc 87062 2014-11-24 11:45:21Z gcosmo $
26 // 27 //
27 //-------------------------------------------- 28 //---------------------------------------------------------------------------
28 // 29 //
29 // ClassName: G4EnergyLossForExtrapolator 30 // ClassName: G4EnergyLossForExtrapolator
30 // 31 //
31 // Description: This class provide calculatio 32 // Description: This class provide calculation of energy loss, fluctuation,
32 // and msc angle 33 // and msc angle
33 // 34 //
34 // Author: 09.12.04 V.Ivanchenko 35 // Author: 09.12.04 V.Ivanchenko
35 // 36 //
36 // Modification: 37 // Modification:
37 // 08-04-05 Rename Propogator -> Extrapolator 38 // 08-04-05 Rename Propogator -> Extrapolator (V.Ivanchenko)
38 // 16-03-06 Add muon tables and fix bug in uni 39 // 16-03-06 Add muon tables and fix bug in units (V.Ivanchenko)
39 // 21-03-06 Add verbosity defined in the const 40 // 21-03-06 Add verbosity defined in the constructor and Initialisation
40 // start only when first public metho 41 // start only when first public method is called (V.Ivanchenko)
41 // 03-05-06 Remove unused pointer G4Material* 42 // 03-05-06 Remove unused pointer G4Material* from number of methods (VI)
42 // 12-05-06 SEt linLossLimit=0.001 (VI) 43 // 12-05-06 SEt linLossLimit=0.001 (VI)
43 // 44 //
44 //-------------------------------------------- 45 //----------------------------------------------------------------------------
45 // 46 //
46 47
47 //....oooOO0OOooo........oooOO0OOooo........oo 48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
48 49
49 #include "G4EnergyLossForExtrapolator.hh" 50 #include "G4EnergyLossForExtrapolator.hh"
50 #include "G4PhysicalConstants.hh" 51 #include "G4PhysicalConstants.hh"
51 #include "G4SystemOfUnits.hh" 52 #include "G4SystemOfUnits.hh"
52 #include "G4ParticleDefinition.hh" 53 #include "G4ParticleDefinition.hh"
53 #include "G4Material.hh" 54 #include "G4Material.hh"
54 #include "G4MaterialCutsCouple.hh" 55 #include "G4MaterialCutsCouple.hh"
55 #include "G4Electron.hh" 56 #include "G4Electron.hh"
56 #include "G4Positron.hh" 57 #include "G4Positron.hh"
57 #include "G4Proton.hh" 58 #include "G4Proton.hh"
58 #include "G4MuonPlus.hh" 59 #include "G4MuonPlus.hh"
59 #include "G4MuonMinus.hh" 60 #include "G4MuonMinus.hh"
60 #include "G4ParticleTable.hh" 61 #include "G4ParticleTable.hh"
61 62
62 //....oooOO0OOooo........oooOO0OOooo........oo 63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63 64
64 #ifdef G4MULTITHREADED << 65 G4TablesForExtrapolator* G4EnergyLossForExtrapolator::tables = 0;
65 G4Mutex G4EnergyLossForExtrapolator::extrMutex <<
66 #endif <<
67 <<
68 G4TablesForExtrapolator* G4EnergyLossForExtrap <<
69 66
70 G4EnergyLossForExtrapolator::G4EnergyLossForEx 67 G4EnergyLossForExtrapolator::G4EnergyLossForExtrapolator(G4int verb)
71 : maxEnergyTransfer(DBL_MAX), verbose(verb) 68 : maxEnergyTransfer(DBL_MAX), verbose(verb)
72 { 69 {
73 emin = 1.*CLHEP::MeV; << 70 currentParticle = 0;
74 emax = 100.*CLHEP::TeV; << 71 currentMaterial = 0;
75 } <<
76 <<
77 //....oooOO0OOooo........oooOO0OOooo........oo <<
78 72
79 G4EnergyLossForExtrapolator::~G4EnergyLossForE << 73 linLossLimit = 0.01;
80 { << 74 emin = 1.*MeV;
81 if(isMaster) { << 75 emax = 10.*TeV;
82 delete tables; << 76 nbins = 70;
83 tables = nullptr; << 77
84 } << 78 nmat = index = 0;
>> 79
>> 80 mass = charge2 = electronDensity = radLength = bg2 = beta2
>> 81 = kineticEnergy = tmax = 0.0;
>> 82 gam = 1.0;
>> 83
>> 84 idxDedxElectron = idxDedxPositron = idxDedxMuon = idxDedxProton
>> 85 = idxRangeElectron = idxRangePositron = idxRangeMuon = idxRangeProton
>> 86 = idxInvRangeElectron = idxInvRangePositron = idxInvRangeMuon
>> 87 = idxInvRangeProton = idxMscElectron = 0;
>> 88
>> 89 electron = positron = proton = muonPlus = muonMinus = 0;
85 } 90 }
86 91
87 //....oooOO0OOooo........oooOO0OOooo........oo 92 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88 93
89 G4double 94 G4double
90 G4EnergyLossForExtrapolator::EnergyAfterStep(G 95 G4EnergyLossForExtrapolator::EnergyAfterStep(G4double kinEnergy,
91 G4double stepLength, 96 G4double stepLength,
92 const G4Material* mat, 97 const G4Material* mat,
93 const G4ParticleDefinition* par 98 const G4ParticleDefinition* part)
94 { 99 {
>> 100 if(0 == nmat) { Initialisation(); }
95 G4double kinEnergyFinal = kinEnergy; 101 G4double kinEnergyFinal = kinEnergy;
96 if(SetupKinematics(part, mat, kinEnergy)) { 102 if(SetupKinematics(part, mat, kinEnergy)) {
97 G4double step = TrueStepLength(kinEnergy,s 103 G4double step = TrueStepLength(kinEnergy,stepLength,mat,part);
98 G4double r = ComputeRange(kinEnergy,part, << 104 G4double r = ComputeRange(kinEnergy,part);
99 if(r <= step) { 105 if(r <= step) {
100 kinEnergyFinal = 0.0; 106 kinEnergyFinal = 0.0;
101 } else if(step < linLossLimit*r) { 107 } else if(step < linLossLimit*r) {
102 kinEnergyFinal -= step*ComputeDEDX(kinEn << 108 kinEnergyFinal -= step*ComputeDEDX(kinEnergy,part);
103 } else { 109 } else {
104 G4double r1 = r - step; 110 G4double r1 = r - step;
105 kinEnergyFinal = ComputeEnergy(r1,part,m << 111 kinEnergyFinal = ComputeEnergy(r1,part);
106 } 112 }
107 } 113 }
108 return kinEnergyFinal; 114 return kinEnergyFinal;
109 } 115 }
110 116
111 //....oooOO0OOooo........oooOO0OOooo........oo 117 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
112 118
113 G4double 119 G4double
114 G4EnergyLossForExtrapolator::EnergyBeforeStep( 120 G4EnergyLossForExtrapolator::EnergyBeforeStep(G4double kinEnergy,
115 G4double stepLength, 121 G4double stepLength,
116 const G4Material* mat, 122 const G4Material* mat,
117 const G4ParticleDefinition* pa 123 const G4ParticleDefinition* part)
118 { 124 {
119 //G4cout << "G4EnergyLossForExtrapolator::En 125 //G4cout << "G4EnergyLossForExtrapolator::EnergyBeforeStep" << G4endl;
>> 126 if(0 == nmat) { Initialisation(); }
120 G4double kinEnergyFinal = kinEnergy; 127 G4double kinEnergyFinal = kinEnergy;
121 128
122 if(SetupKinematics(part, mat, kinEnergy)) { 129 if(SetupKinematics(part, mat, kinEnergy)) {
123 G4double step = TrueStepLength(kinEnergy,s 130 G4double step = TrueStepLength(kinEnergy,stepLength,mat,part);
124 G4double r = ComputeRange(kinEnergy,part,m << 131 G4double r = ComputeRange(kinEnergy,part);
125 132
126 if(step < linLossLimit*r) { 133 if(step < linLossLimit*r) {
127 kinEnergyFinal += step*ComputeDEDX(kinEn << 134 kinEnergyFinal += step*ComputeDEDX(kinEnergy,part);
128 } else { 135 } else {
129 G4double r1 = r + step; 136 G4double r1 = r + step;
130 kinEnergyFinal = ComputeEnergy(r1,part,m << 137 kinEnergyFinal = ComputeEnergy(r1,part);
131 } 138 }
132 } 139 }
133 return kinEnergyFinal; 140 return kinEnergyFinal;
134 } 141 }
135 142
136 //....oooOO0OOooo........oooOO0OOooo........oo 143 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
137 144
138 G4double 145 G4double
139 G4EnergyLossForExtrapolator::TrueStepLength(G4 146 G4EnergyLossForExtrapolator::TrueStepLength(G4double kinEnergy,
140 G4double stepLength, 147 G4double stepLength,
141 const G4Material* mat, 148 const G4Material* mat,
142 const G4ParticleDefinition* part 149 const G4ParticleDefinition* part)
143 { 150 {
144 G4double res = stepLength; 151 G4double res = stepLength;
>> 152 if(0 == nmat) { Initialisation(); }
145 //G4cout << "## G4EnergyLossForExtrapolator: 153 //G4cout << "## G4EnergyLossForExtrapolator::TrueStepLength L= " << res
146 // << " " << part->GetParticleName() << 154 // << " " << part->GetParticleName() << G4endl;
147 if(SetupKinematics(part, mat, kinEnergy)) { 155 if(SetupKinematics(part, mat, kinEnergy)) {
148 if(part == electron || part == positron) { 156 if(part == electron || part == positron) {
149 const G4double x = stepLength* << 157 G4double x = stepLength*
150 ComputeValue(kinEnergy, GetPhysicsTable(fMsc << 158 ComputeValue(kinEnergy, GetPhysicsTable(fMscElectron), idxMscElectron);
151 //G4cout << " x= " << x << G4endl; 159 //G4cout << " x= " << x << G4endl;
152 if(x < 0.2) { res *= (1.0 + 0.5* 160 if(x < 0.2) { res *= (1.0 + 0.5*x + x*x/3.0); }
153 else if(x < 0.9999) { res = -G4Log(1.0 - 161 else if(x < 0.9999) { res = -G4Log(1.0 - x)*stepLength/x; }
154 else { res = ComputeRange(kinEnergy, par << 162 else { res = ComputeRange(kinEnergy,part); }
155 } else { 163 } else {
156 res = ComputeTrueStep(mat,part,kinEnergy 164 res = ComputeTrueStep(mat,part,kinEnergy,stepLength);
157 } 165 }
158 } 166 }
159 return res; 167 return res;
160 } 168 }
161 169
162 //....oooOO0OOooo........oooOO0OOooo........oo 170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
163 171
164 G4bool 172 G4bool
165 G4EnergyLossForExtrapolator::SetupKinematics(c 173 G4EnergyLossForExtrapolator::SetupKinematics(const G4ParticleDefinition* part,
166 const G4Material* mat, 174 const G4Material* mat,
167 G4double kinEnergy) 175 G4double kinEnergy)
168 { 176 {
169 if(mat->GetNumberOfMaterials() != nmat) { In << 177 if(0 == nmat) { Initialisation(); }
170 if(nullptr == part || nullptr == mat || kinE << 178 if(!part || !mat || kinEnergy < keV) { return false; }
171 { return false; } << 179 G4bool flag = false;
172 if(part != currentParticle) { 180 if(part != currentParticle) {
>> 181 flag = true;
173 currentParticle = part; 182 currentParticle = part;
>> 183 mass = part->GetPDGMass();
174 G4double q = part->GetPDGCharge()/eplus; 184 G4double q = part->GetPDGCharge()/eplus;
175 charge2 = q*q; 185 charge2 = q*q;
176 } 186 }
177 if(mat != currentMaterial) { 187 if(mat != currentMaterial) {
178 size_t i = mat->GetIndex(); << 188 G4int i = mat->GetIndex();
179 if(i >= nmat) { 189 if(i >= nmat) {
180 G4cout << "### G4EnergyLossForExtrapolat << 190 G4cout << "### G4EnergyLossForExtrapolator WARNING:index i= "
181 << i << " above number of materials " < << 191 << i << " is out of table - NO extrapolation" << G4endl;
182 return false; <<
183 } else { 192 } else {
>> 193 flag = true;
184 currentMaterial = mat; 194 currentMaterial = mat;
185 electronDensity = mat->GetElectronDensit 195 electronDensity = mat->GetElectronDensity();
186 radLength = mat->GetRadlen(); 196 radLength = mat->GetRadlen();
>> 197 index = i;
187 } 198 }
188 } 199 }
189 if(kinEnergy != kineticEnergy) { << 200 if(flag || kinEnergy != kineticEnergy) {
190 kineticEnergy = kinEnergy; 201 kineticEnergy = kinEnergy;
191 G4double mass = part->GetPDGMass(); <<
192 G4double tau = kinEnergy/mass; 202 G4double tau = kinEnergy/mass;
193 203
194 gam = tau + 1.0; 204 gam = tau + 1.0;
195 bg2 = tau * (tau + 2.0); 205 bg2 = tau * (tau + 2.0);
196 beta2 = bg2/(gam*gam); 206 beta2 = bg2/(gam*gam);
197 tmax = kinEnergy; 207 tmax = kinEnergy;
198 if(part == electron) tmax *= 0.5; 208 if(part == electron) tmax *= 0.5;
199 else if(part != positron) { 209 else if(part != positron) {
200 G4double r = CLHEP::electron_mass_c2/mas << 210 G4double r = electron_mass_c2/mass;
201 tmax = 2.0*bg2*CLHEP::electron_mass_c2/( << 211 tmax = 2.0*bg2*electron_mass_c2/(1.0 + 2.0*gam*r + r*r);
202 } 212 }
203 tmax = std::min(tmax, maxEnergyTransfer); << 213 if(tmax > maxEnergyTransfer) { tmax = maxEnergyTransfer; }
204 } 214 }
205 return true; 215 return true;
206 } 216 }
207 217
208 //....oooOO0OOooo........oooOO0OOooo........oo 218 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
209 219
210 const G4ParticleDefinition* 220 const G4ParticleDefinition*
211 G4EnergyLossForExtrapolator::FindParticle(cons 221 G4EnergyLossForExtrapolator::FindParticle(const G4String& name)
212 { 222 {
213 currentParticle = G4ParticleTable::GetPartic << 223 const G4ParticleDefinition* p = 0;
214 if(nullptr == currentParticle) { << 224 if(name != currentParticleName) {
215 G4cout << "### G4EnergyLossForExtrapolator << 225 p = G4ParticleTable::GetParticleTable()->FindParticle(name);
216 << "FindParticle() fails to find " << nam << 226 if(!p) {
>> 227 G4cout << "### G4EnergyLossForExtrapolator WARNING: "
>> 228 << "FindParticle() fails to find "
>> 229 << name << G4endl;
>> 230 }
>> 231 } else {
>> 232 p = currentParticle;
217 } 233 }
218 return currentParticle; << 234 return p;
219 } 235 }
220 236
221 //....oooOO0OOooo........oooOO0OOooo........oo 237 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
222 238
223 G4double 239 G4double
224 G4EnergyLossForExtrapolator::ComputeDEDX(G4dou 240 G4EnergyLossForExtrapolator::ComputeDEDX(G4double ekin,
225 const G4ParticleDefinition* part, << 241 const G4ParticleDefinition* part)
226 const <<
227 { 242 {
228 if(mat->GetNumberOfMaterials() != nmat) { In <<
229 G4double x = 0.0; 243 G4double x = 0.0;
230 if(part == electron) { 244 if(part == electron) {
231 x = ComputeValue(ekin, GetPhysicsTable(fDe << 245 x = ComputeValue(ekin, GetPhysicsTable(fDedxElectron), idxDedxElectron);
232 } else if(part == positron) { 246 } else if(part == positron) {
233 x = ComputeValue(ekin, GetPhysicsTable(fDe << 247 x = ComputeValue(ekin, GetPhysicsTable(fDedxPositron), idxDedxPositron);
234 } else if(part == muonPlus || part == muonMi 248 } else if(part == muonPlus || part == muonMinus) {
235 x = ComputeValue(ekin, GetPhysicsTable(fDe << 249 x = ComputeValue(ekin, GetPhysicsTable(fDedxMuon), idxDedxMuon);
236 } else { 250 } else {
237 G4double e = ekin*CLHEP::proton_mass_c2/pa << 251 G4double e = ekin*proton_mass_c2/mass;
238 G4double q = part->GetPDGCharge()/CLHEP::e << 252 x = ComputeValue(e, GetPhysicsTable(fDedxProton), idxDedxProton)*charge2;
239 x = ComputeValue(e, GetPhysicsTable(fDedxP <<
240 } 253 }
241 return x; 254 return x;
242 } 255 }
243 256
244 //....oooOO0OOooo........oooOO0OOooo........oo 257 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
245 258
246 G4double 259 G4double
247 G4EnergyLossForExtrapolator::ComputeRange(G4do 260 G4EnergyLossForExtrapolator::ComputeRange(G4double ekin,
248 const G4ParticleDefinition* part, << 261 const G4ParticleDefinition* part)
249 const G4Material* mat) <<
250 { 262 {
251 if(mat->GetNumberOfMaterials() != nmat) { In <<
252 G4double x = 0.0; 263 G4double x = 0.0;
253 if(part == electron) { 264 if(part == electron) {
254 x = ComputeValue(ekin, GetPhysicsTable(fRa << 265 x = ComputeValue(ekin, GetPhysicsTable(fRangeElectron), idxRangeElectron);
255 } else if(part == positron) { 266 } else if(part == positron) {
256 x = ComputeValue(ekin, GetPhysicsTable(fRa << 267 x = ComputeValue(ekin, GetPhysicsTable(fRangePositron), idxRangePositron);
257 } else if(part == muonPlus || part == muonMi 268 } else if(part == muonPlus || part == muonMinus) {
258 x = ComputeValue(ekin, GetPhysicsTable(fRa << 269 x = ComputeValue(ekin, GetPhysicsTable(fRangeMuon), idxRangeMuon);
259 } else { 270 } else {
260 G4double massratio = CLHEP::proton_mass_c2 << 271 G4double massratio = proton_mass_c2/mass;
261 G4double e = ekin*massratio; 272 G4double e = ekin*massratio;
262 G4double q = part->GetPDGCharge()/CLHEP::e << 273 x = ComputeValue(e, GetPhysicsTable(fRangeProton), idxRangeProton)
263 x = ComputeValue(e, GetPhysicsTable(fRange << 274 /(charge2*massratio);
264 /(q*q*massratio); <<
265 } 275 }
266 return x; 276 return x;
267 } 277 }
268 278
269 //....oooOO0OOooo........oooOO0OOooo........oo 279 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
270 280
271 G4double 281 G4double
272 G4EnergyLossForExtrapolator::ComputeEnergy(G4d 282 G4EnergyLossForExtrapolator::ComputeEnergy(G4double range,
273 const G4ParticleDefinition* part, << 283 const G4ParticleDefinition* part)
274 const G4Material* mat) <<
275 { 284 {
276 if(mat->GetNumberOfMaterials() != nmat) { In <<
277 G4double x = 0.0; 285 G4double x = 0.0;
278 if(part == electron) { 286 if(part == electron) {
279 x = ComputeValue(range,GetPhysicsTable(fIn << 287 x = ComputeValue(range, GetPhysicsTable(fInvRangeElectron),
>> 288 idxInvRangeElectron);
280 } else if(part == positron) { 289 } else if(part == positron) {
281 x = ComputeValue(range,GetPhysicsTable(fIn << 290 x = ComputeValue(range, GetPhysicsTable(fInvRangePositron),
>> 291 idxInvRangePositron);
282 } else if(part == muonPlus || part == muonMi 292 } else if(part == muonPlus || part == muonMinus) {
283 x = ComputeValue(range, GetPhysicsTable(fI << 293 x = ComputeValue(range, GetPhysicsTable(fInvRangeMuon), idxInvRangeMuon);
284 } else { 294 } else {
285 G4double massratio = CLHEP::proton_mass_c2 << 295 G4double massratio = proton_mass_c2/mass;
286 G4double q = part->GetPDGCharge()/CLHEP::e << 296 G4double r = range*massratio*charge2;
287 G4double r = range*massratio*q*q; << 297 x = ComputeValue(r, GetPhysicsTable(fInvRangeProton),
288 x = ComputeValue(r, GetPhysicsTable(fInvRa << 298 idxInvRangeProton)/massratio;
289 } 299 }
290 return x; 300 return x;
291 } 301 }
292 302
293 //....oooOO0OOooo........oooOO0OOooo........oo 303 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
294 << 304
295 G4double << 305 void G4EnergyLossForExtrapolator::Initialisation()
296 G4EnergyLossForExtrapolator::EnergyDispersion( <<
297 G4double stepLength, <<
298 const G4Material* mat, <<
299 const G4ParticleDefinition* pa <<
300 { 306 {
301 G4double sig2 = 0.0; << 307 if(verbose>1) {
302 if(SetupKinematics(part, mat, kinEnergy)) { << 308 G4cout << "### G4EnergyLossForExtrapolator::Initialisation" << G4endl;
303 G4double step = ComputeTrueStep(mat,part,k <<
304 sig2 = (1.0/beta2 - 0.5) <<
305 *CLHEP::twopi_mc2_rcl2*tmax*step*electro <<
306 } 309 }
307 return sig2; << 310 currentParticle = 0;
>> 311 currentMaterial = 0;
>> 312 kineticEnergy = 0.0;
>> 313
>> 314 electron = G4Electron::Electron();
>> 315 positron = G4Positron::Positron();
>> 316 proton = G4Proton::Proton();
>> 317 muonPlus = G4MuonPlus::MuonPlus();
>> 318 muonMinus= G4MuonMinus::MuonMinus();
>> 319
>> 320 nmat = G4Material::GetNumberOfMaterials();
>> 321 currentParticleName = "";
>> 322 BuildTables();
308 } 323 }
309 324
310 //....oooOO0OOooo........oooOO0OOooo........oo 325 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
311 326
312 G4double G4EnergyLossForExtrapolator::AverageS << 327 #include "G4AutoLock.hh"
313 G4double kinEnergy, << 328 namespace { G4Mutex G4EnergyLossForExtrapolatorMutex = G4MUTEX_INITIALIZER; }
314 G4double stepLength, << 329
315 const G4Material* mat, << 330 void G4EnergyLossForExtrapolator::BuildTables()
316 const G4ParticleDefinition* part) <<
317 { 331 {
318 G4double theta = 0.0; << 332 G4AutoLock l(&G4EnergyLossForExtrapolatorMutex);
319 if(SetupKinematics(part, mat, kinEnergy)) { << 333
320 G4double t = stepLength/radLength; << 334 if(!tables) {
321 G4double y = std::max(0.001, t); << 335 if(verbose > 0) {
322 theta = 19.23*CLHEP::MeV*std::sqrt(charge2 << 336 G4cout << "### G4EnergyLossForExtrapolator::BuildTables for "
323 /(beta2*gam*part->GetPDGMass()); << 337 << nmat << " materials Nbins= " << nbins
>> 338 << " Emin(MeV)= " << emin << " Emax(MeV)= " << emax << G4endl;
>> 339 }
>> 340 tables = new G4TablesForExtrapolator(verbose, nbins, emin, emax);
324 } 341 }
325 return theta; <<
326 } 342 }
327 343
328 //....oooOO0OOooo........oooOO0OOooo........oo 344 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
329 345
330 void G4EnergyLossForExtrapolator::Initialisati << 346 G4EnergyLossForExtrapolator::~G4EnergyLossForExtrapolator()
331 { 347 {
332 if(verbose>0) { << 348 G4AutoLock l(&G4EnergyLossForExtrapolatorMutex);
333 G4cout << "### G4EnergyLossForExtrapolator << 349 delete tables;
334 << tables << G4endl; << 350 tables = 0;
335 } <<
336 electron = G4Electron::Electron(); <<
337 positron = G4Positron::Positron(); <<
338 proton = G4Proton::Proton(); <<
339 muonPlus = G4MuonPlus::MuonPlus(); <<
340 muonMinus= G4MuonMinus::MuonMinus(); <<
341 <<
342 // initialisation for the 1st run <<
343 if(nullptr == tables) { <<
344 #ifdef G4MULTITHREADED <<
345 G4MUTEXLOCK(&extrMutex); <<
346 if(nullptr == tables) { <<
347 #endif <<
348 isMaster = true; <<
349 tables = new G4TablesForExtrapolator(ver <<
350 tables->Initialisation(); <<
351 nmat = G4Material::GetNumberOfMaterials( <<
352 if(verbose > 0) { <<
353 G4cout << "### G4EnergyLossForExtrapol <<
354 << nmat << " materials Nbins= " <<
355 << nbins << " Emin(MeV)= " << e <<
356 << G4endl; <<
357 } <<
358 #ifdef G4MULTITHREADED <<
359 } <<
360 G4MUTEXUNLOCK(&extrMutex); <<
361 #endif <<
362 } <<
363 <<
364 // initialisation for the next run <<
365 if(isMaster && G4Material::GetNumberOfMateri <<
366 #ifdef G4MULTITHREADED <<
367 G4MUTEXLOCK(&extrMutex); <<
368 #endif <<
369 tables->Initialisation(); <<
370 #ifdef G4MULTITHREADED <<
371 G4MUTEXUNLOCK(&extrMutex); <<
372 #endif <<
373 } <<
374 nmat = G4Material::GetNumberOfMaterials(); <<
375 } 351 }
376 352
377 //....oooOO0OOooo........oooOO0OOooo........oo 353 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 354
378 355