| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
>> 23 // $Id: G4EmCalculator.cc,v 1.19 2005/05/30 08:55:51 vnivanch Exp $
>> 24 // GEANT4 tag $Name: geant4-07-01 $
26 // 25 //
27 // ------------------------------------------- 26 // -------------------------------------------------------------------
28 // 27 //
29 // GEANT4 Class file 28 // GEANT4 Class file
30 // 29 //
31 // 30 //
32 // File name: G4EmCalculator 31 // File name: G4EmCalculator
33 // 32 //
34 // Author: Vladimir Ivanchenko 33 // Author: Vladimir Ivanchenko
35 // 34 //
36 // Creation date: 28.06.2004 35 // Creation date: 28.06.2004
37 // 36 //
>> 37 // Modifications:
>> 38 // 12.09.2004 Add verbosity (V.Ivanchenko)
>> 39 // 17.11.2004 Change signature of methods, add new methods (V.Ivanchenko)
>> 40 // 08.04.2005 Major optimisation of internal interfaces (V.Ivantchenko)
>> 41 // 08.05.2005 Use updated interfaces (V.Ivantchenko)
38 // 42 //
39 // Class Description: V.Ivanchenko & M.Novak << 43 // Class Description:
40 // 44 //
41 // ------------------------------------------- 45 // -------------------------------------------------------------------
42 // 46 //
43 //....oooOO0OOooo........oooOO0OOooo........oo 47 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
44 //....oooOO0OOooo........oooOO0OOooo........oo 48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
45 49
46 #include "G4EmCalculator.hh" 50 #include "G4EmCalculator.hh"
47 #include "G4SystemOfUnits.hh" <<
48 #include "G4LossTableManager.hh" 51 #include "G4LossTableManager.hh"
49 #include "G4EmParameters.hh" <<
50 #include "G4NistManager.hh" <<
51 #include "G4DynamicParticle.hh" <<
52 #include "G4VEmProcess.hh" 52 #include "G4VEmProcess.hh"
53 #include "G4VEnergyLossProcess.hh" 53 #include "G4VEnergyLossProcess.hh"
54 #include "G4VMultipleScattering.hh" 54 #include "G4VMultipleScattering.hh"
55 #include "G4Material.hh" 55 #include "G4Material.hh"
56 #include "G4MaterialCutsCouple.hh" 56 #include "G4MaterialCutsCouple.hh"
57 #include "G4ParticleDefinition.hh" 57 #include "G4ParticleDefinition.hh"
58 #include "G4ParticleTable.hh" 58 #include "G4ParticleTable.hh"
59 #include "G4IonTable.hh" <<
60 #include "G4PhysicsTable.hh" 59 #include "G4PhysicsTable.hh"
61 #include "G4ProductionCutsTable.hh" 60 #include "G4ProductionCutsTable.hh"
62 #include "G4ProcessManager.hh" 61 #include "G4ProcessManager.hh"
63 #include "G4ionEffectiveCharge.hh" 62 #include "G4ionEffectiveCharge.hh"
64 #include "G4RegionStore.hh" 63 #include "G4RegionStore.hh"
65 #include "G4Element.hh" 64 #include "G4Element.hh"
66 #include "G4EmCorrections.hh" 65 #include "G4EmCorrections.hh"
67 #include "G4GenericIon.hh" <<
68 #include "G4ProcessVector.hh" <<
69 #include "G4Gamma.hh" <<
70 #include "G4Electron.hh" <<
71 #include "G4Positron.hh" <<
72 #include "G4EmUtility.hh" <<
73 66
74 //....oooOO0OOooo........oooOO0OOooo........oo 67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
75 68
76 G4EmCalculator::G4EmCalculator() 69 G4EmCalculator::G4EmCalculator()
77 { 70 {
78 manager = G4LossTableManager::Instance(); 71 manager = G4LossTableManager::Instance();
79 nist = G4NistManager::Instance(); <<
80 theParameters = G4EmParameters::Instance(); <<
81 corr = manager->EmCorrections(); 72 corr = manager->EmCorrections();
82 cutenergy[0] = cutenergy[1] = cutenergy[2] = << 73 nLocalMaterials = 0;
83 theGenericIon = G4GenericIon::GenericIon(); << 74 verbose = 0;
84 ionEffCharge = new G4ionEffectiveCharge(); << 75 currentCoupleIndex = 0;
85 dynParticle = new G4DynamicParticle(); << 76 currentCouple = 0;
86 ionTable = G4ParticleTable::GetParticle << 77 currentMaterial = 0;
>> 78 currentParticle = 0;
>> 79 currentLambda = 0;
>> 80 chargeSquare = 1.0;
>> 81 massRatio = 1.0;
>> 82 currentParticleName= "";
>> 83 currentMaterialName= "";
>> 84 ionEffCharge = new G4ionEffectiveCharge();
>> 85 isIon = false;
>> 86 isApplicable = false;
87 } 87 }
88 88
89 //....oooOO0OOooo........oooOO0OOooo........oo 89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
90 90
91 G4EmCalculator::~G4EmCalculator() 91 G4EmCalculator::~G4EmCalculator()
92 { 92 {
93 delete ionEffCharge; 93 delete ionEffCharge;
94 delete dynParticle; << 94 for (G4int i=0; i<nLocalMaterials; i++) {
95 for (G4int i=0; i<nLocalMaterials; ++i) { <<
96 delete localCouples[i]; 95 delete localCouples[i];
97 } 96 }
98 } 97 }
99 98
100 //....oooOO0OOooo........oooOO0OOooo........oo 99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
101 100
102 G4double G4EmCalculator::GetDEDX(G4double kinE << 101 G4double G4EmCalculator::GetDEDX(G4double kinEnergy, const G4ParticleDefinition* p,
103 const G4Parti << 102 const G4Material* mat, const G4Region* region)
104 const G4Mater <<
105 const G4Regio <<
106 { 103 {
107 G4double res = 0.0; 104 G4double res = 0.0;
108 const G4MaterialCutsCouple* couple = FindCou 105 const G4MaterialCutsCouple* couple = FindCouple(mat, region);
109 if(nullptr != couple && UpdateParticle(p, ki << 106 if(couple && UpdateParticle(p, kinEnergy) ) {
110 res = manager->GetDEDX(p, kinEnergy, coupl 107 res = manager->GetDEDX(p, kinEnergy, couple);
111 <<
112 if(isIon) { <<
113 if(FindEmModel(p, currentProcessName, ki <<
114 G4double length = CLHEP::nm; <<
115 G4double eloss = res*length; <<
116 //G4cout << "### GetDEDX: E= " << kinE <<
117 // << " de= " << eloss << G4endl <<
118 dynParticle->SetKineticEnergy(kinEnerg <<
119 currentModel->GetChargeSquareRatio(p, <<
120 currentModel->CorrectionsAlongStep(cou <<
121 res = eloss/length; <<
122 //G4cout << " de1= " << eloss << <<
123 // << " " << p->GetParticleName( <<
124 } <<
125 } <<
126 <<
127 if(verbose>0) { 108 if(verbose>0) {
128 G4cout << "G4EmCalculator::GetDEDX: E(Me << 109 G4cout << "E(MeV)= " << kinEnergy/MeV
129 << " DEDX(MeV/mm)= " << res*mm/Me << 110 << " DEDX(MeV/mm)= " << res*mm/MeV
130 << " DEDX(MeV*cm^2/g)= " << res*g << 111 << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
131 << " " << p->GetParticleName() << 112 << " " << p->GetParticleName()
132 << " in " << mat->GetName() << 113 << " in " << mat->GetName()
133 << " isIon= " << isIon << 114 << G4endl;
134 << G4endl; <<
135 } 115 }
136 } 116 }
137 return res; 117 return res;
138 } 118 }
139 119
140 //....oooOO0OOooo........oooOO0OOooo........oo 120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
141 121
142 G4double G4EmCalculator::GetRangeFromRestricte << 122 G4double G4EmCalculator::GetDEDX(G4double kinEnergy, const G4String& particle,
143 << 123 const G4String& material, const G4String& reg)
144 <<
145 <<
146 { 124 {
147 G4double res = 0.0; << 125 return GetDEDX(kinEnergy,FindParticle(particle),FindMaterial(material),FindRegion(reg));
148 const G4MaterialCutsCouple* couple = FindCou <<
149 if(couple && UpdateParticle(p, kinEnergy)) { <<
150 res = manager->GetRangeFromRestricteDEDX(p <<
151 if(verbose>1) { <<
152 G4cout << " G4EmCalculator::GetRangeFrom <<
153 << kinEnergy/MeV <<
154 << " range(mm)= " << res/mm <<
155 << " " << p->GetParticleName() <<
156 << " in " << mat->GetName() <<
157 << G4endl; <<
158 } <<
159 } <<
160 return res; <<
161 } 126 }
162 127
163 //....oooOO0OOooo........oooOO0OOooo........oo 128 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
164 129
165 G4double G4EmCalculator::GetCSDARange(G4double << 130 G4double G4EmCalculator::GetRange(G4double kinEnergy, const G4ParticleDefinition* p,
166 const G4 << 131 const G4Material* mat, const G4Region* region)
167 const G4 <<
168 const G4 <<
169 { 132 {
170 G4double res = 0.0; 133 G4double res = 0.0;
171 if(!theParameters->BuildCSDARange()) { <<
172 G4ExceptionDescription ed; <<
173 ed << "G4EmCalculator::GetCSDARange: CSDA <<
174 << " use UI command: /process/eLoss/CSD <<
175 G4Exception("G4EmCalculator::GetCSDARange" <<
176 JustWarning, ed); <<
177 return res; <<
178 } <<
179 <<
180 const G4MaterialCutsCouple* couple = FindCou 134 const G4MaterialCutsCouple* couple = FindCouple(mat,region);
181 if(nullptr != couple && UpdateParticle(p, ki << 135 if(couple && UpdateParticle(p, kinEnergy)) {
182 res = manager->GetCSDARange(p, kinEnergy, << 136 res = manager->GetRange(p, kinEnergy, couple);
183 if(verbose>1) { << 137 if(verbose>0) {
184 G4cout << " G4EmCalculator::GetCSDARange << 138 G4cout << "E(MeV)= " << kinEnergy/MeV
185 << " range(mm)= " << res/mm << 139 << " range(mm)= " << res/mm
186 << " " << p->GetParticleName() << 140 << " " << p->GetParticleName()
187 << " in " << mat->GetName() << 141 << " in " << mat->GetName()
188 << G4endl; << 142 << G4endl;
189 } 143 }
190 } 144 }
191 return res; 145 return res;
192 } 146 }
193 147
194 //....oooOO0OOooo........oooOO0OOooo........oo 148 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
195 149
196 G4double G4EmCalculator::GetRange(G4double kin << 150 G4double G4EmCalculator::GetRange(G4double kinEnergy, const G4String& particle,
197 const G4Part << 151 const G4String& material, const G4String& reg)
198 const G4Mate <<
199 const G4Regi <<
200 { 152 {
201 G4double res = 0.0; << 153 return GetRange(kinEnergy,FindParticle(particle),FindMaterial(material),FindRegion(reg));
202 if(theParameters->BuildCSDARange()) { <<
203 res = GetCSDARange(kinEnergy, p, mat, regi <<
204 } else { <<
205 res = GetRangeFromRestricteDEDX(kinEnergy, <<
206 } <<
207 return res; <<
208 } 154 }
209 155
210 //....oooOO0OOooo........oooOO0OOooo........oo 156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
211 157
212 G4double G4EmCalculator::GetKinEnergy(G4double << 158 G4double G4EmCalculator::GetKinEnergy(G4double range, const G4ParticleDefinition* p,
213 const G4 << 159 const G4Material* mat, const G4Region* region)
214 const G4 <<
215 const G4 <<
216 { 160 {
217 G4double res = 0.0; 161 G4double res = 0.0;
218 const G4MaterialCutsCouple* couple = FindCou 162 const G4MaterialCutsCouple* couple = FindCouple(mat,region);
219 if(nullptr != couple && UpdateParticle(p, 1. << 163 if(couple && UpdateParticle(p, 1.0*GeV)) {
220 res = manager->GetEnergy(p, range, couple) 164 res = manager->GetEnergy(p, range, couple);
221 if(verbose>0) { 165 if(verbose>0) {
222 G4cout << "G4EmCalculator::GetKinEnergy: << 166 G4cout << "Range(mm)= " << range/mm
223 << " KinE(MeV)= " << res/MeV << 167 << " KinE(MeV)= " << res/MeV
224 << " " << p->GetParticleName() << 168 << " " << p->GetParticleName()
225 << " in " << mat->GetName() << 169 << " in " << mat->GetName()
226 << G4endl; << 170 << G4endl;
227 } 171 }
228 } 172 }
229 return res; 173 return res;
230 } 174 }
231 175
232 //....oooOO0OOooo........oooOO0OOooo........oo 176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
233 177
>> 178 G4double G4EmCalculator::GetKinEnergy(G4double range, const G4String& particle,
>> 179 const G4String& material, const G4String& reg)
>> 180 {
>> 181 return GetKinEnergy(range,FindParticle(particle),FindMaterial(material),FindRegion(reg));
>> 182 }
>> 183
>> 184 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 185
234 G4double G4EmCalculator::GetCrossSectionPerVol 186 G4double G4EmCalculator::GetCrossSectionPerVolume(G4double kinEnergy,
235 co 187 const G4ParticleDefinition* p,
236 co 188 const G4String& processName,
237 co << 189 const G4Material* mat,
238 co << 190 const G4Region* region)
239 { 191 {
240 G4double res = 0.0; 192 G4double res = 0.0;
241 const G4MaterialCutsCouple* couple = FindCou 193 const G4MaterialCutsCouple* couple = FindCouple(mat,region);
242 194
243 if(nullptr != couple && UpdateParticle(p, ki << 195 if(couple) {
244 if(FindEmModel(p, processName, kinEnergy)) << 196 G4int idx = couple->GetIndex();
245 G4int idx = couple->GetIndex(); << 197 FindLambdaTable(p, processName);
246 G4int procType = -1; << 198 if(currentLambda && UpdateParticle(p, kinEnergy)) {
247 FindLambdaTable(p, processName, kinEnerg << 199 G4bool b;
248 << 200 G4double e = kinEnergy*massRatio;
249 G4VEmProcess* emproc = FindDiscreteProce << 201 res = (((*currentLambda)[idx])->GetValue(e,b))*chargeSquare;
250 if(nullptr != emproc) { << 202 res /= mat->GetTotNbOfAtomsPerVolume();
251 res = emproc->GetCrossSection(kinEnergy, cou <<
252 } else if(currentLambda) { <<
253 // special tables are built for Msc mo <<
254 // procType is set in FindLambdaTable <<
255 if(procType==2) { <<
256 auto mscM = static_cast<G4VMscModel* <<
257 mscM->SetCurrentCouple(couple); <<
258 G4double tr1Mfp = mscM->GetTransport <<
259 if (tr1Mfp<DBL_MAX) { <<
260 res = 1./tr1Mfp; <<
261 } <<
262 } else { <<
263 G4double e = kinEnergy*massRatio; <<
264 res = (((*currentLambda)[idx])->Valu <<
265 } <<
266 } else { <<
267 res = ComputeCrossSectionPerVolume(kin <<
268 } <<
269 if(verbose>0) { 203 if(verbose>0) {
270 G4cout << "G4EmCalculator::GetXSPerVol << 204 G4cout << "E(MeV)= " << kinEnergy/MeV
271 << " cross(cm-1)= " << res*cm << 205 << " cross(cm^2/g)= " << res*gram/cm2
272 << " " << p->GetParticleName( << 206 << " " << p->GetParticleName()
273 << " in " << mat->GetName(); << 207 << " in " << mat->GetName()
274 if(verbose>1) << 208 << G4endl;
275 G4cout << " idx= " << idx << " Esc << 209 }
276 << kinEnergy*massRatio <<
277 << " q2= " << chargeSquare; <<
278 G4cout << G4endl; <<
279 } <<
280 } 210 }
281 } 211 }
282 return res; 212 return res;
283 } 213 }
284 214
285 //....oooOO0OOooo........oooOO0OOooo........oo 215 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
286 216
287 G4double G4EmCalculator::GetShellIonisationCro << 217 G4double G4EmCalculator::GetCrossSectionPerVolume(G4double kinEnergy,
288 const << 218 const G4String& particle,
289 G4int << 219 const G4String& processName,
290 G4Ato << 220 const G4String& material,
291 G4dou << 221 const G4String& reg)
292 { 222 {
293 G4double res = 0.0; << 223 return GetCrossSectionPerVolume(kinEnergy,FindParticle(particle),processName,
294 const G4ParticleDefinition* p = FindParticle << 224 FindMaterial(material),FindRegion(reg));
295 G4VAtomDeexcitation* ad = manager->AtomDeexc << 225 }
296 if(nullptr != p && nullptr != ad) { << 226
297 res = ad->GetShellIonisationCrossSectionPe << 227 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
298 } << 228
>> 229 G4double G4EmCalculator::GetCrossSectionPerAtom(G4double kinEnergy,
>> 230 const G4ParticleDefinition* p,
>> 231 const G4String& processName,
>> 232 const G4Material* mat,
>> 233 const G4Region* region)
>> 234 {
>> 235 G4double res = GetCrossSectionPerVolume(kinEnergy,p,processName,mat,region);
>> 236 if(mat) res /= mat->GetTotNbOfAtomsPerVolume();
299 return res; 237 return res;
300 } 238 }
301 239
302 //....oooOO0OOooo........oooOO0OOooo........oo 240 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
303 241
>> 242 G4double G4EmCalculator::GetCrossSectionPerAtom( G4double kinEnergy,
>> 243 const G4String& particle,
>> 244 const G4String& processName,
>> 245 const G4String& material,
>> 246 const G4String& reg)
>> 247 {
>> 248 return GetCrossSectionPerAtom(kinEnergy,FindParticle(particle),processName,
>> 249 FindMaterial(material),FindRegion(reg));
>> 250 }
>> 251
>> 252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 253
304 G4double G4EmCalculator::GetMeanFreePath(G4dou 254 G4double G4EmCalculator::GetMeanFreePath(G4double kinEnergy,
305 const 255 const G4ParticleDefinition* p,
306 const 256 const G4String& processName,
307 const << 257 const G4Material* mat,
308 const 258 const G4Region* region)
309 { 259 {
310 G4double res = DBL_MAX; 260 G4double res = DBL_MAX;
311 G4double x = GetCrossSectionPerVolume(kinEne 261 G4double x = GetCrossSectionPerVolume(kinEnergy,p, processName, mat,region);
312 if(x > 0.0) { res = 1.0/x; } << 262 if(x > 0.0) res = 1.0/x;
313 if(verbose>1) { 263 if(verbose>1) {
314 G4cout << "G4EmCalculator::GetMeanFreePath << 264 G4cout << "E(MeV)= " << kinEnergy/MeV
315 << " MFP(mm)= " << res/mm << 265 << " MFP(mm)= " << res/mm
316 << " " << p->GetParticleName() << 266 << " " << p->GetParticleName()
317 << " in " << mat->GetName() << 267 << " in " << mat->GetName()
318 << G4endl; << 268 << G4endl;
319 } 269 }
320 return res; 270 return res;
321 } 271 }
322 272
323 //....oooOO0OOooo........oooOO0OOooo........oo 273 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
324 274
>> 275 G4double G4EmCalculator::GetMeanFreePath(G4double kinEnergy,
>> 276 const G4String& particle,
>> 277 const G4String& processName,
>> 278 const G4String& material,
>> 279 const G4String& reg)
>> 280 {
>> 281 return GetMeanFreePath(kinEnergy,FindParticle(particle),processName,
>> 282 FindMaterial(material),FindRegion(reg));
>> 283 }
>> 284
>> 285 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 286
325 void G4EmCalculator::PrintDEDXTable(const G4Pa 287 void G4EmCalculator::PrintDEDXTable(const G4ParticleDefinition* p)
326 { 288 {
327 const G4VEnergyLossProcess* elp = manager->G << 289 const G4VEnergyLossProcess* elp = FindEnergyLossProcess(p);
328 G4cout << "##### DEDX Table for " << p->GetP 290 G4cout << "##### DEDX Table for " << p->GetParticleName() << G4endl;
329 if(nullptr != elp) G4cout << *(elp->DEDXTabl << 291 if(elp) G4cout << *(elp->DEDXTable()) << G4endl;
330 } 292 }
331 293
332 //....oooOO0OOooo........oooOO0OOooo........oo 294 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
333 295
334 void G4EmCalculator::PrintRangeTable(const G4P 296 void G4EmCalculator::PrintRangeTable(const G4ParticleDefinition* p)
335 { 297 {
336 const G4VEnergyLossProcess* elp = manager->G << 298 const G4VEnergyLossProcess* elp = FindEnergyLossProcess(p);
337 G4cout << "##### Range Table for " << p->Get 299 G4cout << "##### Range Table for " << p->GetParticleName() << G4endl;
338 if(nullptr != elp) G4cout << *(elp->RangeTab << 300 if(elp) G4cout << *(elp->RangeTableForLoss()) << G4endl;
339 } 301 }
340 302
341 //....oooOO0OOooo........oooOO0OOooo........oo 303 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
342 304
343 void G4EmCalculator::PrintInverseRangeTable(co 305 void G4EmCalculator::PrintInverseRangeTable(const G4ParticleDefinition* p)
344 { 306 {
345 const G4VEnergyLossProcess* elp = manager->G << 307 const G4VEnergyLossProcess* elp = FindEnergyLossProcess(p);
346 G4cout << "### G4EmCalculator: Inverse Range << 308 G4cout << "### G4EmCalculator: Inverse Range Table for " << p->GetParticleName() << G4endl;
347 << p->GetParticleName() << G4endl; << 309 if(elp) G4cout << *(elp->InverseRangeTable()) << G4endl;
348 if(nullptr != elp) G4cout << *(elp->InverseR <<
349 } 310 }
350 311
351 //....oooOO0OOooo........oooOO0OOooo........oo 312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
352 313
353 G4double G4EmCalculator::ComputeDEDX(G4double 314 G4double G4EmCalculator::ComputeDEDX(G4double kinEnergy,
354 const G4P 315 const G4ParticleDefinition* p,
355 const G4S 316 const G4String& processName,
356 const G4M << 317 const G4Material* mat,
357 G4d 318 G4double cut)
358 { 319 {
359 SetupMaterial(mat); <<
360 G4double res = 0.0; 320 G4double res = 0.0;
361 if(verbose > 1) { 321 if(verbose > 1) {
362 G4cout << "### G4EmCalculator::ComputeDEDX << 322 G4cout << "ComputeDEDX: " << p->GetParticleName()
363 << " in " << currentMaterialName << 323 << " in " << mat->GetName()
364 << " e(MeV)= " << kinEnergy/MeV << << 324 << " e(MeV)= " << kinEnergy/MeV << " cut(MeV)= " << cut/MeV << G4endl;
365 << G4endl; << 325 }
366 } << 326 if(FindEmModel(p, processName, kinEnergy)) {
367 if(UpdateParticle(p, kinEnergy)) { << 327 // G4cout << "currentModel= " << currentModel << G4endl;
368 if(FindEmModel(p, processName, kinEnergy)) << 328 if(UpdateParticle(p, kinEnergy)) {
369 G4double escaled = kinEnergy*massRatio; << 329 G4double escaled = kinEnergy*massRatio;
370 if(nullptr != baseParticle) { << 330 if(baseParticle) {
371 res = currentModel->ComputeDEDXPerVolume(mat << 331 res = currentModel->ComputeDEDXPerVolume(mat, baseParticle, escaled, cut)
372 << 332 * chargeSquare;
373 if(verbose > 1) { << 333 if(verbose > 1) {
374 G4cout << "Particle: " << p->GetParticleNa << 334 G4cout << baseParticle->GetParticleName() << " E(MeV)= " << escaled
375 << " E(MeV)=" << kinEnergy << 335 << " DEDX(MeV/mm)= " << res*mm/MeV
376 << " Base particle: " << baseParticle->Ge << 336 << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
377 << " Escaled(MeV)= " << escaled << 337 << G4endl;
378 << " q2=" << chargeSquare << G4endl; << 338 }
379 } <<
380 } else { 339 } else {
381 res = currentModel->ComputeDEDXPerVolume(mat << 340 res = currentModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cut);
382 if(verbose > 1) { <<
383 G4cout << "Particle: " << p->GetParticleNa <<
384 << " E(MeV)=" << kinEnergy << G4endl; <<
385 } <<
386 } <<
387 if(verbose > 1) { <<
388 G4cout << currentModel->GetName() << ": DEDX <<
389 << " DEDX(MeV*cm^2/g)= " <<
390 << res*gram/(MeV*cm2*mat->GetDensity( <<
391 << G4endl; <<
392 } 341 }
393 // emulate smoothing procedure << 342 if(isIon && currentModel->HighEnergyLimit() > 100.*MeV)
394 if(applySmoothing && nullptr != loweMode << 343 res += corr->HighOrderCorrections(p,mat,kinEnergy);
395 G4double eth = currentModel->LowEnergyLimit( << 344
396 G4double res0 = 0.0; << 345 // emulate boundary region for different parameterisations
397 G4double res1 = 0.0; << 346 G4double eth = currentModel->LowEnergyLimit();
398 if(nullptr != baseParticle) { << 347 if(eth > 0.05*MeV && eth < 10.*MeV && escaled > eth) {
399 res1 = chargeSquare* << 348 G4double res1 = 0.0;
400 currentModel->ComputeDEDXPerVolume(mat, << 349 if(baseParticle) {
401 res0 = chargeSquare* << 350 res1 = currentModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
402 loweModel->ComputeDEDXPerVolume(mat, bas << 351 * chargeSquare;
403 } else { 352 } else {
404 res1 = currentModel->ComputeDEDXPerVolume( 353 res1 = currentModel->ComputeDEDXPerVolume(mat, p, eth, cut);
405 res0 = loweModel->ComputeDEDXPerVolume(mat <<
406 } 354 }
407 if(res1 > 0.0 && escaled > 0.0) { << 355 if(isIon) res1 += corr->HighOrderCorrections(p,mat,eth/massRatio);
408 res *= (1.0 + (res0/res1 - 1.0)*eth/escale << 356 G4double res0 = res1;
409 } << 357 if(FindEmModel(p, processName, eth-1.0*keV)) {
410 if(verbose > 1) { << 358 if(baseParticle) {
411 G4cout << "At boundary energy(MeV)= " << e << 359 res0 = currentModel->ComputeDEDXPerVolume(mat, baseParticle, eth, cut)
412 << " DEDX(MeV/mm)= " << res0*mm/MeV << " << 360 * chargeSquare;
413 << " after correction DEDX(MeV/mm)=" << r << 361 } else {
414 } << 362 res0 = currentModel->ComputeDEDXPerVolume(mat, p, eth, cut);
415 } << 363 }
416 // correction for ions <<
417 if(isIon) { <<
418 const G4double length = CLHEP::nm; <<
419 if(UpdateCouple(mat, cut)) { <<
420 G4double eloss = res*length; <<
421 dynParticle->SetKineticEnergy(kinEnergy); <<
422 currentModel->CorrectionsAlongStep(current <<
423 l <<
424 res = eloss/length; <<
425 <<
426 if(verbose > 1) { <<
427 G4cout << "After Corrections: DEDX(MeV/m <<
428 << " DEDX(MeV*cm^2/g)= " <<
429 << res*gram/(MeV*cm2*mat->GetDensity()) <<
430 } <<
431 } 364 }
>> 365 //G4cout << "eth= " << eth << " escaled= " << escaled << " res0= " << res0 << " res1= "
>> 366 // << res1 << " q2= " << chargeSquare << G4endl;
>> 367 res *= (1.0 + (res0/res1 - 1.0)*eth/escaled);
432 } 368 }
433 if(verbose > 0) { 369 if(verbose > 0) {
434 G4cout << "## E(MeV)= " << kinEnergy/MeV << 370 G4cout << "E(MeV)= " << kinEnergy/MeV
435 << " DEDX(MeV/mm)= " << res*mm/MeV << 371 << " DEDX(MeV/mm)= " << res*mm/MeV
436 << " DEDX(MeV*cm^2/g)= " << res*gram/ << 372 << " DEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
437 << " cut(MeV)= " << cut/MeV << 373 << " cut(MeV)= " << cut/MeV
438 << " " << p->GetParticleName() << 374 << " " << p->GetParticleName()
439 << " in " << currentMaterialName << 375 << " in " << mat->GetName()
440 << " Zi^2= " << chargeSquare << 376 << " Zi^2= " << chargeSquare
441 << " isIon=" << isIon << 377 << G4endl;
442 << G4endl; <<
443 } 378 }
444 } 379 }
445 } 380 }
446 return res; 381 return res;
447 } 382 }
448 383
449 //....oooOO0OOooo........oooOO0OOooo........oo 384 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
450 385
451 G4double G4EmCalculator::ComputeElectronicDEDX << 386 G4double G4EmCalculator::ComputeDEDX(G4double kinEnergy,
452 << 387 const G4String& particle,
453 << 388 const G4String& processName,
454 << 389 const G4String& material,
455 { << 390 G4double cut)
456 SetupMaterial(mat); << 391 {
457 G4double dedx = 0.0; << 392 return ComputeDEDX(kinEnergy,FindParticle(particle),processName,
458 if(UpdateParticle(part, kinEnergy)) { << 393 FindMaterial(material),cut);
459 <<
460 G4LossTableManager* lManager = G4LossTable <<
461 const std::vector<G4VEnergyLossProcess*> v <<
462 lManager->GetEnergyLossProcessVector(); <<
463 std::size_t n = vel.size(); <<
464 <<
465 //G4cout << "ComputeElectronicDEDX for " < <<
466 // << " n= " << n << G4endl; <<
467 <<
468 for(std::size_t i=0; i<n; ++i) { <<
469 if(vel[i]) { <<
470 auto p = static_cast<G4VProcess*>(vel[ <<
471 if(ActiveForParticle(part, p)) { <<
472 //G4cout << "idx= " << i << " " << ( <<
473 // << " " << (vel[i])->Particle()- <<
474 dedx += ComputeDEDX(kinEnergy,part,( <<
475 } <<
476 } <<
477 } <<
478 } <<
479 return dedx; <<
480 } 394 }
481 395
482 //....oooOO0OOooo........oooOO0OOooo........oo 396 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
483 397
484 G4double << 398 G4double G4EmCalculator::ComputeNuclearDEDX(G4double kinEnergy,
485 G4EmCalculator::ComputeDEDXForCutInRange(G4dou << 399 const G4ParticleDefinition* p,
486 const << 400 const G4Material* mat)
487 const <<
488 G4dou <<
489 { 401 {
490 SetupMaterial(mat); <<
491 G4double dedx = 0.0; <<
492 if(UpdateParticle(part, kinEnergy)) { <<
493 402
494 G4LossTableManager* lManager = G4LossTable << 403 G4double res = corr->NuclearDEDX(p, mat, kinEnergy, false);
495 const std::vector<G4VEnergyLossProcess*> v <<
496 lManager->GetEnergyLossProcessVector(); <<
497 std::size_t n = vel.size(); <<
498 <<
499 if(mat != cutMaterial) { <<
500 cutMaterial = mat; <<
501 cutenergy[0] = <<
502 ComputeEnergyCutFromRangeCut(rangecut, <<
503 cutenergy[1] = <<
504 ComputeEnergyCutFromRangeCut(rangecut, <<
505 cutenergy[2] = <<
506 ComputeEnergyCutFromRangeCut(rangecut, <<
507 } <<
508 <<
509 //G4cout << "ComputeElectronicDEDX for " < <<
510 // << " n= " << n << G4endl; <<
511 <<
512 for(std::size_t i=0; i<n; ++i) { <<
513 if(vel[i]) { <<
514 auto p = static_cast<G4VProcess*>(vel[ <<
515 if(ActiveForParticle(part, p)) { <<
516 //G4cout << "idx= " << i << " " << ( <<
517 // << " " << (vel[i])->Particle()-> <<
518 const G4ParticleDefinition* sec = (v <<
519 std::size_t idx = 0; <<
520 if(sec == G4Electron::Electron()) { <<
521 else if(sec == G4Positron::Positron( <<
522 404
523 dedx += ComputeDEDX(kinEnergy,part,( << 405 if(verbose > 1) {
524 mat,cutenergy[id << 406 G4cout << p->GetParticleName() << " E(MeV)= " << kinEnergy/MeV
525 } << 407 << " NuclearDEDX(MeV/mm)= " << res*mm/MeV
526 } << 408 << " NuclearDEDX(MeV*cm^2/g)= " << res*gram/(MeV*cm2*mat->GetDensity())
527 } << 409 << G4endl;
528 } 410 }
529 return dedx; << 411 return res;
530 } <<
531 <<
532 //....oooOO0OOooo........oooOO0OOooo........oo <<
533 <<
534 G4double G4EmCalculator::ComputeTotalDEDX(G4do <<
535 cons <<
536 cons <<
537 G4do <<
538 { <<
539 G4double dedx = ComputeElectronicDEDX(kinEne <<
540 if(mass > 700.*MeV) { dedx += ComputeNuclear <<
541 return dedx; <<
542 } 412 }
543 413
544 //....oooOO0OOooo........oooOO0OOooo........oo 414 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
545 415
546 G4double G4EmCalculator::ComputeNuclearDEDX(G4 416 G4double G4EmCalculator::ComputeNuclearDEDX(G4double kinEnergy,
547 const G4 << 417 const G4String& particle,
548 const G4 << 418 const G4String& material)
549 { 419 {
550 G4double res = 0.0; << 420 return ComputeNuclearDEDX(kinEnergy,FindParticle(particle),FindMaterial(material));
551 G4VEmProcess* nucst = FindDiscreteProcess(p, <<
552 if(nucst) { <<
553 G4VEmModel* mod = nucst->EmModel(); <<
554 if(mod) { <<
555 mod->SetFluctuationFlag(false); <<
556 res = mod->ComputeDEDXPerVolume(mat, p, <<
557 } <<
558 } <<
559 <<
560 if(verbose > 1) { <<
561 G4cout << p->GetParticleName() << " E(MeV <<
562 << " NuclearDEDX(MeV/mm)= " << res* <<
563 << " NuclearDEDX(MeV*cm^2/g)= " <<
564 << res*gram/(MeV*cm2*mat->GetDensit <<
565 << G4endl; <<
566 } <<
567 return res; <<
568 } 421 }
569 422
570 //....oooOO0OOooo........oooOO0OOooo........oo 423 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
571 424
572 G4double G4EmCalculator::ComputeCrossSectionPe 425 G4double G4EmCalculator::ComputeCrossSectionPerVolume(
573 426 G4double kinEnergy,
574 c 427 const G4ParticleDefinition* p,
575 c 428 const G4String& processName,
576 c << 429 const G4Material* mat,
577 430 G4double cut)
578 { 431 {
579 SetupMaterial(mat); <<
580 G4double res = 0.0; 432 G4double res = 0.0;
581 if(UpdateParticle(p, kinEnergy)) { << 433 if(FindEmModel(p, processName, kinEnergy)) {
582 if(FindEmModel(p, processName, kinEnergy)) << 434 UpdateParticle(p, kinEnergy);
583 G4double e = kinEnergy; << 435 G4double e = kinEnergy;
584 G4double aCut = std::max(cut, theParamet << 436 if(baseParticle) {
585 if(baseParticle) { << 437 e *= kinEnergy*massRatio;
586 e *= kinEnergy*massRatio; << 438 res = currentModel->CrossSectionPerVolume(mat, baseParticle, e, cut, e)*chargeSquare;
587 res = currentModel->CrossSectionPerVol << 439 } else {
588 mat, baseParticle, e, aCut, e) * << 440 res = currentModel->CrossSectionPerVolume(mat, p, e, cut, e);
589 } else { << 441 }
590 res = currentModel->CrossSectionPerVol << 442 if(verbose>0) {
591 } << 443 G4cout << "E(MeV)= " << kinEnergy/MeV
592 if(verbose>0) { << 444 << " cross(cm^2/g)= " << res*gram/cm2
593 G4cout << "G4EmCalculator::ComputeXSPe << 445 << " " << p->GetParticleName()
594 << kinEnergy/MeV << 446 << " in " << mat->GetName()
595 << " cross(cm-1)= " << res*cm << 447 << G4endl;
596 << " cut(keV)= " << aCut/keV <<
597 << " " << p->GetParticleName( <<
598 << " in " << mat->GetName() <<
599 << G4endl; <<
600 } <<
601 } 448 }
602 } 449 }
603 return res; 450 return res;
604 } 451 }
605 452
606 //....oooOO0OOooo........oooOO0OOooo........oo 453 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
607 454
608 G4double << 455 G4double G4EmCalculator::ComputeCrossSectionPerVolume(
609 G4EmCalculator::ComputeCrossSectionPerAtom(G4d << 456 G4double kinEnergy,
610 con << 457 const G4String& particle,
611 con << 458 const G4String& processName,
612 G4d << 459 const G4String& material,
613 G4d << 460 G4double cut)
614 { 461 {
615 G4double res = 0.0; << 462 return ComputeCrossSectionPerVolume(kinEnergy,FindParticle(particle),processName,
616 if(UpdateParticle(p, kinEnergy)) { << 463 FindMaterial(material),cut);
617 G4int iz = G4lrint(Z); <<
618 CheckMaterial(iz); <<
619 if(FindEmModel(p, processName, kinEnergy)) <<
620 G4double e = kinEnergy; <<
621 G4double aCut = std::max(cut, theParamet <<
622 if(baseParticle) { <<
623 e *= kinEnergy*massRatio; <<
624 currentModel->InitialiseForElement(bas <<
625 res = currentModel->ComputeCrossSectio <<
626 baseParticle, e, Z, A, aCut) * c <<
627 } else { <<
628 currentModel->InitialiseForElement(p, <<
629 res = currentModel->ComputeCrossSectio <<
630 } <<
631 if(verbose>0) { <<
632 G4cout << "E(MeV)= " << kinEnergy/MeV <<
633 << " cross(barn)= " << res/barn <<
634 << " " << p->GetParticleName( <<
635 << " Z= " << Z << " A= " << A/ <<
636 << " cut(keV)= " << aCut/keV <<
637 << G4endl; <<
638 } <<
639 } <<
640 } <<
641 return res; <<
642 } 464 }
643 465
644 //....oooOO0OOooo........oooOO0OOooo........oo 466 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
645 467
646 G4double << 468 G4double G4EmCalculator::ComputeCrossSectionPerAtom(
647 G4EmCalculator::ComputeCrossSectionPerShell(G4 << 469 G4double kinEnergy,
648 co << 470 const G4ParticleDefinition* p,
649 co << 471 const G4String& processName,
650 G4 << 472 G4double Z, G4double A,
651 G4 << 473 G4double cut)
652 { 474 {
653 G4double res = 0.0; 475 G4double res = 0.0;
654 if(UpdateParticle(p, kinEnergy)) { << 476 if(FindEmModel(p, processName, kinEnergy)) {
655 CheckMaterial(Z); << 477 UpdateParticle(p, kinEnergy);
656 if(FindEmModel(p, processName, kinEnergy)) << 478 G4double e = kinEnergy;
657 G4double e = kinEnergy; << 479 if(baseParticle) {
658 G4double aCut = std::max(cut, theParamet << 480 e *= kinEnergy*massRatio;
659 if(nullptr != baseParticle) { << 481 res = currentModel->ComputeCrossSectionPerAtom(baseParticle, e, Z, A, cut)*chargeSquare;
660 e *= kinEnergy*massRatio; << 482 } else {
661 currentModel->InitialiseForElement(bas << 483 res = currentModel->ComputeCrossSectionPerAtom(p, e, Z, A, cut);
662 res = << 484 }
663 currentModel->ComputeCrossSectionPer << 485 if(verbose>0) {
664 << 486 G4cout << "E(MeV)= " << kinEnergy/MeV
665 } else { << 487 << " cross(barn)= " << res/barn
666 currentModel->InitialiseForElement(p, << 488 << " " << p->GetParticleName()
667 res = currentModel->ComputeCrossSectio << 489 << " Z= " << Z << " A= " << A
668 } << 490 << G4endl;
669 if(verbose>0) { <<
670 G4cout << "E(MeV)= " << kinEnergy/MeV <<
671 << " cross(barn)= " << res/barn <<
672 << " " << p->GetParticleName( <<
673 << " Z= " << Z << " shellIdx= <<
674 << " cut(keV)= " << aCut/keV <<
675 << G4endl; <<
676 } <<
677 } 491 }
678 } 492 }
679 return res; 493 return res;
680 } 494 }
681 495
682 //....oooOO0OOooo........oooOO0OOooo........oo 496 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
683 497
684 G4double << 498 G4double G4EmCalculator::ComputeCrossSectionPerAtom(G4double kinEnergy,
685 G4EmCalculator::ComputeGammaAttenuationLength( << 499 const G4String& particle,
686 << 500 const G4String& processName,
>> 501 const G4Element* elm,
>> 502 G4double cut)
687 { 503 {
688 G4double res = 0.0; << 504 return ComputeCrossSectionPerAtom(kinEnergy,FindParticle(particle),processName,
689 const G4ParticleDefinition* gamma = G4Gamma: << 505 elm->GetZ(),elm->GetN(),cut);
690 res += ComputeCrossSectionPerVolume(kinEnerg <<
691 res += ComputeCrossSectionPerVolume(kinEnerg <<
692 res += ComputeCrossSectionPerVolume(kinEnerg <<
693 res += ComputeCrossSectionPerVolume(kinEnerg <<
694 if(res > 0.0) { res = 1.0/res; } <<
695 return res; <<
696 } <<
697 <<
698 //....oooOO0OOooo........oooOO0OOooo........oo <<
699 <<
700 G4double G4EmCalculator::ComputeShellIonisatio <<
701 const <<
702 G4int <<
703 G4Ato <<
704 G4dou <<
705 const <<
706 { <<
707 G4double res = 0.0; <<
708 const G4ParticleDefinition* p = FindParticle <<
709 G4VAtomDeexcitation* ad = manager->AtomDeexc <<
710 if(p && ad) { <<
711 res = ad->ComputeShellIonisationCrossSecti <<
712 <<
713 } <<
714 return res; <<
715 } 506 }
716 507
717 //....oooOO0OOooo........oooOO0OOooo........oo 508 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
718 509
719 G4double G4EmCalculator::ComputeMeanFreePath(G 510 G4double G4EmCalculator::ComputeMeanFreePath(G4double kinEnergy,
720 c 511 const G4ParticleDefinition* p,
721 c 512 const G4String& processName,
722 c << 513 const G4Material* mat,
723 G << 514 G4double cut)
724 { 515 {
725 G4double mfp = DBL_MAX; 516 G4double mfp = DBL_MAX;
726 G4double x = << 517 G4double x = ComputeCrossSectionPerVolume(kinEnergy, p, processName, mat, cut);
727 ComputeCrossSectionPerVolume(kinEnergy, p, << 518 if(x > 0.0) mfp = 1.0/x;
728 if(x > 0.0) { mfp = 1.0/x; } <<
729 if(verbose>1) { 519 if(verbose>1) {
730 G4cout << "E(MeV)= " << kinEnergy/MeV 520 G4cout << "E(MeV)= " << kinEnergy/MeV
731 << " MFP(mm)= " << mfp/mm << 521 << " MFP(mm)= " << mfp/mm
732 << " " << p->GetParticleName() << 522 << " " << p->GetParticleName()
733 << " in " << mat->GetName() << 523 << " in " << mat->GetName()
734 << G4endl; << 524 << G4endl;
735 } 525 }
736 return mfp; 526 return mfp;
737 } 527 }
738 528
739 //....oooOO0OOooo........oooOO0OOooo........oo 529 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
740 530
741 G4double G4EmCalculator::ComputeEnergyCutFromR << 531 G4double G4EmCalculator::ComputeMeanFreePath(G4double kinEnergy,
742 G4double range, << 532 const G4String& particle,
743 const G4ParticleDefin << 533 const G4String& processName,
744 const G4Material* mat << 534 const G4String& material,
>> 535 G4double cut)
745 { 536 {
746 return G4ProductionCutsTable::GetProductionC << 537 return ComputeMeanFreePath(kinEnergy,FindParticle(particle),processName,
747 ConvertRangeToEnergy(part, mat, range); << 538 FindMaterial(material),cut);
748 } 539 }
749 540
750 //....oooOO0OOooo........oooOO0OOooo........oo 541 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
751 542
752 G4bool G4EmCalculator::UpdateParticle(const G4 << 543 G4bool G4EmCalculator::UpdateParticle(const G4ParticleDefinition* p, G4double kinEnergy)
753 G4double <<
754 { 544 {
755 if(p != currentParticle) { 545 if(p != currentParticle) {
756 <<
757 // new particle <<
758 currentParticle = p; 546 currentParticle = p;
759 dynParticle->SetDefinition(const_cast<G4Pa << 547 baseParticle = 0;
760 dynParticle->SetKineticEnergy(kinEnergy); <<
761 baseParticle = nullptr; <<
762 currentParticleName = p->GetParticleName() 548 currentParticleName = p->GetParticleName();
763 massRatio = 1.0; << 549 }
764 mass = p->GetPDGMass(); << 550 isApplicable = false;
765 chargeSquare = 1.0; << 551 if(p->GetPDGCharge() != 0.0) {
766 currentProcess = manager->GetEnergyLossPr << 552 const G4VEnergyLossProcess* proc = FindEnergyLossProcess(p);
767 currentProcessName = ""; << 553 if(proc) {
768 isIon = false; << 554 isApplicable = true;
769 << 555 baseParticle = proc->BaseParticle();
770 // ionisation process exist << 556 massRatio = 1.0;
771 if(nullptr != currentProcess) { << 557 G4double qbase = 1.0;
772 currentProcessName = currentProcess->Get << 558 G4double q = p->GetPDGCharge()/eplus;
773 baseParticle = currentProcess->BaseParti <<
774 if(currentProcessName == "ionIoni" && p- <<
775 baseParticle = theGenericIon; <<
776 isIon = true; <<
777 } <<
778 559
779 // base particle is used << 560 if(baseParticle) {
780 if(nullptr != baseParticle) { <<
781 massRatio = baseParticle->GetPDGMass() 561 massRatio = baseParticle->GetPDGMass()/p->GetPDGMass();
782 G4double q = p->GetPDGCharge()/basePar << 562 qbase = baseParticle->GetPDGCharge()/eplus;
783 chargeSquare = q*q; << 563 }
784 } << 564
785 } << 565 if(p->GetParticleType() == "nucleus") {
786 } << 566 isIon = true;
787 // Effective charge for ions << 567 chargeSquare =
788 if(isIon && nullptr != currentProcess) { << 568 ionEffCharge->EffectiveChargeSquareRatio(p, currentMaterial, kinEnergy);
789 chargeSquare = << 569 //G4double emax = currentModel->HighEnergyLimit(p);
790 corr->EffectiveChargeSquareRatio(p, curr << 570 //if(emax < 100.*MeV) chargeSquare *= q*q /
791 currentProcess->SetDynamicMassCharge(massR << 571 // ionEffCharge->EffectiveChargeSquareRatio(p, currentMaterial, emax/massRatio);
792 if(verbose>1) { << 572
793 G4cout <<"\n NewIon: massR= "<< massRati << 573 } else {
794 << chargeSquare << " " << currentProce << 574 isIon = false;
>> 575 chargeSquare = q*q;
>> 576 }
>> 577
>> 578 chargeSquare /= (qbase*qbase);
>> 579
795 } 580 }
796 } 581 }
797 return true; << 582 return isApplicable;
798 } 583 }
799 584
800 //....oooOO0OOooo........oooOO0OOooo........oo 585 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
801 586
802 const G4ParticleDefinition* G4EmCalculator::Fi 587 const G4ParticleDefinition* G4EmCalculator::FindParticle(const G4String& name)
803 { 588 {
804 const G4ParticleDefinition* p = nullptr; <<
805 if(name != currentParticleName) { 589 if(name != currentParticleName) {
806 p = G4ParticleTable::GetParticleTable()->F << 590 currentParticle = G4ParticleTable::GetParticleTable()->FindParticle(name);
807 if(nullptr == p) { << 591 if(!currentParticle) {
808 G4cout << "### WARNING: G4EmCalculator:: << 592 G4cout << "### WARNING: G4EmCalculator::FindParticle fails to find " << name << G4endl;
809 << name << G4endl; <<
810 } 593 }
811 } else { << 594 currentParticleName = name;
812 p = currentParticle; <<
813 } 595 }
814 return p; << 596 return currentParticle;
815 } <<
816 <<
817 //....oooOO0OOooo........oooOO0OOooo........oo <<
818 <<
819 const G4ParticleDefinition* G4EmCalculator::Fi <<
820 { <<
821 const G4ParticleDefinition* p = ionTable->Ge <<
822 return p; <<
823 } 597 }
824 598
825 //....oooOO0OOooo........oooOO0OOooo........oo 599 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
826 600
827 const G4Material* G4EmCalculator::FindMaterial 601 const G4Material* G4EmCalculator::FindMaterial(const G4String& name)
828 { 602 {
829 if(name != currentMaterialName) { 603 if(name != currentMaterialName) {
830 SetupMaterial(G4Material::GetMaterial(name << 604 currentMaterial = G4Material::GetMaterial(name);
831 if(nullptr == currentMaterial) { << 605 currentMaterialName = name;
832 G4cout << "### WARNING: G4EmCalculator:: << 606 if(!currentMaterial)
833 << name << G4endl; << 607 G4cout << "### WARNING: G4EmCalculator::FindMaterial fails to find " << name << G4endl;
834 } <<
835 } 608 }
836 return currentMaterial; 609 return currentMaterial;
837 } 610 }
838 611
839 //....oooOO0OOooo........oooOO0OOooo........oo 612 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
840 613
841 const G4Region* G4EmCalculator::FindRegion(con 614 const G4Region* G4EmCalculator::FindRegion(const G4String& reg)
842 { 615 {
843 return G4EmUtility::FindRegion(reg); << 616 return G4RegionStore::GetInstance()->GetRegion(reg);
844 } 617 }
845 618
846 //....oooOO0OOooo........oooOO0OOooo........oo 619 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
847 620
848 const G4MaterialCutsCouple* G4EmCalculator::Fi << 621 const G4MaterialCutsCouple* G4EmCalculator::FindCouple(const G4Material* material,
849 const G4Material* << 622 const G4Region* region)
850 const G4Region* re <<
851 { 623 {
852 const G4MaterialCutsCouple* couple = nullptr << 624 if(!material) return 0;
853 SetupMaterial(material); << 625 currentMaterial = material;
854 if(nullptr != currentMaterial) { << 626 currentMaterialName = material->GetName();
855 // Access to materials << 627 // Access to materials
856 const G4ProductionCutsTable* theCoupleTabl << 628 const G4ProductionCutsTable* theCoupleTable=
857 G4ProductionCutsTable::GetProductionCuts << 629 G4ProductionCutsTable::GetProductionCutsTable();
858 const G4Region* r = region; << 630 const G4Region* r = region;
859 if(nullptr != r) { << 631 if(!r) r = G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld");
860 couple = theCoupleTable->GetMaterialCuts << 632
861 << 633 return theCoupleTable->GetMaterialCutsCouple(material,r->GetProductionCuts());
862 } else { << 634
863 G4RegionStore* store = G4RegionStore::Ge <<
864 std::size_t nr = store->size(); <<
865 if(0 < nr) { <<
866 for(std::size_t i=0; i<nr; ++i) { <<
867 couple = theCoupleTable->GetMaterial <<
868 material, ((*store)[i])->GetProduc <<
869 if(nullptr != couple) { break; } <<
870 } <<
871 } <<
872 } <<
873 } <<
874 if(nullptr == couple) { <<
875 G4ExceptionDescription ed; <<
876 ed << "G4EmCalculator::FindCouple: fail fo <<
877 << currentMaterialName << ">"; <<
878 if(region) { ed << " and region " << regio <<
879 G4Exception("G4EmCalculator::FindCouple", <<
880 FatalException, ed); <<
881 } <<
882 return couple; <<
883 } 635 }
884 636
885 //....oooOO0OOooo........oooOO0OOooo........oo 637 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
886 638
887 G4bool G4EmCalculator::UpdateCouple(const G4Ma 639 G4bool G4EmCalculator::UpdateCouple(const G4Material* material, G4double cut)
888 { 640 {
889 SetupMaterial(material); << 641 if(!material) return false;
890 if(!currentMaterial) { return false; } << 642 currentMaterial = material;
891 for (G4int i=0; i<nLocalMaterials; ++i) { << 643 currentMaterialName = material->GetName();
>> 644 for (G4int i=0; i<nLocalMaterials; i++) {
892 if(material == localMaterials[i] && cut == 645 if(material == localMaterials[i] && cut == localCuts[i]) {
893 currentCouple = localCouples[i]; 646 currentCouple = localCouples[i];
894 currentCoupleIndex = currentCouple->GetI 647 currentCoupleIndex = currentCouple->GetIndex();
895 currentCut = cut; 648 currentCut = cut;
896 return true; 649 return true;
897 } 650 }
898 } 651 }
899 const G4MaterialCutsCouple* cc = new G4Mater 652 const G4MaterialCutsCouple* cc = new G4MaterialCutsCouple(material);
900 localMaterials.push_back(material); 653 localMaterials.push_back(material);
901 localCouples.push_back(cc); 654 localCouples.push_back(cc);
902 localCuts.push_back(cut); 655 localCuts.push_back(cut);
903 ++nLocalMaterials; << 656 nLocalMaterials++;
904 currentCouple = cc; 657 currentCouple = cc;
905 currentCoupleIndex = currentCouple->GetIndex 658 currentCoupleIndex = currentCouple->GetIndex();
906 currentCut = cut; 659 currentCut = cut;
907 return true; 660 return true;
908 } 661 }
909 662
910 //....oooOO0OOooo........oooOO0OOooo........oo 663 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
911 664
912 void G4EmCalculator::FindLambdaTable(const G4P 665 void G4EmCalculator::FindLambdaTable(const G4ParticleDefinition* p,
913 const G4S << 666 const G4String& processName)
914 G4double <<
915 { 667 {
916 // Search for the process 668 // Search for the process
917 if (!currentLambda || p != lambdaParticle || << 669 if (p != currentParticle || processName != currentName) {
918 lambdaName = processName; << 670 currentName = processName;
919 currentLambda = nullptr; << 671 currentLambda = 0;
920 lambdaParticle = p; << 672 G4LossTableManager* lManager = G4LossTableManager::Instance();
921 isApplicable = false; << 673 const std::vector<G4VEnergyLossProcess*> vel = lManager->GetEnergyLossProcessVector();
922 << 674 G4int n = vel.size();
923 const G4ParticleDefinition* part = (isIon) << 675 for(G4int i=0; i<n; i++) {
924 << 676 if((vel[i])->GetProcessName() == currentName && (vel[i])->Particle() == p) {
925 // Search for energy loss process << 677 currentLambda = (vel[i])->LambdaTable();
926 currentName = processName; << 678 isApplicable = true;
927 currentModel = nullptr; << 679 break;
928 loweModel = nullptr; << 680 }
929 << 681 }
930 G4VEnergyLossProcess* elproc = FindEnLossP << 682 if(!currentLambda) {
931 if(nullptr != elproc) { << 683 const std::vector<G4VEmProcess*> vem = lManager->GetEmProcessVector();
932 currentLambda = elproc->LambdaTable(); << 684 G4int n = vem.size();
933 proctype = 0; << 685 for(G4int i=0; i<n; i++) {
934 if(nullptr != currentLambda) { << 686 if((vem[i])->GetProcessName() == currentName && (vem[i])->Particle() == p) {
935 isApplicable = true; << 687 currentLambda = (vem[i])->LambdaTable();
936 if(verbose>1) { << 688 isApplicable = true;
937 G4cout << "G4VEnergyLossProcess is f << 689 break;
938 << G4endl; <<
939 } <<
940 } <<
941 curProcess = elproc; <<
942 return; <<
943 } <<
944 <<
945 // Search for discrete process <<
946 G4VEmProcess* proc = FindDiscreteProcess(p <<
947 if(nullptr != proc) { <<
948 currentLambda = proc->LambdaTable(); <<
949 proctype = 1; <<
950 if(nullptr != currentLambda) { <<
951 isApplicable = true; <<
952 if(verbose>1) { <<
953 G4cout << "G4VEmProcess is found out <<
954 } 690 }
955 } 691 }
956 curProcess = proc; <<
957 return; <<
958 } 692 }
959 << 693 if(!currentLambda) {
960 // Search for msc process << 694 const std::vector<G4VMultipleScattering*> vmsc = lManager->GetMultipleScatteringVector();
961 G4VMultipleScattering* msc = FindMscProces << 695 G4int n = vmsc.size();
962 if(nullptr != msc) { << 696 for(G4int i=0; i<n; i++) {
963 currentModel = msc->SelectModel(kinEnerg << 697 if((vmsc[i])->GetProcessName() == currentName && (vmsc[i])->Particle() == p) {
964 proctype = 2; << 698 currentLambda = (vmsc[i])->LambdaTable();
965 if(nullptr != currentModel) { << 699 isApplicable = true;
966 currentLambda = currentModel->GetCross << 700 break;
967 if(nullptr != currentLambda) { <<
968 isApplicable = true; <<
969 if(verbose>1) { <<
970 G4cout << "G4VMultipleScattering i <<
971 << G4endl; <<
972 } <<
973 } 701 }
974 } 702 }
975 curProcess = msc; <<
976 } 703 }
977 } 704 }
978 } 705 }
979 706
980 //....oooOO0OOooo........oooOO0OOooo........oo 707 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
981 708
982 G4bool G4EmCalculator::FindEmModel(const G4Par 709 G4bool G4EmCalculator::FindEmModel(const G4ParticleDefinition* p,
983 const G4Str 710 const G4String& processName,
984 G4 << 711 G4double kinEnergy)
985 { 712 {
986 isApplicable = false; << 713 G4bool res = false;
987 if(nullptr == p || nullptr == currentMateria <<
988 G4cout << "G4EmCalculator::FindEmModel WAR <<
989 << " or materail defined; particle: <<
990 return isApplicable; <<
991 } <<
992 G4String partname = p->GetParticleName(); <<
993 G4double scaledEnergy = kinEnergy*massRatio; <<
994 const G4ParticleDefinition* part = (isIon) ? <<
995 <<
996 if(verbose > 1) { 714 if(verbose > 1) {
997 G4cout << "## G4EmCalculator::FindEmModel << 715 G4cout << "G4EmCalculator::FindEmModel for " << p->GetParticleName()
998 << " (type= " << p->GetParticleType << 716 << " and " << processName << " at e(MeV)= " << kinEnergy
999 << ") and " << processName << " at <<
1000 << G4endl; 717 << G4endl;
1001 if(p != part) { G4cout << " GenericIon i <<
1002 } 718 }
1003 << 719 // Search for the process
1004 // Search for energy loss process <<
1005 currentName = processName; 720 currentName = processName;
1006 currentModel = nullptr; << 721 currentModel = 0;
1007 loweModel = nullptr; << 722 size_t idx = 0;
1008 std::size_t idx = 0; << 723 G4LossTableManager* lManager = G4LossTableManager::Instance();
1009 << 724 const std::vector<G4VEnergyLossProcess*> vel = lManager->GetEnergyLossProcessVector();
1010 G4VEnergyLossProcess* elproc = FindEnLossPr << 725 G4int n = vel.size();
1011 if(nullptr != elproc) { << 726 for(G4int i=0; i<n; i++) {
1012 currentModel = elproc->SelectModelForMate << 727 if((vel[i])->GetProcessName() == currentName && (vel[i])->Particle() == p) {
1013 currentModel->InitialiseForMaterial(part, << 728 const G4ParticleDefinition* bp = (vel[i])->BaseParticle();
1014 currentModel->SetupForMaterial(part, curr << 729 // G4cout << "i= " << i << " bp= " << bp << G4endl;
1015 G4double eth = currentModel->LowEnergyLim << 730 if(!bp) {
1016 if(eth > 0.0) { << 731 currentModel = (vel[i])->SelectModelForMaterial(kinEnergy, idx);
1017 loweModel = elproc->SelectModelForMater << 732 } else {
1018 if(loweModel == currentModel) { loweMod << 733 G4double e = kinEnergy*bp->GetPDGMass()/p->GetPDGMass();
1019 else { << 734 for(G4int j=0; j<n; j++) {
1020 loweModel->InitialiseForMaterial(part << 735 if((vel[j])->Particle() == bp) {
1021 loweModel->SetupForMaterial(part, cur << 736 currentModel = (vel[j])->SelectModelForMaterial(e, idx);
1022 } << 737 break;
1023 } << 738 }
1024 } << 739 }
1025 <<
1026 // Search for discrete process <<
1027 if(nullptr == currentModel) { <<
1028 G4VEmProcess* proc = FindDiscreteProcess( <<
1029 if(nullptr != proc) { <<
1030 currentModel = proc->SelectModelForMate <<
1031 currentModel->InitialiseForMaterial(par <<
1032 currentModel->SetupForMaterial(part, cu <<
1033 G4double eth = currentModel->LowEnergyL <<
1034 if(eth > 0.0) { <<
1035 loweModel = proc->SelectModelForMater <<
1036 if(loweModel == currentModel) { loweM <<
1037 else { <<
1038 loweModel->InitialiseForMaterial(pa <<
1039 loweModel->SetupForMaterial(part, c <<
1040 } <<
1041 } <<
1042 } <<
1043 } <<
1044 <<
1045 // Search for msc process <<
1046 if(nullptr == currentModel) { <<
1047 G4VMultipleScattering* proc = FindMscProc <<
1048 if(nullptr != proc) { <<
1049 currentModel = proc->SelectModel(kinEne <<
1050 loweModel = nullptr; <<
1051 } <<
1052 } <<
1053 if(nullptr != currentModel) { <<
1054 if(loweModel == currentModel) { loweModel <<
1055 isApplicable = true; <<
1056 currentModel->InitialiseForMaterial(part, <<
1057 if(loweModel) { <<
1058 loweModel->InitialiseForMaterial(part, <<
1059 } <<
1060 if(verbose > 1) { <<
1061 G4cout << " Model <" << currentModel- <<
1062 << "> Emin(MeV)= " << currentMod <<
1063 << " for " << part->GetParticleN <<
1064 if(nullptr != elproc) { <<
1065 G4cout << " and " << elproc->GetProce <<
1066 << G4endl; <<
1067 } <<
1068 if(nullptr != loweModel) { <<
1069 G4cout << " LowEnergy model <" << low <<
1070 } <<
1071 G4cout << G4endl; <<
1072 } <<
1073 } <<
1074 return isApplicable; <<
1075 } <<
1076 <<
1077 //....oooOO0OOooo........oooOO0OOooo........o <<
1078 <<
1079 G4VEnergyLossProcess* <<
1080 G4EmCalculator::FindEnLossProcess(const G4Par <<
1081 const G4Str <<
1082 { <<
1083 G4VEnergyLossProcess* proc = nullptr; <<
1084 const std::vector<G4VEnergyLossProcess*> v <<
1085 manager->GetEnergyLossProcessVector(); <<
1086 std::size_t n = v.size(); <<
1087 for(std::size_t i=0; i<n; ++i) { <<
1088 if((v[i])->GetProcessName() == processNam <<
1089 auto p = static_cast<G4VProcess*>(v[i]) <<
1090 if(ActiveForParticle(part, p)) { <<
1091 proc = v[i]; <<
1092 break; <<
1093 } 740 }
>> 741 break;
1094 } 742 }
1095 } 743 }
1096 return proc; << 744 if(!currentModel) {
1097 } << 745 const std::vector<G4VEmProcess*> vem = lManager->GetEmProcessVector();
1098 << 746 G4int n = vem.size();
1099 //....oooOO0OOooo........oooOO0OOooo........o << 747 for(G4int i=0; i<n; i++) {
1100 << 748 if((vem[i])->GetProcessName() == currentName && (vem[i])->Particle() == p) {
1101 G4VEmProcess* << 749 currentModel = (vem[i])->SelectModelForMaterial(kinEnergy, idx);
1102 G4EmCalculator::FindDiscreteProcess(const G4P << 750 isApplicable = true;
1103 const G4S <<
1104 { <<
1105 G4VEmProcess* proc = nullptr; <<
1106 auto v = manager->GetEmProcessVector(); <<
1107 std::size_t n = v.size(); <<
1108 for(std::size_t i=0; i<n; ++i) { <<
1109 const G4String& pName = v[i]->GetProcessN <<
1110 if(pName == "GammaGeneralProc") { <<
1111 proc = v[i]->GetEmProcess(processName); <<
1112 break; <<
1113 } else if(pName == processName) { <<
1114 const auto p = static_cast<G4VProcess*> <<
1115 if(ActiveForParticle(part, p)) { <<
1116 proc = v[i]; <<
1117 break; 751 break;
1118 } 752 }
1119 } 753 }
1120 } 754 }
1121 return proc; << 755 if(!currentModel) {
1122 } << 756 const std::vector<G4VMultipleScattering*> vmsc = lManager->GetMultipleScatteringVector();
1123 << 757 G4int n = vmsc.size();
1124 //....oooOO0OOooo........oooOO0OOooo........o << 758 for(G4int i=0; i<n; i++) {
1125 << 759 if((vmsc[i])->GetProcessName() == currentName && (vmsc[i])->Particle() == p) {
1126 G4VMultipleScattering* << 760 currentModel = (vmsc[i])->SelectModelForMaterial(kinEnergy, idx);
1127 G4EmCalculator::FindMscProcess(const G4Partic << 761 isApplicable = true;
1128 const G4String <<
1129 { <<
1130 G4VMultipleScattering* proc = nullptr; <<
1131 const std::vector<G4VMultipleScattering*> v <<
1132 manager->GetMultipleScatteringVector(); <<
1133 std::size_t n = v.size(); <<
1134 for(std::size_t i=0; i<n; ++i) { <<
1135 if((v[i])->GetProcessName() == processNam <<
1136 auto p = static_cast<G4VProcess*>(v[i]) <<
1137 if(ActiveForParticle(part, p)) { <<
1138 proc = v[i]; <<
1139 break; 762 break;
1140 } 763 }
1141 } 764 }
1142 } 765 }
1143 return proc; << 766 if(currentModel) res = true;
1144 } <<
1145 <<
1146 //....oooOO0OOooo........oooOO0OOooo........o <<
1147 <<
1148 G4VProcess* G4EmCalculator::FindProcess(const <<
1149 const <<
1150 { <<
1151 G4VProcess* proc = nullptr; <<
1152 const G4ProcessManager* procman = part->Get <<
1153 G4ProcessVector* pv = procman->GetProcessLi <<
1154 G4int nproc = (G4int)pv->size(); <<
1155 for(G4int i=0; i<nproc; ++i) { <<
1156 if(processName == (*pv)[i]->GetProcessNam <<
1157 proc = (*pv)[i]; <<
1158 break; <<
1159 } <<
1160 } <<
1161 return proc; <<
1162 } <<
1163 <<
1164 //....oooOO0OOooo........oooOO0OOooo........o <<
1165 <<
1166 G4bool G4EmCalculator::ActiveForParticle(cons <<
1167 G4VP <<
1168 { <<
1169 G4ProcessManager* pm = part->GetProcessMana <<
1170 G4ProcessVector* pv = pm->GetProcessList(); <<
1171 G4int n = (G4int)pv->size(); <<
1172 G4bool res = false; <<
1173 for(G4int i=0; i<n; ++i) { <<
1174 if((*pv)[i] == proc) { <<
1175 if(pm->GetProcessActivation(i)) { res = <<
1176 break; <<
1177 } <<
1178 } <<
1179 return res; 767 return res;
1180 } 768 }
1181 769
1182 //....oooOO0OOooo........oooOO0OOooo........o 770 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1183 771
1184 void G4EmCalculator::SetupMaterial(const G4Ma << 772 const G4VEnergyLossProcess* G4EmCalculator::FindEnergyLossProcess(
1185 { << 773 const G4ParticleDefinition* p)
1186 if(mat) { <<
1187 currentMaterial = mat; <<
1188 currentMaterialName = mat->GetName(); <<
1189 } else { <<
1190 currentMaterial = nullptr; <<
1191 currentMaterialName = ""; <<
1192 } <<
1193 } <<
1194 <<
1195 //....oooOO0OOooo........oooOO0OOooo........o <<
1196 <<
1197 void G4EmCalculator::SetupMaterial(const G4St <<
1198 { <<
1199 SetupMaterial(nist->FindOrBuildMaterial(mna <<
1200 } <<
1201 <<
1202 //....oooOO0OOooo........oooOO0OOooo........o <<
1203 <<
1204 void G4EmCalculator::CheckMaterial(G4int Z) <<
1205 { 774 {
1206 G4bool isFound = false; << 775 const G4VEnergyLossProcess* elp = 0;
1207 if(nullptr != currentMaterial) { << 776
1208 G4int nn = (G4int)currentMaterial->GetNum << 777 G4LossTableManager* lManager = G4LossTableManager::Instance();
1209 for(G4int i=0; i<nn; ++i) { << 778 const std::vector<G4VEnergyLossProcess*> vel = lManager->GetEnergyLossProcessVector();
1210 if(Z == currentMaterial->GetElement(i)- << 779 G4int n = vel.size();
1211 isFound = true; << 780 for(G4int i=0; i<n; i++) {
1212 break; << 781 if((vel[i])->Particle() == p) {
1213 } << 782 elp = vel[i];
>> 783 break;
1214 } 784 }
1215 } 785 }
1216 if(!isFound) { << 786 return elp;
1217 SetupMaterial(nist->FindOrBuildSimpleMate <<
1218 } <<
1219 } 787 }
1220 788
1221 //....oooOO0OOooo........oooOO0OOooo........o 789 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1222 790
1223 void G4EmCalculator::SetVerbose(G4int verb) 791 void G4EmCalculator::SetVerbose(G4int verb)
1224 { 792 {
1225 verbose = verb; 793 verbose = verb;
1226 } 794 }
1227 795
1228 //....oooOO0OOooo........oooOO0OOooo........o 796 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1229 797
1230 798